The First Three Coefficients in the High Temperature Series

Jul 11, 2013 - The first three coefficients of the high temperature series expansion (HTSE) of the Helmholtz free energy for a number of simple potent...
0 downloads 0 Views 1MB Size
Article pubs.acs.org/JPCB

The First Three Coefficients in the High Temperature Series Expansion of Free Energy for Simple Potential Models with Hard-Sphere Cores and Continuous Tails Shiqi Zhou* School of Physics and Electronics, Central South University, Changsha, Hunan, 410083, People’s Republic of China

J. R. Solana Departamento de Física Aplicada, Universidad de Cantabria, 39005 Santander, España ABSTRACT: The first three coefficients of the high temperature series expansion (HTSE) of the Helmholtz free energy for a number of simple potential models with hard-sphere cores plus continuous tails are obtained for the first time from Monte Carlo simulations. The potential models considered include Square-well, Sutherland, attractive Yukawa, and triangle-well with different potential ranges, as well as a model potential qualitatively resembling the depletion potential in colloidal dispersions. The simulation data are used to evaluate performance of a recent coupling parameter series expansion (CPSE) in calculating for these coefficients, and a traditional macroscopic compressibility approximation (MCA) for the second-order coefficient only. A comprehensive comparison based on these coefficients from the two theoretical approaches and simulations enables one to conclude that (i) unlike one common experience that the widely used MCA usually underestimates the second-order coefficient, the MCA can both overestimate and underestimate the second-order coefficient, and worsens as the range of the potential decreases; and (ii) in contrast, the CPSE not only reproduce the trends in the density dependence of the perturbation coefficients, even the third one, observed in the simulations, but also the agreement is quantitative in most cases, and this clearly highlights the potential of the CPSE in providing accurate estimations for the higher-order coefficients, thus giving rise to an accurate higher-order HTSE. one of us3 has derived a new perturbation scheme based on a coupling parameter series expansion (thereafter abbreviated as CPSE) which, although different from the above-mentioned HTSE, in certain particular situations, can easily calculate the higher-order coefficients of the HTSE.4 The first few coefficients in the HTSE framework had been obtained from MC simulation for some model potentials. Apart from some pioneering works,5 simulation data for the first- and second-order coefficients in the HTSE have been reported for the square-well6 and Sutherland7 potentials with variable range. In a recent work,8 there have been obtained from the MC simulations the third- and fourth-order coefficients for the square-well fluids with variable range. For potentials with continuous tails, calculations of the third-order coefficient from the simulations are much more cumbersome and that of the fourth-order term seems to be at present unfeasible in general. In this work, we determine from the MC simulations the firstto third-order coefficients in the HTSE framework for fluids with different potential models with hard-sphere cores and continuous tails, and the results are used to the analyze the performance of the CPSE in calculating for these coefficients and also that of the MCA in case of the second-order coefficient only. The structure

1. INTRODUCTION Perturbation theories based on the inverse temperature expansion of the free energy constitute one of the most frequently used tools to obtain the equilibrium thermodynamic properties of fluids,1 and are widely used in phase equilibrium calculations over a long period. However, they present some unsolved problems. However, at low enough temperatures, the higher-order perturbation terms of the high temperature series expansion (HTSE) will become significant, and this is mainly because a factor of the nth power of an inverse temperature is associated the nth perturbation term. Consequently, probably these expansions will diverge at low enough temperatures (a reason why they are known as HTSE), which prevents their applications to systems with very short-ranged potentials, that may exist in the fluid phase at very low temperatures. However, in general, these theories provide few satisfactory results for the perturbation coefficients of the HTSE beyond the first-order, as is the case of the BarkerHenderson macroscopic compressibility approximation (MCA) for the second-order coefficient.2 The first few coefficients in the expansion can be obtained from Monte Carlo (MC) simulation. Apart from the usefulness of these calculations to test the accuracy of perturbation theories, this will provide some insight about the rate of convergence of the series and will help to decide whether efforts to theoretically derive accurate expressions for the higher-order coefficients in this expansion are worthy of consideration. However, quite recently © 2013 American Chemical Society

Received: June 1, 2013 Revised: July 10, 2013 Published: July 11, 2013 9305

dx.doi.org/10.1021/jp405406f | J. Phys. Chem. B 2013, 117, 9305−9313

The Journal of Physical Chemistry B

Article

Table 1. Simulation Results of a1 to a3 for the HSSW Fluidsa ρ* λ = 1.2 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 λ = 1.5 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 λ = 2.0 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 λ = 3.0 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90

a1

a2

Table 2. Simulation Results of a1 to a3 for the HSTW Fluidsa ρ*

a3

−0.1700 −0.3799 −0.6378(1) −0.9526(1) −1.3334(1) −1.7888(1) −2.3250(1) −2.9431(1) −3.6350(1)

