Widom Test Particle Insertion Method Research Articles

Molecular Simulations of Hydrogen Sorption in Semicrystalline High-Density Polyethylene: The Impact of the Surface Fraction of Tie-Chains.

The modeling of the barrier properties of semicrystalline polymers has gained interest following the possible application of such materials as protective liners for the safe supply of pressurized hydrogen. The mass transport in such systems is intimately related to the complex intercalation between the crystal and amorphous phases, which was approached in this work through an all-atom representation of high-density polyethylene structures with a tailored fraction of amorphous-crystalline connections (tie-chains). Simulations of the polymer pressure-volume-temperature data and hydrogen sorption were performed by means of molecular dynamics and the Widom test particle insertion method. The discretization of the simulation domains of the semicrystalline structures allowed us to obtain profiles of density, degree of order, and gas solubility. The results indicated that the gas sorption in the crystalline regions is negligible and that the confinement of the amorphous phase between crystals induces a significant increase in density and a drop in the sorption capacity, even in the absence of tie-chains. Adding ties between the crystal and the amorphous phase results in further densification, an increase of the lamella tilt angle, and a decrease in the degree of crystallinity and hydrogen sorption coefficient, in agreement with several literature references.

Rapid screening of gas solubility in ionic liquids using biased particle insertions with pre-sampled liquid trajectories

ABSTRACTWe present an efficient, general-purpose variant of the Widom test particle insertion method for computing chemical potentials of gaseous solutes in fluids or porous solids. The method is implemented in the Monte Carlo molecular simulation engine Cassandra, but receiving phase configurations are independent of this process and may be pre-sampled by other molecular simulation engines such as molecular dynamics codes. Efficiency enhancements present in this method include configurational biasing and accelerated atomic overlap detection. When applied to the estimation of Henry's law constants of atomistic difluoromethane and pentafluoroethane in ionic liquids, the accelerated overlap detection results in a speedup of more than an order of magnitude compared to conventional methods without sacrificing accuracy. We found good agreement between this method and Hamiltonian replica exchange (HREX) for Henry's law constant and absorption isotherm estimation. This embarrassingly parallel method is especially well suited for screening Henry's law constants of many small gases in the same solvents, since a liquid trajectory can be reused for as many solutes as desired.

Thermodynamics of cavity formation in different solvents: Enthalpy, entropy, and the solvophobic effects

For a long time thermodynamics of cavity formation in non-aqueous solvents remained virtually unexplored. We evaluated the entropic and enthalpic contributions into the Gibbs free energy of cavity formation in 22 structurally different organic solvents and water at 298 K. The calculations were carried out using 50–100 ns long molecular dynamics trajectories of solvents and Widom test particle insertion method with up to 1011 insertions per trajectory. It is shown that protic solvents with three-dimensional hydrogen bond networks, namely water, formamide, glycerol, ethylene glycol, propylene glycol, and diethylene glycol, behave differently from other considered solvents. They are characterized with lower values of the entropy of cavity formation than other considered solvents, which also leads to the higher Gibbs free energies. Water has a particularly low entropy as well as the lowest enthalpy of cavity formation at room temperature. We point out the similar tendencies observed for the thermodynamic functions of the cavities and solvation of real molecules in the considered solvents using entropy vs enthalpy and Gibbs free energy vs enthalpy plots. The cavity term governs the difference in solvation properties of protic hydrogen-bonded and aprotic solvents giving rise to the hydrophobic and more general solvophobic effects.

Molecular Simulations and Mechanistic Analysis of the Effect of CO2 Sorption on Thermodynamics, Structure, and Local Dynamics of Molten Atactic Polystyrene.

A simulation strategy encompassing different scales was applied to the systematic study of the effects of CO2 uptake on the properties of atactic polystyrene (aPS) melts. The analysis accounted for the influence of temperature between 450 and 550 K, polymer molecular weights (Mw) between 2100 and 31000 g/mol, and CO2 pressures up to 20 MPa on the volumetric, swelling, structural, and dynamic properties of the polymer as well as on the CO2 solubility and diffusivity by performing molecular dynamics (MD) simulations of the system in a fully atomistic representation. A hierarchical scheme was used for the generation of the higher Mw polymer systems, which consisted of equilibration at a coarse-grained level of representation through efficient connectivity-altering Monte Carlo simulations, and reverse-mapping back to the atomistic representation, obtaining the configurations used for subsequent MD simulations. Sorption isotherms and associated swelling effects were determined by using an iterative procedure that incorporated a series of MD simulations in the NPT ensemble and the Widom test particle insertion method, while CO2 diffusion coefficients were extracted from long MD runs in the NVE ensemble. Solubility and diffusivity compared favorably with experimental results and with predictions of the Sanchez–Lacombe equation of state, which was reparametrized to capture the Mw dependence of polymer properties with greater accuracy. Structural features of the polymer matrix were correctly reproduced by the simulations, and the effects of gas concentration and Mw on structure and local dynamics were thoroughly investigated. In the presence of CO2, a significant acceleration of the segmental dynamics of the polymer occurred, more pronouncedly at low Mw. The speed-up effect caused by the swelling agent was not limited to the chain ends but affected the whole chain in a similar fashion.

