4656
J. Phys. Chem. B 2008, 112, 4656-4661
Simulation of Chemical Potentials and Phase Equilibria in Two- and Three-Dimensional Square-Well Fluids: Finite Size Effects† Horst L. Vo1 rtler,*,‡ Katja Scha1 fer,‡ and William R. Smith§ Molecular Dynamics and Computer Simulation Research Group, Institute of Theoretical Physics, UniVersity of Leipzig, Postfach 100920, 04009 Leipzig, Germany, and Faculty of Science, UniVersity of Ontario Institute of Technology, 2000 Simcoe Street N, Oshawa ON L1H 7K4, Canada ReceiVed: May 15, 2007; In Final Form: January 28, 2008
We study the simulation cell size dependence of chemical potential isotherms in subcritical square-well fluids by means of series of canonical ensemble Monte Carlo simulations with increasing numbers of particles, for both three-dimensional bulk systems and two-dimensional planar layers, using Widom-like particle insertion methods. By estimating the corresponding vapor/liquid coexistence densities using a Maxwell-like equal area rule for the subcritical chemical potential isotherms, we are able to study the influence of system size not only on chemical potentials but also on the coexistence properties. The chemical potential versus density isotherms show van der Waals-like loops in the subcritical vapor/liquid coexistence range that exhibit distinct finite size effects for both two- and three-dimensional fluids. Generally, in agreement with recent findings for related studies of Lennard-Jones fluids, the loops shrink with increasing number of particles. In contrast to the subcritical isotherms themselves, the equilibrium vapor/liquid densities show only a weak system size dependence and agree quantitatively with the best-known literature values for three-dimensional fluids. This allows our approach to be used to accurately predict the phase coexistence properties. Our resulting phase equilibrium results for two-dimensional square-well fluids are new. Knowledge concerning finite size effects of square-well systems is important not only for the simulation of thermodynamic properties of simple fluids, but also for the simulation of models of more complex fluids (such as aqueous or polymer fluids) involving square-well interactions.
1. Introduction The study of the influence of the system size in computer simulations of fluid many particle systems, which originates from the approximation of a macroscopic system by a small periodically replicated simulation sample cell consisting typically only of several hundreds or thousands of particles, is one of the most important problems in molecular simulation methodology. The periodic boundary conditions result in the generation of an infinitely extended system with an imprinted periodicity on the order of the length of the simulation cell, where the differences between the properties of the macroscopic original system and the properties of the simulated periodic system vanish with increasing size of the simulated sample cell. Such simulated infinite systems with periodic continuation are well defined in the thermodynamic limit (V f ∞ and N f ∞) and therefore allow a direct comparison with statisticalmechanical theories. Whereas for the simulation of thermodynamic properties of single fluid phases with short-ranged forces it is well-known that systems with particles numbers on the order of some hundred to thousand match the macroscopic properties,1-3 the knowledge for phase equilibria properties is much less satisfactory4 and requires further studies. Recently we studied chemical potential versus density isotherms for simple square-well fluids in the bulk and under † Originally submitted for the “Keith E. Gubbins Festschrift”, published as the November 1, 2007, issue of J. Phys. Chem. C (Vol. 111, No. 43). * Corresponding author. Phone: +49-341-235 2435. Fax: +49-341235 2307. E-mail:
[email protected]. ‡ University of Leipzig. § University of Ontario Institute of Technology.
inhomogeneous conditions5 as well as primitive water models, where hydrogen bonding is described by square-well interactions.6,7 We estimated the chemical potentials using canonical Monte Carlo (MC) simulations combined with particle insertion methods. From the subcritical chemical potential versus density isotherms we estimated the corresponding vapor/liquid coexistence densities and chemical potentials by means of a Maxwelllike equal area rule.5 For bulk square-well fluids, we found good agreement between the vapor/liquid coexistence data and the literature values, even for moderate numbers of particles (N ≈ 200-500), but a systematic study of finite size effects has not been performed. In general, a study of finite size effects requires a series of simulations for the same system with increasing numbers of particles in the simulation cell, up to particle numbers that allow extrapolation to the infinite system (finite size scaling). For complex fluid models, such simulations are rather time-consuming, and therefore a detailed study of the system size dependence over a wide range of particle numbers is usually restricted to simplified models. In general, a fluid at subcritical temperatures in the density range between the coexistence vapor and liquid density has to be considered as an inhomogeneous system that exhibits an equilibrium between liquid domains and the surrounding gas. In the framework of Ising-like lattice gas models in the grand canonical ensemble, Furukawa and Binder8 showed some time ago that, in inhomogeneous mixed state systems with periodic boundary conditions, van der Waals-like loops in the chemical potential versus density isotherms have to be expected, which will disappear with increasing system size according to a power law in the linear system size.9 In our studies, however, a direct
10.1021/jp073726r CCC: $40.75 © 2008 American Chemical Society Published on Web 03/22/2008
Chemical Potential Isotherms in Square-Well Fluids
Figure 1. Chemical potential versus density isotherms of the bulk square-well fluid with N ) 1000 for the number of subcritical temperatures ranging from T* ) 0.95 to T* ) 1.2. For comparison, some results for N ) 216 are included; T* ) 1.2: × - simulations, 2 - simulations (N ) 216); T* ) 1.15: f- simulations; T* ) 1.10: + - simulations; T* ) 1.05: b - simulations; T* ) 1.0: 9 simulations; T* ) 0.95: ( - simulations, O - simulations (N ) 216). Solid lines are numerical fits to eq 4.
application of these scaling laws is not possible, since we work on canonical conditions (NVT ensemble), where we have to set the number of particles and to adjust the system size in order to fix the density; i.e., in contrast to µVT conditions, we cannot simply scale our system size in terms of the linear dimension. This is quite the same situation as in Gibbs ensemble simulations of phase equilibria, which usually are performed in the NVT ensemble.4 Recently MacDowell et al.10,11 reported distinct finite size effects on chemical potential versus density isotherms of Lennard-Jones fluids in the subcritical vapor/liquid transition range. The present paper investigates such effects for squarewell systems in the bulk as well as in two-dimensional layers. The main goal of the paper is to examine to what extent finite size effects affect phase equilibrium properties. We study the square-well fluid mainly for three reasons: (i) It is one of the simplest possible molecular fluid models that represents the realistic fluid phase equilibria of real simple liquids. (ii) It is a basic intermolecular potential model that is involved in molecular models of several classes of more complex molecular liquids such as polymers12 or aqueous liquids.13 (iii) Two-dimensional square-well monolayers are closely related to the two-dimensional Ising model, one of the most studied many-particle lattice models, where analytical results are available and which plays a crucial role not only for a microscopic understanding of phase transitions of lattice systems but also for fluid (off-lattice) systems, taking into account the representation of a fluid in terms of a lattice-gas model, which is isomorphic to the original Ising model.14 We compare our results for the phase equilibria of squarewell fluids with the most accurate literature data for such fluids15 and discuss the finite size effects in the context of the abovementioned findings for grand canonical MC simulations on bulk Lennard-Jones systems.10,11 The paper is organized as follows. In the next section we describe the fluid models and the simulation methods used, including some technical details. In the third and fourth sections we present our results for finite size effects on chemical potentials and phase equilibria for both bulk systems and planar
J. Phys. Chem. B, Vol. 112, No. 15, 2008 4657
Figure 2. Chemical potential versus density isotherms for different numbers of particles of the bulk square-well fluid at T* ) 1.05; N ) 64: f - simulations; N ) 216: b - simulations; N ) 1000: ( simulations; N ) 3375: × - simulations. Solid lines are numerical fits to eq 4.
Figure 3. Chemical potential versus the number of particles of a bulk square-well fluid at T* ) 1.05, for selected densities in the stable gas and liquid phase and in the liquid-gas density gap; coexistence chemical potential: × - simulations, solid line - del Rio et al.;15 stable liquid (Fσ3 ) 0.65): b - simulations; stable gas (Fσ3 ) 0.03): f simulations; unstable fluid (Fσ3 ) 0.274): 2 - simulations; unstable fluid (Fσ3 ) 0.391): 1 - simulations.
monolayers. This is followed by a discussion and comparison of our bulk and monolayer results with existing literature results. Finally, as a summary of our studies, we draw some general conclusions concerning the system size dependence of chemical potentials and fluid phase equilibria. 2. Systems under Consideration and Simulation Methods We consider a fluid with intermolecular potential model modeled by square-well potential in its usual form:
{
∞ for r e σ u(r) ) - for σ < r e λσ 0 for r > λσ
(1)
where σ is the diameter of the particles, and and λ describe the depth and the width of the potential well, respectively. Throughout the paper we set λ ) 1.5 and use reduced units
r* ) r/σ, u* ) u/, T* ) kT/
(2)
4658 J. Phys. Chem. B, Vol. 112, No. 15, 2008
Vo¨rtler et al.
TABLE 1: Vapor/Liquid Coexistence Densities and Chemical Potentials for Three-Dimensional Square-Well Fluids, Obtained by the Maxwell Equal Area Rule from Chemical Potential versus Density Isotherms, and Comparison with Literature Dataa T*
Fvσ3 N ) 1000
Fvσ3 N ) 216
0.95 1.00 1.05 1.10 1.15 1.20
0.0202 0.0303 0.0445 0.0636 0.0954 0.155
0.0203
a
Fvσ3 a
0.0424
0.0307 0.0445
0.144
0.1545
Flσ3 N ) 1000
F lσ 3 N ) 216
0.669 0.646 0.617 0.583 0.539 0.468
0.668
Flσ3 a
0.618
0.6443 0.6177
0.479
0.4515
µco/kT N ) 1000
µco/kT N ) 216
-4.186 -3.887 -3.628 -3.399 -3.196 -3.017
-4.197
µco/kT a
-3.634
-3.877 -3.619
-3.022
-3.017
See ref 15.
for the distance and for the energy and the temperature, respectively. In this paper we study both normal three-dimensional bulk square-well fluids and planar monolayers (i.e., two-dimensional arrangements) of square-well particles. In order to study finite size effects of the chemical potential in both the bulk system and the monolayer, we perform standard canonical ensemble MC simulations1,4 using periodic boundary conditions in three and two dimensions, respectively. We generate the corresponding Markov chains of configurations and measure the chemical potential by means of (virtual) test-particle insertions. To examine the efficiency of the insertions, we have analyzed several variants of virtual particle insertion methods. Starting with the well-known Widom method16 particularly, we considered the “unbonded” particle insertion of Tripathi and Chapman,17 including some recent extensions,7,18 and the scaledparticle insertion attributed to Labik and Smith.19 Since the densities in the gas/liquid transition range of interest are not extremely high, we found that the standard Widom test particle method sufficed for our calculations. According to Widom,16,20 the excess chemical potential over the ideal gas value for both homogeneous and inhomogeneous fluids with an intermolecular potential UN(r bN) may be written in the form
βµex ) -ln(〈exp[-β∆WN+1,N]〉N)rN+1
(3)
if WN(r bN) ) UN(r bN) + Uext bN) is the total potential energy of N (r the fluid in an external potential Uext bN). The canonical MC N (r average 〈exp[-β∆WN+1,N]〉N of the difference of the Boltzmann factors between systems with N and N + 1 particles in eq 3 is measured by virtual insertions of test particles in the N particle system at random positions in the usual way.4 In order to achieve reproducible results for subcritical chemical potential versus density isotherms, extremely long Markov chains were required. Particularly, in the range between the equilibrium vapor and liquid densities, where van der Waalslike loops occur, we typically generated 106-108 MC moves per particle and performed about the same number of virtual particle insertion attempts to measure the chemical potential. The number of MC moves per particle required to reach a given precision (variance) was found to increase significantly with the number of particles in the simulation cell and with decreasing temperature. The precisions of the simulated chemical potentials were analyzed by means of the block average method,21 which permits the estimation of the variance of highly correlated data. Typically, the relative error of the chemical potential is on the order of (1-2) × 10-3. To accelerate the simulations we developed a parallel MPI code that allows the calculation of a complete chemical potential isotherm in a single run on a multiprocessor cluster.
The coexistence densities were calculated by means of a Maxwell-like equal area integration, as discussed in refs 5 and 6 from accurate chemical potential µ versus F isotherms. For the thermodynamic integration, we used smoothed µ(F) isotherms obtained by fitting the simulation results to analytical polynomials of the form k
βµ )
bi(Fσ3)i-2 ∑ i)0
(4)
with adjustable parameters bi and k ) 11. 3. Finite Size Effects on Subcritical Bulk Square-Well Fluids 3.1. Chemical Potentials. We have simulated chemical potential versus density isotherms of three-dimensional bulk square-well fluids for a number of temperatures in the subcritical vapor/liquid coexistence range. Figure 1 shows typical van der Waals loops, which indicate phase coexistence between a liquid and a gas phase, for systems of 1000 particles in a temperature range between T* ) 0.95 and T* ) 1.2. For comparison, the results for systems of 216 particles are included in Figure 1 for selected temperatures. The temperature range investigated is below the well-known reduced critical temperature of T/c ) 1.22 for the square-well fluid, which has been estimated by several authors using grand canonical simulations combined with finite size scaling techniques22 or by phase equilibrium simulations combined with scaling assumptions (Wegner expansions,23,24 see also section 4.2) near the critical point.5,15,25 Comparing the 216-particle curves with the 1000-particle curves, we observe a significant increase in the size effects with decreasing temperature. For selected subcritical isotherms, we studied the size dependence in detail by a series of simulations with different numbers of particles, ranging from N ) 64 up to N ) 3375. In Figure 2, we present the results for isotherms at a reduced temperature of T* ) 1.05 and particle numbers N ) {64, 216, 1000, 3375}. Significant finite size effects are observed in the van der Waals loops. The size of the loops decreases with increasing number of particles N, and all isotherms intersect with high accuracy in a single density point. Furthermore, the number of MC steps per particle required to reach comparable accuracy increases strongly with increasing N. All these findings are compatible with the expected general properties of the chemical potential in the gas/liquid transition range8,9 and in qualitative agreement with the recent results of grand canonical MC simulations10,11 for Lennard-Jones fluids in the slightly subcritical temperature range. Figure 3 shows the simulation results of the chemical potential versus the number of particles for the same temperature at selected densities comprising the stable gas and liquid phase and the liquid-gas density gap. We again find strong finite size
Chemical Potential Isotherms in Square-Well Fluids
Figure 4. Vapor/liquid coexistence densities versus the number of particles for the bulk square-well fluid at T* ) 1.05; liquid: f simulation, solid line - literature;15 gas: × - simulation, solid line literature.15
effects on the chemical potential in the liquid-gas density gap, while we observe, in agreement with the literature,3 only weak effects in the stable gas and liquid phases. 3.2. Fluid Phase Equilibria. As mentioned previously, we have calculated the liquid and vapor coexistence densities and the corresponding equilibrium chemical potentials from the subcritical chemical potential isotherms by means of a Maxwelllike equal area rule. In Table 1 we present our results for the N ) 1000 particle and N ) 216 particle systems and compare them with the results of del Rio et al.,15 which show an accurate overall description of the coexistence data. Even the 216-particle results are close to the literature values. In our previous paper5 we estimated the critical data from the phase equilibria results obtained by µ(F) simulations of systems with a moderate number of particles (200-500) and found good agreement with literature data. The liquid-gas coexistence chemical potential (included in Figure 3) also shows, similarly to the stable gas and liquid chemical potentials, only a weak size dependence, and reaches with increasing number of particles the best known literature value,15 which is also included in Figure 3. In Figure 4 we present the vapor/liquid coexistence densities as a function of the number of particles for the bulk squarewell fluid at T* ) 1.05. We find almost insignificant finite size effects on the coexistence densities. The numerical values agree very well with the best known literature data for the phase equilibrium of the square-well fluid (del Rio et al.15), which are included as the horizontal lines in the figure. Only for the smallest number of particles studied (N ) 64) is a small, but significant, deviation found. 4. Finite Size Effects on Subcritical Two-Dimensional Square-Well Fluids 4.1. Chemical Potentials. Analogously to the bulk squarewell system, we have simulated chemical potential isotherms of fluid square-well monolayers (two-dimensional square-well fluids) for several subcritical temperatures. In Figure 5 we have plotted some typical critical/subcritical isotherms for a system with 200 particles in the range T* ) 0.54 to T* ) 0.58. For selected isotherms we have studied a wide range of particle numbers (ranging from N ) 100 to N ) 3200) in the subcritical vapor/liquid two-phase region. The results for particle numbers N ) {100, 200, 800, 1600, 3200} at a reduced temperature of T* ) 0.55 are presented in Figure 7. We observe
J. Phys. Chem. B, Vol. 112, No. 15, 2008 4659
Figure 5. Subcritical chemical potential versus density isotherms for the two-dimensional square-well fluid (N ) 200); T* ) 0.58: × simulation; T* ) 0.565: b - simulation; T* ) 0.55: 1 - simulation; T* ) 0.54: 2 - simulation. Solid lines are numerical fits to eq 4.
Figure 6. Temperature-density vapor-liquid coexistence curve for the two-dimensional square-well fluid (N ) 200); × - simulation results + Maxwell integration; 0 - results of independent Gibbs simulations; solid line - numerical fit to eq 5 for the full temperature range studied (T* ) 0.5-0.58); dotted line - fit only for the restricted temperature range (T* ) 0.55-0.58) close to the critical point; b - critical point estimated from the full temperature range; f - critical point estimated from the restricted temperature range). The rectilinear line (Fl + Fv)/2 is also shown for both the simulation results (symbols) and the fits (lines).
van der Waals-like loops with distinct finite size effects decreasing with increasing particle numbers, qualitatively similar to the results for the bulk fluid. However, the loops are less symmetric and more scattered than in case of the bulk fluid. For the largest system size with N ) 3200 particles, the loop nearly disappears, and a horizontal line, which represents the gas/liquid coexistence equilibrium chemical potential, connects the stable gas and liquid phases (compare also the next section and Figure 9). In the gas/liquid density gap, for large systems (N ) 1600 and N ) 3200), extremely long MC chains (up to 2 × 108 test particle insertions/per particle) were necessary to reach stable (equilibrated) values of the chemical potential. 4.2. Fluid Phase Equilibria. Similarly to the case of the bulk fluid, we present here the liquid/vapor coexistence densities of the monolayer (two-dimensional fluid). Somewhat surprisingly, for this simple system we found that no accurate literature data for phase equilibria and critical data are available. We estimated
4660 J. Phys. Chem. B, Vol. 112, No. 15, 2008
Figure 7. Chemical potential versus density isotherms for various numbers of particles for the two-dimensional square-well fluid at T* ) 0.55; N ) 100: × - simulations; N ) 200: 2 - simulations; N ) 800: b - simulations; N ) 1600: f - simulations; N ) 3200: ( simulations. Solid lines are numerical fits to eq 4.
Figure 8. Coexistence densities versus number of particles for the two-dimensional square-well fluid, at T* ) 0.55; liquid: × simulations, 1 - Gibbs ensemble simulation; vapor: f - simulations, 2 - Gibbs ensemble simulation.
from our chemical potential isotherms for a small system (N ) 200) the vapor/liquid equilibrium densities Fl/v in the same way as in case of the bulk square-well fluid5 using the equal area rule. The results are shown in Figure 6. Since there are no accurate comparison data in the literature, we performed additional Gibbs ensemble simulations for the coexistence densities18 as an independent check of our simulations. The results are included in Figure 6 and Figure 8 and found to be consistent with our results obtained by the Maxwell integration of the chemical potentials. Also, we included in Figure 6 the critical density Fc and the critical temperature T/c obtained by means of a nonlinear regression procedure to fit the coexistence densities to a scaling expression
Fl/v ) Fc + C2(1 - T*/T/c ) ( 0.5B0(1 - T*/T/c )β; T < Tc (5) based on a Wegner expansion,23 where, following the common practice, higher order terms were neglected.24,25 The constants C2, B0, and the critical exponent β are treated as adjustable parameters. The best fit extended over the com-
Vo¨rtler et al.
Figure 9. Estimation of the coexistence chemical potential of the twodimensional square well fluid for the largest system studied, with N ) 3200 particles at T* ) 0.55: × - simulations, solid line - linear fit.
Figure 10. Coexistence chemical potential versus number of particles of the two-dimensional square-well fluid at T* ) 0.55; b are the simulation results, and the solid line shows a fit to an exponential rise to a maximum (infinite system). Error bars for the simulations are included.
plete temperature range studied (T* ) 0.5-0.58) results in18
Fc ≈ 0.37; T/c ≈ 0.58; β ≈ 0.167 Calculating the critical data from a restricted temperature range near the critical point (T* ) 0.55-0.58), the critical temperature and density are only slightly changed, but we obtain as an effective critical exponent β ≈ 0.127, which is very close to the analytical value β ) 1/8 of a corresponding two-dimensional Ising model (lattice gas).14 The estimated critical temperature T/c ≈ 0.58 is also not far from to the exact Ising model value of T* ) 0.567. Figure 8 shows, analogously to Figure 3, the coexistence densities of the monolayer fluid versus the number of particles at the subcritical temperature T* ) 0.55. Since the isotherm for the largest system with N ) 3200 is nearly a horizontal line (compare Figure 7), the estimation of coexistence densities by means of the Maxwell equal area rule was not possible. The finite size effects for the monolayer are found to be significantly larger than for the bulk case. If we extrapolate the vapor coexistence density Fvσ3 and the liquid coexistence density
Chemical Potential Isotherms in Square-Well Fluids Flσ3 to the infinite system by means of exponential fits, we obtain
Fvσ3 ≈ 0.095; Flσ3 ≈ 0.61 Figure 9 illustrates the estimation of the coexistence chemical potential for the largest number of particles (N ) 3200) studied, for which the van der Waals loop in the µ(F) isotherm nearly disappears, and approaches with good approximation a horizontal line (with some scatter), which represents the coexistence value of the chemical potential. Fitting the chemical potential to a linear function results in an equilibrium value of µ/kT ≈ -3.734. The shrinking of the µ(F) loop indicates the approach to the infinite system. In Figure 10 we present the simulation results for the coexistence chemical potential versus the number of particles, which follows within the expected numerical uncertainties an exponential rise to a maximum µ/kT ≈ -3.737(4), which may be interpreted as the value of the chemical potential of an infinite extended system, i.e., for N f ∞. 5. Discussion and Conclusions For both bulk square-well fluids and two-dimensional layers we find significant finite size effects in the subcritical µ(F) isotherms, where the size of the observed van der Waals-like loops decreases with increasing particle numbers N. Bulk square-well fluids show distinct finite size effects for the chemical potential in the liquid-gas density gap, but weak effects in the stable gas and liquid phases. The influence of the system size on the coexistence properties of (three-dimensional) bulk fluids is also weak. The results for the coexistence densities and chemical potentials obtained by our Maxwell-like area rule are in excellent agreement with literature values (del Rio et al.15). We conclude from these findings that, for an accurate simulation of phase equilibria of (square-well-like) hard-core fluids, the use of systems with particle numbers in the order of N ≈ 1000 should be sufficient. Even for system sizes with a few hundred particles, the estimated phase equilibria data should be considered to be good approximations of the infinite system properties. This is in agreement with our previous results for square-well systems,5 and justifies the study of more complex aqueous fluids with a moderate number of particles (N ) 200500).6,7 Comparing our bulk results with recent findings of finite size effects on µ(F) van der Waals loops in bulk Lennard-Jones fluids (MacDowell et al.10,11), we observe quite similar results for square-well fluids, including an increasing scatter of the data with increasing system size, which may indicate droplet formation. Although we did not study in detail the range of lower subcritical temperatures, where MacDowell et al. found sharp discontinuities connected with droplet formation, our results for the lowest temperatures seem to be consistent with such findings (compare the relatively large scatter in our data for T* ) 0.95 and N ) 1000 included in Figure 1). Our results for fluid square-well monolayers are new, and they show significantly larger size effects in both the chemical potentials and the derived phase equilibria data. In contrast to grand canonical simulations, we performed NVT simulations with a preset number of particles, in connec-
J. Phys. Chem. B, Vol. 112, No. 15, 2008 4661 tion with virtual particle insertion to measure the chemical potential, which are conceptually much simpler, and provide chemical potential finite size effects in the NVT ensemble in terms of the particle numbers, which has to be distinguished from the grand canonical results, where the size effects are given in terms of the linear extension of the system. The simulations become extremely time-consuming with increasing number of particles and decreasing temperatures, requiring high performance (parallel) computing facilities for their implementation. The simulated chemical potential versus density isotherms represent important quasi-experimental reference data for a detailed theoretical study of the properties of the chemical potential, particularly in the subcritical vapor/liquid-phase transition range. We currently use the simulated chemical potentials as a benchmark for analytical statistical mechanical theories of the square-well fluid, such as perturbed-hard-sphere equations of state and perturbation theories, where we study in detail the physical meaning of the van der Waals loops in the context of statistical mechanics and simulations. The results of these studies will be presented in a subsequent paper. Acknowledgment. The authors thank Professor Kurt Binder, Mainz, for helpful discussions, which brought to our attention recent results on finite size effects. The financial support of research stays of H.L.V. and K.S. at UOIT Oshawa under the Natural Sciences and Engineering Research Council of Canada Grant OPG1041, and the use of the computer resources of the SHARCNET high-performance computing facilities (http:// www.sharcnet.ca) are gratefully acknowledged. References and Notes (1) Allen, M. P.; Tildesley, D. J. Computer Simulations of Liquids; Clarendon: Oxford, 1987. (2) Siepmann, J. I.; McDonald, I. R.; Frenkel, D. Finite-size corrections to the chemical potential. J. Phys.: Condens. Matter 1992, 4, 679. (3) Kolafa, J.; Vo¨rtler, H. L.; Aim, K.; Nezbeda, I. Mol. Simul. 1993, 11, 305. (4) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications; Academic Press: San Diego, CA, 2002. (5) Vo¨rtler, H. L.; Smith, W. R. J. Chem. Phys. 2000, 112, 5168. (6) Vo¨rtler, H. L.; Kettler, M. Chem. Phys. Lett. 2003, 377, 557. (7) Vo¨rtler, H. L.; Kettler, M. Mol. Phys. 2006, 104, 233. (8) Furukawa, H.; Binder, K. Phys. ReV. A 1982, 26, 556. (9) Kaski, K.; Binder, K.; Gunton, J. D. Phys. ReV. B 1984, 29, 3996. (10) MacDowell, L. G.; Virnau, P.; Mu¨ller, M.; Binder, K. J. Chem. Phys. 2004, 120, 5293. (11) MacDowell, L. G.; Chen, V. K.; Errington, J. R. J. Chem. Phys. 2006, 125, 034705. (12) Yethiraj, A.; Hall, C. K. J. Chem. Phys. 1991, 95, 1999. (13) Nezbeda, I.; Slovak, J. Mol. Phys. 1997, 90, 353. (14) Chandler, D. Introduction to Modern Statistical Mechanics; Oxford University Press: Oxford, 1987. (15) del Rio, F.; Avalos, E.; Espindola, R.; Rull, L. F.; Jackson, G.; Lago, S. Mol. Phys. 2002, 100, 2531. (16) Widom, B. J. Chem. Phys. 1963, 39, 2808. (17) Tripathi, S.; Chapman, W. G. Mol. Phys. 2003, 101, 1199. (18) Scha¨fer, K. Diploma Thesis, Universita¨t Leipzig, 2006. (19) Labik, S.; Malijevsky, A.; Kao, R.; Smith, W. R.; del Rio, F. Mol. Phys. 1999, 96, 849. (20) Widom, B. J. Stat. Phys. 1978, 19, 563. (21) Flyvbjerg, H.; Peterson, H. G. J. Chem. Phys. 1989, 91, 461. (22) Orkoulas, G.; Panagiotopoulos, A. Z. J. Chem. Phys. 1999, 110, 1581. (23) Wegner, F. J. Phys. ReV. 1972, 5B, 4529. (24) Green, D. G.; Jackson, G. J. Chem. Phys. 1994, 101, 3190. (25) Vega, L.; de Miguel, E.; Rull, L. F.; Jackson, G.; McLure, I. A. J. Chem. Phys. 1992, 96, 2296.