−0.0742(1) −0.1416(2) −0.1972(5) −0.2374(6) −0.2579(6) −0.2600(8) −0.2482(7) −0.2294(5) −0.2134(6)

−0.0204(2) −0.031(2) −0.032(2) −0.028(4) −0.028(4) −0.010(5) −0.014(5) −0.006(5) 0.000(3)

−0.5343 −1.1422 −1.8192(1) −2.5545(1) −3.3293(2) −4.1168(1) −4.8823(2) −5.5856(3) −6.1847(4)

−0.1944(4) −0.2922(4) −0.3215(7) −0.3173(5) −0.3093(9) −0.308(1) −0.307(9) −0.286(1) −0.240(2)

−0.056(3) −0.076(3) −0.074(5) −0.054(6) −0.029(4) −0.010(8) −0.005(8) 0.00(1) 0.02(1)

−1.5083(1) −3.0777(1) −4.6716(1) −6.2574(1) −7.8158(1) −9.3509(1) −10.9007(3) −12.5420(5) −14.377(1)

−0.447(1) −0.551(2) −0.547(1) −0.537(1) −0.552(2) −0.597(1) −0.674(4) −0.758(4) −0.804(5)

−0.237(8) −0.285(6) −0.18(1) −0.08(1) −0.05(1) −0.06(2) −0.05(2) −0.08(5) −0.01(8)

−5.4890(3) −11.0427(3) −16.6374(3) −22.2519(2) −27.8596(3) −33.4412(3) −39.0238(5) −44.7253(8) −50.732(1)

−1.425(3) −1.596(3) −1.492(5) −1.404(4) −1.398(5) −1.484(6) −1.654(6) −1.86(1) −2.06(3)

−2.70(7) −3.18(6) −2.1(1) −1.09(6) −0.43(7) −0.26(5) −0.2(1) −0.3(2) 0.1(2)

δ = 1.2 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 δ = 1.4 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 δ = 1.6 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 δ = 1.8 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 δ = 2.0 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90

a

The numbers between parentheses are the statistical uncertainties in the last decimal place.

of the present work is organized as follows: in Section 2, the potential models considered and simulation details are described, and the procedure used to calculate the coefficients of the HTSE in the framework of the CPSE is detailed; the results are collected and analyzed in Section 3; our conclusions are summarized in Section 4.

2. MODELS AND METHODS The model potentials considered include a hard sphere (HS) + square well (HSSW) potential, a HS + attractive Yukawa (HSAY) potential, a HS + Sutherland (HSS) potential, a HS + Triangle-well (HSTW) potential, and a HS + oscillatory tail (HSOT) potential, and their mathematical expressions are described as follows: ⎧∞ , r ≤ σ ⎪ u(r ) = ⎨− ε , σ < r ≤ λσ ⎪ ⎩ 0, r > λσ

a1

a2

a3

−0.0806 −0.1820 −0.3096 −0.4699(1) −0.6710 −0.9226(1) −1.2363(1) −1.6252(2) −2.1032(3)

−0.0234 −0.0470 −0.0693(2) −0.0887(2) −0.1029(3) −0.1103(3) −0.1086(5) −0.0973(6) −0.0774(6)

−0.0050(1) −0.0086(4) −0.0103(3) −0.0111(6) −0.010(1) −0.009(2) −0.007(2) −0.006(1) −0.004(4)

−0.1795 −0.3968 −0.6573 −0.9669(1) −1.3302(1) −1.7507(1) −2.2294(1) −2.7633(1) −3.3451(1)

−0.0458(1) −0.0794(1) −0.0985(3) −0.1027(3) −0.0933(3) −0.0748(2) −0.0529(2) −0.0338(2) −0.0209(2)

−0.0092(2) −0.0134(7) −0.0144(6) −0.0147(8) −0.014(1) −0.011(1) −0.0072(8) −0.0041(6) −0.0013(4)

−0.2988 −0.6470 −1.0468(1) −1.4983(1) −1.9987(1) −2.5412(1) −3.1156(1) −3.7067(1) −4.2963(1)

−0.0670(1) −0.1010(1) −0.1065(3) −0.0922(4) −0.0684(3) −0.0451(2) −0.0278(1) −0.0177(1) −0.0120(1)

−0.0141(6) −0.0204(8) −0.022(1) −0.0203(9) −0.0153(8) −0.0085(5) −0.0033(3) −0.0006(3) 0.0003(2)

−0.4405 −0.9368 −1.4854(1) −2.0793(1) −2.7081(1) −3.3577(1) −4.0121 −4.6546(1) −5.2704(1)

−0.0877(2) −0.1166(2) −0.1067(3) −0.0797(3) −0.0515(2) −0.0307(1) −0.0184(1) −0.0122(1) −0.0095(1)

−0.0212(9) −0.032(1) −0.031(1) −0.0239(4) −0.0135(7) −0.0053(4) −0.0012(2) 0.0000(1) 0.0001(1)

