Pressure derivatives in the classical molecular-dynamics ensemble

Classical statistical mechanics in the μVL and μpR ensembles

Molecular expressions for thermodynamic properties and derivatives of the entropy up to third order in the adiabatic grand-isochoric μVL and adiabatic grand-isobaric μpR ensembles are systematically derived using the methodology developed by Lustig for the microcanonical and canonical ensembles [J. Chem. Phys. 100, 3048 (1994)10.1063/1.466446; Mol. Phys. 110, 3041 (2012)10.1080/00268976.2012.695032]. They are expressed by phase-space functions, which represent derivatives of the entropy with respect to the chemical potential, the volume, and the Hill energy L in the μVL ensemble and with respect to the chemical potential, the pressure, and the Ray energy R in the μpR ensemble. The derived expressions are validated for both ensembles by Monte Carlo simulations for the simple Lennard-Jones model fluid at three selected state points.

Rigorous expressions for thermodynamic properties in the NpH ensemble

Molecular expressions for thermodynamic properties of fluids and derivatives of the entropy up to third order in the isoenthalpic-isobaric ensemble are derived by using the methodology developed by Lustig for the microcanonical and canonical ensembles [J. Chem. Phys. 100, 3048 (1994)JCPSA60021-960610.1063/1.466446; Mol. Phys. 110, 3041 (2012)MOPHAM0026-897610.1080/00268976.2012.695032]. They are expressed in a systematic way by phase-space functions, which represent derivatives of the phase-space volume with respect to enthalpy and pressure. The expressions for thermodynamic properties contain only ensemble averages of combinations of the kinetic energy and volume of the system. Thus, the calculation of thermodynamic properties in the isoenthalpic-isobaric ensemble does not require volume derivatives of the potential energy. This is particularly advantageous in Monte Carlo simulations when the interactions between molecules are described by very accurate ab initio pair and nonadditive three-body potentials. The derived expressions are validated by Monte Carlo simulations for the simple Lennard-Jones model fluid as a test case.

Classical statistical mechanics in the grand canonical ensemble

The methodology developed by Lustig for calculating thermodynamic properties in the microcanonical and canonical ensembles [J. Chem. Phys. 100, 3048 (1994)JCPSA60021-960610.1063/1.466446; Mol. Phys. 110, 3041 (2012)MOPHAM0026-897610.1080/00268976.2012.695032] is applied to derive rigorous expressions for thermodynamic properties of fluids in the grand canonical ensemble. All properties are expressed by phase-space functions, which are related to derivatives of the grand canonical potential with respect to the three independent variables of the ensemble: temperature, volume, and chemical potential. The phase-space functions contain ensemble averages of combinations of the number of particles, potential energy, and derivatives of the potential energy with respect to volume. In addition, expressions for the phase-space functions for temperature-dependent potentials are provided, which are required to account for quantum corrections semiclassically in classical simulations. Using the Lennard-Jones model fluid as a test case, the derived expressions are validated by Monte Carlo simulations. In contrast to expressions for the thermal expansion coefficient, the isothermal compressibility, and the thermal pressure coefficient from the literature, our expressions yield more reliable results for these properties, which agree well with a recent accurate equation of state for the Lennard-Jones model fluid. Moreover, they become equivalent to the corresponding expressions in the canonical ensemble in the thermodynamic limit.

Liquid dibromomethane under pressure: a computational study

Molecular dynamics simulations have been performed on liquid dibromomethane at thermodynamic states corresponding to temperature in the range 268-328 K and pressure varying from 1 bar to 3000 bar. The interaction model is a simple effective two-body pair potential with atom-atom Coulomb and Lennard-Jones interactions and molecules are rigid. Thermodynamic properties have been studied, including the isobaric thermal expansion coefficient, the isothermal compressibility, the heat capacities and the speed of sound. The simulation results exhibit a crossing of the isotherms of the isobaric thermal expansion coefficient at about 800 bar in very good agreement with the prediction of an isothermal fluctuation equation of state predicting such a crossing in the pressure range 650-900 bar, though experimental results up to 1000 bar do not find any crossing.

Effect of the range of particle cohesion on the phase behavior and thermodynamic properties of fluids

