Journal of Environmental Quality 31:470-476 (2002)
© 2002 American Society of Agronomy, Crop Science Society of America, and Soil Science Society of America
Article
SYMPOSIUM PAPERS
Flow and Diffusion Measurements in Natural Porous Media Using Magnetic Resonance Imaging
Thomas Baumann*,a,
Rainer Petschb,
Gunther Feslb and
Reinhard Niessnera
a Institute for Hydrochemistry, Technical University of Munich, Marchioninistr. 17, D-81377 München, Germany
b Dep. of Neuroradiology, Ludwig-Maximilians-University Munich, Marchioninistr. 15, D-81377 München, Germany
* Corresponding author (thomas.baumann{at}ch.tum.de)
Received for publication June 2, 2000.
 |
ABSTRACT
|
|---|
Flow and diffusion of water in natural porous media, quartz sand, and calcareous gravel were measured using a 1.5-T clinical magnetic resonance tomograph. The spatial resolution of the dynamic measurements was 1.32 x 1.32 x 5 mm3, and the time between two cross-sectional measurements was approximately 10 s. The measured coefficients of molecular diffusion for water were in good agreement with theoretical data. Flow was measured without any tracer at velocities between 0.15 and 6.67 mm/s. The results, based on a calibration within one part of the column, were in good agreement with data obtained from a tracer experiment and from a numerical model. It was possible to measure the flow velocity in larger pores and preferential flow paths directly. The results of the flow measurements in smaller pores reflected the mean velocity within that volume element. In that case the obtained values were close to the average linear velocity. Since the time resolution is high a monitoring of flow processes is possible. The pore space was imaged with a spatial resolution of 0.5 x 0.5 x 0.5 mm3. Here, the porosity of pores that are larger than 0.2 mm can be measured directly; for smaller pores a calibration is necessary.
Abbreviations: EPI, echoplanar imaging sequence • MR, magnetic resonance • MRI, magnetic resonance imaging • RF, radio frequency • S/N, signal to noise ratio • TE, echo time • TR, repetition time
|
INTRODUCTION
|
|---|
THE FLOW VELOCITY in natural porous media is of high relevance to contaminant and colloidal transport in subsurface environments. The transport of dissolved contaminants is subject to sorption and desorption processes with a defined sorption kinetic (Wu and Gschwend, 1986; Song et al., 1994), and the flow velocity at the pore scale determines the extent of dynamic sorption processes. The spatial distribution of dissolved contaminants is again a function of the flow velocity on a pore scale (Corapcioglu et al., 1997). The different flow velocities and different flow paths are the reason for dispersion. In a three-phase system, consisting of fluid, colloids, and stationary matrix, the colloids may facilitate contaminant transport (McCarthy and Zachara, 1989). Apart from surface interactions, the transport of the colloids is depending on the pore space itself (pore-size distribution, connectivity) and on the flow velocity within the pores (Harvey and Garabedian, 1991; Rehmann et al., 1999). Preferential flow paths (e.g., wormholes, fractures, or sediments with a higher flow velocity) contribute to the transport behavior of dissolved and colloidally bound contaminants.
From the physical description of the transport processes it becomes quite clear that the average linear velocities, as obtained from, for example, tracer experiments, are an approximation. The linear distance between the injection and monitoring locations is always the shortest possible distance and does not account for the actual flow path, nor for the different flow velocities in larger and smaller pores, or for the parabolic flow profile in a pore. These factors are summarized as mechanical dispersion. Still, mechanical dispersion is calculated without knowledge of the length of the actual flow path. Thus, the flow velocities calculated from tracer experiments are always slower than the actual flow velocities at the pore scale. The differences become relevant, if a numerical model includes the sorption kinetic measured in laboratory experiments or derived from thermodynamic databases, or if filtration theory is applied for colloid transport phenomena. In those cases the knowledge of the flow velocity at the pore scale would be useful to scale laboratory experiments to real-world problems. There is also need for (near-) real-time measurements of flow and diffusion and the spatial evolution of the pore space to investigate, for example, the obstruction of pore spaces by colloids (Veerapaneni and Wiesner, 1997) or the development of preferential flow paths in the time domain.
Among the techniques for a spatial description of the pore space, thin sectioning and light microscopy are quite common (Vogel, 1997; Vogel and Roth, 1998). The use of X-ray computer tomography (CT) is reported for the characterization of the pore space of soils (Reinken et al., 1995) and glass beads (Jasti et al., 1993). In the context of oil exploration, CT was used to determine the pore space (Akin et al., 1996) or to characterize the swelling behavior of shales (Onaisi et al., 1993). Computer tomography has a good spatial resolution, down to 5 µm in high-resolution mode (Reinken et al., 1995). However, CT measurements of flow processes have several drawbacks. Primarily, CT images the sediment matrix, but the fluid itself yields almost no signal and contrast-enhancing fluids may influence the flow by changes in fluid density or fluid viscosity. Second, CT is not able to monitor flow phenomena directly; instead, the flow measurements are based on an evaluation of time series images. Finally, the sensitivity of CT for flow phenomena in small pores is low, therefore only preferential flow paths can be assessed.
To overcome these limitations another technique was established to measure the flow within porous media directly. 1H–nuclear magnetic resonance (NMR) tomography, also known as magnetic resonance tomography (MRT) or magnetic resonance imaging (MRI), is described in principle by Becker (1980), Bovey (1988), and Skoog and Leary (1992). Magnetic resonance imaging has become a versatile technique in the field of medical imaging (Morneburg, 1995). Besides medical imaging, MRI has been used for the monitoring of exchange processes in plants (Rokitta et al., 1999). Applications of NMR microscopy are given by Callaghan (1991). Nuclear magenetic resonance tracers (Cu, Ni) were used to visualize flow in packed glass beads (Oswald et al., 1997). As these metals underlie sorption and desorption a noninvasive monitoring of contaminant transport becomes possible. Magnetic resonance imaging is sensitive to the 1H protons. Since water contains many protons, MRI seems ideally suited to image the distribution of water within porous media directly. There also have been attempts to monitor the movement of oil and hydrocarbons using magnetic resonance imaging (Chudek and Reeves, 1998).
In this study, flow and diffusion processes in natural porous media are measured using a clinical 1.5-T magnetic resonance scanner with a time resolution of less than 10 s and without a tracer.
|
MATERIALS AND METHODS
|
|---|
Materials
The column, constructed from polymethyl methacrylate, was 20 cm in length and 10 cm in diameter. All tubes and fittings, except the pump tubings (Tygon; Ismatec, Wertheim, Germany) were made of PTFE. Seven distinct layers of natural sediments were packed carefully into the column. The bottom-most layer was composed of medium gravel (Sediment Type D), the uppermost layer was made of coarse sand (Sediment Type C). At the center of the medium sand layer (Type A) a core of coarse sand was installed to obtain different flow velocities within that layer. Figure 1
shows the orientation, distribution, and thickness of the various sediment layers within the column.
Sediments A to C were provided by E+M Bohr GmbH (Hof, Germany), while Sediment D was extracted from the Munich gravel plain. The properties of the sediments are given in Table 1. The sediments were classified according to DIN 18123 (1983), while hydraulic conductivities were determined according to DIN 18130 T 1 (1989). The porosity of each layer was measured gravimetrically, and the grain density 
s was determined according to DIN 18124 (1989).
A tracer experiment was performed with NO-3. After the column was saturated with a 10 mmol/kg NaCl solution, flow was established at a flow rate of 0.047 mL/s from the bottom to the top of the column with a low-pulsation peristaltic pump (Ismatec). Then, 4 mL of a 1.35 g/L NaNO3 solution were injected as a tracer pulse and measured with a DU 650 UV-Vis spectrometer (Beckman, Mülheim, Germany) at a wavelength of 220 nm.
The breakthrough curve (BTC) showed a bimodal development that is caused by different flow velocities in the center and in the outer shell of Layer 6. A numerical two-dimensional, radially symmetric, ground water flow model was calibrated with the hydraulic properties given in Table 1 and data from the tracer test to obtain the flow velocities in each layer. During the magnetic resonance (MR) experiments flow was established between two chambers with constant head and different water levels 
h (see Fig. 1).
Magnetic Resonance Imaging
Under idealized conditions the signal intensity in a MR experiment (i.e., the proton spin density) is proportional to the number of protons in the voxel, which is the smallest three-dimensional grid element that can be measured individually. The porosity ne of a given slice can then be calculated from the proton spin density I with Eq. [1]:
 | [1] |
