III.A General Geometric Characterization Theories

A general geometric characterization of porous media should satisfy the following requirements:

It should be well defined in terms of geometric quantities.
It should involve only parameters which are directly observable or measurable in an experiment independent of the phenomenon of interest.
It should not require the specification of too many parameters. The required independent experiments should be simple and economical to carry out. What is economical depends on the available data processing technology. With current data processing technology a characterization requiring more than $10^{{12}}$ numbers must be considered uneconomical.
The characterization should be usable in exact or approximate solutions of the equations of motion governing the phenomenon of interest.

The following sections discuss methods based on porosities $\overline{\phi}$ , correlation functions $S_{n}({\bf r}),C_{n}({\bf r})$ , local porosity distributions $\mu(\phi)$ , pore size distributions $\Pi(r)$ and capacities $T(\mathbb{K})$ . Table I collects the advantages and disadvantages of these methods according to the specified criteria.

Table I: Advantages and disadvantages of different geometric characterization methods for porous media. A geometrical characterization is called economical if it requires the specification of less than $10^{{12}}$ numbers. $\overline{\phi},...$ stands for porosity and other numbers defined in section III.A.1. The correlation functions $S_{n}({\bf r})$ and $C_{n}({\bf r})$ are discussed in section III.A.2. The pore size distribution of mercury porosimetry $\Pi(r)$ is defined in section III.A.3. The local porosity distributions $\mu(\phi)$ are discussed in section III.A.5. The capacities $T(\mathbb{K})$ are defined in III.A.6.

Characterization	well defined	predictive	economical	easily usable
$\overline{\phi},...$	yes	yes	yes	yes
$S_{2}({\bf r}),...$	yes	yes	yes	yes
$S_{n}({\bf r}),...(n\geq 3)$	yes	yes	no	yes
$\Pi(r)$	no	no	yes	yes
$\mu(\phi)$	yes	yes	yes	yes
$T(\mathbb{K})$	yes	no	no	no

III.A.1 Porosity and Other Numbers

III.A.1.a Porosity

The porosity of a porous medium is its most important geometrical property. Most physical properties are influenced by the porosity.

The porosity $\phi(\mathbb{S})$ of a two component porous medium $\mathbb{S}=\mathbb{P}\cup\mathbb{M}$ consisting of a pore space $\mathbb{P}$ (component one) and a matrix space $\mathbb{M}$ (component two) is defined as the ratio

$\phi(\mathbb{S})=\frac{V_{3}(\mathbb{P})}{V_{3}(\mathbb{S})}$

(3.1)

which gives the volume fraction of pore space. Here $V_{3}(\mathbb{P})$ denotes the volume of the pore space defined in (2.3) and $V_{3}(\mathbb{S})$ is the total sample volume. In the following the shorthand notation $V(\mathbb{G})=V_{3}(\mathbb{G})$ will often be employed.

The definition (3.1) is readily extended to stochastic porous media. In that case $V(\mathbb{P})$ and $\phi$ are random variables. If the medium is stationary then one finds using (2.3) and (2.11)

$\displaystyle\left\langle\phi\right\rangle$	$\displaystyle=$	$\displaystyle\frac{\left\langle V(\mathbb{P})\right\rangle}{V(\mathbb{S})}$	(3.2)
	$\displaystyle=$	$\displaystyle\frac{1}{V(\mathbb{S})}\left\langle\int _{\mathbb{S}}\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r})d^{3}{\bf r}\right\rangle$
	$\displaystyle=$	$\displaystyle\frac{1}{V(\mathbb{S})}\int _{\mathbb{S}}\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r})\right\rangle d^{3}{\bf r}$
	$\displaystyle=$	$\displaystyle\frac{1}{V(\mathbb{S})}\int _{\mathbb{S}}\mbox{Pr}\{{\bf r}\in\mathbb{P}\} d^{3}{\bf r}$
	$\displaystyle=$	$\displaystyle\mbox{Pr}\{{\bf r}_{0}\in\mathbb{P}\}=\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{0})\right\rangle$

where the last line holds only if the medium is stationary. ${\bf r}_{0}$ in the last line is an arbitrary point. Although the use of the expectation value $\left\langle...\right\rangle$ from (2.11) requires an underlying discretization a continuous notation was used to indicate that the result holds also in the continuous case. If the stochastic porous medium is not only stationary but also mixing or ergodic, and if it can be thought of as being infinitely extended, then the limit

$\overline{\phi}=\lim _{{R(\mathbb{S})\rightarrow\infty}}\phi(\mathbb{S})=\left\langle\phi\right\rangle$

(3.3)

exists and equals $\left\langle\phi\right\rangle$ . Here the diameter $R(\mathbb{G})$ of a set $\mathbb{G}$ is defined as $R(\mathbb{G})=\sup\{|{\bf r}_{1}-{\bf r}_{2}|:{\bf r}_{1},{\bf r}_{2}\in\mathbb{G}\}$ as the supremum of the distance between pairs of points. The notation $\overline{\phi}=\overline{\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r})}$ indicates a spatial average while $\left\langle\phi\right\rangle=\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r})\right\rangle$ is a configurational average.

Equation (3.3) represents always an idealization. Geological porous media for example are often heterogeneous on all scales [5]. This means that their composition or volume fraction $\phi(R(\mathbb{S}))$ does not approach a limit for $R(\mathbb{S})\rightarrow\infty$ . Equation (3.3) assumes the existence of a length scale beyond which fluctuations of the porosity decrease. This scale is used traditionally to define so called “representative elementary volumes” [76, 5]. The problem of macroscopic heterogeneity is related to the remarks in the disccusion of stationarity in section II.B.2. It will be taken up again in section III.A.5 below.

The definition of porosity in (3.1) gives the so called total porosity which has to be distinguished from the open porosity or effective porosity. Open porosity is the ratio of accessible pore volume to total volume. Accessible means connected to the surface of the sample.

The porosity of a simple porous medium $\mathbb{S}=\mathbb{P}\cup\mathbb{M}$ is related to the bulk density $\rho\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{S}}}$ the density of the matrix material $\rho\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{M}}}$ and the density of the pore space material $\rho\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}$ through

$\phi=\frac{\rho\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{M}}}-\rho\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{S}}}}{\rho\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{M}}}-\rho\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}}.$

(3.4)

Therefore porosity is conveniently determined from measuring densities using liquid buoyancy or gas expansion porosimetry [3, 1, 77, 2]. Other methods of measuring porosity include small angle neutron, small angle X-ray scattering and quantitative image analysis for total porosity [2, 77, 78, 43, 44]. Open porosity may be obtained from Xylene and water impregnation, liquid metal impregnation, Nitrogen adsorption and air or Helium penetration [77, 44].

Porosity in rocks originates as primary porosity during sedimentation or organogenesis and as secondary porosity at later stages of the geological development [1]. In sedimentary rocks the porosity is further classified as intergranular porosity between grains, intragranular or intercrystalline porosity within grains, fracture porosity caused by mechanical or chemical processes, and cavernous porosity caused by organisms or chemical processes.

III.A.1.b Specific Internal Surface Area

Similar to the porosity the specific internal surface area is an important geometric characteristic of porous media. In fact, a porous medium may be loosely defined as a medium with a large “surface to volume” ratio. The specific internal surface area is a quantitative measure for the surface to volume ratio. Often this ratio is so large that it has been idealized as infinite [78, 43, 79, 80, 81, 82, 83, 84, 85] and the application of fractal concepts has found much recent attention [58, 86, 87, 88, 84, 42, 89, 90, 91]

The specific internal surface $S$ of a two component porous medium is defined as

$S=\frac{V_{2}(\partial\mathbb{P})}{V_{3}(\mathbb{S})}$

(3.5)

where $V_{2}(\partial\mathbb{P})$ is the surface area, defined in eq. (2.3), of the boundary set $\partial\mathbb{P}$ . The surface area $V_{2}(\partial\mathbb{P})$ exists only if the internal surface or interface $\partial\mathbb{P}$ fulfills suitable smoothness requirements. Fractal surfaces would have $V_{2}(\partial\mathbb{P})=\infty$ and in such cases it is necessary to replace the Lebesgue measure in (2.3) with the Hausdorff measure or another suitable measure of the “size” of $\partial\mathbb{P}$ [58, 59, 60, 61].

