Calculation of surface tension via area sampling

We examine the performance of several molecular simulation techniques aimed at evaluation of the surface tension through its thermodynamic definition. For all methods explored, the surface tension is calculated by approximating the change in Helmholtz free energy associated with a change in interfacial area through simulation of a liquid slab at constant particle number, volume, and temperature. The methods explored fall within three general classes: free-energy perturbation, the Bennett acceptance-ratio scheme, and the expanded ensemble technique. Calculations are performed for both the truncated Lennard-Jones and square-well fluids at select temperatures spaced along their respective liquid-vapor saturation lines. Overall, we find that Bennett and expanded ensemble approaches provide the best combination of accuracy and precision. All of the methods, when applied using sufficiently small area perturbation, generate equivalent results for the Lennard-Jones fluid. However, single-stage free-energy-perturbation methods and the closely related test-area technique recently introduced by Gloor et al. [J. Chem. Phys. 123, 134703 (2005)] generate surface tension values for the square-well fluid that are not consistent with those obtained from the more robust expanded ensemble and Bennett approaches, regardless of the size of the area perturbation. Single-stage perturbation methods fail also for the Lennard-Jones system when applied using large area perturbations. Here an analysis of phase-space overlap produces a quantitative explanation of the observed inaccuracy and shows that the satisfactory results obtained in these cases from the test-area method arise from a cancellation of errors that cannot be expected in general. We also briefly analyze the variation in method performance with respect to the adjustable parameters inherent to the techniques.     

Test-area simulation method for the direct determination of the interfacial tension of systems with continuous or discontinuous potentials

A novel test-area (TA) technique for the direct simulation of the interfacial tension of systems interacting through arbitrary intermolecular potentials is presented in this paper. The most commonly used method invokes the mechanical relation for the interfacial tension in terms of the tangential and normal components of the pressure tensor relative to the interface (the relation of Kirkwood and Buff [J. Chem. Phys. 17, 338 (1949)]). For particles interacting through discontinuous intermolecular potentials (e.g., hard-core fluids) this involves the determination of delta functions which are impractical to evaluate, particularly in the case of nonspherical molecules. By contrast we employ a thermodynamic route to determine the surface tension from a free-energy perturbation due to a test change in the surface area. There are important distinctions between our test-area approach and the computation of a free-energy difference of two (or more) systems with different interfacial areas (the method of Bennett [J. Comput. Phys. 22, 245 (1976)]), which can also be used to determine the surface tension. In order to demonstrate the adequacy of the method, the surface tension computed from test-area Monte Carlo (TAMC) simulations are compared with the data obtained with other techniques (e.g., mechanical and free-energy differences) for the vapor-liquid interface of Lennard-Jones and square-well fluids; the latter corresponds to a discontinuous potential which is difficult to treat with standard methods. Our thermodynamic test-area approach offers advantages over existing techniques of computational efficiency, ease of implementation, and generality. The TA method can easily be implemented within either Monte Carlo (TAMC) or molecular-dynamics (TAMD) algorithms for different types of interfaces (vapor-liquid, liquid-liquid, fluid-solid, etc.) of pure systems and mixtures consisting of complex polyatomic molecules.     

Finite-size effects in dissipative particle dynamics simulations

We have performed dissipative particle dynamics (DPD) simulations to evaluate the effect that finite size of transversal area has on stress anisotropy and interfacial tension. The simulations were carried out in one phase and two phases in parallelepiped cells. In one-phase simulations there is no finite-size effect on stress anisotropy when the simulation is performed using repulsive forces. However, an oscillatory function of stress anisotropy is found for attractive-repulsive interactions. In the case of liquid-liquid interfaces with repulsive interaction between molecules, there is only a small effect of surface area on interfacial tension when the simulations are performed using the Monte Carlo method at constant temperature and normal pressure. An important but artificial finite-size effect of interfacial area on surface tension is found in simulations in the canonical ensemble. Reliable results of interfacial tension from DPD simulations can be obtained using small systems, less than 2000 particles, when they interact exclusively with repulsive forces.     

On interfacial tension calculation from the test-area methodology in the grand canonical ensemble