Thermodynamics and simulation of hard-sphere fluid and solid: Kinetic Monte Carlo method versus standard Metropolis scheme.

The paper aims at a comparison of techniques based on the kinetic Monte Carlo (kMC) and the conventional Metropolis Monte Carlo (MC) methods as applied to the hard-sphere (HS) fluid and solid. In the case of the kMC, an alternative representation of the chemical potential is explored [E. A. Ustinov and D. D. Do, J. Colloid Interface Sci. 366, 216 (2012)], which does not require any external procedure like the Widom test particle insertion method. A direct evaluation of the chemical potential of the fluid and solid without thermodynamic integration is achieved by molecular simulation in an elongated box with an external potential imposed on the system in order to reduce the particle density in the vicinity of the box ends. The existence of rarefied zones allows one to determine the chemical potential of the crystalline phase and substantially increases its accuracy for the disordered dense phase in the central zone of the simulation box. This method is applicable to both the Metropolis MC and the kMC, but in the latter case, the chemical potential is determined with higher accuracy at the same conditions and the number of MC steps. Thermodynamic functions of the disordered fluid and crystalline face-centered cubic (FCC) phase for the hard-sphere system have been evaluated with the kinetic MC and the standard MC coupled with the Widom procedure over a wide range of density. The melting transition parameters have been determined by the point of intersection of the pressure-chemical potential curves for the disordered HS fluid and FCC crystal using the Gibbs-Duhem equation as a constraint. A detailed thermodynamic analysis of the hard-sphere fluid has provided a rigorous verification of the approach, which can be extended to more complex systems.

Liquid-crystal phase equilibria of Lennard-Jones chains

ABSTRACTLiquid-crystal phase equilibria of Lennard-Jones chain fluids and the solubility of a Lennard-Jones gas in the coexisting phases are calculated from Monte Carlo simulations. Direct phase equilibria calculations are performed using an expanded formulation of the Gibbs ensemble. Monomer densities, order parameters, and equilibrium pressures are reported for the coexisting isotropic and nematic phases of: (1) linear Lennard-Jones chains, (2) a partially-flexible Lennard-Jones chain, and (3) a binary mixture of linear Lennard-Jones chains. The effect of chain length is determined by calculating the isotropic-nematic coexistence of linear Lennard-Jones chain fluids made of 8, 10, and 12 segments (8-, 10-, 12-mer). The effect of molecular flexibility on the isotropic-nematic equilibrium is studied for a Lennard-Jones 10-mer chain fluid with one freely-jointed segment at the end of the chain. An isotropic-nematic phase split and fractionation are reported for a binary mixture of linear 7-mer and 12-mer chains. Simulation results are compared with theoretical results as obtained from a recently developed analytical equation of state based on perturbation theory. Excellent agreement between theory and simulations is observed. The solubility of a monomer Lennard-Jones gas in the coexisting isotropic and nematic phases is estimated using the Widom test-particle insertion method. A linear relationship between solubility difference and density difference at isotropic-nematic coexistence is observed. It is shown that gas solubility is independent of the nematic ordering of the fluid, at constant temperature and density conditions.

Molecular Simulation Studies of the Diffusion of Methane, Ethane, Propane, and Propylene in ZIF-8

ZIF-8 is a strong candidate for propane/propylene separation, which is regarded as one of the most industrially demanding. Molecular simulation of this separation must account for the flexibility of the structure, which enables the adsorption and diffusion of molecules with kinetic diameter larger than the apertures of the pores. Moreover, this simulation requires modeling subtle changes since the strong sieving effect upon the mixture depends on the very small differences between propane and propylene molecular sizes (∼0.2 A). In this work, a new force-field for the ZIF-8 structure has been developed from DFT calculations in simplified structures. The new parameter set reproduces structural properties in very good agreement with the experimental measurements reported in literature. Molecular dynamics simulations and the Widom test particle insertion method were then employed for the calculation of diffusivities, activation energies and adsorption properties of propane and propylene. The results are in ag...

