Computation of some thermodynamic properties of helium–neon and helium–krypton fluid mixtures using molecular dynamics simulation

Fluid Phase Equilibria 291 (2010) 117–124

Contents lists available at ScienceDirect

Fluid Phase Equilibria journal homepage: www.elsevier.com/locate/fluid

Computation of some thermodynamic properties of helium–neon and helium–krypton fluid mixtures using molecular dynamics simulation Mohsen Abbaspour a , Elaheh K. Goharshadi b,∗ , Majid Namayandeh Jorabchi b a b

Department of Chemistry, Sabzevar Tarbiat Moallem University, Sabzevar, Iran Department of Chemistry, Ferdowsi University of Mashhad, Mashhad 91779, Iran

a r t i c l e

i n f o

Article history: Received 31 October 2009 Received in revised form 15 December 2009 Accepted 19 December 2009 Available online 4 January 2010 Keywords: Potential energy function Molecular dynamics simulation Quantum corrections Three-body interactions

a b s t r a c t We have performed molecular dynamics simulations to obtain internal energy and pressure of helium–neon and helium–krypton mixtures at different densities using accurate recently two-body ab initio potentials supplemented by quantum corrections following the Feynman–Hibbs approach. The significance of this work is that the three-body expression of Wang and Sadus [22] was used to improve prediction of the pressures and internal energies of helium + krypton and helium + neon mixtures without requiring an expensive three-body calculation. Our results show a good agreement with the corresponding experimental data. © 2009 Elsevier B.V. All rights reserved.

1. Introduction Simulation of binary mixtures has been lagging behind that of pure fluids [1]. There exist very few simulation data for the properties of fluid mixtures since the experimental data for thermodynamics of fluid mixtures are scarce or nonexistent. It is well known that quantum-mechanical effects in bulk properties of dilute gases can be important when dealing with systems of small reduced masses at low temperatures [2]. The rare gases, except helium and neon, show the classical behavior in their liquid states. The quantum effects of helium give rise to superfluidity. Although neon does not display such dramatic effects, the quantum effects cannot be ignored [3,4]. Two approaches have been proposed to consider the quantum effects, the so-called Wigner–Kirkwood (WK) [5–7] and Feynman–Hibbs (FH) [8] potentials. The WK potential arises from an expansion in powers of  of the partition function which has been used in the literature to estimate the quantum corrections on different properties. The FH potential used in the present work is based on the Feynman–Hibbs variational estimate of the quantum partition function [8] which depends on temperature and is easy to implement in a standard molecular dynamics (MD) or Monte Carlo (MC) simulation code [3]. It is well known that a simple pair potential, though giving the essential features of the thermodynamic properties for noble

∗ Corresponding author. Tel.: +98 511 8797022; fax: +98 511 8796416. E-mail address: [email protected] (E.K. Goharshadi). 0378-3812/$ – see front matter © 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.fluid.2009.12.024

gases and their mixtures, is not sufficient for a quantitative description and many-body effects most notably three-body interactions have to be taken into account. In fact, the structural and thermodynamic properties of fluids are influenced significantly by three-body interactions. These important three-body effects have previously remained undetected because earlier work was confined to effective potentials. Study of three-body forces began in 1943 when Axilrod and Teller [9] used third-order perturbation theory to calculate the three-body triple–dipole dispersion energy for atoms with a spherical charge distribution. Nasrabad et al. [10] predicted the thermophysical properties of neon, argon, and the binary mixtures of neon–argon and argon–krypton by Monte Carlo simulation using ab initio potentials enhanced by Axilord–Teller triple–dipole interactions. Calculation of three-body interactions is computationally prohibitive for routine applications. In the worst case, the evaluation of three-body interactions for N atoms involves N3 calculations compared with N2 calculations for two-body potentials. To remedy this problem, Marcelli and Sadus [11] proposed a simple formula to obtain the three-body contribution to the potential energy from the two-body interaction potential energy. Goharshadi and Abbaspour [12] performed molecular dynamics simulation to obtain the internal energy and pressure of argon, krypton, and xenon using the two-body HFD-like potential plus the three-body interactions using the Marcelli and Sadus expression [11]. Their results were in a good overall agreement with the experiment.

118

M. Abbaspour et al. / Fluid Phase Equilibria 291 (2010) 117–124

Table 1 Intermolecular potential parameters for different like and unlike binary interactions (Eqs. (1)–(4)). Parameter A (J) a1 (m−1 ) a−1 (m) a2 (m−2 ) a−2 (m2 ) d1 * d2 d3 * b (m−1 ) C6 (J m6 ) C8 (J m8 ) C10 (J m10 ) C12 (J m12 ) C14 (J m14 ) C16 (J m16 ) f2n *(r)

He–He [16]

Ne–Ne [17] −17

