This series of papers is devoted to identifying and explaining the properties of strongly correlating liquids, i.e., liquids with more than 90% correlation between their virial W and potential energy U fluctuations in the NVT ensemble. Paper IV [N. Gnan et al, J. Chem. Phys.131, 234504 (2009) https://doi.org/10.1063/1.3265957] showed that strongly correlating liquids have “isomorphs,” which are curves in the phase diagram along which structure, dynamics, and some thermodynamic properties are invariant in reduced units. In the present paper, using the fact that reduced-unit radial distribution functions are isomorph invariant, we derive an expression for the shapes of isomorphs in the WU phase diagram of generalized Lennard-Jones systems of one or more types of particles. The isomorph shape depends only on the Lennard-Jones exponents; thus all isomorphs of standard Lennard-Jones systems (with exponents 12 and 6) can be scaled onto a single curve. Two applications are given. One tests the prediction that the solid-liquid coexistence curve follows an isomorph by comparing to recent simulations by Ahmed and Sadus [J. Chem. Phys.131, 174504 (2009)] https://doi.org/10.1063/1.3253686. Excellent agreement is found on the liquid side of the coexistence curve, whereas the agreement is less convincing on the solid side. A second application is the derivation of an approximate equation of state for generalized Lennard-Jones systems by combining the isomorph theory with the Rosenfeld-Tarazona expression for the temperature dependence of the potential energy on isochores. It is shown that the new equation of state agrees well with simulations.

This is the final paper in a series1–4 investigating the properties of strongly correlating liquids,5 i.e., liquids that have strong correlations between their constant-volume equilibrium fluctuations of potential energy, U(t), and virial6,7

$W(t)\equiv -1/3 \sum _i {\bf r}_i \cdot {\bf \nabla }_{{\bf r}_i} U({\bf r}_1, ..., {\bf r}_N)$
W(t)1/3iri·riU(r1,...,rN) where ri is the position of particle i at time t. As is well known, the average virial W gives the configurational contribution to the pressure:

\begin{equation}pV \,=\, Nk_BT+W\,.\end{equation}
pV=NkBT+W.
(1)

Letting Δ denote instantaneous deviations from equilibrium mean values, the WU correlation is quantified by the correlation coefficient R (with

$\left< ... \right>$
... denoting equilibrium average):

\begin{equation}R \,=\, \frac{\langle \Delta W\Delta U\rangle }{\sqrt{\langle (\Delta W)^2\rangle \langle (\Delta U)^2\rangle }}\,.\end{equation}
R=ΔWΔU(ΔW)2(ΔU)2.
(2)

Perfect correlation gives R = 1. As a pragmatic definition we have chosen “strongly correlating liquids” to designate liquids that have R ⩾ 0.9 in the NVT ensemble (constant volume, V, and temperature, T).

Strongly correlating liquids have simpler physics than liquids in general, an observation that has particular significance for the highly viscous phase.8–18 Thus it has been shown that strongly correlating viscous liquids to a good approximation have all eight frequency-dependent thermoviscoelastic response functions19–21 given in terms of just one22 (i.e., are single-parameter liquids in the sense of having dynamic Prigogine-Defay ratio19 close to unity2,20,22). Strongly correlating viscous liquids moreover obey density scaling23–27 to a good approximation, i.e., their dimensionless relaxation time

$\tilde{\tau }\equiv \tau \rho ^{1/3}\sqrt{k_BT/m}$
τ̃τρ1/3kBT/m (where m is the average particle mass) depends on density ρ = N/V and temperature as
$\tilde{\tau }= F(\rho ^\gamma /T)$
τ̃=F(ργ/T)
.28–30 

Paper I1 presented computer simulations of 13 different systems, showing that van der Waals type liquids are strongly correlating, whereas hydrogen-bonding liquids like methanol or water are not. Strongly correlating liquids include,1,2,5,22,29 for instance, the standard Lennard-Jones (LJ) liquid, the Kob-Andersen binary LJ (KABLJ) mixture, an asymmetric rigid-bond dumbbell model, a seven-site united-atom toluene model, and the Lewis-Wahnström OTP model.

Paper II2 analyzed the cause of WU correlations with a focus on the LJ potential. The strong correlations were related to the well-known fact that an inverse power-law (IPL) pair potential, v(r)∝rn where r is the distance between two particles,31–42 implies perfect WU correlation,2,5

\begin{equation}\Delta W(t)\, =\, \gamma \Delta U(t)\,\end{equation}
ΔW(t)=γΔU(t)
(3)

with γ = n/3. Around the potential energy minimum, the LJ potential is well described by an “extended” inverse power-law potential (eIPL),2 vLJ(r)≅Arn + B + Cr. At constant volume the linear term contributes little to the virial and potential-energy fluctuations: when one nearest-neighbor distance increases, others decrease in such a way that the sum is almost constant. Thus systems interacting via the LJ potential inherit strong WU correlations from an underlying inverse power-law—they have a “hidden scale invariance.”3,28

Paper III3 gave further numerical evidence for the explanation for strong WU correlations presented in Paper II, and theoretical results were given on the statistical mechanics and thermodynamics of the hidden scale invariance that characterizes strongly correlating liquids. It was also shown that strong virial potential-energy correlations are present even in out-of-equilibrium situations—the hidden scale invariance is a property of the potential energy surface, not just of the equilibrium states.

Paper IV4 introduced the concept of “isomorphs” in the phase diagram of a strongly correlating liquid. Starting from a single assumption a number of isomorph invariants were derived. In particular, structure and dynamics were shown to be invariant on isomorphs when reduced units are used.

In the present paper further simulation results supporting the isomorph predictions are presented for systems interacting via a multicomponent generalized LJ potential:

\begin{equation}v_{ij}(r_{ij})=v_{ij}^{(m)}(r_{ij})+v_{ij}^{(n)}(r_{ij}),\,\,\end{equation}
vij(rij)=vij(m)(rij)+vij(n)(rij),
(4)

where

$v_{ij}^{(k)}(r_{ij})$
vij(k)(rij) is an IPL potential involving the two particles i and j:

\begin{equation}v_{ij}^{(k)}(r_{ij})\equiv \varepsilon ^{(k)}_{ij}\big({\sigma ^{(k)}_{ij}}/{r_{ij}}\big)^k.\end{equation}
vij(k)(rij)ɛij(k)σij(k)/rijk.
(5)

Section II briefly reviews the isomorph concept. For systems interacting via a generalized LJ potential, a prediction for the shape of the isomorphs in the WU phase diagram is derived in Sec. III and demonstrated to fit well to simulation results.43 Interestingly, the isomorph shape depends only on the exponents m and n. Thus, e.g., all 12-6 LJ systems have the same isomorphs in the WU phase diagram. In Sec. IV we present two applications of the theory. One tests the Paper IV prediction that solid-liquid coexistence curves are isomorphs. The second application gives an approximate equation of state for systems interacting via generalized LJ potentials; this is arrived at by combining the present theory with Rosenfeld and Tarazona's expression for the isochoric temperature dependence of the potential energy.

We term a microscopic configuration “physically relevant” if its influence on the thermodynamics and dynamics of the system is not a priori negligible (Paper IV). For instance, any configuration with very strong particle overlap is physically irrelevant. Note, however, that some unlikely configurations like transition states are physically relevant.

Two state points (1) and (2) with temperatures T1 and T2 and densities ρ1 and ρ2, respectively, are defined to be isomorphic (Paper IV) if they obey the following: any two physically relevant configurations of state points (1) and (2),

$\big({\bf r}_1^{(1)},\ldots , {\bf r}_N^{(1)}\big)$
r1(1),...,rN(1) and
$\big({\bf r}_1^{(2)},\ldots , {\bf r}_N^{(2)}\big)$
r1(2),...,rN(2)
, which trivially scale into one another,

\begin{equation}\rho _1^{1/3}\mathbf {r}_i^{(1)} \,=\,\rho _2^{1/3}\mathbf {r}_i^{(2)}\,\, (i=1,...,N)\,,\end{equation}
ρ11/3ri(1)=ρ21/3ri(2)(i=1,...,N),
(6)

have proportional configurational Boltzmann statistical weights:

\begin{equation}e^{-U({\bf r}_1^{(1)},\ldots , {\bf r}_N^{(1)})/k_BT_1}\, =\, C_{12}e^{-U({\bf r}_1^{(2)},\ldots , {\bf r}_N^{(2)})/k_BT_2}\,.\end{equation}
eU(r1(1),...,rN(1))/kBT1=C12eU(r1(2),...,rN(2))/kBT2.
(7)

Here U(r1, …, rN) is the potential energy function and it is understood that the constant C12 depends only on the state points (1) and (2).

The property of being isomorphic defines a mathematical equivalence relation on the set of state points. The corresponding equivalence classes are smooth curves in the phase diagram termed isomorphs or isomorphic curves.

Equation (7) implies identity of the normalized Boltzmann probabilities for scaled physically relevant configurations of isomorphic state points. As detailed in Paper IV this identity implies that several quantities are invariant along an isomorph. Examples include the configurational entropy, the isochoric specific heat, all multiparticle distribution functions in reduced units (in particular, the radial distribution function(s)), reduced-unit dynamics (both for Newtonian and stochastic dynamics), normalized reduced-time autocorrelation functions, reduced-unit transport coefficients, etc. It was further shown in Paper IV that a jump from equilibrium at one state point to an isomorphic state point brings the system instantaneously to equilibrium (see also Ref. 44).

IPL potentials are Euler homogeneous functions, i.e., obey Ur1, …, λrN) = λnU(r1, …, rN). Such purely repulsive potentials – which do not describe real intermolecular interactions – are the only potentials that have 100% correlation between virial and potential-energy fluctuations; this reflects the well-known IPL virial-theorem identity W = (n/3)U. Likewise, only IPL liquids have exact isomorphs; they satisfy Eq. (7) with C12 = 1 for all configurations of two state points obeying