The specific internal surface is a characteristic inverse length giving the surface to volume ratio of a porous medium. Typical values for unconsolidated sand are $2\times 10^{{4}}$ m ${}^{{-1}}$ , and range from $10^{{5}}$ m ${}^{{-1}}$ to $10^{{7}}$ m ${}^{{-1}}$ for sandstones [92, 3]. A piece of sandstone measuring $10$ cm on each side and having a specific internal surface of $10^{{7}}$ m ${}^{{-1}}$ contains the same area as a sports arena of dimensions $100$ m $\times 100$ m. This illustrates the importance of surface effects for all physical properties of porous media.

Specific internal surface area can be measured by similar techniques as porosity. Some commonly employed methods are given in Figure 1 together with their ranges of applicability. Particularly important methods are based on physisorption isotherms [93, 94]. The interpretation of the BET-method [93] is restricted to certain types of isotherms, and its interpretation requires considerable care. In particular, if micropores are present these will be filled spontaneously and application of the BET-analysis will lead to wrong results [44]. Other methods to determine $S$ measure the two point correlation function. As discussed further in the next section the specific internal surface area can for statistically homogeneous media be deduced from the slope of the correlation function at the origin [95].

III.A.2 Correlation Functions

Porosity and specific internal surface area are merely two numbers characterizing the geometric properties of a porous medium. Obviously these two numbers are not sufficient for a full statistical characterization of the system. A full characterization can be given in terms of multipoint correlation functions [96, 97, 98, 99, 100, 101, 102, 103, 104, 6, 105, 106, 107, 108, 109, 110, 111, 112, 8].

The average porosity $\left\langle\phi\right\rangle$ of a stationary two component porous medium is given by equation (3.2) as

$\left\langle\phi\right\rangle=\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{0})\right\rangle=\mbox{Pr}\{{\bf r}_{0}\in\mathbb{P}\}$

(3.6)

in terms of the expectation value of the random variable $\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{0})$ taking the value $1$ if the point ${\bf r}_{0}$ lies in the pore space and $0$ if not. This is an example of a so called one-point function. An example of a two-point function is the covariance function $C_{2}({\bf r}_{0},{\bf r})$ defined as the covariance of two random variables $\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r})$ and $\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{0})$ at two points ${\bf r}_{0}$ and ${\bf r}$ ,

$C_{2}({\bf r}_{0},{\bf r})=\left\langle\left[\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{0})-\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{0})\right\rangle\right]\left[\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r})-\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r})\right\rangle\right]\right\rangle.$

(3.7)

For a stationary medium the covariance function depends only on the difference ${\bf r}-{\bf r}_{0}$ which allows to set ${\bf r}_{0}={\bf 0}$ without loss of generality. This gives

$C_{2}({\bf r})=\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf 0})\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r})\right\rangle-\left\langle\phi\right\rangle^{2}.$

(3.8)

Because $\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}^{2}({\bf r})=\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r})$ it follows that $C_{2}({\bf 0})=\left\langle\phi\right\rangle(1-\left\langle\phi\right\rangle)$ . The correlation coefficient of two random variables $X$ and $Y$ is in general defined as the ratio of the covariance cov $(X,Y)$ to the two standard deviations of $X$ and $Y$ [73, 74]. It varies between $1$ and $-1$ corresponding to complete correlation or anticorrelation. The covariance function is often normalized analogous to the correlation coefficient by division with $C_{2}({\bf 0})$ to obtain the two-point correlation function

$G_{2}({\bf r})=\frac{C_{2}({\bf r})}{C_{2}({\bf 0})}=\frac{C_{2}({\bf r})}{\left\langle\phi\right\rangle(1-\left\langle\phi\right\rangle)}.$

(3.9)

An illustration of a two point correlation function can be seen in Figure 14. The porosity in (3.6) is an example of a moment function. The general $n$ -th moment function is defined as

$S_{n}({\bf r}_{1},...,{\bf r}_{n})=\left\langle\prod _{{i=1}}^{n}\left[\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{i})\right]\right\rangle$

(3.10)

where the average is defined in eq. (2.11) with respect to the probability density of microstructures given in eq. (2.10). The covariance function in (3.7) or (3.8) is an example of a cumulant function (also known as Ursell or cluster functions in statistical mechanics). The n-th cumulant function is defined as

$C_{n}({\bf r}_{1},...,{\bf r}_{n})=\left\langle\prod _{{i=1}}^{n}\left[\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{i})-\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{i})\right\rangle\right]\right\rangle=\left\langle\prod _{{i=1}}^{n}\left[\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{i})-\left\langle\phi\right\rangle\right]\right\rangle$

(3.11)

where the second equality assumes stationarity. The cumulant functions are related to the moment functions. For $n=1,2,3$ one has the relations

$\displaystyle C_{1}({\bf r}_{1})$	$\displaystyle=$	$\displaystyle 0$	(3.12)
$\displaystyle C_{2}({\bf r}_{1},{\bf r}_{2})$	$\displaystyle=$	$\displaystyle S_{2}({\bf r}_{1},{\bf r}_{2})-S_{1}({\bf r}_{1})S_{1}({\bf r}_{2})$	(3.13)
$\displaystyle C_{3}({\bf r}_{1},{\bf r}_{2},{\bf r}_{3})$	$\displaystyle=$	$\displaystyle S_{3}({\bf r}_{1},{\bf r}_{2},{\bf r}_{3})-S_{2}({\bf r}_{1},{\bf r}_{2})S_{1}({\bf r}_{3})-S_{2}({\bf r}_{1},{\bf r}_{3})S_{1}({\bf r}_{2})$	(3.14)
	$\displaystyle-$	$\displaystyle S_{2}({\bf r}_{2},{\bf r}_{3})S_{1}({\bf r}_{1})+2S_{1}({\bf r}_{1})S_{1}({\bf r}_{2})S_{1}({\bf r}_{3})$	(3.14)

The analogous moment functions may be defined for the matrix space $\mathbb{M}$ by replacing $\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}$ in all formulas with $\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{M}}}$ , and they have been called $n$ -point matrix probability functions [105, 107] or simply correlation functions [111] From (3.10),(2.10) and (2.11) the probabilistic meaning of the moment functions is found as

$S_{n}({\bf r}_{1},...,{\bf r}_{n})=\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{1})\cdot...\cdot\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{n})\right\rangle=\mbox{Pr}\{{\bf r}_{1}\in\mathbb{P},...,{\bf r}_{n}\in\mathbb{P}\}.$

(3.15)

Therefore $S_{n}({\bf r}_{1},...,{\bf r}_{n})$ is the probability that all the points ${\bf r}_{1},...,{\bf r}_{n}$ fall into the pore space.

The case $n=2$ of the second moment is of particular interest. If the pore space is stationary and isotropic then $S_{2}({\bf r})=S_{2}(|{\bf r}|)$ and one has $S_{2}(0)=\left\langle\phi\right\rangle$ . If the porous medium is also mixing then $S_{2}(\infty)=\left\langle\phi\right\rangle^{2}$ . If the pore space $\mathbb{P}$ is three dimensional, and does not contain flat twodimensional surfaces of zero thickness then its derivative at the origin is related to the specific internal surface area $S$ through

$S_{2}^{\prime}(0)=\left.\frac{dS_{2}(r)}{dr}\right|_{{r=0}}=-\frac{S}{4}.$

(3.16)

In two dimensions an analogous formula holds in which $S$ is replaced with a “specific internal length” and the denominator $4$ is replaced with $\pi$ .

The practical measurement of two point correlation functions is based on Minkowski addition and subtraction of sets [10, 37]. The Minkowski addition of two sets $\mathbb{A}$ and $\mathbb{B}$ in $\mathbb{R}^{d}$ is defined as the set

$\mathbb{A}\oplus\mathbb{B}=\{{\bf q}+{\bf r}:{\bf q}\in\mathbb{A},{\bf r}\in\mathbb{B}\}.$

(3.17)

Note that $\mathbb{A}\oplus\{{\bf r}\}=\mathbb{A}+{\bf r}$ is the translation defined in equation (2.17). Therefore $\mathbb{A}\oplus\mathbb{B}=\bigcup _{{{\bf q}\in\mathbb{B}}}(\mathbb{A}+{\bf q})=\bigcup _{{{\bf r}\in\mathbb{A}}}(\mathbb{B}+{\bf r})$ is the union of the translates $(\mathbb{B}+{\bf r})$ as ${\bf r}$ runs through $\mathbb{A}$ . The dual operation to Minkowski addition is Minkowski subtraction defined as

$\mathbb{A}\ominus\mathbb{B}=\bigcap _{{{\bf q}\in\mathbb{B}}}(\mathbb{A}+{\bf q})=C(C\mathbb{A}\oplus\mathbb{B})$