4.2382 × 10 −3.8120 × 1010 −2.2841 × 10−11 −1.6414 × 1019 3.8793 × 10−22 1.6712 × 10−3 (m−1 ) 1.7664 × 100 3.8480 × 1010 1.3987 × 10−79 3.7845 × 10−99 1.3788 × 10−118 6.8634 × 10−138 4.5001 × 10−157 3.7513 × 10−176 1

Parameter −16

5.5611 × 10 −4.2870 × 1010 −5.3462 × 10−11 −3.3384 × 1018 5.0177 × 10−21 C8 * 3.3694 × 1010

Kr–Kr [18]

Parameter −21

2.7296 × 10 3.5710 × 10−10 −15.2071 4071.7500 2.5721

ε (J)  (m) A* ˛* C6 * 6.3004 C10 *

4.9244 × 10 6.0813 × 10−79 2.2755 × 10−98 1.0909 × 10−117 6.6998 × 10−137 5.2718 × 10−156 5.3145 × 10−175

The values of coefficients A, a1 , a2 , a−1 , a−2 , b, d1 *, d2 , d3 *, f2n *(r), and C6 to C16 have been given in Table 1. For neon–neon interaction, we have used the interatomic potential energy of Hellmann et al. [17] determined from quantummechanical ab initio calculations and described with the following analytical representation: U(r) = A exp(a1 r + a2 r 2 + a−1 r −1 + a−2 r −2 ) −

8  C2n n=3

2. Theory 2.1. Intermolecular potentials For the site–site interactions of helium–neon and helium–krypton mixtures, we have used the following binary interactions. For helium–helium interaction, we have used a new interaction potential of Hellmann et al. [16] obtained from highly accurate ab initio calculations at a large number of helium–helium separations given by: U(r) = A exp(a1 r + a2 r 2 + a−1 r −1 + a−2 r −2 + d1∗ sin(d2 r + d3∗ )) −

n=3



r 2n

r

2n k  (br)

× 1 − exp(−br)

k=0

k!



(1)

2n k  (br)

1 − exp(−br)

k=0



(2)

k!

The values of coefficients A, a1 , a2 , a−1 , a−2 , b, and C6 to C16 have been given in Table 1. For krypton–krypton interaction, we have used the potential of Goharshadi et al. [18] which has been obtained from the inversion of the viscosity collision integrals at zero pressure and fitted to a HFD-like potential. It has the following form: U(r) = A∗ exp(˛∗ y) − ε



C6∗ y6

+

C8∗ y8

+

∗ C10



(3)

y10

In this equation, ε is the well depth of the potential and y = r/ ( is the distance at which the intermolecular potential is zero). The values of parameters , ε, A*, ˛*, and C6 * to C10 * have been given in Table 1. For helium–neon and helium–krypton interactions, we have used the interaction potential energies of Partridge et al. [19] obtained using high-quality ab initio calculations with the following analytical representation: r U(r) = exp a − 



3.0850 × 10−11 5.7070 1.2529 × 10−77

2.1759 × 10 5.0010 5.8205 × 10−79

 (m) a* C6 (J m6 )

−8.8520



C2n ∗ f2n (r) 2n

He–Kr [19] −11

10

Wang and Sadus [13] used the Gibbs–Duhem Monte Carlo simulations for predicting the vapor-liquid phase coexistence of binary argon–krypton mixtures. They employed accurate two-body potentials in addition to contributions of three-body interactions from the Marcelli and Sadus formula [11]. They showed that the addition of three-body interactions provides a near perfect agreement with experiment for the vapor branch of the existence curve while simultaneously improving the agreement with experiment on the liquid branch. Recently, Wang and Sadus [14] reported molecular dynamics simulation data for two-body and three-body interactions in noble gases at densities covering the gas, liquid, and solid phases. They proposed a new expression for incorporating three-body interactions which is valid for both the normal fluid region and solid phases of noble gases and possibly other atoms. Moosavi and Goharshadi [15] performed molecular dynamics (MD) to compute the pressure and internal energy of binary mixtures of argon–neon, argon–krypton, and argon–xenon at different temperatures and compositions using two-body Hartree–Fock dispersion-like (HFD-like), total (two-body + three-body) HFD-like, and Lennard–Jones (LJ) potentials. Their results were in good overall agreement with the experimental data. The purpose of the present paper is to perform MD simulations to obtain the internal energy and pressure of helium–neon and helium–krypton mixtures at different densities using twobody potential supplemented by quantum corrections following the Feynman–Hibbs approach. The contribution of three-body interactions using the Wang and Sadus relationship [22] between two-body and three-body interactions has been also incorporated in the simulations.

8 

He–Ne [19]

∗

8  



−

2n  r  k (r/)

1 − exp −

n=3



k=0

k!



C2n r 2n

(4)

The values of coefficients a*, , and C6 for He–Ne and He–Kr have been given in Table 1. The values of the parameters C8 and C10 are approximated using the relations: 3 C8 = 2 C6 C10 C6

   ra4 r4  2 +  b2 and ra