We propose the extension of the test-area methodology, originally proposed to evaluate the surface tension of planar fluid-fluid interfaces along a computer simulation in the canonical ensemble, to deal with the solid-fluid interfacial tension of systems adsorbed on slitlike pores using the grand canonical ensemble. In order to check the adequacy of the proposed extension, we apply the method for determining the density profiles and interfacial tension of spherical molecules adsorbed in slitlike pore with different pore sizes and solid-fluid dispersive energy parameters along the same simulation. We also calculate the solid-fluid interfacial tension using the original test-area method in the canonical ensemble. Agreement between the results obtained from both methods indicate that both methods are fully equivalent. The advantage of the new methodology is that allows to calculate simultaneously the density profiles and the amount of molecules adsorbed onto a slitlike pore, as well as the solid-fluid interfacial tension. This ensures that the chemical potential at which all properties are evaluated during the simulation is exactly the same since simulations can be performed in the grand canonical ensemble, mimicking the conditions at which the adsorption experiments are most usually carried out in the laboratory.     

Vapor-liquid interfacial properties of rigid-linear Lennard-Jones chains

We have obtained the interfacial properties of short rigid-linear chains formed from tangentially bonded Lennard-Jones monomeric units from direct simulation of the vapour-liquid interface. The full long-range tails of the potential are accounted for by means of an improved version of the inhomogeneous long-range corrections of Janec?ek [J. Phys. Chem. B 110, 6264-6269 (2006)] proposed recently by MacDowell and Blas [J. Chem. Phys. 131, 074705 (2009)] valid for spherical as well as for rigid and flexible molecular systems. Three different model systems comprising of 3, 4, and 5 monomers per molecule are considered. The simulations are performed in the canonical ensemble, and the vapor-liquid interfacial tension is evaluated using the test-area method. In addition to the surface tension, we also obtain density profiles, coexistence densities, critical temperature and density, and interfacial thickness as functions of temperature, paying particular attention to the effect of the chain length and rigidity on these properties. According to our results, the main effect of increasing the chain length (at fixed temperature) is to sharpen the vapor-liquid interface and to increase the width of the biphasic coexistence region. As a result, the interfacial thickness decreases and the surface tension increases as the molecular chains get longer. The surface tension has been scaled by critical properties and represented as a function of the difference between coexistence densities relative to the critical density.     

A unified methodological framework for the simulation of nonisothermal ensembles

A general framework is developed for the simulation of nonisothermal statistical-mechanical ensembles. This framework is intended to synthesize the formulation of advanced Monte Carlo simulation methods such as multihistogram reweighting, replica-exchange methods, and expanded ensemble techniques so that they can be applied to different nonisothermal ensembles. Using Lennard-Jones systems as test cases, novel implementations of these methods are demonstrated with different ensembles including the microcanonical, isobaric-isoenthalpic, and isobaric-semigrand ensembles. In particular, it is shown that the use of multiensemble methods allows the efficient simulation of microcanonical density of states, entropies, vapor-liquid and solid-liquid equilibrium for pure component systems, and fluid-phase coexistence for binary mixtures. In these applications, comparisons are also presented that highlight the advantages of the proposed multiensemble implementations over alternative methods used before.     

Isothermalisobaric molecular dynamics simulation of diatomic liquids and their mixtures

Isothermalisobaric ensemble molecular dynamics simulations have been performed for systems of two-center Lennard-Jones molecules for pure fluids as well as binary mixtures. The results obtained from various ensemble averages have been compared for pure fluid simulations of Lennard-Jones model diatomic fluids. The excess enthalpy and excess volume of three equimolar mixtures (argonnitrogen, argonoxygen, and nitrogenoxygen) have been calculated and compared with values obtained from previous NVT molecular dynamics and perturbation theory. Pair distribution functions obtained from various methods are compared for the equimolar mixture of nitrogen and oxygen and used to study the effect of attractive forces on the local structure. For four other systems (argonethane, methaneethane, carbon disulfidecarbon tetrachloride, and carbon disulfidetetrachloroethylene), excess enthalpies and excess volumes from NPT simulations are used to test the ability of perturbation theory to predict these properties and are also compared with experimental data. The comparison of simulation and experiment clearly shows the need to improve the available parameters for cross interactions in binary mixtures.     

Monte Carlo simulation methods for computing the wetting and drying properties of model systems

