Density functional theory of the structure of magnesium-doped helium nanodroplets

We have studied the structure of He droplets doped with magnesium atoms using density functional theory. We have found that the solvation properties of this system strongly depend on the size of the He droplet. For small drops, Mg resides in a deep surface state, whereas for large-size drops it is fully solvated but radially delocalized in their interior. We have studied the 3s3pP1←3sS0 transition of the dopant, and have compared our results with experimental data from laser-induced fluorescence LIF . Line-broadening effects due to the coupling of dynamical deformations of the surrounding helium with the dipole excitation of the impurity are explicitly taken into account. We show that the Mg radial delocalization inside large droplets may help reconcile the apparently contradictory solvation properties of magnesium as provided by LIF and electronimpact ionization experiments. The structure of He drops doped with two magnesium atoms is also studied and used to interpret the results of resonant two-photon-ionization R2PI and LIF experiments. We have found that the two solvated Mg atoms do not easily merge into a dimer, but rather form a weakly bound state due to the presence of an energy barrier caused by the helium environment that keeps them some 9.5 Å apart, preventing the formation of the Mg2 cluster. From this observation, we suggest that Mg atoms in He drops may form, under suitable conditions, a soft “foamlike” aggregate rather than coalesce into a compact metallic cluster. Our findings are in qualitative agreement with recent R2PI experimental evidence. We predict that, contrarily, Mg atoms adsorbed in He droplets do not form such metastable aggregates.


I. INTRODUCTION
The study of magnesium atoms attached to helium drops has unraveled an interesting and somewhat unexpected solvation behavior as a function of the number ͑N͒ of helium atoms in the drop.Diffusion Monte Carlo ͑DMC͒ calculations 1 carried out for small drops containing a number of helium atoms up to N = 50 indicate that a Mg atom is not fully solvated in drops of sizes below N ϳ 30.More recent quantum Monte Carlo ͑MC͒ calculations 2,3 suggest a surface Mg state for 4 He clusters with up to ϳ200 atoms.Density functional theory ͑DFT͒ calculations 4 for 3 He N and 4 He N nanodroplets with N Ն 300 doped with alkaline-earth atoms have shown that Mg atoms are solvated in their interior, in agreement with the analysis of laser-induced fluorescence ͑LIF͒ ͑Ref.5͒ and resonant two-photon-ionization ͑R2PI͒ experiments. 6LIF experiments on the absorption and emission spectra of Mg atoms in liquid 3 He and 4 He have been reported and successfully analyzed within a vibrating atomic bubble model, where full solvation of the impurity atom is assumed. 7,8A more recent experiment, 9 in which electronimpact ionization data from Mg-doped 4 He drops with about 10 4 atoms seem to indicate that magnesium is instead at the surface of the droplet, is in disagreement with the abovementioned LIF and R2PI experiments.
There is some ambiguity associated with the notion of solvation in a helium droplet.For not too small droplets, one may consider that Mg is fully solvated when its position inside the droplet is such that its solvation energy or atomic shift does not appreciably differ from its asymptotic value in bulk liquid helium, as both quantities approach such limit fairly alongside each other.However, for very small drops the energy or atomic shifts of an impurity atom at the center of the drop may still differ appreciably from the bulk liquid values because there is not enough helium to saturate these quantities.This is the case, e.g., of Mg@ 4 He 50 studied in Ref. 1.
Mella et al. 1 discussed how the solvation properties of magnesium are affected by the number of helium atoms in small 4 He drops.Since DMC calculations cannot be extended to the very large drops involved in LIF experiments, 5 they could not carry out a detailed comparison between their calculated atomic shifts and the experiments.They also pointed out the sensitivity of the Mg solvation properties on apparently small differences in the He-Mg pair potentials available in the literature.
The aim of the present work is to re-examine the solvation of magnesium in 4 He nanodroplets from the DFT perspective, [10][11][12][13][14] extending our previous calculations 4 down to drops with N ϳ 50 atoms and improving the DFT approach 15 by fully taking into account the coupling of the dipole excitation of the impurity with the dynamical deformations of the helium around the Mg atom.Our results confirm that full solvation of a Mg atom in 4 He drops requires a minimum number of helium atoms, and disclose some unusual results for small drops.We calculate the absorption spectrum of a Mg atom attached to small and large drops, finding good agreement with the experiments for the latter.We discuss in a qualitative way the effect of the impurity angular momentum on the electron-impact ionization yield and on the absorption spectrum.We also address the struc-ture of a two-magnesium-doped drop; the results are used to discuss the scenario proposed by Przystawik et al. 6 to interpret their experimental results on R2PI.The experimental status of doped helium drops is reviewed in, e.g., Refs.16-19.This work is organized as follows.In Sec.II we discuss the results we have obtained for the structure of Mg@ 4 He N drops.The method we have employed to obtain the atomic shifts is presented in Sec.III, and it is applied to the case of Mg-doped 4 He droplets in Sec.IV.In Sec.V we study how the presence of a second magnesium atom in the droplet may alter both the helium drop structure and the calculated atomic shift.A summary is presented in Sec.VI, and some technical details of our calculations are described in Appendixes A and C.

