text:
# A Method to Extract Potentials from the Temperature Dependence of Langmuir Constants for Clathrate-Hydrates

## 1 Introduction

A central mission of chemical physics is to determine intermolecular interactions from experimental phase equilibrium data. This is generally a very difficult task because macroscopic equilibrium constants reflect averaging over vast numbers of poorly understood microscopic degrees of freedom. Since the advent of the computer, empirically guided numerical fitting has become the standard method to obtain intermolecular potentials. An appealing and often overlooked alternative, however, is to solve inverse problems based on simple statistical mechanical models.

A well-known class of inverse problems is to determine the density of states $`n(E)`$ from the partition function $`Q(\beta )`$, or one of its derivatives, as a function of inverse temperature $`\beta =1/kT`$. For classical Maxwell-Boltzmann statistics, these quantities are simply related by an Laplace transform

$$Q(\beta )=_0^{\mathrm{}}e^{\beta \epsilon }n(\epsilon )𝑑\epsilon $$
(1)
which is straightforward to invert. Analogous quantum-mechanical problems have also been solved using Laplace (or Mellin) transform methods. For example, in Bose-Einstein statistics,

$$Q(\beta )=_0^{\mathrm{}}\frac{n(\epsilon )d\epsilon }{e^{\beta \epsilon }+1},$$
(2)
the phonon density of states can be obtained from the specific heat of a crystal , and the area distribution of a blackbody radiator can be obtained from its power spectrum . Similarly, in Fermi-Dirac statistics,

$$Q(\beta )=_0^{\mathrm{}}\frac{n(\epsilon )d\epsilon }{e^{\beta (\epsilon \epsilon _F)}1},$$
(3)
the band structure of a doped semiconductor can be obtained from the temperature dependence of its carrier density . In each of these examples, it is possible to extract the microscopic density of states from a temperature-dependent thermodynamic quantity because each is a function of just one variable.

In the more complicated case of chemical systems, thermodynamic quantities are related to classical configurational integrals,

$$Z(\beta )=_Ve^{\beta \mathrm{\Phi }(\stackrel{}{r})}𝑑\stackrel{}{r},$$
(4)
where the total energy $`\epsilon `$ in Eq. (1) is replaced by the intermolecular potential $`\mathrm{\Phi }(\stackrel{}{r})`$ , and the integral extends over the interaction volume $`V`$. Unfortunately, these integrals are multi-dimensional, which leads to under-determined inverse problems. Perhaps as a result, statistical inversion methods have apparently not been developed for intermolecular potentials in chemical systems, and instead empirical fitting has been used exclusively to describe phase equilibrium data.

In this article, we show that within the common spherical-cell approximation the intermolecular potential is completely determined by the temperature dependence of the Langmuir constant. In this case, the linear integral equation (1) is replaced by a non-linear equation of the form

$$Z(\beta )=4\pi _0^{\mathrm{}}e^{\beta w(r)}r^2𝑑r$$
(5)
where $`w(r)`$ is a spherically averaged “cell potential” . As shown below in the important case of Langmuir constants for clathrate-hydrates, this simple inverse problem can be solved exactly without resorting to numerical fitting schemes.

Before proceeding, we mention some related ideas in the recent literature of solid state physics. As with chemical systems, empirical fitting is also the standard approach to derive interatomic potentials for metals and semiconductors. Since the pioneering work of Carlsson, Gelatt, and Ehrenreich in 1980 , however, exact inversion methods have been developed to obtain potentials from cohesive energy curves . These theoretical advances, discussed briefly in Appendix A, have recently led to improvements in the modeling of silicon, beyond what had previously been obtained by empirical fitting alone . Inspired by such developments in solid state physics, here we seek similar insights into clathrate-hydrate intermolecular forces, albeit using a very different statistical mechanical formalism.

## 2  The Inverse Problem for Clathrate-Hydrates

### 2.1  The Statistical Theory of van der Waals and Platteeuw

Clathrate-hydrates exist throughout nature and are potentially very useful technological materials . For example, existing methane hydrates are believed to hold much more energy than any fossil fuel in use today. Carbon dioxide hydrates are being considered as effective materials for the sequestration and/or storage of CO2. In spite of their great importance, however, the theory of clathrate-hydrate phase behavior is not very well developed, still relying for the most part on the ad hoc empirical fitting of experimental data. Therefore, we have chosen to develop our statistical inversion method in the specific context of clathrate-hydrate chemistry.

Since being introduced in 1959, the statistical thermodynamical model of van der Waals and Platteeuw (vdWP) has been used almost exclusively to model the phase behavior of clathrate-hydrates, usually together with a spherical cell (SC) model for the interaction potential between the enclathrated or “guest” molecule and the cage of the clathrate-hydrate. The SC model was also introduced by vdWP, inspired by an analogous approximation made by Lennard-Jones and Devonshire in the case of liquids .

In the general formulation of vdWP , the chemical potential difference between an empty, unstable hydrate structure with no guest molecules, labeled MT, and the stable hydrate, labeled H, is related to the so-called Langmuir hydrate constant $`C_{Ji}`$ and the fugacity of the guest molecule $`\widehat{f}_J`$

$$\mathrm{\Delta }\mu ^{MTH}=kT\underset{i}{}\nu _i\mathrm{ln}(1+\underset{J}{}C_{Ji}\widehat{f}_J)$$
(6)
where $`i`$ designates the type of cage, $`\nu _i`$ the number of cages of type $`i`$ per water molecule and $`J`$ the type of guest molecule. In practice, experimental phase equilibria data is used to determine $`\mathrm{\Delta }\mu ^{MTH}`$.

The connection with intermolecular forces within vdWP theory is made by expressing the Langmuir hydrate constant as the configurational integral $`Z_{Ji}`$ divided by $`kT`$, which is written explicitly as an integral over the volume $`V`$

$$C_{Ji}(T)=\frac{1}{8\pi ^2kT}_Ve^{\mathrm{\Phi }(r,\theta ,\varphi ,\alpha ,\xi ,\gamma )/kT}r^2sin\theta sin\xi 𝑑r𝑑\theta 𝑑\varphi 𝑑\alpha 𝑑\xi 𝑑\gamma $$
(7)
where $`\mathrm{\Phi }(r,\theta ,\varphi ,\alpha ,\xi ,\gamma )`$ is the general six dimensional form of the interaction potential between the guest molecule at spherical coordinates $`(r,\theta ,\varphi )`$ oriented with Euler angles $`(\alpha ,\xi ,\gamma )`$ with respect to all of the water molecules in the clathrate-hydrate.

In the SC approximation, which is made without any careful mathematical justification, the intermolecular potential $`\mathrm{\Phi }`$ is replaced by a spherically averaged cell potential $`w(r)`$, which reduces the Langmuir constant formula (7) to a single, radial integration

$$C_{Ji}(T)=\frac{4\pi }{kT}_0^Re^{w(r)/kT}r^2𝑑r$$
(8)
where the cutoff distance $`R`$ is often arbitrarily taken as the radius of the cage. (As shown below, the exact value rarely matters because temperatures are typically so low that the high energy portion of the cage $`rR`$ makes a negligible contribution to the integral.) Although the SC approximation may appear to be a drastic simplification, it is nevertheless very useful for theoretical studies of intermolecular forces based on Langmuir constant measurements.

### 2.2  Numerical Fitting Schemes

Before this work, the functional form of the cell potential $`w(r)`$ has always been obtained by first choosing a model interaction potential between the guest molecule in a cage and each nearest neighbor water molecule essentially ad hoc, and then performing the spherical average from (7) to (8) analytically. The most common potential form in use today is the Kihara potential, which is simply a shifted Lennard-Jones potential with a hard-core. Using the Kihara potential and spherically averaging the interaction energy, typically over the first-shell only, yields the following functional form for $`w(r)`$:

$$w(r)=2zϵ[\frac{\sigma ^{12}}{R^{11}r}(\delta _{10}+\frac{a}{R}\delta _{11})\frac{\sigma ^6}{R^5r}(\delta _4+\frac{a}{R}\delta _5)]$$
(9)
where

$$\delta _N=\frac{1}{N}[(1\frac{r}{R}\frac{a}{R})^N(1+\frac{r}{R}\frac{a}{R})^N]$$
(10)
and $`z`$ is the coordination number, $`R`$ is the radius of the cage, and $`\sigma `$, $`ϵ`$, and $`a`$ are the Kihara parameters. As a result of the averaging process leading from (7) to (8), the functional form of $`w(r)`$ is fairly complicated, and the parameters $`ϵ`$ and $`\sigma `$ are generally determined by fitting monovarient equilibrium temperature-pressure data numerically .

There are several serious drawbacks to this ubiquitous numerical fitting procedure, which suggest that the Kihara parameters lack any physical significance: $`(i)`$ The Kihara parameters are not unique, and many different sets can fit the experimental data well; $`(ii)`$ the Kihara parameters found by fitting Langmuir curves do not match those of found by fitting other experimental data, such as the second virial coefficient or the gas viscosity ; and $`(iii)`$ comparisons of Langmuir constants found via the SC approximation (8) and via explicit multi-dimensional quadrature (7) show that the two can differ by over 12 orders of magnitude  (which results from the exponentially strong sensitivity of the Langmuir constant to changes in the cell potential). These problems call into question the validity of using the Kihara potential as the basis for the empirical fitting, and even the use of the SC approximation itself.

### 2.3  Inversion of Langmuir Curves

It would clearly be preferable to extract more reliable information about the interatomic forces in clathrate-hydrates directly from experimental data without any ad hoc assumptions about their functional form. In principle, such an approach is possible for clathrate-hydrates which contain a single type of guest molecule occupying only one type of cage. In this case, each of the sums in Eq. (6) contains only one term, and by using an equation of state to compute the fugacity $`\widehat{f}_J`$, the Langmuir constant $`C_{Ji}`$ can be determined directly from experimental phase equilibria data. Typical data sets obtained in this manner are shown in Fig. 1 for Structure I ethane and cyclopropane clathrate-hydrates .

Because the full potential $`\mathrm{\Phi }`$ in (7) is multi-dimensional (while the Langmuir constant only depends on a single parameter $`T`$), the general vdWP theory is too complex to pose a well-defined inverse problem for the interatomic forces. The SC approximation, on the other hand, introduces a very convenient theoretical construct, the spherically averaged potential $`w(r)`$, which has the same dimensionality as the Langmuir curve $`C_{Ji}(T)`$ of a single type of guest molecule occupying a single type of cage. (Since we consider only this case, we drop the subscripts $`Ji`$ hereafter.) Although one can question the accuracy of the SC approximation, its simplicity at least allows precise connections to be made between the Langmuir curve and the cell potential.

As an appealing alternative to empirical fitting, in this article we view Eq. (8) as an integral equation to be solved analytically for $`w(r)`$, given a particular Langmuir curve $`C(T)`$. Letting $`\beta =1/kT`$, we rewrite (8) as

$$C(\beta )=4\pi \beta _0^{\mathrm{}}e^{\beta w(r)}r^2𝑑r,$$
(11)
where we have also set the upper limit of integration to $`R=\mathrm{}`$, which introduces negligible errors due to the very low temperatures (large $`\beta `$) accessible in experiments. (This will be justified a posteriori with a precise definition of “low” temperatures below.) Note that Eq. (11) has the form of Eq. (5) with $`Z(\beta )=\beta C(\beta )`$.

### 2.4  Application to Experimental Data

