Skip to content

Theory: Ornstein-Zernike equation and closure relations

previous (Fourier, Fourier-Bessel, and Sine transformations)
theory index
next (pyOZ algorithm)

This section contains some information about the Ornstein-Zernike equation and closure relations supported in pyOZ. It doesn’t, however, contain any derivation of the equation – such a discussion would be very much beyond the scope of this simple introduction.

Ornstein-Zernike equation

Let’s first define two important correlation functions-total correlation function h(\mathbf{r}_{ij}) and pair correlation function g(\mathbf{r}_{ij}). The pair correlation function shows us the probability of finding a particle of type j in a certain distance of the particle i divided by the probability of such an arrangement in the case where both i and j behave as an ideal gas. In other words, these functions show the influence of particle i on a particle j in a distance defined by the vector \mathbf{r}_{ij}.

The relationship between h and g is very simple.

h(\mathbf{r_{ij}}) = g(\mathbf{r_{ij}})-1

From definition, pair correlation function can be written as

\displaystyle g(\mathbf{r}_{ij}) = e^{-\beta W(\mathbf{r}_{ij})}

where W is the total interaction potential (potential of mean force). The way how the total potential is treated depends on the approximation used (see HNC and PY below).

The influence of particle i on particle j can be split into two contributions-direct, defined by direct correlation function c(\mathbf{r_{ij}}), and indirect. The indirect influence can be described by a direct influence of particle i on particle n, which in turn directly and indirectly influences the particle j.

This arrangement is described by the Ornstein-Zernike equation

\displaystyle h(\mathbf{r}_{ij}) = c(\mathbf{r}_{ij}) + \sum\limits_{n=1}^N \rho_n \int\limits_{-\infty}^{\infty} c(\mathbf{r}_{in}) h(\mathbf{r}_{nj}) d\mathbf{r}_{in}

where the sum goes over all components present and ρ is the number density of a given component. This equation assumes that the origin is set at the position of the particle i. This is the reason why the infinitesimal volume element contains the combination in

Furthermore, the direct interaction of i with j is described using vector \mathbf{r}_{ij}. The direct interaction with n then corresponds to vector \mathbf{r}_{in}, and iteraction between n and j to \mathbf{r}_{nj}. Simple vector algebra then yields

\mathbf{r}_{ij} = \mathbf{r}_{in} + \mathbf{r}_{nj}
\mathbf{r}_{nj} = \mathbf{r}_{ij}-\mathbf{r}_{in}

And the Ornstein-Zernike equation can be written in the following form

\displaystyle h(\mathbf{r}_{ij}) = c(\mathbf{r}_{ij}) + \sum\limits_{n=1}^N \rho_n \int\limits_{-\infty}^{\infty} c_{in}(\mathbf{r}_{in}) h_{nj}(\mathbf{r}_{ij}-\mathbf{r}_{in}) d\mathbf{r}_{in}

The integral is exactly the definition of convolution, therefore,

\displaystyle h(\mathbf{r}_{ij}) = c(\mathbf{r}_{ij}) + \sum\limits_{n=1}^N \rho_n (c_{in}* h_{nj})(\mathbf{r}_{ij})

This equation is extremely complicated to solve. However, as described in the theory of Fourier transformation, Fourier transformation of this function replaces the convolution by a product of Fourier transforms of individual functions.

\displaystyle \^h(\mathbf{k}_{ij}) = \^c(\mathbf{k}_{ij}) + \sum_{n=1}^N \rho_n \^c(\mathbf{k}_{in}) \^h(\mathbf{k}_{nj})

‘Hat’-ed letters denote Fourier-transformed correlation functions. The equation in Fourier space is easier to solve (at least numerically).

Closure relation

To be really able to solve this problem, however, we need one more equation-the so called closure relation (recall that we have two unknowns-both correlation functions, but only one equation at the moment). There exists an (exact) general closure relation, that can be derived using the density functional theory of liquids and graph theory, and which reads (from now on we drop the vector notation and the subscripts, although they are still operational)

g(r) = e^{-\beta U(r) + h(r)-c(r) + b(r)} = e^{-\beta U(r) + \Gamma(r) + b(r)}

where β is the reciprocal of the Boltzmann factor kT, U(r) is the pair potential acting between the particles, Γ(r)= h(r)-c(r) is called the indirect correlation function (Gamma function), and b(r) is term arising from the graph theory (bridged graphs), which is not known exactly. The expressions for h(r) and c(r) follow:

h(r) = e^{-\beta U(r) + \Gamma(r) + b(r)}-1

c(r) = e^{-\beta U(r) + \Gamma(r) + b(r)}-1-\Gamma(r)

For the closure relation, several approximations exist. In the Percus-Yewick (PY) and Hypernetted Chain (HNC) closures, the bridged graphs are neglected and some further simplifications are done.

In the Percus-Yewick closure, we expand the exponential of the Γ(r) using the Taylor expansion and take just the first two terms. Then

\displaystyle c(r) = e^{-\beta U(r) + \Gamma(r)}-\Gamma(r)-1 = e^{-\beta U(r)}e^{\Gamma(r)}-\Gamma(r)-1 =\\= e^{-\beta U(r)}(1 + \Gamma(r))-\Gamma(r)-1 = e^{-\beta U(r)}(1+\Gamma(r))-g(r) + c(r)

After rearrangement, we obtain an expression for the pair correlation function

\displaystyle g(r) = e^{-\beta U(r)}(1+\Gamma(r))

The PY closure (approximation) works well for noncharged systems. For the Lennard-Jones potential, it even provides analytical solution of the Ornstein-Zernike equation.

For the Hypernetted Chain closure (HNC), just neglect of bridged graphs and some rearrangement of terms is necessary.

\displaystyle c(r) = e^{-\beta U(r) + \Gamma(r)}-\Gamma(r)-1 = e^{-\beta U(r) + \Gamma(r)}-h(r) + c(r)-1 =\\= e^{-\beta U(r) + \Gamma(r)}-g(r) + c(r)

\displaystyle g(r) = e^{-\beta U(r) + \Gamma(r)}

HNC closure performs better than PY (it’s closer to the exact general closure, after all), at least for charged systems.

It’s possible to obtain these expressions also in other way. Consider the following expression for the total correlation function

c(r) = g(r)-g_{\mathrm{indirect}}(r) = e^{-\beta W(r)}-e^{-\beta(W(r)-U(r)) U(r)}

W(r) is the total interaction in terms of the potential of mean force (PMF) and U(r) is the pair potential. The second term, therefore, corresponds to the direct interaction between the particles. The exponential term with the PMF is indeed equal to the pair correlation function, therefore, for the PY closure

\displaystyle c(r) = g(r)-e^{\beta U(r)} g(r) = h(r) + 1-e^{\beta U(r)} g(r)

g(r) = e^{-\beta U(r)} [h(r)-c(r) + 1] = e^{-\beta U(r)} [\Gamma(r) + 1]

The HNC closure uses the Taylor expansion of the second exponential term

\displaystyle c(r) = g(r)-1 + \beta W(r)-\beta U(r) = g(r)-1-\beta U(r)-\ln g(r) =\\= h(r)-\beta U(r)-\ln g(r)

\displaystyle \ln g(r) = h(r)-c(r)-\beta U(r)

and, therefore,

\displaystyle g(r) = e^{-\beta U(r) + h(r)-c(r)} = e^{-\beta U(r) + \Gamma}

As can be seen, these expressions for both closures agree with those derived previously.

previous (Fourier, Fourier-Bessel, and Sine transformations)
theory index
next (pyOZ algorithm)