|
|
||||||||
a Institute of Terrestrial Ecology, Swiss Federal Institute of Technology, ETH Zurich, Switzerland
b Centre for Non-destructive Testing, Swiss Federal Laboratories for Materials Testing and Research, EMPA, Switzerland
c Spallation Neutron Source Division, Paul Scherrer Institute, PSI, Switzerland
d Institut für Computeranwendungen im Bauingenieurwesen, TU Braunschweig, Germany
5 GKSS-Research Centre, Geesthacht, Germany
* Corresponding author (peter.lehmann{at}env.ethz.ch)
Received 10 December 2004.
| ABSTRACT |
|---|
|
|
|---|
Abbreviations: EPC, Euler-Poincaré characteristic HASYLAB, Hamburger Synchrotron Laboratories, Germany SLS, Swiss Light Source, Paul Scherrer Institute, Switzerland
| INTRODUCTION |
|---|
|
|
|---|
To understand transport processes in porous media, we need to determine the pore structure with a voxel size of a few micrometers. To assess the size of the pores and the pore-size distribution, indirect methods are often used: the solid phase is complementary to the pore space, and some soil properties can be estimated using the measured particle-size distribution (Arya and Paris, 1981; Arya et al., 1999). However, the structure of the porous medium does not depend on the size of the particles alone. Structure is the result of soil formation processes. A living organism, for example, may form a burrow or may produce organic substances that act as cementing agents forming soil aggregates. Such structures cannot be derived or deduced from the sizes of the particles. Another indirect approach to determine the pore sizes is the measurement of the water retention curve. To drain a porous medium, a pressure difference is applied between soil water and air. The water is retained in pores that are small enough to withstand the applied suction. However, this interpretation is partially misleading because a large pore surrounded by small pores will not drain until the suction is large enough to drain the small pores. Therefore, with the applied pressure difference, only the size of the smallest pores will be detected. Vogel and Roth (1998) used a network of capillaries and investigated the effect of network topology on soil water characteristic and hydraulic conductivity to show the importance of the connectivity. Lehmann (2002) showed that the interpretation of the water retention curve as pore-size distribution to predict the hydraulic conductivity holds for soils with a wide pore- and particle-size distribution. But the prediction based on the derivation of the pore-size distribution from the water retention curve failed for soils with a narrow particle-size distribution (sand). For media with a narrow size distribution, the soil water characteristic is strongly affected by the connectivity. Also, the wetting process is affected by network topology because a tiny pore surrounded by large pores will not be refilled before the neighboring large pores are saturated. The soil water retention function may be hysteretic because small structures dominate the drainage and large pores control the wetting process (Mualem, 1974; Stauffer, 1996; Lehmann et al., 1998).
The size and shape of pores can be measured directly by thin-cut analysis with a pixel size of a few micrometers. While it is possible to deduce the connectivity of the pore space using pairs of parallel images (DeHoff et al., 1971), so-called dissectors (Sterio, 1984; Vogel and Kretzschmar, 1996), the resolution in the third dimension of consecutive thin-cuts is limited. To minimize the distance between parallel images, the sample must be impregnated with a resin. This method is time consuming, destructive and limited in the resolution of the third dimension.
To measure the three-dimensional pore network nondestructively with the same resolution in all three spatial dimensions, we must use tomographic techniques. The optimal technique depends on the size of the sample, the required resolution, and the most relevant structures of the sample. So, it makes a difference whether we are interested in the density distribution of the whole sample, the distribution of the water phase alone, or the spatial distribution of certain constituents. X-rays were used to scan the density distribution in soil columns (Hopmans et al., 1994; Tidwell et al., 2000) and to map large pores that are relevant for fast flow processes (Heijs et al., 1995; Vanderborght et al., 2002). Normally, the beam is generated with a device known as an X-ray tube, where a metallic target is bombarded with electrons from a cathode and X-rays are released. Alternatively, very intense X-rays can be emitted by very fast electrons traveling through a magnetic field. To distinguish between this synchrotron radiation and X-rays from tubes, we use the terms X-rays from tubes and X-rays from synchrotrons. The advantage of the synchrotron radiation is the wide range of the wavelength, ranging from ultraviolet light to hard X-rays, and the very intense photon flux density. To scan samples with a resolution of 1 µm or less, a high photon flux per micrometer is needed to limit the exposure time. A comparison between tomography with micro-X-ray and synchrotron radiation was given by Bernhardt et al. (2004).
To analyze fine details of the pore structure of small samples, X-rays from synchrotrons were used to measure the distribution of pore sizes in sandstone (Lindquist et al., 2000) and the porosity in chalk stone (Betson et al., 2004). To determine variations in water saturation, neutron tomography can be applied (Solymar et al., 2003). Also the nuclear magnetic resonance technique is sensitive to the amount of water in the porous media and can be used to measure the size distribution of water filled pores (Bird et al., 2005). In this study, we focus on tomography with X-rays from tubes and synchrotrons and thermal neutrons to measure the size and the arrangement of pores.
In soil physics, we want to understand and predict flow phenomena and transport processes for partially water saturated conditions. Our ultimate goal is the prediction of transport processes for a given soil type based on the information about its structure, that is the spatial arrangement of the pores and particles. In this study, we measured the structure of various porous media and calculated the water and air distribution. The manuscript has three key aspects:
| MEASUREMENTS WITH COMPUTER TOMOGRAPHY |
|---|
|
|
|---|
|
Tomography of Sand Material with X-rays from Tubes
The porosity and particle-size distribution of the sand was determined in the laboratory. The measured porosity was 0.37 m3 m3. The masses of the following particle-size fractions were measured with sieves: 80 to 125, 125 to 250, 250 to 400, 400 to 625, 625 to 800, 800 to 1000, and 1000 to 1250 µm. The mean radius of the measured particle-size distribution was 263 µm. The pore structure was investigated using a micro-X-ray scanner. A sand column 2.5 cm in diameter and 0.9 cm in height was irradiated with a cone beam of energy 140 keV maximum, filtered with 1-mm aluminum. The voxel size of the resulting image was 60 µm. The pore structure of the sand sample after segmentation is shown in Fig. 1A. A threshold was applied to classify the voxels as pores or particles with the constraint of the measured porosity. To scan a large (5.3-cm diam., 7-cm height) column of the same sand material, we used an industrial tomographic scanner with a 450-kV X-ray tube (focal spot size of 2.5 mm). The beam is collimated to a fan of 2-mm height. The height of the scanned sector was 5.3 cm, corresponding to 750 scans with a spacing of 70 µm. The resolution was 70 by 70 by 210 µm3, because 210 µm was the minimum opening of the collimator to receive a signal above the noise level. To enhance the resolution in the z-direction, the vertical increments between two measured planes were 70 µm. Thus, the resulting image had a grid spacing of 70 by 70 by 70 µm3, but the gray value of a voxel depended on the material density in a cuboid of size 70 by 70 by 210 µm3. Before we tried to calculate the material density in a cube of size 70 by 70 by 70 µm3, the gray-scale value images had to be filtered due to a bad signal/noise ratio. Each voxel is affected by the material in the plane below and in the plane above. We faced the problem that for the first and the last of the scanned planes, the measurements were influenced by neighboring planes whose properties are unknown. For the reconstruction of N planes, a system with N equations with N + 2 unknowns must be solved. For each position x, y of the image, a sequence of gray values as a function of the plane in height z was obtained with tomography. This sequence is denoted as g(z). From g(z), we want to derive the unknown gray values in cubic voxels, v(z). The measured gray value may thus be described as g(z) = h(1)v(z 1) + h(0)v(z) + h(1)v(z + 1), with a function h that defines the contribution of the planes to the measured gray value. In case of industrial X-ray tomography, the contribution of the three planes to the detected signal was equal; therefore, h(1) = h(0) = h(1) = 1.0. In mathematical terms, this is a convolution that can be described as a multiplication after the Fourier transformation of the functions, G(
) = V(
) H(
), with the space frequency
.
The goal was to find an optimized response function h* with its Fourier transform H*(
), then to calculate the inverse Fourier transform of G(
)/H*(
) and to use this function as an estimate for the unknown v(z). The response function h*(
), where
is the shift with respect to the position z, is non-zero for h*(0), h*(1), and h*(N 1). The argument N 1 equals a shift of 1. We tested different functions h*(
), transformed them into the Fourier domain, divided G(
) with H*(
), and calculated the inverse transformation of G(
)/H*(
). With an adequate choice of H*(
), the effect of the unknown boundary values can be filtered out, and a valuable estimate for the unknown function v(z) can be calculated. Lehmann et al. (2001) estimated the function h*(
) by analyzing a sequence of random numbers v(z). Here, we used numerically generated sand media with the same particle-size distribution and porosity as the real one. This procedure is specified in the section "Sand Material with Wide Particle-Size Distribution." Cross sections through the images are shown in Fig. 2
.
|
Tomography of Sand Material using Synchrotron Radiation
The synchrotron X-rays of an electron storage ring (Swiss Light Source [SLS], Paul Scherrer Institute, Switzerland) and a positron storage ring (Hamburger Synchrotron Laboratories, HASYLAB, Germany) were used to investigate the pore space of sand samples in detail. Tomographic scans of a sand sample (1.5-cm diam.) containing 300- to 900-µm sand particles were performed at beamline W2 of the HASYLAB at the Deutsches Elektronen-Synchrotron DESY using a photon energy of about 42 keV and a voxel size of about 11.5 µm. At SLS a sand sample (0.5-cm diam.) containing 100- to 200-µm particles was scanned with a monochromatic energy of 16 keV, with a resulting voxel size of 3.5 µm. To classify the voxels as pores or solid phase, two thresholds were applied. Voxels with gray values smaller than the first were designated as pores, and voxels with gray values higher than the second threshold were designated as solid phase. In a second step, voxels with a gray value between the two thresholds were classified as a pore if one of the six neighboring voxels was already designated as pore. This second step was repeated until no further voxels with intermediate gray values were detected in the neighborhood of a pore. The remaining gray voxels were designated as solids. A section of the sample measured at SLS is presented in Fig. 1B.
| DESCRIPTION OF SIMULATED MEDIA |
|---|
|
|
|---|
Numerical Generation of Porous Media
Medium with Particles of Different Sizes
The particle diameters ranged from 80 to 1250 µm with a mean radius of 263 µm. The mass of seven size fractions was determined in the laboratory by sieving. We summarized these mass fractions to get the measured cumulative particle size distribution function (CDF). We chose a lognormal distribution m(r) as a model for the particle size distribution as a function of the particle radius r. The cumulative distribution of this analytical function m(r) (i.e., an error function with the logarithm of the particle size in the argument) was fitted to the measured CDF. The fitted lognormal distribution m(r) was divided by the mass of a spherical particle
4/3
r3 to obtain the number of particles per particle size p(r), with the constant density of the particles
. The function p(r) was integrated over the particle size r and normalized to obtain values between 0 and 1.0. This integral, in the following denoted as P(r), determines the probability of finding a particle of size equal to or smaller than r. To determine the size of the modeled particles, a random number between 0 and 1 was chosen. This random number n corresponds to P(r) for a certain particle radius r and a spherical particle with radius r was generated. This procedure was repeated until the volume of the generated particles was equal to the measured volume fraction of the solid phase (1
) L3, with the measured porosity
= 0.37 m3 m3 and the size of the modeled system L. To model the media analyzed with micro-X-ray tomography, 15 139 spherical particles were randomly arranged in a cube of size L = 6.25 mm. In the case of the industrial X-ray tomography, 99 870 spheres were distributed in a cube of size L = 12.5 mm. In a second step of the numerical packing procedure, the particles were distributed randomly in space, starting with the largest particles. If one particle overlapped another, then a new random position was chosen. In literature, such an approach is denoted as a random-parking algorithm (Cooper, 1988).
Medium with Single Particle Size
In the case of a single particle size, media with small porosities with nonoverlapping particles cannot be generated with a parking algorithm. The porosity of the glass bead packing was 0.40 m3 m3. The code we used was based on the description of the algorithm given by Jodrey and Tory (1981). The spherical particles were randomly distributed and the overlapping of the particles was allowed in a first step. In a next step, all overlapping particles were rearranged by moving the particles in a direction opposite to the center of the nearest neighbor. This procedure of redistribution was repeated until a system with nonoverlapping particles was generated. To optimize the redistribution of overlapping particles, the size of the spheres was slightly enlarged for the initial distribution of overlapping spheres. The initial size was 100.2% of the true size. When the redistribution of the particles became less effective with respect to overlapping, the size of the particles was decreased. To model the tomography of the glass bead columns, 9168 spheres of 1-mm radius were distributed in a cubic box of size 40 mm.
With this method, packings with higher densities can be generated compared with the algorithms described in the previous subsection. However, for a porosity of 0.37 m3 m3 and the wide particle size distribution of the sand, the random parking algorithm was more efficient.
Analyzed Properties of Simulated Media
We mapped the generated structure onto a grid. The grid spacing corresponded to the voxel size of the resulting image. For each voxel of the medium we computed the volume fraction of the solid phase and determined a gray-scale value between 0 (i.e., voxel contains no solid) and 255 (i.e., voxel contains no voids). Then, a threshold was applied to obtain an image of the pores with the measured porosity. An example of the resulting density image (gray) and the corresponding image of the pore space (black and white) is given in Fig. 3
. We calculated the density map and the image of the pore space and determined the following quantities.
|
Minkowski Functionals
To determine the morphology of a structure, the four Minkowski functionals of the solid phase were calculated. They include the volume, surface, integral mean curvature (or mean breadth), and the Euler characteristic. The Minkowski functionals characterize the morphology of disordered media (Arns et al., 2002). The mean breadth is proportional to the integral mean curvature and depends on the curvature of the object. If the center of the curvature is within the analyzed object, the curvature is positive. In the opposite case, the value is negative. For isolated convex objects, the mean breadth corresponds to the mean diameter. The integral mean curvature equals the mean breadth multiplied with 2
. For a complex network containing convex structures, the mean breadth becomes negative. The Euler characteristic, also known as the Euler-Poincaré characteristic (EPC), quantifies the connectivity of a structure. It is negative for highly connected networks with redundant connections and positive for a set of isolated elements. The EPC for the solid phase is high for spheres that are completely separated and pores that are embedded within the solid phase. In the case of an imaged structure, the structure consists of cubic voxels. Each voxel consist of one cubic body, six faces, 12 edges, and eight vertices. The Minkowski functionals of an imaged structure can be obtained by counting the total number of cubes nc, faces nf, edges ne, and vertices nv according to the following formula (Michielsen and De Raedt, 2001):
If two neighboring voxels have a face, an edge, or a vertex in common, this element is counted only once. For the surface area (and the mean breadth) of spheres this formula is not exact because it is impossible to describe the smooth surface of spherical solid particles by counting the quadratic faces of the voxels. Each voxel of the solid phase has six faces, and each of the six faces is the boundary between two voxels. If one of the neighbored voxels belongs to the pore space, this face with an area of
L2 is an element of the surface of the solid phase, where
L is the length of a voxel. A more accurate estimation of the area was computed with a method described in Ohser and Mücklich (2002, p. 114121). For each voxel of the solid phase at position x, y, z they included 26 neighbors at position x +
x, y +
y, z +
z, with the values 1, 0, or 1 for
x,
y, and
z. If the voxel at a position x +
x, y +
y, z +
z belongs to the pore space, an interface of wi is added to the total interface of the solid phase. The coefficients wi depend on the distance between the two voxels and are equal to 0.183
L2 for w1 with
x2 +
y2 +
z2 = 1, 0.105
L2 for w2 with
x2 +
y2 +
z2 = 2 and 0.081
L2 for w3 with
x2 +
y2 +
z2 = 3. For the mean breadth, the method is similar but with
L instead of
L2 and with the coefficients divided by 2. Because we generated the media with spherical particles, we will use this more accurate estimation.
Pore-Size Distribution
For each voxel determined as pore, the radius of a sphere that fits entirely into the pore space was calculated. The pore voxel is the center of the sphere. First, the radius of this sphere was attributed to the pore voxel. This result equals the distance transform (Soille, 2003), determining the minimum distance of each pore voxel from the solid wall. However, a pore voxel can be an element of a larger sphere. In such a case, the radius of the larger sphere is assigned to this voxel. The fraction of pore voxels that belongs to the different sphere sizes was computed. Another method to compute the pore size distribution was applied by Vogel (1997). He calculated the pore size distribution with a morphological filtering denoted as opening (Serra, 1982) without the preceding determination of the distance transform.
| GEOMETRICAL PROPERTIES AS A FUNCTION OF RESOLUTION |
|---|
|
|
|---|
Sand Material with Wide Particle-Size Distribution
The density distribution in the generated media was calculated for voxel sizes ranging from 6.25 to 125 µm, corresponding to a grid with 10003 to 503 voxels. To compare the resolutions for media with different particle sizes, we define a dimensionless resolution R by dividing the voxel size with the mean particle radius.
The fraction of the solid phase was calculated for each voxel resulting in a gray value between 0 and 255. The image was segmented into a map of pore space and solid phase with the required porosity of 0.37 m3 m3. Cross sections through the density map and the pore space image are shown in Fig. 4 for different resolutions.
|
Next, we analyzed the Minkowski functionals for different resolutions. Due to the condition that the fraction of pore space must correspond to the measured porosity, the volume of the solids is constant. The surface of the spheres in the image is compared with the analytically determined surface of all spheres. The surface increases with decreasing voxel size (Fig. 5A ) and converges to the true value. In Fig. 5B, the change of surface area with voxel size is given. For voxel sizes <30 µm (R = 0.11), the increase of surface with decreasing voxel size becomes less negative, indicating the convergence to the true value. The calculated value of the mean breadth was divided by the analytically determined value of the spheres (Fig. 6A ). At high resolutions, the mean breadth converges to the true value. The mean breadth becomes positive for resolutions better than 60 µm (R = 0.23). At lower resolutions, the solid phase is a complex network with inclusions of pore voxels. For high resolutions, the solid phase consists of separated convex spheres. For the latter case, positive breadth values are characteristic, and negative numbers result in the first case. The Euler characteristic (EPC) divided by the number of particles was computed for the solid phase (Fig. 6B). At poor resolutions, the spheres are connected, as indicated by the negative EPC. With increasing resolution, the spheres are separated and the EPC becomes positive for voxel sizes <20 µm (R = 0.08). For the best resolutions, the EPC converges to the value 1.0, with all spheres separated.
|
|
30 µm. Hence, the fraction of large pores is underestimated with the voxel size of the measurement, and the mapped pore-size distribution deviates from the true function. The true function (i.e., the pore-size distribution mapped with the maximum resolution) is characterized by a small volume fraction of small and large pores and large volume fraction of intermediate pore sizes. Therefore, the second derivative of the cumulative pore-size distribution must change from positive to negative values with increasing pore size. This criterion is fulfilled for resolutions
30 µm.
|
|
)/H*(
), denoted as v*(z), can be compared with v(z). We fitted the function H*(
) to v(z) and applied the same procedure with the same function H*(
) for the measured gray values to estimate the density distribution in cubic voxels. To optimize the function H*(
), we applied thresholds to v(z) and v*(z) to obtain segmented images of the pore space with the same porosity as the measured medium. Then we counted the number of voxels with different phases in the two computed media. The smallest error was found for the values 0.70, 0.67, and 0.60 for h*(1), h*(0), and h*(1). The error in the segmented image does not depend on the absolute values but only on the ratios between h*(1), h*(0), and h*(1). With higher absolute values for h*, the range of the gray values becomes larger. The results for the generated and the scanned media are shown in Fig. 2.
Glass Bead Medium with Single Particle Size
To investigate numerically the effect of resolution on geometric properties of a sample with particles of equal size, 9168 spheres with a radius of 1 mm were distributed in a 40-mm cube. The fraction of glass in each voxel was calculated and expressed in terms of a gray-scale value between 0 and 255 for voxel sizes ranging from 40 to 500 µm, corresponding to images of the size 10003 and 803 voxels. Cross sections through the images are shown in Fig. 4.
The fraction of voxels determined as pores with non-zero gray values decreased with increasing resolution from 92% (voxel size 500 µm) to 11% (voxel size 40 µm). For the voxel size corresponding to the resolution obtained with thermal neutrons (167 µm, R = 0.17), 56% of the voxels determined as pores had a gray value of 0 (Fig. 3). For voxel sizes <200 µm (R = 0.20) and 70 µm (R = 0.07), the fraction of pore voxels with non-zero density values was <50% and 20%, respectively.
The calculated area divided by the analytically determined value converges to 1 with decreasing voxel size. For large voxel sizes, neighboring spheres overlapped in the image, and the surface of the intersection will be neglected. Therefore, the surface decreases with increasing voxel size (Fig. 5A). The change of surface with decreasing voxel size became constant for voxel sizes <150 µm (R = 0.15). This change of surface with voxel size differs from the function calculated for the sand. This is caused by the slightly different numerical packing procedure. In the case of the glass beads, almost all spheres are in contact with other spheres. Even for highly resolved images, the spheres are not completely separated, and the area of the intersection is neglected. In the case of the spheres of the sand model, fewer spheres are in contact with other particles. The most spheres become isolated for small voxel sizes, and the surface area is not affected by a further increase of resolution. The mean breadth of the solid phase (Fig. 6A) is positive for voxel size <250 µm (R = 0.25). Compared with the sand model, the values for high resolutions are smaller because the spheres touch each other, and these connections decrease the curvature. Then the Euler characteristic of the solid phase was determined (Fig. 6B). With resolutions better than 90 µm (R = 0.09), the spheres become separated and the EPC was positive. For high resolutions, the EPC divided by the total number of particles was higher than the expected value of 1. This deviation was caused by the inclusion of pore voxels in the vicinity of touching spheres. As mentioned above, pores enclosed in the solid phase increase the EPC. In the case of the packing algorithm used to model the glass beads, the most spheres touch other spheres because they were initially overlapping. In the discretized image, the surface of the spheres is rough, and it is possible that a pore voxel can be embedded between the voxels of two spheres that are very close. In the case of the packing algorithm used to simulate the sand media, there is only a small chance that a sphere touches another one and fewer voxels of the pore space are enclosed.
The pore-size distribution as a function of voxel size is shown in Fig. 7B. With the resolution of the tomography with thermal neutrons (167 µm), the shape of the curve corresponds to that mapped with the highest resolution. The second derivative of the cumulative pore-size distribution measured with a voxel size <200 µm is positive for small pores. The mean pore radius computed with the highest resolution is 319 µm (R = 0.32). Thus, the voxel size of the tomographic measurement was about 50% of the mean pore radius. In the case of the sand, the voxel size was 150% of the mean pore radius.
To summarize, the properties of the mapped pore space turn toward the true values for voxel sizes between 10 and 20% of the particle radius. The resolution required to map the pore space of the glass bead mixture must be smaller than 100 to 200 µm, or less than 10 or 20% of the mean particle radius (Table 1). With the resolution of the tomography with thermal neutrons (167 µm) this condition was partially fulfilled.
| APPLICATION: FLUID PHASE DISTRIBUTION |
|---|
|
|
|---|
Measured Fluid Phase Distribution in Glass Bead Column Scanned with Thermal Neutrons
For thermal neutrons, the attenuation is dominated by the presence of water, so this method can be applied to measure the water distribution. Comparing the images of a dry and a wet sample, the differences due to the presence of water can be quantified. Often, heavy water (D2O) is used instead of normal water, because the attenuation of the neutrons is smaller for heavy water. The intensity of the beam attenuated by a system filled with normal water could be too small to be used for the reconstruction of the structure. We scanned a glass bead column (4-cm diam.) containing 3-mm glass beads. The bottom of the column was filled with heavy water. The reconstruction of the system and the boundary between the water saturated and dry zone is shown in Fig. 1D. For X-rays, the attenuation is density dependent and is smaller for water than for quartz sand. For this reason, the differences between the wet and a dry sample images are small, and it is difficult to distinguish between the two phases. To enhance the contrast, a potassium solution could be used instead of normal water (Hirono et al., 2003).
Fluid Phase Distribution Calculated with Lattice-Boltzmann Approach
A sample of diameter 25 mm with sand grains ranging from 80 to 1250 µm was scanned with X-rays from tubes with a resolution of 60 µm. For a reconstructed sand cube of length 3 mm (Fig. 8A
), the water and air distribution was calculated using a LatticeBoltzmann approach. The method describes the dynamics of fluid particle distributions on a regular grid. Within each time step, particle distributions are modified by a simplified collision process and propagate to a finite set of neighboring nodes. The definition of the collision rules ensures that the dynamics of the moments of the particle distributions are equivalent to the dynamics of the variables of Navier-Stokes and the continuity equation up to second order. Thus, the Lattice-Boltzmann approach can be utilized to obtain approximate solutions of the Navier-Stokes equations. Details of this method are given, for example, in Wolf-Gladrow (2000) and Krafczyk (2001). In addition, the interface between the two-fluid phases can be computed with an extended approach (Tölke et al., 2002; Sukop and Or, 2004). For that purpose, the dynamics of two probability distributions representing the air and the water phase are computed.
|
The simulation of the full two-phase dynamics in a porous medium is a demanding task in terms of computational resources. To calculate the fluid phase distribution in larger systems, we used another approach that we denote as a "morphological pore network model." We introduce this approach in the following subsection.
Fluid Phase Distribution Calculated with Morphological Pore Network Model
In network models, the pore space consists of elements with well-defined size. Most often, the pore space is described as a composite of pore bodies connected with pore throats. Size, shape, and arrangement of these elements can be adapted to mimic the fluid distribution within a porous media. In contrast to such a conceptual representation of the porous medium, the measured image of the pore structure is used in case of morphological pore networks. The "true" pore space is subdivided into elements of well-defined size. In the approach used in this study, each pore voxel is an element of the pore space and its size is the radius of a sphere as explained. The size of each pore voxel will be needed to calculate the water and air distribution. Because the whole geometrical information is included in this network approach, we denote this type as morphological network. A similar approach to map the size of the pores and to calculate the fluid phase distribution in numerically generated sphere packings was applied by Hilpert and Miller (2001). Here, we use the morphological network to calculate the water and air distribution in sand media measured with X-rays from synchrotrons. We calculated the air and water distribution for a coarse sand with particle sizes ranging from 300 to 900 µm scanned at HASYLAB with voxel size 11.5 µm and for a fine sand with particles ranging from 100 to 200 µm scanned with a voxel size of 3.5 µm at SLS. The voxel sizes correspond to 4% of the mean particle radius (75 and 304 µm). So, the requirement of a voxel size smaller than 10% of the mean particle radius was fulfilled. In a numerical drainage experiment, we calculated the water and air distribution for a suction applied at the lower boundary. We used the three-dimensional segmented image and computed the size of each pore voxel as explained above. In Fig. 9A
and 9B, the calculated pore sizes for the voxels of the pore space are shown for two scanned sand samples. For the sand shown in Fig. 9A, a small section is enlarged to explain the determination of the pore size in more detail. First, the distance from the next solid was calculated for each voxel of the pore space (Fig. 9A1). This distance corresponds to the radius of a sphere within the pore space. A voxel can be an element of different spheres and the radius of the largest sphere determined the size of the voxel. This is shown in Fig. 9A2.
|
|
| SUMMARY AND CONCLUSIONS |
|---|
|
|
|---|
Based on the analysis of numerically generated porous media, the voxel size must be smaller than 10 to 20% of the mean particle radius to reveal the geometrical properties of the pore space. This criterion is based on the investigation of the effect of resolution on the pore-size distribution, the surface of the particles, the breadth and the Euler characteristic of the solid phase. For a voxel size of 10 to 20% of the mean particle radius, these geometric properties turn toward true values. In addition, more than 50% of the voxels denoted as pores are entirely within the pore space for this resolution. The requirement of a voxel size smaller than 10 to 20% of the mean particle radius was not fulfilled for the tomography of a sand sample with particle diameters ranging from 80 to 1250 µm using X-rays from tubes. The resolution of the tomography divided by the mean particle size (263 µm) was 0.23. For the tomography with thermal neutrons, the resolution divided by the particle radius of 1000 µm was 0.17, and the requirements were partially fulfilled. For the tomography of two sand samples with particles of diameter 100 to 200 and 300 to 900 µm using synchrotron radiation, the resolution was about 4% of the mean particle radius (75 and 304 µm), and all requirements were fulfillednote that this analysis is valid for the quantification of pores between particles. In the case of worm burrows, the pore structure is independent of the size of the particles, and other models of structure generation must be applied to analyze such structures and to define similar requirements.
We used the images of dry sand samples to calculate the water and air distribution for different boundary conditions. With the images obtained by X-rays from synchrotrons, we could predict the drainage of the sand sample in absolute terms with a morphological pore network model. In the case of the sand scanned with X-rays from tubes, the hysteresis phenomenon could be simulated with a Lattice Boltzmann approach.
By combining the appropriate technique of tomography with the numerical tools to model multiphase flow, new insights in the processes at the pore scale and more accurate predictions of the hydraulic properties of soil materials can be expected.
| REFERENCES |
|---|
|
|
|---|
This article has been cited by other articles:
![]() |
A. Carminati, D. Vetterlein, U. Weller, H.-J. Vogel, and S. E. Oswald When Roots Lose Contact Vadose Zone J., August 11, 2009; 8(3): 805 - 809. [Abstract] [Full Text] [PDF] |
||||
![]() |
S. E. Oswald, M. Menon, A. Carminati, P. Vontobel, E. Lehmann, and R. Schulin Quantitative Imaging of Infiltration, Root Growth, and Root Water Uptake via Neutron Radiography Vadose Zone J., August 13, 2008; 7(3): 1035 - 1047. [Abstract] [Full Text] [PDF] |
||||
![]() |
T. Sander and H. H. Gerke Noncontact Shrinkage Curve Determination for Soil Clods and Aggregates by Three-Dimensional Optical Scanning Soil Sci. Soc. Am. J., August 9, 2007; 71(5): 1448 - 1454. [Abstract] [Full Text] [PDF] |
||||
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
| HOME | HELP | FEEDBACK | SUBSCRIPTIONS | ARCHIVE | SEARCH | TABLE OF CONTENTS |
| The SCI Journals | Agronomy Journal | Crop Science | |||
| Journal of Natural Resources and Life Sciences Education |
Soil Science Society of America Journal | ||||
| Journal of Plant Registrations | Journal of Environmental Quality |
The Plant Genome | |||