Many-body correlations from integral geometry

TL;DR

This paper develops a morphometric approach based on integral geometry to predict many-body correlations in hard sphere liquids, extending scaled particle theory and validating with simulations.

Contribution

It introduces a new morphometric theory for many-body correlations in hard sphere liquids, generalizing scaled particle theory and demonstrating its effectiveness through numerical tests.

Findings

New morphometric theory for many-body correlations

Effective in predicting local structures in liquids

Validated against simulation data

Abstract

In a recent letter we presented a framework for predicting the concentrations of many-particle local structures inside the bulk liquid as a route to assessing changes in the liquid approaching dynamical arrest. Central to this framework was the morphometric approach, a synthesis of integral geometry and liquid state theory, which has traditionally been derived from fundamental measure theory. We present the morphometric approach in a new context as a generalisation of scaled particle theory, and derive several morphometric theories for hard spheres of fundamental and practical interest. Our central result is a new theory which is particularly suited to the treatment of many-body correlation functions in the hard sphere liquid, which we demonstrate by numerical tests against simulation.

Equations

$$⟨\cdots ⟩=\frac{1}{\Xi }\sum _{N=0}^{\infty }\frac{z^N}{N!}\int \left(\cdots \right)e^{-\beta U_N}d{\mathbf{r}}^N,$$

$$\mathrm{Prob}\left[\text{any }n\text{ particles in volume }d{\mathbf{r}}^n\right]:={\rho }^{(n)}({\mathbf{r}}^n)d{\mathbf{r}}^n.$$

$${\rho }^{(n)}({\mathbf{r}}^n)=\frac{1}{\Xi }\sum _{N=n}^{\infty }\frac{z^N}{(N-n)!}\int e^{-\beta U_N}d{\mathbf{r}}^{(N-n)}.$$

$${\rho }^{(n)}({\mathbf{r}}^n)$$

$$=z^n{e}^{-\beta U_n}⟨e^{-\beta U_{n\leftrightarrow N}}⟩$$

$${\rho }^{(n)}({\mathbf{r}}^n)=z^n{e}^{-\beta (U_n+\Omega -{\Omega }_{\mathrm{hom}})}.$$

$${\rho }^{(n)}({\mathbf{r}}^n)={\rho }^n{e}^{-\beta (U_n+\Omega -{\Omega }_{\mathrm{hom}}-n{\mu }^{\mathrm{ex}})}.$$

$$g^{(n)}({\mathbf{r}}^n):=\frac{{\rho }^{(n)}({\mathbf{r}}^n)}{{\rho }^n}=e^{-\beta (U_n+\Delta \Omega -n{\mu }^{\mathrm{ex}})}$$

$$\begin{array}{rl}{\varphi }^{(n)}({\mathbf{r}}^n)& :=-k_{B}T\ln g^{(n)}({\mathbf{r}}^n)\\ & =U_n+\Delta \Omega -n{\mu }^{\mathrm{ex}},\end{array}$$

$${\Omega }_{\mathrm{hom}}=-pV.$$

$$\Omega =-p{V}_{\mathcal{R}}+{\Omega }_{\mathrm{ex}}[\partial \mathcal{L}],$$

$$\Delta \Omega :=\Omega -{\Omega }_{\mathrm{hom}}=p{V}_{\mathcal{L}}+{\Omega }_{\mathrm{ex}}[\partial \mathcal{L}].$$

$${\Omega }_{\mathrm{ex}}[\partial \mathcal{L}]=\gamma [\partial \mathcal{L}]A_{\partial \mathcal{L}}$$

$$\Delta \Omega [\mathcal{L}]=p{V}_{\mathcal{L}}+\gamma [\partial \mathcal{L}]A_{\partial \mathcal{L}}.$$

$$\Delta \Omega (R)=p\frac{4\pi R^3}{3}+a_{2}4\pi R^2+a_{1}4\pi R+a_{0}4\pi ,$$

$$\Delta \Omega [K]=pV[K]+a_{2}A[K]+a_{1}C[K]+a_{0}X[K],$$

$$A[K]$$

$$C[K]$$

$$X[K]$$

$$\Delta \Omega [K]=p{V}_{+}[K]+a_2^{\prime }{A}_{+}[K]+a_1^{\prime }{C}_{+}[K]+a_0^{\prime }{X}_{+}[K],$$

$$X_{+}[K]$$

$$C_{+}[K]$$

$$A_{+}[K]$$

$$V_{+}[K]$$

$$a_0^{\prime }$$

$$a_1^{\prime }$$

$$a_2^{\prime }$$

$$\beta a_0^{SPT/PY}$$

$$\beta a_1^{SPT/PY}$$

$$\beta a_2^{SPT/PY}$$

$$\frac{\beta p^{SPT/PY}}{\rho }$$

$$\beta a_0^{SPT/CS}$$

$$\beta a_1^{SPT/CS}$$

$$\beta a_2^{SPT/CS}$$

$$\frac{\beta p^{SPT/CS}}{\rho }$$