II. STRUCTURE AND ENERGETICS OF Mg-DOPED HELIUM NANODROPLETS WITHIN DENSITY FUNCTIONAL THEORY
In this section we recall the essentials of our method and refer the interested reader to Ref. 15 for the details.Our starting point is the Orsay-Trento density functional, 10 together with the Mg-He adiabatic ground-state potential X 1 ⌺ of Ref. 20, here denoted as V Mg-He .To check the sensitivity of our results to the details of different available pair potentials describing the Mg-He interaction, we also use the slightly less attractive potential computed in Ref. 21.For the sake of comparison, we plot both potentials in Fig. 1.Despite the apparently minor differences between these two potential curves, they cause very different solvation properties of Mg in small 4 He drops, as we will show in the following.
The energy of the Mg-helium system E͓⌿ , ⌽͔ is written as a functional of the Mg wave function ⌽͑r͒ and the 4 He order parameter ⌿͑r͒ = ͱ ͑r͒, where ͑r͒ is the 4 He atomic density.Minimizing E under the constraints of a given N and a normalized Mg wave function, one obtains an Euler-Lagrange equation for ⌿ coupled to a Schrödinger equation for ⌽ that we have self-consistently solved as outlined in Ref. 15 in order to determine the ground state ͑gs͒ of the drop-impurity complex as well as the helium chemical potential and lowest Mg eigenenergy .
To address the solvation of the Mg atom, we have found it convenient to minimize E subjected to the additional constraint of a fixed distance Z 0 between the centers of mass of the helium moiety and of the impurity atom, which, due to the symmetry of the problem, can both be taken on the z axis.This is done following a method borrowed from nuclear physics-similar to that used to describe the fission of rotating 3 He drops. 22We minimize the expression E + C ͓Z − Z 0 ͔ 2 / 2, where Z is the distance between the impurity and the center of mass of the helium droplet.C is an arbitrary constant.The value C ϳ 1000 K Å −2 has been used in our calculations, which ensures that the desired Z 0 value is obtained within a 0.1% accuracy.
We have first solved the mentioned coupled equations for C = 0 and several N values, namely N = 30, 50, 100, 200, 300, 500, 1000, and 2000.This yields the ground state of the drop-impurity complex, and will allow us to study the atomic shift for selected cluster sizes.Figure 2 shows the energy of a magnesium atom in a drop, defined as We compare in Fig. 2 the S N ͑Mg͒ values here obtained with those of Ref. 4, where the Mg atom was treated as an infinitely massive particle-i.e., as a fixed external field acting on the 4  The DMC energy found for the same system is ϳ −168.2K, 1 whereas our DFT result is ϳ−157.0K.The "asymptotic" DMC value for S 200 ͑Mg͒, −33.1 K, compares well with the DFT value, still far from the limit value corresponding to a very large helium drop ͑see Fig. 2͒.The solvation properties of the Mg atom are determined by a delicate balance between the different energy terms- surface, bulk, and helium-impurity-whose contribution is hard to disentangle and depends on the number of atoms in the drop, as shown by DMC and DFT calculations.To gain more insight into the solvation process of magnesium in small 4 He droplets, we have computed the energy of the doped droplet as a function of the impurity position Z 0 .
The bottom panel in Fig. 3 shows E͑Z 0 ͒ for N = 50, as computed using the two different Mg-He pair interactions shown in Fig. 1 ͑the discussion of the top panel is postponed until Sec.IV A͒.It can be seen that E͑Z 0 ͒ in both cases displays two local minima.In one case ͑i.e., for the slightly less attractive potential in Fig. 1͒ a "surface" state for the Mg atom is energetically preferred, while in the other case the impurity prefers to sit in the interior to the droplet ͑although not exactly at its center͒.These results are in agreement with the DMC calculations in Ref. 1, where it was shown that a bimodal distributions for the Mg radial probability density function with respect to the center of mass of the helium moiety appears for N Յ 30. 1 More recent DMC calculations carried out up to N ϳ 200 drops 3 seem to point out that Mg is always in a surface state, although somewhat beneath the drop surface.Our DFT calculations yield that Mg is already solvated for N = 200.
For both pair potentials, the bottom panel of Fig. 3 shows that the two local minima in E͑Z 0 ͒ are separated by an energy barrier of about 1 K height, allowing the impurity to temporarily visit, even at the experimental temperature T ϳ 0.4 K, [16][17][18][19] the less energetically favored site.This causes changes in the total energy of the system by less than 1%, but has a large effect on the value of the atomic shift.We will address this important issue and its consequences on the computed spectral properties of Mg@ He N in Secs.IV and V.
Figure 4 shows the helium configurations for the four stationary points displayed in Fig. 3, namely, those corresponding to Z 0 = 0, 2, 4, and 6 Å.Although the preference for the surface or the solvated state depends on the He-Mg pair potential used in the calculations, similar to the case of other alkali atoms, 4 this does not seem to be the case for the stationary points at Z 0 = 2 and 4 Å that are present in both curves.We have compared in Fig. 5 the density profiles along the z axis for the Z 0 =2͑4͒ Å configuration with the profile of the pure 4 He drop, finding that the appearance of these stationary points is related to the position of the density peak in the first helium solvation shell with respect to a maximum ͑minimum͒ of the density of the pure drop.We see that a minimum ͑maximum͒ in the energy is associated with a constructive ͑destructive͒ interference in the oscillation pattern of the He density.We are lead to conclude that the interplay between the density oscillations already present in pure drops and the solvation shells generated by the impurity plays an important role in the solvation properties of the Mg atom.This effect is also present, although to a lesser extent, in Ca-doped 4 He nanodroplets. 15ventually, for larger drops Mg becomes fully solvated.We have found that this is the case whichever of these two  ͑Color online͒ Helium density plots of the Mg@ 4 He 50 droplet on the y = 0 plane obtained using the Mg-He potential of Ref. 20.From top to bottom and left to right, the Z 0 values correspond to the stationary points displayed in Fig. 3, namely, 0, 2, 4, and 6 Å.The brighter regions are the higher-density ones.potentials we use ͑see, for instance, the bottom panel of Fig. 6͒.For this reason, the results we discuss in the following have been obtained with the Mg-He pair potential of Ref. 20, unless differently stated.
When the Mg atom is fully solvated, e.g., for N = 1000, we have found that E͑Z 0 ͒ grows monotonously as Z 0 increases ͑i.e., as the Mg atom approaches the droplet surface͒.This is shown in the bottom panel of Fig. 6.To better understand how E͑Z 0 ͒ depends on N, we have plotted in the top panel of Fig. 6 the energies of the N = 1000 and 2000 doped drops, referred to their equilibrium values, as functions of the distance from the dividing surface, i.e., the radius R 1/2 at which the density of the pure drop equals b / 2, b being the liquid density value ͑R 1/2 = r 0 N 1/3 , with r 0 = 2.22 Å͒.These radii are 22.2 and 28.0 Å for N = 1000 and 2000, respectively.
The top panel of Fig. 6 shows that most of the change in ⌬E takes place in the 15 Å outer region of the drop, irrespective of its size.As a consequence of the flatness of E͑Z 0 ͒, the magnesium atom is very delocalized in the radial direction even in relatively small drops.This delocalization might affect the absorption spectrum of the attached Mg atom, and also must be explicitly considered in the interpretation of the electron-impact yield experiments. 9 We address this issue in Sec.IV.

III. EXCITATION SPECTRUM OF A Mg ATOM ATTACHED TO A 4 He DROP
To determine the absorption spectrum of a Mg atom embedded in a helium drop, we have used the method of Lax, 23 together with the diatomics-in-molecules approach. 24As in Sec.II, we recall here the essentials of the method and refer again the interested reader to Ref. 15 for the details.
To obtain the line shape, one has to carry out an average on the possible initial states of the system that may be thermally populated.This is likely not needed for helium drops, as their temperature is much smaller than the vibrational excitation energies of the Mg atom in the mean field created by the helium environment.In small helium droplets thermal effects can show up in the Mg absorption spectrum due to the high mobility of the atom.For large drops, thermal motion plays a minor role, as the Mg atom hardly gets close enough to the drop surface to have some effect on the line shape.In this case, however, dynamical deformations of the  cavity around the impurity may be relevant 25,26 and these degrees of freedom have to be taken into account.
The line shape I͑͒ for electronic transitions from the gs to the excited state ͑ex͒ in a condensed phase system is evaluated using the Born-Oppenheimer approximation and the Franck-Condon principle.If ͉⌿ gs ͘ is the ground state of the system, projecting on eigenstates of the orbital angular momentum of the excited electron ͉m͘, one obtains where ប X gs and X gs ͑r , ͓␣͔͒ are the energy and the wave function of the rovibrational ground state of the frozen heliumimpurity system, and H m ex ͑r , ͓␣͔͒ is the rovibrational excited Hamiltonian with potential energy V m ex ͑r , ͓␣͔͒ determined by the electronic energy eigenvalue, as obtained for a p ← s transition.At this point, we have introduced the variables ͓␣͔ to represent the degrees of freedom needed to describe possible deformations of the system, corresponding to the zeropoint oscillations of the helium bubble around the impurity.If this effect is neglected, the deformation parameters ͓␣͔ are dropped and the ground-state wave function X gs coincides with the Mg wave function ⌽ found by solving the equation of Schrödinger we referred to in Sec.II.
Treating as classical the relevant excited states for the transition, which is justified if they have large quantum numbers, one obtains the expression where ⍀ m ͑͒ is the surface defined by the equation + X gs − V m ex ͑r , ͓␣͔͒ / ប = 0.If the atom is in bulk liquid helium or at the center of the drop, the problem has spherical symmetry and the above equation can be further simplified. 15he potential-energy surfaces ͑PESs͒ V m ex ͑r , ͓␣͔͒ needed to carry out the calculation of the atomic shifts have been obtained in the pairwise sum approximation, using the V ⌸ ͑r͒ and V ⌺ ͑r͒ Mg-He adiabatic potentials from Ref. 1.In Cartesian coordinates and assuming that the He-impurity spinorbit interaction is negligible for magnesium, the eigenvalues of the symmetric matrix are the V m ex ͑r , ͓␣͔͒ potentials which define the PESs as a function of the distance between the center of mass of the droplet and that of the impurity.

͑5͒
In this case, two of the PESs are degenerate, namely, 1 ͑r , ͓␣͔͒ = 2 ͑r , ͓␣͔͒ 3 ͑r , ͓␣͔͒. 27This holds true for r 0, and it is relevant when we take into account the delocalization of the impurity inside the bubble due to its quantum motion.Otherwise, since at r = 0 all the i coincide, they are threefold degenerate.