In our analytical approach, some straight-forward fitting of the raw experimental data is needed to construct the function $`C(\beta )`$, but after that, the “inversion” process leading to $`w(r)`$ is exact. For example, typical sets of experimental data are well described by a van’t Hoff temperature dependence

$$C(\beta )=C_oe^{m\beta }$$
(12)
as shown in Fig. 2 for ethane and cyclopropane clathrate hydrates , and the constant $`m`$ is generally positive. (Note that the exponential dependence which we call ”van’t Hoff dependence” in this paper can be expected based on quite general thermodynamic considerations .)

Aided by the analysis, the quality and functional form of these fits are discussed below in section 7. In order to allow for deviations from the dominant van’t Hoff behavior, however, in this article we consider the more general form

$$C(\beta )=\beta F(\beta )e^{m\beta }$$
(13)
where $`m`$ is a constant defined by

$$m=\underset{\beta \mathrm{}}{lim}\mathrm{log}C(\beta )/\beta $$
(14)
whenever this limit exists and is finite, i.e. when the prefactor $`F(\beta )`$ in Eq. (13) is dominated by the exponential term at low temperatures. We exclude the possibility of hyper-exponential behavior at low temperatures, $`\mathrm{log}C/\beta \mathrm{}`$, which is not physically meaningful, as explained below. The set of possible prefactors includes power-laws, $`F(\beta )=\beta ^\mu `$, as well as various rational functions.

The rest of article is organized as follows. In section 3 we discuss various necessary and sufficient conditions for the existence of physically reasonable solutions, and we also derive the asymptotics of the Langmuir curve at low temperature from the behavior of the cell potential at near its minimum. In section 4, we perform the analysis in the general case, and in section 5, we discuss the specific case of van’t Hoff dependence (12), which leads to a cubic solution as well as various unphysical solutions involving cusps. In section 6, we derive analytical solutions for several different temperature dependences, which reveal the significance of possible deviations from van’t Hoff behavior for the form of the potential $`w(r)`$. The theoretical curves are compared with the experimental data in section 7, and the physical conclusions of the analysis are summarized in section 8. Relevant mathematical theorems are proved in the Appendix B.

## 3 General Analysis of the Inverse Problem

### 3.1 Necessary Conditions for the Existence of Solutions

On physical grounds, it expected that the cell potential $`w(r)`$ is continuous (at least piecewise) and has a finite minimum at $`r_o0`$ somewhere inside the clathrate cage, $`w(r)w(r_o)=w_o`$. We also allow the possibility that $`w(r)`$ is infinite for certain values of $`r`$ (e.g. outside a “hard-wall radius”) by simply omitted such values from the integration in Eq. (11). As proved in the Appendix B, these simple physical requirements suffice to imply the asymptotic relation (14), where $`m=w_o`$. They also place important constraints on the prefactor $`F(\beta )`$ defined in (13), which must be

* analytic in the half plane $`\text{Re}\beta >c`$ and
* real, positive and non-increasing for $`\beta >c`$ on the real axis

where $`c0`$ is a real number. (Note that we view the inverse temperature $`\beta `$ as a complex variable, for reasons soon to become clear.) Moreover, if the set

$$S=\{r0|w_o<w(r)<\mathrm{}\}$$
(15)
has nonzero measure, then $`F(\beta )`$ is strictly decreasing on the positive real axis. It is straightforward to generalize these rigorous results to the multi-dimensional integral of vdWP theory, Eq. (7), without making the spherical cell approximation, as described in the Appendix B, but hereafter we discuss only the spherically averaged integral equation, Eq. (11), because it makes possible an exact inversion.

In this section, we give simple arguments to explain the results proved in Theorem 1 of the Appendix B. First, we consider the illustrative example of a constant cell potential with a hard wall at $`r=r_{hw}>0`$,

$$w(r)=\{\begin{array}{cc}w_o\hfill & \text{ if }0r<r_{hw}\hfill \\ \mathrm{}\hfill & \text{ if }r>r_{hw}\hfill \end{array}$$
(16)
which satisfies the assumptions stated above. The integral (11) is easily performed in this case to yield

$$C(\beta )=\frac{4}{3}\pi r_{hw}^3\beta e^{w_o\beta }$$
(17)
which implies $`m=w_o`$, the well depth, and $`F(\beta )=\frac{4}{3}\pi r_{hw}^3`$, the volume of negative energy. Consistent with the general results above, $`F(\beta )`$ is constant in this case since $`S=\mathrm{}`$. For continuous potentials $`w(r)`$, however, the prefactor $`F(\beta )`$ must be strictly decreasing because $`S\mathrm{}`$.

The integral equation (11) can be simplified by a change of variables from radius to volume. In terms of a shifted cell potential versus volume

$$u(x)=w(r)w_o\text{ where }x=\frac{4\pi }{3}r^3$$
(18)
the integral equation is reduced to the form

$$C(\beta )=\beta F(\beta )e^{w(r_o)\beta }$$
(19)
where

$$F(\beta )=_0^{\mathrm{}}e^{\beta u(x)}𝑑x.$$
(20)
Since $`u(x)0`$ by construction, the function $`F(\beta )`$ is clearly non-increasing. In the Appendix B, it is proved that if the potential varies continuously near its minimum (in a very general sense), then $`F(\beta )`$ does not decay exponentially. Since $`F(\beta )`$ also positive and bounded above, we conclude

$$\underset{\beta \mathrm{}}{lim}\mathrm{log}C(\beta )/\beta =w_o$$
(21)
which implies $`m=w_o`$. Therefore, the slope of the Langmuir curve on a “van’t Hoff plot” ($`\mathrm{log}C`$ versus $`\beta =1/kT`$) in the low temperature limit is equal to (minus) the minimum energy of the cell potential. Since it is generally observed that $`m`$ is positive, as in the case of ethane and cyclopropane clathrate-hydrates shown in Figs. 12, the cell potential must be attractive, $`w_o=m<0`$, which simply indicates that the total internal energy is lowered by the introduction of guest molecules into clathrate-hydrates.

The fact that $`F(\beta )`$ must be nondecreasing has important consequences for the existence of solutions which are piecewise continuous and bounded below. For example, consider the class of Langmuir curves of the form

$$C(\beta )=\beta ^\nu e^{m\beta }$$
(22)
which is useful in fitting experimental data (see below). We have already addressed the borderline case, $`\nu =1`$, in which a discontinuous hard-wall solution is possible. Although it is not obvious a priori, there are no solutions to the inverse problem if $`\nu >1`$, since in that case $`F(\beta )=\beta ^{\nu 1}`$ is increasing. On the other hand, if $`\nu <1`$, then well-behaved continuous solutions are possible, because $`F(\beta )`$ is strictly decreasing.

### 3.2 Low Temperature Asymptotics of the Langmuir Curve

¿From Eq. (21), we see that the minimum energy $`w_o=w(r_o)`$ determines the leading order asymptotics of $`C(\beta )`$ in the low temperature limit. More generally, one would expect that $`C(\beta )`$ at low temperatures is completely determined by shape of the cell potential at low energies, close to its minimum. Using standard methods for the asymptotic expansion of Laplace integrals , it is straightforward to provide a mathematical basis for this intuition. For simplicity, here we consider the usual case of a parabolic minimum

$$w(r)=w_o+\frac{1}{2}k(rr_o)^2+O((rr_o)^3)$$
(23)
for some constants $`k>0`$ and $`r_o0`$, although below we will derive many exact solutions with non-parabolic minima.

Due to the factor of $`r^2`$ appearing in the integrand in Eq. (11), the two cases of a non-central or central minimum, $`r_o>0`$ and $`r_o=0`$, respectively, must be treated separately. Physically, this qualitative difference between central and non-central-wells is due to the spherical averaging process going from Eq. (7) to Eq. (11): A central well in $`w(r)`$ corresponds to a unique local minimum of the multidimensional potential $`\mathrm{\Phi }`$, but a non-central-well in $`w(r)`$ corresponds to a nonlocal minimum of $`\mathrm{\Phi }`$ which is smeared across a sphere of radius $`r_o`$.

Beginning with non-central-well case, $`r>0`$, we have the following asymptotics as $`\text{Re}\beta \mathrm{}`$:

$`C(\beta )`$ $``$ $`4\pi r_o^2\beta {\displaystyle _{r_oϵ}^{r_o+ϵ}}e^{\beta (w_o+\frac{1}{2}k(rr_o)^2)}𝑑r`$ (24)
$``$ $`4\pi r_o^2\left({\displaystyle \frac{2\beta }{k}}\right)^{1/2}e^{\beta w_o}{\displaystyle _{\mathrm{}}^{\mathrm{}}}e^{t^2}𝑑t`$
$`=`$ $`4\pi r_o^2\left({\displaystyle \frac{2\pi \beta }{k}}\right)^{1/2}e^{\beta w_o}`$
which is the usual leading order term in the expansion of a Laplace integral . Therefore, the experimental signature of a non-central-well is a Langmuir constant which behaves at low temperatures like

$$C(\beta )C_oe^{m\beta }\beta ^{1/2}\text{ as }\beta \mathrm{}.$$
(25)
Comparing (24) and (25), we can identify the well depth $`w_o=m`$, consistent with the general arguments above, but it is impossible to determine independently the location of $`r_o`$ and the curvature $`k`$ of the minimum. Instead, any $`r_o`$ and $`k`$ satisfying $`4\pi r_o^2\sqrt{2\pi /k}=C_o`$ would exactly reproduce the same large-$`\beta `$ asymptotics of the Langmuir curve (as would a completely different central-well solution described in section 6). This degeneracy of non-central-well solutions revealed in the low temperature asymptotics is actually characteristic of all non-central-well solutions, as explained below.

In the central-well case, $`r_o=0`$, the asymptotics must be carried out more carefully because the leading term derived in (24) vanishes:

$`C(\beta )`$ $``$ $`4\pi \beta {\displaystyle _0^ϵ}e^{\beta (w_o+\frac{1}{2}kr^2)}r^2𝑑r`$ (26)
$``$ $`4\pi \left({\displaystyle \frac{2}{k^3\beta }}\right)^{1/2}e^{\beta w_o}{\displaystyle _0^{\mathrm{}}}t^{1/2}e^t𝑑t`$
$`=`$ $`\left({\displaystyle \frac{2\pi }{k}}\right)^{3/2}{\displaystyle \frac{e^{\beta w_o}}{\beta ^{1/2}}}.`$
The experimental signature of a parabolic central well in the Langmuir curve,

$$C(\beta )C_oe^{m\beta }\beta ^{1/2}\text{ as }\beta \mathrm{},$$
(27)
is qualitatively different from (25), which provides an unambiguous way to separate the two cases using low temperature measurements. Moreover, unlike the non-central-well case, the curvature $`k=2\pi C_o^{2/3}`$ of a parabolic central minimum is uniquely determined by the low temperature asymptotics of the Langmuir curve. Consistent with asymptotic results, we shall see in section 4 that central-well solutions to the inverse problem are unique, while non-central-well solutions are not.

### 3.3 Sufficient Conditions for the Existence of Solutions

The primary difficulty in solving Eq. (11) lies in its being a nonlinear integral equation of the “first kind” for which no general theory of the existence and uniqueness of solutions exists . In the linear case, however, there is a special class of first-kind equations which can be solved using Laplace, Fourier or Mellin transforms, namely integral equations of the additive or multiplicative convolution type

$`\mathrm{\Psi }(x)`$ $`=`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}K(xy)\psi (y)𝑑y`$ (28)
$`\text{or }\mathrm{\Psi }(x)`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}K(xy)\psi (y)𝑑y,`$ (29)
respectively, where $`\psi (x)`$ is the unknown function and $`\mathrm{\Psi }(x)`$ is given. Integral equations of the multiplicative form (29) often arise in statistical mechanics as explained in the Introduction.