Isotropic-nematic phase equilibria of hard-sphere chain fluids-Pure components and binary mixtures.

The isotropic-nematic phase equilibria of linear hard-sphere chains and binary mixtures of them are obtained from Monte Carlo simulations. In addition, the infinite dilution solubility of hard spheres in the coexisting isotropic and nematic phases is determined. Phase equilibria calculations are performed in an expanded formulation of the Gibbs ensemble. This method allows us to carry out an extensive simulation study on the phase equilibria of pure linear chains with a length of 7 to 20 beads (7-mer to 20-mer), and binary mixtures of an 8-mer with a 14-, a 16-, and a 19-mer. The effect of molecular flexibility on the isotropic-nematic phase equilibria is assessed on the 8-mer+19-mer mixture by allowing one and two fully flexible beads at the end of the longest molecule. Results for binary mixtures are compared with the theoretical predictions of van Westen et al. [J. Chem. Phys. 140, 034504 (2014)]. Excellent agreement between theory and simulations is observed. The infinite dilution solubility of hard spheres in the hard-sphere fluids is obtained by the Widom test-particle insertion method. As in our previous work, on pure linear hard-sphere chains [B. Oyarzún, T. van Westen, and T. J. H. Vlugt, J. Chem. Phys. 138, 204905 (2013)], a linear relationship between relative infinite dilution solubility (relative to that of hard spheres in a hard-sphere fluid) and packing fraction is found. It is observed that binary mixtures greatly increase the solubility difference between coexisting isotropic and nematic phases compared to pure components.

The phase behavior of linear and partially flexible hard-sphere chain fluids and the solubility of hard spheres in hard-sphere chain fluids

The liquid crystal phase behavior of linear and partially flexible hard-sphere chain fluids and the solubility of hard spheres in hard-sphere chain fluids are studied by constant pressure Monte Carlo simulations. An extensive study on the phase behavior of linear fluids with a length of 7, 8, 9, 10, 11, 12, 13, 14, 15, and 20 beads is carried out. The phase behavior of partially flexible fluids with a total length of 8, 10, 14, and 15 beads and with different lengths for the linear part is also determined. A precise description of the reduced pressure and of the packing fraction change at the isotropic-nematic coexistence was achieved by performing long simulation runs. For linear fluids, a maximum in the isotropic to nematic packing fraction change is observed for a chain length of 15 beads. The infinite dilution solubility of hard spheres in linear and partially flexible hard-sphere chain fluids is calculated by the Widom test-particle insertion method. To identify the effect of chain connectivity and molecular anisotropy on free volume, solubility is expressed relative to that of hard spheres in a hard sphere fluid at same packing fraction as relative Henry's law constants. A linear relationship between relative Henry's law constants and packing fraction is observed for all linear fluids. Furthermore, this linearity is independent of liquid crystal ordering and seems to be independent of chain length for linear chains of 10 beads and longer. The same linear relationship was observed for the solubility of hard spheres in nematic forming partially flexible fluids for packing fractions up to a value slightly higher than the nematic packing fraction at the isotropic-nematic coexistence. At higher packing fractions, the small flexibility of these fluids seems to improve solubility in comparison with the linear fluids.

Calculating thermodynamic factors of ternary and multicomponent mixtures using the Permuted Widom test particle insertion method

We present an approach to calculate thermodynamic factors for multicomponent liquid mixtures from molecular simulations by the Permuted Widom test particle insertion method. The thermodynamic factors relate experimentally obtained Fick diffusion coefficients and the Maxwell–Stefan diffusion coefficients computed from molecular simulations/theory. In our earlier work, we presented the Permuted Widom test particle insertion method (Balaji et al. in Mol Phys 2012) for calculating thermodynamic factors for binary systems from a single simulation. In this paper, we have extended this approach to ternary systems consisting of Lennard–Jones particles. The thermodynamic factors computed from the simulations and from numerically differentiating the activity coefficients obtained from the conventional Widom single test particle insertion method were compared. An excellent agreement was observed.

A direct method for calculating thermodynamic factors for liquid mixtures using the Permuted Widom test particle insertion method