IV. RESULTS FOR THE ABSORPTION SPECTRUM OF MAGNESIUM ATOMS
We consider first the ͓␣͔ = 0 case ͑i.e., no zero-point deformations of the helium cavity hosting the impurity͒.The general situation, in particular the homogeneous width calculation, is presented later on in this section.
If ͓␣͔ = 0, the model is expected to yield at most the energies of the atomic transitions, but not the line shapes since the impurity-droplet excitation interactions as well as inho-mogeneous broadening effects are not taken into account. 28,29We discuss here two illustrative examples, namely, the atomic shifts of magnesium in N = 50 and 1000 nanodroplets.The homogeneous width is calculated in Sec.IV B for the N = 1000 droplet.
For Mg@ 4 He 50 , the calculated shifts at Z 0 =0 ͑spherical configuration͒, 2, and 6 Å are 500, 450, and 281 cm −1 , respectively ͑281.3, 281.7, and 283.0 nm wavelengths͒.No experimental information is available for such a small drop.Contrarily, the absorption spectrum of the 3s3p 1 P 1 ← 3s 2 1 S 0 transition of Mg atoms attached to large helium drops has been measured, 5 displaying a broad peak strongly blueshifted from its position in the gas phase.This spectrum is remarkably close to the one obtained by Moriwaki and Morita 7 in bulk liquid helium, and hence it has been concluded that Mg is in the interior of the 4 He droplet.We want to point out that, while the absorption line in liquid helium was attributed to a single broad peak of energy of 281.5 nm ͑35 524 cm −1 ͒, i.e., a shift of 474 cm −1 , 7 in large drops a similar line profile 30 was fitted by two Gaussians centered at 35 358 and 35 693 cm −1 ͑282.8 and 280.2 nm wavelengths͒, respectively, i.e., shifted 307 and 642 cm −1 from the gasphase line. 5The origin of the two peaks was attributed to the splitting of the degenerate ⌸ state by dynamical quadrupole deformations of the cavity surrounding the dopant since this argument had qualitatively explained similar doubly shaped D 2 excitation spectra of Rb and Cs atoms in liquid 4 He due to a quadrupole oscillation of the helium bubble ͑dynamic Jahn-Teller effect͒. 26IF experiments on the heavier alkaline earths Ca and Sr in large 4 He droplets 31 have disclosed the existence of strongly blueshifted broad peaks with no apparent structure, although it cannot be discarded that this broad line could be a superposition of unresolved peaks.The same happens for Ba. 32The surface location of Ca, Sr, and Ba in these drops has been further confirmed by DFT calculations. 4,15It is also interesting to recall that LIF experiments on Ca atoms in liquid 4 He and 3 He have found a broad line in the region of the 4s4p 1 P 1 ← 4s 21 S 0 transition with no apparent splitting, 33 contrary to the case of Mg.
Since Mg is fully solvated in the N = 1000 drop, the calculated atomic shift ⌬ may be sensibly compared with the experimental data where drops with N in the 10 3 -10 4 range are studied.We have obtained ⌬ = 659.0cm −1 ͑280.0 nm wavelength͒; this peak nearly corresponds to the Gaussian that takes most of the intensity of the absorption line ͑about 87%͒. 5We have carried out a detailed analysis for this drop, determining the equilibrium structure of Mg@ 4 He 1000 as a function of Z 0 , and have used it to evaluate ⌬.The results are displayed in Table I, showing the actual sensitivity of the absorption spectrum to the Mg-atom environment.
The impurity-drop excitations will determine the homogeneous width of the spectral line, and the population of excited states may be relevant given the limit temperature attained by the droplets. 17,18,35In this context, the relevant excitation modes of the helium bubble are radial oscillations of monopole type ͑breathing modes͒ and multipole shape oscillations about the equilibrium configuration, as well as displacements of the helium bubble inside the droplet.We will address these issues in Secs.IV and IV B.

A. Thermal motion and angular momentum effects
To describe the displacement of the helium bubble inside the droplet, we have fitted the E͑Z 0 ͒ curve of the Mg@ 4 He 1000 system to a parabola, and have obtained the excitation energy ប for this three-dimensional ͑3D͒ isotropic harmonic oscillator.The hydrodynamic mass of the impurity atom has been estimated by its bulk liquid helium value, M ‫ء‬ ϳ 40 a.u., obtained by the method outlined in Appendix A, Eq. ͑A23͒.We find ប = 0.1 K, indicating that thermal motion, i.e., the population of the excited states, of the "mean field" E͑Z 0 ͒ is important at the experimental temperature T = 0.4 K and may produce observable effects in the absorption spectrum and the electron-impact ionization yield.
To describe in more detail the delocalization of the Mg atom inside the drop, we have used an effective Hamiltonian where we interpret Z 0 as the radial distance R between the impurity and center of mass of the helium moiety, and E͑Z 0 ͒ as the "potential energy" V͑R͒ associated with this new degree of freedom of the impurity in the drop.Namely, where M ‫ء‬ is the Mg hydrodynamic mass.In the canonical ensemble, the total probability distribution W as a function of R can be written as where the partition function is defined as The radial probability density w is This expression has been evaluated for N = 50 and 1000 by solving the Schrödinger equation for Hamiltonian ͑6͒ to obtain the orbitals nᐉ and eigenenergies E nᐉ .For larger N values, we have used the semiclassical approximation where the effective potential is Integrating p R in phase space, we obtain for the probability density

͑10͒
with the normalization Lacking a better choice, we have weighted any possible angular momentum value with a Boltzmann energy factor.It has been shown 36 that some of the angular momentum deposited in the droplet during the pick-up process may be kept in the impurity atom, resulting in a different angular momentum distribution from the Boltzmann one.This could yield that some Mg atoms are actually closer to the drop surface.
If the angular momentum associated with the motion of the magnesium atom-whose "radius" is ϳ5 Å ͑see Fig. 1͒-is such that Mg can be some 10 Å beneath the drop surface, the shift of the absorption line would be hardly distinguishable from that of the totally solvated case-as seen in Table I.At the same time, the electron-energy dependence of the Mg + yield observed in electron-impact ionization experiments 9 ͑and which was considered as evidence of a surface location of Mg atoms on 4 He droplets͒ could indeed be due to the Penning ionization of the impurity in a collision with a metastable He ‫ء‬ atom that occupies a surface bubble state in the drop, instead of being due to the transfer of a positive hole ͑He + ͒ to the Mg atom, which is the primary ionization mechanism when the impurity is very attractive and resides in the deep bulk of the droplet.
We show in the top panel of Fig. 3 the probability densities w͑R͒ at T = 0.4 K corresponding to the configurations displayed in the bottom panel. 37Similarly, the top panel of Fig. 6 shows that for Mg@ 4 He 1000 and Mg@ 4 He 2000 , if thermal motion is taken into account and the impurity retains some of the pick-up angular momentum, the maximum density probability of Mg is at ϳ15 Å beneath the drop surface in both cases.To obtain it, we have taken for M ‫ء‬ the bulk value of 40 a.u.As seen from Table I, the absorption line shift changes by a small 2% with respect to the R = 0 configuration.The values of the angular momentum corresponding to these maximum density probabilities are ͗L 2 ͘ 1/2 ϳ 9ប and ϳ10ប, respectively. 38For an N = 10 000 drop, 9 whose radius is R 1/2 = 47.8Å, we have extrapolated inward the ⌬E curves of the calculated N = 1000 and 2000 drops, and have obtained from it the probability distribution displayed in Fig. 6.Its maximum is at ϳ18.5 Å beneath the surface, with ͗L 2 ͘ 1/2 ϳ 28ប.
Finally, using the effective potential in Eq. ͑9͒, we have determined that for N = 1000 ͑2000͒, a Mg atom with ᐉ ϳ 24ប ͑44ប͒ is in an "equilibrium position" some 10 Å beneath the drop surface.For the N = 10 000 drop, this value is ᐉ ϳ 148ប.These values look reasonable, and Mg atoms holding this angular momentum or larger might thus be the origin of the primary electron-collision ionization yield by the Penning process, without questioning the conclusion drawn from LIF experiments that magnesium is fully solvated in 4 He drops.