$$\begin{array}{rl}{\varphi }^{(2)}(r):=& -k_{B}T\ln g^{(2)}(r)\\ =& pV(r)+a_{2}A(r)+a_{1}C(r)+a_{0}X(r)-2{\mu }^{\mathrm{ex}}[p].\end{array}$$

$$g^{(2)}(\sigma )=\frac{3}{2\pi {\sigma }^3\rho }\left(\frac{\beta p}{\rho }-1\right).$$

Full text

I Introduction

Since the beginnings of modern liquid state theory Kirkwood (1935), the hard sphere liquid has remained the archetypal model for atomic systems and soft matter. The dynamics of the system at high density in the metastable regime above the freezing transition are hotly debated, despite relentless study. Proposed mechanisms for dynamical phenomena all loosely fall under the broad umbrella of many-body correlations; nucleation occurs via crystal seed formation Sear (2007), and to explain dynamic arrest approaching the glass transition thermodynamic theories invoke cooperatively rearranging regions Lubchenko and Wolynes (2007) or elastic soft modes Brito and Wyart (2009) while kinetic theories posit the existence of dynamical defects Chandler and Garrahan (2010). In a recent letter Robinson et al. (2019) we proposed a framework for treating many-body correlations, and developed an operational scheme for predicting the populations and dynamics of local structural motifs within a uniform liquid. Central to this is the use of the morphometric approach.

The morphometric approach provides an efficient means of treating the thermodynamics of a bulk liquid without fully determining its equilibrium density profile König et al. (2004); Roth et al. (2006); Hansen-Goos et al. (2007); Robinson et al. (2019). Detailed investigations have shown that it is highly accurate in the hard sphere liquid regime Oettel et al. (2009); Ashton et al. (2011); Laird et al. (2012); Blokhuis (2013); Urrutia (2014); Hansen-Goos (2014), so we can expect an accurate treatment where the bulk system provides background depletion interactions while its detailed microstates remain unimportant. This feature makes it ideally suited for many-body correlations if we can identify relevant dynamical degrees of freedom. While existing morphometric theories have been proven accurate in the liquid regime, we require a theory which works in the supercooled regime. Here we derive such a theory using scaled particle theory (SPT).

SPT determines bulk properties from consideration of a spherical solute of varying radius. It remains one of most enduring theories of simple liquids; though 60 years old as of this year Reiss et al. (1959), aspects of this approach remain in modern theories. This is particularly true for hard spheres where SPT has been unified with the Percus-Yevick integral equation solution Wertheim (1963), another old theory, in the form of fundamental measure theory (FMT) Rosenfeld (1989). Though originally a theory of single-component hard spheres Reiss et al. (1959), SPT has been extended to other potentials Reiss et al. (1960); Helfand et al. (1960); Reiss and Mayer (1961) and shapes Gibbons (1969, 1970), mixtures Lebowitz et al. (1965), dimers Stillinger et al. (2006); Chatterjee et al. (2006) and disks Helfand et al. (1961); Martin et al. (2018); Hansen-Goos (2019). Morphological thermodynamics can be seen as a modern generalisation of SPT for a wide class of physically relevant geometries. Its basis in integral geometry replaces the semi-empirical approach of classical SPT with clearly defined postulates. In this work we present the morphometric approach in the context of SPT and derive a new theory suitable for high densities above freezing. In additional appendices, we show that minor modifications of our arguments can be used to derive previous theories: the classical SPT coefficients, and the White Bear II morphometric coefficients of Ref. Hansen-Goos and Roth (2006a).

In section II we show how one can map the problem of treating many-body correlations onto a solvent-solute problem. We spend the rest of the paper discussing the solvation problem through the lens of SPT. We introduce the morphometric approach as a useful generalisation of SPT in section III, and derive a theory well-suited for treating many-body correlations using scaled particle arguments. In section IV we numerically test these theories’ two– and three–body correlation functions to demonstrate their effectiveness in treating correlation functions.

II Solvation expression for many-body correlations

II.1 Correlations in terms of the insertion cost

We will show that correlations of $n$ particles at positions $\boldsymbol{\mathbf{r}}^{n}:=\{\boldsymbol{\mathbf{r}}_{1},\cdots,\boldsymbol{\mathbf{r}}_{n}\}$ can be expressed in terms of the free energy cost of inserting them at $\boldsymbol{\mathbf{r}}^{n}$, by generalising the potential distribution theorem Widom (1963, 1982) to many particles. The classical approach, also known as Widom’s insertion method, expresses the (excess) chemical potential $\mu^{\mathrm{ex}}$ of a single-component system as the free energy cost of inserting an additional particle. See Ref. Rowlinson and Widom (2002) and references therein for a detailed review of this classical approach. Our generalisation results in a potential of mean force for interactions between the $n$ particles, which is formally identical to the chemical potential of a solute; this latter form is particularly suitable for geometric approximation schemes.

We consider a bulk liquid (the solvent) of $N$ particles with interaction potential energy $U_{N}$. Integrating over all solvent arrangements in the absence of any external field gives the (grand-canonical) average

[TABLE]