rb

   6  4  4  ra6 rb r 21 ra  2  b2 = 2 + 2 + 2 2 ra

rb

5

ra

rb

(5)

which can be derived from the general formulation of Starkshall and Gordon [20] for the interaction of atoms a and b. The quantities rn  are the radial expectation values for the atoms and are tabulated

M. Abbaspour et al. / Fluid Phase Equilibria 291 (2010) 117–124 Table 2 The values coefficients.

of

three-body

 (J m9 )a

He–He–He Ne–Ne–Ne Kr–Kr–Kr

2.0939 × 10−110 1.7047 × 10−109 2.2039 × 10−107

a

where  is the nonadditive coefficient,  = N/V is the number density,  is the distance at which the intermolecular potential is zero, and ε is the well depth of the potential. To apply Eq. (7) for mixtures, Wang and Sadus [22] proposed a simple one-fluid approximation:

nonadditivity

Interaction

U3B = −

 C 3 2n C2n−2

0.85U2B

Ref. [10].

C2n−4

6 ε11 11



(6) UT = U2B

1−

Wang and Sadus [14] showed that there is a simple and accurate relationship between the two-body (U2B ) and three-body (U3B ) potential energies of a fluid as: U3B = −

0.85U2B ε 6

(x12 111 + x1 x2 112 + x1 x2 221 + x22 222 )

(8)

where xi is the molar fraction of component i and ijk is the threebody potential for the three different components i, j, and k. The overall interaction potential, two-body plus three-body, (UT ) is as follows:

by Desclaux [21]. Higher order (2n > 10) dispersion coefficients are obtained from recursion relations [19]: C2n+2 =

119

0.85  6 ε11 11



(x12 111 + x1 x2 112 + x1 x2 221 + x22 222 ) (9)

The significance of this equation is that it allows us to use twobody potentials to predict accurately the properties of real fluids without incurring the computational cost of three-body calcula-

(7)

Table 3 Our simulated values of pressure and internal energy of He–Kr mixture using different potentials at T = 223.15 K.  (103 kg/m3 )

P (MPa)

E (kJ/mol)

U2B

UQ2B

UT

UQT

U2B

UQ2B

UT

UQT

xHe = 0.2 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.60 0.70 0.80 0.90 1.00

1.35 2.61 3.81 4.92 5.98 6.88 7.77 8.67 8.93 9.68 11.10 11.60 13.40 13.80 15.20

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.03 0.07 0.14 0.23 0.35 0.53 0.63 0.56 0.92 1.04 0.88 1.44 2.20 2.18 2.26

1.35 2.61 3.82 4.94 5.90 6.82 7.81 8.49 9.20 9.49 11.40 12.20 14.10 14.60 15.20

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.03 0.07 0.16 0.24 0.32 0.57 0.79 0.44 0.58 1.08 1.00 1.32 2.10 1.86 2.23

1.34 2.61 3.85 4.97 6.07 6.99 8.15 8.99 9.88 10.68 12.79 13.81 15.83 18.98 20.72

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.03 0.06 0.17 0.26 0.26 0.44 0.74 0.53 0.74 0.80 0.92 1.36 1.71 1.86 2.51

1.34 2.57 3.79 4.87 5.86 6.74 7.81 8.60 9.29 9.63 11.49 12.71 14.23 14.88 17.02

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.03 0.07 0.16 0.24 0.31 0.49 0.75 0.60 0.45 0.91 1.03 1.08 1.58 2.55 2.63

2.69 2.56 2.41 2.29 2.13 2.00 1.85 1.74 1.59 1.40 1.10 0.78 0.54 0.23 −0.02

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.04 0.04 0.01 0.03 0.04 0.04 0.03 0.04 0.10 0.05 0.04 0.08 0.04 0.04

2.69 2.55 2.41 2.29 2.13 2.00 1.85 1.74 1.59 1.39 1.10 0.79 0.56 0.25 −0.03

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.04 0.04 0.02 0.02 0.04 0.07 0.04 0.02 0.09 0.04 0.03 0.04 0.06 0.05

2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.77 2.78 2.78 2.78

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.03 0.04 0.09 0.03 0.02 0.06 0.02 0.03 0.05 0.03 0.07 0.05 0.05 0.04

2.69 2.56 2.41 2.29 2.13 2.02 1.87 1.73 1.60 1.42 1.14 0.82 0.59 0.30 0.08

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.04 0.05 0.01 0.03 0.02 0.06 0.03 0.06 0.07 0.03 0.04 0.05 0.05 0.05

xHe = 0.4 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.60 0.70 0.80 0.90 1.00

1.78 3.51 5.16 6.86 8.43 9.91 11.20 13.10 15.00 16.30 19.20 22.70 26.00 30.00 35.90

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.03 0.10 0.20 0.28 0.31 0.54 0.48 0.46 1.41 0.92 2.03 1.63 2.25 2.04 3.48

