The delocalization of SF6 in helium clusters

2 September 1994

ZLSEVIER

CHEMICAL PHYSICS LETTUIS

Chemical Physics Letters 227 (1994) 143-148

The delocalization of SF6 in helium clusters * E. Krotscheck avb,S.A. Chin b l

Institutefor Theoretical Physics, University of California, Santa Barbara, CA 93106, USA b Department of Physics, Texas A&M University, College Station, TX 77843. USA Received 28 April 1994; in fmal form 29 June 1994

Variational calculations of the ground state structure and the excitation spectrum of 4He clusters doped with an SFs molecule show that, for clusters size larger than approximately 100 atoms, a spherically symmetric configuration with the impurity molecule at the center of the cluster is unstable. The instability is indicated by the ‘softening’ of the droplet’s collective dipole excitation. This instability strongly suggests that the SF, molecule is delocalized inside larger helium droplets and may be found near a non-isotropic environment such as the surface.

Based on the observed splitting of the infrared absorption spectrum of SF6 in helium clusters as reported by Goyal et al. [ 1,2], one can conclude that the SF6 molecule is not rigidly located at the center of the cluster, but must be allowed to be near a non-isotropic environment such as the surface. In this Letter, we provide the first theoretical confirmation of this delocalization effect. Our analysis is based on an unbiased variational calculation of the ground state and the excitation spectrum of such doped helium clusters. By studying the excitation energies, we can indeed conclude that a spherically symmetric contiguration of helium atoms, with the SF6 molecule rigidly fixed at their center, is dynamically unstable against infinitesimal relative displacement for sufficiently large cluster sizes with Z 100 particles. Our work is thus supportive of an explanation for the observed SF6 spectrum splitting based on impurity de* Work supported, in part, by the National Science Foundation grantsNo. PHY92-13502 (to SAC), PHY91-08066 (to EK), and PHY89-04035 to the Institute for Theoretical Physics (ITP) in Santa Barbara.

localization. However, the issue of whether delocalization itself can completely account for the amount of observed splitting, which requires detail knowledge of the impurity distribution after delocalization, is beyond the scope of the present Letter, and must be addressed by other means. Technically, our calculation is based on a ground state calculation of the many-boson system with an approximate wavefunction of the Jastrow-Feenberg form

+

C d-i, rj)+

iej

C

i-cjck

u3(ri,

rj,

rkk) >

(1)

and uses integral equation techniques to calculate the physical quantities of interest. The one-body function uI (r) impresses the non-uniform structure, the two-body function Uz(ri, rj) takes care of the shortrange repulsion, and the three-body function u3( ri, r, rk) describes triplet correlations. The correlations are assumed to be the most general ones permitted by the

0009-2614/94/$07.00 0 1994 Elsevier Science B.V. All rights reserved SSDIOOO9-2614(94)00807-8

144

E. Krotscheck, S.A. Chin /Chemical PhysicsLetters 227 (1994) 143-148

geometry, for example, the pair correlation function ul( rj, ri) depends on the distance between these two particles, and on the distance of both particles from the center of the droplet, Essential to the method is that the correlations are determined by minimizing the ground state energy E,, with respect to the correlation functions uZ(r, , .... rn )