where N is the number of voxels and I0 is the intensity of a voxel that is completely filled with water. Within narrow pores the spin–spin relaxation times (T2) are influenced by matrix effects (e.g., pendular water on the grain surfaces) and the relation between proton spin density and porosity becomes nonlinear. In that case a calibration procedure is necessary to measure the porosity.
The principles of the measurement sequence used for flow measurements are given by Stejskal and Tanner (1965), based on the classical spin echo sequence. At the beginning of the measurement the spins are oriented along the main magnetic field (z axis). With a 90° radio frequency (RF) pulse the spins are oriented along the y axis. Assuming a constant gradient along the x axis the spins will start to dephase along that gradient. With a second 180° "refocusing pulse" at t = TE/2 (where TE is echo time) the magnetic moments are rotated over the x axis while continuing to move along the x gradient. At t = TE the spins are in phase and oriented along the negative y axis, so the sum of their magnetization vectors reaches a maximum. The signal induced in the receiver by this magnetization is called "spin echo."
It is important to understand that the spin echo as described above requires a static system for the refocusing pulse to be effective. That means that the location of the spins has to be constant in the time domain. Diffusion causes a random displacement of the spins. As a consequence at t = TE spins are not completely in phase and aligned along the negative y axis. Since the diffusion is random, only a decrease of the signal intensity of the spin echo but no clear deviation from the y axis can be observed. Flow without diffusion causes a uniform lateral displacement of the spins. At t = TE the spins will all be in phase, but the direction of the phase is different from the negative y axis. This phase difference can be used for the measurement of flow processes. In natural systems both flow and diffusion are present and have to be included in the evaluation.
Conventional flow sensitive sequences used in medical imaging have a sensitivity for flow velocities of 10 mm/s or faster, reflecting the physiological range of flowing blood or cerebrospinal fluid. However, the flow velocities within porous media are usually in a range of 0.2 mm/s or lower. To achieve the phase shifts needed to measure the flow of water in porous media, the gradients in such sequences have to be applied very intensively over a long time period. In that case the sequences are sensitive for both flowing water and diffusion processes in the water itself. In clinical applications of MRI similar sequences have been used to examine stroke patients (Neumann-Haefelin et al., 1999).
The difficulty is that both flow and diffusion will alter the MR signal. Therefore, a formula has to be found that describes both. First, the signal intensity I of a MR diffusion experiment is described by:
 | [2] |