B. Homogeneous width from shape deformations of the helium bubble
We have shown that for large drops, the Mg atom is fully solvated and its thermal motion only produces small changes in the absorption shift.This allows us to decouple the effect of the translational motion of the helium bubble on the absorption line from that of its shape fluctuation.Moreover, we can address shape fluctuations in the much simpler spherically symmetric ground state, when magnesium is located at the center of the drop.
To quantify the effect of these fluctuations, we have first used the spherical cap model 39 to estimate the excitation energies of the helium bubble around the impurity in liquid helium.To this end, we have fitted the V Mg-He potential to a Lennard-Jones potential with depth = 7 K at a minimum distance r min = 5 Å.Minimizing the total energy within this model yields a configuration with an equilibrium radius of R 0 = 0.97ϫ 2 −1/6 r min , which we approximate by R 0 =2 −1/6 r min to obtain the excitation energies.
Deformations of the 4 He around the Mg atom are modeled as 40,41 where R is the radius of the bubble cavity hosting the solvated Mg atom, ⍀ represents the solid angle variables ͑ , ͒, Z ͑⍀͒ is a real spherical harmonic, and ␣ is the amplitude of the -multipole deformation.The condition R ␣ = R 0 ͑1 − 1 4 ͚ ͉␣ ͉ 2 ͒ ensures the conservation of the number of particles up to second order in ␣ .The dipole mode amplitude ␣ 1 is absent since, for an incompressible fluid, it corresponds to a translation of the bubble, and this has been considered in Sec.IV A.
If S͓R͑⍀͔͒ represents the bubble surface and represents the surface tension, the energy for a large drop can be written as where C 0 =8R 0 2 ͑1+12 A ͒ and C = R 0 2 ͓͑ −1͒͑ +2͒ +6 A ͔ are the stiffness parameters, and A = b 2 −1/6 r min −1 is the impurity-He solvation parameter. 39The mass parameters are B 0 =4 b m He R 0 5 and B = b m He R 0 5 / ͑ +1͒, 42 and the excitation energies are determined from ប = ប ͱ C / B , yielding ប 0 = 10.2K for the breathing mode and ប 2 = 9.8 K for the quadrupole mode.Given the droplet temperature of 0.4 K, we conclude that only the ground state is populated.The mean amplitude of the shape oscillations is estimated from the variance ␥ = ប 1/2 ͑B C ͒ −1/4 / 2, giving ␥ 0 = 0.03 and ␥ 2 = 0.15.This model thus yields that the bubble can experience monopole oscillations of ϳ3% amplitude and quadrupole deformations of ␥ 2 ͱ 3 / 4 ϳ 8% ampli- tude.Amplitudes of this order have been determined within the atomic bubble model for cesium atoms in liquid helium. 26Since their effect in the absorption spectrum is expected to be relevant, we have undertaken a more refined calculation within DFT taking Mg@ 4 He 1000 as a case study.
For helium droplets, we have described bubble deformations in a way similar to that in Refs.26 and 40-42, namely, if 0 ͑r͒ is the helium spherical ground-state density, deformations are introduced as ͑r , t͒ = 0 ͓R͑r , t͔͒K −1 , with where the real spherical harmonics are normalized as 2+1 ␦ Ј ␦ Ј for convenience, and the normalization K = N −1 ͐d 3 r 0 ͓R͑r , t͔͒ ensures particle number conservation.If the Mg wave function follows adiabatically the helium density deformation, it can be shown that to second order in ␣ , the total energy of the system can be written as

͑14͒
where E gs is the ground-state energy, M ‫ء‬ is the hydrodynamic mass associated with the mode, and E ͑2͒ is the second derivative of the total energy with respect to ␣ .This equation represents the Hamiltonian of a set of uncoupled harmonic oscillators, whose quantization yields a ground state to whose energy each mode contributes with ‫ء‬ , and a ground-state wave function

͑15͒
where ␣ 2 ϵ͚ =− ␣ 2 .Details are given in Appendix A. We have computed the hydrodynamic masses assuming that the drop is large enough to use Eq.͑A23͒.We have obtained M 0 ‫ء‬ = 15.0mHe + 0.28m Mg ϳ 66.7 a.u. and M 2 ‫ء‬ = 1.9mHe + 0.56m Mg ϳ 21.0 a.u.In actual calculations, instead of using Eq.͑A10͒, the energies E 0 ͑2͒ and E 2 ͑2͒ have been numerically obtained by computing the total energy of the system for different small values of ␣ 0 and ␣ 2 .This has been carried out by numerically introducing the desired deforma-tion parameter into the ground-state density and renormalizing it, solving next the Schrödinger equation for the Mg atom to determine the ground state of the impurity, and computing the total energy of the system.Fitting these curves to a parabola, we have obtained E 0 ͑2͒ = 49.7 K Å −2 and E 2 ͑2͒ = 16.8K Å −2 .We have then calculated the ground-state energies ប / 2 and deformation mean amplitudes ␥ , obtaining ប 0 / 2 = 10.6 K and ប 2 / 2 = 6.3 K, with mean amplitudes ␥ 0 = 0.18 Å ͑ϳ3.7%͒ and ␥ 2 = 0.42 Å ͑ϳ8.5%͒.
To quantitatively determine the effect of these deformations on the absorption spectrum, we have developed Eq.͑3͒ to first order in the deformation parameters, and have explicitly shown that to this order, only the breathing and quadrupole modes affect the dipole absorption spectrum.The details are given in Appendix B, where we show that the breathing mode affects the shift and shape of the line, whereas quadrupole modes only affect the shift.
Consequently, we restrict in Eq. ͑3͒ the deformation parameters needed to properly describe the homogeneous broadening of the absorption dipole line, namely, d͓␣͔ → d␣ 0 d 5 ␣ 2 and X gs ͑r , ͓␣͔͒ → ͑␣ 0 , ␣ 2 ͒⌽͑r , ␣ 0 ͒, ⌽͑r , ␣ 0 ͒ being the wave function of the Mg atom for a given ␣ 0 value.We compute the spectrum as where ex ͑r , ␣ 0 , ͓␣ 2 ͔͒ are the eigenvalues m ͑r , ͓␣͔͒ of the excited potential matrix in Eq. ͑C4͒, and r m ͑͒ is the root of the equation + X gs Expression ͑16͒ has been integrated using a Monte Carlo method.We have sorted M =10 6 sets of values i ͖ using the square of the wave function in Eq. ͑15͒ as probability density.Next, for each set we have found the eigenvalues V m ex ͑r , ␣ 0 i , ͓␣ 2 ͔ i ͒ of the U i,j matrix that define the potential-energy surfaces and have used a trapezoidal rule to evaluate the integrals using a discretized representation of the delta function. 15Finally, we have obtained the spectrum as Figure 7 shows the absorption spectrum of one Mg atom attached to 4 He 1000 in the vicinity of the 3s3p 1 P 1 ← 3s 2 1 S 0 transition when homogeneous broadening is considered.We have decomposed the absorption line into its three components, the higher-frequency component being the ⌺ one.The starred vertical line represents the gas-phase transition, and the experimental curve, adapted from Ref. 5, has been vertically offset for clarity.Also shown is the absorption spectrum obtained by neglecting homogeneous broadening ͑hatched region͒.This figure shows that the both the energy and width of the absorption peak that takes most of the experimental intensity are correctly described by our calculations.