1.78 3.52 5.12 6.81 8.42 9.97 11.30 13.20 14.50 16.30 19.00 23.30 27.60 31.10 38.10

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.03 0.10 0.19 0.21 0.43 0.54 0.55 0.90 0.81 1.18 1.06 1.55 2.89 2.13 2.80

1.79 3.54 5.25 6.92 8.54 10.22 11.91 13.70 15.59 16.69 20.78 25.18 29.50 34.81 41.43

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.03 0.08 0.20 0.27 0.24 0.59 0.42 0.76 0.82 0.84 1.22 1.66 2.32 2.49 3.11

1.78 3.50 5.11 6.87 8.36 9.88 11.51 13.20 14.79 16.09 18.99 22.88 27.89 32.61 36.53

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.03 0.10 0.20 0.17 0.37 0.50 0.60 0.71 0.92 0.86 1.41 2.03 2.03 2.08 2.43

2.70 2.61 2.50 2.39 2.26 2.15 2.02 1.91 1.81 1.64 1.44 1.26 0.98 0.78 0.57

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.02 0.03 0.03 0.02 0.03 0.01 0.03 0.05 0.04 0.08 0.04 0.08 0.03 0.04

2.70 2.61 2.50 2.38 2.24 2.16 2.02 1.92 1.80 1.67 1.42 1.26 1.01 0.79 0.59

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.03 0.03 0.02 0.03 0.04 0.04 0.03 0.05 0.04 0.05 0.03 0.06 0.05 0.05

2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.77 2.77 2.78 2.77

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.02 0.02 0.03 0.02 0.03 0.02 0.03 0.02 0.02 0.06 0.04 0.04 0.07 0.05

2.70 2.61 2.50 2.39 2.26 2.16 2.04 1.92 1.82 1.69 1.46 1.27 1.07 0.83 0.66

± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.03 0.04 0.02 0.03 0.03 0.03 0.02 0.02 0.04 0.06 0.06 0.05 0.05 0.05

xHe = 0.7 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.60 0.70 0.80 0.90

3.38 6.84 10.30 13.90 17.80 21.80 26.50 31.00 36.00 41.10 53.20 67.40 85.60 109.00

± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.05 0.14 0.35 0.37 0.55 0.87 1.18 1.40 1.57 2.00 2.19 2.49 4.13 4.16

3.38 6.86 10.40 14.00 18.00 22.40 26.10 30.70 36.00 41.40 54.70 67.80 87.60 110.00

± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.05 0.13 0.27 0.44 0.72 0.69 0.95 0.88 0.99 1.73 1.84 2.74 3.08 5.65

3.40 6.87 10.49 14.08 18.17 22.46 26.75 31.44 36.24 41.93 54.92 71.01 89.30 112.89

± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.06 0.14 0.30 0.35 0.45 0.77 1.00 1.52 1.29 1.66 1.96 3.70 2.95 3.75

3.39 6.85 10.39 13.98 17.97 21.66 26.25 30.64 36.14 41.73 52.72 68.11 86.90 109.89

± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.05 0.12 0.26 0.56 0.71 0.72 1.13 1.19 1.16 1.73 2.21 3.88 4.07 5.37

2.72 2.65 2.58 2.52 2.45 2.38 2.31 2.24 2.18 2.09 1.96 1.87 1.74 1.68

± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.02 0.01 0.02 0.03 0.01 0.05 0.03 0.04 0.06 0.05 0.01 0.03 0.04

2.72 2.65 2.57 2.50 2.45 2.39 2.29 2.23 2.18 2.12 1.96 1.87 1.78 1.69

± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.01 0.02 0.03 0.03 0.03 0.03 0.03 0.04 0.04 0.03 0.01 0.06

2.78 2.78 2.78 2.77 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.78 2.77 2.78

± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.02 0.03 0.03 0.02 0.02 0.05 0.02 0.03 0.04 0.05 0.04 0.06

2.72 2.65 2.57 2.52 2.45 2.40 2.30 2.24 2.19 2.13 2.00 1.90 1.81 1.72

± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.02 0.04 0.01 0.02 0.03 0.04 0.02 0.03 0.04 0.07 0.05 0.08

120

M. Abbaspour et al. / Fluid Phase Equilibria 291 (2010) 117–124

For the energy-like parameters :

tions. In this work, we have used this equation in our simulations along with the two-body potentials. We have also applied the following formula proposed by Smit et al. [23] for calculating the configuration pressure when the total potential (two-body plus three-body interaction potential) is used. Pconf = −



122 =

1 dU2B (rij ) rij 3V drij

 

−

6 ε11 11

i
2/92 U2B 6 ε11 11

(x12 111 +x1 x2 112 +x1 x2 221 +x22 222 ) ×

dU2B (rij ) drij

rij

 (x12 111

+ x1 x2 112 + x1 x2 221 +

(11)

3

1 22

(12)

112 =

21 + 2 3

(13)

122 =

1 + 22 3

(14)

(10)

x22 222 )