$\rho _1^\gamma /T_1=\rho _2^\gamma /T_2$
ρ1γ/T1=ρ2γ/T2 where γ = n/3.

Appendix A of Paper IV showed that strongly correlating liquids generally have isomorphs to a good approximation. Consequently, to a good approximation these liquids inherit a number of the exact invariants that IPL liquids have along their curves given by ργ/T = Const., namely, all IPL invariants that do not rely on the IPL identity C12 = 1. Examples include the above mentioned thermodynamic, static, and dynamic isomorph invariants. IPL invariants that are not inherited by strongly correlating liquids in general include, e.g., the Helmholtz free energy over temperature, the potential energy over temperature, and the reduced-unit compressibility.

The fact that strongly correlating liquids have isomorphs to a good approximation should be understood as follows. Any liquid has curves in its phase diagram of constant excess entropy, curves of constant isochoric specific heat, curves of constant reduced relaxation time, curves of constant reduced-unit transport coefficients, etc. For a strongly correlating liquid these curves are almost identical, and they define its (approximate) isomorphs. Points on these isomorphs have approximately identical structure and approximately identical dynamics as probed, e.g., by normalized reduced-time-autocorrelation functions.

For clarity of presentation below we shall not repeatedly mention that the existence of isomorphs is not exact for generalized LJ systems. Accordingly, we shall speak of isomorphs as a unique, well-defined concept. In other words: all properties of isomorphs derived below in generalized LJ systems would be exact if these systems did have exact isomorphs, e.g., fulfilled Eq. (7); since this is not the case, the properties are in principle – as well as practice – approximate.

The structure in reduced units (

$\tilde{r} \equiv \rho ^{1/3}r$
r̃ρ1/3r⁠) is predicted to be invariant along an isomorph (Paper IV). The excess entropy, Sex = SSideal, depends only on structure and was shown in Paper IV also to be invariant along an isomorph (“excess” refers to the quantity in question in excess to that of an ideal gas with same density and temperature). In the below simulations we generate an isomorph as a set of state points in the phase diagram with constant Sex. To change density and temperature keeping Sex constant, the following identity is used (Paper IV):

\begin{eqnarray}\gamma &\equiv & -\left( \frac{\partial \ln T}{\partial \ln V}\right)_{\rm S_{ex}} = \frac{V \left( \frac{\partial S_{\rm ex}}{\partial V}\right)_T}{T \left( \frac{\partial S_{\rm ex}}{\partial T}\right)_V} = \frac{V\beta _{\rm v}^{\rm ex}}{C_{\rm v}^{\rm ex}} \,,\end{eqnarray}
γlnTlnVS ex =VS ex VTTS ex TV=Vβv ex Cv ex ,
(8)

where

$\beta _{\rm v}^{\rm ex} \equiv ( {\partial p_{\rm ex}}/{\partial T}\big)_V = \frac{1}{V}( {\partial W}/{\partial T})_V$
βv ex (p ex /T)V=1V(W/T)V and a Maxwell relation was applied. Utilizing the standard fluctuation formulae this leads (Paper IV) to:

\begin{eqnarray}\gamma &=& \left( \frac{\partial \ln T}{\partial \ln \rho }\right)_{\rm S_{ex}} = \frac{\left\langle \Delta U \Delta W \right\rangle }{\left\langle (\Delta U)^2\right\rangle }.\end{eqnarray}
γ=lnTlnρS ex =ΔUΔW(ΔU)2.
(9)

The procedure applied to generate isomorphs can be summarized as: (1) an equilibrium NVT simulation was performed at a given state point; (2) γ was calculated from the fluctuations using Eq. (9); (3) the density ρ was changed slightly (of order 1%), and γ was used to calculate the corresponding change in temperature in order to keep

$\rm S_{ex}$
S ex constant (Eq. (9)); 4) a simulation at the new state point was performed, and the procedure was repeated.

Figure 1(a) demonstrates the isomorph invariance of the radial distribution function for the large (A) particles of the KABLJ mixture.45 For comparison, Fig. 1(b) shows the same on an isotherm with a similar (but smaller) density change. Clearly, structure is to a good approximation invariant on the isomorph, whereas this is not the case on the isotherm.

Taking the logarithm of Eq. (7) shows that the potential energy surface in reduced units (

$\tilde{U} \equiv U/k_BT$
ŨU/kBT⁠) is the same for two isomorphic state points, except for an additive constant and scaling of the coordinates. This constant does not influence the dynamics, and consequently the dynamics in reduced units is the same for two isomorphic state points. Figure 2(a) demonstrates the isomorph invariance of the incoherent intermediate scattering function of the A particles of the KABLJ mixture. For comparison, Fig. 2(b) shows the same property on the T = 0.55 isotherm. Like the structure, the dynamics is also invariant to a good approximation on the isomorph—both regarding average relaxation time and the shape of the relaxation function; this is far from the case on the corresponding isotherm.

The exponent γ calculated from Eq. (9) generally depends on the state point. Figure 3 shows (as red squares) γ along the isomorph of Figs. 1(a) and 2(a). Temperature changes by almost a factor of two, but γ changes less than 5%. Paper IV argued that for a system with good isomorphs γ depends to a good approximation only on density. This is supported by Fig. 3 in which γ on an isochore is shown to vary only ∼1% when temperature is changed by a factor of two (black circles). Note, however, that γ eventually approaches m/3 = 4, but slowly so: at T = 10 we find γ = 4.6.

FIG. 3.

γ calculated from equilibrium fluctuations via Eq. (9) for the KABLJ system. Red squares: state points belonging to the isomorph presented in Figs. 1(a) and 2(a). Black circles: isochoric state points. Straight lines are linear regression fits. The standard errors on the slopes are 0.06 (isomorph) and 0.04 (isochore).