V. TWO MAGNESIUM ATOMS ATTACHED TO A 4 He DROP
The attachment of magnesium atoms in 4 He droplets has been recently addressed using resonant two-photon ionization. 6In particular, the authors of Ref. 6 obtained the absorption spectrum for drops doped with different selected numbers of Mg atoms.From their measurements it appears that two main features contribute to the observed line shapes, one peaked at about 279 nm and another at about 282 nm.This is in agreement with the results of Refs. 5 and 7 ͑we recall that, actually, the two peaks were not resolved by the authors of the bulk liquid experiment͒.The structure at 282 nm, however, only appears if the droplet contains more than one Mg atom.Thus the two-peak structure cannot be due to the splitting of the absorption line due to dynamical quadrupole deformations of the helium bubble around the impurity, as previously believed.We have indeed shown in Sec.IV that this coupling only produces a broad peak, in good agreement with the results of Ref. 6 for helium drops containing just one Mg atom.
Another interesting observation reported in Ref. 6 is that their experimental results for multiatom-doped 4 He droplets are not consistent with the formation of compact metallic Mg clusters inside the 4 He droplet.The magnesium atoms in the droplet appear instead to be relatively isolated from each other, showing only a weak interaction and leading to the 282 nm shift in the observation.
To confirm this scenario and find an explanation for the origin of the low-energy component in the absorption peak, we have carried out DFT calculations to determine the structure of a two-magnesium-doped 4 He drop.Our goal is to verify whether the helium density oscillation around a magnesium atom may result in an energy barrier preventing the Mg atoms from merging into a Mg 2 dimer, as suggested by Przystawik et al. 6 To obtain the structure of two Mg atoms in a 4 He drop, we have minimized the energy of the system written as

͑19͒
where V Mg-Mg ͉͑r − rЈ͉͒ is the Mg-Mg pair potential in Ref. 43, ⌽ i ͑r͒ is the wave function of the ith Mg atom and ⌿͑r͒ = ͱ ͑r͒ represents again the 4 He order parameter, where ͑r͒ is the 4 He atomic density.In this expression, E͑͒ is the 4 He potential-energy density. 10here are at least two additional effects which are not considered when modeling the Mg-Mg interaction via the pair potential in vacuum, as implied in the above expression.The first is due to three-body ͑and higher͒ correlation effects involving the 4 He atoms surrounding the Mg pair: these should exert an additional, albeit small, screening effect due to He polarization, which is expected to reduce the absolute value of the dispersion coefficients in the long-range part of the Mg-Mg interaction.The second is a possible reduction in  4 He 1000 in the vicinity of the 3s3p 1 P 1 ← 3s 2 1 S 0 transition.The line has been decomposed into its two ⌸ and one ⌺ components; the former is the higher-frequency transition.The starred vertical line represents the gas-phase transition.The experimental curve, adapted for Ref. 5, has been vertically offset for clarity.Also shown is the absorption spectrum obtained neglecting homogeneous broadening ͑hatched region͒.
the Mg-atom polarizability due to the presence of the surrounding 4 He cavity which, because of the repulsive character of the electron-He interaction, should make the electronic distribution of the impurity atom slightly "stiffer," thus reducing further the values of the dispersion coefficients in the Mg-Mg pair interactions.Although in principle these effects might reduce the net interaction between a pair of Mg atoms embedded in liquid 4 He, in practice in the present system they are indeed very small.The correction to the leading term of the long-range dispersion interaction, −C 6 / r 6 , due to three-body correlation effects can be written to first order 44 as −C 6 ͑1−2n␣ / 3͒ / r 6 , ␣ being the static polarizability of the host fluid ͑␣ He = 1.39a 0 3 ͒.Such correction is on the order of only 1% in our case.To estimate the change in the Mg atomic polarizability due to the surrounding He, we computed, using ab initio pseudopotential calculations, the ͑static͒ polarizability of a Mg atom in the presence of an effective ͑mainly repulsive͒ potential acting on the Mg valence electrons due to the presence of the surrounding He.The effective interaction is derived from the equilibrium shape of the 4 He bubble hosting the Mg atom, as predicted by our DFT calculations, and assuming a ͑local͒ electron-He density-dependent interaction which was proposed by Cheng et al. 45 We find a very small change in the static atomic polarizability ␣ of Mg.Assuming that, roughly, C 6 ϰ ␣ 2 , we find a reduction in the C 6 coefficient of about 1 -2 %.
The minimization of the total-energy functional written above under the constraint of a given number of helium atoms and normalized Mg ground-state wave functions should in principle yield the equilibrium configuration of the system.In practice, depending on the initial configuration, we have found several local minima, whose origin is again the "interference" of the He solvation shells around the Mg atoms.We have found three such metastable configurations for ͑Mg+ Mg͒ @ 4 He 1000 if we start the minimization procedure with one Mg atom near the center of the droplet, and the other placed off center, at some distance from the first.They are displayed in Fig. 8.The energy difference between the innermost ͑Mg-Mg distance d = 9.3 Å͒ and the outermost ͑d = 18.5 Å͒ configurations is 12.5 K.The energy of the d = 9.3 Å configuration is sensibly that of the configuration specularly symmetric about the z = 0 plane ͑−5581.4K͒ also shown in the figure.It is worth noticing that, since R 1/2 = 22.2 Å and the radius of Mg is ϳ5 Å, only the upper left corner configuration has the Mg impurity in a surface state.
Figure 9 shows two density profiles of the symmetric configuration obtained along the z axis ͑solid line͒ and the x or y axis ͑dashed line͒.It shows a relatively high-density helium ring around the two-bubble waist, clearly visible in Fig. 8, where the local density is almost three times the bulk liquid 4 He density, and that prevents the collapse of the two-bubble configuration.The Mg wave functions are peaked at ϳ Ϯ 4.75 Å and very narrow.This justifies a posteriori the assumptions we have made to write the total energy of the system, Eq. ͑19͒.
Our results confirm the existence of the energy barrier suggested in Ref. 6, which prevents the two Mg atoms from coming closer than some 9 Å, and thus hinders, at least temporarily, the formation of the Mg 2 dimer.This barrier is shown in Fig. 10 as a function of the Mg-Mg distance d, which is kept fixed in a constrained energy minimization.Note that the energy of the Mg+ Mg system increases as d does because the two Mg atoms in a drop form a state more bound than that of a pure drop with the impurities well apart.
Since the height of the barrier is larger than the experimental temperatures, two solvated Mg atoms will not easily merge into a dimer, but rather form a metastable weakly bound state.Based on these findings, we suggest that several Mg atoms solvated inside 4 He drops might form a sparse, weakly bound "foamlike" aggregate rather than coalesce into a more tightly bound metallic cluster.Partial coagulation of impurities was already invoked by Lewerenz et al. 46 to ex- plain their experimental findings for the successive capture of foreign atoms and molecules in helium clusters ͑see also Ref. 47 for the case of bulk liquid helium͒.Very recently, a kind of "quantum gel" has been predicted to be formed in 4 He drops doped with neon atoms. 48Although some degree of mutual isolation between foreign atoms is expected in the case of strongly attractive impurities ͑such as those studied in the two cases mentioned above͒, where they are kept apart by the presence a solid 4 He layer coating the impurity, 47 our calculations show that this effect is possible even for relatively weakly attractive impurities such as Mg, where such solidlike 4 He layer is absent.
One may estimate the mean life of the metastable state as