(3.18)

where $C\mathbb{A}$ denotes the complement of $\mathbb{A}$ . With these definitions the two-point function is given as

$S_{2}({\bf r})=\mbox{Pr}\{{\bf 0}\in\mathbb{P},{\bf r}\in\mathbb{P}\}=\mbox{Pr}\{{\bf 0}\in\mathbb{P}\ominus(-\mathbb{B})\}$

(3.19)

where $\mathbb{B}=\{{\bf 0},{\bf r}\}$ is the set consisting of the origin and the point ${\bf r}$ . This formula is the basis for the statistical estimation of $S_{2}({\bf r})$ and $G_{2}({\bf r})$ in image analysers from the area of the “eroded” set $\mathbb{P}\ominus(-\mathbb{B})$ . The operation $\mathbb{P}\rightarrow\mathbb{P}\ominus(-\mathbb{B})$ is called erosion of the set $\mathbb{P}$ with the set $\mathbb{B}$ and it has been used in methods to define pore size distributions which will be discussed in the next section. An example for the erosion of a pore space image by a set $\{{\bf 0},{\bf r}\}$ is shown in Figure 9.

Figure 9: Erosion of a pore space image by the set $\relax\{\relax{\bf 0},\relax{\bf r}\relax\}$ consisting of the origin $\relax{\bf 0}$ and a vector $\relax{\bf r}$ . The vector - $\relax{\bf r}$ is shown as an arrow in the image. The overlap of the two images is rendered darker than the rest.

The original image (shown Figure 12) is obtained from a cross section micrograph of a Savonnier oolithic sandstone. Two copies of the image are displaced relative to each other by a vector ${\bf r}$ as indicated in Figure 9. The two images are rendered in grey, and their intersection is coloured black. The area of the intersection is an estimate for $S_{2}({\bf r})$ .

The main advantage of the correlation function method for characterizing porous media is that it provides a set of well defined functions of increasing complexity for the geometrical description. In practice one truncates the hierarchy of correlation functions at the two-point functions. While this provides much more information about the geometry than the porosity and specific surface area alone, many important properties of the medium (such as its connectivity) are buried in higher order functions. ² (This is a footnote:) ² Two points are called connected if there exists a path between them which lies completely inside the pore space. Therefore the probabilistic description of connectedness properties requires multipoint correlation functions involving all the points which make up the path. Depending on the required accuracy a simple two point function for a three dimensional stationary but anisotropic two component medium could be specified by $10^{6}$ to $10^{9}$ data points which would be economical according to the criterion adopted previously. An $n$ point function with the same accuracy would require $10^{{6(n-1)}}$ to $10^{{9(n-1)}}$ data points. Specifying five or higher point functions quickly becomes just as impractical as specifying a given geometry completely.

III.A.3 “Pore Size” Distributions

In certain porous materials such as wood (see Figure 5) it is natural to identify cylindrically shaped pores and to represent their disorder through a distribution of pore diameters. In other media such as systems with cavernous or oomoldic porosity it is possible to identify roughly convex pore bodies analogous to convex sand grains dispersed in a uniform background. If the radius $R$ of the cylindrical capillaries or spherical pore bodies in such media is randomly distributed then the pore size distribution function $\Pi(r)$ can be defined as

$\Pi(r)=\mbox{Pr}\{ R\leq r\}$

(3.20)

giving the probability that the random radius $R$ of the cylinders or spheres is smaller than $r$ . For general porous microstructures, however, it is difficult to define “pores” or “pore bodies”, and the concept of pore size distribution remains ill defined.

Nevertheless many authors have introduced a variety of well defined probability distributions of length for arbitrary media which intended to overcome the stated difficulty [3, 2, 113, 114, 115, 116, 117, 118, 119]. The concept of pore size distributions enjoys continued popularity in most fields dealing with porous materials. Recent examples can be found in chromatography [120, 121], membranes [122, 123, 124], polymers [125], ceramics [126, 127, 128, 129], silica gels [130, 131, 132], porous carbon [133, 134], cements [135, 136, 137], rocks and soil science [138, 139, 140, 141, 142, 143], fuel research [144], separation and adhesion technology [30, 145] or food engineering [146]. The main reasons for this popularity are adsorption measurements [147, 30, 148] and mercury porosimetry [149, 150, 151, 152].

III.A.3.a Mercury Porosimetry

The “pore size distribution” of mercury porosimetry is not a geometric but a physical characteristic of a porous medium. Mercury porosimetry is a transport and relaxation phenomenon [43, 153], and its discussion would find a more appropriate place in chapter V below. On the other hand “pore size distributions” are routinely measured in practice using mercury porosimetry, and many readers will expect its discussion in a section on pore size distributions. Therefore pore size distributions from mercury porosimetry are discussed already here together with other definitions of this important concept.

Mercury porosimetry is based on the fact that mercury is a strongly nonwetting liquid on most substrates, and that it has a high surface tension. To measure the “pore size distribution” a porous sample $\mathbb{S}$ with pore space $\mathbb{P}$ is evacuated inside a pycnometer pressure chamber at elevated temperatures and low pressures [43]. Subsequently the sample is immersed into mercury and an external pressure is applied. As the pressure is increased mercury is injected into the pore space occupying a subset $\mathbb{P}_{{\mbox{\tiny Hg}}}(P)$ of the pore space which depends on the applied external pressure $P$ . The experimenter records the injected volume of mercury $V(\mathbb{P}_{{\mbox{\tiny Hg}}}(P))$ as a function of the applied external pressure. If the volume of the pore space $V(\mathbb{P})$ is known independently then this gives the saturation $S_{{\mbox{\tiny Hg}}}(P)=V(\mathbb{P}_{{\mbox{\tiny Hg}}}(P))/V(\mathbb{P})$ as a function of pressure. The cumulative “pore size” distribution function $\Pi _{{\mbox{\tiny Hg}}}(r)=\mbox{Pr}\{ R\leq r\}$ of mercury porosimetry is now defined by

$\Pi _{{\mbox{\tiny Hg}}}(r)=1-S_{{\mbox{\tiny Hg}}}\left(\frac{2\sigma _{{\mbox{\tiny Hg}}}\cos\theta}{r}\right).$

(3.21)

For rocks a contact angle $\theta _{{\mbox{\tiny Hg}}}\approx 135^{\circ}$ and surface tension with vacuum of $\sigma _{{\mbox{\tiny Hg}}}\approx 0.48$ Nm ${}^{{-1}}$ are commonly used [1, 43, 153]. The definition of $\Pi _{{\mbox{\tiny Hg}}}(r)$ is based on the equation

$P_{{cap}}=\frac{2\sigma _{{\mbox{\tiny Hg}}}\cos\theta}{r}$

(3.22)

for the capillary pressure $P_{{cap}}$ which expresses the force balance in a single cylindrcal capillary tube. Equation (3.21) follows from (3.22) if it is assumed that the saturation history $S_{{\mbox{\tiny Hg}}}(P)$ is identical to that obtained from the so called capillary tube model discussed in section III.B.1 below. The capillary tube model is a hypothetical porous medium consisting of parallel nonintersecting cylindrical capillaries of random diameter.

The fact that $\Pi _{{\mbox{\tiny Hg}}}(r)$ is not a geometrical quantity but a capillary pressure function is obvious from its definition. It depends on physical properties such as the nature of the injected fluid or wetting properties of the walls. For a suitable choice of tube diameter the function $S_{{\mbox{\tiny Hg}}}(P)$ of pressure could equally well be translated into a distribution of the wetting angles $\theta$ . $S_{{\mbox{\tiny Hg}}}(p)$ shows hysteresis implying that the pore size distribution $\Pi _{{\mbox{\tiny Hg}}}(r)$ is process dependent.

Although $\Pi _{{\mbox{\tiny Hg}}}(r)$ is not a geometrical quantity it contains much useful information about the microstructure of the porous sample. An example for the information obtained from mercury porosimetry is shown in Figure 10 together with an image of the rock for which it was measured [1].³ (This is a footnote:) ³1 MPa = 10 bar = 10 ${}^{7}$ dyn/cm ${}^{2}$ = 9.869 atm = 145.04 psi The rock is an example for a medium with hollow pores. The correct interpretation of the saturation history $S_{{\mbox{\tiny Hg}}}(p)$ obtained from mercury porosimetry continues to be an active research topic [154, 155, 156, 157, 158, 159, 160].