FIG. 3.

γ calculated from equilibrium fluctuations via Eq. (9) for the KABLJ system. Red squares: state points belonging to the isomorph presented in Figs. 1(a) and 2(a). Black circles: isochoric state points. Straight lines are linear regression fits. The standard errors on the slopes are 0.06 (isomorph) and 0.04 (isochore).

Close modal

It follows from Eq. (8) that the slope of isochores in the WU phase diagram is given by γ:

\begin{eqnarray}\gamma &=& \left( \frac{\partial W}{\partial U}\right)_{\rm V}.\end{eqnarray}
γ=WUV.
(10)

The insignificant change of γ along isochores (black circles in Fig. 3) means that these to a good approximation are straight lines in the WU phase diagram. This is illustrated in Fig. 4 for four different binary generalized LJ mixtures.

FIG. 4.

Isochores for four different binary generalized LJ mixtures. 12-6 KABLJ is the Kob-Andersen binary 80:20 mixture45 using the standard 12-6 LJ potential. 18-6 KABLJ is the same system with the generalized LJ potential (m = 18, n = 6, Eq. (4)) and parameters chosen to place the minimum of the potential at the same distance and depth as in the 12-6 KABLJ system. Likewise, 18-9 KABLJ is defined by using the exponents m = 18 and n = 9. Finally 12-6 WBLJ is the Wahnström 50:50 mixture46 using the standard 12-6 LJ potential. Straight lines are linear regression fits. The inset shows a zoom on the data for the two 12-6 LJ systems. For each system the simulated temperatures range from where non-exponential relaxation sets in (τα ≈ 100) to a viscous state with two-step relaxation (τα ≈ 103 − 105).

FIG. 4.

Isochores for four different binary generalized LJ mixtures. 12-6 KABLJ is the Kob-Andersen binary 80:20 mixture45 using the standard 12-6 LJ potential. 18-6 KABLJ is the same system with the generalized LJ potential (m = 18, n = 6, Eq. (4)) and parameters chosen to place the minimum of the potential at the same distance and depth as in the 12-6 KABLJ system. Likewise, 18-9 KABLJ is defined by using the exponents m = 18 and n = 9. Finally 12-6 WBLJ is the Wahnström 50:50 mixture46 using the standard 12-6 LJ potential. Straight lines are linear regression fits. The inset shows a zoom on the data for the two 12-6 LJ systems. For each system the simulated temperatures range from where non-exponential relaxation sets in (τα ≈ 100) to a viscous state with two-step relaxation (τα ≈ 103 − 105).

Close modal

What is the shape of an isomorph in the WU phase diagram? This section answers this question in generality for the multi-component generalized LJ potential. Recall that this potential (Eq. (4)) is the sum of two IPL’s. Correspondingly, the potential energy and virial can be expressed as sums of two IPL terms:

\begin{eqnarray}U &=& U_m + U_n, \qquad U_k \equiv \left< \sum _{i>j}v_{ij}^{(k)}(r_{ij})\right>.\end{eqnarray}
U=Um+Un,Uki>jvij(k)(rij).
(11)

For pair interactions the virial is given by6 

\begin{eqnarray}W \equiv -\frac{1}{3}\left< \sum _{i>j}r_{ij}v_{ij}^{^{\prime }}(r_{ij})\right>,\end{eqnarray}
W13i>jrijvij(rij),
(12)

where the prime denotes the derivative with respect to rij. From this we get

\begin{eqnarray}W &=& \frac{m}{3}U_m + \frac{n}{3}U_n .\end{eqnarray}
W=m3Um+n3Un.
(13)

For any point in the WU phase diagram we can solve Eqs. (11) and (13) for (Um, Un):

\begin{eqnarray}U_m &=& \frac{ 3W-nU}{m-n} ,\end{eqnarray}
Um=3WnUmn,
(14)
\begin{eqnarray}U_n &=& \frac{-3W+mU}{m-n} .\end{eqnarray}
Un=3W+mUmn.
(15)

These equations allow one to determine which regions of the WU phase diagram are accessible to the system. Adopting the convention

$m>n$
m>n⁠, we restrict ourselves to the case
$\epsilon _{ij}^{(m)}>0$
εij(m)>0
in order to have a repulsive core for all interactions. This implies Um > 0 for all configurations, and therefore also for all thermodynamic averages. Combining this with Eq. (14) implies
$W>\frac{n}{3}U$
W>n3U
. If all interactions have an attractive part
$(\epsilon _{ij}^{(n)}<0)$
(εij(n)<0)
we have Un < 0, which in combination with Eq. (15) implies
$W>\frac{m}{3}U$
W>m3U
. As an example, the region of the WU phase diagram accessible to the standard 12-6 LJ potential is given by W > 4U and W > 2U (the later inequality being relevant when U ⩽ 0). See Fig. 5 for an illustration.

FIG. 5.

Schematic drawing of the WU phase diagram for systems interacting via the standard 12-6 LJ potential. Included in this figure are also two isomorphs, plotted using Eq. (28). The accessible points are those above both dashed lines (see text). The part of the isomorphs between the dashed lines in the W, U > 0 quadrant are thus not accessible to standard 12-6 LJ systems (they are accessible to the purely repulsive 12-6 LJ potential with

$\epsilon _{ij}^{(6)}>0$
εij(6)>0⁠).

FIG. 5.

Schematic drawing of the WU phase diagram for systems interacting via the standard 12-6 LJ potential. Included in this figure are also two isomorphs, plotted using Eq. (28). The accessible points are those above both dashed lines (see text). The part of the isomorphs between the dashed lines in the W, U > 0 quadrant are thus not accessible to standard 12-6 LJ systems (they are accessible to the purely repulsive 12-6 LJ potential with

$\epsilon _{ij}^{(6)}>0$
εij(6)>0⁠).

Close modal

Along an isomorph the structure is invariant in reduced units (Fig. 1). From this it follows that Ukk/3 is invariant along an isomorph (

$\tilde{r} \equiv \rho ^{1/3}r$
r̃ρ1/3r⁠):

\begin{eqnarray}\frac{ U_k }{\rho ^{k/3}}& =& \left< \sum _{i>j}\varepsilon ^{(k)}_{ij}\left(\frac{\sigma ^{(k)}_{ij}}{\tilde{r}_{ij}}\right)^{k}\right>\end{eqnarray}
Ukρk/3=i>jɛij(k)σij(k)r̃ijk
(16)
\begin{eqnarray}&=& \sum _{i>j}\varepsilon ^{(k)}_{ij} \left({\sigma ^{(k)}_{ij}}\right)^k \left<{\tilde{r}_{ij}}^{-k}\right>\,.\end{eqnarray}
=i>jɛij(k)σij(k)kr̃ijk.
(17)

If we let “*” denote a reference point, (W*, U*), we can use Eqs. (14) and (15) to get

$(U_m^*, U_n^*)$
(Um*,Un*) and find for other state points on the same isomorph (where
$\tilde{\rho }\equiv \rho /\rho _*$
ρ̃ρ/ρ*
):

\begin{equation}U_k=\left(\frac{\rho }{\rho _*}\right)^{k/3}U_k^*=\tilde{\rho }^{k/3}U_k^*.\end{equation}
Uk=ρρ*k/3Uk*=ρ̃k/3Uk*.
(18)

Combining Eq. (18) with Eqs. (11) and (13) we obtain a parametric description of an isomorph in the WU phase diagram:

\begin{eqnarray}U &=& \tilde{\rho }^{m/3}U_m^* + \tilde{\rho }^{n/3}U_n^*,\end{eqnarray}
U=ρ̃m/3Um*+ρ̃n/3Un*,
(19)
\begin{eqnarray}W &=& \frac{m}{3}\tilde{\rho }^{m/3}U_m^* +\frac{n}{3}\tilde{\rho }^{n/3}U_n^* \,.\end{eqnarray}
W=m3ρ̃m/3Um*+n3ρ̃n/3Un*.
(20)

Once the numbers

$U_m^*$
Um* and
$U_n^*$
Un*
have been determined at the reference state point using Eqs. (14) and (15), the entire isomorph to which this state point belongs is thus traced out in the WU phase diagram using Eqs. (19) and (20). Alternatively, given a collection of isomorphic state point generated, e.g., as described in Sec. II C, the predicted relation between U, W, and
$\tilde{\rho }$
ρ̃
can be tested by linear regression:

\begin{eqnarray}\frac{U}{\tilde{\rho }^{n/3}} &=& \tilde{\rho }^{(m-n)/3}U_m^* + U_n^*,\end{eqnarray}
Uρ̃n/3=ρ̃(mn)/3Um*+Un*,
(21)
\begin{eqnarray}\frac{W}{\tilde{\rho }^{n/3}} &=& \frac{m}{3}\tilde{\rho }^{(m-n)/3}U_m^* + \frac{n}{3}U_n^* \,.\end{eqnarray}
Wρ̃n/3=m3ρ̃(mn)/3Um*+n3Un*.
(22)

This is done in Fig. 6 for the 12-6 KABLJ isomorph presented in Figs. 1 and 2. The potential energies were fitted by linear regression by means of Eq. (21). Subsequently, the virials were compared to Eq. (22) using the parameters estimated from the potential energies (Eq. (21)), i.e., without performing a new fit.

FIG. 1.

(a) Radial distribution function (AA particles) for three isomorphic state points of the Kob-Andersen binary LJ mixture45 (KABLJ). (b) Radial distribution function along an isotherm of the KABLJ mixture. A similar picture applies to the AB and BB radial distribution functions, though they are slightly less isomorph invariant (Fig. 2 in Paper IV).

FIG. 1.

(a) Radial distribution function (AA particles) for three isomorphic state points of the Kob-Andersen binary LJ mixture45 (KABLJ). (b) Radial distribution function along an isotherm of the KABLJ mixture. A similar picture applies to the AB and BB radial distribution functions, though they are slightly less isomorph invariant (Fig. 2 in Paper IV).

Close modal
FIG. 6.

Test of a 12-6 KABLJ isomorph (same as in Figs. 1 and 2). Top panel: U fitted to Eq. (21) using ρ* = 1 (i.e.,

$\tilde{\rho }= \rho$
ρ̃=ρ⁠). Bottom panel: W compared to Eq. (22), using the parameters
$U_m^*$
Um*
and
$U_n^*$
Un*
found in the top panel.

FIG. 6.

Test of a 12-6 KABLJ isomorph (same as in Figs. 1 and 2). Top panel: U fitted to Eq. (21) using ρ* = 1 (i.e.,

$\tilde{\rho }= \rho$
ρ̃=ρ⁠). Bottom panel: W compared to Eq. (22), using the parameters
$U_m^*$
Um*
and
$U_n^*$
Un*
found in the top panel.

Close modal
FIG. 2.

(a) Incoherent intermediate scattering function of the large (A) particles for isomorphic state points (density variation 11.6%) of the KABLJ mixture; time is given in reduced units:

$\tilde{t} \equiv t \rho ^{1/3}\sqrt{k_BT/m}$
t̃tρ1/3kBT/m and the q-vector is kept constant in reduced units. (b) Same as in (a), but along an isotherm (density variation 9.5%).

FIG. 2.

(a) Incoherent intermediate scattering function of the large (A) particles for isomorphic state points (density variation 11.6%) of the KABLJ mixture; time is given in reduced units:

$\tilde{t} \equiv t \rho ^{1/3}\sqrt{k_BT/m}$
t̃tρ1/3kBT/m and the q-vector is kept constant in reduced units. (b) Same as in (a), but along an isotherm (density variation 9.5%).

Close modal

We note here two important consequences of the above. (1) If the reference state point (U*, W*) generates the isomorph

$(U(\tilde{\rho }),W(\tilde{\rho }))$
(U(ρ̃),W(ρ̃))⁠, then the reference state point (aU*, aW*) generates the isomorph
$(aU(\tilde{\rho }),aW(\tilde{\rho }))$
(aU(ρ̃),aW(ρ̃))
. This means that all isomorphs for a given system have the same shape in the UW diagram and can be scaled onto one another. (2) The shape of the isomorphs depends only on the exponents, m and n. This means that all 12-6 LJ systems have the same isomorphs in the WU phase diagram – they may be scaled onto a single “master isomorph” (see also Eq. (27) below).

Figure 7 shows four isomorphs in the WU phase diagram for three different 12-6 LJ systems.47 Figure 8 shows the same isomorphs scaled onto the 12-6 master isomorph, Figure 9 shows an isomorph for the 18-6 KABLJ system and compares it to a 12-6 LJ isomorph. The agreement between simulation results and predicted shapes is quite good.

FIG. 7.

Four different 12-6 LJ isomorphs in the WU phase diagram. The points give simulation results, the full curves the predictions of Eqs. (19) and (20). Two isomorphs belong to the KABLJ system, the third results from simulations of the Wahnström binary LJ liquid46 (WBLJ), and the fourth (SCLJ) is the single component LJ liquid.

FIG. 7.

Four different 12-6 LJ isomorphs in the WU phase diagram. The points give simulation results, the full curves the predictions of Eqs. (19) and (20). Two isomorphs belong to the KABLJ system, the third results from simulations of the Wahnström binary LJ liquid46 (WBLJ), and the fourth (SCLJ) is the single component LJ liquid.

Close modal
FIG. 9.

Isomorph of the 18-6 KABLJ liquid compared to the theoretical prediction (Eqs. (19) and (20)). For reference, the isomorph prediction for 12-6 LJ systems is included.

FIG. 9.

Isomorph of the 18-6 KABLJ liquid compared to the theoretical prediction (Eqs. (19) and (20)). For reference, the isomorph prediction for 12-6 LJ systems is included.

Close modal
FIG. 8.

The “master isomorph”: collapsing the four isomorphs of Fig. 7 by scaling with

$W_0^*$
W0* defined as the virial on each isomorph when U = 0. The points give simulation results, the dashed curve the prediction of Eq. (28).

FIG. 8.

The “master isomorph”: collapsing the four isomorphs of Fig. 7 by scaling with

$W_0^*$
W0* defined as the virial on each isomorph when U = 0. The points give simulation results, the dashed curve the prediction of Eq. (28).

Close modal

Combining Eq. (18) with Eqs. (14) and (15) we get

\begin{eqnarray}\tilde{\rho }^{m/3} = \frac{U_m}{U_m^*} = \frac{3W-nU}{3W_*-nU_*},\end{eqnarray}
ρ̃m/3=UmUm*=3WnU3W*nU*,
(23)
\begin{eqnarray}\tilde{\rho }^{n/3} = \frac{U_n}{U_n^*} = \frac{3W-mU}{3W_*-mU_*}.\end{eqnarray}
ρ̃n/3=UnUn*=3WmU3W*mU*.
(24)

Eliminating density and rearranging, we find an invariant for the isomorph:

\begin{eqnarray}\frac{ \left( {W -\frac{m}{3}U }\right)^m }{ \left( {W -\frac{n}{3}U } \right)^n } = \frac{ \left( {W_*-\frac{m}{3}U_*}\right)^m }{ \left( {W_*-\frac{n}{3}U_*} \right)^n }.\end{eqnarray}
Wm3UmWn3Un=W*m3U*mW*n3U*n.
(25)

Choosing the reference point with zero potential energy