‫ء‬
UЉ͑d eq ͒ exp͓⌬U/͑k B T͔͒, ͑20͒ where Mg ‫ء‬ = M Mg ‫ء‬ / 2 ϳ 20 a.u. is the hydrodynamic reduced mass of the Mg+ Mg system and ⌬U is the barrier height.From Fig. 10 we have that UЉ͑d eq ͒ϳ40 K Å −2 .This yields a mean life of a few nanoseconds, which is about 5-6 orders of magnitude smaller than the time needed for its experimental detection. 6,49The mean life becomes increasingly large as the relative angular momentum L deposited into the two Mg system increases.Writing one obtains the L-dependent energy barriers displayed in Fig. 10.For L = 30 we have ϳ 0.6 s, and for L = 40, ϳ 0.1 ms.Thus, there is an angular momentum window that may yield mean lives compatible with the experimental findings.Increasing L much further would produce too a distant Mg+ Mg system which would correspond to two independent Mg impurities in a drop.
To check whether this foamlike structure of the Mg aggregate also appears in 3 He drops, we have carried out calculations for ͑Mg+ Mg͒ @ 3 He 1000 using the same density functional as in Ref. 50, and the method presented in this section.Figure 11 shows the energy of the ͑Mg +Mg͒ @ 3 He 1000 complex as a function of d.For distances smaller than some 8.6 Å, we have found that the system has a tendency to collapse into a dimer-physically unreachable from our starting point, Eq. ͑19͒.We are led to conclude that there is no barrier in the case of liquid 3 He.The configuration corresponding to the closest d we have calculated is shown in Fig. 12.
We are now in the position to determine the effect of these weakly bound systems on the LIF and R2PI experiments on 4 He droplets containing more than one Mg atom.Notice that the bottom right panel of Fig. 8 shows that the helium bubbles have a nonzero static quadrupole moment, whereas they are spherically symmetric for a single Mg atom in the drop.It is precisely the existence of this static quadrupole moment that causes an additional separation between the ⌺ and ⌸ spectral components in the absorption spectra, which results, as a consequence of the broadening of each line, in a double-peak structure of the computed spectra, in semiquantitative agreement with the experimental data.Details of our calculation are given in Appendix C; see Eq. ͑C6͒. Figure 13 shows such two-peak structure corresponding to the specularly symmetric configuration displayed in Fig. 8, and indicates that the 282 nm structure observed in the experiments may be attributed to the distortion produced by neighbor Mg bubbles.These bubbles contribute incoherently to the absorption spectrum, and the relative intensities of the 282 and 279 nm peaks might reflect the different populations of drops doped with one and two Mg atoms since those hosting a compact cluster ͑dimer, trimer, etc.͒ would not yield an absorption signal in the neighborhood of the monomer 3s3p 1 P 1 ← 3s 2 1 S 0 transition.Notice that a large static distor- .͑Color online͒ Energy ͑K͒ of the Mg+ Mg@ 4 He 1000 system as a function of the Mg-Mg distance ͑Å͒.The energies have been referred to that of the metastable equilibrium configuration ͑local minimum͒ at L = 0.The lines have been obtained by a cubic spline of the actual calculations.From bottom to top, the curves correspond to L =0-50ប in 10ប steps.tion of the helium bubble could also arise if the Mg atom were in a shallow dimple at the drop surface, but Fig. 6 discards this possibility.
We finally note that we have not considered in our work another source for an additional splitting of the spectral lines of a Mg atom in the field produced by a neighboring one, i.e., the resonant dipole-dipole interaction occurring during the electronic excitation.This effect could in principle lead to an additional ͑but probably small, compared with the effect discussed here͒ splitting of the calculated lines.An accurate determination of the dipole moment is required for a proper inclusion of this effect, which is beyond the scope of the present paper.

VI. SUMMARY
We have obtained, within DFT, the structure of 4 He droplets doped with Mg atoms and have discussed in detail the magnesium solvation properties.In agreement with previous DMC calculations, 1,3 we have found that Mg is not fully solvated in small 4 He drops, whereas it becomes fully solvated in a large droplet.
As a consequence of its interaction with the helium environment, it turns out that magnesium is radially quite delocalized inside the droplets.This large delocalization provides a way to reconcile two contradictory results on the solution of one Mg atom in a 4 He drop, namely, center localization ͑LIF and R2PI experiments 5,6 ͒ and surface localization ͑electron-impact ionization experiments 9 ͒.
We have calculated the absorption spectrum of magnesium in the vicinity of the 3s3p 1 P 1 ← 3s 2 1 S 0 transition.For the large Mg@ 4 He 1000 droplet, where Mg is fully solvated, we reproduce the more intense component of the absorption line found by LIF and R2PI experiments in large drops and in liquid helium.This agreement is only achieved when homogeneous broadening due to the coupling of the dipole excitation with the quadrupole deformations of the helium bubble is fully taken into account.This coupling is naturally included in quantum Monte Carlo simulations of the absorption spectrum, 51,52 whereby one takes advantage of the inherent fluctuations present in these simulations.These fluctuations are the full quantal equivalent of the dynamical distortions of the helium bubble we have introduced for the description of homogeneous broadening.An alternative method to include shape fluctuations within DFT has been proposed and applied to the case of Cs in bulk liquid helium. 53It would be interesting to adapt this method to the drop geometry since it is not simple to handle dynamical bubble distortions in a nonspherical environment or in 3 He drops.
To explain the origin of the low-energy peak in the absorption line and confirm the likely existence of soft, foam-   4 He 1000 in the vicinity of the 3s3p 1 P 1 ← 3s 21 S 0 transition in the distorted environment created by the presence of another Mg atom.The line has been decomposed into its two ⌸ and one ⌺ components; the former one is the higher-frequency transition.The starred vertical line represents the gas-phase transition.The experimental curve, adapted for Ref. 5, has been vertically offset for clarity.Also shown is the absorption spectrum obtained neglecting homogeneous broadening ͑hatched region͒.like structure of Mg aggregates in 4 He drops as proposed by Przystawik et al., 6 we have addressed the properties of two Mg atoms in 4 He 1000 and have found that, indeed, Mg atoms are kept apart by the presence of helium atoms that prevent the formation of a compact Mg cluster.We have estimated that the height of the energy barrier for the formation of the Mg dimer in 4 He drops is ϳ2 -3 K, which should be enough, at the droplet experimental temperature of 0.4 K, to guarantee a relatively long lifetime to these weakly bound Mg aggregates.We predict that, contrarily, Mg atoms adsorbed in 3 He droplets do not form such metastable states.
The presence of neighboring Mg atoms in these structures induces a static quadrupole deformation in the helium bubble accommodating a Mg atom.As a consequence, the dipole absorption line around the 3s3p 1 P 1 ← 3s 2 1 S 0 transition splits.We attribute to this static quadrupole moment the origin of the low-energy peak in the absorption line, and confirm the suggestion made by the Rostock group 6 that the splitting of the absorption line, rather than being due to a dynamical ͑Jahn-Teller͒ deformation of the helium bubble, is due to the presence of more than one magnesium atom in the same droplet.
Our previous study on Ca-doped helium drops 15 and the present work show that, within the diatomics-in-molecules approach and provided the impurity-helium pair potentials are accurately determined, DFT is able to quantitatively address the dipole absorption of dopants in 4 He drops whose size is relevant for the analysis of current experiments.Clearly, the influence of an inaccurate pair potential on the solvation properties of a given impurity shows up not only within the DFT approach but also within any other method.A clear example is presented in Ref. 1, where the structure of small drops is addressed within quantum Monte Carlo.Another key issue is how the choice of the density functional may affect the final results.This is as relevant as the choice of the He-He potential within any Monte Carlo approach.In this respect, the parallel evolutions of both DFT and quantum MC approaches have come to establish the Aziz potential 54 and the DF of Ref. 10 as the tools of reference for addressing helium systems.
Finally, we want to point out that, while we have a consistent scenario that explains the results of LIF and R2PI experiments, the understanding of the electron-impact ionization experiment reported in Ref. 9 still requires further analysis.Indeed, since Mg atoms may be in the bulk of the drop or just beneath the drop surface, the experimental ion yield curve should reflect both possibilities, whereas apparently it does not ͑see Fig. 2 of Ref. 9͒.One possible explanation may be that for electron-impact experiments, a N = 10 000 drop is still small, so that the impurity is always close enough to the drop surface to make the Penning ionization process prevail in the direct formation of a He + ion.