with partition function $\Xi:=e^{-\beta\Omega_{\mathrm{hom}}}$, where $\Omega_{\mathrm{hom}}=-pV$ is the usual homogeneous grand potential. The activity is $z=\exp{\beta\mu}/\Lambda^{d}$ in terms of the (total) chemical potential $\mu$ and the thermal de Broglie wavelength $\Lambda$.

Descriptions of many-body correlations naturally employ the $n$-particle density $\rho^{(n)}$, defined as

[TABLE]

The $n$-density can be obtained by integrating the full (configurational) probability distribution over the remaining degrees of freedom. For the single-component system this yields Hansen and McDonald (2013)

[TABLE]

Changing the summation limits $N\rightarrow N+n$ we obtain

[TABLE]

where in the latter step we decomposed the total potential $U_{N+n}$ into purely local and solvent terms, i.e. $U_{N+n}=U_{n}+U_{N}+U_{n\leftrightarrow N}$, where $U_{\alpha}$ for $\alpha\in\{n,N\}$ indicates the internal interactions between particles in component $\alpha$. The “interspecies” interactions are contained within $U_{n\leftrightarrow N}$ which acts as an external field for the solvent. Thus, (2) becomes

[TABLE]

where $\Omega$ is the grand potential of the solvent in the presence of the $n$-particle inhomogeneity. Splitting the chemical potential into its ideal and excess parts so that $\beta\mu=\ln{\Lambda^{d}\rho}+\beta\mu^{\mathrm{ex}}$ gives

[TABLE]

The $n$-particle distribution functions are then determined from Hansen and McDonald (2013)

[TABLE]

where $\Delta\Omega:=\Omega-\Omega_{\mathrm{hom}}$ is the reversible (free energy) cost of inserting the particles at fixed position $\boldsymbol{\mathbf{r}}^{n}$, or equivalently describes the average depletion interactions between mobile particles. For $n=1$ we have $\Delta\Omega=\mu^{\mathrm{ex}}$ and this is identical to the potential distribution theorem of Widom Widom (1963, 1982). The distribution functions are written in terms of the potential

[TABLE]

which we call the generalised potential of mean force. For the case $n=2$ this reduces to the usual potential of mean force in the liquid state literature Hansen and McDonald (2013).

This completes our proof that the correlations can be transformed to a potential, and we can proceed with a geometrical construction for $\Delta\Omega$.

II.2 Representing the insertion cost as a solvation problem

For systems with excluded volume interactions, we can divide the space into a local component $\mathcal{L}\subset\mathbb{R}^{d}$ of volume $V_{\mathcal{L}}$ inaccessible to solvent degrees of freedom and the remaining space $\mathcal{R}=\mathbb{R}^{d}\setminus\mathcal{L}$ of volume $V_{\mathcal{R}}$ filled by the rest of the liquid (Fig. 1). The total volume is $V=V_{\mathcal{L}}+V_{\mathcal{R}}$ so the homogeneous grand potential is

[TABLE]

After inserting the inhomogeneity the total volume accessible to the rest of the liquid will be reduced by $V_{\mathcal{L}}$, so the grand potential becomes

[TABLE]

where $\Omega_{\mathrm{ex}}$ is an excess term brought about by the introduction of a dividing surface $\partial\mathcal{L}$ between the two liquid components. Subtracting these two expressions gives

[TABLE]

This dividing surface has area $A_{\partial\mathcal{L}}$, creating a surface tension $\gamma$ so we can write the excess term as

[TABLE]

which is a formal definition of surface tension and depends on the choice of dividing surface (see two examples in Fig. 1b). We know from density functional theory Evans (1979) that the excess free energy is a functional of the density profile, which will in turn depend on the shape of the boundary; we write $\gamma=\gamma[\partial\mathcal{L}]$ to indicate this functional dependence on the surface shape. The solvation form of the inhomogeneous grand potential term in (4) is then

[TABLE]

The problem of determining the $n$-particle distributions has been reduced to a solvation problem: we must find the surface tension between a solute (the specific local arrangement) and a solvent (the rest of the liquid). We will use the solute–solvent terminology, but one could also think of local–bulk nomenclature.

III Obtaining a morphological theory for many-body correlations

We will consider a single-component hard sphere fluid, for particles of diameter $\sigma$ and bulk volume fraction $\eta$. Using the correspondence between many-body correlation functions and chemical potentials, we require an approximate model for solvation and a choice of surface in (5) to evaluate $\Delta\Omega$ in (4). We introduce our central approximation in section III.1 and our choice of surface in III.2. Then, we show that previous theories fail to produce accurate correlation functions at high densities in III.3 and derive a new theory to rectify this in III.4.

III.1 Our central approximation: the morphometric/scaled particle ansatz

Our key approximation, the morphometric approach, can be understood as a generalisation of scaled particle theory. In every formulation of scaled particle theory one considers a hard spherical solute of radius $R$. In most approaches, the cost $\Delta\Omega$ is assumed to have an analytic expansion in powers of the radius; in classical approaches this was simply postulated, however we will be able provide proper justification below through geometric arguments. Recognising that terms scaling faster than $R^{3}$ must be zero for it to remain well-defined in the limit of large solutes leads to the third-order polynomial Reiss et al. (1959)