Figure 10: a) Crossed nicols image of Dolomite oncolithe. The image width corresponds to 2.3 mm. The porosity is 11%, the permeability 11 md. b) Capillary presssure curve for the rock shown in the image. The right axis shows the ”pore size” scale. [Reproduced with permission from [1]. Copyright Springer Verlag 1982.]

III.A.3.b Random Point Methods

Several authors [3, 113, 116] suggest to define the “pore size” by first choosing a point ${\bf r}\in\mathbb{P}$ at random in the pore space, then to choose a compact set $\mathbb{K}$ containing ${\bf r}$ , and finally to enlarge $\mathbb{K}$ until it first intersects the matrix space $\mathbb{M}$ . In its simplest version [3, 116] the set $\mathbb{K}$ is chosen as a small sphere $\mathbb{B}({\bf r},\epsilon)$ of radius $\epsilon$ (see (2.5) above). Then the pore size distribution $\Pi _{{rp}}(x)$ of the random point method is defined as the distribution function $\Pi _{{rp}}(x)=\mbox{Pr}\{ R({\bf r})\leq x\}$ of the random variable $R({\bf r})$ defined as

$R({\bf r})=\inf\{\epsilon:\mathbb{B}({\bf r},\epsilon)\cap\mathbb{M}\neq\emptyset\}.$

(3.23)

Here $R$ is a random variable because ${\bf r}\in\mathbb{P}$ is chosen at random.

In a more sophisticated version of the same idea the set $\mathbb{K}$ is chosen as a small coordinate cross whose axes are then increased independently until they first touch the matrix space. [113, 2]. This gives direction dependent pore size distributions.

The main weakness of such a definition is that it is imprecise. This becomes apparent from the fact that the randomness of the pore sizes $R({\bf r})$ does not arise from the irregularities of the pore space, but from the random placement of ${\bf r}$ . Consider a regular pore space consisting of nonintersecting spheres of equal radii $R_{0}$ centered at the vertices of a simple (hyper)cubic lattice. Assuming that the points ${\bf r}$ are chosen at random with a uniform distribution it follows that the pore size distribution $\Pi _{{rp}}$ is not given by a $\delta$ -function at $R_{0}$ but instead as a uniform distribution on the interval $[0,R_{0}]$ . More dramatically, exactly the same pore size distribution is obtained for every pore space made from nonintersecting spheres of equal radii, no matter whether they are placed randomly or not. In addition $\Pi _{{rp}}(x)$ can be changed arbitrarily by changing the distribution function governing the random placement of ${\bf r}$ .

A more precise formulation of the same idea is possible in terms of conditional probabilities. In section III.A.4 it will be seen that the vague ideas underlying random point methods are given a precise definition in the form of contact and chord length distributions [10, 37]

III.A.3.c Erosion Methods

Another approach to the definition of a geometrical pore size distribution [117, 118, 119] is borrowed from the erosion operation in image processing [161, 162]. Erosion is defined in terms of Minkowski addition and subtraction of sets introduced above in (3.17) and (3.18). The erosion of a set $\mathbb{A}$ by a set $\mathbb{B}$ is defined as the map $\mathbb{A}\rightarrow\mathbb{A}\ominus(-\mathbb{B})=\bigcap _{{{\bf q}\in\mathbb{B}}}(\mathbb{A}-{\bf q})$ . The erosion was illustrated in Figure 9 for a pore space image and a set $\mathbb{B}=\{{\bf 0},{\bf r}\}$ .

The method for locating “pore chambers”, “pore channels” and “pore throats” suggested in [117] is based on eroding the matrix space $\mathbb{M}$ of a two component porous medium. A ball $\mathbb{B}({\bf 0},\epsilon)$ is chosen as the structuring element. The erosion $\mathbb{M}\ominus(-\mathbb{B}({\bf 0},\epsilon))$ shrinks the matrix space $\mathbb{M}$ . The erosion operation is repeated until the matrix space decomposes into disconnected fragments. Continuing the erosion the pieces may either fragment again or become convex. If a piece becomes convex it is called a “grain”. The centroid of the conves grain is called a “grain center”. Reversing the erosion process allows to locate the point of first/last contact of two fragments. Connecting neighbouring grain centers by a path through their last contact point produces a network model of the grain space. Having defined a network of grain centers and last contact points the authors of [117], and their followers [119, 118], suggest to erect contact “surfaces” in each contact point. A contact plane is defined as a “minimum area cross section” of $\mathbb{M}$ . Subsequently a ball is placed at each grain center and continually enlarged. When the enlargement encounters a surface plane the ball is truncated at the surface plane and only its non-truncated pieces continue to grow until the sample space is completely filled with the inflated grains. The intersection points of three or more planes in the resulting tesselation of space are defined to be pore chambers. The intersection lines of two planes are called pore channels, and “minimal-area cross sections” of the pore space $\mathbb{P}$ along the pore channels are called pore throats. The pore throats are not unique, and sensitive to details of the local geometry. The pore chambers and the pore channels will in general not lie in the pore space.

A drawback of this procedure is that it is less unique than it seems at first sight. The network constructed from eroding the pore space is not unique because the erosion operation involves the set $\mathbb{B}$ as a structuring element, and hence there are infinitely many erosions possible. The resulting grain network depends on the choice of the set $\mathbb{B}$ , a fact which is not discussed in [117, 119]. Figure 11 shows an example where erosion with a sphere produces three grains while erosion with an ellipsoid produces only two grains. The original set consists of the grey and black region, the eroded part is coloured grey, and the residual set is coloured black. The theoretically described procedure for determining the network was not carried out in practice [117]. Instead “subjective human preprocessing” ([117],p.4158) was used to determine the network.

Figure 11: Erosion $\mathbb{A}\oplus(-\mathbb{B})$ of a set $\mathbb{A}$ by a set $\mathbb{B}$ for the cases where $\mathbb{B}$ is a sphere or an ellipsoid. For a sphere the eroded set (shown in black) has three components while for an ellipsoid it consists of only two connected components.

III.A.3.d Hydraulic Radius Method

The hydraulic radius method [114, 115, 2] for determining pore size distributions is based on the idea of “symbolically closing pore throats”. The definition of pore throats is given in terms of “cross sections” of the pore space. A cross section $\mathbb{C}$ could be defined as the intersection of a plane $\mathbb{E}({\bf q},{\bf r}_{0})$ , characterized by its unit normal ${\bf q}\in\partial\mathbb{B}({\bf 0},1)$ and a point ${\bf r}_{0}$ in the plane, with the pore space $\mathbb{P}$ and some suitable set $\mathbb{G}({\bf q},{\bf r}_{0})$ which represents the region of interest and could depend on the choice of plane. In symbols $\mathbb{C}({\bf q},{\bf r}_{0})=\mathbb{P}\cap\mathbb{E}({\bf q},{\bf r}_{0})\cap\mathbb{G}({\bf q},{\bf r}_{0})$ . A pore throat containing the point ${\bf r}_{0}$ is then defined as

$\mathbb{C}^{*}({\bf r}_{0})=\min _{{{\bf q}\in\partial\mathbb{B}({\bf 0},1)}}\frac{V_{2}(\mathbb{C}({\bf q},{\bf r}_{0}))}{V_{1}(\partial\mathbb{C}({\bf q},{\bf r}_{0}))}$

(3.24)

where the minimum is taken over the unit sphere of orientations of the planes. A pore throat is now defined as a local minimum of the function $\mathbb{C}^{*}({\bf r}_{0})$ as ${\bf r}_{0}$ is varied over the pore space. This ideal definition has in practice been replaced with a subjective choice of orientations based on the assumption of isotropy [114, 115, 2]. After constructing all pore throats of a medium the pore space becomes divided into separate compartments called “pore bodies” whose “size” can then be measured by a suitable measure such as the volume to the power $1/3$ .

The definition of pore throats in hydraulic radius methods is very sensitive to surface roughness. This is readily seen from an idealized spherical pore with a few spikes. Another problem as remarked in [115], page 586, is that “the size of a pore body is not readily related in a unique manner to any measurable physical quantity”.

III.A.4 Contact and Chord Length Distributions

Chord length distributions [163, 164, 165, 166] are special cases of so called contact distributions [10, 37]. Consider the random matrix space $\mathbb{M}$ of a two component stochastic porous medium and choose a compact set $\mathbb{K}$ containing the origin ${\bf 0}$ . Then the contact distribution is defined as the conditional probability