12 2

and for :

 2/9

3





i
+



112 =

In these Eqs. (1) stands for helium and 2 stands for neon or krypton. Table 2 shows the values of the nonadditivity coefficients for the pure components [10]. In order to incorporate the quantum effects, we have also used the Feynman–Hibbs potential [8] with two-body and total poten-

where the angle bracket represent ensemble averages. In a binary mixture, when the three-body potential is incorporated into a simulation program, the cross intermolecular parameters can be estimated from the values of pure components.

Table 4 Our simulated values of pressure and internal energy of He–Ne mixture using different potentials at T = 173.15 K.  (103 kg/m3 )

P (MPa)

E (kJ/mol)

U2B

UQ2B

UT

UQT

U2B

UQ2B

UT

UQT

xHe = 0.2 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.53 0.57 0.60 0.65 0.67 0.70 0.73 0.75

4.25 8.48 12.80 17.20 22.50 26.90 32.50 39.10 45.10 53.10 57.00 65.00 71.00 80.50 85.50 94.60 103.50 109.00

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.06 0.21 0.4 0.68 0.91 0.84 1.52 1.24 2.37 2.24 3.921 3.68 3.77 3.85 4.10 3.82 6.23 5.395

4.26 8.51 12.80 17.20 21.70 27.10 32.00 39.20 44.50 54.20 58.10 66.00 71.30 83.80 89.00 98.50 106.00 112.00

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.08 0.21 0.3 0.49 0.75 0.88 1.09 1.48 2.56 2.49 2.69 3.39 3.66 3.82 4.84 3.95 5.21 5.15

4.27 8.57 13.09 17.69 22.59 27.78 33.88 39.58 46.08 54.08 60.50 68.10 74.58 87.00 91.00 98.78 109.00 117.00

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.05 0.22 0.41 0.50 0.95 1.24 1.63 1.32 2.10 1.80 2.94 3.86 4.13 3.88 3.50 4.20 4.26 5.66

4.26 8.52 12.89 16.99 22.19 27.18 32.48 38.68 44.88 53.88 58.00 66.20 73.78 83.20 88.40 97.88 106.00 112.00

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.08 0.19 0.26 0.62 0.75 1.14 1.88 1.64 1.45 1.21 2.76 3.57 2.86 5.50 4.16 3.83 4.34 4.50

2.07 1.98 1.89 1.79 1.71 1.62 1.54 1.45 1.37 1.29 1.23 1.18 1.14 1.06 1.03 0.98 0.95 0.93

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.02 0.01 0.03 0.02 0.03 0.02 0.02 0.04 0.08 0.06 0.05 0.05 0.05 0.04 0.08 0.05

2.07 1.97 1.89 1.79 1.71 1.62 1.51 1.44 1.35 1.29 1.24 1.19 1.13 1.07 1.04 1.00 0.96 0.93

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.01 0.02 0.01 0.01 0.02 0.02 0.02 0.06 0.04 0.05 0.02 0.04 0.07 0.06 0.06 0.05

2.15 2.15 2.15 2.15 2.16 2.15 2.16 2.15 2.15 2.15 2.16 2.15 2.16 2.16 2.15 2.16 2.16 2.16

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.02 0.01 0.04 0.02 0.03 0.02 0.03 0.03 0.04 0.07 0.06 0.06 0.04 0.04 0.06 0.06

2.07 1.97 1.88 1.80 1.71 1.62 1.53 1.45 1.38 1.29 1.25 1.20 1.15 1.08 1.05 1.01 0.98 0.95

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.01 0.02 0.02 0.01 0.03 0.02 0.02 0.04 0.04 0.05 0.02 0.08 0.05 0.03 0.04 0.04

xHe = 0.4 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.43 0.45 0.47 0.50 0.53 0.55 0.57 0.60

5.33 10.80 16.70 22.60 30.20 37.40 46.70 55.70 62.80 67.10 72.00 80.30 89.80 95.70 102.00 111.00

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.10 0.22 0.44 0.85 0.86 1.14 1.50 2.31 2.41 1.94 3.50 3.24 2.95 3.53 3.76 4.10

5.36 10.80 16.60 22.90 30.60 38.30 45.70 56.50 63.30 67.80 73.30 82.20 88.30 97.10 105.00 116.00

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.08 0.27 0.43 0.73 0.94 1.31 1.53 2.57 2.36 2.25 2.93 2.07 3.83 3.73 5.74 3.97

5.36 10.89 16.89 23.29 30.28 38.68 46.98 56.98 64.40 66.98 74.30 81.98 91.70 100.00 106.00 116.97

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.08 0.19 0.45 0.81 0.91 1.44 1.56 1.89 2.05 2.93 2.25 2.61 3.54 3.86 3.43 3.62