−0.6072 −1.2723(1) −1.9858(1) −2.7351(1) −3.5057(1) −4.2831(1) −5.0563 −5.8206(1) −6.5825(1)

−0.1088(3) −0.1302(3) −0.1072(3) −0.0727(3) −0.0438(3) −0.0257(1) −0.0160(1) −0.0116(1) −0.0095(1)

−0.033(1) −0.048(1) −0.0431(7) −0.0276(6) −0.0133(7) −0.0049(2) −0.0016(2) −0.0006(2) −0.0003(1)

a

The numbers between parentheses are the statistical uncertainties in the last decimal place.

for the HSSW potential, where λ is the potential width in units of the diameter σ of the HS; ⎧∞ r ≤ σ ⎪ u(r ) = ⎨ ε (r /σ − δ)/(δ − 1), σ < r ≤ δσ ⎪ ⎩ 0 r > δσ

(1) 9306

(2)

dx.doi.org/10.1021/jp405406f | J. Phys. Chem. B 2013, 117, 9305−9313

The Journal of Physical Chemistry B

Article

Table 3. Simulation Results of a1 to a3 for the HSAY Fluidsa ρ* κ=3 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 κ=5 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 κ = 10 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90

a1

a2

Table 4. Simulation Results of a1 to a3 for the HSS Fluidsa ρ*

a3

−0.2975 −0.6329 −1.0071 −1.4209(1) −1.8745(1) −2.3678(1) −2.9002(1) −3.4709(1) −4.0787(1)

−0.0401(1) −0.0596(1) −0.0640(2) −0.0591(2) −0.0490(1) −0.0374(2) −0.0265(2) −0.0177(1) −0.0111(1)

−0.0085(2) −0.0129(5) −0.0135(4) −0.0125(3) −0.0096(6) −0.0063(4) −0.0038(4) −0.0021(5) −0.0010(5)

−0.1643 −0.3581 −0.5849 −0.8484(1) −1.1521(1) −1.4996(1) −1.8944(1) −2.3398(1) −2.8387(2)

−0.0273(1) −0.0464 −0.0573(2) −0.0610(2) −0.0585(1) −0.0518(2) −0.0426(2) −0.0327(2) −0.0233(2)

−0.0052(1) −0.0080(3) −0.0090(2) −0.0093(4) −0.0081(7) −0.0064(7) −0.0048(5) −0.0033(4) −0.0019(4)

−0.0768 −0.1719 −0.2892 −0.4335(1) −0.6108 −0.8278(1) −1.0927(1) −1.4153(2) −1.8069(3)

−0.0153 −0.0299 −0.0429(1) −0.0536(1) −0.0613(1) −0.0655(2) −0.0656(3) −0.0620(3) −0.0547(4)

−0.0029(1) −0.0050(2) −0.0063(2) −0.0071(4) −0.0068(6) −0.0063(8) −0.0057(8) −0.005(1) −0.003(1)

n=6 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 n=9 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 n = 18 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90

a

The numbers between parentheses are the statistical uncertainties in the last decimal place.

for the HSTW potential, where the parameter δ is used to adjust the potential range; ⎧∞ , r ≤ σ ⎪ ⎪ u(r ) = ⎨ −κ * (r − σ )/ σ ⎪− εσ e ,r>σ ⎪ ⎩ r

a3

−0.2207 −0.4724 −0.7576 −1.0788(1) −1.4386(1) −1.8399(1) −2.2853(1) −2.7778(1) −3.3202(2)

−0.0287(1) −0.0465 −0.0552(2) −0.0568(2) −0.0532(1) −0.0464(2) −0.0379(2) −0.0293(2) −0.0212(2)

−0.0056(1) −0.0086(3) −0.0094(2) −0.0094(4) −0.0079(6) −0.0059(6) −0.0043(5) −0.0028(3) −0.0016(6)

−0.1148 −0.2525 −0.4169 −0.6123(1) −0.8436 −1.1162(1) −1.4363(1) −1.8111(2) −2.2483(2) −0.0774 −0.1726 −0.2894 −0.4325(1) −0.6072 −0.8203(1) −1.0794(1) −1.3942(2) −1.7760(3)

−0.0193 −0.0352 −0.0472(1) −0.0551(1) −0.0585(1) −0.0581(2) −0.0541(2) −0.0476(3) −0.0392(2) −0.0145 −0.0281 −0.0403(1) −0.0504(1) −0.0578(1) −0.0621(1) −0.0629(3) −0.0603(3) −0.0543(4)

−0.0036(1) −0.0059(2) −0.0069(2) −0.0075(4) −0.0069(6) −0.0060(8) −0.0051(7) −0.0040(7) −0.0026(8) −0.0027(1) −0.0047(2) −0.0059(2) −0.0067(3) −0.0065(5) −0.0060(7) −0.0054(7) −0.001(7) −0.001(2)