Understanding mass transport in liquids by mutual diffusion is an important topic for many applications in chemical engineering. The reason for this is that diffusion is often the rate limiting step in chemical reactors and separators. In multicomponent liquid mixtures, transport diffusion can be described by both generalized Fick's law and the Maxwell–Stefan theory. The Maxwell–Stefan and Fick approaches in an n-component system are related by the so-called thermodynamic factor [R. Taylor and H.A. Kooijman, Chem. Eng. Commun, 102, 87 (1991)]. As Fick diffusivities can be measured in experiments and Maxwell–Stefan diffusivities can be obtained from molecular simulations/theory, the thermodynamic factors bridge the gap between experiments and molecular simulations/theory. It is therefore desirable to be able to compute thermodynamic factors from molecular simulations. Unfortunately, presently used simulation techniques for computing thermodynamic factors are inefficient and often require numerical differentiation of simulation results. In this work, we propose a modified version of the Widom test-particle method to compute thermodynamic factors from a single simulation. This method is found to be more efficient than the conventional Widom test particle insertion method combined with numerical differentiation of simulation results. The approach is tested for binary systems consisting of Lennard-Jones particles. The thermodynamic factors computed from the simulation and from numerically differentiating the activity coefficients obtained from the conventional Widom test particle insertion method are in excellent agreement.

A Molecular Dynamics Technique to Extract Forces in Soft Matter Systems Under Compression With Constant Solvent Chemical Potential.

Molecular dynamics simulations of opposing polymer brushes at varying surface separation distances were performed to develop a method for conducting a static compression of soft matter. As all separation distances were represented by independent simulations, the proper solvent density for every level of compression needed to be determined to acquire realistic data. This was accomplished by maintaining a constant solvent chemical potential for each separation distance. In doing so, each independent simulation is equilibrated with all others, reproducing conditions encountered experimentally in force spectroscopy measurements. Chemical potential was determined using the Widom test particle insertion method. Force information was extracted from pressure profiles, such that unphysical forces occurring within the surface layers were not accounted for in the calculation. Each individual simulation was a canonical ensemble molecular dynamics simulation, but taken together they approximate a grand canonical ensemble for the solvent particles by holding their chemical potential constant.

On the Behavior of Solutions of Xenon in Liquid n-Alkanes: Solubility of Xenon in n-Pentane and n-Hexane

The solubility of xenon in liquid n-pentane and n-hexane has been studied experimentally, theoretically, and by computer simulation. Measurements of the solubility are reported for xenon + n-pentane as a function of temperature from 254 to 305 K. The uncertainty in the experimental data is less than 0.15%. The thermodynamic functions of solvation such as the standard Gibbs energy, enthalpy, and entropy of solvation have been calculated from Henry's law coefficients for xenon + n-pentane solutions and also for xenon + n-hexane, which were reported in previous work. The results provide a further example of the similarity between the xenon + n-alkane interaction and the n-alkane + n-alkane interactions. Using the SAFT-VR approach we were able to quantitatively predict the experimental solubility for xenon in n-pentane and semiquantitatively that of xenon in n-hexane using simple Lorentz-Berthelot combining rules to describe the unlikely interaction. Henry's constants at infinite dilution for xenon + n-pentane and xenon + n-hexane were also calculated by Monte Carlo simulation using a united atom force field to describe the n-alkane and the Widom test particle insertion method.

Molecular Dynamics Simulation of Sorption of Gases in Polystyrene

A computational method for the calculation of the solubility of gases, including argon, hydrogen, nitrogen, carbon dioxide, methane, and propane, in polystyrene over a wide range of temperatures and pressures is described. The excess chemical potentials and the partial molar volumes of gases in polystyrene were calculated using Widom's test-particle insertion method. Using the calculated chemical potentials and partial molar volumes of the sorbed gases at a fixed temperature, the chemical potentials in the polymer phase were expanded in terms of the pressure. Performing a grand-canonical ensemble molecular dynamics simulation of the gas phase, we calculated the phase coexistence point using a recent method by Vrabec and Hasse [Mol. Phys. 2002, 100, 3375−3383]. The results on the calculated solubility coefficients and solubilities over a wide range of temperatures and pressures are compared with experimental data.

Investigation of the Salting Out of Methane from Aqueous Electrolyte Solutions Using Computer Simulations