Molecular simulations are performed for the (m + 1, m) potential to systematically investigate the effect of changing the range of particle cohesion on both vapor-liquid equilibria and thermodynamic properties of fluids. The results are reported for m = 4-11, which represent a progressive narrowing of the potential energy well. The conventional Lennard-Jones potential is used as a reference point for normal fluid behavior. Small values of m result in a broadening of the phase envelope compared with the Lennard-Jones potential, whereas a contraction is observed in other cases. The critical properties are reported, and a relationship between the critical temperature and the Boyle temperature is determined. The low values of the critical compressibility factor when m < 6 reflect the behavior observed for real fluids such as n-alkanes. The results for supercritical thermodynamic properties are much more varied. Properties such as pressure, potential energy, isochoric thermal pressure coefficient, and thermal expansion coefficient vary consistently with m, whereas other properties such as the Joule-Thomson coefficient exhibit much more nuanced behavior. Maximum and minimum values are reported for both the isochoric heat capacity and isothermal compressibility. A minimum in the speed of sound is also observed.

Thermodynamic properties and anomalous behavior of double-Gaussian core model potential fluids

The structural, thermodynamic, and vapor-liquid equilibria properties of the double-Gaussian core model (DGM) potential are studied via molecular simulation. Results are presented for the pressure (p), potential energy (U), isochoric and isobaric heat capacities (C_{V,p}), isothermal compressibility (β_{T}), isochoric thermal pressure coefficient (γ_{V}), thermal expansion coefficient (α_{p}), speed of sound (ω_{0}), and the Joule-Thomson coefficient (μ_{JT}), which are compared with simulations for the Gaussian core model (GCM) potential. A feature of the simulations is that both the GCM and DGM potentials reproduce many of the anomalous properties of water, such as a maximum density, γ_{V}<0, maximum values for both α_{p} and β_{T}, and minimum values in both C_{p} and ω_{0}. The presence of attractive interaction enhances the anomalies and also yields some additional features such as a more structured vapor phase and Joule-Thomson inversion.

Molecular simulation of orthobaric isochoric heat capacities near the critical point

A molecular simulation strategy is investigated for detecting the divergence of the isochoric heat capacity (C_{V}) on the vapor and liquid coexistence branches of a fluid near the critical point. The procedure is applied to the empirical Lennard-Jones potential and accurate state-of-the-art ab initio two-body and two-body + three-body potentials for argon. Simulations with the Lennard-Jones potential predict the divergence of C_{V}, and the phenomenon is also observed for both two-body and two-body + three body potentials. The potentials also correctly predict the crossover between vapor and liquid C_{V} values and the subcritical liquid C_{V} minimum, which marks the commencement C_{V} divergence. The effect of three-body interactions is to delay the onset of divergence to higher subcritical temperatures.

Molecular hydrodynamics: Vortex formation and sound wave propagation

In the present study, quantitative feasibility tests of the hydrodynamic description of a two-dimensional fluid at the molecular level are performed, both with respect to length and time scales. Using high-resolution fluid velocity data obtained from extensive molecular dynamics simulations, we computed the transverse and longitudinal components of the velocity field by the Helmholtz decomposition and compared them with those obtained from the linearized Navier-Stokes (LNS) equations with time-dependent transport coefficients. By investigating the vortex dynamics and the sound wave propagation in terms of these field components, we confirm the validity of the LNS description for times comparable to or larger than several mean collision times. The LNS description still reproduces the transverse velocity field accurately at smaller times, but it fails to predict characteristic patterns of molecular origin visible in the longitudinal velocity field. Based on these observations, we validate the main assumptions of the mode-coupling approach. The assumption that the velocity autocorrelation function can be expressed in terms of the fluid velocity field and the tagged particle distribution is found to be remarkably accurate even for times comparable to or smaller than the mean collision time. This suggests that the hydrodynamic-mode description remains valid down to the molecular scale.

Galilean-invariant Nosé-Hoover-type thermostats