$(U_*,W_*) = (0, W^*_{0})$
(U*,W*)=(0,W0*) (this reference point is guaranteed to exist if all interactions have attraction), we get:

\begin{eqnarray}\frac{ \left( {W -\frac{m}{3}U }\right)^m }{ \left( {W -\frac{n}{3}U } \right)^n } = \left(W^*_{0}\right)^{m-n}.\end{eqnarray}
Wm3UmWn3Un=W0*mn.
(26)

$W^*_{0}$
W0* is a unique number identifying the isomorph. Equation (26) implies that the isomorph is the solution to

\begin{eqnarray}\left( { \frac{W}{W^*_{0}} -\frac{m}{3} \frac{U}{W^*_{0}} }\right)^{m/n} = \left( { \frac{W}{W^*_{0}} -\frac{n}{3} \frac{U}{W^*_{0}} } \right).\end{eqnarray}
WW0*m3UW0*m/n=WW0*n3UW0*.
(27)

Here it can be seen directly that – as argued above – for a given set of exponents (m, n) all isomorphs have the same shape, suggesting the name “master isomorph”;

$W^*_{0}$
W0* is the parameter that determines the scale of each isomorph. For m = 2n (e.g., the standard 12-6 LJ potential) the solution to this (quadratic) equation is:

\begin{equation}2\frac{W}{W^*_{0}}={1 + 4\frac{n}{3}\frac{U}{W^*_{0}}\pm \sqrt{1 + 4\frac{n}{3}\frac{U}{W^*_{0}}}}{}.\end{equation}
2WW0*=1+4n3UW0*±1+4n3UW0*.
(28)

This has real solution(s) whenever

${U}/{W^*_{0}}\ge -\frac{3}{4n}\,$
U/W0*34n⁠, where equality corresponds to W = 0. Figure 5 plots Eq. (28) with n = 6 for two values of
$W^*_{0}$
W0*
.

In Paper IV it was argued briefly that Eq. (7) predicts that solid-liquid coexistence curves are isomorphs. In more detail the argument goes as follows. First one notes that an isomorph cannot cross a solid-liquid coexistence curve: recall that the isomorph concept relates to equilibrium ensemble probabilities. For two isomorphic state points all pairs of microscopic configurations related by Eq. (6) have identical Boltzmann probabilities according to Eq. (7). Consequently, the two isomorphic state points must also have the same macroscopic phase behavior. In other words, if an isomorph were to cross the coexistence curve, part of the isomorph would have liquid microstates as the most likely and part of it would have mixed crystalline-liquid microstates as the most likely (having in mind a density-temperature phase diagram). This would violate the proportionality of Boltzmann statistical weights for scaled states, the property that defines an isomorph.

Given that an isomorph cannot cross the solid-liquid coexistence curve, the next point is to note that—because all state points belong to some isomorph—a curve infinitesimally close to the coexistence curve must be an isomorph. This argument applies in the W, U phase diagram, as well as in the density-temperature or the pressure-temperature phase diagrams. In the first two cases there is a multitude of isomorphs in the coexistence region; for instance, one particular isomorph will be characterized by 53% liquid and 47% crystal, etc. In the pressure-temperature phase diagram these “coexistence” isomorphs collapse into one, of course.

The prediction that the melting curve, both slightly to the liquid side and slightly to the solid side, are isomorphs has a number of consequences, and it sheds new light on some well-known phenomenological melting rules (as well as exceptions to these rules for non-strongly correlating liquids, see Paper IV). Thus, for instance, the observation that the Lindemann melting criterion is pressure independent for a given liquid follows from the isomorph property of the melting line (whereas the existence of a universal Lindemann criterion does not follow). Indeed, as is easy to show from paper IV the isomorph theory implies Ross’ “generalized Lindemann melting law” from 1969,48 according to which the reduced configurational partition function of the crystalline phase

$Q^*=\int d\tilde{\bf R} \exp [-\beta (U(\tilde{\bf R})-U(0))]$
Q*=dR̃exp[β(U(R̃)U(0))] is invariant along the melting curve (here, U(0) is the potential energy of the crystal with all atoms located at their lattice sites, and the tilde denotes reduced coordinates). Further well-known melting rules, which the isomorph theory explains, include (Paper IV): the reduced viscosity, the reduced surface tension, the reduced diffusion constant, and the reduced heat conductivity are all invariant in the liquid phase along the melting curve. Likewise, the reduced-unit static structure factor is invariant along the melting curve. These rules have been confirmed for several liquids in both experiment and simulation (see the references of Paper IV).

Like other isomorph predictions following from Eq. (7), the prediction that solid-liquid coexistence curves are isomorphs is only exactly obeyed for systems that obey Eq. (7) exactly, i.e., perfectly correlating systems. Indeed, it is straightforward to show that in IPL systems the phase behavior depends only on ρn/3/T, and thus the co-existence curve is an exact isomorph. The interesting question now is: how well do strongly correlating systems, e.g., generalized LJ systems, follow the prediction? Note that we do not expect the liquid-gas co-existence curves to be isomorphs (and indeed they are not); the gas-phase is not strongly correlating and thus Eq. (7) cannot be fulfilled in the relevant part of configuration space. In contrast, at solid-liquid coexistence both phases involved are strongly correlating for LJ systems (Paper I). Below we briefly test to which degree the isomorph prediction is fulfilled.

Ahmed and Sadus used in a recent paper49 a novel numerical method for determining solid-liquid coexistence. They reported results for the generalized LJ potential with the repulsive exponent m varying from 12 to 7 and fixed attractive exponent n = 6 (Ahmed and Sadus termed the repulsive exponent n, not m; we keep, however, the above notation). Figure 10 shows the state points reported by Ahmed and Sadus for the liquid side of the solid-liquid coexistence curve compared to the isomorphic predictions as expressed in Eqs. (21) and (22). The right panel shows that the parameters found by fitting to Eqs. (21) and (22), respectively, are consistent. Figure 11 shows the same data plotting potential energy (left panel) and virial (right panel) versus density. Figure 12 shows the 12-6 LJ solid-liquid phase behavior in the WU-diagram and in the ρT-diagram (inset), including here also data from Mastny and de Pablo.50 

FIG. 10.

Test of isomorph predictions for the liquid side of the solid-liquid coexistence curve of the single component generalized LJ potential (n = 6). Left panel: Eq. (21), middle panel: Eq. (22); data points are simulation results from Ahmed and Sadus,49 straight lines are linear regression fits. The right panel plots the fit parameters from Eqs. (21) and (22) against each other, demonstrating that the fits to energies and virials, respectively, are consistent.

FIG. 10.

Test of isomorph predictions for the liquid side of the solid-liquid coexistence curve of the single component generalized LJ potential (n = 6). Left panel: Eq. (21), middle panel: Eq. (22); data points are simulation results from Ahmed and Sadus,49 straight lines are linear regression fits. The right panel plots the fit parameters from Eqs. (21) and (22) against each other, demonstrating that the fits to energies and virials, respectively, are consistent.

Close modal
FIG. 12.

Solid-liquid phase behavior for the single-component 12-6 LJ system. Main figure: filled symbols mark the coexistence region in the WU diagram. Blue circles are from Ref. 49 and red squares are from Ref. 50. Full curves are the isomorph prediction (Eq. (28)) with

$W_0^*/N=49.8$
W0*/N=49.8 and 56.8, respectively. Open circles indicate results from simulations in the solid and liquid phase, respectively. Both phases are strongly correlating with R ⩾ 0.99. The inset shows the solid-liquid phase behavior in a ρ, T-diagram with logarithmic axes. Dashed red curves are the isomorph prediction ργ/T = const using γ found using Eq. (9) from the simulations of the two phases (open circles) giving γ = 5.24 and 5.18, respectively. Full black curves are free fits giving γ = 5.13 and 6.49, respectively.