$\displaystyle\Pi _{\mathbb{K}}(x)$	$\displaystyle=$	$\displaystyle 1-\mbox{Pr}\{{\bf 0}\notin\mathbb{M}\oplus(-x\mathbb{K}))\|{\bf 0}\notin\mathbb{M}\}$		(3.25)
	$\displaystyle=$	$\displaystyle 1-\frac{\mbox{Pr}\{{\bf 0}\notin\mathbb{M}\oplus(-x\mathbb{K})\}}{1-\phi}$		(3.25)

for $x\geq 0$ . Here $\phi$ denotes the bulk porosity as usual. Two special choices of the compact set $\mathbb{K}$ are of particular importance. These are the unit sphere $\mathbb{K}=\mathbb{B}({\bf 0},1)\subset\mathbb{R}^{3}$ and the the unit interval $\mathbb{K}=[0,1]\subset\mathbb{R}$ .

For a unit sphere $\mathbb{K}=\mathbb{B}({\bf 0},1)$ the quantity $1-\Pi _{{\mathbb{B}({\bf 0},1)}}(x)$ is the conditional probability that a fixed point in the pore space is the center of a sphere of radius $x$ contained completely in the pore space, under the condition that the chosen point does not belong to $\mathbb{M}$ . If $\mathbb{M}$ is isotropic and its boundary is sufficiently smooth then the specific internal surface $S$ can be obtained from the derivative at the origin as [37]

$S=\phi\left.\frac{d}{dx}\Pi _{{\mathbb{B}({\bf 0},1)}}(x)\right|_{{x=0}}.$

(3.26)

The spherical contact distribution $\Pi _{{\mathbb{B}({\bf 0},1)}}(x)$ provides a more precise formulation of the random point generation methods for pore size distributions [3, 116] discussed in subsection III.A.3.b.

For the unit interval $\mathbb{K}=[0,1]$ the contact distribution $\Pi _{{[0,1]}}(x)$ is related to the chord length distribution $\Pi _{{cl}}(x)$ giving the probability that an interval in the intersection of $\mathbb{P}$ with a straight line containing the unit interval has a length smaller than $x$ . This provides a more precise formulation of the random point generation ideas in [113, 2]. The relation between the contact distribution and the chord length distribution is given by the equation

$\Pi _{{[0,1]}}(x)=1-\frac{\int _{x}^{\infty}(y-x)d\Pi _{{cl}}(y)}{\int _{0}^{\infty}yd\Pi _{{cl}}(y)}$

(3.27)

where the denominator on the right hand side gives the mean chord length $\overline{\ell}=\int _{0}^{\infty}yd\Pi _{{cl}}(y)$ . The mean chord length is related to the specific internal surface area through $\overline{\ell}=4\phi/S$ . In section III.A.2 it was mentioned that the specific internal surface area can be obtained from the two point function (see 3.16). Therefore also the mean chord length can be related to the correlation function through

$\overline{\ell}=-\frac{\phi}{S_{2}^{\prime}(0)}$

(3.28)

Along these lines it has been suggested in [167] that the full chord length distribution can be obtained directly from small angle scattering experiments.

Contact and chord length distributions provide much more geometrical information about the porous medium than the porosity and specific surface area, and are at the same time not as unnecessarily detailed as the complete specification of the deterministic or stochastic geometry. Depending on the demands on accuracy a contact or chord length distribution may be specified by 10 to 1000 numbers irrespective of the microscopic resolution. This should be compared with $10^{{24}}$ or $2^{{10^{{24}}}}$ numbers for a full deterministic or stochastic characterization.

III.A.5 Local Geometry Distributions

III.A.5.a Local Porosity Distribution

Local porosity distributions, or more generally local geometry distributions, provide a well defined general geometric characterization of stochastic porous media. [168, 169, 170, 171, 172, 173, 174, 175]. Local porosity distributions were mainly developed as an alternative to pore size distributions (see section III.A.3). They are intimately related with the theory of finite size scaling in statistical physics [176, 177, 64, 178]. Although fluctuations in the porosities have been frequently discussed [179, 10, 37, 2, 5, 180, 181, 76, 182], the concept of local porosity distributions and its relation with correlation functions was developed only recently [168, 169, 170, 171, 172, 173, 174, 175]. More applications are being developed [183, 184].

Local porosity distributions can be defined for deterministic as well as for stochastic porous media. For a single deterministic porous medium consider a partitioning $\mathcal{K}=\{\mathbb{K}_{1},...,\mathbb{K}_{M}\}$ of the sample space $\mathbb{S}$ into $M$ mutually disjoint subsets, called measurement cells $\mathbb{K}_{j}$ . Thus $\bigcup _{{j=1}}^{M}\mathbb{K}_{j}=\mathbb{S}$ and $\mathbb{K}_{i}\cap\mathbb{K}_{j}=\emptyset$ if $i\neq j$ . A particular partitioning was used in the orginal publications [178, 169, 170, 171] where the $\mathbb{K}_{j}$ are unit cells centered at the vertices of a Bravais lattice superimposed on $\mathbb{S}$ . This has the convenient feature that the $\mathbb{K}_{j}$ are translated copies of one and the same set, and they all have the same shape. An example is illustrated in Figure 12 showing a quadratic lattice as the measurement grid in two dimensions superposed on a thin section of an oolithic sandstone.

Figure 12: Measurement lattice of squares superimposed on the discretized thin section image of a sandstone. The pore space $\mathbb{P}$ is coloured black, the matrix space is rendered white.

The local porosity inside a measurement cell $\mathbb{K}_{j}$ is defined as

$\phi(\mathbb{K}_{j})=\frac{V(\mathbb{P}\cap\mathbb{K}_{j})}{V(\mathbb{K}_{j})}=\frac{1}{M_{j}}\sum _{{{\bf r}_{i}\in\mathbb{K}_{j}}}\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{i})$

(3.29)

where the second equality applies in case of discretized space and $M_{j}$ denotes the number of volume elements or voxels in $\mathbb{K}_{j}$ . Thus the empirical one cell local porosity density function is defined as

$\widetilde{\mu}(\phi;\mathcal{K})=\frac{1}{M}\sum _{{j=1}}^{M}\;\delta(\phi-\phi(\mathbb{K}_{j}))$

(3.30)

where $\delta(x)$ is the Dirac $\delta$ -distribution. Obviously the distribution depends on the choice of partitioning the sample space. Two extreme partitions are of immediate interest. The first arises from setting $M=N$ and thus each $\mathbb{K}_{j}$ contains only one individual volume element $\mathbb{K}_{j}=\{{\bf r}_{j}\}$ with $1\leq j\leq M=N$ . In this case $\phi(\mathbb{K}_{j})=0$ or $\phi(\mathbb{K}_{j})=1$ depending on whether the volume element falls into matrix space (0), or pore space (1). This gives immediately

$\widetilde{\mu}(\phi;\{\{{\bf r}_{1}\},...,\{{\bf r}_{N}\}\})=\phi(\mathbb{S})\delta(\phi-1)+(1-\phi(\mathbb{S}))\delta(\phi)$

(3.31)

where $\phi(\mathbb{S})$ is the total porosity. The other extreme arises for $M=1$ and thus $\mathbb{K}_{1}=\mathbb{S}$ the measurment cell coincides with the sample space. In this case obviously

$\widetilde{\mu}(\phi;\{\mathbb{S}\})=\delta(\phi-\phi(\mathbb{S})).$

(3.32)

Note that in both extreme cases the local porosity density is completely determined by the total porosity $\phi(\mathbb{S})$ , which equals $\phi(\mathbb{S})=\overline{\phi}$ if the sample is sufficiently large and mixing or ergodicity (3.3) holds.

For a stochastic porous medium the one cell local porosity density function is defined for each measurement cell as

$\mu(\phi;\mathbb{K}_{j})=\left\langle\delta(\phi-\phi(\mathbb{K}_{j}))\right\rangle$

(3.33)

where $\mathbb{K}_{j}\in\mathcal{K}$ is an element of the partitioning of the sample space. For the finest partition with $M=N$ and $\mathbb{K}_{j}=\{{\bf r}_{j}\}$ one finds now using eq. (3.2)

$\displaystyle\mu(\phi;\{{\bf r}_{j}\})$	$\displaystyle=$	$\displaystyle\mbox{Pr}\{ X({\bf r}_{j})=1\}\delta(\phi-1)+\mbox{Pr}\{ X({\bf r}_{j})=0\}\delta(\phi)$		(3.34)
	$\displaystyle=$	$\displaystyle\left\langle\phi\right\rangle\delta(\phi-1)+(1-\left\langle\phi\right\rangle)\delta(\phi)$		(3.34)