[TABLE]

where we identified the largest power with the work term $pV$ from comparison with (5), and $\{a_{0},a_{1},a_{2}\}$ are thermodynamic coefficients describing the subleading corrections. We have chosen to introduce factors of $4\pi$ in front of the subleading terms to lead into the generalisation beyond spherical geometries. For a general solute $K\subset\mathbb{R}^{3}$ we then write the morphometric insertion cost as

[TABLE]

where $C$ and $X$ are the integrated mean and Gaussian curvatures. All of these functionals act on $K$ but the latter three can also be expressed as surface integrals, as in

[TABLE]

where $\kappa$ is the curvature tensor for the surface $\partial K$. For a spherical solute these reduce to the values given in (6), so this represents a proper generalisation of SPT for more general geometries. We give a brief justification of the ansatz (7) using integral geometric arguments in appendix A.

The key advantage of a geometric expansion of the free energy is that the role of thermodynamics and geometry are kept separate. Thermodynamics only enters through the coefficients $\{p,a_{2},a_{1},a_{0}\}$, so they can be determined in simple geometries to obtain a general theory. As a linear theory, only four (independent) equations are required to fix these coefficients; with many thermodynamic relations to choose from this approximate theory is overconstrained in general. We must use physical intuition to choose suitable equations, after which the accuracy of the resulting coefficients can be assessed. After determining these coefficients all the complexity of computing $\Delta\Omega$ is reduced to measuring the geometric quantities $\{V,A,C,X\}$ of the specific solute. However, we must first specify a choice of surface $\partial\mathcal{L}$ in (5) before we can proceed.

III.2 Choice of dividing surface

All coefficients we give are for the molecular geometry bounded by the molecular surface ($\partial\mathcal{L}_{1}$ in Fig. 1b), the surface where interactions occur between the solute and a test particle representing the remaining liquid. However, it is usually more convenient to do calculations with the excluded geometry: the space inaccessible to the centre of a test particle bounded by the solvent accessible surface ($\partial\mathcal{L}_{2}$ in Fig. 1b). Note that there is also an infinite family of well-defined parallel surfaces between these two extremes, but they are not widely used in practice so we will not consider them Oettel et al. (2009). The choice of dividing surface will change the surface tension, and thus requires new coefficients $\{a_{0}^{\prime},a_{1}^{\prime},a_{2}^{\prime}\}$ i.e.

[TABLE]

where the excluded geometry terms transform via the canonical relations Hansen-Goos and Roth (2006a); Oettel et al. (2009); Santaló (2004); Klain and Rota (1997)

[TABLE]

It is straightforward to transform between these two conventions via Hansen-Goos and Roth (2006a); Oettel et al. (2009)

[TABLE]

The resulting $\Delta\Omega$ will be identical whichever surface is chosen, except when there is a topological change in the molecular surface marking the breakdown of the theory; this is discussed in detail in Ref. Oettel et al. (2009).

III.3 Failure of previous morphometric theories in treating correlations

Having specified the surface, we can examine the self-consistency of correlation functions determined through previously known morphological theories. We briefly state the main theories below, then proceed to show how they produce inaccurate correlation functions at high densities. This underscores the need for a more accurate theory, and the specific inconsistency we highlight in this section will be used to construct one in the next section.

With either the scaled particle or morphometric ansatzes, (6) or (7), a specific theory comprises the set of coefficients $\{p,a_{2},a_{1},a_{0}\}$. In appendix B we summarise the classical scaled particle arguments of Refs. Reiss et al. (1959); Lebowitz et al. (1965) using modern notation, which produce coefficients

[TABLE]

In this classical approach, the Percus-Yevick (PY) equation of state emerges as an output of the theory. More recently, morphometric theories have been obtained as the bulk limit of FMT, with the hitherto most successful theory determined in Ref. Hansen-Goos and Roth (2006a) as

[TABLE]

obtained from a functional constructed to impose the Carnahan-Starling (CS) equation of state (13d). The latter equation of state is known to be highly accurate across the whole stable liquid regime, and even at the high density limits accessible to simulation in the supercooled regime Berthier et al. (2016). The same equations are also obtained in the bulk limit of the functional of Ref. Santos (2012), which similarly imposes the CS pressure but is slightly more self-consistent. Curiously, we can make a minor modification to SPT arguments to impose the CS equation of state as an input to obtain the above coefficients without invoking FMT (details in appendix C). We thus label this theory as SPT/CS.

To demonstrate the inaccuracy of the correlation functions produced by these known theories using (3), we consider what happens to the pair correlation at high densities. The potential of mean force (4) for non-overlapping spheres with the morphometric ansatz (7) is written

[TABLE]

As a self-consistency test, we will compare this explicit result at contact against the exact value of $g^{(2)}(\sigma)$ predicted by the virial theorem as Hansen and McDonald (2013)

[TABLE]