FIG. 12.

Solid-liquid phase behavior for the single-component 12-6 LJ system. Main figure: filled symbols mark the coexistence region in the WU diagram. Blue circles are from Ref. 49 and red squares are from Ref. 50. Full curves are the isomorph prediction (Eq. (28)) with

$W_0^*/N=49.8$
W0*/N=49.8 and 56.8, respectively. Open circles indicate results from simulations in the solid and liquid phase, respectively. Both phases are strongly correlating with R ⩾ 0.99. The inset shows the solid-liquid phase behavior in a ρ, T-diagram with logarithmic axes. Dashed red curves are the isomorph prediction ργ/T = const using γ found using Eq. (9) from the simulations of the two phases (open circles) giving γ = 5.24 and 5.18, respectively. Full black curves are free fits giving γ = 5.13 and 6.49, respectively.

Close modal
FIG. 11.

Same data as in Fig. 10, plotting potential energy versus density (left panel) and virial versus density (right panel). Dashed curves: isomorph predictions (Eqs. (19) and (20)), using the parameters determined in Fig. 10.

FIG. 11.

Same data as in Fig. 10, plotting potential energy versus density (left panel) and virial versus density (right panel). Dashed curves: isomorph predictions (Eqs. (19) and (20)), using the parameters determined in Fig. 10.

Close modal

The prediction that solid-liquid coexistence curves are isomorphs agrees well on the liquid side with the simulation results of Ahmed and Sadus49 as well as with the results for 12-6 LJ from Mastny and de Pablo.50 On the solid side the agreement is less convincing (see Fig. 12); in particular the logarithmic derivative in the ρT diagram (slope in the inset) differs by 25%. At present we do not have an explanation for this. It might be related to inter-particle distances longer than the first peak of the radial distribution function playing a more important role for the phase behavior on the solid side51 than on the liquid side.

Equations (19) and (20) give the potential energy and the virial as functions of density along an isomorph. The temperature variation along the isomorph can be found by integrating Eq. (9). Since γ generally changes only little (Fig. 3), in many situations it is a good approximation to assume it constant, leading to the relation

$T=\tilde{\rho }^\gamma T_*$
T=ρ̃γT* along an isomorph, where T* is the temperature at the reference point. Taken together, Eqs. (19), (20), and (9) thus imply that from a single reference point (ρ*, T*, U*, W*) we have a prediction for (ρ, T, U, W) on the isomorph to which the reference point belongs.

Figure 4 demonstrated that in the WU phase diagram, isochores to a good approximation are straight lines with slope γ, i.e., that one can write W(ρ, T) = W0(ρ) + γ(ρ)U(ρ, T). The only ingredient missing to generate a full equation of state expressed as U = U(ρ, T) and W = W(ρ, T) is a relation between temperature and potential energy or virial on an isochore.

Rosenfeld and Tarazona52 derived from density functional theory an expression for the potential energy on an isochore:

\begin{eqnarray}U(\rho , T) = U_0(\rho ) + \alpha (\rho ) T^{3/5}.\end{eqnarray}
U(ρ,T)=U0(ρ)+α(ρ)T3/5.
(29)

Rosenfeld and Tarazona noted that the expression “… provides a good estimate of the thermal energy (...) near freezing densities only for predominant repulsive interactions,” and confirmed their expression by simulations of IPL potentials with different exponents. Equation (29) has since been shown, however, to be an excellent approximation also for several models with attraction, including the KABLJ mixture.53–55 This was regarded as a bit of a mystery, but can now be understood as a consequence of the fact that dynamics and heat capacity of strongly correlating liquids are well reproduced by a corresponding IPL system54 – a consequence of the hidden scale invariance of strongly correlating liquids discussed briefly in the introduction (see Papers I–III).

Figure 13 tests the Rosenfeld-Tarazona prediction for the systems investigated in Fig. 4. The prediction works very well, but less so for the supercooled Wahnström BLJ system, for which the formation of extended “Frank-Kasper” clusters at low temperatures affects the temperature dependence of the potential energy.56 

FIG. 13.

Test of Rosenfeld-Tarazona scaling (Eq. (29)) plotting potential energy versus T3/5 on isochores for the systems investigated in Fig. 4.

FIG. 13.

Test of Rosenfeld-Tarazona scaling (Eq. (29)) plotting potential energy versus T3/5 on isochores for the systems investigated in Fig. 4.

Close modal

Combining Eq. (29) with our results on isomorphs and isochores we can construct an approximate equation of state for generalized LJ systems. The idea is to map any state point (ρ, T, U, W) to the corresponding isomorphic state point (ρ*, T*, U*, W*) on a reference isochore where it is implicitly understood that the parameters (W0, γ, U0, α) are evaluated at the reference isochore ρ*:

\begin{eqnarray}T_*(\rho ,T) &=& T \tilde{\rho }^{-\gamma }, \tilde{\rho }\equiv \rho /\rho _*,\end{eqnarray}\vspace*{-19pt}
T*(ρ,T)=Tρ̃γ,ρ̃ρ/ρ*,
(30)
\begin{eqnarray}U_*(\rho ,T) &=& U_0 + \alpha (T \tilde{\rho }^{-\gamma })^{3/5},\end{eqnarray}
U*(ρ,T)=U0+α(Tρ̃γ)3/5,
(31)
\begin{eqnarray}W_*(\rho ,T) &=& W_0 + \gamma U_*(\rho ,T).\end{eqnarray}
W*(ρ,T)=W0+γU*(ρ,T).
(32)

From (U*, W*) we can calculate the contribution to the potential energy from the two IPL terms of the potential (compare Eqs. (14) and (15)):

\begin{eqnarray}U_{m}^*(\rho ,T) &=& \frac{ 3W_0 - (n - 3\gamma )\left( U_0 + \alpha (T \tilde{\rho }^{-\gamma })^{3/5}\right)}{m-n}, \nonumber\\\end{eqnarray}
Um*(ρ,T)=3W0(n3γ)U0+α(Tρ̃γ)3/5mn,
(33)
\begin{eqnarray}U_{n}^*(\rho ,T) &=& \frac{-3W_0 + (m - 3\gamma )\left( U_0 + \alpha (T \tilde{\rho }^{-\gamma })^{3/5}\right) }{m-n}. \nonumber\\\end{eqnarray}
Un*(ρ,T)=3W0+(m3γ)U0+α(Tρ̃γ)3/5mn.
(34)

Finally, we use isomorphic scaling (Eqs. (19) and (20)) to go back to the (ρ, T, U, W) state point:

\begin{eqnarray}\!\!\!\!\!\!\!\!\!\!\!\!\!U(\rho ,T) &=&\tilde{\rho }^{m/3} U_{m}^* (\rho ,T) + \tilde{\rho }^{n/3}U_{n}^*(\rho ,T),\end{eqnarray}
U(ρ,T)=ρ̃m/3Um*(ρ,T)+ρ̃n/3Un*(ρ,T),
(35)
\begin{eqnarray}\!\!\!\!\!\!\!\!\!\!\!\!\!W(\rho ,T) &=&\frac{m}{3}\tilde{\rho }^{m/3} U_{m}^*(\rho ,T) + \frac{n}{3}\tilde{\rho }^{n/3}U_{n}^*(\rho ,T).\end{eqnarray}
W(ρ,T)=m3ρ̃m/3Um*(ρ,T)+n3ρ̃n/3Un*(ρ,T).
(36)

The equation of state given by Eqs. (33) to (36) contains four parameters (W0, γ, U0, α) evaluated at a reference density ρ*. The parameters can be calculated from