independent of $j$ . If mixing (3.3) holds then $\phi(\mathbb{S})=\left\langle\phi\right\rangle=\overline{\phi}$ if the sample becomes sufficiently large, and the result becomes identical to equation (3.31) for deterministic media. In the other extreme of the coarsest partition one finds

$\mu(\phi;\mathbb{S})=\left\langle\delta(\phi-\phi(\mathbb{S}))\right\rangle$

(3.35)

which may in general differ from (3.32) even if the sample becomes sufficiently large, and mixing holds. This is an important observation because it emphasizes the necessity to consider more carefully the infinite volume limit $\mathbb{S}\rightarrow\mathbb{R}^{d}$ .

If a large deterministic porous medium is just a realization of a stochastic medium obeying the mixing property, and if the $\mathbb{K}_{j}$ are chosen such that the random variables $\phi(\mathbb{K}_{j})$ are independent, then Gliwenkos theorem [185] of mathematical statistics guarantees that the empirical one cell distribution approaches $\mu(\phi;\mathbb{K})$ in the limit $M\rightarrow\infty$ . In symbols

$\lim _{{M\rightarrow\infty}}\widetilde{\mu}(\phi;\mathcal{K})=\mu(\phi;\mathbb{K}_{j})$

(3.36)

where the right hand side is independent of the choice of $\mathbb{K}_{j}$ . Therefore $\mu$ and $\widetilde{\mu}$ are identified in the following. This identification emerges also from considering average and variance as shown next.

Define the average local porosity $\overline{\phi}=\int _{0}^{1}\phi\mu(\phi)d\phi$ as the first moment of the local porosity distribution. For a stationary (=homogeneous) porous medium the definitions (3.33) and (3.29) immediately yield

$\displaystyle\overline{\phi(\mathbb{K}_{j})}$	$\displaystyle=$	$\displaystyle\int _{0}^{1}\phi\mu(\phi;\mathbb{K}_{j})d\phi$	(3.37)
	$\displaystyle=$	$\displaystyle\left\langle\phi(\mathbb{K}_{j})\right\rangle$
	$\displaystyle=$	$\displaystyle\frac{1}{M_{j}}\sum _{{{\bf r}_{i}\in\mathbb{K}_{j}}}\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{i})\right\rangle$
	$\displaystyle=$	$\displaystyle\left\langle\phi\right\rangle$

where $\mathbb{K}_{j}$ is a measurement cell. Similarly, the variance of local porosities reads

$\displaystyle\overline{(\phi(\mathbb{K}_{j})-\overline{\phi(\mathbb{K}_{j})})^{2}}$	$\displaystyle=$	$\displaystyle\int _{0}^{1}(\phi-\overline{\phi(\mathbb{K}_{j})})^{2}\mu(\phi;\mathbb{K}_{j})d\phi$	(3.38)
	$\displaystyle=$	$\displaystyle\int _{0}^{1}\phi^{2}\left\langle\delta(\phi-\phi(\mathbb{K}_{j}))\right\rangle d\phi\;-\;\overline{\phi}^{2}$
	$\displaystyle=$	$\displaystyle\left\langle\phi(\mathbb{K}_{j})^{2}\right\rangle-\left\langle\phi(\mathbb{K}_{j})\right\rangle^{2}$
	$\displaystyle=$	$\displaystyle\frac{1}{M_{j}^{2}}\left\langle\right[M_{j}\phi-\sum _{{i=1}}^{{M_{j}}}\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{i})\left]{}^{2}\right\rangle$
	$\displaystyle=$	$\displaystyle\frac{\left\langle\phi\right\rangle}{M_{j}}(1-\left\langle\phi\right\rangle)+\frac{2}{M_{j}^{2}}\sum _{{{i\neq k}\atop{{\bf r}_{i},{\bf r}_{k}\in\mathbb{K}_{j}}}}C_{2}({\bf r}_{i}-{\bf r}_{k}).$

This result is important for two reasons. Firstly it relates local porosity distributions to correlation functions discussed above in section III.A.2. Secondly it shows that the variance depends inversely on the volume $M_{j}$ of the measurement cell. This observation reconciles eq. (3.35) for stochastic media with eq. (3.32) for deterministic media because it shows that $\mu(\phi;\mathbb{K}_{j})$ approaches a degenerate $\delta$ -distribution in the limit $M_{j}\rightarrow\infty$ . Together with (3.37) and $M_{j}=N$ and $M=1$ for $\mathbb{K}_{j}=\mathbb{S}$ this shows that (3.35) and (3.32) become equivalent.

The $n$ -cell local porosity density function $\mu _{n}(\phi _{1},...,\phi _{n};\mathbb{K}_{1},...,\mathbb{K}_{n})$ is the probability density to find local porosity $\phi _{1}$ in measurement cell $\mathbb{K}_{1}$ , $\phi _{2}$ in $\mathbb{K}_{2}$ and so on until $n$ . Formally it is defined by generalizing (3.33) to read

$\mu _{n}(\phi _{1},...,\phi _{n};\mathbb{K}_{1},...,\mathbb{K}_{n})=\left\langle\delta(\phi _{1}-\phi(\mathbb{K}_{1}))\cdot...\delta(\phi _{n}-\phi(\mathbb{K}_{n}))\right\rangle$

(3.39)

where the sets $\mathbb{K}_{1},...,\mathbb{K}_{n}\in\mathcal{K}$ are a subset of measurement cells in the partition $\mathcal{K}$ . Note that the $n$ -cell functions $\mu _{n}$ are only defined for $M>n$ , i.e. if there are sufficient number of cells. In particular for the extreme case $M=1$ only the one-cell function is defined. In the extreme case $M=N$ of highest resolution ( $\mathbb{K}_{j}=\{{\bf r}_{j}\}$ ) the moments of the $n$ -cell local porosity distribution reproduce the moment functions (3.10) of the correlation function approach (section III.A.2) as

$\displaystyle\overline{\phi _{{i_{1}}}\cdot...\cdot\phi _{{i_{m}}}}$	$\displaystyle=$	$\displaystyle\int _{0}^{1}...\int _{0}^{1}\phi _{{i_{1}}}\cdot...\cdot\phi _{{i_{m}}}\mu _{n}(\phi _{1},...,\phi _{n};\{{\bf r}_{1}\},...,\{{\bf r}_{n}\})d\phi _{1}...d\phi _{n}$	(3.40)
	$\displaystyle=$	$\displaystyle\left\langle\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{{i_{1}}})\cdot...\cdot\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\mathbb{P}}}({\bf r}_{{i_{m}}})\right\rangle$
	$\displaystyle=$	$\displaystyle S_{m}({\bf r}_{{i_{1}}},...,{\bf r}_{{i_{m}}})$

where $m\leq n$ must be fulfilled. This provides a connection between local porosity approaches and correlation function approaches. Simultaneously it shows that the local porosity distributions are significantly more general than correlation functions. Already the one-cell function $\mu$ contains more information than the two-point correlation functions $S_{2}$ or $C_{2}$ as demonstrated on test images in [171].

The most important practical aspect of local porosity distributions is that they are easily measurable in an independent experiment. Experimental determinations of local porosity distributions have been reported in [186, 173, 174, 175]. They were obtained from two dimensional sections through a sample of sintered glass beads. An example from [175] is shown in Figure 13.

Figure 13: Pore space image of sintered $250~\mu$ m glass beads with total porosity of 10.7%. The pore space is rendered in black, the matrix space in white. [175]

The porousmedium was made from sintering glass beads of roughly $250~\mu$ m diameter. Scanning electron migrographs obtained from the specimen were then digitized with a spatial resolution of $4.1~\mu$ m per pixel. The pore space is represented black in Figure 13 and the total porosity is $\overline{\phi}=10.7\%$ . The pixel-pixel correlation function of the pore space image is shown in Figure 14.

Figure 14: Pixel-pixel porosity correlation function $G_{2}(r)$ for the image shown in Figure 13. The distance $r$ is measured in pixels corresponding to 4.1 $\mu$ m per pixel [175].

The local porosity distribution was then measured using a square lattice with lattice constant $L$ . The results are shown in Figure 15 for measurement cell sizes of $L=40,50,60$ and $L=80$ pixels.

Figure 15: Local porosity probability density function $\mu(\phi;L)$ measured from square shaped measurement cells from the image shown in Figure 13. The dash-dotted line corresponds to $L=40$ pixels, the solid line to $L=50$ , the dashed line to $L=60$ and the dotted line to $L=80$ pixels. One pixels corresponds to $4.1~\mu$ m. [175]