5.34 10.89 16.79 22.79 30.08 37.28 46.38 56.18 63.10 68.88 72.60 79.08 90.30 97.90 104.00 114.97

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.09 0.23 0.50 0.70 1.01 1.51 1.44 1.49 1.77 2.52 2.02 2.87 2.93 3.53 4.34 4.06

2.08 2.01 1.93 1.86 1.79 1.72 1.66 1.59 1.56 1.53 1.50 1.47 1.44 1.42 1.39 1.36

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.01 0.01 0.02 0.01 0.01 0.02 0.01 0.01 0.04 0.03 0.03 0.03 0.02 0.02

2.08 2.01 1.93 1.86 1.80 1.72 1.66 1.59 1.57 1.53 1.51 1.48 1.43 1.42 1.41 1.38

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.01 0.02 0.03 0.01 0.02 0.04 0.03 0.02 0.02 0.02 0.02 0.03 0.02 0.02

2.15 2.16 2.15 2.16 2.16 2.16 2.16 2.15 2.15 2.16 2.16 2.15 2.16 2.15 2.16 2.16

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.01 0.02 0.02 0.03 0.02 0.02 0.02 0.04 0.02 0.03 0.02 0.02 0.03 0.03

2.08 2.01 1.93 1.86 1.79 1.73 1.66 1.59 1.56 1.55 1.51 1.47 1.45 1.43 1.41 1.40

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.01 0.02 0.02 0.02 0.02 0.02 0.01 0.02 0.01 0.03 0.02 0.03 0.05 0.04

xHe = 0.7 0.05 0.10 0.15 0.20 0.25 0.27 0.30 0.33 0.35 0.37 0.40 0.43

8.51 18.20 28.70 41.60 56.40 62.50 73.70 85.90 94.40 104.00 120.00 139.00

± ± ± ± ± ± ± ± ± ± ± ±

0.17 0.36 0.76 1.29 1.79 1.92 2.57 2.46 2.86 2.92 3.26 3.78

8.52 18.10 29.00 41.80 56.60 62.80 75.50 86.90 97.80 107.00 125.00 143.00

± ± ± ± ± ± ± ± ± ± ± ±

0.14 0.46 1.00 1.14 1.85 1.89 2.16 3.26 3.14 3.20 3.72 4.09

8.57 18.09 28.99 41.49 56.99 64.10 75.28 86.60 95.48 106.00 122.98 141.00

± ± ± ± ± ± ± ± ± ± ± ±

0.18 0.40 0.71 1.21 1.57 2.30 2.46 2.91 3.54 2.41 3.17 4.36

8.51 17.99 28.89 41.59 56.69 63.70 75.68 88.70 97.88 107.00 122.98 144.00

± ± ± ± ± ± ± ± ± ± ± ±

0.15 0.42 0.73 1.05 1.55 2.25 1.88 2.55 3.72 3.36 3.14 3.55

2.10 2.05 1.98 1.95 1.90 1.88 1.85 1.85 1.83 1.81 1.80 1.80

± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.02 0.01 0.01 0.01 0.02 0.01 0.03 0.03 0.02 0.03

2.10 2.05 2.00 1.95 1.90 1.88 1.87 1.84 1.84 1.82 1.82 1.80

± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.02 0.02 0.01 0.04 0.02 0.03 0.03 0.03 0.02 0.03 0.02

2.16 2.15 2.15 2.15 2.15 2.16 2.16 2.16 2.15 2.15 2.16 2.16

± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.02 0.02 0.02 0.02 0.03 0.02 0.04 0.02 0.03 0.05

2.10 2.05 1.99 1.95 1.90 1.89 1.88 1.86 1.85 1.83 1.82 1.81

± ± ± ± ± ± ± ± ± ± ± ±

0.01 0.01 0.02 0.02 0.02 0.03 0.02 0.02 0.04 0.03 0.03 0.02

M. Abbaspour et al. / Fluid Phase Equilibria 291 (2010) 117–124

Fig. 1. Comparison between experimental data (Eq. (17)) and our simulated results of pressure for He–Kr mixture using the two-body potential (U2B ) and final potential (UQT ) at 223.15 K at three mole fractions: (a) x(He) = 0.2, (b) x(He) = 0.4, (c) x(He) = 0.7. Table 5 Virial coefficients of He–Ne and He–Kr mixtures at 173.15 and 223.15 K, respectively. Coefficient

He–Nea

He–Kra

B11 (10−6 m3 /mol) B12 (10−6 m3 /mol) B22 (10−6 m3 /mol) C111 (10−6 m3 /mol)2 C222 (10−6 m3 /mol)2

11.91 10.91 6.45 194.40 199.00

11.59 19.74 −93.05 176.00 3011.00

a

Ref. [25].

121

Fig. 2. Comparison between experimental data (Eq. (17)) and our simulated results of pressure for He–Ne mixture using the two-body potential (U2B ) and final potential (UQT ) at 173.15 K at three mole fractions: (a) x(He) = 0.2, (b) x(He) = 0.4, (c) x(He) = 0.7.