To evaluate (14) we need to calculate the size measures for the two particle solute resembling a “dumbbell”. It is easier to calculate the excluded volume geometry, after which we can obtain the molecular volumes using the canonical relations (10). The excluded volume consists of the union of two balls of radius $\sigma$ separated by a distance $r$. The geometric properties at contact are then Oettel et al. (2009)

[TABLE]

Transforming to the parallel molecular surface using the inverse transformation of (10) gives the solute parameters as

[TABLE]

Fig. 2 shows the contact value $g^{(2)}(\sigma)$ from inserting the geometric parameters above into (14), and the quasi-exact result of (15) assuming the CS equation of state (13d). We find that both the known morphometric theories (12) and (13) are reasonably accurate until around the freezing density, above which contact correlations spuriously decay. Thus, a new theory is needed to treat correlations at high densities; in the next section we will construct one which satisfies (15) by construction.

III.4 Obtaining the new theory by self-consistency of the contact value of $g^{(2)}(r)$ with the virial theorem

Our goal is to develop a morphometric theory which produces accurate correlation functions $g^{(n)}$. As described at the end of the last section, the correlation functions produced by an SPT approach are inaccurate at high densities. We will correct the spurious decay of the contact value of the pair distribution function $g^{(2)}(r)$ at high densities by building this into the theory explicitly, with the aim of producing more accurate correlation functions. A working understanding of scaled particle arguments is necessary to follow the details of this derivation, which we lay out in appendices B and C.

Inserting the volumes at contact (16) into (14) and applying the virial theorem (15) for the contact value of $g^{(2)}$ gives the final expression

[TABLE]

We will use this last expression instead of the contact theorem (31) in order to obtain new coefficients. Together (24a), (25) and (17) solve to give coefficients:

[TABLE]

We refer to coefficients obtained this way for the CS pressure (13d) as virial/CS, but we will not give them explicitly. Unlike the WBII coefficients above these are new. The pair correlation produced by these coefficients (black line in Fig. 2) is self-consistent with CS at contact by construction.

IV Numerical results

We apply the thermodynamic coefficients determined in previous sections for a system of hard spheres to obtain two– and three–body distribution functions using the generalised potential of mean force (4) with the morphometric approach (7), and compare these against molecular dynamics simulations. For the analytics we determine the input geometric quantities $\{V,A,C,X\}$ using the algorithms of Refs. Mecke et al. (1994); Klenin et al. (2011). For the simulations we performed event-driven molecular dynamics of $N=1372$ monodisperse hard spheres using the DynamO software package Bannerman et al. (2010). We measure the pair and triplet distribution functions $g^{(2)}$ and $g^{(3)}$ for simulations at $\eta=0.45$. For simulations above freezing $\eta\simeq 0.494$ we used a 5-component equimolar distribution with $\sim 8\%$ polydispersity.

For $g^{(2)}$ shown in Fig. 3 we find the virial/CS theory outperforms the SPT/CS theory even away from contact. The agreement with the molecular dynamics simulations is excellent, until $r\gtrsim\sqrt{3}\sigma$ where the solute boundary self-intersects marking the end of the theory’s regime of validity. Geometrically, the regime $r<\sqrt{3}\sigma$ is the regime where the canonical relations (10) apply so the thermodynamics is independent of the choice of surface definition. Physically, for $r>\sqrt{3}\sigma$ interactions between solvent particles can occur through the solute, and these correlations are not captured by the theory. More discussion of this breakdown can be found in Ref. Oettel et al. (2009). Only the contact value was fixed, so accuracy for $r>\sigma$ was not guaranteed; the accuracy is a welcome bonus. We can quantify this accuracy through the integrated value

[TABLE]

shown in the top panel of Fig. 4 where we take $\delta=0.4\sigma$. We find this integrated quantity is within $10\%$ accuracy across the liquid regime for all three theories, with the new theory performing substantially better overall. Despite the improved accuracy, the errors begin to increase in magnitude at the end of the liquid regime so we expect them to become significant with very deep supercooling.

Next we compare the theories’ predicted surface tension against simulation data. The surface tension at a planar wall is simply $a_{2}$ because it conjugates with the area. In Ref. Davidchack et al. (2015) a highly accurate $a_{2}$ was measured for hard spheres through extensive simulation, which was parameterised by the following expression

[TABLE]

Comparing this highly accurate expression against the values predicted from the morphometric coefficients, we find the virial/CS surface tension is less accurate than the SPT/CS prediction (Fig. 4 bottom panel) despite its superior correlation functions at high densities. Moreover, we find that at low densities the new theory is less accurate than classical SPT/PY theory. This discrepancy occurs because both SPT/PY and SPT/CS feature the correct low density asymptotics of $a_{2}\sim\mathcal{O}(\eta)$, which is imposed through the radial derivative of $\Delta\Omega(r)$ in the point solute limit (24b). This suggests that the new virial/CS theory sacrifices asymptotic accuracy at low densities, for more self-consistency of the surface tension at moderate to high densities. One of the great strengths of the SPT/CS theory is its accuracy in the planar limit Hansen-Goos and Roth (2006a), and so SPT/CS coefficients may give more accurate grand potentials (and thus correlations) for large solutes where the surface becomes approximately planar.