ACKNOWLEDGMENTS
We would like to thank Massimo Mella and Fausto Cargnoni for useful comments and for providing us with their Mg-He excited pair potentials, and Marius Lewerenz, Vitaly Kresin, Karl-Heinz Meiwes-Broer, Josef Tigges-bäumker, Andreas Przystawik, Carlo Callegari, and Kevin Lehmann for useful comments and discussions.This work was performed under Grant No. FIS2008-00421 from DGI, Spain ͑FEDER͒ and Grant No. 2005SGR00343 from Generalitat de Catalunya, and under the HPC-EUROPA project ͑Project No. RII3-CT-2003-506079͒, with the support of the European Community-Research Infrastructure Action under the FP6 "Structuring the European Research Area" Programme.

APPENDIX A
In this appendix we obtain the energy of the doped drop up to second order in the deformation parameters and the hydrodynamic mass of the helium bubble.To this end, the helium order parameter and Mg wave function are expressed as ⌿͑r , t͒ = ͱ ͑r , t͒ exp͓i that allow to identify S͑r , t͒ and ͑r , t͒ as velocity field potentials, and the first two terms in Eq. ͑A1͒ as a collective kinetic energy, whose density we denote as t͓ , S , ͉⌽͉ , ͔. Thus, In the adiabatic approximation, the dynamics of the system requires the following steps: ͑i͒ introduce a set of collective variables ͑or deformation parameters͒ ͓␣͑t͔͒ that define the helium density, ͑r , t͒ = ͑r , ͓␣͑t͔͒͒; ͑ii͒ for each helium configuration defined by ͓␣͔, solve the timeindependent Schrödinger equation obeyed by ͉⌽͑r , ͓␣͔͉͒; ͑iii͒ obtain the potential surface V͓ , ͉⌽͉͔ by computing the static energy for each configuration; ͑iv͒ determine the velocity field potentials S͑r , ͓␣͑t͔͒͒ and ͑r , ͓␣͑t͔͒͒ by solving the continuity equations; ͑v͒ compute the collective kinetic energy to obtain the hydrodynamic mass; and ͑vi͒ solve the equation of motion associated with the effective Hamiltonian written as a function the deformation parameters.
We aim to describe harmonic deformations of a spherical helium bubble created by an impurity in the ground state, and have to determine the helium density ͑r , ͓␣͑t͔͒͒ resulting from a change in the radial distance to the center of the spherical bubble induced by the ͓␣͑t͔͒ parameters: where the = 1 deformation is now introduced to allow for displacements of the bubble.We recall that the real spherical harmonics have been normalized as where from now on, the prime will denote the derivative of the function with respect to its argument.

Impurity wave function
To first order, the wave function ͉⌽͑r , ͓␣͔͉͒ is written as The amplitudes ⌽ ͑1͒ ͑r͒ are determined in first-order perturbation theory from the multipole expansion of the impurityhelium pair potential which defines U ͑1͒ ͑r͒.We obtain which shows that, actually, ⌽ ͑1͒ is independent.Once we have obtained the wave function, we can compute the energy surface V͓ , ͉⌽͉͔.Since we describe deformations around a spherically symmetric ground state, the first-order term vanishes, and the derivative g ϵ ‫ץ‬ 2 v / ‫␣ץ‬ ‫ץ‬ ␣ Ј Ј is also spherically symmetric.We can evaluate the second-order contribution to the collective potential energy as

͑A10͒
which defines E ͑2͒ .The parameters ͓␣͔ are the dynamical variables that describe the evolution of the system.

Velocity field potentials
Introducing the expansion S͑r , t͒ϵ͚ ␣ ˙ ͑t͒S ˜͑r͒Z ͑r ˆ͒, where the dot denotes the time derivative, the continuity equation for the liquid helium is, to first order,

͑A12͒
When r → ϱ, the density vanishes for a drop and approaches b for the liquid.In the latter case, Eq. ͑A12͒ reduces to the radial part of the Laplace equation whose general solution is We have solved Eq. ͑A12͒ adapting the method proposed in Ref. 55.Let r i be the first point where 0 ͑r͒ is significantly different from zero ͓ 0 ͑r͒ = 0 for r Յ r i ͔.At this point, Eq. ͑A12͒ implies that which determines the boundary condition at r i : For the liquid, the other boundary condition is that when r → ϱ, the solution behaves as in Eq. ͑A14͒ with A = 0 to have a physically acceptable solution.If the bubble has a sharp surface of radius r i and the liquid is uniform, S ˜͑r͒ is completely determined by Eq. ͑A14͒ and the velocity field potential that fulfills Eq. ͑A16͒ corresponds to the coefficients To find the velocity field potential in a large drop, we have defined a radial distance r b , far from the bubble and from the drop surface at R 1/2 , around which one may consider that the density is that of the liquid.Starting from r = r b with the liquid solution fixed by the coefficients given in Eq. ͑A17͒, we have integrated inward Eq. ͑A12͒, finding the solutions S ˜ inh ͑r͒ and S ˜ h ͑r͒ that correspond, respectively, to the nonhomogeneous and to the homogeneous differential equations that result by setting to zero the left-hand side of Eq. ͑A12͒.The general solution that satisfy the boundary condition Eq. ͑A16͒ is obtained as The field is analogously obtained after introducing the expansion ͑r , t͒ = ͚ ␣ ˙ ͑t͒ ˜͑r͒Z ͑r ˆ͒.We have assumed that the wave function of the impurity in the ground state is a Gaussian ⌽ 0 ͑r͒ = A exp͑−␤r 2 ͒ whose shape has been determined by fitting it to the actual wave function, and have introduced a cutoff distance r g such that safely ⌽ 0 ͑r͒ =0 if r Ն r g .We then obtain the following differential equation to determine ˜:

͑A20͒
The appropriate boundary conditions are and the solution is analytical but very involved ͑we have found it by using the MATHEMATICA package͒.The dipole mode is the only exception; if = 1 we have ˜1 =−r.

Kinetic energy
Once the velocity fields S͑r , t͒ and ͑r , t͒ have been determined, the collective kinetic energy can be easily calculated to second order in the collective parameters: which defines the hydrodynamic mass M ‫ء‬ for each mode.Using the continuity equations and the Gauss theorem, one can find an alternative expression for M ‫ء‬ M ‫ء‬ = 4 2 + 1 ͭ m He͵ drr 2 0 Ј͑r͒S ˜͑r͒ + m Mg͵ drr 2 2⌽ 0 ͑r͒⌽ 0 Ј͑r͒ ˜͑r͒ ͮ .͑A23͒ We have checked that both expressions yield the same values for M ‫ء‬ , which constitutes a test on the numerical accuracy of the method.It is easy to see from Eq. ͑A22͒ that in bulk liquid helium, the = 1 hydrodynamic mass coincides with that given in Ref. 55.Using ˜1 =−r, it is also easy to check from the above expressions that the impurity contribution to the = 1 hydrodynamic mass is just the bare mass of the Mg atom.The sum of Eqs.͑A10͒ and ͑A22͒ represents the Hamiltonian of a set of uncoupled harmonic oscillators whose frequency only depends on .

͑B2͒
where ⑀ i ͑r͒ = ‫ץ‬V i ex ͑r͒ / ‫␣ץ‬ ͉ ͓␣͔ Ն2 =0 .The first integral is zero due to the orthogonality of the spherical harmonics.To evaluate the second integral, we expand ⑀ i ͑r͒ as a power series of r.Taken into account that the PESs have a stationary point at r = 0 due to the spherical geometry-the first-order term is zero-we can safely stop the expansion at the zeroth-order term since the wave function ⌽͑r , ␣ 0 ͒ is very narrow, ͵ dr͉r⌽͑r,␣ 0 ͉͒ 2 ␦Ј͓ − ͑␣ 0 ͒ − V i ex ͑r,␣ 0 ͒/ប͔ͩ 1 ប ⑀ i ͑0͒␣ ͪ.

