Computation of phase equilibria

The computation of phase equilibria using computer simulation can follow a number of different strategies. Here we will focus mainly on first-order transitions in fluid phases, usually liquid-vapour equilibria. Thermodynamic equilibrium implies, for two phases Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \alpha } and Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \beta } :

• equal temperatures: Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle T_{\alpha} = T_{\beta} }
• equal pressures: Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle p_{\alpha} = p_{\beta} }
• equal chemical potentials: Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \mu_{\alpha} = \mu_{\beta} }

Independent simulations for each phase at fixed temperature in the canonical ensemble

Simulations can be carried out using either the Monte Carlo or the molecular dynamics technique. Assuming that one has some knowledge on the phase diagram of the system, one can try the following recipe:

1. Fix a temperature and a number of particles
2. Perform a limited number of simulations in the low density region (where the gas phase density is expected to be)
3. Perform a limited number of simulations in the moderate to high density region (where the liquid phase should appear)
4. In these simulations we can compute for each density (at fixed temperature) the values of the pressure and the chemical potentials (for instance using the Widom test-particle method)

A quick 'first guess' method

Using the previously obtained results the following, somewhat unsophisticated, procedure can be used to obtain a first inspection of the possible phase equilibrium:

1. Fit the simulation results for each branch by using appropriate functional forms: Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. \mu_{v}(\rho) \right. ; p_v(\rho);\mu_l(\rho); p_l(\rho) }
2. Use the fits to build, for each phase, a table with three entries: Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \rho, p, \mu } , then plot for both tables Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \mu } as a function of Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle p } and check to see if the two lines intersect.
3. The crossing point provides (to within statistical uncertainty, the errors due to finite size effects, etc.) the coexistence conditions.

Improving the 'first guess' method

It can be useful to take into account classical thermodynamics to improve the previous analysis. This is because is is not unusual have large uncertainties in the results for the properties. The basic idea is to use thermodynamic consistency requirements to improve the analysis.

Methodology in the NpT ensemble

Low temperature: Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. T << T_c \right. }

For temperatures well below the critical point, provided that the calculation of the chemical potential of the liquid phase using Widom test-particle method gives precise results, the following strategy can be used to obtain a 'quick' result:

1. Perform an $NpT$ simulation of the liquid phase at zero pressure, i.e. Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle p \simeq 0 }
2. Arrive at an initial estimate, Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \mu^{(1)} } for the coexistence value of the chemical potential by computing, in the liquid phase: Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. \mu^{(1)} = \mu_l (N,T,p=0) \right. }
3. Make a first estimate of the coexistence pressure, Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle p^{(1)} } , by computing, either via simulation or via the virial coefficients of the gas phase, the pressure at which the gas phase fulfils: Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. \mu_g(N,T, p^{(1)} ) = \mu^{(1)} \right. }
4. Refine the results, if required, by performing a simulation of the liquid phase at Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. p^{(1)} \right. } , or use estimates of Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left( \partial V / \partial p \right)_{N,T,p=0} } (from the initial simulation) and the gas equation of state data to correct the initial estimates of pressure and chemical potential at coexistence.

Note that this method works only if the liquid phase remains metastable at zero pressure.

Weak first order transitions

There are situations where other strategies based on NpT ensemble simulations can be used. Similar approaches have also been applied in the Grand Canonical ensemble. In practice, the free energy barrier between the two phases has to be low enough to allow the simulated system to cross from one phase to the other when it is being simulated. Somehow, the situation is now the opposite to that described in Section 2.1. Taking into the account the classical partition function in the NpT ensemble it can be written:

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. P(V|N,p,T) \propto \exp \left[ - \frac{ p V}{k_B T} \right] \exp \left[ - \frac{A(N,V,T)}{k_B T} \right] \right. } ,

where Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. P(V|N,p,T) \right. } stands for the probability of the volume Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. V \right. } at given conditions Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. N,P,T \right. } . Let Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle p_{eq} } the pressure at which the phase transition occurs. In such a case the following scenario is expected for Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. P(V|N,p,T) \right. } :

• Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. P(V|N,p_{eq},T) \right. } has two maxima, corresponding to the liquid and vapor pure phases, with Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. P(V_v|N,p_{eq},T) = P(V_l|N,p_{eq},T) = P_{v/l} \right. }

• The probability of a given intermediate volume at Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. p_{eq} \right. } can be estimated (from macroscopic arguments) as:

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. P(V|N,p_{eq},T) \simeq P_{v/l} \times \exp \left[ - \frac{ \gamma(T) \mathcal A }{k_B T } \right] \right. } ,

where $\left.\gamma (T)\right.$ is the surface tension of the vapor-liquid interface, and Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. {\mathcal A} \right. } is the surface area, which depends on the thermodynamic variables Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. (N,V,T) \right. } and the geometry of the simulation box. For small values of the surface tension, small system sizes and good simulation algorithms it could be possible for pressures close to Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle p_{eq} } to sample in a simulation the whole region of densities between vapour and liquid densities. If such is the case the phase equilibria conditions can be computed by reweighing techniques applied on the volume histograms of the simulation.