−0.0467 −0.1055 −0.1796 −0.2730 −0.3910 −0.5402(1) −0.7290(1) −0.9685(2) −1.2729(3)

−0.0095 −0.0196 −0.0300(1) −0.0404(1) −0.0503(1) −0.0590(1) −0.0658(2) −0.0701(3) −0.0706(5)

−0.0018 −0.0034(1) −0.0045(1) −0.0055(2) −0.0058(4) −0.0059(5) −0.0055(6) −0.005(1) −0.004(2)

The numbers between parentheses are the statistical uncertainties in the last decimal place.

(3)

⎧∞ , r / σ ≤ 1 ⎪ ⎪ ⎪(εb + 1)(r /σ − 1)/(b − 1) − 1, b ≥ r /σ ≥ 1 u(r )/ε = ⎨ ⎪ ⎡ ⎛r ⎞⎤ ⎛ r ⎞n ⎞⎤ ⎡ ⎛ r ⎪ ⎪ εbexp⎢⎣ − α⎜⎝ σ − b⎟⎠⎥⎦cos⎢⎣ks⎜⎝ σ − b⎟⎠⎥⎦ /⎜⎝ bσ ⎟⎠ , b < r /σ ⎩ (6)

(4)

for the HSOT potential, which is introduced9 as a simple model capable of reproducing the main features of the depletion potential in colloidal suspensions as well as the effective interionic pair potential in some alkali metals and alloys, and combines the potential parameters εb, α, b, ks, and n to adjust the shape and range of the model potential. One points out that throughout eqs 1−6, ε is an energy parameter of the potential models.

for HSS potential, where n is a positive integral number that determines the effective range of the potential; ⎧∞ , r / σ ≤ 1 ⎪ ⎪ ⎪ uts(r ) = ⎨ u(r ) − u(rcσ ), 1 < r /σ ≤ rc ⎪ ⎪ ⎪ 0, r < r /σ ⎩ c

a2

a

for the HSAY potential, where the parameter κ* determines the effective range of the potential; ⎧∞ , r ≤ σ ⎪ ⎪ u(r ) = ⎨ n ⎪− ε⎛⎜ σ ⎞⎟ , r > σ ⎪ ⎝ ⎠ ⎩ r

a1

In the HTSE, the excess Helmholtz free energy per particle fex of a fluid with an intermolecular potential u(r) can be expressed

(5) 9307

dx.doi.org/10.1021/jp405406f | J. Phys. Chem. B 2013, 117, 9305−9313

The Journal of Physical Chemistry B

Article

Table 5. Simulation Results of a1 to a3 for the HSOT Fluidsa ρ* set 1 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 set 2 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 set 3 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90

a1

a2

a3

0.0032 0.0014 −0.0073 −0.0259 −0.0592 −0.1157(1) −0.2089(1) −0.3604(2) −0.6003(5)

−0.0259(1) −0.0570(1) −0.0938(2) −0.1373(2) −0.1877(2) −0.2447(5) −0.3035(8) −0.357(1) −0.394(2)

−0.0017(1) −0.0034(4) −0.0044(9) −0.005(1) −0.003(2) −0.002(5) 0.000(4) 0.012(9) 0.01(3)

0.0689 0.1442 0.2280(1) 0.3230 0.4313(1) 0.5532(1) 0.6834(1) 0.8062(2) 0.8893(5)

−0.0636(1) −0.1370(2) −0.2207(4) −0.3153(4) −0.4205(6) −0.5378(7) −0.6683(7) −0.806(2) −0.940(4)

0.0045(5) 0.008(1) 0.012(2) 0.014(4) 0.014(6) 0.008(8) 0.01(2) 0.01(5) −0.03(6)

−0.0239 −0.0530 −0.0886 −0.1329 −0.1892 −0.2624 −0.3605(1) −0.4967(1) −0.6932(3)

−0.0144 −0.0320(1) −0.0533(1) −0.0791(1) −0.1101(1) −0.1471(2) −0.1901(3) −0.2375(6) −0.2870(9)

−0.0017 −0.0034(1) −0.0049(3) −0.0062(5) −0.0078(9) −0.007(1) −0.007(1) −0.004(5) 0.001(7)

a3 = −

n=1

an(ρ) T *n

1 N

a2 = −

∑ ⟨Ni⟩0 u1*(ri) i

1 2N

* * ∑ {⟨NN i j⟩0 − ⟨Ni⟩0 ⟨Nj⟩0 }u1 (ri)u1 (rj) ij

a1 = −

⟨M ⟩0 N

a2 = −

1 1 (M − ⟨M ⟩0 )2 2N

a3 = −

1 1 (M − ⟨M ⟩0 )3 6N

(11)

0

(12)

0

(13)