Our goal was to develop a theory capable of treating correlations at the many-body level, so we now examine three-body correlation functions. Triplet geometries are characterised by a triangle of side lengths $r,s,t$ so $g^{(3)}=g^{(3)}(r,s,t)$. We also compare the morphometric theories against the Kirkwood approximation Kirkwood (1935) i.e.

[TABLE]

where we take the values of $g^{(2)}$ from the virial/CS theory because of its already demonstrated accuracy at the two body level. Comparison of the morphometric correlation functions, and the Kirkwood closure, against molecular dynamics are shown in Fig. 5. The virial/CS closure most closely matches the simulations at high densities, suggesting the theory is suitable for modeling complex many-particle local structures Robinson et al. (2019). For comparison we also include the tabulated values of Ref. Müller and Gubbins (1993) where $g^{(3)}$ is used to treat polyatomic molecules 111We believe there is a misprint in Eq. (A1) of Ref. Müller and Gubbins (1993), which should read $(1-\eta)^{3}$ in the denominator. We refitted their source data Attard and Stell (1992) with this corrected form to be sure.; our theory is marginally more accurate, and more importantly it provides a recipe for treating the higher-order correlation functions.

To quantify accuracy at the three-body level we consider the concentration of triangles with side lengths $r,s,t\in[\sigma,\sigma+\delta]$ in the bulk liquid, from (1) we find this as Krumhansl and Wang (1972)

[TABLE]

Comparison with molecular dynamics simulations in Fig. 6 shows similar levels of accuracy for small $\delta$, though the performance decreases as it is increased above the first minimum of the $g^{(2)}(r)$; this is not surprising as our virial closure only enforces accuracy approaching contact. Notably, the Kirkwood approximation (21) performs surprisingly well at the three-body level in both of these tests.

V Discussion and summary

We have presented the morphometric approach as a generalisation of SPT, thus placing the scaled particle ansatz on more precise and physically motivated assumptions i.e. those underlying the theorems of integral geometry. Using the scaled particle approach we have systematically derived a new theory capable of accurately calculating many-body correlations in the hard sphere liquid; we recently used this to accurately treat local structures in Ref. Robinson et al. (2019). Our scaled particle formalism is flexible enough to derive all known morphometric theories without invoking fundamental measure theory.

In principle this approach could be extended to simple liquids where the interaction potential can be approximated as a perturbation around a hard core e.g. the Lennard-Jones potential. However, as we exploited features of the hard sphere interaction potential to achieve closed form expressions for the thermodynamic coefficients, more realistic interaction potentials would likely require numerical expressions. Additionally, attractions can introduce non-analytic behaviour from wetting/drying transitions which would not be accounted for in our theory Evans et al. (2003, 2004).

By making the underlying assumptions explicit we can better understand the limits of the theory: any deviation from the morphometric/SPT ansatz must be due to a violation of translation/rotation invariance, additivity or continuity. The fact that these theories are very accurate for hard spheres suggests that the assumptions are only weakly violated for this system. While translational/rotational invariance and continuity are physically plausible conditions on $\Delta\Omega$, additivity is a very strong assumption. In particular, we expect significant deviations from additivity where the liquid develops a static length scale exceeding the size of the solute König et al. (2004). As such, we expect the validity of the morphometric approach to require the solute to be larger than the point-to-set length Montanari and Semerjian (2006), which acts as an upper bound for all structural length scales Yaida et al. (2016). The morphometric ansatz must break down approaching a critical point, so it cannot be used to obtain asymptotics in the event of a thermodynamic glass transition.

Finally, we remark that while it is tempting to call the treatment of bulk degrees of freedom with the morphometric approach mean-field, this is not a completely accurate characterisation. Mean-field theories typically become formally exact in the limit of infinite spatial dimensions, where the thermodynamic role of fluctuations disappears. By contrast, the morphometric approach (and related theories like SPT and FMT) become formally exact in the one-dimensional limit of hard rods. Though this theory does not explicitly describe fluctuations, they are built into the choice of thermodynamic coefficients entering the theory. In this sense it is more accurate to describe the morphometric approach (and related theories) as an excluded volume theory, or as a free volume theory because the thermodynamics only shows divergent behaviour as $\eta\to 1$.

Appendix A Justification of the scaled particle ansatz

We now give a brief justification of our main approximation (7), in particular why there are only four terms in the expansion. Radius is the only natural parameter for a sphere, however for more general geometries there might be arbitrarily many parameters so one may wonder if they should be included in a general geometric expansion. Nevertheless, there are compelling arguments from integral geometry König et al. (2004) to only retain the four terms listed which we will summarise below.

The basis of the morphometric approach is that the functionals $\{V,A,C,X\}$ are normalisations of the so-called intrinsic volumes. These play a central role in integral geometry as the only physically meaningful size measures in the sense that they:

Are invariant with respect to translations and rotations. 2. 2.