tials (UQ2B and UQT ):



U  (r) ˇh ¯  U2B (r) + 2 2B UQ2B (r) = U2B (r) + 24 r 2



2 U  (r) ˇh ¯ UQT (r) = UT (r) + UT (r) + 2 T 24 r

 (15)

 (16)

122

M. Abbaspour et al. / Fluid Phase Equilibria 291 (2010) 117–124

where ˇ = 1/kB T,  is the reduced mass and  = h/2 (h is the Planck’s constant). The prime and the double prime are the first and second r derivatives, respectively. As Eqs. (15) and (16) show, the Feynman–Hibbs quantum correction term depends on mass and temperature.

2.2. Simulation details The MD simulations using MOLDY software [24] have been performed for a system of 1000 atoms in a cubic box and the conventional periodic boundary condition has been applied. We have used the NVT ensemble using a Nose–Hoover thermostat for atoms interacting via the two-body potential, U2B , (Eqs. (1)–(4)), quantum two-body potential, UQ2B , (Eqs. (1)–(4) with Eq. (15)), total (two-body plus three-body) potential, UT , (Eqs. (1)–(4) and (9)), and quantum total potential (final potential), UQT , (Eqs. (1)–(4), (16) and (9)). The size of time steps, t, and the number of time steps, nt , and the cut-off radius, rc , have been chosen as 0.001 ps, 500,000, and 3, respectively. The long-range correction terms have been evaluated using MOLDY to recover the contribution of the longrange cut-off of the intermolecular potential on the pressure and energy.

3. Results and discussion We have performed MD simulations to obtain the pressure and internal energy of helium–krypton and helium–neon mixtures using two-body potential (Eqs. (1)–(4)). We have also incorporated the quantum corrections, the FH potential (Eqs. (15) and (16)). The simple expression of Wang and Sadus [22] has also been used for incorporating the three-body interactions in our simulations with the two-body potential interactions (Eq. (9)). Our results for pressure and internal energy of He–Kr mixture at 223.15 K and He–Ne mixture at 173.15 K in the NVT ensemble and three different molar fractions and different densities using different potentials have been given in Tables 3 and 4, respectively. In these tables 2B, Q2B, T, and QT denote two-body, quantum twobody (FH), total (two-body plus three-body), and quantum total potential (final potential), respectively. As Tables 3 and 4 show, the quantum effects using the FH potential (UQ2B ) tend to increase the pressure results for two-body potential (U2B ) especially at higher densities. The similar results were also achieved by Tchouar et al. [3,4]. They computed the properties of neon and helium using the Lennard–Jones potential with the FH quantum correction. The data in Tables 3 and 4 also indicate that the three-body interactions via the expression of Wang and Sadus [22] included in the total potential (UT ) tend to increase the pressure and energy results of He–Kr and He–Ne mixtures for two-body potential (U2B ) especially at higher densities. The similar results were also obtained by Moosavi and Goharshadi [15] for some noble gas mixtures. These tables also show that the increasing of the simulation data for the mixtures due to using the three-body effects of Wang and Sadus is bigger than those of FH quantum effects. Tables 3 and 4 also show that the total corrections (quantum and three-body potential) are more important at high densities rather than low densities. To the best of our knowledge, there is very little corresponding experimental data for the pressure of these mixtures to which we can compare our results. In order to compare our simulated results with the experimental data, we have used the virial expansion for low and moderate pressures (up to 120 MPa) and truncated after

Fig. 3. Comparison between our simulated results of energy for He–Kr mixture using the two-body potential (U2B ) with the final potential (UQT ) at 223.15 K at three mole fractions: (a) x(He) = 0.2, (b) x(He) = 0.4, (c) x(He) = 0.7.

the third term [25]: PV = 1 + B + C2 RT

(17)

with B = x12 B11 + 2x1 x2 B12 + x22 B22

(18)

C = x13 C111 + 3x12 x2 C112 + 3x22 x1 C122 + x23 C222

(19)

where xi is the mole fraction of component in the mixture and Bii and Ciii are the second and third virial coefficients, respectively, for

M. Abbaspour et al. / Fluid Phase Equilibria 291 (2010) 117–124

123