where M is the number of pairs separated a distance x ≤ λ. The expressions 8−13 can be evaluated from simulation in the HS reference system, as stated before. Such calculations have been carried out by the above-mentioned authors, among others, for several potential models, although mainly for the HSSW potential. In the case of potentials with continuous tails, the simulations have limited up to now to the first- and secondorder coefficients, and for a reduced number of potential models, because the calculation of the third order term is computationally very demanding and that of the fourth-order term seems to be unfeasible at present. The calculation of the latter from simulation is more affordable for the HSSW potential, and so Alder et al.5 determined in this way the thirdand fourth-order coefficients for λ = 1.5 and more recently ́ Espindola-Heredia et al.8 carried out this kind of simulations for 1.1 ≤ λ ≤ 3.0. According to these results, the second-order coefficient is typically 1 order of magnitude smaller than that of the first order. In contrast, the behavior of the third- and fourth-order coefficients is more complex as revealed by the simulation data ́ reported by Espindola-Heredia et al.:8 at low densities their contributions are more important than those of the second-order coefficient for high values of λ, but become negligible at high densities. In this work, we have carried out simulations in the NVT ensemble of the reference HS fluid to obtain the values of a1 to a3 from the expressions 8−13 for the potential models described in Section 1 with different potential parameters. We considered systems with N = 500 particles, initially placed in an FCC configuration, 5 × 104 cycles for equilibration and 106 cycles for production, each cycle consisting in N attempted particle moves. Averages over ten independent runs for each density were performed and the statistical error was estimated as the standard deviation. We took Δr = 0.02σ for the HSSW, HSTW, and HSOT potentials and Δr = 0.04σ for the HSAY and HSS potentials. A cutoff distance rc = 4.0 was settled for the HSAY and HSS potentials without shifting, and a cutoff distance rc = 3.0 was settled for the HSOT potentials with shifting. We have found that the results are only slightly

(7)

In this expansion, fex−ref is the excess free energy per particle of a reference fluid, the HS fluid for the potential models considered here, β = 1/kBT, and an is the coefficient of the nth perturbation term, and is related to fluctuations of the energy due to the perturbation calculated in the canonical ensemble of the reference system. The first few of these coefficients can be obtained by means of simulations in the reference system from the fluctuations in the perturbation energy. Thus, the first three terms, in which we are interested here, are given by the following:10,5 a1 =

(10)

where Ni is the number of molecular distances in the interval (ri,ri+1), with Δr = ri+1−ri ≪ σ, i = 0,1,..., angular brackets indicate averages, subscript 0 means that the averages are performed in the reference system, and u1* = u1/ε, where u1(r) is the perturbation part of the potential and equal to u(r) subtracted by the HS potential. For the particular case of a system with the HSSW potential, u1* (r) = −1 within the potential well 1 ≤ x ≤ λ, where x = r/σ, and the expressions of these perturbation coefficients take the simple form5

as an expansion in power series of the inverse of the reduced temperature T* = kBT/ε, where T is the absolute temperature, and kB is the Boltzmann constant, in the form4,10 ∝

ijk

}

The numbers between parentheses are the statistical uncertainties in the last decimal place.



