- Email: [email protected]
Path integral simulation of positronium
Nuclear Instruments and Methods in Physics Research B 192 (2002) 176–179 www.elsevier.com/locate/nimb
Path integral simulation of positronium Bruce N. Miller a
a,*
, Terrence Reese
b
Department of Physics and Astronomy, Texas Christian University, Fort Worth, TX 76129, USA Department of Physics, Southern University and A&M College, Baton Rouge, LA 70813, USA
b
Abstract We use path integral Monte Carlo to investigate the dependence of the positronium pick-off annihilation rate on density in Xenon on two isotherms above the critical temperature. Using a simple interaction potential, reasonable agreement with experiment is obtained over a wide range of density. Surprisingly, we find that positronium remains localized at twice the critical point density. Ó 2002 Published by Elsevier Science B.V. PACS: 71.60.þz; 78.70.Bj; 31.15.Kb; 36.10.Dr; 82.30.Gg Keywords: Positronium; Path integral; Monte Carlo
1. Introduction The generic model consisting of a low mass particle in thermal equilibrium with a simple fluid encompasses a range of interesting phenomena, including electron transport and positron and positronium decay [1–4]. It has been known for some time that gross changes in the electron mobility [5] or the positron lifetime [6] can accompany relatively small changes in fluid density in a neighborhood of the liquid–vapor critical point. Theorists have modeled these systems using variants of mean field theory [7,8,13] and Feynman’s polaron approximation [9]. Experimental measurements of the ortho-positronium pick-off decay rate are of special interest because they provide direct information concerning the local environment of the equilibrated quantum particle (qp) [10]. The pur-
pose of this study is to investigate the properties of positronium in thermal equilibrium in Xenon fluid. We are particularly interested in the dependence of the ortho-positronium pick-off annihilation rate on the fluid density and temperature which, in turn, depends directly on the positronium–atom correlation in position [10]. The former has been measured experimentally [11,12] while we have shown that both the eþ and Ps annihilation rates can be computed from first principles using path integral Monte Carlo (PIMC) [10,14,15] At low fluid density, the path integral approach can be used to evaluate terms in a density expansion of the annihilation rate [16]. In the following we describe our recent study of the super-critical Ps–Xenon system. 2. Theory 2.1. Adiabatic approximation
*
Corresponding author. Tel.: +1-817-257-7123; fax: +1-817257-7742. E-mail address: [email protected] (B.N. Miller).
The system consists of positronium in thermal equilibrium with fluid Xenon above the critical
0168-583X/02/$ - see front matter Ó 2002 Published by Elsevier Science B.V. PII: S 0 1 6 8 - 5 8 3 X ( 0 2 ) 0 0 8 6 4 - 9
B.N. Miller, T. Reese / Nucl. Instr. and Meth. in Phys. Res. B 192 (2002) 176–179
temperature. The Ps–Xe interaction, which is dominated by exchange forces [7], is represented by a hard sphere potential, and the Xe–Xe interaction is represented by a Lennard-Jones (6,12) potential. Assume the atoms obey classical mechanics and Ps obeys quantum mechanics. Then the partition function of the hybrid system consisting of N atoms and a single qp is [10,14,15] Z b Z ¼ ð1=K3N N !Þ dRebU ðRÞ Trfeb H g; ð1Þ where U is the total interatomic potential energy, R is the set of atomic positions Ri , K is the thermal wavelength of an atom, and the Ps hamiltonian is b ¼b H p 2 =2m þ
N X
177
varying points along the Feynman chain, represented by the imaginary time indices t and t0 . Important structural information concerning the qp and the local properties of the fluid in its vicinity is provided by the qp–atom radial distribution function * + p1 N X X gfp ðrÞ ¼ dðRj ri rÞ=p =q: j¼1
i¼0
The decay rate is an observable and can be represented by the operator b k¼
N X
f ðRj rþ Þ;
ð4Þ
j¼1
V ðjRj rjÞ;
ð2Þ
j¼1
where V is the Ps–atom pair potential. The mean value of a physical observable O is then Z b ð3Þ hOi ¼ dRebU ðRÞ Trfeb H Og=ZðK3N N !Þ:
2.2. Path integral Monte Carlo For the purposes of computing statistical averages, Feynman showed that the qp can be represented by a chain of p harmonically coupled ‘‘pseudo-particles’’ with a temperature dependent force constant [17]. Each pseudo-particle also interacts with each classical atom through the reduced potential V =p. This construction transforms the calculation of the equilibrium properties of a single qp into that of p classical pseudo-particles. Statistical averages are performed by averaging over many realizations of the chain and atom configurations. The size of the chain is increased until convergence is obtained. In order to improve the convergence, we used the image approximation for the qp–atom hard sphere potential [18]. The details are discussed in one of our earlier papers [14]. 2.3. Important quantities Information concerning the degree of localization of the qp in the host fluid is provided by the mean square separation of pseudo-particles at
where f ðRj rþ Þ is the local electron density contributed by the jth atom at the positron position rþ . The mean Ps decay rate can be directly expressed in terms of gfp [10], D E Z b b k ¼ dRebU ðRÞ Trfeb H b k g=ðZK3N N !Þ ¼ ½q=ð8pa30 Þ Z Z dR0 dx f ðjR0 x=2jÞex=a0 gfp ðR0 Þ: ð5Þ With the additional assumption that the electronic charge is concentrated at the nucleus, we obtain Z D E 0 3 b k ¼ q=ðpa0 Þ dR0 e2R =a0 gfp ðR0 Þ; ð6Þ for the mean decay rate.
3. Results We performed simulations over a wide range of fluid densities at 300 and 340 K, where experimental measurements of the o-Ps decay rate in fluid Xenon have been reported [11]. The calculations were made using eight density values (q ¼ 0:017, 0.035, 0.088, 0.17, 0.3, 0.35, 0.4 and 0.7) on each isotherm. Here the scaled density q ¼ qr3 where r is the Lennard-Jones radius. In these units the density of the critical point is about 0.35.
178
B.N. Miller, T. Reese / Nucl. Instr. and Meth. in Phys. Res. B 192 (2002) 176–179
Fig. 1. RMSD at T ¼ 300 K versus (t t0 ) at three density values. The extent of localization of the Ps increases with density.
D E 2 1=2 Fig. 1 provides a plot of jrðtÞ r0 ðt0 Þj , the root mean square displacement (RMSD) versus pseudo-particle index ðt t0 Þ along the Feynman chain, at T ¼ 300 K for three density values and the free particle. The top curve is the plot of the free chain (q ¼ 0:0). As the density is increased, the curve becomes flatter, indicating that the spread of the polymer is becoming smaller, and the Ps more localized in the fluid. Self-trapping appears complete at q ¼ 0:7: Fig. 2 presents a plot of the radial distribution function versus position for T ¼ 300 K at three different densities. The local fluid distortion induced by the qp is greatest at the critical density,
Fig. 3. Comparison of the experimental and theoretical pick-off decay rate at T ¼ 300 K. The straight line is the linear extrapolation of the low density dependence.
while the qp confinement is greatest at twice this density. Fig. 3 presents a plot of the decay rate versus q for T ¼ 300 K. It displays both the theoretical computations and the experimental results. The straight line is the linear extrapolation of the low density dependence of the decay rate. The theoretical curve is close to the experimental result. Similar results were obtained on the 340 K isotherm [19].
4. Conclusions
Fig. 2. Ps–atom radial distribution function at three densities versus position at T ¼ 300 K.
The results of the RMSD computations show that as the density is increased, the qp becomes more confined. At the highest density, the distance between almost all of the pseudo-particles in the polymer representation is the same, indicating a completely self-trapped state. In the plots of the polymer-fluid radial distribution functions it is seen that the greatest number of fluid molecules are expelled from the vicinity of the Ps–atom at the critical density. This results in the sharply reduced probability of annihilation observed in experiments. Mean field theory suggests that well above the critical density the decay rate returns to the predicted linear density dependence due to collapse of the ‘‘bubble’’ [8,13]. Here we see for the
B.N. Miller, T. Reese / Nucl. Instr. and Meth. in Phys. Res. B 192 (2002) 176–179
first time that no such thing occurs. At the highest density value (q ¼ 0:7), it is noted that the large local density of fluid molecules at the edge of the bubble results in the increased probability of positron annihilation. Compared with earlier PIMC computations, the plots of the decay rate versus density now show curvature, sharing some of the attributes of the experimental curves. Both at 300 and 340 K they also show the rapid, non-linear, increase above the critical density seen in other experiments. Unfortunately there isn’t sufficient data available for Xenon to compare with the simulations at these higher densities. However, on the 300 K isotherm near the critical density, the theoretical curve doesn’t show a plateau and lies above the experimental results above the critical density. The 300 and 340 K theoretical isotherms have similar shapes: In particular, they both show an anomalous rise (bump) above the critical density which is not present in the experimental data. It will be interesting to determine if this is a consequence of either the simplified interaction potential or atomic electron density used to perform the simulations.
179
References [1] J. Hernandez, Rev. Mod. Phys. 63 (1991) 675. [2] D. Chandler, K. Leung, Annu. Rev. Phys. Chem. 45 (1994) 557. [3] G.N. Chuev, Izv. Akad. Nauk., Ser. Fiz. 61 (1997) 1770. [4] B.J. Berne, D. Thirumali, Ann. Rev. Phys. 37 (1986) 401. [5] A.F. Borghesani, M. Santini, Phys. Rev. A 45 (1992) 8803. [6] P. Hautojarvi, K. Rytsola, P. Tuovinen, P. Jauho, Phys. Lett. A 57 (1976) 175. [7] I.T. Iakubov, A.G. Khrapak, Rep. Prog. Phys. 45 (1982) 697. [8] B.N. Miller, T. Reese, Phys. Rev. A 39 (1989) 4735. [9] A.L. Nichols III, D. Chandler, Y. Singh, D.M. Richardson, J. Chem. Phys. 81 (1984) 5109. [10] B.N. Miller, Y. Fan, Phys. Rev. A 42 (1990) 2228. [11] M. Tuomisaari, K. Rytsola, P. Hautojarvi, J. Phys. B 21 (1988) 3917. [12] T.J. Murphy, C.M. Surko, J. Phys. B 23 (1990) L727. [13] T. Reese, B.N. Miller, Phys. Rev. A 42 (1990) 6068. [14] T.L. Reese, B.N. Miller, Phys. Rev. E 47 (1993) 2581. [15] G.A. Worrell, B.N. Miller, Phys. Rev. A 46 (1992) 3380. [16] B.N. Miller, T.L. Reese, G. Worrell, Phys. Rev. E 47 (1993) 4083. [17] R.P. Feynman, Statistical Mechanics, Addison-Wesley, Reading, MA, 1972. [18] P.A. Whitlock, M.H. Kalos, J. Chem. Phys. 30 (1978) 361. [19] T.L. Reese, B.N. Miller, Phys. Rev. E 64 (2001) 061201.
Path integral simulation of positronium Bruce N. Miller a
a,*
, Terrence Reese
b
Department of Physics and Astronomy, Texas Christian University, Fort Worth, TX 76129, USA Department of Physics, Southern University and A&M College, Baton Rouge, LA 70813, USA
b
Abstract We use path integral Monte Carlo to investigate the dependence of the positronium pick-off annihilation rate on density in Xenon on two isotherms above the critical temperature. Using a simple interaction potential, reasonable agreement with experiment is obtained over a wide range of density. Surprisingly, we find that positronium remains localized at twice the critical point density. Ó 2002 Published by Elsevier Science B.V. PACS: 71.60.þz; 78.70.Bj; 31.15.Kb; 36.10.Dr; 82.30.Gg Keywords: Positronium; Path integral; Monte Carlo
1. Introduction The generic model consisting of a low mass particle in thermal equilibrium with a simple fluid encompasses a range of interesting phenomena, including electron transport and positron and positronium decay [1–4]. It has been known for some time that gross changes in the electron mobility [5] or the positron lifetime [6] can accompany relatively small changes in fluid density in a neighborhood of the liquid–vapor critical point. Theorists have modeled these systems using variants of mean field theory [7,8,13] and Feynman’s polaron approximation [9]. Experimental measurements of the ortho-positronium pick-off decay rate are of special interest because they provide direct information concerning the local environment of the equilibrated quantum particle (qp) [10]. The pur-
pose of this study is to investigate the properties of positronium in thermal equilibrium in Xenon fluid. We are particularly interested in the dependence of the ortho-positronium pick-off annihilation rate on the fluid density and temperature which, in turn, depends directly on the positronium–atom correlation in position [10]. The former has been measured experimentally [11,12] while we have shown that both the eþ and Ps annihilation rates can be computed from first principles using path integral Monte Carlo (PIMC) [10,14,15] At low fluid density, the path integral approach can be used to evaluate terms in a density expansion of the annihilation rate [16]. In the following we describe our recent study of the super-critical Ps–Xenon system. 2. Theory 2.1. Adiabatic approximation
*
Corresponding author. Tel.: +1-817-257-7123; fax: +1-817257-7742. E-mail address: [email protected] (B.N. Miller).
The system consists of positronium in thermal equilibrium with fluid Xenon above the critical
0168-583X/02/$ - see front matter Ó 2002 Published by Elsevier Science B.V. PII: S 0 1 6 8 - 5 8 3 X ( 0 2 ) 0 0 8 6 4 - 9
B.N. Miller, T. Reese / Nucl. Instr. and Meth. in Phys. Res. B 192 (2002) 176–179
temperature. The Ps–Xe interaction, which is dominated by exchange forces [7], is represented by a hard sphere potential, and the Xe–Xe interaction is represented by a Lennard-Jones (6,12) potential. Assume the atoms obey classical mechanics and Ps obeys quantum mechanics. Then the partition function of the hybrid system consisting of N atoms and a single qp is [10,14,15] Z b Z ¼ ð1=K3N N !Þ dRebU ðRÞ Trfeb H g; ð1Þ where U is the total interatomic potential energy, R is the set of atomic positions Ri , K is the thermal wavelength of an atom, and the Ps hamiltonian is b ¼b H p 2 =2m þ
N X
177
varying points along the Feynman chain, represented by the imaginary time indices t and t0 . Important structural information concerning the qp and the local properties of the fluid in its vicinity is provided by the qp–atom radial distribution function * + p1 N X X gfp ðrÞ ¼ dðRj ri rÞ=p =q: j¼1
i¼0
The decay rate is an observable and can be represented by the operator b k¼
N X
f ðRj rþ Þ;
ð4Þ
j¼1
V ðjRj rjÞ;
ð2Þ
j¼1
where V is the Ps–atom pair potential. The mean value of a physical observable O is then Z b ð3Þ hOi ¼ dRebU ðRÞ Trfeb H Og=ZðK3N N !Þ:
2.2. Path integral Monte Carlo For the purposes of computing statistical averages, Feynman showed that the qp can be represented by a chain of p harmonically coupled ‘‘pseudo-particles’’ with a temperature dependent force constant [17]. Each pseudo-particle also interacts with each classical atom through the reduced potential V =p. This construction transforms the calculation of the equilibrium properties of a single qp into that of p classical pseudo-particles. Statistical averages are performed by averaging over many realizations of the chain and atom configurations. The size of the chain is increased until convergence is obtained. In order to improve the convergence, we used the image approximation for the qp–atom hard sphere potential [18]. The details are discussed in one of our earlier papers [14]. 2.3. Important quantities Information concerning the degree of localization of the qp in the host fluid is provided by the mean square separation of pseudo-particles at
where f ðRj rþ Þ is the local electron density contributed by the jth atom at the positron position rþ . The mean Ps decay rate can be directly expressed in terms of gfp [10], D E Z b b k ¼ dRebU ðRÞ Trfeb H b k g=ðZK3N N !Þ ¼ ½q=ð8pa30 Þ Z Z dR0 dx f ðjR0 x=2jÞex=a0 gfp ðR0 Þ: ð5Þ With the additional assumption that the electronic charge is concentrated at the nucleus, we obtain Z D E 0 3 b k ¼ q=ðpa0 Þ dR0 e2R =a0 gfp ðR0 Þ; ð6Þ for the mean decay rate.
3. Results We performed simulations over a wide range of fluid densities at 300 and 340 K, where experimental measurements of the o-Ps decay rate in fluid Xenon have been reported [11]. The calculations were made using eight density values (q ¼ 0:017, 0.035, 0.088, 0.17, 0.3, 0.35, 0.4 and 0.7) on each isotherm. Here the scaled density q ¼ qr3 where r is the Lennard-Jones radius. In these units the density of the critical point is about 0.35.
178
B.N. Miller, T. Reese / Nucl. Instr. and Meth. in Phys. Res. B 192 (2002) 176–179
Fig. 1. RMSD at T ¼ 300 K versus (t t0 ) at three density values. The extent of localization of the Ps increases with density.
D E 2 1=2 Fig. 1 provides a plot of jrðtÞ r0 ðt0 Þj , the root mean square displacement (RMSD) versus pseudo-particle index ðt t0 Þ along the Feynman chain, at T ¼ 300 K for three density values and the free particle. The top curve is the plot of the free chain (q ¼ 0:0). As the density is increased, the curve becomes flatter, indicating that the spread of the polymer is becoming smaller, and the Ps more localized in the fluid. Self-trapping appears complete at q ¼ 0:7: Fig. 2 presents a plot of the radial distribution function versus position for T ¼ 300 K at three different densities. The local fluid distortion induced by the qp is greatest at the critical density,
Fig. 3. Comparison of the experimental and theoretical pick-off decay rate at T ¼ 300 K. The straight line is the linear extrapolation of the low density dependence.
while the qp confinement is greatest at twice this density. Fig. 3 presents a plot of the decay rate versus q for T ¼ 300 K. It displays both the theoretical computations and the experimental results. The straight line is the linear extrapolation of the low density dependence of the decay rate. The theoretical curve is close to the experimental result. Similar results were obtained on the 340 K isotherm [19].
4. Conclusions
Fig. 2. Ps–atom radial distribution function at three densities versus position at T ¼ 300 K.
The results of the RMSD computations show that as the density is increased, the qp becomes more confined. At the highest density, the distance between almost all of the pseudo-particles in the polymer representation is the same, indicating a completely self-trapped state. In the plots of the polymer-fluid radial distribution functions it is seen that the greatest number of fluid molecules are expelled from the vicinity of the Ps–atom at the critical density. This results in the sharply reduced probability of annihilation observed in experiments. Mean field theory suggests that well above the critical density the decay rate returns to the predicted linear density dependence due to collapse of the ‘‘bubble’’ [8,13]. Here we see for the
B.N. Miller, T. Reese / Nucl. Instr. and Meth. in Phys. Res. B 192 (2002) 176–179
first time that no such thing occurs. At the highest density value (q ¼ 0:7), it is noted that the large local density of fluid molecules at the edge of the bubble results in the increased probability of positron annihilation. Compared with earlier PIMC computations, the plots of the decay rate versus density now show curvature, sharing some of the attributes of the experimental curves. Both at 300 and 340 K they also show the rapid, non-linear, increase above the critical density seen in other experiments. Unfortunately there isn’t sufficient data available for Xenon to compare with the simulations at these higher densities. However, on the 300 K isotherm near the critical density, the theoretical curve doesn’t show a plateau and lies above the experimental results above the critical density. The 300 and 340 K theoretical isotherms have similar shapes: In particular, they both show an anomalous rise (bump) above the critical density which is not present in the experimental data. It will be interesting to determine if this is a consequence of either the simplified interaction potential or atomic electron density used to perform the simulations.
179
References [1] J. Hernandez, Rev. Mod. Phys. 63 (1991) 675. [2] D. Chandler, K. Leung, Annu. Rev. Phys. Chem. 45 (1994) 557. [3] G.N. Chuev, Izv. Akad. Nauk., Ser. Fiz. 61 (1997) 1770. [4] B.J. Berne, D. Thirumali, Ann. Rev. Phys. 37 (1986) 401. [5] A.F. Borghesani, M. Santini, Phys. Rev. A 45 (1992) 8803. [6] P. Hautojarvi, K. Rytsola, P. Tuovinen, P. Jauho, Phys. Lett. A 57 (1976) 175. [7] I.T. Iakubov, A.G. Khrapak, Rep. Prog. Phys. 45 (1982) 697. [8] B.N. Miller, T. Reese, Phys. Rev. A 39 (1989) 4735. [9] A.L. Nichols III, D. Chandler, Y. Singh, D.M. Richardson, J. Chem. Phys. 81 (1984) 5109. [10] B.N. Miller, Y. Fan, Phys. Rev. A 42 (1990) 2228. [11] M. Tuomisaari, K. Rytsola, P. Hautojarvi, J. Phys. B 21 (1988) 3917. [12] T.J. Murphy, C.M. Surko, J. Phys. B 23 (1990) L727. [13] T. Reese, B.N. Miller, Phys. Rev. A 42 (1990) 6068. [14] T.L. Reese, B.N. Miller, Phys. Rev. E 47 (1993) 2581. [15] G.A. Worrell, B.N. Miller, Phys. Rev. A 46 (1992) 3380. [16] B.N. Miller, T.L. Reese, G. Worrell, Phys. Rev. E 47 (1993) 4083. [17] R.P. Feynman, Statistical Mechanics, Addison-Wesley, Reading, MA, 1972. [18] P.A. Whitlock, M.H. Kalos, J. Chem. Phys. 30 (1978) 361. [19] T.L. Reese, B.N. Miller, Phys. Rev. E 64 (2001) 061201.