$U,W,C_{\rm V}^{\rm ex}$
U,W,CV ex ⁠, and
$\beta _{\rm V}^{\rm ex}$
βV ex
at a single state point. A more convenient way of estimating the parameters might be: for a collection of (ρ, T, U, W) state points, at least two of which are non-isomorphic, use isomorphic scaling (Eqs. (19) and (20)) to get the corresponding (U*, W*) at the chosen reference density and fit them to Eq. (32) to determine W0 and γ. Next, using the fitted value for γ, fit to Eq. (31) in order to find U0 and α.

To put the new equation of state to a test, the procedure described above was used to determine the four parameters (W0, γ, U0, α) for the 12-6 KABLJ system using merely two state points. The resulting equation of state is shown as full curves in Fig. 14, where it is compared to simulation results shown as data points. The agreement is good. The largest deviations are seen when changing density away from the reference density, most evident for the isotherm in the WU diagram and the isomorph in the pT diagram. Better agreement is expected if the density dependence of γ is taken into account, but such a correction will not be attempted here. To summarize, combining the isomorph theory with two further approximations, a state-point independent exponent γ and the Rosenfeld-Tarazona expression, a realistic equation of state with four parameters is arrived at for generalized LJ systems.

FIG. 14.

Left panel: test of the equation of state in the WU phase diagram using ρ = 1.200 as reference density. Curves and lines: theoretical predictions (Eqs. (33) to (36)) for three isochores, an isomorph, and an isotherm. Data points: simulations of the 12-6 KABLJ mixture. The parameters of the equation of state were found from the two state points indicated by circles, giving (W0/N, γ, U0/N, α/N) = (38.4, 5.17, −8.63, 2.61). By comparison, the corresponding parameters found by fitting to the ρ = 1.200 isochore in Figs. 4 and 13 are (38.3, 5.16, −8.63, 2.61), illustrating the robustness of the fitted equation of state. Applying the fitting procedure to the simulated state points on the T = 0.550 isotherm gives the parameters (38.5, 5.19, −8.63, 2.61). Right panel: same data as above in a pT diagram.

FIG. 14.

Left panel: test of the equation of state in the WU phase diagram using ρ = 1.200 as reference density. Curves and lines: theoretical predictions (Eqs. (33) to (36)) for three isochores, an isomorph, and an isotherm. Data points: simulations of the 12-6 KABLJ mixture. The parameters of the equation of state were found from the two state points indicated by circles, giving (W0/N, γ, U0/N, α/N) = (38.4, 5.17, −8.63, 2.61). By comparison, the corresponding parameters found by fitting to the ρ = 1.200 isochore in Figs. 4 and 13 are (38.3, 5.16, −8.63, 2.61), illustrating the robustness of the fitted equation of state. Applying the fitting procedure to the simulated state points on the T = 0.550 isotherm gives the parameters (38.5, 5.19, −8.63, 2.61). Right panel: same data as above in a pT diagram.

Close modal

Paper IV demonstrated that strongly correlating liquids have isomorphs, i.e., curves in the phase diagram along which structure and dynamics to a good approximation are invariant in reduced units. The main new results in this paper are predictions specific to generalized Lennard-Jones systems and supporting numerical evidence.

Starting from the invariance of structure along isomorphs, a prediction for the shape of these in the WU phase diagram was presented and shown to agree well with simulation results. It was shown that for a given system, the isomorphs all have the same shape in the WU phase diagram; they are simply scaled versions of each other. Furthermore, the isomorph shape depends only on the exponents m and n, i.e., all systems with, say, 12-6 LJ interactions have the same isomorphs in the WU phase diagram. What differs from system to system is the values of the density and temperature along an isomorph.

Throughout Papers I-V the WU phase diagram was referred to extensively. This particular phase diagram has not been used much before, but it is evident that for strongly correlating liquids several features stand out when presented in the WU phase diagram. Prominent examples are that isochores are straight lines (with slope γ) and the existence of the “master” isomorph discussed above for generalized LJ systems.

The class of strongly correlating liquids includes, we believe, all van der Waals liquids and metallic liquids (but excludes hydrogen-bonding, covalent, as well as strongly ionic liquids). More simulations, as well as experiments, are needed to substantiate that the class of strongly correlating liquids is this large, but preliminary simulations of simple molecular models like28 the asymmetric dumbbell and the Lewis-Wahnstrom OTP three-site model are encouraging. In particular, the isomorph properties of more complex molecule models need to be simulated, as well as theoretically contemplated; the scaling properties of molecules with fixed bond lengths are not trivial. It should be noted that even a quite complex model system like a phospholipid membrane can be strongly correlating with respect to its slow degrees of freedom,57 the interactions of which are dominated by van der Waals forces. Another interesting question is to study in detail the properties of liquids where the exponent γ varies significantly throughout the phase diagram, for instance, the Weeks-Chandler-Andersen cutoff variant of the generalized binary LJ liquid that was recently investigated by Coslovich and Roland,58 as well as Berthier and Tarjus.59 Another system that deserves further investigation is the 200-100 LJ system60 which we find to be an exception to the general observation that generalized LJ systems are strongly correlating—presumably because of its extremely narrow range of interactions.

This series of papers investigated in depth the properties of strongly correlating liquids. We believe to have demonstrated that strongly correlating liquids are simpler than liquids in general. Hopefully this insight will turn out to be useful in furthering the understanding of the properties of liquids in general.

U.R.P. is supported by the Danish Council for Independent Research/Natural Sciences. The center for viscous liquid dynamics “Glass and Time” is sponsored by the Danish National Research Foundation (DNRF).