A new pairwise Nosé-Hoover type thermostat for molecular dynamics (MD) simulations which is similar in construction to the pair-velocity thermostat of Allen and Schmid, [Mol. Simul. 33, 21 (2007)] (AS) but is based on the configurational thermostat is proposed and tested. Both thermostats generate the canonical velocity distribution, are Galilean invariant, and conserve linear and angular momentum. The unique feature of the pairwise thermostats is an unconditional conservation of the total angular momentum, which is important for thermalizing isolated systems and those nonequilibrium bulk systems manifesting local rotating currents. These thermostats were benchmarked against the corresponding Nosé-Hoover (NH) and Braga-Travis prescriptions, being based on the kinetic and configurational definitions of temperature, respectively. Some differences between the shear-rate-dependent shear viscosity from Sllod nonequilibrium MD are observed at high shear rates using the different thermostats. The thermostats based on the configurational temperature produced very similar monotically decaying shear viscosity (shear thinning) with increasing shear rate, while the NH method showed discontinuous shear thinning into a string phase, and the AS method produced a continuous increase of viscosity (shear thickening), after a shear thinning region at lower shear rates. Both pairwise additive thermostats are neither purely kinetic nor configurational in definition, and possible directions for further improvement in certain aspects are discussed.

Thermodynamic properties and entropy scaling law for diffusivity in soft spheres

The purely repulsive soft-sphere system, where the interaction potential is inversely proportional to the pair separation raised to the power n, is considered. The Laplace transform technique is used to derive its thermodynamic properties in terms of the potential energy and its density derivative obtained from molecular dynamics simulations. The derived expressions provide an analytic framework with which to explore soft-sphere thermodynamics across the whole softness-density fluid domain. The trends in the isochoric and isobaric heat capacity, thermal expansion coefficient, isothermal and adiabatic bulk moduli, Grüneisen parameter, isothermal pressure, and the Joule-Thomson coefficient as a function of fluid density and potential softness are described using these formulas supplemented by the simulation-derived equation of state. At low densities a minimum in the isobaric heat capacity with density is found, which is a new feature for a purely repulsive pair interaction. The hard-sphere and n = 3 limits are obtained, and the low density limit specified analytically for any n is discussed. The softness dependence of calculated quantities indicates freezing criteria based on features of the radial distribution function or derived functions of it are not expected to be universal. A new and accurate formula linking the self-diffusion coefficient to the excess entropy for the entire fluid softness-density domain is proposed, which incorporates the kinetic theory solution for the low density limit and an entropy-dependent function in an exponential form. The thermodynamic properties (or their derivatives), structural quantities, and diffusion coefficient indicate that three regions specified by a convex, concave, and intermediate density dependence can be expected as a function of n, with a narrow transition region within the range 5 < n < 8.

Riemannian geometry study of vapor-liquid phase equilibria and supercritical behavior of the Lennard-Jones fluid

The behavior of thermodynamic response functions and the thermodynamic scalar curvature in the supercritical region have been studied for a Lennard-Jones fluid based on a revised modified Benedict-Webb-Rubin equation of state. Response function extrema are sometimes used to estimate the Widom line, which is characterized by the maxima of the correlation lengths. We calculated the Widom line for the Lennard-Jones fluid without using any response function extrema. Since the volume of the correlation length is proportional to the Riemannian thermodynamic scalar curvature, the locus of the Widom line follows the slope of maximum curvature. We show that the slope of the Widom line follows the slope of the isobaric heat capacity maximum only in the close vicinity of the critical point and that, therefore, the use of response function extrema in this context is problematic. Furthermore, we constructed the vapor-liquid coexistence line for the Lennard-Jones fluid using the fact that the correlation length, and therefore the thermodynamic scalar curvature, must be equal in the two coexisting phases. We compared the resulting phase envelope with those from simulation data where multiple histogram reweighting was used and found striking agreement between the two methods.

Periodic boundary condition induced breakdown of the equipartition principle and other kinetic effects of finite sample size in classical hard-sphere molecular dynamics simulation

We examine consequences of the non-Boltzmann nature of probability distributions for one-particle kinetic energy, momentum, and velocity for finite systems of classical hard spheres with constant total energy and nonidentical masses. By comparing two cases, reflecting walls (NVE or microcanonical ensemble) and periodic boundaries (NVEPG or molecular dynamics ensemble), we describe three consequences of the center-of-mass constraint in periodic boundary conditions: the equipartition theorem no longer holds for unequal masses, the ratio of the average relative velocity to the average velocity is increased by a factor of [N/(N-1)]1/2, and the ratio of average collision energy to average kinetic energy is increased by a factor of N/(N-1). Simulations in one, two, and three dimensions confirm the analytic results for arbitrary dimension.