The experimental values of parameters B11 , B12 , B22 , C111 , and C222 for He–Kr mixture at 223.15 K and He–Ne mixture at 173.15 K presented in Table 5 [25]. We have compared our simulated pressure for He–Kr mixture at 223.15 K and He–Ne mixture at 173.15 K using NVT ensemble at three different mole fractions and different densities using the two-body potential (U2B ) and final potential (UQT ) with the experiment (Eq. (17)) in Figs. 1 and 2, respectively. As these figures show our results are in good agreement with the experiment and the improvement in agreement is mostly at concentrations high in Kr (i.e., x(He) = 0.2 in Fig. 1) and high densities. This makes sense because three-body effects would be expected to make a greater contribution in the larger Kr atoms compared with He. Similarly, three-body effects become more important at high densities rather than low densities. We have also compared our simulated internal energy results for He–Kr mixture at 223.15 K and He–Ne mixture at 173.15 K and three molar fractions and different densities using the two-body potential (U2B ) with final potential (UQT ) in Figs. 3 and 4, respectively. To the best of our knowledge, there is no corresponding experimental data for the internal energy of these mixtures to which we can compare our results. Figs. 3 and 4 show that there is small deviation between the energy results using the U2B and UQT at low densities but at higher densities the UQT increases the internal energy values. The reason is that the role of the intermolecular interactions in internal energy becomes more important at higher densities and hence different potentials give different values. 4. Conclusions We have performed molecular dynamics simulations to obtain internal energy and pressure of helium–krypton and helium–neon mixtures at different densities using two-body potentials supplemented by quantum corrections following the Feynman–Hibbs approach. The significance of this work is that the threebody expression of Wang and Sadus [22] can be used to improve the prediction of the pressures of helium + krypton mixtures without requiring an expensive three-body calculation. Our results show that including the FH quantum and Wang and Sadus three-body corrections to the two-body potential increases the values of internal energy and pressure and gives better agreement with the experiment relative to the two-body potential (U2B ) especially at high densities. References

Fig. 4. Comparison between our simulated results of energy for He–Ne mixture using the two-body potential (U2B ) with the final potential (UQT ) at 173.15 K at three mole fractions: (a) x(He) = 0.2, (b) x(He) = 0.4, (c) x(He) = 0.7.

a pure component i. In these equations, 1 stands for helium and 2 stands for neon or krypton. The parameters C112 and C122 computed using:



C112 =

3

2 C C111 222

(20)

2 C111 C222

(21)

 C122 =

3

[1] M.F. Pas, B.J. Zwolinski, Mol. Phys. 73 (1991) 483–494. [2] R.J. Munn, F.J. Smith, E.A. Mason, L. Monchick, J. Chem. Phys. 42 (1965) 537–539. [3] N. Tchouar, F. Ould-Kaddour, D.J. Levesque, J. Chem. Phys. 121 (2004) 7326–7331. [4] N. Tchouar, M. Benyettou, F.J. Ould-Kaddour, J. Mol. Liq. 122 (2005) 69–73. [5] E.P. Wigner, Phys. Rev. 40 (1932) 749–759. [6] J.G. Kirkwood, Phys. Rev. 44 (1933) 31–37. [7] J.P. Hansen, J.J. Weis, Phys. Rev. 188 (1969) 314–318. [8] R.P. Feynman, A. Hibbs, Quantum Mechanics and Path Integral, McGraw-Hill, New York, 1965. [9] B.M. Axilrod, E. Teller, J. Chem. Phys. 11 (1943) 299–300. [10] A.E. Nasrabad, R. Laghaei, U.K. Deiters, J. Chem. Phys. 121 (2004) 6423– 6434. [11] G. Marcelli, R.J. Sadus, J. Chem. Phys. 112 (2000) 6382–6385. [12] E.K. Goharshadi, M. Abbaspour, J. Chem. Theory Comput. 2 (2006) 920– 926. [13] L. Wang, R.J. Sadus, J. Chem. Phys. 125 (2006) 74503. [14] L. Wang, R.J. Sadus, J. Chem. Phys. 125 (2006) 144509. [15] M. Moosavi, E.K. Goharshadi, Fluid Phase Equilib. 274 (2008) 51–58. [16] R. Hellmann, E. Bich, E. Vogel, Mol. Phys. 105 (2007) 3013–3023. [17] R. Hellmann, E. Bich, E. Vogel, Mol. Phys. 106 (2008) 133–140. [18] E.K. Goharshadi, M.J. Alahmadi, B. Najafi, Can. J. Chem. 81 (2003) 866–871. [19] H. Partridge, J.R. Stallcop, E. Levin, J. Chem. Phys. 115 (2001) 6471–6488.

124 [20] [21] [22] [23]

M. Abbaspour et al. / Fluid Phase Equilibria 291 (2010) 117–124 G. Starkshall, R.G. Gordon, J. Chem. Phys. 56 (1972) 2801–2806. J.P. Desclaux, At. Data Nucl. Data Tables 12 (1973) 311–406. L. Wang, R.J. Sadus, J. Chem. Phys. 125 (2006) 074503. B. Smit, T. Hauschild, M. Prausnitz, Mol. Phys. 77 (1992) 1021–1031.

[24] The MOLDY program was coded by K. Refson and can be downloaded from the internet at http://www.earth.ox.ac.uk/%7Ekeith/moldy.html. [25] J.H. Dymond, E.B. Smith, The Virial Coefficients of Pure Gases and Mixtures, Oxford University Press, Oxford, 1980.