Measurements of the local porosity distribution from three dimensional pore space images are currently carried out [183]. The pore space is reconstructed from serial sections which is a standard, but costly, technique to obtain threedimensional pore space representations [162, 187, 117, 113, 188, 114, 2, 183]. The advent of synchrotron microtomography [119, 4] and laser scanning confocal microscopy [188a] promises to reduce the cost and effort.

III.A.5.b Local Geometry Entropies

The linear extension $L$ of the measurement cells is the length scale at which the pore space geometry is described. The length $L$ can be taken as the side length of a hypercubic measurement cell, or more generally as the diameter $R(\mathbb{K})=\sup\{|{\bf r}_{1}-{\bf r}_{2}|:{\bf r}_{1},{\bf r}_{2}\in\mathbb{K}\}$ of a cell $\mathbb{K}$ defined as the supremum of the distance between pairs of points. As $L\rightarrow\infty$ the local porosity distribution approaches a $\delta$ -distribution concentrated at $\overline{\phi}$ according to (3.32) and (3.35). For $l\rightarrow 0$ on the other hand it approaches two $\delta$ -distributions concentrated at $0$ and $1$ according to (3.31) and (3.34). In both limits the local porosity distribution contains only the bulk porosity $\overline{\phi}=\left\langle\phi\right\rangle=\phi(\mathbb{S})$ as a geometric parameter. At intermediate scales the distribution contains additional information, such as the variance of the porosity fluctuations. This suggests to search for an intermediate scale $L^{*}$ which provides an optimal description. Several criteria for determining $L^{*}$ were discussed in [171]. One interesting possibility is to optimize an information measure or entropy associated with $\mu(\phi;L)$ , and to define the entropy function

$I(L)=\int _{0}^{1}\mu(\phi;L)\log\mu(\phi;L)\; d\phi$

(3.41)

relative to the conventional a priori uniform distribution. The so called entropy length is then obtained from the extremality condition

$\left.\frac{d\mu(\phi;L)}{dL}\right|_{{L=L^{*}}}=0.$

(3.42)

That the entropy length $L^{*}$ exists and is well defined was first demonstrated in [171] using synthetic computer generated images. Figure 16 shows the function $I(L)$ calculated for the image displayed in Figure 13. A clear minimum appears at $L\approx 40$ pixels corresponding to $164~\mu$ m.

Figure 16: Entropy function as a function of the linear size of the measurement cell in pixels for the local porosity probability density functions obtained from the image shown in Figure 13. One pixel corresponds to $4.1~\mu$ m.[175]

A similar entropy analysis was recently discussed in [189] for thin film morphologies of composites. The entropy function $H^{*}(L)$ defined in [189] is not equivalent to $I(L)$ above because it is not additive. Nevertheless it was recently found [184] that $H^{*}(L)$ and $I(L)$ give the same value of $L^{*}$ .

III.A.5.c Local Specific Internal Surface Distributions

Local specific internal surface area distributions are a natural generalization of local porosity ditributions which was first introduced in [171] in the study of fluid transport in porous media. Define the local specific internal surface area in a cell $\mathbb{K}_{j}$ as

$S(\mathbb{K}_{j})=\frac{V_{2}(\partial\mathbb{P}\cap\mathbb{K}_{j})}{V_{3}(\mathbb{K}_{j})}=\frac{1}{M_{j}}\sum _{{{\bf r}_{i}\in\mathbb{K}_{j}}}\chi\rule[-4.3pt]{0.0pt}{8.6pt}_{{\partial\mathbb{P}}}({\bf r}_{i})$

(3.43)

which is analogous to equation (3.29). Generalizing equations (3.30) and (3.33) the local specific internal surface area probability density is defined as

$\mu(S;\mathbb{K}_{j})=\left\langle\delta(S-S(\mathbb{K}_{j}))\right\rangle.$

(3.44)

in analogy with equation (3.33). The joint probability density $\mu(\phi,S;\mathbb{K}_{j})$ to find a local porosity $\phi$ and local specific internal surface area in the range $\phi$ to $\phi+d\phi$ and $S$ to $S+dS$ will be called local geometry distribution and it is defined as the probability density

$\mu(\phi,S;\mathbb{K}_{j})=\left\langle\delta(\phi-\phi(\mathbb{K}_{j}))\delta(S-S(\mathbb{K}_{j}))\right\rangle.$

(3.45)

The average specific internal surface area in a measurement region $\mathbb{K}$ is then obtained from the local geometry distribution as

$\overline{S(\mathbb{K})}=\int _{0}^{\infty}\int _{0}^{\infty}S\mu(\phi,S;\mathbb{K})\; d\phi dS$

(3.46)

and it represents an important local length scale.

Of course local geometry distributions can be extended to include other well defined geometric characteristics such as mean curvature or topological invariants. The definition of the generalized local geometry distribution is then obtained by generalizing (3.45).

III.A.5.d Local Percolation Probability

In addition to the local porosity distributions and local specific internal surface area distributions it is necessary to characterize the geometrical connectivity properties of a porous medium. This is important for discussing transport properties which depend critically on the connectedness of the pore space, but are less sensitive to its overall porosity or specific internal surface.

Two points inside the pore space $\mathbb{P}$ of a two component porous medium are called connected if there exists a path contained entirely within the pore space which connects the two points. Using this connectivity criterion a cubic measurement cell $\mathbb{K}$ is called percolating if there exist two points on opposite surfaces of the cell which are connected to each other. The local percolation probability $\lambda(\phi,S;\mathbb{K})$ is defined as the probability to find a percolating geometry in measurement cells $\mathbb{K}$ whose local porosity is $\phi$ and whose local specific internal surface area is $S$ . In practice the estimator for $\lambda(\phi,S;\mathbb{K})$ is the fraction of percolating measurement cells which have the prescribed values of $\phi$ and $S$ .

The average local percolation probability defines an important global geometric characteristic

$p(\mathbb{K})=\int _{0}^{\infty}\int _{0}^{1}\lambda(\phi,S,;\mathbb{K})\mu(\phi,S;\mathbb{K})\; d\phi dS$

(3.47)

which gives the total fraction of percolating local geometries. It will be seen in section V.B.4 and V.C.6 that $p(\mathbb{K})$ is an important control parameter for the connectivity of the porous medium.

III.A.5.e Large Scale Local Porosity Distributions

This section reviews the application of recent results in statistical physics [190, 63, 64, 65, 178, 66, 67, 68, 69] to the problem of describing the macroscopic heterogeneity on all scales. The original definition (3.33) of the local porosity distributions depends upon the size and shape of the measurement cells, respectively on the partitioning $\mathcal{K}$ of the sample space. This dependence on the choice of a test set or “structuring element” is characteristic for many methods of mathematical morphology [10, 37, 71], and many of those discussed in sections III.A.3 and III.A.4 above. On the other hand subsection III.A.5.a has shown that in the limit of large measurement cells $R(\mathbb{K})\rightarrow\infty$ the form of the local porosity distribution $\mu(\phi;\mathbb{K})$ becomes independent of $\mathbb{K}$ and approaches one and the same universal limit given by $\delta(\phi-\overline{\phi})$ . This behaviour is an expression of the central limit theorem. Local porosity distributions have support in the unit interval, hence their second moment is always finite, and average local porosities must become sharp in the limit. It will be seen now that this behaviour is indeed characteristic for macroscopically homogeneous porous media, while other limiting distributions may arise for macroscopically heterogeneous media.

Consider a convex measurement cell $\mathbb{K}$ of volume $V_{0}$ . Let $b$ be the random scale factor at which the pore space volume $V(\mathbb{P}\cap b\mathbb{K})$ of the inflated measurement cell $b\mathbb{K}$ first exceeds $V_{0}$ , i.e. define $b$ as

$b=\inf\{ c\geq 1:V(\mathbb{P}\cap c\mathbb{K})\geq V_{0}\}.$

(3.48)

Consider $N$ mutually disjoint measurement cells $\mathbb{K}_{i}$ $(1\leq i\leq N)$ all having the same volume $V(\mathbb{K}_{i})=V_{0}$ . Let $b_{i}$ denote the scale factors associated with the cells, and let $V_{i}=V(\mathbb{P}\cap b_{i}\mathbb{K}_{i})$ be the $N$ values of the pore space volumes. If the medium is homogeneous there exists a finite correlation length beyond which fluctuations decrease. Then the nonoverlapping cells $\mathbb{K}_{i}$ can be chosen such that the inflated cells remain nonoverlapping, and such that they are separated more than the correlation length. Then the $N$ local porosities $\phi(\mathbb{K}_{i})=V_{0}/V_{i}$ are uncorrelated random variables. For macroscopically heterogeneous media the correlation length may be infinite, and thus it is necessary to consider the limit $V_{0}\rightarrow\infty$ of infinitely large cells to obtain uncorrelated porosities. In [190, 63, 64, 65, 178] the resulting ensemble limit $N,V_{0}\rightarrow\infty$ has been defined and studied in detail. In the present context the ensemble limit can be used to study the limiting distribution of the N-cell porosity