∑ {⟨NNN i j k⟩0 − 3⟨NN i j⟩0 ⟨Nk⟩0

+ 2⟨Ni⟩0 ⟨Nj⟩0 ⟨Nk⟩0 u1*(ri)u1*(rj)u1*(rk)

a

βfex = βfex − ref +

1 6N

(8)

(9) 9308

dx.doi.org/10.1021/jp405406f | J. Phys. Chem. B 2013, 117, 9305−9313

The Journal of Physical Chemistry B

Article

Figure 1. HTSE coefficients a1 to a3 for fluids with the HSSW potential (1) of different potential ranges. Points are the simulation data of this work and the dashed curves are the theoretical results from the CPSE. The continuous curves for a2 are the results of the BH MCA.

fluid, the two quantities are obtained from the Verlet-Weis correction12 of the Percus−Yevick solution,13 and a Carnahan− Starling (CS) equation of state.14 The CPSE can be presented as follows:1,3,4,9,11,15,16

sensitive to the size of Δr within the above-mentioned ranges and that with the truncation at rc = 4.0, there are not significant differences between the shifted and unshifted potentials. Besides the MC simulations, there are theoretical approaches allowing to calculate for these coefficients. An extensively used scheme for estimation of a2 is the MCA proposed by Barker and Henderson;2 although it is known for a long time that the MCA significantly underestimates the magnitude of a2, a recent study indicates11 that it also overestimates the magnitude of a2 for some cases, and a general conclusion is that the MCA unfortunately fails to give even a qualitative description of a2 for a lot of moderate circumstances. However, in this work, the MCA is still used for the aim of comparison with the CPSE. The MCA reads as follows: a 2 = πρ



(15)

where, βfper − n = (1/n!)2πρ



drr 2βu1(r )

∂(n − 1)g imag (r , ζ , ρ , T ) ∂ζ (n − 1) ζ=0

(16)

As in the HTSE, the HS fluid is chosen as the reference fluid, then, fex−ref is the excess Helmholtz free energy per particle of the HS fluid, and is acquired by integrating the HS CS equation of state as done in the practices of the CPSE.1,3,4,9,11,15,16 gimag(r,ζ, ρ,T) is the rdf of an imaginary bulk fluid with a pair potential u(r;ζ) = uHS(r)+ζu1(r), ∂(n−1)gimag/∂ζ(n−1)|ξ=0 is (n−1)−order derivative of gimag w.r.t. the coupling parameter ζ



ref

∑ βfper− n n=1



∫ drr 2(βu1(r))2 gref (r , ρ) 1β ⎝ ∂∂Pρ ⎠ ⎜



βfex − Zhou = βfex − ref +

(14)

where gref(r,ρ) is radial distribution function (rdf) of the reference fluid, and (1/β)((∂ρ)/(∂P))ref is inverse of a reduced compressibility of the reference fluid (P denoting pressure, and ρ number density); given that the HS fluid is used as reference 9309

dx.doi.org/10.1021/jp405406f | J. Phys. Chem. B 2013, 117, 9305−9313

The Journal of Physical Chemistry B

Article

Figure 2. HTSE coefficients a1 to a3 for fluids with the HSTW potential (2) of different potential ranges. Points are the simulation data of this work and the dashed curves are the theoretical results from the CPSE. The continuous curves for a2 are the results of the BH MCA.

evaluated at ζ = 0, ∂0gimag(r,ζ, ρ,T)/∂ζ0|ζ=0=gimag(r,0, ρ,T) is actually the rdf of the HS fluid. For the numerical details of calculating ∂(n−1)gimag∂ζ(n−1)|ξ=0, the interested readers can have a reference to the original papers1,3,4,9,11,15,16 relevant to the CPSE, and not repeated here. It’s important to point out that Zhou4 has proven that if T*=1, one has the following: an = βfper − n |T *= 1 , n ≥ 1

The results of a1 and a2 for the HSSW fluids have very small statistical uncertainty and are in generally excellent agreement with the data reported in ref 8, with the exception of a2 for λ = 3.0, for which the absolute values of the present data are greater for ρ* ≤ 0.5. These differences are attributed to the fact that we have used a greater number of particles in the simulations as well as a considerable number of independent runs. We first analyze the behavior pertaining to the HSSW potential. The situation for a3 is more complex and will be discussed with details. The statistical error increases considerably, especially for values of a3 close to zero, and the error rises with the number of particles used in the simulation, a fact that was also noted in ref 8. The relative error is more notable for λ ≤ 1.5, because the values of a3 are very small. This may have a considerable importance for λ = 1.2, as in this case the fluid phase extends to reduced temperatures well below T* = 1 and, from eq 7, the term corresponding to a3 contributes to fex with a factor 1/T*3. For λ = 2.0 and especially for λ = 3.0, the absolute values of a3 are considerably high for reduced densities ρ* ≤ 0.4 and the relative errors correspondingly small.

(17)

Equation 17 supplies a mathematical relation between the HTSE and CPSE perturbation frameworks, and is the basis of using the CPSE to calculate an, the most essential input information in implementing the HTSE.

3. RESULTS AND DISCUSSIONS The MC results are reported in Tables 1−5 for the potential models cited above with different sets of parameters and reduced densities ρ* = ρσ3 ranging from 0.1 to 0.9. Some comments are relevant. 9310

dx.doi.org/10.1021/jp405406f | J. Phys. Chem. B 2013, 117, 9305−9313

The Journal of Physical Chemistry B

Article

Figure 3. HTSE coefficients a1 to a3 for fluids with the HSAY potential (3) of different potential ranges. Points are the simulation data of this work and the dashed curves are the theoretical results from the CPSE. The continuous curves for a2 are the results of the BH MCA.

Figure 4. HTSE coefficients a1 to a3 for fluids with the HSS potential (4) of different potential ranges. Points are the simulation data of this work and the dashed curves are the theoretical results from the CPSE. The continuous curves for a2 are the results of the BH MCA.

However, the importance of a3 will be nearly negligible, as the reduced temperatures for the fluid phase are usually much greater than 1. As an example, the critical temperature for λ = 3.0 is Tc* = 9.87.17 In any case, our results for a3 are in quite reasonable agreement with those reported in ref 8. Concerning the remaining potential models considered, one can see in the tables that the statistical errors in a1 and a2 are negligible. The absolute values of a3 are generally small and decrease, while the relative errors increase with decreasing the effective potential range, with the highest statistical errors occurring at high densities where, for the potential with very small effective potential ranges, it becomes enormous as compared with the values of a3. In spite of this, the values of a3 listed in the tables are considered reasonably accurate as they are averaged values over ten independent runs, as said before. The results for a1 to a3 obtained from the CPSE are compared in Figures 1−5 with the simulation data reported in Tables 1−5. The results provided for a2 by the BH MCA2 are also included for comparison. For a1, there is little to say, because the CPSE provides for this coefficient the same expression as the HTSE,

which is known to be accurate for this purpose. Concerning a2, it is clearly seen in these figures that the CPSE provides much better predictions for a2 than the BH MCA. The BarkerHenderson local compressibility approximation (LCA)2 is slightly more accurate than the MCA, but the difference is small6,18 and so it has not been included in the figure. It is particularly noteworthy the excellent performance of the CPSE for a2 in the HSSW potential, for which the theory reproduces the simulation data nearly exactly for large potential widths and continues to be very accurate for relatively small width λ = 1.2. For the HSTW, HSAY, and HSS potentials, the theory is also remarkably accurate for a2, although some deviations with respect to the simulation data appear at high densities. For the HSOT potential, again the CPSE is in perfect agreement with simulations for a2. As for a3, Figure 1 shows that the CPSE is again impressively accurate for the HSSW potential, even for the shortest potential width considered. It is also very accurate for the HSTW potential with range δ ≥ 1.6, and somewhat less accurate for shorter potential ranges at moderate to high densities, as well as for the HSAY and HSS potentials considered in this work. 9311

dx.doi.org/10.1021/jp405406f | J. Phys. Chem. B 2013, 117, 9305−9313

The Journal of Physical Chemistry B

Article

Figure 5. HTSE coefficients a1 to a3 for fluids with the HSOT potential (5) of different potential ranges. Points are the simulation data of this work and the dashed curves are the theoretical results from the CPSE. The continuous curves for a2 are the results of the BH MCA.

approaches: one the well-known BH MCA and the other a recently proposed CPSE theorized by one of the present authors, are used to calculate these coefficients, and the relevant results are compared with those from the MC simulations. A comprehensive comparison based on these coefficients from the two routes leads to several natural conclusions: (1) First of all, some comments about the simulations are timely. Whereas simulations for the first- and secondorder coefficients of the HTSE have been reported in the literature for several potential models, and even up to the fourth-order coefficient for some discrete potential models, to the best of our knowledge, this is the first time that are reported simulations for the third-order coefficient for potential models with continuous tails. This is due to the fact that the simulation of the third-order coefficient for continuous potentials is computationally very demanding, it is generally very small in magnitude, and the relative statistical errors involved are considerable at high densities. In spite of this, the third- and higher-order terms may yield a non-negligible contribution to the thermodynamic properties at very low temperatures, at which systems with very short-ranged attractive interactions may remain in the fluid state. As there are some perturbation theories, noteworthy the CPSE, capable of providing a number of perturbation coefficients beyond second-order, the availability of simulation data for them to test the performance of theories would be desirable. Unfortunately, going beyond the third-order term by computer simulation seems to be unfeasible at present, because the computational costs would become prohibitive and, moreover, the statistical errors, which quickly increase with the order of the perturbation term beyond the second one, would become enormous.

For the HSOT potential, the CPSE results for a3, not only are in general excellent agreement with simulations, but also seems to be even more reliable than the latter, on accounting for the considerable uncertainty in the simulation data at high densities. In contrast, as the MCA provides poor results for the second-order term, it is not expected that it will work better for the third- and higher-order terms. As a matter of fact, Praestgaard and Toxvaerd19 resummed the perturbation expansion in the MCA approximation and reported expressions for the resulting free energy and pressure. Their results reveal that the resummed expansion considerably underestimates the contribution of the terms beyond the first-order one. Also Barker and Henderson10 discussed the performance of the MCA for the second- and third-order terms and, in particular, they concluded that the theory underestimates the magnitude the third-order term. Therefore, again it is not worth going further with this approach. From the preceding analysis, it seems clear that the CPSE provides very accurate predictions of the first three terms of the HTSE for a wide variety of hard-core potential models either with discrete or continuous potential tails. It is worth recalling here that the CPSE is not an inverse temperature expansion, although on a special occasion can provide the perturbation coefficients of the latter, and that it can be carried out beyond the third-order term. In contrast, performing computer simulations for perturbation coefficients higher than the third one for potentials with continuous tails seems to be unfeasible at present.

4. CONCLUSIONS To summarize, in the present study, systematic MC simulation data are provided for the first three coefficients of the HTSE for fluids with several potential models consisting of hard-spherical cores and continuous tails; in addition to this, two theoretical 9312

dx.doi.org/10.1021/jp405406f | J. Phys. Chem. B 2013, 117, 9305−9313

The Journal of Physical Chemistry B

Article

R. T.; Ullmann, G. M. J. Phys. Chem. B 2011, 115, 507. Ferrando, N.; de Hemptinne, J. C.; Mougin, P.; Passarello, J. P. J. Phys. Chem. B 2012, 116, 367. (2) Barker, J. A.; Henderson, D. J. Chem. Phys. 1967, 47, 2856−2861. (3) Zhou, S. Phys. Rev. E 2006, 74, 031119. (4) Zhou, S. AIP Adv. 2011, 1, 040703. Zwanzig, R. W. J. Chem. Phys. 1954, 22, 1420. (5) Barker, J. A.; Henderson, D. Proceedings of the Fourth Symposium on Thermophysical Properties; Moszynsky, J. R., Ed.; American Society of Mechanical Engineers: New York, 1968; pp 30−36. Alder, B. J.; Young, D. A.; Mark, M. A. J. Chem. Phys. 1972, 56, 3013. (6) Largo, J.; Solana, J. R. Mol. Simul. 2003, 29, 363. Largo, J.; Solana, J. R. J. Phys. Chem. B 2004, 108, 10062. (7) Díez, A.; Largo, J.; Solana, J. R. J. Chem. Phys. 2006, 125, 074509. (8) Espíndola-Heredia, R.; del Río, F.; Malijevsky, A. J. Chem. Phys. 2009, 130, 024509. (9) Zhou, S.; Solana, J. R. Phys. Rev. E 2008, 78, 021503. (10) Barker, J. A.; Henderson, D. Rev. Mod. Phys. 1976, 48, 587. (11) Zhou, S.; Solana, J. R. J. Chem. Phys. 2009, 131, 134106. (12) Verlet, L.; Weis, J.-J. Phys. Rev. A 1972, 5, 939. (13) Thiele, E. J. Chem. Phys. 1963, 39, 474. Wertheim, M. S. Phys. Rev. Lett. 1963, 19, 321. (14) Carnhan, N. F.; Starling, K. E J. Chem. Phys. 1969, 51, 635. (15) Zhou, S. J. Phys. Chem. B 2007, 111, 10736. Zhou, S.; Solana, J. R. Phys. Chem. Chem. Phys. 2009, 11, 11528. (16) Zhou, S. Phys. Rev. E 2008, 77, 041110. (17) Orkoulas, G.; Panagiotopoulos, A. Z. J. Chem. Phys. 1999, 110, 1581. (18) Largo, J.; Solana, J. R. Mol. Phys. 2003, 101, 2981. (19) Praestgaard, E.; Toxvaerd, S. J. Chem. Phys. 1969, 51, 1895. Toxvaerd, S.; Praestgaard, E. J. Chem. Phys. 1970, 53, 2389.

(2) The widely used BH MCA only provides very rough estimations of the second-order coefficient for all models considered. To be specific, we have found both overestimation and underestimation of the second-order coefficient by the MCA, and this finding differs from one common misperception that the MCA usually underestimates the second-order coefficient. Another fact that emerges from the comparison performed in this work is that the performance of the MCA worsens as the range of the potential decreases, and this may contribute to explain why the second-order HTSE based on the MCA worsens as the potential range is shortened. (3) In contrast, the CPSE shows excellent performance in calculating higher-order coefficients of the HTSE. For the cases in which the statistical errors of the simulation data are small enough, and thus clear comparisons are possible, the CPSE results are in excellent agreement with the simulation ones for all of the models considered; this confirms the validity of the ML HS bridge function approximation used in calculating the derivatives of the imaginary fluid rdf w.r.t., the coupling parameter, ζ, which contributes to the excellent performance of the CPSE especially evident in the case of the HSOT potential. The only weak point here is that the CPSE results deviate appreciably in some cases from the corresponding simulation results at high densities and this might point out the need of further improvement on the HS bridge function approximation; another possible explanation is that characteristics of solid states more or less occur in the liquid states at high densities, and they can be accounted for by simulations, but goes beyond the role of the liquid OZ theory; the observed deviations at high densities for the third-order coefficient also may come from the fact that the simulation errors are very large, and so a reliable comparisons are not possible. In any case, the CPSE results not only reproduce the trends in the density dependence of the perturbation coefficients, even the third one, observed in the simulations, but also the agreement is quantitative in most cases. This clearly highlights the potential of the CPSE in providing accurate estimations for the higherorder coefficients of the HTSE beyond the third-order one, although the reliability of these theoretical estimates cannot be tested at present because of the above-mentioned limitations of the simulations, thus giving rise to an accurate higher-order HTSE.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS J.R.S. acknowledges financial support from the Spanish Ministerio de Ciencia e Innovación (MICINN) under Grant No. FIS2009-09616. This project is supported by the National Natural Science Foundation of China (No. 21173271).



REFERENCES

(1) Zhou, S.; Solana, J. R. Chem. Rev. 2009, 109, 2829. Beierlein, F. R.; Michel, J.; Essex, J. W. J. Phys. Chem. B 2011, 115, 4911. Ullmann, 9313

dx.doi.org/10.1021/jp405406f | J. Phys. Chem. B 2013, 117, 9305−9313