Simple reweighing of the volume probability distribution

Suppose that a precise Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. NpT \right. } simulation has been carried out at pressure $\left.p_{0}\right.$. From that simulation a volume probability distribution, Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. P_0(V|N,p_0,T) \right. } has been computed. It is possible to use this function to estimate the distributions of volume for pressure values close to Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. p_0 \right. } ;

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. P(V|N,p,T) \propto \exp \left[ - \frac{ (p-p_0) V }{k_B T } \right] P_0(V|N,p_0,T) \right. }

Using this reweighing procedure we can estimate the value of Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. p_{eq} \right. } ; that is: the value of pressure for which the volume distribution presents two peaks of equal height. The analysis of the form of the distributions at the equilibrium conditions for different system sizes (i.e. in the current case for different values of $\left.N\right.$) can be useful to distinguish between continuous and first order phase transitions.

See also

Surface tension

Van der Waals loops in the canonical ensemble

It is possible to compute the liquid-vapour equilibrium without explicit calculations of the chemical potential (or the pressure) by performing a number of simulations sampling appropriately the vapour, liquid, and intermediate regions. As an example, consider a simple fluid at a given subcritical temperature (Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. T < T_c \right. } ). We can perform a number of simulations for a given number of particles, Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. N \right. } and different densities:

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \rho_i = i \times \Delta \rho; \; \; \; \; i=1, 2, \cdots, m }

In these simulations, we can compute the pressure (or the chemical potential) and fit the result to an appropriate equation. With such an equation of state the phase equilibria can be estimated. If two phase equilibria exists, a loop in the representation of Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. p = p (\rho) \right. } (or $\left.\mu =\mu (\rho )\right.$) should appear.

• Computing Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. \mu(\rho) \right. } from the equation of state given as Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. p(\rho) \right. } :

For fixed temperature and number of particles:

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. d A = - p d V \right. }

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. d A = \frac{ p N }{\rho^2} d \rho ; \; \; \; \; d a = \frac{ p }{\rho^2} d \rho \right. }

where $\left.A\right.$ is the Helmholtz energy function and Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. a = A/N \right. } , integrating:

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. a(\rho) = a(\rho_0) + \int_{\rho_0}^{\rho} \frac{ p(\rho') }{(\rho')^2} d \rho' \right. }

On the other hand

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. A = - p V + \mu N ; \; \; \; a = - p/ \rho + \mu \right. }

therefore:

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. \mu(\rho) - p(\rho)/\rho = \mu(\rho_0) - p(\rho_0)/\rho_0 + \int_{\rho_0}^{\rho} \frac{ p(\rho') }{(\rho')^2} d \rho' \right. }

$\left.\mu (\rho )=\mu (\rho _{0})+p(\rho )/\rho -p(\rho _{0})/\rho _{0}+\int _{\rho _{0}}^{\rho }{\frac {p(\rho ')}{(\rho ')^{2}}}d\rho '\right.$

A similar procedure can be built up to compute Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. p(\rho) \right. } from Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. \mu(\rho) \right. } . Once Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. p(\rho) \right. } and Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \left. \mu(\rho) \right. } are known it is straightforward to compute the coexistence point.

Practical details

Some precautions should be taken if this procedure is used:

• The precision of the simulation results in the two phase region will be poor (so, large simulations are required to have a good estimation of the equation of state)
• The simulation results in the two phase region will depend dramatically on the system size (calculations with different number of particles become convenient to check the quality of the phase equilibria results)

Direct simulation of the two phase system

Direct simulation of the two phase system was first implemented by Abraham . It has since been applied to the Lennard-Jones model and water .

Gibbs ensemble Monte Carlo for one component systems

The Gibbs ensemble Monte Carlo method is often considered as a smart variation of the standard canonical ensemble procedure (See ). The simulation is, therefore, carried out at constant volume, temperature and number of particles. The whole system is divided into two non-interacting parts, each one has its own simulation box with its own periodic boundary conditions. This separation of the two phases into different boxes is in order to suppress any influence due to interfacial effects. The two subsystems can interchange volume and particles. The rules for these interchanges are built up so as to guarantee conditions of both chemical and mechanical equilibrium between the two phases. If the overall conditions are of phase separation, it is expected that two phases will appear in different simulation boxes.

External links

• Gibbs ensemble Monte Carlo code on the Panagiotopoulos Group Homepage

Mixtures

Symmetric mixtures

Examples of symmetric mixtures can be found both in lattice of continuous model. The Ising model can be viewed as mixture of two different chemical species which de-mix at low temperature. The symmetry in the interactions can be exploited to simplify the calculation of phase diagrams.

See also

• Constrained cell method
• Gibbs-Duhem integration
• Phase switch Monte Carlo