1.
N. P.
Bailey
,
U. R.
Pedersen
,
N.
Gnan
,
T. B.
Schrøder
, and
J. C.
Dyre
,
J. Chem. Phys.
129
,
184507
(
2008
) (Paper I).
2.
N. P.
Bailey
,
U. R.
Pedersen
,
N.
Gnan
,
T. B.
Schrøder
, and
J. C.
Dyre
,
J. Chem. Phys.
129
,
184508
(
2008
) (Paper II).
3.
T. B.
Schrøder
,
N. P.
Bailey
,
U. R.
Pedersen
,
N.
Gnan
, and
J. C.
Dyre
,
J. Chem. Phys.
131
,
234503
(
2009
) (Paper III).
4.
N.
Gnan
,
T. B.
Schrøder
,
U. R.
Pedersen
,
N. P.
Bailey
, and
J. C.
Dyre
,
J. Chem. Phys.
131
,
234504
(
2009
) (Paper IV).
5.
U. R.
Pedersen
,
N. P.
Bailey
,
T. B.
Schrøder
, and
J. C.
Dyre
,
Phys. Rev. Lett.
100
,
015701
(
2008
).
6.
M. P.
Allen
and
D. J.
Tildesley
,
Computer Simulation of Liquids
(
Oxford University Press
,
New York
,
1987
).
7.
J. P.
Hansen
and
J. R.
McDonald
,
Theory of Simple Liquids
, 3ed ed. (
Academic
,
New York
,
2005
).
8.
W.
Kauzmann
,
Chem. Rev.
43
,
219
(
1948
).
9.
G.
Harrison
,
The Dynamic Properties of Supercooled Liquids
(
Academic
,
New York
,
1976
).
10.
S.
Brawer
,
Relaxation in Viscous Liquids and Glasses
(
American Ceramic Society
,
Columbus, OH
,
1985
).
11.
I.
Gutzow
and
J.
Schmelzer
,
The Vitreous State: Thermodynamics, Structure, Rheology, and Crystallization
(
Springer
,
Berlin
,
1995
).
12.
M. D.
Ediger
,
C. A.
Angell
, and
S. R.
Nagel
,
J. Phys. Chem.
100
,
13200
(
1996
).
13.
C. A.
Angell
,
K. L.
Ngai
,
G. B.
McKenna
,
P. F.
McMillan
, and
S. W.
Martin
,
J. Appl. Phys.
88
,
3113
(
2000
).
14.
C.
Alba-Simionesco
,
C. R. Acad. Sci., Ser. IV
2
,
203
(
2001
).
15.
P. G.
Debenedetti
and
F. H.
Stillinger
,
Nature (London)
410
,
259
(
2001
).
16.
K.
Binder
and
W.
Kob
,
Glassy Materials and Disordered Solids: An Introduction to their Statistical Mechanics
(
World Scientific
,
Singapore
,
2005
).
17.
F.
Sciortino
,
J. Stat. Mech.
P05015
(
2005
).
18.
J. C.
Dyre
,
Rev. Mod. Phys.
78
,
953
(
2006
).
19.
N. L.
Ellegaard
,
T.
Christensen
,
P. V.
Christiansen
,
N. B.
Olsen
,
U. R.
Pedersen
,
T. B.
Schrøder
, and
J. C.
Dyre
,
J. Chem. Phys.
126
,
074502
(
2007
).
20.
N. P.
Bailey
,
T.
Christensen
,
B.
Jakobsen
,
K.
Niss
,
N. B.
Olsen
,
U. R.
Pedersen
,
T. B.
Schrøder
, and
J. C.
Dyre
,
J. Phys.: Condens. Matter
20
,
244113
(
2008
).
21.
T.
Christensen
and
J. C.
Dyre
,
Phys. Rev. E
78
,
021501
(
2008
).
22.
U. R.
Pedersen
,
T.
Christensen
,
T. B.
Schrøder
, and
J. C.
Dyre
,
Phys. Rev. E
77
,
011201
(
2008
).
23.
A.
Tölle
,
Rep. Prog. Phys.
64
,
1473
(
2001
).
24.
C.
Dreyfus
,
A.
Aouadi
,
J.
Gapinski
,
M.
Matos-Lopes
,
W.
Steffen
,
A.
Patkowski
, and
R. M.
Pick
,
Phys. Rev. E
68
,
011204
(
2003
).
25.
C.
Alba-Simionesco
,
A.
Cailliaux
,
A.
Alegria
, and
G.
Tarjus
,
Europhys. Lett.
68
,
58
(
2004
).
26.
R.
Casalini
and
C. M.
Roland
,
Phys. Rev. E
69
,
062501
(
2004
).
27.
C. M.
Roland
,
S.
Hensel-Bielowka
,
M.
Paluch
, and
R.
Casalini
,
Rep. Prog. Phys.
68
,
1405
(
2005
).
28.
T. B.
Schrøder
,
U. R.
Pedersen
,
N. P.
Bailey
,
S.
Toxvaerd
, and
J. C.
Dyre
,
Phys. Rev. E
80
,
041502
(
2009
).
29.
D.
Coslovich
and
C. M.
Roland
,
J. Phys. Chem. B
112
,
1329
(
2008
).
30.
D.
Coslovich
and
C. M.
Roland
,
J. Chem. Phys.
130
,
014508
(
2009
).
31.
W. G.
Hoover
,
D. A.
Young
, and
E.
Grover
,
J. Chem. Phys.
56
,
2207
(
1972
).
32.
Y.
Hiwatari
,
H.
Matsuda
,
T.
Ogawa
,
N.
Ogita
, and
A.
Ueda
,
Prog. Theor. Phys.
52
,
1105
(
1974
).
33.
L. V.
Woodcock
,
Phys. Rev. Lett.
54
,
1513
(
1985
).
34.
J. L.
Barrat
,
J. P.
Hansen
,
G.
Pastore
, and
E. M.
Waisman
,
J. Chem. Phys.
86
,
6360
(
1987
).
35.
P. G.
Debenedetti
,
F. H.
Stillinger
,
T. M.
Truskett
, and
C. J.
Roberts
,
J. Phys. Chem. B
103
,
7390
(
1999
).
36.
E.
La Nave
,
F.
Sciortino
,
P.
Tartaglia
,
M. S.
Shell
, and
P. G.
Debenedetti
,
Phys. Rev. E
68
,
032103
(
2003
).
37.
R. J.
Speedy
,
J. Phys.: Condens. Matter
15
,
S1243
(
2003
).
38.
M. S.
Shell
,
P. G.
Debenedetti
,
E.
La Nave
, and
F.
Sciortino
,
J. Chem. Phys.
118
,
8821
(
2003
).
39.
G.
Rickayzen
and
D. M.
Heyes
,
Phys. Rev. E
71
,
061204
(
2005
).
40.
A. C.
Branka
and
D. M.
Heyes
,
Phys. Rev. E
74
,
031202
(
2006
).
41.
R.
Casalini
,
U.
Mohanty
, and
C. M.
Roland
,
J. Chem. Phys.
125
,
014505
(
2006
).
42.
D. M.
Heyes
and
A. C.
Branka
,
Phys. Chem. Chem. Phys.
9
,
5570
(
2007
).
43.
Unless otherwise stated, simulations were performed using a newly developed molecular dynamics code optimized for NVIDIA graphics cards. The code is available as open source at http://rumd.org. Potentials were cut and shifted at 2.5σαβ. Reported potentials and virials do not include contributions beyond the cut-off.
44.
N.
Gnan
,
C.
Maggi
,
T. B.
Schrøder
, and
J. C.
Dyre
,
Phys. Rev. Lett.
104
,
125902
(
2010
).
45.
W.
Kob
and
H. C.
Andersen
,
Phys. Rev. Lett.
73
,
1376
(
1994
).
46.
G.
Wahnström
,
Phys. Rev. A
44
,
3752
(
1991
).
47.
The simulations reported in Figs. 7 and 8 were performed using Gromacs,
H. J. C.
Berendsen
,
D.
van der Spoel
, and
R.
van Drunen
,
Comp. Phys. Comm.
91
,
43
(
1995
);
E.
Lindahl
,
B.
Hess
, and
D.
van der Spoel
,
J. Mol. Mod.
7
,
306
(
2001
).
48.
M.
Ross
,
Phys. Rev.
184
,
233
(
1969
).
49.
A.
Ahmed
and
R. J.
Sadus
,
J. Chem. Phys.
131
,
174504
(
2009
).
50.
E. A.
Mastny
and
J. J.
de Pablo
,
J. Chem. Phys.
127
,
104504
(
2007
).
51.
F. H.
Stillinger
,
J. Chem. Phys.
115
,
5208
(
2001
).
52.
Y.
Rosenfeld
and
P.
Tarazona
,
Mol. Phys.
95
,
141
(
1998
).
53.
F.
Sciortino
,
W.
Kob
, and
P.
Tartaglia
,
J. Phys.: Condens. Matter
12
,
6525
(
2000
).
54.
U. R.
Pedersen
,
T. B.
Schrøder
, and
J. C.
Dyre
,
Phys. Rev. Lett.
105
,
157801
(
2010
).
55.
B.
Coluzzi
,
G.
Parisi
, and
P.
Verrocchio
,
J. Chem. Phys.
112
,
2933
(
2000
).
56.
U. R.
Pedersen
,
T. B.
Schrøder
,
J. C.
Dyre
, and
P.
Harrowell
,
Phys. Rev. Lett.
104
,
105701
(
2010
).
57.
U. R.
Pedersen
,
G. H.
Peters
,
T. B.
Schrøder
, and
J. C.
Dyre
,
J. Phys. Chem. B
114
,
2124
(
2010
).
58.
D.
Coslovich
and
C. M.
Roland
,
J. Chem. Phys.
131
,
151103
(
2009
).
59.
L.
Berthier
and
G.
Tarjus
, e-print arXiv:1103.0432 (
2011
).
60.
L.
Angelani
,
G.
Foffi
,
F.
Sciortino
, and
P.
Tartaglia
,
J. Phys.: Condens. Matter
17
,
L113
(
2005
).