Although our nonlinear, first-kind equation (11) is not of the convolution type because the unknown function $`w(r)`$ appears in the exponent, it does somewhat resemble a Laplace transform. This connection is more obvious in the alternative formulation (20) relating $`F(\beta )`$ and $`u(x)`$, which is equivalent to the original equation (11) according to the analysis above. In the next section, it is shown that physically reasonable solutions exist whenever $`F(\beta )`$ has an inverse Laplace transform $`f(y)`$ which is positive, nondecreasing and non-constant for $`y>0`$. (Sufficient conditions on $`F(\beta )`$ to ensure these properties of $`f(y)`$ are given in Appendix B.)

## 4 Analytical Solutions for Arbitrary Langmuir Curves

### 4.1 The Unique Central-Well Solution

It is tempting to change variables $`y=u(x)`$ in the integral (20) to reduce it to a Laplace transform, but care must be taken since $`u^1(y)`$ may not be single-valued. This leads us to treat solutions which are monotonic separately from those from those which are not, an important distinction fore-shadowed by the asymptotic analysis above. As a natural first case, we seek differentiable solutions $`u(x)`$ which are strictly increasing without bound ($`u(\mathrm{})=\mathrm{}`$) from a central minimum ($`u(0)=0`$). Such “central-well solutions” correspond to cell potentials $`w(r)`$ which are strictly increasing from a finite minimum $`w(0)=w_o`$ at the center of the cage. We proceed by considering the inverse cell potential $`v(y)=u^1(y)`$ with units of volume as a function of energy, which is single-valued and strictly increasing with $`v(0)=0`$, as shown in Fig. 3(a). With the substitution $`y=u(x)`$, Eq. (20) is reduced to Laplace’s integral equation  for the unknown function $`v^{}(y)`$

$$F(\beta )=_0^{\mathrm{}}\frac{e^{\beta y}dy}{u^{}(u^1(y))}=_0^{\mathrm{}}e^{\beta y}v^{}(y)𝑑y.$$
(30)
Upon taking inverse Laplace transforms, we arrive at a differential equation for $`v(y)`$

$$v^{}(y)=f(y)$$
(31)
whose unique solution is

$$v(y)=_0^yf(y)𝑑y$$
(32)
using the boundary condition $`v(0)=0`$. According to (31), the continuity of $`f(y)`$ for $`y>0`$ (which is not assumed) would the guarantee differentiability of $`v(y)`$ for $`y>0`$, and hence of $`u(x)`$ for $`x>0`$.

Equivalently, we can also simplify (30) with an integration by parts

$$F(\beta )=\beta _0^{\mathrm{}}e^{\beta y}v(y)𝑑y.$$
(33)
Therefore, the inverse cell potential is given by

$$v(y)=g(y)$$
(34)
where $`g(y)`$ is the inverse Laplace transform of the function

$$G(\beta )=\frac{F(\beta )}{\beta }=\frac{C(\beta )e^{\beta w_o}}{\beta ^2}.$$
(35)
The cell potential $`u(x)`$ is determined implicitly by the algebraic equation

$$g(u)=x.$$
(36)
Returning to the original variables, we have a general expression for $`w(r)`$ in the central-well case

$$w(r)=w_o+g^1\left(\frac{4}{3}\pi r^3\right).$$
(37)
This equation uniquely determines the central-well potential that exactly reproduces any admissible Langmuir curve.

### 4.2 Non-Central-Well Solutions

The simplest kind of non-central-well solution is the central-well (37) shifted by a “hard-core” radius $`r_{hc}>0`$