in an unrestricted function space. In practice, one can include up to triplet correlations in this scheme, which up to now has always been found to be sufficient. The method has matured in recent years to a level were the quality of its predictions are comparable to those of Green functions and diffusion Monte Carlo techniques. In terms of computational effort, the method is an order of magnitude greater than density functional techniques but an order or two less than Monte Carlo methods. For the problem at hand, one of the most important and desirable feature of the theory is that it has a ‘built-in’ indicator for deciding whether an assumed system configuration is stable against infinitesimal changes. This is a unique feature, and it is worth contrasting it with alternative theoretical approaches to the same problem: Variational Monte Carlo calculations start from some assumed analytic wavefunctions similar to ( I), and perform the integrations by Monte Carlo methods. This precludes the solution of the Euler equations (2) in an unrestricted function space. Quite sophisticated variational wavefunctions have been constructed [ 31 that are capable of producing a shell-like distribution of helium atoms around an impurity, but these are strongly dependent on one’s physical intuition (or prejudice) for the system under consideration. Green function and diffusion Monte Carlo techniques can in principle overcome this problem, but in practice these methods may be biased by an ill-chosen importance sampling function or require an excessively long convergence time to ‘find the correct geometry. Similar statements hold for density functional calculations, which are gaining popularity in the investigation of cluster properties [4-61. In these ap proaches, one assumes a certain geometry, and can a posteriori determine whether the assumed geometry is stable, for example by studying the excited states of the system [ 4,5 1. A revealing example contrasting

the different philosophies underlying microscopic and density-functional methods is clusters doped with alkali metal atoms, which were recently studied using a non-local density functional theory [ 61. Such calculations are impossible within the optimized variational scheme since spherically symmetric configurations of helium atoms with an alkali metal atom in the center are unstable and do not occur in nature. In this Letter, we report on the structure and excitation calculations of helium droplets doped with a single SF6 molecule. The overall quality of our theoretical prediction has been checked in the case of pure droplets where the identical computational procedure was used to compute the ground state energy, one and two-body densities, collective excitations and dynamic structure functions for droplet size up to N= 1000 [ 7 1. When compared to our earlier DMC calculations [ 81 for cluster size N< 112, excellent agreements were found for the droplet energetics and dynamic structure functions, only the droplet central densities ( ~50.020 Am3) were slightly lower than the DMC results ( x 0.022 A-‘). Additional comparison with available diffusion Monte Carlo calculations on doped helium droplets [ 31 will be given below. The technical aspects of the optimization of the correlation functions through the Euler equations (2) have been described elsewhere (see, for example, Ref. [ 91) and will not be repeated here. We treat the impurity particle as an object with infinite mass located at the coordinate origin and use the spherically symmetric approximate potential by Pack et al. [ lo]. The assumption of spherical symmetry significantly simplifies the computational requirement of our method. It is known that deviations from a spherically symmetric potential can cause notable changes in the energetics and structure of small clusters [ 3 1, but are expected to have little effect on the question of stability for large droplets considered here. For clusters doped with an SF, molecule, Monte Carlo calculations have been carried out by Bamett and Whaley [ 31. Fig. 1 shows our optimized variational ground state energies as compared to their VMC and DMC results. Also for comparison, we show our variational results for the case of pure droplets. On the energy scale of the graph, our variational energies are indistinguishable from our earlier DMC results [ 81 at N= 20, 40, 70, and 112. The best global fit to our energy per particle for pure clusters with 70 or more particles is

E. Krotscheck, S.A. Chin / Chemical Physics Letters 22 7 (1994) I43- 148 0

145

r

-5 -10

_ ? “S

8 -15

0.06

is -20 -2s -30 d 1000

VMC o

240

DMC A

112

40

70

20

N on N-II3 scale

1

Fig. 1. The ground state energy of pure (dashed line with solid diamonds) and SF,-doped helium clusters (solid line with solid squares) is shown as a function of N-‘13. The &shed line is our tit to the energies (3) the DMC data of Ref. [ 81 are indistinguishable from our variational calculation. Also shown are the variational (open diamonds) and DMC (open triangles) data of Bamett and Whaley [ 31.

E(N) -=-7.21+17.71N-1’3 N

-5.95N-2’3K.

(3)

Clusters with 20 and 40 particles were excluded from the tit since these are surface dominated clusters whose central densities have not yet reached the bulk value. The fit Eq. (3) extrapolates to a bulk energy per particle of - 7.21 K and a surface energy of 0~0.272 K AB2, which compares favorably with the experimental value [ 111 of 6% 0.256 K Am2, For the case of SF,, our energies, through having the correct trend and are lower than Bamett and Whaley’s [ 31 VMC results, are still somewhat higher than their DMC energies. This is not unexpected, since although our treatment of the many-body correlations has been demonstrated to be accurate around the normal bulk density [ 121 there is no a priori reason to expect that this same accuracy can be maintained up to four times (see the next paragraph) its intended density. Our density profiles for doped droplets are shown in Fig. 2. Because of the strong SF6 attraction, the first shell of helium atoms surrounding the impurity is compressed to a density nearly four times that of the equilibrium density of bulk 4He. In comparison with the DMC results of Bamett and Whaley [ 3 1, our results for the densities are only slightly lower, resulting in a slightly more rapid filing of the ‘shells’ around

(4

Fig. 2. The variational ground state density of SFc,-doped‘He clusters with N=40, 70, and 90 particles (solid lines) and their corresponding lowest dipole excitation function f(r) (dashed lines). Also shown are the DMC densities of Bamett and Whaley [ 31 for N=39 and 69 (dotted lines); errorbarsare omitted for the sake of clarity.

the impurity molecule. The significance of the dashed curves will be discussed below. Having established the validity of our ground state calculations, let us now turn to the stability problem and the computation of excited states. As a part of the algorithm for optimizing the correlations demanded by the extremum condition (2)) one calculates the dynamic response function, which is given by the usual RPA relation X(ri, r2; 0) x vph(r3,

=xo(~

r4)x(r4,

r2; WI+

I

r2; w),

d3r3 d3r4zoh

r3; 0)

(4)

where x0 (r,, r2; co) is the response function of a ‘noninteracting’ system defined by the one-body Hamiltonian

(5)

In Eq. (4), the local particle-hole interaction is defined as the second variational derivative of the correlation energy functional with respect to the density

VP&, r2) =

S2E Wc)

&P(r2) *

(7)

146

E. Krotscheck, S.A. Chin /Chemical Physics Letters 227 (1994) 143-148

function, the energy difference iio~ E,-EO determined by Eq. ( 10) would be an upper bound to the exact excitation energy. The lowest possible upper bound for the excitation energy is obtained by minimizing the energy difference fiw with respect to the excitation functionf( r). This leads to the eigenvalue equation

The dynamic response function (4) is solved iteratively by noting that the static structure function is given by

(8)

X --03

In the iterative process, the particle-hole interaction can be expressed in terms of the bare interaction, the one-body density, and the static structure function, and the process is repeated until convergence is reached. From this procedure, it is evident that the iterations can be carried out only if all collective excitations are stable. A different route to calculate excitations, which leads to exactly the same result, is to generalize the Feynman theory of excitations [ 131 to non-uniform geometries [ 8,141. We write a trial function of the excited state as Ul,(r, 9***,rN) =F(5, = i~lfw

(11) with the coordinate space representation of the static structure function

.... rN) %U,(c, .... rN)

%I(5 3.... fN) .

(9)

If !&(r,, ...) rN) is the exact ground state wavefunction, or if it is a Jastrow-Feenberg wavefunction with optimized pair correlations, the energy difference between the ground state and the trial excited state can be expressed as

+j d3rd3r’f(r)

[p2(r, r7 -pi (r)p, (r’) V(r)}-’ (10)

where&=( YF]H] YF)/( ‘v,] YF) istheenergycorresponding to the wavefunction (9)) E,,is the ground state energy, and p2 (r, r’) is the pair density. In the above equations, we have left out the center-of-mass corrections that have to be applied in the pure droplet cases. If Yo(r,, ...) rN) were the exact ground state wave-

,

Eq. ( 11) is the desired generalization of the bulk Feynman dispersion relation [ 13 1. For the interpretation of our results, it is useful to consider the Feynman excitation functionf(r) since it gives information on the change of the microscopic wave function. For example, f(r) =const. corresponds to a renormalization of the wavefunction (i.e. no motion of the particle at all), whereas f(r) =c-r corresponds to a rigid translation of the droplet. Both cases will be relevant below. We have calculated the excitation energies of a close sequence of SF6 doped helium droplets from cluster size N=20 to N= 112 and two larger droplets with N= 150 and N= 240. In the course of these calculations, we found a characteristic dependence of the low-lying excitations on the droplet’s ‘shell structure’ surrounding the impurity atom. Similar structures have been found in liquid films, where the dependence of the excitation energy on the amount of 4He can be strong enough to drive a series of phase transitions [ 9,151. Such transitions, corresponding to a spontaneous breaking of the spherical symmetry in small droplets, with the impurity atom remaining at the droplet’s center-of-mass, were not found in the present case. What we did discover, however, was a softening of the dipoleexcitationaround droplet sizes of 85-90 particles. This excitation corresponds to a relative oscillation between the helium droplet and the impurity atom. An instability of this mode means that the spherically symmetric configuration of helium atoms around the impurity is unstable against

E. Krotscheck S.A. Chin /Chemical PhysicsLetters 227 (1994) 143-148

infinitesimal displacements. While such an instability gives no information about the energetically most favorable configuration (i.e. the exact manner in which the impurity delocalizes), it does imply that the assumed spherically symmetric configuration with the impurity atom at the center of the droplet is not the energetically lowest state. Our results for the dipole excitation energy of SF6-doped droplets are shown in Fig. 3. As the droplet increases in size, the dipole excitation energy plunges toward zero at around N= 85-90. For additional checks on the onset of this instability, to be discussed below, the size dependence of the quadrupole energy for both doped and pure clusters are also shown. If SF6 indeed delocalizes for NZ 85-90, how were the larger droplet energies shown in Fig. 1 obtained? According to our earlier argument, with the assumed geometry, one should not be able to obtain solution beyond Nk 85-90. The answer is that up to particle numbers of N= 112, the calculation was continued by simply omitting the unstable mode from the energy integration in Eq. (8). While the process is not entirely justifiable, it had little effect on other calculated quantities such as the ground state energy, densities or excitation energies at higher angular momentum. The only noticeable effect is the slight ‘upward’ curvature of the quadrupole energy at 2.5

0.01,’1000

/

,s'

I 240

112

70

40

20

N on N-“3 scale

Fig. 3. The dipole excitation energies of SF,doped helium drop lets as a function of cluster size N (solid line with solid diamonds). Also shown are the quadrupole excitation energies (dashed line with triangles), the quadrupole energy ofpure drop lets with the same ground state theory (dashed line), the density functional result of Barranco and Hemandez [ 5 ] (dotted line) and the liquid drop limit (dash-dotted line).

147

N> 100 particles, as seen in Fig. 3. However, even when the stability of the droplet was forced by ignoring the unstable dipole excitation, the quadrupole excitation eventually became unstable at cluster size NZ 150. To obtain solutions for even larger N, the Euler equation (2) was discretized and solved on a mesh that was slightly ‘too small’ for the size of the droplet. Effectively, the droplet was solved in a spherical cavity with hard walls. This again has little effect on the ground state energy, but increases all excitation energies such that the instability is avoided. Since this procedure is totally ad hoc with questionable physical relevance, it clearly demonstrates that reasonably looking energies need not reflect reasonable underlying physics. The physical nature of the ‘soft’ modes is most clearly displayed by considering the Feynman excitation functionf( r) in the dipole channel. As pointed out earlier, a constant value means that the system is rigid, whereas a linear function indicates translation. Back in Fig. 2, three such excitation functions were shown. One clearly sees that the first shell of atoms around the impurity is rigid and does not participate in the excitation. The second shell shows a modest amount of translation, and the rest shows only translational motion. What appears to be happening physically is that at least one shell of 4He atoms is rigidly connected to the impurity molecule. This shell of helium atoms essentially shields other helium atoms from the impurity. The impurity, plus this shell of highly compressed helium atoms at about four times their equilibrium density can be thought of another ‘effective impurity molecule’ with reduced impurityhelium interaction. As illustrated in Fig. 2, these atoms do not participate in the motion of the rest of cluster. Moreover, being at a very high density, the compressed helium atoms will repel others. This results in the prominent ‘dip’ of the density after the first shell. But it apparently also has the effect of shielding the impurity such that the attraction between the compound impurity and the rest of the cluster is not strong enough to keep this compound at the center. We have carried out a number of different tests to assure that our basic result on the softening of the dipole mode is not an artefact of the approximations made when solving the optimization problem (2). It is a well-known feature of variational theory that the

148

E. Krotscheck,S.A. Chin /Chemical PhysicsLetters227 (1994) 143-148

particle-hole interaction derived from the Euler equation for the pair correlations ( Eq. ( 2 ) for n = 2 ) , which appears in the RPA equation (4) is identical with the definition (7) only if no other approximations were made. This is never the case, it is therefore important to verify that the instability is not driven by the approximations. In order to verify that, we have also calculated the excitation energies for pure droplets using exactly the same method. We did not find any such instability, moreover we found [ 7 ] excellent agreement between the features of the excitation spectrum and the surface energy obtained from the mass-formula fit (3). It is noteworthy that there appears to be no convergence of the doped quadrupole energy towards it’s pure cluster counterpart with increasing particle number. Finally, we note that Barranco and Hemandez [ 51 have found a similar lack of convergence between the quadrupole mode of the pure cluster and a xenon-doped cluster. In this case, the quadrupole instability occurs at somewhat larger particle numbers [ 16 1. Bamett and Whaley’s DMC calculation shows no indications of instability similar to the kind discussed here. For their largest droplet size of N= 499, their SF6 remains localized to within z 0.5 8, of the droplet’s center. We believe that this is due to their choice of the helium-impurity correlation function which unduly constrains the impurity to the cluster center [ 171. Removing this constraint, the average separation between SF, and the cluster center in a VMC calculation [ 17 ] can be as large as 3 8, for cluster size N=240. We are currently redoing the DMC calculation with such an improved correlation function. A more recent path integral Monte Carlo (PIMC ) calculation by Ceperley [ 18 ] also reaches a similar conclusion. For a neon impurity in a cluster of 78 helium atoms, his calculation shows that the center-of-mass of the helium atoms is on the average about 2 A away from the impurity atom. We have described in this Letter a structural instability of doped quantum liquid drops which indicates that the energetically most favourable position of a heavy impurity atom or molecule in a helium cluster is not necessarily at it’s center. Our result is in general supportive of a qualitative interpretation of infrared spectroscopy data based on impurity delocalization [ 1,2]. We are presently studying the structure, stabil-

ity and excitation spectra of helium clusters doped with rare gas impurities. Since these impurities are lighter with weaker attractions, one would expect their dipole instabilities to occur at smaller cluster sizes. While such systems are more difficult to study experimentally, they are more accessible to exact methods of ground state calculations.

EK thanks the ITP for warm hospitality and support during the spring of 1994. The paper profited from discussions with C.E. Campbell, D. Ceperley, W. Kohn, M. Saarela, and B. Whaley. We thank R.N. Barnett, M. Barranco, D. Ceperley and K-B. Whaley for communication of unpublished material. The computations have been carried out on an IBM RS/ 6000-590 which was provided by IBM for benchmarking purposes. We would like to thank IBM for providing us with this opportunity.

References [ I] S. Goyal, D.L. Schutt and G. Stoles, Phys. Rev. Letters 69 (1992) 933. [2] S. Goyal, D.L. Schutt and G. Stoles, J. Phys. Chem. 97 (1993) 2236. [ 31 R.N. Bamett and KB. Whaley, J. Chem. Phys. 99 (1993) 9730. [4] M. Casas and S. Stringari, J. Low. Temp. Phys. 79 ( 1990) 135. [ 51 M. Barranco and E.S. Hemandez, Phys. Rev. B, in press. [ 61 F. Dalfovo, Z. Physik D 29 ( 1994) 61. [ 71 E. Rrotscheck and S.A. Chin, submitted for publication. [S] S.A. Chin and E. Rrotscheck, Phys. Rev. B 45 (1992) 852. [ 91 B.E. Clements, J.L. Epstein, E. Rrotscheck and M. Saarela, Phys. Rev. B 48 (1993) 7450. [lo] R.T Pack, E. Piper, G.A. Pfeffer and J.P. Toennies, J. Chem. Phys. 80 ( 1984) 4940. [ 111 M. Iino, M. Suzuki and A.J. Ikhushima, J. Low Temp. Phys 61 (1985) 155. [ 121 E. Rrotscheck and M. Saarela, Phys. Rept. 232 ( 1993) 1. [ 13 ] R.P. Feynman, Phys. Rev. 94 ( 1954) 262. [ 141 C.C. Chang and M. Cohen, Phys. Rev. B 11 (1975) 1059. [ 151 B.E. Clements, E. Rrotscheck and H.J. Lauter, Phys. Rev. Letters 70 (1993) 1287. [ 161 M. Barranco, private communication. [ 171 S.A. Chin and E. Rrotscheck, talk presented at the Second Workshop on Quantum Liquid Clusters, June 12-l 5,1994, S&loss Ringberg, Bavaria, Germany. [ 18 ] D. Ceperley, private communication.