We calculate the excess chemical potential of methane in aqueous electrolyte solutions of NaCl using Monte Carlo computer simulations. In a recent work [Docherty et al. J. Chem. Phys. 2006, 125, 074510], we presented a new potential model for methane in water which is capable of describing accurately the excess chemical potential of methane in pure water over a range of temperatures, a quantity that can be related to the solubility and which is commonly used to study the hydrophobic effect. Here, we use the same potential model for the water-methane interactions and investigate the effect of added salt on the chemical potential of methane in the solution. The methane molecules are modeled as single Lennard-Jones (LJ) interaction sites, and the water molecules are modeled with the TIP4P/2005 model. A correcting factor of chi = 1.07 for the energetic Berthelot (geometric) combining rule of the methane-water interaction is also used, which mimics the polarization of methane in water. We consider NaCl as the salt and treat the ions with the Smith and Dang model (i.e., as charged LJ interaction sites). Ion-water, ion-ion, and ion-methane interactions are treated using Lorentz-Berthelot combining rules. In addition, the Coulombic potential is used to model charge-charge interactions which are calculated using the Ewald sum. We have carried out isobaric-isothermal (NpT) simulations to determine the equilibrium densities of the solutions. The simulation data is in excellent agreement with experimental densities of aqueous NaCl solutions of different concentration. Hydration numbers are also obtained and found to be in agreement with reported data. Canonical (NVT) simulations at the averaged densities are then performed using the Widom test-particle insertion method to obtain the excess chemical potential of methane in the saline solutions. An increase in the chemical potential of methane, corresponding to a salting out effect, is observed when salt is added to the solution. We investigate different concentrations and ion sizes. An overprediction of the salting out effect as compared with experimental data is observed, which we believe is due to the polarizing effect of the ions in the solution, which is not taken into account by the model. We also find a direct correlation between the increase in the chemical potential and the packing fraction of the solution and argue that the main cause of the observed salting out effect (as represented by an increase in the excess chemical potential) is the increase in the packing fraction of the solutions due to the added salt. Together, with this, we put forward an argument toward explaining the anomalous Hofmeister effect of Li(+).

Atomistic Simulation of Poly(dimethylsiloxane): Force Field Development, Structure, and Thermodynamic Properties of Polymer Melt and Solubility of n-Alkanes, n-Perfluoroalkanes, and Noble and Light Gases

Poly(dimethylsiloxane) (PDMS) is a widely used polymer for a number of industrial applications. In order that PDMS is selected for a specific application, accurate knowledge of its physical properties is necessary. Physical properties can be either measured or calculated based on reliable suitable methods. Molecular simulation using realistic models is a powerful tool for the elucidation of microscopic structure of polymers and the subsequent estimation of macroscopic physical properties. In this work, a force field is developed for the prediction of thermodynamic and structure properties of PDMS melts. Force field development was based on existing force fields for PDMS together with fitting to experimental thermodynamic data at ambient pressure. Extensive NPT molecular dynamics (MD) simulations were performed at different temperature and pressure values. In all cases, good agreement was obtained between literature experimental data and model predictions for the melt density. Calculations are reported also for the solubility parameter of the polymer melt at different temperatures. Furthermore, radial distribution functions for the intra- and intermolecular interactions are presented and shown to be in good agreement with previous literature work. The new force field is used subsequently for the calculation of solubility of 17 different compounds in PDMS using the Widom test particle insertion method. The solubility of n-alkanes from methane to n-hexane at 300 and 450 K and different pressures was calculated. In addition, solubility calculations for n-perfluoroalkanes at 300 and 450 K and for noble and light gases at 300, 375, and 450 K and ambient pressure were performed. Model predictions are in very good agreement with experimental data, in all cases. The infinite dilution solubility coefficient is shown to increase with temperature for very light gases and decrease for the heavier ones.

Analysis of Henry’s constant for carbon dioxide in water via Monte Carlo simulation

We present calculations of the Henry constant for carbon dioxide in water by Monte Carlo simulations over a broad range of temperatures, from 0 ∘C to the critical temperature of water. A range of intermolecular potential models is examined for each species. Carbon dioxide is modelled by two three-site (EPM2, and Errington and Panagiotopoulos) potentials and water is modelled by four three-site (SPC, SPC/E, MSPC/E, and Errington and Panagiotopoulos) potentials, by the four-site TIP4P potential and by the five-site TIP5P potential. Henry’s constant is computed via the Widom test-particle insertion method and by means of a staged free-energy perturbation method. The performance of the various potential models with respect to the accuracy of their prediction of the Henry constant is discussed. The staged free-energy perturbation method, employed at several representative temperatures, allows further analysis of the Henry constant with respect to the free energy of cavity formation for hosting the CO 2 solute molecule in the H 2O solvent and the free energy of interactions between the CO 2 solute molecule and the H 2O solvent. We found that all CO 2/H 2O models predicted a qualitatively correct temperature dependence of the Henry constant but only the Errington and Panagiotopoulos CO 2/H 2O model gave values for the Henry constant in reasonable agreement with experimental data.