͑B3͒
This equation may be interpreted as the expansion to first order of a shift in , so it can be written as Thus the eigenvalues ⑀ i ͑0͒ are related to the diagonalization of the expansion of U ij ͑r , ͓␣͔͒ defined in Eq. ͑4͒.Writing this matrix equation as a function of the real spherical harmonics This shows that, to first order, only quadrupolar deformations are coupled to the dipole electronic transition, and that its effect is a shift of the spectral line, as shown in Eq. ͑B4͒.At variance with the approach of Ref. 26, where the above matrix is approximated by its diagonal expression, implying that only the = 0 and 2 components of the quadrupole deformation are considered, our approach incorporates all five components.The relation between the eigenvalues i ͑r , ͓␣͔͒ of U͑r , ͓␣͔͒ and the coefficients ⑀ i ͑0͒ is Finally, the total spectrum is written as I͑͒ ϰ 4 ͵ d␣ 0 d 5 ␣ 2 ͉͑␣ 0 ,␣ 2 ͉͒ 2 I͑,͓␣͔͒ Ӎ 4 ͵ d␣ 0 d 5 ␣ 2 ͉͑␣ 0 ,␣ 2 ͉͒ 2 ͚ m ͵ dr͉r⌽͑r,␣ 0 ͉͒ 2 ϫ␦ ͫ + X gs ͑␣ 0 ͒ − 1 ប i ͑r,͓␣͔͒ ͬ .

͑B8͒ APPENDIX C
In this appendix we consider that the doped drop is cylindrically symmetric.We expand the helium density and ground-state wave function into spherical harmonics with =0, 4␤ ͐d⍀͑r͒Z 0 ͑r ˆ͒, with analogous definitions for ⌽ ˜ and ␥ .If ␤ 0 ӷ ␤ Ͼ0 ͑⇒␥ 0 ӷ ␥ Ͼ0 ͒, we can compute the line shape to first order in ␤ Ͼ0 ͑␥ Ͼ0 ͒; in analogy with Eq. ͑B2͒ we write with ⑀ 0 i defined now by the potential matrix

͑C4͒
where we have defined ␣ ˜20 ϵ ␣ 20 + ␤ 2 C, with These equations show that the computation of the line shape for this geometry is as in the spherical case but with a shift in the ␣ 20 parameter.The dipole absorption spectrum is finally obtained as He drop.The neglect of the quantum kinetic energy of the impurity overestimates the impurity solvation energies by quite a large amount, about 19.7 K for N = 50 and 18.8 K for N = 2000.The value we have found for S 50 ͑Mg͒, −18.4 K, compares well with the DMC result of Ref. 3 ͑ −21 K͒, showing that DFT performs quite well for small clusters, 14 far from the regime for which it was parametrized.

50 FIG. 3 .
FIG. 3. ͑Color online͒ Bottom panel: Total energy ͑K͒ of Mg@ 4 He 50 as a function of Z 0 ͑Å͒ obtained using the Mg-He potentials of Ref. 20 ͑squares͒ and of Ref. 21 ͑circles͒.The energies are referred to their equilibrium values, −157.0 and −153.8K, respectively.The vertical line locates the drop surface at R 1/2 = r 0 N 1/3 , with r 0 = 2.22 Å.The horizontal line has been drawn 0.4 K above the equilibrium energy.Top panel: Probability densities for the configurations displayed in the bottom panel; the single-peak distribution corresponds to the Mg-He interaction of Ref. 21.
FIG.4.͑Color online͒ Helium density plots of the Mg@ 4 He 50 droplet on the y = 0 plane obtained using the Mg-He potential of Ref.20.From top to bottom and left to right, the Z 0 values correspond to the stationary points displayed in Fig.3, namely, 0, 2, 4, and 6 Å.The brighter regions are the higher-density ones.

FIG. 5 .
FIG. 5. ͑Color online͒ Helium density profiles of the Mg@ 4 He 50 droplet along the z axis obtained using the Mg-He potential of Ref. 20.The Z 0 value is indicated in each panel.Dashed lines, pure drops; solid lines, doped drops.

FIG. 6 .
FIG. 6. ͑Color online͒ Bottom panel: Total energy ͑K͒ of Mg@ 4 He 1000 as a function of Z 0 ͑Å͒ obtained using the Mg-He potentials of Ref. 20 ͑squares͒ and of Ref. 21 ͑circles͒.The energies are referred to their equilibrium values, −5482.0 and −5476.8K, respectively.The vertical lines roughly delimit the drop surface region, conventionally defined as the radial distance between the points where the density equals 0.1 b and 0.9 b , with b = 0.0218 Å −3 as the bulk liquid density.The horizontal line has been drawn 0.4 K above the equilibrium energy.Top panel: Total energy ͑K͒ of Mg@ 4 He N with N = 1000 and 2000 ͑vertical right scale͒ as a function of Z 0 ͑Å͒ obtained using the Mg-He potential of Ref. 20.The energies are referred to their equilibrium values, −5482.0 and −11 629.8 K, respectively, and the distances ͑horizontal scale͒ are referred to the R 1/2 radius.Also shown are the corresponding probability densities ͑vertical left scale͒: dotted ͑dot-dashed͒ line, N = 1000 ͑2000͒.The solid line represents the probability density of the N = 10 000 drop.

FIG. 7 .
FIG.7.͑Color online͒ Total absorption spectrum of one Mg atom attached to4 He 1000 in the vicinity of the 3s3p 1 P 1 ← 3s 2 1 S 0 transition.The line has been decomposed into its two ⌸ and one ⌺ components; the former is the higher-frequency transition.The starred vertical line represents the gas-phase transition.The experimental curve, adapted for Ref. 5, has been vertically offset for clarity.Also shown is the absorption spectrum obtained neglecting homogeneous broadening ͑hatched region͒.

50 FIG. 10
FIG.10.͑Color online͒ Energy ͑K͒ of the Mg+ Mg@ 4 He 1000 system as a function of the Mg-Mg distance ͑Å͒.The energies have been referred to that of the metastable equilibrium configuration ͑local minimum͒ at L = 0.The lines have been obtained by a cubic spline of the actual calculations.From bottom to top, the curves correspond to L =0-50ប in 10ប steps.

FIG. 11 .
FIG.11.͑Color online͒ Energy ͑K͒ of the Mg+ Mg@ 3 He 1000 system as a function of the Mg-Mg distance ͑Å͒.The line has been drawn to guide the eyes.

FIG. 12 .
FIG.12.͑Color online͒ Top panel: Helium density profiles of the ͑Mg+ Mg͒ @ 3 He 1000 complex at d = 8.6 Å along the z axis ͑solid line͒ and the x or y axis ͑dashed line͒.Bottom panel: Equidensity lines corresponding to the same configuration.The brighter regions are the higher-density ones.

FIG. 13 .
FIG.13.͑Color online͒ Total absorption spectrum of one Mg atom attached to4 He 1000 in the vicinity of the 3s3p 1 P 1 ← 3s 21 S 0 transition in the distorted environment created by the presence of another Mg atom.The line has been decomposed into its two ⌸ and one ⌺ components; the former one is the higher-frequency transition.The starred vertical line represents the gas-phase transition.The experimental curve, adapted for Ref. 5, has been vertically offset for clarity.Also shown is the absorption spectrum obtained neglecting homogeneous broadening ͑hatched region͒.

TABLE I .
Atomic shift ⌬ of Mg@ 4 He 1000 ͑R 1/2 = 22.2 Å͒ as a function of the average distance between the magnesium atom and the center of mass of the 4 He 1000 moiety.Also indicated is the corresponding wavelength .The value of the transition energy in the gas phase is 35051 cm −1 ͑Ref.34͒.