$$w(r)=\{\begin{array}{cc}\mathrm{}\hfill & \text{ if }0r<r_{hc}\hfill \\ w_o+g^1\left[\frac{4}{3}\pi (r^3r_{hc}^3)\right]\hfill & \text{ if }rr_{hc}\hfill \end{array}$$
(38)
which exemplifies a peculiar general property of our integral equation: An arbitrary hard core can be added to any solution. Note that, if $`u(x)`$ is any solution of the rescaled equation (20), then so is

$$\stackrel{~}{u}(x)=\{\begin{array}{cc}\mathrm{}\hfill & \text{ if }0x<x_{hc}\hfill \\ u(xx_{hc})\hfill & \text{ if }xx_{hc}\hfill \end{array}$$
(39)
for any hard-core volume $`x_{hc}0`$. The proof is simple:

$$_0^{\mathrm{}}e^{\beta \stackrel{~}{u}(x)}𝑑x=_{x_{hc}}^{\mathrm{}}e^{\beta u(xx_{hc})}𝑑x=_0^{\mathrm{}}e^{\beta u(x)}𝑑x=F(\beta ).$$
(40)
Physically, a hard core for the cell potential could represent the presence of a second guest molecule (in a spherically symmetric model) in the same clathrate-hydrate cage. Alternatively, a hard-core could represent a water molecule (again in a spherically symmetric model) at the node of several adjacent clathrate cages, in which case the cell potential actually describes the “super-cage” surrounding the central water molecule.

The arbitrary hard-core just described only hints at the vast multiplicity of non-monotonic solutions to the integral equation (20), which is a common characteristic of first-kind equations . Next, we consider the general case of a non-central-well, shown in Fig. 3(b), which includes (38) as a special case. To be precise, we seek continuous solutions $`u(x)`$ on an interval $`x_1<x<x_2`$ composed of a non-increasing function $`u_{}(x)`$ and a nondecreasing function $`u_{}(x)`$ which are piecewise differentiable and non-negative. We also allow for a possible hard-core in the central region $`x<x_1`$ as well as a “hard wall” beyond the clathrate cage boundary $`x>x_2`$. The general form of such a non-central-well solution is

$$u(x)=\{\begin{array}{cc}\mathrm{}\hfill & \text{ if }0<x<x_1\hfill \\ u_{}(x)\hfill & \text{ if }x_1<xx_o\hfill \\ u_+(x)\hfill & \text{ if }x_ox<x_2\hfill \\ \mathrm{}\hfill & \text{ if }x<x_2\hfill \end{array}$$
(41)
where $`u_{}(x_o)=u_+(x_o)`$. We do not assume $`u_+^{}(x_o)=u_{}^{}(x_o)`$, which would imply differentiability at the minimum $`u^{}(x_o)=0`$, although we do not rule out this case either. Instead, we allow for the mathematical possibility of a discontinuous first derivative at $`x_o`$, i.e. a “cusp” at the well position, at least for the moment.

As before, it is convenient to express the solution (41) in terms of two differentiable functions $`v_{}(y)=u_{}^1(y)`$ and $`v_+(y)=u_+^1(y)`$ which describe the multi-valued inverse cell potential. Note that $`v_{}(\mathrm{})=x_1`$, $`v_+(\mathrm{})=x_2`$ and $`v_{}(0)=v_+(0)=x_o`$. In terms of the inverse cell potentials, the integral equation (20) takes the form

$`F(\beta )`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}e^{\beta u(x)}𝑑x`$ (42)
$`=`$ $`{\displaystyle _{x_1}^{x_o}}e^{\beta u_{}(x)}𝑑x+{\displaystyle _{x_o}^{x_2}}e^{\beta u_+(x)}𝑑x`$
$`=`$ $`{\displaystyle _0^{\mathrm{}}}e^{\beta y}\left[v_+^{}(y)v_{}^{}(y)\right]𝑑y`$
which implies

$$v_+^{}(y)v_{}^{}(y)=f(y).$$
(43)
In this case, the continuity of $`f(y)`$ would only guarantee the differentiability of the difference $`v_+(y)v_{}(y)`$, but not of the individual functions $`v_+(y)`$ and $`v_{}(y)`$. Integrating (42) by parts before taking the inverse transform yields a general expression for the solution

$$v_+(y)v_{}(y)=g(y)$$
(44)
where again $`g(y)`$ is the inverse Laplace transform of $`F(\beta )/\beta `$. Unfortunately, we have two unknown functions and only one equation, so the set of non-central-well solutions is infinite.

The scaled Langmuir curve $`F(\beta )`$ uniquely determines only $`v_+(y)v_{}(y)`$, the volume difference as a function of energy between the two branches of the scaled cell potential $`u(x)`$, but not the branches $`v_+(y)`$ and $`v_{}(y)`$ themselves. An infinite variety of non-central-well solutions, which exactly reproduce the same Langmuir curve as the central well solution, can be easily generated by choosing any non-increasing, non-negative, piecewise differentiable function $`v_{}(y)`$ such that the function $`v_+(y)`$ defined by (44) is nondecreasing. Even the position of the well $`v_{}(0)=x_o`$ can be chosen arbitrarily.

For example, one such family of solutions with a central “soft-core” ($`x_1=0`$) is given by

$$u(x)=\{\begin{array}{cc}u_{}(x)\hfill & \text{ if }0xx_o\hfill \\ u_+(x)\hfill & \text{ if }xx_o\hfill \end{array}$$
(45)
where

$$v_{}(y)=u_{}^1(y)=\{\begin{array}{cc}x_oay^b\hfill & \text{ if }0yy_c\hfill \\ 0\hfill & \text{ if }yy_c\hfill \end{array}$$
(46)
and

$$v_+(y)=u_+^1(y)=v_{}(y)+g(y),$$
(47)
for any $`a,b>0`$ and $`x_o0`$. (In the limit $`a0`$, we recover the unique central-well solution.) Note that $`y_c=u(0)=(x_o/a)^{1/b}`$ is the height of the central maximum of the potential. These solutions, all derived from a single Langmuir curve, exist whenever $`g(y)`$ increases quickly enough that $`v_+(y)`$ is nondecreasing, which is guaranteed if $`g^{}(y)aby^{b1}`$ for $`0<y<y_c`$.

Another family of non-central-well solutions with a soft-core can be constructed with the choice

$$v_{}(y)=\{\begin{array}{cc}x_o\lambda g(y)\hfill & \text{ if }0yy_c\hfill \\ 0\hfill & \text{ if }yy_c\hfill \end{array}$$
(48)
for any $`0<\lambda <1`$ and $`x_o0`$, where $`y_c=g^1(x_o/\lambda )`$. In this case, the cell potential is easily expressed in terms of $`g^1(x)`$ as

$$u(x)=\{\begin{array}{cc}g^1\left(\frac{x_ox}{\lambda }\right)\hfill & \text{ if }0xx_o\hfill \\ g^1\left(\frac{x_ox}{1\lambda }\right)\hfill & \text{ if }x_oxx_c\hfill \\ g^1(x)\hfill & \text{ if }xx_c\hfill \end{array}$$
(49)
where $`x_c=v_+(y_c)=x_o/\lambda `$. This class of solutions exists whenever $`g(y)`$ is nondecreasing (or $`f(y)0`$). If $`y_c=\mathrm{}`$, then there is a hard core $`u(0)=\mathrm{}`$. Otherwise, if there is a soft core $`u(0)=y_c<\mathrm{}`$, then there is typically a cusp (discontinuous derivative) at $`x_c`$, as explained below.

As demonstrated by the preceding examples, it is simple to generate an enormous variety of non-central-well solutions, with an arbitrarily shaped soft or hard core, and an arbitrary position of the minimum. In spite of the multiplicity of non-central-well solutions, however, our analysis of the inverse problem at least determines $`v_+(y)v_{}(y)`$ uniquely from any experimental Langmuir curve. This important analytical constraint is not satisfied by empirical fitting procedures.

### 4.3 Soft Cores and Outer Cusps

Non-central-well solutions with a soft-core satisfy $`v_{}(y)=0`$ for $`yy_c>0`$, as in the examples above. In such cases, $`v_+(y)=g(y)`$ for $`yy_c`$ regardless of whether or not there is an outer hard wall, which implies that $`u(x)=g^1(x)`$ for $`xx_c`$, where $`x_c=g^1(y_c)`$. If $`f(y)`$ is continuous for $`y>0`$, then, unless $`u_{}(x)`$ has an “inverted cusp” at the origin ($`v_{}^{}(y_c)=0`$ and $`u_{}^{}(0^+)=\mathrm{}`$), any non-central-well solution $`u(x)`$ with a soft-core must have a cusp at $`x=x_c`$, as in the examples above. This “outer cusp” in $`u(x)`$ could only be avoided if $`f(y)`$ itself has a cusp at $`y_c`$ which would allow $`v_+(y)`$ to be continuous. However, an inverted cusp in $`u(x)`$ at the origin does not necessarily imply an cusp in $`w(r)`$ at the origin due the transformation $`x=\frac{4}{3}\pi r^3`$. For example, if $`v_{}(y)(y_cy)^{3/2}`$ as $`yy_c`$, or $`u(x)y_cx^{2/3}`$ as $`x0`$, then $`w(r)`$ would have a physically reasonable, parabolic soft core $`w(r)w_o+y_c(4\pi /3)^{2/3}r^2`$ as $`r0`$. Nevertheless, even in such cases, if $`f(y)`$ were continuous for all $`y>0`$, then both $`u(x)`$ and $`w(r)`$ would have unphysical second-derivative discontinuities at $`x=x_c`$ related to the soft core. In general, a continuously differentiable, non-central-well solution with a central soft core could only arise if $`f(y)`$ were discontinuous at some $`y_c>0`$, and such discontinuities are generally not present.

### 4.4 Cusps at a Non-Central Minimum

As mentioned above, the behavior of the cell potential near its minimum (whether central or not) is determined by the behavior of the Langmuir curve at low temperature, or equivalently, at large inverse temperature, $`\beta =1/T`$. The Laplace transform formalism makes this connection transparent and mathematically rigorous. The asymptotic behavior of $`G(\beta )`$ as $`\text{Re}\beta \mathrm{}`$ is related to the asymptotics of the inverse transform $`g(y)`$ as $`y0`$, which in turn governs the local shape of the energy minimum through Eq. (34) for a central well or Eq. (44) for a non-central well. The leading order asymptotics has already been computed above for parabolic minima, but the general solutions above show how various non-local properties of the potential are related to finite temperature features of the Langmuir curve. Here, we comment on a subtle difference in differentiability between central and non-central-wells, related to the small-$`y`$ behavior of $`f(y)`$.

For typical sets of experimental data, including the van’t Hoff form (12), the prefactor $`F(\beta )`$ has a bounded inverse Laplace transform in the neighborhood of the origin

$$\underset{y0}{lim}f(y)=f(0)<\mathrm{}.$$
(50)
This generally implies the existence of a cusp at a non-central minimum of $`u(x)`$, which is signified by a nonzero right and/or left derivative. When $`u(x)`$ is differentiable at its minimum, it satisfies $`u_{}^{}(x_o^{})=u_+^{}(x_o^+)=0`$. In the central-well case $`x_o=0`$, the existence of a cusp in $`u(x)`$ follows from (31)

$$u^{}(0^+)=1/v^{}(0^+)=1/f(0)>0,$$
(51)
but this does not imply a cusp in the unscaled potential $`w(r)`$ as long as $`f(0)>0`$ because in that case

$$w(r)w_o=u(4\pi r^3/3)(4\pi /3f(0))r^3\text{ as }r0^+.$$
(52)
In the non-central-well case, however, the bounded inverse transform (50) implies a cusp at the minimum because, $`v_+^{}(0)v_{}^{}(0)=f(0)<\mathrm{}`$ from (44) along with $`v_+^{}(0)0`$ and $`v_{}^{}(0)0`$ implies that $`v_+^{}(0)<\mathrm{}`$ and/or $`v_{}^{}(0)>\mathrm{}`$ which in turn implies $`u_+^{}(0)>0`$ and/or $`u_{}^{}(0)<0`$. Unlike the central-well case, however, a cusp in $`u(x)`$ at the non-central minimum $`x_o>0`$ implies a cusp at the corresponding non-central minimum of $`w(r)`$. Therefore, we conclude that whenever (50) holds, the only physically reasonable solution is the central-well solution (37).

### 4.5 Asymptotics at High Energy and Temperature

The high energy behavior of the cell potential is related to (but not completely determined by) the high temperature asymptotics of the Langmuir hydrate constant, through the function $`g(y)`$. For example, the cell potential would have a hard wall at $`x_2<\mathrm{}`$, if and only if $`g(y)`$ were unbounded

$$\underset{y\mathrm{}}{lim}g(y)=\mathrm{}.$$
(53)
Since the empirical modeling of Langmuir curves using Kihara potentials assumes an outer hard-core, Eq. (53) could be used to test the suitability of using the Kihara potential form, although experimental data is often not available at sufficiently high temperatures to make a fully adequate comparison (see below). Whenever (53) holds, the non-central-well solutions $`u(x)`$ are also universally asymptotic to the central-well solution

$$u(x)g^1(x)$$
(54)
at large volumes $`x\mathrm{}`$. This follows from (44) and the fact that $`v_{}(y)`$ is bounded, which implies $`v_+(y)g(y)`$. The exact inversions performed in section 6 provide further insight into the relationship between small $`\beta `$ asymptotics of the Langmuir constant and high energy behavior of the cell potential.

## 5 Langmuir Curves with van’t Hoff Temperature Dependence

Experimental Langmuir hydrate-constant curves $`C(\beta )`$ are well fit by an ideal van’t Hoff temperature dependence (12), demonstrated by straight lines on Arrhenius log-linear plots

$$\mathrm{log}C=m\beta +\mathrm{log}C_o$$
(55)
as shown in Figs. 1 and 2 for ethane ($`C_o=4.733\times 10^7`$ atm-1, $`m=9.4236`$ kcal/mol) and cyclopropane ($`C_o=1.9041\times 10^7`$ atm-1, $`m=10.5939`$ kcal/mol) clathrate-hydrates . This data is analyzed carefully in section 7, where alternative functional forms are considered. In the ideal van’t Hoff case, we have $`F(\beta )=C_o/\beta `$ and $`G(\beta )=C_o/\beta ^2`$. The inverse Laplace transforms of these functions are simply $`f(y)=C_oH(y)`$ and $`g(y)=C_oyH(y)`$, respectively, where $`H(y)`$ is the Heaviside step function.

We begin by discussing the unique central-well solution, which is illustrated by the solid line in Fig. 4 for the case of ethane. The central-well solution is linear in volume $`u(x)=g(x)=C_oyH(y)`$, and cubic in radius

$$w(r)=\frac{4\pi r^3}{3C_o}m.$$
(56)
A curious feature of this exact solution is that it has a vanishing “elastic constant”, $`w^{\prime \prime }(0)=0`$, a somewhat unphysical property which we address again in section 7.

The simple form of (56) makes it very appealing as a means of interpreting experimental data with van’t Hoff temperature dependence. We have already noted that the slope of a van’t Hoff (Fig. 2) plot of the Langmuir constant is equal to the well depth $`m=w_o`$, but now we see that the $`y`$-intercept $`\mathrm{log}C_o`$ is related to the well-size, e.g. measured by the volume of negative energy $`mC_o`$. This volume corresponds to a spherical radius of

$$r_s=\left(\frac{3mC_o}{4\pi }\right)^{1/3}$$
(57)
which is $`0.4180`$ Å for ethane and $`0.3208`$ Å for cyclopropane. Since the van der Waals radius of ethane is less than that of cyclopropane, it makes physical sense that $`r_s^{\mathrm{ethane}}>r_s^{\mathrm{cyclopropropane}}`$. Moreover, these volumes fall within the ranges determined from two different experimental modeling approaches: Using the radius of the water cage from x-ray scattering experiments , Lennard-Jones potentials from gas viscosity data give 0.79 Å for ethane and 0.61 Å for cyclopropane , while computations with van der Waals radii give 0.18 Å for ethane and 0.03 Å for cyclopropane .

There are infinitely many non-central-well solutions reproducing van’t Hoff temperature dependence, but each of them has unphysical cusps (discontinuous derivatives). There will always be a cusp at the minimum of the potential, since $`f(y)`$ satisfies the general condition (50). For example, the central-well solution can be shifted by an arbitrary hard-core radius $`r_o0`$

$$w(r)=\{\begin{array}{cc}\mathrm{}\hfill & \text{ if }0r<r_o\hfill \\ \frac{4\pi (r^3r_o^3)}{3C_o}m\hfill & \text{ if }rr_o\hfill \end{array}$$
(58)
In the case of a soft core, there must be a second cusp in the outer branch of the potential at the same energy as the inner core due to the continuity of $`f(y)`$, as explained above. This is illustrated by the following piecewise cubic family of soft-core solutions of the general form (49):

$$w(r)=\{\begin{array}{cc}\frac{8\pi |r_o^3r^3|}{3C_o}m\hfill & \text{ if }0r2^{1/3}r_o\hfill \\ \frac{4\pi r^3}{3C_o}m\hfill & \text{ if }r_o2^{1/3}r_o\hfill \end{array}$$
(59)
which are shown in Fig. 4 in the case of ethane guest molecules. An infinite variety of other piecewise differentiable solutions exactly reproducing van’t Hoff dependence of the Langmuir curve could easily be generated, as described above, but each would have unphysical cusps.

Previous studies involving ad hoc fitting of Kihara potentials have reported non-central-wells , but these empirical fits may be only approximating various exact, cusp-like, non-central-well solutions, such as those described above. Moreover, given that the central-well solution (56) can perfectly reproduce the experimental data, it is clear that the results obtained by fitting Kihara potentials to Langmuir curves are simply artifacts of the ad hoc functional form, without any physical significance. Kihara fits also assume a hard wall at the boundary of the clathrate cage (by construction), whereas all of the exact analytical solutions (both central and non-central-wells) have the asymptotic dependence

$$w(r)\frac{4\pi r^3}{3C_o}$$
(60)
as $`r\mathrm{}`$ according to (54). Any deviation from the cubic shape at large radii, such as a hard wall, would be indicated by a deviation from van’t Hoff behavior at high temperatures, but such data would be difficult to attain in experiments (see below).

The preceding analysis shows that the only physical information contained in a Langmuir curve with van’t Hoff temperature dependence is the depth $`w_o`$ and the effective radius $`r_s`$ of the spherically averaged cell potential, which takes the unique form (56) in the central-well case. In hindsight, the simple two-parameter form of the potential is not surprising since a van’t Hoff dependence is described by only two parameters, $`m`$ and $`C_o`$. It is clearly inappropriate to fit more complicated ad hoc functional forms, such as Eq. (9) derived from the Kihara potential, since they contain extraneous fitting parameters and do not reproduce the precise shape of any exact solution.

## 6 Analysis of Possible Deviations from van’t Hoff Behavior

### 6.1 Dimensionless Formulation

The general analysis above makes it possible to predict analytically the significance of possible deviations from van’t Hoff temperature dependence, which could be present in the experimental data (see below). We have already discussed the experimental signatures of various low and high energy features of the cell potential in the asymptotics of the Langmuir curve. In this section, we derive exact solutions for Langmuir curves of the form (13) where $`F(\beta )`$ is a rational function. Such cases correspond to logarithmic corrections of linear behavior on a van’t Hoff plot of the Langmuir curve, which are small enough over the accessible temperature range to be of experimental relevance, in spite of the dominant van’t Hoff behavior seen in the data.

Fitting to the dominant van’t Hoff behavior (55) introduces natural scales for energy, $`m`$, and pressure, $`C_o^1`$, so it is convenient and enlightening to introduce dimensionless variables. With the definitions

$$\stackrel{~}{\beta }=m\beta ,\stackrel{~}{C}(\stackrel{~}{\beta })=C(\stackrel{~}{\beta }/m)/C_o,\text{ and }\stackrel{~}{F}(\stackrel{~}{\beta })=F(\stackrel{~}{\beta }/m)/mC_o,$$
(61)
the Langmuir curve can be expressed in the dimensionless form

$$\stackrel{~}{C}(\stackrel{~}{\beta })=\stackrel{~}{\beta }\stackrel{~}{F}(\stackrel{~}{\beta })e^{\stackrel{~}{\beta }}.$$
(62)
For consistency with these definitions, the other energy-related functions in the analysis are nondimensionalized as follows

$$\stackrel{~}{G}(\stackrel{~}{\beta })=G(\stackrel{~}{\beta }/m)/m^2C_o,\stackrel{~}{y}=y/m,\stackrel{~}{f}(\stackrel{~}{y})=f(m\stackrel{~}{y})/C_o,\stackrel{~}{g}(\stackrel{~}{y})=g(m\stackrel{~}{y})/mC_o,$$
(63)
where $`\stackrel{~}{f}(\stackrel{~}{y})`$ and $`\stackrel{~}{g}(\stackrel{~}{y})`$ are the inverse Laplace transforms of $`\stackrel{~}{F}(\stackrel{~}{\beta })`$ and $`\stackrel{~}{G}(\stackrel{~}{\beta })`$, respectively. The natural scales for energy and pressure also imply natural scales for volume, $`mC_o`$, and distance, $`r_s`$, as described in the previous section, which motivates the following definitions of the dimensionless cell potential versus volume

$$\stackrel{~}{x}=x/mC_o,\stackrel{~}{u}(\stackrel{~}{x})=u(mC_o\stackrel{~}{x})/m$$
(64)
and radius

$$\stackrel{~}{r}=r/r_s,\stackrel{~}{w}(\stackrel{~}{r})=w(r_s\stackrel{~}{r})/m.$$
(65)
Note that $`\stackrel{~}{x}=\stackrel{~}{r}^3`$. With these definitions, the central-well solution takes the simple form,

$$\stackrel{~}{u}(\stackrel{~}{x})=\stackrel{~}{g}^1(\stackrel{~}{x})$$
(66)
in terms of the dimensionless volume, or

$$\stackrel{~}{w}(\stackrel{~}{r})=1+\stackrel{~}{g}^1(\stackrel{~}{r}^3)$$
(67)
in terms of the dimensionless radius. We now consider various prefactors $`\stackrel{~}{F}(\stackrel{~}{\beta })`$ which encode valuable information about the energy landscape in various regions of the clathrate cage.

### 6.2 The Interior of the Clathrate Cage

#### 6.2.1 Power-Law Prefactors

The simplest possible correction to van’t Hoff behavior involves a power-law prefactor

$$\stackrel{~}{F}(\stackrel{~}{\beta })=\stackrel{~}{\beta }^\mu \text{ for any }\mu >0,$$
(68)
which corresponds to a logarithmic correction on a van’t Hoff plot of the Langmuir constant

$$\mathrm{log}\stackrel{~}{C}=\stackrel{~}{\beta }+(1\mu )\mathrm{log}(\stackrel{~}{\beta })$$
(69)
as shown in Fig. 5(a). In this case, we have

$$\stackrel{~}{f}(\stackrel{~}{y})=\stackrel{~}{y}^{\mu 1}H(\stackrel{~}{y})/\mathrm{\Gamma }(\mu )\text{ and }\stackrel{~}{g}(\stackrel{~}{y})=\stackrel{~}{y}^\mu H(\stackrel{~}{y})/\mathrm{\Gamma }(\mu +1),$$
(70)
where $`\mathrm{\Gamma }(z)`$ is the gamma function. In general, power-law prefactors at low temperatures signify an energy minimum with a simple polynomial shape.

#### 6.2.2  The Central-Well Solution

The unique central-well solution is also a simple power law

$$\stackrel{~}{u}(\stackrel{~}{x})=\left[\mathrm{\Gamma }(\mu +1)\stackrel{~}{x}\right]^{1/\mu }$$
(71)
or equivalently

$$\stackrel{~}{w}(\stackrel{~}{r})=1+\mathrm{\Gamma }(\mu +1)^{1/\mu }\stackrel{~}{r}^{3/\mu }.$$
(72)
The cubic van’t Hoff behavior is recovered in the case $`\mu =1`$, as is the (asymptotic) parabolic behavior from (23) and (27) in the case $`\mu =3/2`$. Because $`\stackrel{~}{w}(\stackrel{~}{r})+1\stackrel{~}{r}^{3/\mu }`$, a power-law correction to van’t Hoff behavior with a positive exponent ($`\mu <1`$) corresponds one which is “wider” than a cubic, while a negative exponent ($`\mu >1`$) corresponds to a potential which is “more narrow” than a cubic, as shown in Fig. 5(b). On physical grounds, the smooth polynomial behavior described by (72) is always to be expected near the minimum energy of the cell potential. Therefore, the power-law correction to van’t Hoff behavior (69) has greatest relevance for low temperature measurements in the range $`\stackrel{~}{\beta }1`$, from which it determines interatomic forces in the interior of the clathrate cage at low energies $`|\stackrel{~}{w}(\stackrel{~}{r})|1`$.

#### 6.2.3  Non-Central-Well Solutions

As described above, there are infinitely many non-central-well solutions. One family of solutions of the form (49) with $`\lambda =1/2`$ is given by

$$\stackrel{~}{w}(\stackrel{~}{r})+1=\{\begin{array}{cc}\left[2\mathrm{\Gamma }(\mu +1)|\stackrel{~}{r}^3\stackrel{~}{r}_o^3|\right]^{1/\mu }\hfill & \text{ if }0\stackrel{~}{r}2^{1/3}\stackrel{~}{r}_o\hfill \\ \left[\mathrm{\Gamma }(\mu +1)\stackrel{~}{r}^3\right]^{1/\mu }\hfill & \text{ if }\stackrel{~}{r}2^{1/3}\stackrel{~}{r}_o\hfill \end{array}$$
(73)
where $`\stackrel{~}{r}_o=r_o/r_s`$ is arbitrary, as shown in Fig. 5(c) for the case $`\stackrel{~}{r}_o=0.65`$. These solutions are unphysical since they all have cusps at $`\stackrel{~}{r}=2^{1/3}\stackrel{~}{r}_o`$ near the outer wall of the cage. However, they can still have reasonable behavior near the minimum at $`\stackrel{~}{r}_o`$ for certain values of $`\mu `$, which could have experimental relevance for low temperature measurements. Near the minimum, the exact solutions (73) have the asymptotic form

$$\stackrel{~}{w}(\stackrel{~}{r})1+\left[6\mathrm{\Gamma }(\mu +1)\stackrel{~}{r}_o^2|\stackrel{~}{r}\stackrel{~}{r}_o|\right]^{1/\mu }\text{ as }\stackrel{~}{r}\stackrel{~}{r}_o,$$
(74)
which is cusp-like for $`\mu >1/2`$, but differentiable for $`0<\mu 1/2`$. For example, the non-central-well has a parabolic shape in the case $`\mu =1/2`$, which agrees with the asymptotic analysis in Eqs. (23)–(25) when the units are restored, and it has a cubic shape when $`\mu =1/3`$. On the other hand, in the central-well case $`\mu =3/2`$ and $`\mu =1`$ correspond to parabolic and cubic minima, respectively. Therefore, this example nicely illustrates the difference between the low-energy asymptotics of central and non-central-wells described above in section 3, which would be useful in interpreting any experimental Langmuir constant data showing deviations from van’t Hoff behavior.

### 6.3 The Outer Wall of the Clathrate Cage

#### 6.3.1  Rational Function Prefactors

The behavior of the Langmuir curve in the high temperature region $`\stackrel{~}{\beta }=O(1)`$ is directly linked to properties of the outer wall of the clathrate cage, described by the cell potential at high energies $`\stackrel{~}{w}(\stackrel{~}{r})+1=O(1)`$. Although this region of the Langmuir curve does not appear to be accessible in experiments (see below), in this section we derive exact solutions possessing different kinds of outer walls, whose faint signature might someday be observed in experiments at moderate temperatures. In order to isolate possible effects of the outer wall, we consider Langmuir curves which are exactly asymptotic to the usual van’t Hoff behavior at low temperatures with small logarithmic corrections (on a van’t Hoff plot) at moderate temperatures. These constraints suggest choosing rational functions for $`\stackrel{~}{F}(\stackrel{~}{\beta })`$ such that $`\stackrel{~}{F}(\stackrel{~}{\beta })1/\stackrel{~}{\beta }`$ as $`\beta \mathrm{}`$.

#### 6.3.2 Central Wells with Hard Walls

We begin by considering a “shifted power-law” prefactor

$$\stackrel{~}{F}(\stackrel{~}{\beta })=1/(\stackrel{~}{\beta }+\alpha ),\text{ for any }\alpha >0$$
(75)
which corresponds to a shifted logarithmic deviation from van’t Hoff behavior,

$$\mathrm{log}\stackrel{~}{C}=\stackrel{~}{\beta }\mathrm{log}(1+\alpha /\stackrel{~}{\beta }).$$
(76)
As shown in Fig. 6(a), this suppresses the Langmuir constant at high temperatures, which intuitively should be connected with an enhancement of the strength of the outer wall compared to the cubic van’t Hoff solution. Taking inverse Laplace transforms we have

$$\stackrel{~}{f}(\stackrel{~}{y})=e^{\alpha \stackrel{~}{y}}H(\stackrel{~}{y})\text{ and }\stackrel{~}{g}(\stackrel{~}{y})=(1e^{\alpha \stackrel{~}{y}})H(\stackrel{~}{y})/\alpha ,$$
(77)
and indeed, since $`\stackrel{~}{g}(\stackrel{~}{y})`$ is bounded, all solutions must have a hard wall regardless of whether or not the well is central, as described above. For example, the unique central-well solution is

$$\stackrel{~}{w}(\stackrel{~}{r})=1\mathrm{log}(1\alpha \stackrel{~}{r}^3)/\alpha \text{ for }0\stackrel{~}{r}<\alpha ^{1/3}$$
(78)
which has an outer hard wall at $`\stackrel{~}{r}=\alpha ^{1/3}`$, as shown in Fig. 6(b). The solution is also asymptotic to the cubic van’t Hoff solution at small radii $`\stackrel{~}{r}\alpha ^{1/3}`$. Therefore, in the limit $`\alpha 0`$, the radius of the outer hard wall diverges, and the solution reduces to the cubic shape as the deviation from van’t Hoff behavior is moved to increasingly large temperatures. Since empirical fitting with Kihara potential forms arbitrarily assumes an outer hard wall, this example provides analytical insight into the nature of the approximation at moderate to high temperatures, where the Langmuir constant should be suppressed according to (76).

#### 6.3.3  Central Wells with Soft Walls

Next we consider the opposite case of a Langmuir constant which is enhanced at high temperatures compared to van’t Hoff behavior, which intuitively should indicate the presence of a “soft wall”, rising much less steeply than a cubic function. An convenient choice is

$$\stackrel{~}{F}(\stackrel{~}{\beta })=\stackrel{~}{\beta }/(\stackrel{~}{\beta }^2\gamma ^2),\text{ for any }\gamma >0.$$
(79)
which is analytic except for poles at $`\beta =\pm \gamma `$ on the real axis. Although this function diverges at $`\beta =\gamma `$ due to the overly soft outer wall, the corresponding Langmuir curve

$$\mathrm{log}\stackrel{~}{C}=\stackrel{~}{\beta }\mathrm{log}\left[1(\gamma /\stackrel{~}{\beta })^2\right]$$
(80)
shown in Fig. 6(a) could have experimental relevance at moderate temperatures $`\stackrel{~}{\beta }\gamma `$, if $`\gamma `$ were sufficiently small. In this case, we have

$$\stackrel{~}{f}(\stackrel{~}{y})=\mathrm{cosh}(\gamma \stackrel{~}{y})H(\stackrel{~}{y})\text{ and }\stackrel{~}{g}(\stackrel{~}{y})=\mathrm{sinh}(\gamma \stackrel{~}{y})H(\stackrel{~}{y})/\gamma $$
(81)
which yields the central-well solution

$$\stackrel{~}{w}(\stackrel{~}{r})=1+\mathrm{sinh}^1(\gamma \stackrel{~}{r}^3)/\gamma .$$
(82)
As shown in Fig. 6(b), this function follows the van’t Hoff cubic at small radii $`\stackrel{~}{r}\gamma ^{1/3}`$ but “softens” to a logarithmic dependence for large radii $`\stackrel{~}{r}\gamma ^{1/3}`$.

## 7 Interpretation of Experimental Data

We begin by fitting Langmuir curves, computed from experimental phase equilibria data, an equation of state, and reference thermodynamic properties  for ethane and cyclopropane clathrate-hydrates to the van’t Hoff equation

$$\mathrm{log}C=m\beta +b$$
(83)
using least-squares linear regression. This leads to rather accurate results, as indicated by the small uncertainties in the parameters displayed in Table 1 (63% confidence intervals corresponding to much less than one percent error). The high quality of the regression of $`\mathrm{log}C`$ on $`\beta `$ is further indicated by correlation coefficients very close to unity, $`0.99650`$ and $`0.99998`$ for the ethane and cyclopropane data, respectively. Using the fitted values for $`m`$ and $`C_o=e^b`$, the data for the two clathrate-hydrates can be combined into a single plot in terms of the dimensionless variables $`\stackrel{~}{C}`$ and $`\stackrel{~}{\beta }`$, as shown in in Fig. 7, which further demonstrates the common linear dependence.

Converting the experimental data to dimensionless variables also reveals that the measurements correspond to extremely “low temperatures”. This is indicated by large values of $`\stackrel{~}{\beta }=m/kT`$ in the range of 16 to 24, which imply that $`kT`$ is less than $`6\%`$ of the well depth $`m`$. As such, physical intuition tells us that the experiments can probe the cell potential only very close to its minimum. This intuition is firmly supported by the asymptotic analysis above, which (converted to dimensionless variables) links the asymptotics of the Langmuir constant for $`\stackrel{~}{\beta }1`$ to that of the cell potential for $`|\stackrel{~}{r}\stackrel{~}{r}_o|1`$. In this light, it is clear that any features of the cell potential other than the local shape of its minimum, which are determined by empirical fitting, e.g. using Eq. (9) based on the Kihara potential, are simply artifacts of an ad hoc functional form, devoid of any physical significance.

Since the shape of the potential very close to its minimum should always be well approximated by a polynomial (the leading term in its Taylor expansion), the analysis above implies that only simple power-law prefactors to van’t Hoff behavior should be considered in fitting low temperature data. Therefore, we refit the experimental data, allowing for a logarithmic correction,

$$\mathrm{log}C=m\beta +b+\nu \mathrm{log}(\beta )$$
(84)
as in Eq. (69). The results are shown in Table 1, and the best-fit functions are displayed in dimensionless form in Fig. 7. In the case of ethane, the best-fit value of $`\mu =1\nu `$ corresponds to a roughly linear central-well solution $`\stackrel{~}{w}\stackrel{~}{r}^{0.9}`$ or a cusp-like non-central-well solution $`\stackrel{~}{w}|\stackrel{~}{t}\stackrel{~}{r}_o|^{0.3}`$. Although these solutions are not physically reasonable, perhaps the qualitative increase in $`\mu `$ compared with ideal van’t Hoff behavior ($`\mu =1`$) is indicative of a parabolic central well ($`\mu =3/2`$). In the case of cyclopropane, we have $`\mu =1.4\pm 0.9`$, which violates the general condition $`\mu 0`$ needed for the existence of solutions to the inverse problem. If this fit were deemed reliable, then the basic postulate of vdWP theory, Eq. (7), would be directly contradicted, with or without the spherical cell approximation (see the Appendix B). It is perhaps more likely that the trend of decreasing $`\mu <1`$ could indicate a non-central parabolic minimum in the spherically averaged cell potential ($`\mu =1/2`$).

Although it appears there may be systematic deviations from ideal van’t Hoff behavior in the experimental data for ethane and cyclopropane, $`\nu 0`$ or $`\mu 1`$, the results are statistically ambiguous. For both types of guest molecules, adding the third degree of freedom $`\nu `$ substantially degrades the accuracy of the two linear parameters $`m`$ and $`b`$, with errors increased by several hundred percent. Moreover, the uncertainty in $`\nu `$ is comparable to its best-fit value. Therefore, it seems that we cannot trust the results with $`\nu 0`$, and, by the principle of Occam’s razor, we are left with the more parsimonious two-parameter fit to van’t Hoff behavior, which after all is quite good, and its associated simple cubic, central-well solution.

On the other hand, there are different two-parameter fits, motivated by the inversion theory, which can describe the experimental data equally well, but which are somewhat more appealing than the cubic solution in that they possess a non-vanishing elastic constant (second spatial derivative of the energy). For example, the fits can be done using (84) with the parameter $`\nu `$ fixed at either $`1/2`$ or $`1/2`$, corresponding to either a non-central or central, parabolic minimum, respectively. The results shown in Table 1 reveal that these physically significant changes in the functional form have little effect on the van’t Hoff parameters $`m`$ and $`b=\mathrm{log}C_o`$.

The difficulty with the experimental data as a starting point for inversion is its limited range in $`\stackrel{~}{\beta }`$ of roughly one decade, which makes it nearly impossible to detect corrections proportional to $`\mathrm{log}\beta `$ related to different polynomial shapes of the minimum. It would be very useful to extend the range of the data, using the analytical predictions to interpret the results. In general, it is notoriously difficult to determine power-law prefactors multiplying a dominant exponential dependence, but at least the present analysis provides important guidance regarding the appropriate fitting functions, which could not be obtained by ad hoc numerical fitting. Moreover, the clear physical meaning of the dominant van’t Hoff parameters elucidated by the analysis also makes them much more suitable to describe experimental data than the artificial Kihara potential parameters.

## 8 Summary

In this article, we have shown that spherically averaged intermolecular potentials can be determined analytically from the temperature dependence of Langmuir constants. Starting from the statistical theory of van der Waals and Platteeuw, the method has been developed for the case of clathrate-hydrates which contain a single type of guest molecule occupying a single type of cage. Finally, the method has been applied to experimental data for ethane and cyclopropane clathrate-hydrates. Various conclusions of the analysis are summarized below.

General Theoretical Conclusions

* Physically reasonable intermolecular potentials (which are piecewise continuous and bounded below) exist only if the Langmuir curve has a dominant exponential (van’t Hoff) dependence at low temperatures, $`lim_\beta \mathrm{}\mathrm{log}C/\beta =m`$, with a prefactor $`F(\beta )=C(\beta )e^{m\beta }/\beta `$ which is smooth and non-increasing.
* The slope $`m`$ of an experimental “van’t Hoff plot” of $`\mathrm{log}C`$ versus inverse temperature $`\beta `$ is precisely equal to the well depth, i.e. (minus) the minimum of the potential. This is true not only for the spherically averaged cell potential, $`\mathrm{min}w(r)=m`$, but also for the exact multi-dimensional potential, $`\mathrm{min}\mathrm{\Phi }(\stackrel{}{r})=m`$.
* For any physically reasonable Langmuir curve, the unique central-well potential can be determined from Eq. (37).
* There also exist infinitely many non-central-well solutions of the form (41), constrained only to satisfy Eq. (44). Several classes of such solutions with a central “soft core” (a finite maximum at the center of the cage) are described explicitly in Eqs. (45)–(49).
* Each one of the multitude of non-central-well solutions with a soft-core typically possesses unphysical cusps (slope discontinuities), while the unique central-well solution is a well-behaved analytic function.
* For ideal van’t Hoff temperature dependence, $`C(\beta )=C_oe^{m\beta }`$, the central-well solution is a simple cubic given by Eq. (56). The attractive region of the potential has depth $`m`$, volume $`mC_o`$, and radius $`r_s=(3mC_o/4\pi )^{1/3}`$. Each non-central-well solution for van’t Hoff dependence has two unphysical cusps, one at the minimum.
* The experimental signature of a parabolic, non-central-well is a Langmuir curve that behaves like $`C(\beta )C_oe^{m\beta }\beta ^{1/2}`$ at low temperatures ($`\beta \mathrm{}`$), while a parabolic central well corresponds to $`C(\beta )C_oe^{m\beta }\beta ^{1/2}`$.
* If there is a pure power-law prefactor multiplying van’t Hoff behavior $`C(\beta )=C_o(m\beta )^{1\mu }e^{m\beta }`$ with $`\mu >0`$, the central-well solution is also a power-law (72). For certain values of the prefactor exponent $`0<\mu 1/2`$, there are also non-central-well solutions with differentiable minima such as (73), although such solutions still possess cusps at higher energies.
* Rational function prefactors multiplying van’t Hoff behavior, such as (76) or (80), are associated with non-cubic behavior at the outer wall of the cage, such as a “hard wall” (78) or a “soft wall” (82), respectively.

Conclusions for Clathrate-Hydrates

* Since Langmuir constants must increase with temperature, $`m>0`$, on the basis of general thermodynamical arguments, the intermolecular potential must be attractive (with a region of negative energy).
* The depth $`w_o`$ and radius $`r_s`$ of the attractive region of the cell potential can be estimated directly from experimental data using the simple formulae $`w_o=m`$ and $`r_s=(3mC_o/4\pi )^{1/3}`$ without any numerical fitting. The resulting values for ethane and cyclopropane hydrates are consistent with typical estimates obtained by other means.
* The experimental Langmuir constant data for ethane and cyclopropane clathrate-hydrates is very well fit by an ideal van’t Hoff dependence, which corresponds to a cubic central well,

$$w(r)=\frac{4\pi r^3}{3C_o}m$$

as given by Eq. (56). However, the data is also equally consistent with a central parabolic well,

$$w(r)=\frac{\pi r^2}{C_o^{2/3}}m,$$

or various non-central (spherically averaged) parabolic wells. The range of temperatures is insufficient to distinguish between these cases.
* Experimental data tends to be taken at very “low” temperatures, $`kTm`$, which means that only the region of the potential very close to the minimum $`|rr_o|r_s`$ is probed. Therefore, only simple polynomial functions are to be expected, and fitting to more complicated functional forms, such as the Kihara potential, has little physical significance.
* In practical applications to clathrate-hydrates, the full power of our analysis could be exploited by measuring Langmuir hydrate constants over a broader range of temperatures than has previously been done.
* The availability of the inversion method obviates the need for empirical fitting procedures , at least for single-component hydrates in which guest molecules occupy only one type of cage. Moreover, the method also allows a systematic analysis of empirical functional forms, such as the Kihara potential, which cannot be expected to have much predictive power beyond the data sets used in parameter fitting.
* The general method of “exact inversion” developed here could also be applied to other multi-phase chemical systems, including guest-molecule adsorption at solid surfaces or in bulk liquids.

## Acknowledgments

We would like to thank Z. Cao for help with the experimental figures, J. W. Tester for comments on the manuscript, and H. Cheng for useful discussions. This work was supported in part by the Idaho National Engineering and Environmental Laboratory.

## Appendix A: Inversion of Cohesive Energy Curves for Solids

The basic idea of obtaining interatomic potentials by “exact inversion” has also recently been pursued in solid state physics (albeit based on a very different mathematical formalism having nothing to do with statistical mechanics). The inversion approach was pioneered by Carlsson, Gelatt and Ehrenreich in 1980 in the case of pair potentials for crystalline metals . These authors had the following insight: Assuming that the total (zero-temperature) cohesive energy $`E(x)`$ of a crystal with nearest neighbor distance $`x`$ can be expressed as a lattice sum over all pairs of atoms $`(i,j)`$

$$E(x)=\underset{ij}{}\varphi (xs_{ij})$$
(85)
where $`s_{ij}`$ are normalized atomic separation distances, then a unique pair potential $`\varphi (r)`$ can be derived which exactly reproduces the cohesive energy curve $`E(x)`$. (Lattice sums also appear in some clathrate-hydrates models , but to our knowledge they have never been used as the basis for an inversion procedure.)

The mathematical theory for the inversion of cohesive energy curves has been developed considerably in recent years and applied to wide variety of solids . The extension of the inversion formalism to semiconductors has required solving a nonlinear generalization of Eq. (85) representing many-body angle-dependent interactions :

$$F(x)=\underset{ijk}{}g(xs_{ij})g(xs_{ik})h(\theta _{ijk})$$
(86)
where $`F(x)`$ is the many-body energy, $`h(\theta _{ijk})`$ is the energy of the angle between two covalent bonds $`\stackrel{}{r}_{ij}=x\stackrel{}{s}_{ij}`$ and $`\stackrel{}{r}_{ik}=x\stackrel{}{s}_{ik}`$, and $`g(r)`$ is a radial function which sets the range of the interaction. In general, $`\varphi (r)`$, $`g(r)`$, and $`h(\theta )`$ can be systematically obtained from a set of multiple cohesive energy curves for the same material . The angular interaction can also be obtained directly from cohesive energy curves for non-isotopic strains .

## Appendix B: Mathematical Theorems

The first theorem provides necessary conditions on the Langmuir $`C(\beta )`$ so that the cell potential $`w(r)`$ is bounded below and continuous. It also interprets the slope of a van’t Hoff plot of the Langmuir curve in the low temperature limit as the well depth, under very general conditions. As pointed out in the main text, it is convenient to view the inverse temperature $`\beta `$ as a complex variable.

###### Theorem 1

Let $`w(r)`$ be real and continuous (except at possibly a finite number of discontinuities) for $`r0`$ with a finite minimum, $`w(r)w_o=w(r_o)>\mathrm{}`$ for some $`r_o0`$, and suppose that the integral

$$C(\beta )=4\pi \beta _0^{\mathrm{}}e^{\beta w(r)}r^2𝑑r$$
(87)
converges for some $`\beta =c`$ on the real axis. Then

$$C(\beta )=\beta F(\beta )e^{w_o\beta }$$
(88)
where the complex function $`F(\beta )`$ is

* real, positive and non-increasing on the real axis for $`\beta >c`$ and
* analytic in the half plane $`\text{Re}\beta >c`$.

If, in addition, the set $`S_ϵ=\{r0|w_o<w(r)<w_o+ϵ\}`$ has nonzero measure for some $`ϵ=ϵ_o>0`$, then $`F(\beta )`$ is strictly decreasing on the positive real axis (for $`\beta >c`$). Moreover, if $`S_ϵ`$ has finite, nonzero measure for every $`0<ϵ<ϵ_o`$, then

$$\underset{\beta \mathrm{}}{lim}\mathrm{log}C(\beta )/\beta =w_o$$
(89)
where the limit is taken on the real axis.

Proof: Define a shifted cell potential versus volume, $`u\left(\frac{4\pi }{3}r^3\right)=w(r)w_o`$. Substituting $`u(x)`$ for $`w(r)`$ reduces Eq. (87) to Eq. (88), where

$$F(\beta )=_0^{\mathrm{}}e^{\beta u(x)}𝑑x.$$
(90)
Since $`u(x)0`$ is real, the function $`F(\beta )`$ is real and positive for all real $`\beta `$ for which the integral converges. Moreover, for any complex $`\beta `$ and $`\beta ^{}`$ with $`\text{Re}\beta >\text{Re}\beta ^{}>c`$, we have the bound

$$|F(\beta )|_0^{\mathrm{}}e^{\text{Re}\beta u(x)}𝑑x_0^{\mathrm{}}e^{\text{Re}\beta ^{}u(x)}𝑑xF(c)<\mathrm{}$$
(91)
which establishes that the defining integral (90) converges in the right half plane $`\text{Re}\beta c`$ and is non-increasing on the real axis, thus completing the proof of (i).

Next let $`w(r)`$ be larger than its minimum value (but finite), $`w(r_o)<w(r)<\mathrm{}`$, on a set $`S_{\mathrm{}}`$ of nonzero measure, so that $`0<u(x)<\mathrm{}`$ for the corresponding set of volumes. Then for every $`\beta >\beta ^{}>c`$ on the real axis we have

$$_S_{\mathrm{}}e^{\beta u(x)}𝑑x<_S_{\mathrm{}}e^{\beta ^{}u(x)}𝑑x.$$
(92)
On the complement $`S_{\mathrm{}}^c=(0,\mathrm{})S_{\mathrm{}}`$, either $`u(x)=0`$ or $`u(x)=\mathrm{}`$, which implies

$$_{S_{\mathrm{}}^c}e^{\beta u(x)}𝑑x=_{S_{\mathrm{}}^c}e^{\beta ^{}u(x)}𝑑x.$$
(93)
¿From Eqs. (92)–(93) we conclude that $`F(\beta )`$ is strictly decreasing on the real axis.

Next we establish the low-temperature limit (89). Given $`0<ϵ<ϵ_o`$, we have the following lower bound for any $`\beta >c`$ on the real axis:

$$e^{ϵ\beta }F(\beta )=_0^{\mathrm{}}e^{\beta [u(x)ϵ]}𝑑xe^{ϵ\beta /2}_{S_{ϵ/2}}𝑑x+_{S_{ϵ/2}^c}e^{\beta [u(x)ϵ]}𝑑xe^{ϵ\beta /2}_{S_{ϵ/2}}𝑑x.$$
(94)
Combining this with the upper bound, $`F(\beta )F(c)<\mathrm{}`$, we obtain

$$e^{ϵ/2}M_ϵF(\beta )F(c)$$
(95)
where $`M_ϵ=_{S_{ϵ/2}}𝑑x`$ is a finite, nonzero constant (because $`S_{ϵ/2}`$ is assumed to have finite, nonzero measure). Substituting Eq. (88) in Eq. (95), we arrive at

$$\mathrm{log}\beta +\mathrm{log}M_ϵϵ\beta /2\beta w_o\mathrm{log}C(\beta )\mathrm{log}\beta +\mathrm{log}F(c)\beta w_o$$
(96)
which yields

$$w_oϵ/2\underset{\beta 0}{lim}\mathrm{log}C(\beta )/\beta w_o.$$
(97)
The desired result is obtained in the limit $`ϵ0`$.

Finally, we establish the analyticity of $`F(\beta )`$ in the open half plane $`\text{Re}\beta >c`$ by showing that its derivative exists and is given explicitly by

$$F^{}(\beta )=_0^{\mathrm{}}e^{\beta u(x)}u(x)𝑑x.$$
(98)
This requires justifying the passing a derivative inside the integral (20), which we have just shown to converge for $`\text{Re}\beta c`$. Using a classical theorem of analysis , it suffices to show that the integral in (98) converges uniformly for $`\text{Re}\beta >c+ϵ`$ for every $`ϵ>0`$ because the integrand is a continuous function of $`\beta `$ and $`x`$. (The possibility of a finite number of discontinuities in $`u(x)`$ is easily handled by expressing (98) as finite sum of integrals with continuous integrands.) It is a simple calculus exercise to show that $`te^t<1/e`$, and hence

$$te^{(c+ϵ)t}\frac{e^{ct}}{eϵ}$$
(99)
for all real $`t0`$. This allows us to derive a bound on the “tail” of the integral (98):

$$|_X^{\mathrm{}}e^{\beta u(x)}u(x)𝑑x|_X^{\mathrm{}}e^{\text{Re}\beta u(x)}u(x)𝑑x_X^{\mathrm{}}e^{(c+ϵ)u(x)}u(x)𝑑x<\frac{1}{eϵ}_X^{\mathrm{}}e^{cu(x)}𝑑x$$
(100)
which is independent of $`\beta `$. This uniform bound vanishes in the limit $`X\mathrm{}`$ because it is proportional to the tail of the convergent integral defining $`F(\beta )`$, which completes the proof. $`\mathrm{}`$

The proof of Theorem 1 does not depend in any way on the dimensionality of the integral and thus can be trivially extended to the general multi-dimensional case of vdWP theory without the spherical cell approximation.

###### Theorem 2

Let $`\mathrm{\Phi }(r,\theta ,\varphi ,\alpha ,\xi ,\gamma )\mathrm{\Phi }(r_o,\theta _o,\varphi _o,\alpha _o,\xi _o,\gamma _o)=w_o`$ be real and continuous, and suppose that the integral

$$C(\beta )=\frac{\beta }{8\pi ^2}_Ve^{\beta \mathrm{\Phi }(r,\theta ,\varphi ,\alpha ,\xi ,\gamma )}r^2sin\theta sin\xi 𝑑r𝑑\theta 𝑑\varphi 𝑑\alpha 𝑑\xi 𝑑\gamma $$
(101)
converges for some $`\beta =c`$ (real). Then all the conclusions of Theorem 1 hold.

The six-dimensional integral (101) of Theorem 2 does not present a well-posed inverse problem for the intermolecular potential $`\mathrm{\Phi }`$. However, the spherically averaged integral equation (87) of Theorem 1 can be solved for the cell potential $`w(r)`$ for a broad class of Langmuir curves $`C(\beta )`$ specified in the following theorem. The proof is spread throughout section 4 of the main text.

###### Theorem 3

If the inverse Laplace transform $`f(y)`$ of $`F(\beta )`$ exists and is nondecreasing and non-constant for $`y>0`$, then there exist a unique central-well solution ($`r_o=0`$) and infinitely many non-central-well solutions ($`0<r_o<\mathrm{}`$) to the inverse problem (11). If $`f(y)`$ is also continuous, then the central-well solution is the only continuously differentiable solution.

Finally, we state sufficient assumptions on $`F(\beta )`$ to guarantee the assumed properties of $`f(y)`$. In light of the necessary condition that $`F(\beta )`$ be analytic the right half plane $`\text{Re}\beta >c`$, the defining contour integral for $`f(y)`$

$$f(y)=\frac{1}{2\pi i}_{c^{}i\mathrm{}}^{c^{}+i\mathrm{}}e^{\beta y}F(\beta )𝑑\beta $$
(102)
must converge for any $`c^{}>c`$. By closing the contour in the left half plane, it can be shown that a sufficient (but not necessary) condition to ensure the assumed properties of $`f(y)`$ is that $`F(\beta )`$ decay in the left half plane ($`lim_\rho \mathrm{}|F(\rho e^{i\theta })|=0`$ for $`\pi /2\theta 3\pi /2`$) and have isolated singularities only on the negative real axis or at the origin with positive real residues. The particular examples of $`F(\beta )`$ considered in section 6 satisfy these conditions, but the weaker assumptions above regarding $`f(y)`$ suffice for the general derivation in section 4.


format: text

text:
# Enhanced interactions, mergers and starbursts in dense cluster environments

## 1. Introduction: sub-cluster merging and cluster disk galaxy evolution

Both optical studies of sky-projected cluster galaxy distributions (e.g. Baier 1979; Geller & Beers 1982) and studies of cluster X-ray morphologies (e.g. Mohr, Fabricant & Geller 1993; Buote & Tsai 1996) reveal the presence of significant substructure in clusters of galaxies. Subgroup velocity statistics (e.g. Dressler & Shectman 1988; Bird 1994; Girardi et al. 1997; Solanes, Salvador-Solé & González-Casado 1999) confirm these results, and are in substantial agreement in showing that some 30–40% of clusters have statistically significant substructure.

Numerical simulations also show that subclustering is a typical property of cluster formation in hierarchical theories of structure formation. Such subclustering is expected for a broad range of cosmologies, and its predicted extent is in agreement with that observed (e.g. Knebe & Müller 2000). Its prevalance indicates that clusters are young objects, since such substructure is expected to be smoothed out over a few cluster crossing times. Moreover the existing substructure indicates that recent major mergers have occurred in many clusters.

The evolution of cluster substructure is potentially of great importance for its influence on the star formation histories and morphological evolution of cluster galaxies (e.g Tully & Shaya 1984; Lavery & Henry 1988; Gnedin 1999; Bekki 1999). However, until very recently this potential influence has largely been neglected (cf. Bekki 1999). Of particular interest, is its relevance to the evolution of disk galaxies in clusters.

As is well known, Butcher & Oemler (1978) found that the cores of intermediate redshift clusters contain a higher fraction of blue, star-forming galaxies than similar environments at the present. These galaxies have been shown to be normal spiral or irregular galaxies, a fraction of which are obviously interacting or disturbed (e.g. Dressler et al. 1994; Oemler, Dressler & Butcher 1997; Smail et al. 1997). In rich clusters they constitute up to 50% of the population, but by the present epoch have been depleted by a factor of two and replaced by a corresponding increase in S0s (e.g. Dressler et al. 1997). Wang & Ulmer (1997) have found a correlation between the fraction of blue galaxies in intermediate redshift clusters and the cluster X-ray ellipticities. They interpret the latter as an age sequence, and conclude that the younger clusters have a higher fraction of blue, star-forming galaxies. As the authors note, this interpretation is in accord with a hierarchical clustering scenario, in which a typical distant cluster is assembled from smaller units which tend to contain more gas-rich, late type galaxies. The relative depletion of such star-forming galaxies in older, more relaxed clusters, suggests that the process of sub-cluster merging plays a crucial role in the transformation of an initial population of cluster spirals to a subsequent S0 population.

What effect is the sub-cluster merging process likely to have on a population of gas-rich, star-forming disk galaxies? Bekki (1999) using numerical simulation of a merger between a small group of galaxies and a cluster, shows that the time-dependent tidal gravitational field of the merger gives strong nonaxisymmetric perturbations to the disk galaxies in the group, inducing simultaneous moderate secondary starbursts in the central regions of the galaxies.

Similarly, Gnedin (1999) using self-consistent cluster simulations, demonstrates that for substantial infall of additional material into a cluster, the time-varying gravitational potential causes a sequence of tidal shocks on individual galaxies over a wide area of the cluster. These shocks enhance galaxy–galaxy interactions as well as amplifying galaxy merger rates. Again the effect is also likely to induce simultaneous starbursts in gas-rich infalling disk galaxies over a wide area of the cluster. Moreover, the tidal shocks are shown to typically thicken the disks of large spirals by a factor of two, making it unlikely that spiral structure and gaseous shocks in the inner regions of the galaxies will form. A significant amount of the dark matter halo is stripped, possibly removing the resevoir of fresh gas which maintains the star formation activity (cf. Larson, Tinsley & Caldwell 1980). While gas in the outer regions of the disk can be stripped by the ram pressure of the ambient medium, interpenetrating collisions with neighbouring galaxies could remove all their gas. The effect of all these tidal transformations is the conversion of spirals to S0s.

## 2. Observational evidence for tidally-induced starbursts in cluster galaxies

The theoretical considerations discussed above suggest that (moderate) tidally-induced starbursts should be widespread in cluster gas-rich disk galaxies in clusters which show evidence of on-going or recent merger activity. Together with enhanced starburst activity, enhanced galaxy–galaxy interactions and mergers are also expected. In fact, there is considerable observational evidence to support this scenario, both from studies of early-type and late-type galaxies in nearby clusters, and from evidence of poststarburst galaxies in intermediate redshift clusters. This evidence is now briefly reviewed.

Caldwell et al. (1993) have surveyed early-type galaxies in the Coma cluster and found an unusually high number with spectra that reflect recent enhanced star formation activity in a substructure of the cluster. As noted by Bekki (1999), this may readily be explained as due principally to tidal gravitational effects of a group-cluster merger with rather large relative velocity. Caldwell & Rose (1997) also found a larger fraction of poststarburst galaxies in nearby clusters with obvious double structure. In a similar fashion, Bekki notes that these results can also be explained by the tidal effects of cluster mergers.

Moss & Whittle (2000) have undertaken an H$`\alpha `$ survey of a magnitude-limited complete sample of 320 galaxies of types Sa and later within 1.5$`h^1`$ Mpc of the centres of 8 nearby Abell clusters. Some 116 of the sample were detected in emission. A subset of the emission-line galaxies are found to have tidally-induced circumnuclear starburst emission, as evidenced from the particular H$`\alpha `$ morphology, and its strong correlation with a disturbed galaxy morphology. Starburst emission is found to be most prevalent in the richest clusters in the survey. The percentage of the total sample of late-type galaxies with starburst emission is $``$ 40% for the two richest clusters, Coma and Abell 1367. In contrast, the corresponding percentage for field galaxies in the sample is only a few percent. Moreover X-ray morphologies and temperature structures for Coma and Abell 1367 indicate that both these clusters are recent post-merger systems (cf. Donnelly et al. 1998; Honda et al. 1996). Again, a plausible explanation for the enhanced starburst activity seen in the clusters, is tidal effects of subcluster merging.

Moss & Whittle also identify a class of late-type galaxies in their sample which are likely to be on-going mergers, in which the companion is indistinguishable from its merger partner. These on-going mergers are most prevalent in the Coma cluster. They are likely to represent a later stage of close double, interacting systems with tidally induced emission which comprise a significant fraction ($``$ 14%) of the galaxies they detect with starburst emission.

For intermediate redshift clusters, there is also evidence of an enhancement of starburst emission. Galaxies with spectra showing strong Balmer lines in absorption, indicative of a poststarburst phase, are an order of magnitude more frequent in the cluster environment as compared to the high redshift field (e.g. Dressler et al. 1999; Poggianti et al. 1999). A surprising fraction of the star-forming galaxies at intermediate redshift show signs of morphological disturbance, of which many are suggestive of merging systems, virtually always involving a disk galaxy (cf. Ellis 1999). This suggests a tidal origin for the starburst phase, although no clear difference in the incidence of disturbed morphology between the high redshift field and cluster environments has been established (cf. Dressler et al. 1999).

## 3. Conclusion

Strong theoretical and observational data suggest that sub-cluster merging plays a crucial role in the transformation of spirals to S0s in the process of cluster formation, due to the strong tidal forces on the disk galaxy members of an accreting sub-group. Theoretical studies predict simultaneous moderate starbursts in disk galaxies of the accreting sub-group due to tidal forces. This prediction is supported by observational data for both nearby and intermediate redshift clusters. Tidal shocks due to the time-varying potential associated with sub-cluster merging are expected to effect morphological transformation of star-forming disk galaxies to S0s. Thus these can explain the transformation of the majority of the spiral population in younger, less relaxed clusters to a population of S0s in older, more relaxed clusters as an integral part of the process of cluster formation.

## References

Baier, F. 1979, Astron. Nachr., 300, 85

Bekki, K. 1999, ApJ, 510, L15

Bird, C.M. 1994, ApJ, 107, 1637

Buote, D.A. & Tsai, J.C. 1996, ApJ, 458, 27

Butcher, H. & Oemler, A. 1978, ApJ, 219, 18

Caldwell, N. & Rose, J. 1997, AJ, 113, 492

Caldwell, N., Rose, J., Sharples, R.M, Ellis, R.S. & Bower, R.G. 1993, AJ, 106, 473

Donnelly, R.H. et al. 1998, ApJ, 500, 138

Dressler, A. & Shectman, S.A. 1988, AJ, 95, 985

Dressler, A., Oemler, A., Butcher, H.R. & Gunn, J.E. 1994, ApJ, 430, 107

Dressler, A. et al. 1997, ApJ, 490, 577

Dressler, A. et al. 1999, ApJS, 122, 51

Ellis, R.S. 1999, in IAU Symp. 183, The Formation and Evolution of Galaxies, ed. K. Sato, (Dordrecht: Boston), 107

Geller, M.J., & Beers, T.C. 1982, PASP, 94, 421

Girardi, M. et al. 1997, ApJ, 482, 41

Honda, H. et al. 1996, ApJ, 473, L71

Gnedin, O.Y. 1999, Ph.D. Thesis, Princeton Univ.

Knebe, A. & Müller, V. 2000, 354, 761

Larson, R.B., Tinsley, B.M. & Caldwell, C.N. 1980, ApJ, 237, 692

Lavery, R.J. & Henry, J.P. 1988, ApJ, 330, 596

Mohr, J.J., Fabricant, D.G. & Geller, M.J. 1993, ApJ, 413, 492

Moss, C., & Whittle, M. 2000, MNRAS, 317, 667

Oemler, A., Dressler, A. & Butcher, H.R. 1997, ApJ, 474, 561

Poggianti et al. 1999, ApJ, 518, 576

Smail, I., Dressler, A., Couch, W.J., Ellis, R.S., Oemler, A., Butcher H. & Sharples, R.M. 1997, ApJS, 110, 213

Solanes, J.M., Salvador-Solé, E. & González-Casado, G. 1999, A&A, 343, 733

Tully, R.B. & Shaya, E.J. 1984, ApJ, 281, 31

Wang, Q.D. & Ulmer, M.P. 1997, MNRAS, 292, 920