Calculation of the chemical potential of chain molecules using the staged particle deletion scheme

A formulation is presented for the calculation of the chemical potential of chain molecules from molecular simulation based on the staged deletion of a test molecule. This formulation is an extension of a formulation presented recently [Boulougouris et al., Mol. Phys. 96, 905 (1999)] for the case of spherical molecules. An analytical method for the calculation of the volume accessible to a hard-sphere chain molecule is used together with the excluded volume map sampling technique. The new method is applied to a fluid of tangent sphere dimers and to ethane under various conditions. NPT, NVT, and Gibbs ensemble Monte Carlo simulation results are presented. Extensive comparison is made against calculations with the Widom test particle insertion method. In all cases, the new method results in considerable savings in CPU time.

A new Monte Carlo simulation approach for the prediction of sorption equilibria of oligomers in polymer melts: Solubility of long alkanes in linear polyethylene

This work is concerned with the study of the solubility of small alkanes and ethylene oligomers in polyethylene melts of prescribed molecular weight distribution. The infinite dilution case is examined first, using Widom test particle insertion method. Simulation estimates of Henry’s constant for alkanes from C4 to C8 in a melt with average chain length C78 at T=450 K and P=1 atm are compared to older simulation results and to experimental data and are found to be adequately accurate. Cases where the oligomer is dissolved in the melt at finite concentrations are of greater practical interest. To address these cases, two different methods are implemented. The first one employs Monte Carlo (MC) simulations in the [f1NpnPTμ*] statistical ensemble, in which conventional insertion and deletion moves are used to equilibrate the concentration of the small molecules in the PE melt. It is observed that this technique leads to satisfactory results only for short alkanes, of length up to C10. Configurationally-biased insertions and deletions of alkane molecules in dense polymer phases become extremely difficult and time-consuming as the length of the inserted molecule is increased, rendering equilibration within reasonable CPU time practically impossible. To study the solubility of longer alkanes and oligomers in a polyethylene melt two novel MC moves are introduced: scission of a polymer chain to create an alkane molecule and a shorter chain, and fusion of a polymer chain and an alkane to create a longer polymer. For implementing these two new moves, a new statistical ensemble, the [f 1′Npn0PTμ*], is formulated. The new moves lead to extremely fast equilibration of the concentration of alkane molecules in the polymer melt and allow studying the solubility of long oligomers in the polymer matrix at high pressures. Their advantage lies exactly in the fact that they generate alkane molecules from already existing polymer chains, thereby eliminating the need for insertion and deletion moves. Results are presented from both MC methods for the sorption isotherms of C5, C10, and C20 in molten linear polyethylene and their dependence on the average chain length of the polymer matrix. The structure of the polymer–alkane mixtures and the volume changes (swelling) of the polyethylene matrix upon sorption are also examined. Predictions obtained from the atomistic simulations are compared extensively to available experimental data and to calculations based on the Flory–Huggins theory and the SAFT equation of state (EoS). The agreement between simulation and experiment is seen to be very good, demonstrating the value of the new simulation scheme designed and implemented in the course of this work.

Computer simulation studies of a square-well fluid in a slit pore. Spreading pressure and vapor–liquid phase equilibria using the virtual-parameter-variation method

We study model square-well fluids with well-width parameter λ=1.5 confined to hard planar slits. We derive a general computer simulation method for numerically calculating an arbitrary first derivative of the canonical ensemble partition function with respect to a simulation parameter, which we call the virtual-parameter-variation method. Two special cases of this approach are the Widom test-particle insertion method for calculating the excess chemical potential, and a method for calculating the pressure due to Eppenga and Frenkel [Mol. Phys. 52, 52, 1303 (1984)]. We use this approach to calculate the volume derivative parallel to the slit walls of the Helmholtz free energy in an (N,V,T) Monte Carlo simulation, and show that this spreading pressure is numerically consistent with the thermodynamic pressure obtained by integration of the Gibbs–Duhem equation using the simulated chemical potentials of the confined fluid as a function of density. We obtain new simulation results for the spreading pressure and the phase equilibrium properties of the confined square-well fluid, and we also estimate its critical point properties, observing a decrease of the critical temperature in comparison to the bulk fluid.