We introduce general Monte Carlo simulation methods for determining the wetting and drying properties of model systems. We employ an interface-potential-based approach in which the interfacial properties of a system are related to the surface excess free energy of a thin fluid film in contact with a surface. Two versions of this approach are explored: a "spreading" method focused on the growth of a thin liquid film from a surface in a mother vapor and a "drying" method focused on the growth of a thin vapor film from a surface in a mother liquid. The former provides a direct measure of the spreading coefficient while the latter provides an analogous drying coefficient. When coupled with an independent measure of the liquid-vapor surface tension, these coefficients enable one to compute the contact angle. We also show how one can combine information gathered from application of the spreading and drying methods at a common state point to obtain direct measures of the contact angle and liquid-vapor surface tension. The computational strategies introduced here are applied to two model systems. One includes a monatomic Lennard-Jones fluid that interacts with a structureless substrate via a long-ranged substrate potential. The second model contains a monatomic Lennard-Jones fluid that interacts with an atomistically detailed substrate via a short-ranged potential. Expanded ensemble techniques are coupled with the interface potential approach to compile the temperature- and substrate strength-dependence of various interfacial properties for these systems. Overall, we find that the approach pursued here provides an efficient and precise means to calculate the wetting and drying properties of model systems.     

气体混合物各组分化学势的Monte Carlo模拟

在则系综测试粒子Ｍｏｎｔｅ Ｃａｒｌｏ （ＧＣＭＣ）方法模拟常温下空气（以氮气为代表）及其污染物微量有机物（以苯为例）的混合物中各组分的化学势。模拟中,氮气和苯分子采用ＬＪ球型分子势能模型,采用Ｍｅｔｒｏｐｏｌｉｓ抽样及周期边界条件。通过模拟并拟合得到了３００．２Ｋ、苯的摩尔分数为０．０６２５,氮气及苯化学势与压力的关联式,以用于狭缝碳孔中该混合物体系的选择性吸附。     

Discontinuous molecular dynamics for rigid bodies: applications

Event-driven molecular dynamics simulations are carried out on two rigid-body systems which differ in the symmetry of their molecular mass distributions. First, simulations of methane in which the molecules interact via discontinuous potentials are compared with simulations in which the molecules interact through standard continuous Lennard-Jones potentials. It is shown that under similar conditions of temperature and pressure, the rigid discontinuous molecular dynamics method reproduces the essential dynamical and structural features found in continuous-potential simulations at both gas and liquid densities. Moreover, the discontinuous molecular dynamics approach is demonstrated to be between 3 and 100 times more efficient than the standard molecular dynamics method depending on the specific conditions of the simulation. The rigid discontinuous molecular dynamics method is also applied to a discontinuous-potential model of a liquid composed of rigid benzene molecules, and equilibrium and dynamical properties are shown to be in qualitative agreement with more detailed continuous-potential models of benzene. The few qualitative differences in the angular dynamics of the two models are related to the relatively crude treatment of variations in the discontinuous repulsive interactions as one benzene molecule rotates by another.     

Solid-liquid phase coexistence of the Lennard-Jones system through phase-switch Monte Carlo simulation

The phase-switch Monte Carlo method of Wilding and Bruce [Phys. Rev. Lett. 85, 5138 (2000)] is extended to enable calculation of solid-liquid phase coexistence for soft potentials. The method directly accesses coexistence information about a system while avoiding simulation of the interfacial region. Order parameters are introduced that allow one to define a path that connects liquid and crystalline phases. Transition matrix methods are employed to bias the sampling such that both phases are sampled in a rapid and efficient manner. Coexistence properties are determined through an analysis of specific volume probability distributions, which are generated naturally during a biased simulation. The approach is demonstrated with the Lennard-Jones system. Finite-size effects are examined and compared to those for the hard sphere system. In addition, two techniques are considered for accounting for long-range interactions. The methodology presented here is general and therefore provides a basis for its application to other soft systems.     

Investigation of the interfacial tension of complex coacervates using field-theoretic simulations