$\phi(\mathbb{K}_{1},...,\mathbb{K}_{N})=\frac{NV_{0}}{V_{1}+...+V_{N}}=N\left[\sum _{{i=1}}^{N}\left(\phi(\mathbb{K}_{i})\right)^{{-1}}\right]^{{-1}}$

(3.49)

obtained from the $N$ measurements. In the ensemble limit the $N$ porosities $\phi(\mathbb{K}_{i})$ become independent but ill defined, and this suggests to consider instead the limiting behaviour of the renormalized sums of positive random variables

$\varphi _{N}=\frac{1}{D_{N}}\left(\sum _{{i=1}}^{N}\left(\phi(\mathbb{K}_{i})\right)^{{-1}}-C_{N}\right)$

(3.50)

where $D_{N}>0$ and $C_{N}$ are renormalization constants. Note that $\phi(\mathbb{K}_{i})^{{-1}}\geq 1$ for all $i$ .

If the sequence of distribution functions of the random variables $\varphi _{N}$ converges in the ensemble limit $N\rightarrow\infty$ and $V_{0}\rightarrow\infty$ then the limiting distribution is given by a stable law [74]. The existence of this limit is an indication of fractional stationarity [63, 64, 65, 178, 68, 69]. The limiting probability density function for the variables $\varphi _{N}$ is obtained along the lines of [64] as

$h(\varphi;\varpi,C,D)=\frac{1}{\varpi(\varphi-C)}H^{{10}}_{{11}}\left(\frac{D^{{1/\varpi}}}{\varphi-C}\left|\begin{array}[]{cc}(0,1)\\ (0,1/\varpi)\end{array}\right.\right)$

(3.51)

where $\varphi\geq C$ and the parameters obey the restrictions $0<\varpi<1$ , $D\geq 0$ and $C\geq 1$ . The function $H^{{10}}_{{11}}$ appearing on the left hand side is a generalized hypergeometric function which can be defined through a Mellin-Barnes contour integral [191]. The limiting local porosity density is obtained as the distribution of the random variable $\phi=1/\varphi$ , and reads

$\mu(\phi;\varpi,C,D)=\frac{1}{\varpi\phi(1-C\phi)}H^{{10}}_{{11}}\left(\frac{D^{{1/\varpi}}\phi}{1-C\phi}\left|\begin{array}[]{cc}(0,1)\\ (0,1/\varpi)\end{array}\right.\right)$

(3.52)

for $0\leq\phi\leq 1/C$ and $\mu(\phi)=0$ for $\phi>1/C$ . The universal limiting local porosity distributions $\mu(\phi;\varpi,C,D)$ depend on only three parameters, and are independent of the diameter or size of the measurement cells $\mathbb{K}_{i}$ . It is plausible that the limiting distributions will also be independent of the shape of the measurement cells, at least for the classes and sequences of convex measurement sets usually employed in studying the thermodynamic limit [192].

Because the limiting distributions have their support in the interval $0\leq\phi\leq 1/C$ all its moment exist, and one has $\overline{\phi^{n}}=\overline{\phi^{n}}(\varpi,C,D)$ . If the moments for $n=1,2,3$ can be inverted then the parameters $\varpi,C,D$ can be written as

$\varpi=\varpi(\overline{\phi},\overline{\phi^{2}},\overline{\phi^{3}}),C=C(\overline{\phi},\overline{\phi^{2}},\overline{\phi^{3}})D=D(\overline{\phi},\overline{\phi^{2}},\overline{\phi^{3}})$

(3.53)

in terms of the first three integer moments. Figure 17 displays the form of $\mu(\phi;1/2,2,D)$ for various values of $D$ .

Figure 17: Universal limiting local porosity density (3.52) for $\varpi=0.5,C=2$ and $D=0.05,0.5,1.0,2.0,5.0,10.0$ . For $D\rightarrow 0$ the density function is concentrated at $\phi=1/C$ , and it vanishes generally for $\phi>1/C$ .

In the limit of small porosities $\phi\rightarrow 0$ the result (3.52) behaves as a power law

$\mu(\phi)\propto\phi^{{\varpi-1}}.$

(3.54)

Within local porosity theory this behaviour can give rise to scaling laws in transport and relaxation properties of porous media [170]. The importance of universal limiting local geometry distributions arises from the fact that there exists a class of limit laws which remains broad even after taking the macroscopic limit $N\rightarrow\infty,V_{0}\rightarrow\infty$ . This is the signature of macroscopic heterogeneity, and it occurs for $\varpi<1$ . Macroscopically homogeneous systems, corresponding to $\varpi=1$ , converge instead towards a $\delta$ -distribution concentrated at the bulk porosity $\overline{\phi}$ .

III.A.6 Capacities

While local porosity distributions (in their one cell form) give a useful practical characterization of stochastic porous media they do not characterize the medium completely. A complete characterization of a stochastic medium is given by the so called Choquet capacities [72, 10]. Although this characterization is very important for theoretical and conceptual purposes it is not practical because it requires to specify the set of “all” compact subsets (see discussion in section II.B.2).

Consider the pore space $\mathbb{P}$ of a stochastic two component porous medium as a random set. Let $\mathcal{F}$ denote the family of all closed sets and $\mathcal{K}$ the set of all compact sets as in section II.B.2 above. For any $\mathbb{K}\in\mathcal{K}$ the “hitting function” or capacity functional is defined as

$T(\mathbb{K})=\mbox{Pr}\{\mathcal{F}_{\mathbb{K}}\}=\mbox{Pr}\{\mathbb{P}\cap\mathbb{K}\neq\emptyset\}$

(3.55)

where Pr is the probability law governing the random set $\mathbb{P}$ . Then $T$ has the following properties: (i) $T(\emptyset)=0$ and $0\leq T(\mathbb{K})\leq 1$ . (ii) If $\mathbb{K}_{1},...,\mathbb{K}_{n},...$ is a sequence of compact subsset then $\mathbb{K}_{1}\downarrow\mathbb{K}$ implies $T(\mathbb{K}_{n})\downarrow\mathbb{K}$ . (iii) For all $n$ the numbers $S_{n}(\mathbb{K};\mathbb{K}_{1},...,\mathbb{K}_{n})\geq 0$ are nonnegative where the $S_{n}$ are determined by the recursion relation

$S_{n}(\mathbb{K};\mathbb{K}_{1},...,\mathbb{K}_{n})=S_{{n-1}}(\mathbb{K};\mathbb{K}_{1},...,\mathbb{K}_{{n-1}})-S_{{n-1}}(\mathbb{K}\cup\mathbb{K}_{n};\mathbb{K}_{1},...,\mathbb{K}_{{n-1}})$

(3.56)

and $S_{0}(\mathbb{K})=1-T(\mathbb{K})$ . The number $S_{n}(\mathbb{K};\mathbb{K}_{1},...,\mathbb{K}_{n})$ gives the probability that $\mathbb{P}\cap\mathbb{K}$ is empty, but $\mathbb{P}\cap\mathbb{K}_{i}$ is not empty for all $i$ . The functional $T(\mathbb{K})$ is called an alternating Choquet capacity of infinite order[10, 71].

Choquet’s theorem says that the converse is also true. Explicitly, if $T$ is a functional on $\mathcal{K}$ , then there exists a necessarily unique distribution Pr in $\mathfrak{F}$ with

$\mbox{Pr}(\mathcal{F}_{\mathbb{K}})=T(\mathbb{K})$

(3.57)

for all $\mathbb{K}\in\mathcal{K}$ if and only if $T$ is an alternating Choquet capacity of infinite order. This theorem shows that capacity functionals play the same defining role for random sets in continuous space as do the numbers $\mu(x_{1},...,x_{N})$ in (2.10) in the discrete case.

The main problem with this theoretically important result is that the family of hitting sets $\mathcal{F}_{\mathbb{K}}$ is much too large for both practical and theoretical purposes.

Author:	R. Hilfer
originally published in:	Advances in Chemical Physics 92, 299-424 (1996)
originally submitted:	May 9th, 1995