Increase additively, i.e. they transform under combination of subsystems via the inclusion/exclusion relation e.g.

[TABLE]

and similar expressions for $A$, $C$, and $X$. 3. 3.

Are continuous (specifically with respect to the Hausdorff metric). Loosely speaking, this means that the size measures converge as the object is approximated by increasingly finely meshed polyhedra excluding e.g. fractal geometries. As a simple intuitive example, the measurement of a length will converge continuously to some number as one uses rulers with progressively finer distance markings.

More details on properties of intrinsic volumes can be found in standard texts, e.g. Refs. Santaló (2004); Klain and Rota (1997).

The central assertion of the morphometric approach is that the insertion cost $\Delta\Omega$ exactly possesses the properties above, providing the connection between geometry and thermodynamics König et al. (2004). A classic theorem of integral geometry due to Hadwiger Hadwiger (1957) states that the intrinsic volumes are the only class of functionals with the properties listed above; a corollary of this is that they form a linear vector space for any functional possessing these properties. The morphometric form (7) then follows. In addition to providing a more general ansatz than SPT, this approach lays out its underlying assumptions explicitly eschewing the ad-hoc way in which the original SPT ansatz (6) was obtained. Moreover, classical SPT assumes hard spheres from the outset while our generalisation based on integral geometry is more flexible, allowing for generalisations to mixtures, more realistic pair potentials and non-spherical particles without compromising its assumptions.

The morphometric approach is certainly an approximation, as the insertion cost will not rigorously possess the three properties above in reality. Notably, in SPT $\Delta\Omega$ is known to contain singularities in its high order derivatives with solute radius Reiss et al. (1959); these non-analytic terms result from violations of the additivity assumption. Nevertheless, the approximation is accurate in hard spheres Oettel et al. (2009); Ashton et al. (2011); Laird et al. (2012); Blokhuis (2013); Urrutia (2014); Hansen-Goos (2014) so these violations should be small.

Appendix B Classical scaled particle relations

Following the protocol of scaled particle theories, we consider the insertion of a hard spherical solute of radius $R$ into the liquid. Assuming the morphometric form for the insertion cost returns us to the ansatz (6). Below we give the exact thermodynamic relations for hard spheres which produce the classical SPT coefficients.

It is possible to consider the insertion of a solute with a negative radius: the hard core interaction between the two particles only occurs when the solute is ‘inside’ a solvent particle. In this limit the insertion cost can be determined exactly as Reiss et al. (1959)

[TABLE]

for $-\frac{\sigma}{2}\leq R\leq 0$. It may appear concerning that this result does not possess the morphometric form (7); however, this does not discount the validity of the morphometric approach as the nonphysical geometry violates the continuity assumption (section III) because it cannot be approximated by polyhedra. This places the result for $R<0$ outside the theory’s stated regime of validity, however $\Delta\Omega$ is continuous up to its second derivative across $R=0$ with a discontinuity in its third derivative Reiss et al. (1959). In the limit $R\to 0$ the expression above corresponds to the cost of inserting a hard point giving

[TABLE]

Note that (24a) can also be justified by considering that the probability of a randomly selected position in space being empty is simply the free volume $1-\eta$.

Together applying (24) to (6) fixes the coefficients $\{a_{0},a_{1},a_{2}\}$, so the theory requires an additional thermodynamic relation to determine the pressure. When $R=\frac{\sigma}{2}$ the solute is equivalent to the solvent particles themselves and we recover $\Delta\Omega=\mu^{\mathrm{ex}}$, so from (6) we have

[TABLE]

Combining this expression with the thermodynamic relation (27) gives a differential equation for $\beta p$ whose solution gives the classical SPT coefficients for hard spheres (12). The equation of state (12d) is equivalent to the one obtained through the solution of the Percus-Yevick (PY) integral equation Wertheim (1963); these two routes have been unified within FMT Rosenfeld (1989).

Appendix C First generalisation: SPT with an empirical equation of state

In the classical SPT approach described in the previous section, the SPT/PY equation of state emerges as an output of the theory. Taking inspiration from the White Bear free energy functional Roth et al. (2002), we reformulate the SPT argument so that the equation of state is an input to the theory. In so doing we aim to construct a theory from a more accurate equation of state, with the trade-off being that we must sacrifice some self-consistency. The main equation of state we impose is the CS relation (13d) Carnahan and Starling (1969). This ultimately results in a theory previously known as a limit of a free energy functional Hansen-Goos and Roth (2006a), but through simpler arguments. We extend these arguments in the main text to arrive at a new theory capable of accurately treating correlation functions.

A crucial component of scaled particle approaches is thermodynamic consistency of the (osmotic) pressure via

[TABLE]

where $f^{\mathrm{ex}}=F^{\mathrm{ex}}/V$ is the (excess) free energy density. The form most useful for a single-component system comes from taking the derivative with respect to density, and noting the definition of the excess chemical potential

[TABLE]

giving

[TABLE]