Complex coacervation, a liquid-liquid phase separation that occurs when two oppositely charged polyelectrolytes are mixed in a solution, has the potential to be exploited for many emerging applications including wet adhesives and drug delivery vehicles. The ultra-low interfacial tension of coacervate systems against water is critical for such applications, and it would be advantageous if molecular models could be used to characterize how various system properties (e.g., salt concentration) affect the interfacial tension. In this article we use field-theoretic simulations to characterize the interfacial tension between a complex coacervate and its supernatant. After demonstrating that our model is free of ultraviolet divergences (calculated properties converge as the collocation grid is refined), we develop two methods for calculating the interfacial tension from field-theoretic simulations. One method relies on the mechanical interpretation of the interfacial tension as the interfacial pressure, and the second method estimates the change in free energy as the area between the two phases is changed. These are the first calculations of the interfacial tension from full field-theoretic simulation of which we are aware, and both the magnitude and scaling behaviors of our calculated interfacial tension agree with recent experiments.     

Size dependence, stability, and the transition to buckling in model reverse bilayers

Molecular dynamics simulations of a model bilayer made of surfactant dimers in a Lennard-Jones solvent are reported for three sizes of the systems up to an area of 100sigma x 100sigma and for a large interval of the specific areas: from hole formation under tension deep into the floppy state of a buckling compressed bilayer. The transition to the floppy state appears quite abrupt and discontinuous; in the floppy state the lateral tension is negative and scales with size while vanishing from below. The structure factor was also determined for all three sizes and all areas; for most part the apparent tension is larger than the lateral tension whereas the apparent rigidity constant--always positive--is low in the floppy state and increasing in the tensioned state. Both do not scale visibly with size. The replacement of the 1q(2) capillary-wave divergence by another pole is accounted for and explained.     

Phase equilibria and interfacial tension of fluids confined in narrow pores

Correlation between phase behaviors of a Lennard-Jones fluid in and outside a pore is examined over wide thermodynamic conditions by grand canonical Monte Carlo simulations. A pressure tensor component of the confined fluid, a variable controllable in simulation but usually uncontrollable in experiment, is related with the pressure of a bulk homogeneous system in equilibrium with the confined system. Effects of the pore dimensionality, size, and attractive potential on the correlations between thermodynamic properties of the confined and bulk systems are clarified. A fluid-wall interfacial tension defined as an excess grand potential is evaluated as a function of the pore size. It is found that the tension decreases linearly with the inverse of the pore diameter or width.     

Surface tension of a Lennard-Jones liquid under supersaturation

A formally exact Kirkwood-Buff virial formula for the surface tension of a supersaturated interface is derived. A modified Gibbs ensemble method is given that allows the creation of interacting supersaturated phases of equal chemical potential, and which enables the Kirkwood-Buff formula to be applied. The methods are tested by Monte Carlo simulation of a supersaturated Lennard-Jones fluid with a planar liquid-vapor interface. The Kirkwood-Buff results for the supersaturated surface tension are found to be in reasonable agreement with new results obtained here using the recently developed, formally exact, ghost interface method, [M. P. Moody and P. Attard, J. Chem. Phys., 2002, 117, 6705]. The surface tension is obtained as a function of supersaturation at four temperatures, and it is found to decrease with increasing supersaturation, and to vanish at the vapor spinodal. The relevance of the present results to the nucleation of droplets in a supersaturated vapor is discussed.     

Stress anisotropy induced by periodic boundary conditions

Finite size effects due to periodic boundary conditions are investigated using computer simulations in the canonical ensemble. We study liquids with densities corresponding to typical liquid coexistence densities, and temperatures between the triple and critical points. The components of the pressure tensor are computed in order to analyze the finite size effects arising from the size and geometry of the simulation box. Two different box geometries are considered: cubic and parallelepiped. As expected the pressure tensor is isotropic in cubic boxes, but it becomes anisotropic for small noncubic boxes. We argue this is the origin of the anomalous behavior observed recently in the computation of the surface tension of liquid-vapor interfaces. Otherwise, we find that the bulk pressure is sensitive to the box geometry when small simulation boxes are considered. These observations are general and independent of the model liquid considered. We report results for liquids interacting through short range forces, square well and Lennard-Jones, and also long range Coulombic interactions. The effect that small surface areas have on the surface tension is discussed, and some preliminary results at the liquid vapor-interface for the square well potential are given.     

On the long-range corrections to computer simulation results for the Lennard-Jones vapor-liquid interface