where 
equals 2
0B0, 0 is the Larmor frequency (42.58 x 109 Hz), B0 is the magnetic field strength, G is the amplitude of the gradient (i.e., the RF pulse), T is the duration assuming a rectangular shape, 
is the time between the 90° gradient and the 180° gradient, and D is the coefficient of molecular diffusion. The factor b is specific to the sequences used, so that only I0 and D are unknown. Second, the phase shift 
caused by a flow velocity v is:
 | [3] |
If the phase shift << /2 (i.e., the flow velocities are low as in the column experiment), both effects can be combined in one formula for the signal intensity:
 | [4] |
where b is the constant specific to the sequence, and a is a regression factor, which in the case of this study was determined with a least-squares fit of the measured data from Layer 3 (arbitrarily chosen) at five different flow velocities to the data obtained from the tracer experiment and the numerical model.
The measurements were done with a clinical MR tomograph (Magnetom Vision, Siemens Medical Engineering, Erlangen, Germany) with B0 = 1.5 T and an inner diameter of 0.8 m. A combined transmitting and receiving coil (CP Head, Siemens) was used. The spatial pore structure and the porosity of the sediments were measured using a standard Siemens three-dimensional high-resolution T2-weighted sequence (CISS, see Nitz, 1999), with TE = 5.9 ms and TR (repetition time) = 12.25 s. The column was measured in three series of 150 slices each (thickness = 0.5 mm) and an overlap of approximately 15 mm. This overlap is necessary for a subsequent correction of the intensities. The data matrix for the column was then 256 x 256 x 450 voxels. Pictures were received as ACR-NEMA files and evaluated with common statistical and imaging packages.
To calculate the coefficient of molecular diffusion from Eq. [2] the same slice was measured with three different amplitudes of the gradient (b = 50, 500, 100), with TR = 0.6 s, TE = 69 ms, and slice thickness remaining unchanged.
Two T2-weighted echoplanar sequences (EPI0, EPIdiff) were used for the combined flow and diffusion experiments. The two sequences differ only in their b values: EPI0 has a b value < 50 and is less diffusion sensitive, while the b value of the EPIdiff is approximately 500. With both sequences the whole slice is measured with a single HF excitation pulse followed by a single 180° RF pulse. All other necessary echos are achieved with a series of inversed read-out gradients (TE = 69 ms, TR = 0.6 s, inversion time = 170 ms, number of averages = 1). The data matrix was 256 x 256, the voxel size 1.32 x 1.32 x 5 mm3.
|
RESULTS AND DISCUSSION
|
|---|
Spatial Description of the Pore Space
Representative slices are given in Fig. 2
. Lighter grayscales indicate a higher proton spin density. The minimal voxel size that could be obtained with the current setup was 0.5 x 0.5 x 0.5 mm3 with a signal to noise ratio (S/N) between 15 (medium sand) and 50 (medium gravel). The spatial resolution is sufficient to map the pore structure of the medium and fine gravel (Layers 1, 2, and 4) accurately. Also, larger structural elements like the different sediments in Layer 6 can be seen clearly. The medium sand in Layer 6 shows a homogeneous proton spin density with no structural information.
If the pore diameter is less than 0.5 mm a direct mapping is not possible. Up to now this is partly an artificial limitation, because for medical applications the resolution was sufficient, so there was no need for further enhancements. Based on the measured S/N a voxel size of 0.2 x 0.2 x 0.2 mm3 seems possible without a change of the setup. The S/N can be enhanced by increasing the number of averages (thus increasing the acquisition time), modifying the sequences, or using more sensitive receiving coils.
A three-dimensional rendering of the pore space in the medium gravel is given in Fig. 3
. This data may serve as input data for numerical flow models and as a benchmark for numerical models that generate the structure of porous media. Additionally, the propagation of single particles through the pore space can be simulated.
Diffusion Measurements
The measured coefficients of molecular diffusion for water are given in Table 2. Only those voxel elements with a water content of at least 10% were taken into account. The values of the diffusion coefficients are well in the range given in literature (Atkins, 1982). It can be seen that the diffusion coefficients will be underestimated if the S/N is low (see Layer 4). The standard deviation of the diffusion coefficients within one slice is less than 1%. If the voxels with less than 10% water content are included in the calculation the standard deviation is 50 to 100%. The distribution of the diffusion coefficients is homogeneous within the measured slices. The difference between the coarse sand and the mediums and in Layer 6 is not significant, which indicates that the diffusion in the medium sand is not restricted at the time scale of the sequence.
The spatial resolution of the diffusion measurements is limited to the voxel size of 1.32 x 1.32 x 5 mm3. With that resolution only larger flow structures (e.g., preferential flow paths or fissures) can be measured individually. Usually a voxel contains several pores and the measured D is actually a mean value for all pore spaces. Further work will therefore focus on an enhancement of the spatial resolution.
Flow Measurements
After measuring the diffusion coefficients, flow at five different velocities from 0.15 mm/s to approximately 6.67 mm/s was measured. The relation between signal intensity and flow rate for the EPI0 sequence is given in Fig. 4
. The slope of the regression lines (in logarithmic scaling) is proportional to the flow velocities. The regression lines for the EPIdiff sequence are very similar, although the regression parameters differ. For both sequences the regression coefficients in all layers were R2 > 0.98. The standard deviation of the regression factor a was less than 5%. The signal in the center of Layer 6 at a volume flow of 3.33 mL/s (corresponding to a flow velocity of about 6.67 mm/s) is not significantly different from the background (6 ± 3 arbitrary units [A.U.]). This flow velocity marks the upper limit of flow velocities that can be measured with the described setup and sequences.
In Fig. 5
the mean flow velocities measured with MRI are plotted against the velocities obtained from the numerical model together with the regression lines and the confidence intervals (95% level) for the regression. The measurements with MRI are in good agreement to the modeled velocities. In general, the EPI0 sequence underestimates the flow velocity. Both sequences are not capable of measuring the highest flow velocity (Layer 6, center) correctly. The more pronounced deviation of the EPI0 values may be attributable to the relatively low flow velocity for that sequence.
View larger version (29K):
[in this window]
[in a new window]
|
Fig. 5. Measured vs. modeled velocities. Confidence intervals for the regression lines are given at the 95% level.
|
|
The relative frequencies of the mean voxel velocities (mean flow velocity within one voxel) are given in Fig. 6
at a flow rate of 1.97 mL/s. In the medium gravel (Layer 1), the voxel size is small against the grain size. Therefore, many voxels are located completely within the sediment matrix where the flow velocity is zero. The histogram is significantly different from the other layers in the sense that there is a linear decrease of the relative frequency for higher flow velocities. At Layer 2 (fine gravel) some voxels can be considered completely filled with sediment. The relative frequency is almost normally distributed with a maximum at 0.50 to 0.83 mm/s. This shift to higher voxel velocities is continued in the coarse sand (Layer 3). With decreasing grain and pore size more pore spaces are averaged within one voxel. Hence, the measured voxel velocity distribution decreases and approaches the mean velocity in the pores. The different flow velocities in the coarse sand and the medium sand of Layer 6 lead to a wide velocity distribution. In the coarse sand core velocities up to 5 mm/s were measured at a flow rate of 1.97 mL/s. The velocities in the medium sand were up to 0.67 mm/s.
The problems of calculating the flow velocities at a voxel scale are comparable with those of the diffusion measurements. At the current spatial resolution the velocities in preferential flow paths and pore spaces within the medium gravel can be measured individually. In other cases regions of different flow velocities can be mapped, such as the center in Layer 6 (Fig. 7)
, some dead end pore spaces at the top and bottom of the column, or a flow heterogeneity at the column inlet.
Given the voxel velocities and the pore space, the volume flux at distinct locations can be calculated. The volume fluxes at vertical cross-sections and at different flow rates are given in Table 3. The fluxes were calculated using the velocities obtained with the EPIdiff sequence. The agreement between the mean calculated volume flow and the measured flux is satisfactory. The standard deviations are in the range of 25%. At low velocities the standard deviation increases. In Layer 6 about 54 ± 20% of the volume passes through the coarse sand core, which has a ratio of 25% of the cross-sectional area. The ratio calculated with the hydraulic conductivity values given in Table 1 is 37%.
Further development will focus on an enhancement of the spatial resolution, while keeping the almost real-time measurement conditions. Currently, the time between two measurements of the same slice is <10 s. An increase of the acquisition time, which would have been the easy way to increase S/N, does not seem appropriate because of possible changes of the flow during acquisition.
Overall, MRI may contribute to improved soil column studies by providing reliable data on the pore structure and the flow processes without the need for a tracer and without destroying soil structure. Column studies focusing on preferential flow or colloidal transport may be improved significantly when performed under MRI control. The data provided may be fed into numerical models that attempt to link physical and chemical processes on the pore scale to contaminant transport. We must admit that the sediments used in this study had rather large pores. But even if the pores are smaller than the voxel size, the benefit of measuring a mean velocity within the pores at a high spatial resolution, compared with an average linear velocity at a much larger scale, remains. Numerical methods with a more empirical background will benefit from benchmark data for validation.
The technique presented may also contribute to studies on the details of unsaturated flow and multiphase flow in subsurface environments. Other applications include the simulation of remediation techniques, ground water flow to wells, or the development of fissures and diffusive transport processes in mineral liners at landfill sites.
|
REFERENCES
|
|---|
- Akin, S., M.R.B. Demiral, and E. Okandan. 1996. A novel method of porosity measurement utilizing computerized tomography. In Situ 20:347–365.
- Atkins, P.W. 1982. Physical chemistry. Oxford Univ. Press, Oxford.
- Becker, E.D. 1980. High resolution NMR. Academic Press, New York.
- Bovey, F.A. 1988. Nuclear magnetic resonance spectroscopy. Academic Press, New York.
- Callaghan, P.T. 1991. Principles of nuclear magnetic resonance microscopy. Claredon Press, Oxford.
- Chudek, J.A., and A.D. Reeves. 1998. An application of nuclear magnetic resonance imaging to study migration rates of oil-related residues in estuarine sediments. Biodegradation 9:443–449.[Web of Science]
- Corapcioglu, M.Y., S. Chowdhury, and S.E. Roosevelt. 1997. Micromodel visualization and quantification of solute transport in porous media. Water Resour. Res. 33:2547–2558.
- DIN 18123. 1983. Bestimmung der Korngrößssenverteilung. Beuth, Berlin/Köln.
- DIN 18124. 1989. Bestimmung der Korndichte-Kapillarpyknometer-Weithalspyknometer. Beuth, Berlin/Köln.
- DIN 18130 T 1. 1989. Bestimmung des Wasserdurchlässigkeitsbeiwerts. Beuth, Berlin/Köln.
- Harvey, R.W., and S.P. Garabedian. 1991. Use of colloid filtration theory in modeling movement of bacteria through a contaminated sandy aquifer. Environ. Sci. Technol. 25:178–185.
- Jasti, J.K., G. Jesion, and L. Feldkamp. 1993. Microscopic imaging of porous media with x-ray computer tomography. SPE Form. Eval. 8:189–193.
- McCarthy, J.F., and J.M. Zachara. 1989. Subsurface transport of contaminants. Environ. Sci. Technol. 23:496–502.
- Morneburg, H. 1995. Bildgebende Systeme für die medizinische Diagnostik. Wiley-VCH, Weinheim, Germany.
- Neumann-Haefelin, T., H.J. Wittsack, F. Wenserski, M. Siebler, R.J. Seitz, U. Modder, and H.J. Freund. 1999. Diffusion- and perfusion-weighted MRI—The DWI/PWI mismatch region in acute stroke. Stroke 30:1591–1597.[Abstract/Free Full Text]
- Nitz, W.R. 1999. MR imaging acronyms and clinical applications. Eur. Radiol. 9:979–997.[Medline]
- Onaisi, A., A. Audibert, M.T. Bieber, L. Bailey, J. Denis, and P.S. Hammond. 1993. X-ray tomography vizualization and mechanical modeling of swelling shale around the wellbore. J. Petrol. Sci. Eng. 9:313–329.
- Oswald, S., W. Kinzelbach, A. Greiner, and G. Brix. 1997. Observation of flow and transport processes in artificial porous media via magnetic resonance imaging. Geoderma 80:417–429.
- Rehmann, L.L., C. Welty, and R.W. Harvey. 1999. Stochastic analysis of virus transport in aquifers. Water Resour. Res. 35:1987–2006.
- Reinken, G., F. Fuehr, N. Zadgorsky, H. Halling, and O. Schult. 1995. Application of a novel high-resolution tomographic method for non-destructive characterization of pore structure in soil cores. BCPC Monogr. 62:39–44.
- Rokitta, M., U. Zimmermann, and A. Haase. 1999. Fast NMR flow measurements in plants using flash imaging. J. Magn. Reson. 137:29–32.[Medline]
- Skoog, D.A., and J.J. Leary. 1992. Principles of instrumental analysis. 4th ed. Saunders College Publ., Orlando, FL.
- Song, L., P.R. Johnson, and M. Elimelech. 1994. Kinetics of colloid deposition onto heterogenously charged surfaces in porous media. Environ. Sci. Technol. 28:1164–1171.
- Stejskal, E.O., and J.E. Tanner. 1965. Spin diffusion measurements: Spin echoes in the presence of a time-dependent field gradient. J. Chem. Phys. 42:288.
- Veerapaneni, S., and M.R. Wiesner. 1997. Deposit morphology and head loss development in porous media. Environ. Sci. Technol. 31:2738–2744.
- Vogel, H.J. 1997. Morphological determination of pore connectivity as a function of pore size using serial sections. Eur. J. Soil Sci. 48:365–377.
- Vogel, H.J., and K. Roth. 1998. A new approach for determining effective soil hydraulic functions. Eur. J. Soil Sci. 49:547–556.
- Wu, S.C., and P.M. Gschwend. 1986. Sorption kinetics of hydrophobic organic compounds to natural sediments and soils. Environ. Sci. Technol. 20:717–725.[Web of Science]
TOP
ABSTRACT
INTRODUCTION