Note that the consistency relation (26) provides the route to generalising all of our arguments to arbitrary mixtures Rosenfeld (1989); Sollich et al. (2001); Santos (2012). The free energy remains well-defined, even for polydisperse mixtures, if the composition dependence enters only through a set of weighted moments of the density $\{\xi_{k}\}$. Then (26) becomes

[TABLE]

For weighted densities consistent with an SPT (or FMT) approach, it has been shown Ref. Santos (2012) that the thermodynamic coefficients for mixtures are determined once an equation of state for the single-component system is known.

We can relate the radial derivative of $\Delta\Omega(R)$ to the solvent density at contact; by connecting this to the virial theorem we can obtain a new thermodynamic relation. Following Ref. Bryk et al. (2003) we take the normal derivative of $\Omega$ with respect to $R$, and noting that $\Delta\Omega(R)=\Omega(R)-\Omega_{\mathrm{hom}}$ gives

[TABLE]

where $\rho_{0}$ is the equilibrium density profile and $\phi_{\mathrm{ext}}$ is the external potential (i.e. the potential of the solute). In equilibrium $\Omega$ is minimised so

[TABLE]

and the first integral in (29) vanishes. As the solute is hard, the external potential and its derivative are zero everywhere except at a distance $\frac{\sigma}{2}$ from the surface where both $\rho_{0}$ and $\phi_{\mathrm{ext}}$ are discontinuous. We consider its Boltzmann weight, i.e.

[TABLE]

Taking the (distributional) derivative of both sides gives

[TABLE]

Inserting this expression into (29) and using the fact that $\rho(\boldsymbol{\mathbf{r}})e^{\beta\phi_{\mathrm{ext}}(\boldsymbol{\mathbf{r}})}$ is a continuous function (cf. Ref. Hansen and McDonald (2013)) gives the contact theorem

[TABLE]

and the contact density in this inhomogeneous system is $\rho(\sigma;\phi_{\mathrm{ext}}):=\rho^{(2)}(\sigma)/\rho=\rho\,g^{(2)}(\sigma)$ recalling the definition of $\rho^{(n)}$ for the homogeneous system (1) (cf. Ref. Percus (1962)), so we have

[TABLE]

So inserting the SPT ansatz (6) gives

[TABLE]

Inserting the virial theorem (15) into the right-hand side of (30) yields the final expression:

[TABLE]

This relation is satisfied by the coefficients (12), which is surprising given that it was obtained from a completely different thermodynamic route and the ansatz (6) is inexact. Nonetheless, this self-consistency is a testament to the effectiveness of SPT and related approaches.

Since the pressure is now a known input, the excess chemical potential can be determined by integrating (27) i.e.

[TABLE]

To keep the expressions simple we will not evaluate the chemical potential until the very end, but it should be recognised as a known variable wherever it appears. With the pressure fixed we have three free parameters in the theory $\{a_{0},a_{1},a_{2}\}$; we must thus choose three out of the five available thermodynamic relations in (24), (25) and (31) to satisfy. Therefore, we must lose consistency with two of these relations to obtain a more accurate theory for practical applications.

To set the correct energy scale we choose to fix $\Delta\Omega(R=0)$ and $\Delta\Omega(R=\sigma/2)$ through equations (24a) and (25) using the chemical potential determined above in (32). This in turn imposes the consistency of the osmotic pressure (26). For the final equation we choose to set the contact value of $g^{(2)}$ through (31) which better represents solutes of interest than the two relations for point geometries at $R=0$. Solving these three equations gives the generalised SPT coefficients

[TABLE]

It can be verified that inserting the Percus-Yevick equation of state (12d) into these expressions yields the previously obtained coefficients (12), as expected. Inserting the CS equation of state we obtain (13) which are identical to the coefficients derived from the White Bear II (WBII) free energy functional of Ref. Hansen-Goos and Roth (2006a), although this is only clear after transforming to the excluded geometry through the canonical relations (11). Remarkably, we have obtained these coefficients through a route completely different from their original derivation.

In Ref. Hansen-Goos and Roth (2006a) the coefficients were determined within FMT by taking the limit of a binary mixture where one component is infinitely dilute. Here we completely avoided FMT, in favour of geometrical arguments similar to the classical SPT approach outlined in the previous section. This suggests that this generalised scaled particle argument is built into the structure of the WBII functional of Ref. Hansen-Goos and Roth (2006a); this is not an obvious fact as the derivation of this functional did not explicitly involve these arguments. Rather, the WBII functional was constructed based on a novel extension of the CS equation to mixtures by requiring self-consistency of the pressure in (28) Hansen-Goos and Roth (2006b). We imposed this relation by setting the chemical potential in (25) using the chemical potential obtained via (27). It is unclear to us how our final choice of using (25) instead of one of the two relations at the origin, i.e. (24b) or (24c), is built into the WBII functional.

Acknowledgements.

We are grateful to Bob Evans who encouraged writing up the results of section III, which eventually became this manuscript. JFR, FT and CPR acknowledge the European Research Council under the FP7 / ERC Grant Agreement No. 617266 “NANOPRS”. CPR would like to acknowledge the Royal Society for financial support.