The long-range corrections (LRCs) to the configurational energy have been taken into consideration in the Monte Carlo simulation of the vapor-liquid interface for a pure Lennard-Jones (LJ) fluid. The simulated bulk densities agree satisfactorily with those obtained from the Gibbs ensemble method, and the simulated surface tension values agree reasonably well with those reported in the literature for a larger number of molecules and a larger cut-off distance. To compare the influence of the potential forms on the simulation results, a truncated LJ potential, and a shifted and truncated LJ potential have been examined. Although the bulk densities and surface tensions calculated for different model fluids are strongly affected by the LRC, the different potentials essentially lead to similar density values and similar surface tension values when the respective calculated values are compared on the basis of a reduced temperature scale.     

Oscillatory surface tension due to finite-size effects

The simulation results of surface tension at the liquid-vapor interface are presented for fluids interacting with Lennard Jones and square-well potentials. From the simulation of liquids we have reported [M. González-Melchor et al., J. Chem. Phys. 122, 4503 (2005)] that the components of pressure tensor in parallelepiped boxes are not the same when periodic boundary conditions and small transversal areas are used. This fact creates an artificial oscillatory stress anisotropy in the system with even negative values. By doing direct simulations of interfaces we show in this work that surface tension has also an oscillatory decay at small surface areas; this behavior is opposite to the monotonic decay reported previously for the Lennard Jones fluid. It is shown that for small surface areas, the surface tension of the square-well potential artificially takes negative values and even increases with temperature. The calculated surface tension using a direct simulation of interfaces might have two contributions: one from finite-size effects of interfacial areas due to box geometry and another from the interface. Thus, it is difficult to evaluate the true surface tension of an interface when small surface areas are used. Care has to be taken to use the direct simulation method of interfaces to evaluate the predicted surface tension as a function of interfacial area from capillary-wave theory. The oscillations of surface tension decay faster at temperatures close to the critical point. It is also discussed that a surface area does not show any important effect on coexisting densities, making this method reliable to calculate bulk coexisting properties using small systems.     

Molecular simulation of the phase behavior of fluids and fluid mixtures using the synthetic method

A new molecular simulation procedure is reported for determining the phase behavior of fluids and fluid mixtures, which closely follows the experimental synthetic method. The simulation procedure can be implemented using Monte Calro or molecular dynamics in either the microcanonical or canonical statistical ensembles. Microcanonical molecular dynamics simulations are reported for the phase behavior of both the pure Lennard-Jones fluid and a Lennard-Jones mixture. The vapor pressures for the pure fluid are in good agreement with Monte Carlo Gibbs ensemble and Gibbs-Duhem calculations. The Lennard-Jones mixture is composed of equal size particles, with dissimilar energy parameters (?(2)∕?(1) = 1∕2, ?(12)∕?(1) = 1∕2). The binary Lennard-Jones mixture exhibits liquid-liquid equilibria at high pressures and the simulation procedure allows us to estimate the coordinates of the high-pressure branch of the critical curve.     

Study of Winsor I to Winsor II transitions in a lattice model

Experiments show that with increasing temperature, microemulsion systems undergo Winsor transitions. The transitions occur from Winsor I (oil droplets in water media) to Winsor II (water droplets in oil media) via Winsor III (bicontinuous phase) with an increase in the temperature. In this paper, it has been shown, for the first time, how one can study the qualitative effects of temperature, head, tail, and oil chain lengths, on these transitions. Simple cubic lattice with excluded volume and periodic boundary conditions is used to mimic the box of the simulation as a bulk of solution. The simulations have been done using the standard traditional Metropolis algorithm in the canonical ensemble (N, V, T). Configurational bias Monte Carlo and reptation moves are used with an equal probability to relax the systems. A very simple interaction model, i.e., the repulsions of water (or heads of surfactants) with oil (or tails of surfactants), is used due to the main characteristic of oil-water mixtures or amphiphilic molecule that is the hydrophobicity. The interfacial tension between oil and water (gammaow) is related to the averaged total energy of the lattice. The model shows that the Winsor III has a minimum interfacial tension (gammaow) similar to experimental results. Changing the phase structure from Winsor III to Winsor I (or Winsor II), increases the interfacial tension which is in agreement with experiments. To relate interfacial tension with the interaction parameter, the simple theory of Bragg-Williams has been used. All of the results such as the effects of oil chain length, head and tail beads number are all similar to the experimental results. Using the Davies method for calculating hydrophilic-lypophilic Balance (HLB), similar to the experimental results, Winsor III phase is formed at HLB value nearly to 10.