A multiple-point statistics algorithm for 3D pore space reconstruction from 2D images

Fluid flow behavior in a porous medium is a function of the geometry and topology of its pore space. The construction of a three dimensional pore space model of a porous medium is therefore an important first step in characterizing the medium and predicting its flow properties. A stochastic technique for reconstruction of the 3D pore structure of unstructured random porous media from a 2D thin section training image is presented. The proposed technique relies on successive 2D multiple point statistics simulations coupled to a multi-scale conditioning data extraction procedure. The Single Normal Equation Simulation Algorithm (SNESIM), originally developed as a tool for reproduction of long-range, curvilinear features of geological structures, serves as the simulation engine. Various validating criteria such as marginal distributions of pore and grain, directional variograms, multiple-point connectivity curves, single phase effective permeability and two phase relative permeability calculations are used to analyze the results. The method is tested on a sample of Berea sandstone for which a 3D micro-CT scanning image is available. The results confirm that the equi-probable 3D realizations obtained preserve the typical patterns of the pore space that exist in thin sections, reproduce the long-range connectivities, capture the characteristics of anisotropy in both horizontal and vertical directions and have single and two phase flow characteristics consistent with those of the measured 3D micro-CT image.     

On the relationship between indicators of geostatistical,flow and transport connectivity

Connectivity of high conductivity (K) paths is important because it can lead to channeling, i.e. flow along preferential paths, which can reduce travel times very significantly. Nevertheless, limited effort has been devoted to defining the concept quantitatively. We propose and evaluate nine indicators of connectivity. Three account for the presence of flow connectivity, that is, the flow rate increase caused by preferential flow paths. Two account for the presence of transport connectivity, that is, the existence of fast paths allowing early solute arrival. The remaining four are statistical indicators based on two- and multiple-point statistics. We test these indicators on heterogeneous conductivity fields with different visual connectivity. The indicators of flow connectivity and one of the transport connectivity indicators succeed in identifying the increased presence of connected high-K features. The two-point statistical indicators fail to do so. The directional multi-point statistical indicator performs better. None of the statistical indicators correlate with the flow and transport indicators. We find only weak dependence between the flow and transport indicators. Our results suggest that transport connectivity is much less sensitive to barriers, which control flow connectivity. Instead, transport connectivity appears to be controlled by the existence of narrow, possibly discontinuous high-K paths. Therefore, we conclude that connectivity is a process-dependent concept.     

Multiscale flow and transport model in three-dimensional fractal porous media

This paper proposes a multiscale flow and transport model which can be used in three-dimensional fractal random fields. The fractal random field effectively describes a field with a high degree of variability to satisfy the one-point statistics of Levy-stable distribution and the two-point statistics of fractional Levy motion (fLm). To overcome the difficulty of using infinite variance of Levy-stable distribution and to provide the physical meaning of a finite domain in real space, truncated power variograms are utilized for the fLm fields. The fLm model is general in the sense that both stationary and commonly used fractional Brownian motion (fBm) models are its special cases. When the upper cutoff of the truncated power variogram is close to the lower cutoff, the stationary model is well approximated. The commonly used fBm model is recovered when the Levy index of fLm is 2. Flow and solute transport were analyzed using the first-order perturbation method. Mean velocity, velocity covariance, and effective hydraulic conductivity in a three-dimensional fractal random field were derived. Analytical results for particle displacement covariance and macrodispersion coefficients are also presented. The results show that the plume in an fLm field moves slower at early time and has more significant long-tailing behavior at late time than in fBm or stationary exponential fields. The proposed fractal transport model has broader applications than those of stationary and fBm models. Flow and solute transport can be simulated for various scenarios by adjusting the Levy index and cutoffs of fLm to yield more accurate modeling results.     

Impact of geometrical properties on permeability and fluid phase distribution in porous media

To predict fluid phase distribution in porous media, the effect of geometric properties on flow processes must be understood. In this study, we analyze the effect of volume, surface, curvature and connectivity (the four Minkowski functionals) on the hydraulic conductivity and the water retention curve. For that purpose, we generated 12 artificial structures with 800³ voxels (the units of a 3D image) and compared them with a scanned sand sample of the same size. The structures were generated with a Boolean model based on a random distribution of overlapping ellipsoids whose size and shape were chosen to fulfill the criteria of the measured functionals. The pore structure of sand material was mapped with X-rays from synchrotrons.     

Reconstruction of porous media using multiple-point statistics with data conditioning

The reconstruction of porous media is of great importance in predicting fluid transport properties, which are widely used in various fields such as catalysis, oil recovery, medicine and aging of building materials. The real three-dimensional structural data of porous media are helpful to describe the irregular topologic structures of porous media. By using multiple-point statistics (MPS) to extract the characteristics of real porous media acquired from micro computed tomography (micro-CT) scanning, the probabilities of structural characteristics of pore spaces are obtained first, and then reproduced in the reconstructed regions. One solution to overcome the anisotropy of training images is to use real 3D volume data as a training image (TI). The CPU cost and memory burden brought up by 3D simulations can be reduced greatly by selecting the optimal multiple-grid template size that is determined by the entropy of a TI. Moreover, both soft data and hard data are integrated in MPS simulation to improve the accuracy of reconstructed images. The variograms and permeabilities, computed by lattice Boltzmann method, of the reconstructed images and the target image obtained from real volume data are compared, showing that the structural characteristics of reconstructed porous media using our method are similar to those of real volume data.     

A lattice Boltzmann study on the impact of the geometrical properties of porous media on the steady state relative permeabilities on two-phase immiscible flows

In the current paper, the effect of the geometrical characteristics of 2-D porous media on the relative permeability in immiscible two-phase flows is studied. The generation of the different artificial porous media is performed using a Boolean model based on a random distribution of overlapping circles/ellipses, the size and shape of which are chosen to satisfy the specific Minkowski functionals (i.e. volume fraction, solid line contour length, connectivity). The study aims to identify how each different Minkowski functional affects the relative permeability of each phase at various saturations of the non-wetting phase. A 2-D multi-relaxation time (MRT) lattice Boltzmann model (LBM) that can handle high density ratios is employed in the simulation. The relationship between the driving forces G and the relative permeabilities of the two phases for every artificial structure is quantified. It is found that for high non-wetting phase saturations (fully connected flow), a non-linear relationship exists between the non-wetting phase flow rate and the driving force, whilst this relationship becomes linear at higher magnitudes of the latter. The force magnitude required to approach the linear region is highly influenced by the pore size distribution and the connectivity of the solid phase. For lower non-wetting phase saturation values, its relative permeability in the linear regime decreases as the fraction of small pores in the structure increases and the non-wetting phase flow becomes disconnected. A strong influence of the solid phase connectivity is also observed.     

Generation of three dimensional flow fields for statistically anisotropic heterogeneous porous media

A methodology for generating three dimensional (3D) flow fields for statistically anisotropic heterogeneous porous media is presented and demonstrated. The simulated flow fields are shown to exhibit the input spatial correlation structure and observe mass continuity. Sample flow fields are presented in the form of cross sectional slices of the 3D formation. These cross sections demonstrate visually the characteristics of subsurface flow. The method was found to be faster than traditional techniques in terms of its computational requirements. Given this method, it is possible to generate the large number of realizations of a velocity field necessary to compute high order statistics in transport problems.     

Properties and limitations of sequential indicator simulation

The sequential indicator algorithm is a widespread geostatistical simulation technique that relies on indicator (co)kriging and is applicable to a wide range of datasets. However, such algorithm comes up against several limitations that are often misunderstood. This work aims at highlighting these limitations, by examining what are the conditions for the realizations to reproduce the input parameters (indicator means and correlograms) and what happens with the other parameters (other two-point or multiple-point statistics). Several types of random functions are contemplated, namely: the mosaic model, random sets, models defined by multiple indicators and isofactorial models. In each case, the conditions for the sequential algorithm to honor the model parameters are sought after. Concurrently, the properties of the multivariate distributions are identified and some conceptual impediments are emphasized. In particular, the prior multiple-point statistics are shown to depend on external factors such as the total number of simulated nodes and the number and locations of the samples. As a consequence, common applications such as a flow simulation or a change of support on the realizations may lead to hazardous interpretations.     

Topological groundwater hydrodynamics

Topological groundwater hydrodynamics is an emerging subdiscipline in the mechanics of fluids in porous media whose objective is to investigate the invariant geometric properties of subsurface flow and transport processes. In this paper, the topological characteristics of groundwater flows governed by the Darcy law are studied. It is demonstrated that: (i) the topological constraint of zero helicity density during flow is equivalent to the Darcy law; (ii) both steady and unsteady groundwater flows through aquifers whose hydraulic conductivity is an arbitrary scalar function of position and time are confined to surfaces on which the streamlines of the flow are geodesic curves; (iii) the surfaces to which the flows are confined either are flat or are tori; and (iv) chaotic streamlines are not possible for these flows, implying that they are inherently poorly mixing in advective solute transport.     

Pattern-based conditional simulation with a raster path: a few techniques to make it more efficient

Mass transport is known to depend on heterogeneity in geological formations. This entails geological bodies with complex geometries. The major interest of multiple-point simulation is its ability to reproduce such geological features through the use of a training image. The idea behind the training image is to describe a geological concept with the expected geological architecture. Its structural content is then used to infer multiple-point statistics. This yields a database with a variety of possible patterns or events. In this paper, we present a hybrid algorithm combining geostatistical multiplepoint and texture synthesis techniques for simulating geological reservoir models constrained to hard data. The proposed algorithm is a two steps process, involving first analysis with the building of an organized database from the training image content, and second synthesis with the simulation of a realization. Various tests are performed to investigate the potential of the algorithm in terms of computation time and ability to properly reproduce the shapes and connectivity features of the objects represented in the training image. We also propose a few improvements to make the algorithm more efficient. Last, six examples are presented based upon different kinds of training images depicting large-scale channelized and fractured media as well as fine-scale porous media.     

The effects of topographic depressions on multiscale overland flow connectivity: A high‐resolution spatiotemporal pattern analysis approach based on connectivity statistics

In watershed modelling, the traditional practice of arbitrarily filling topographic depressions in digital elevation models has raised concerns. Advanced high‐resolution remote sensing techniques, including airborne scanning laser altimetry, can identify naturally occurring depressions that impact overland flow. In this study, we used an ensemble physical and statistical modelling approach, including a 2D hydraulic model and two‐point connectivity statistics, to quantify the effects of depressions on high‐resolution overland flow patterns across spatial scales and their temporal variations in single storm events. Computations for both models were implemented using graphic processing unit‐accelerated computing. The changes in connectivity statistics for overland flow patterns between airborne scanning laser altimetry‐derived digital elevation models with (original) and without (filled) depressions were used to represent the shifts of overland flow response to depressions. The results show that depressions can either decrease or increase (to a lesser degree and shorter duration) the probability that any two points (grid locations) are hydraulically connected by overland flow pathways. We used macro‐connectivity states (Φ) as a watershed‐specific indicator to describe the spatiotemporal thresholds of connectivity variability caused by depressions. Four states of Φ are identified in a studied watershed, and each state represents different magnitudes of connectivity and connectivity changes (caused by depressions). The magnitude of connectivity variability corresponds to the states of Φ, which depend on the topological relationship between depressions, the rising/recession limb, and the total rainfall amount in a storm event. In addition, spatial distributions of connectivity variability correlate with the density of depression locations and their physical structures, which cause changes in streamflow discharge magnitude. Therefore, this study suggests that depressions are “nontrivial” in watershed modelling, and their impacts on overland flow should not be neglected. Connectivity statistics at different spatial scales and time points within a watershed provide new insights for characterizing the distributed and accumulated effects of depressions on overland flow.     

Effective dispersion in a chemically heterogeneous medium under temporally fluctuating flow conditions

We investigate effective solute transport in a chemically heterogeneous medium subject to temporal fluctuations of the flow conditions. Focusing on spatial variations in the equilibrium adsorption properties, the corresponding fluctuating retardation factor is modeled as a stationary random space function. The temporal variability of the flow is represented by a stationary temporal random process. Solute spreading is quantified by effective dispersion coefficients, which are derived from the ensemble average of the second centered moments of the normalized solute distribution in a single disorder realization. Using first-order expansions in the variances of the respective random fields, we derive explicit compact expressions for the time behavior of the disorder induced contributions to the effective dispersion coefficients. Focusing on the contributions due to chemical heterogeneity and temporal fluctuations, we find enhanced transverse spreading characterized by a transverse effective dispersion coefficient that, in contrast to transport in steady flow fields, evolves to a disorder-induced macroscopic value (i.e., independent of local dispersion). At the same time, the asymptotic longitudinal dispersion coefficient can decrease. Under certain conditions the contribution to the longitudinal effective dispersion coefficient shows superdiffusive behavior, similar to that observed for transport in s stratified porous medium, before it decreases to its asymptotic value. The presented compact and easy to use expressions for the longitudinal and transverse effective dispersion coefficients can be used for the quantification of effective spreading and mixing in the context of the groundwater remediation based on hydraulic manipulation and for the effective modeling of reactive transport in heterogeneous media in general.     

Generalized Radial Flow in Synthetic Flow Systems

Traditional analysis methods used to determine hydraulic properties from pumping tests work well in many porous media aquifers, but they often do not work in heterogeneous and fractured‐rock aquifers, producing non‐plausible and erroneous results. The generalized radial flow model developed by Barker (1988) can reveal information about heterogeneity characteristics and aquifer geometry from pumping test data by way of a flow dimension parameter. The physical meaning of non‐integer flow dimensions has long been a subject of debate and research. We focus on understanding and interpreting non‐radial flow through high permeability conduits within fractured aquifers. We develop and simulate flow within idealized non‐radial flow conduits and expand on this concept by simulating pumping in non‐fractal random fields with specific properties that mimic persistent sub‐radial flow responses. Our results demonstrate that non‐integer flow dimensions can arise from non‐fractal geometries within aquifers. We expand on these geometric concepts and successfully simulate pumping in random fields that mimic well‐test responses seen in the Culebra Dolomite above the Waste Isolation Pilot Plant.     

Stochastic analysis of immiscible displacement of the fluids with arbitrary viscosities and its dependence on support scale of hydrological data

Stochastic analysis is commonly used to address uncertainty in the modeling of flow and transport in porous media. In the stochastic approach, the properties of porous media are treated as random functions with statistics obtained from field measurements. Several studies indicate that hydrological properties depend on the scale of measurements or support scales, but most stochastic analysis does not address the effects of support scale on stochastic predictions of subsurface processes. In this work we propose a new approach to study the scale dependence of stochastic predictions. We present a stochastic analysis of immiscible fluid–fluid displacement in randomly heterogeneous porous media. While existing solutions are applicable only to systems in which the viscosity of one phase is negligible compare with the viscosity of the other (water–air systems for example), our solutions can be applied to the immiscible displacement of fluids having arbitrarily viscosities such as NAPL–water and water–oil. Treating intrinsic permeability as a random field with statistics dependant on the permeability support scale (scale of measurements) we obtained, for one-dimensional systems, analytical solutions for the first moments characterizing unbiased predictions (estimates) of system variables, such as the pressure and fluid–fluid interface position, and we also obtained second moments, which characterize the uncertainties associated with such predictions. Next we obtained empirically scale dependent exponential correlation function of the intrinsic permeability that allowed us to study solutions of stochastic equations as a function of the support scale. We found that the first and second moments converge to asymptotic values as the support scale decreases. In our examples, the statistical moments reached asymptotic values for support scale that were approximately 1/10000 of the flow domain size. We show that analytical moment solutions compare well with the results of Monte Carlo simulations for moderately heterogeneous porous media, and that they can be used to study the effects of heterogeneity on the dynamics and stability of immiscible flow.     

Impact of fractional probability distributions on statistics of hydraulic conductivity,dynamics of groundwater flow and solute transport at a low-permeability site

The fractional Brownian motion (fBm) and fractional Lévy motion (fLm) can easily describe the geometry and the statistical structure of hydraulic conductivity (K) for real-world. However, the fBm and fLm models have not been systematically evaluated when building the K field for a low-permeability site. In this study, both the fBm and fLm are used to simulate the low-K field at NingCheGu (NCG), Tianjin, China. Groundwater flow and solute transport are then computed using MODFLOW and MT3DMS, respectively, and the influence of the fBm/fLm models for K on groundwater flow and solute transport is discussed. Results show that the fLm fits better the statistics of the low-K medium than fBm, and the random logarithmic K (LnK) field generated by fLm is more stable because the resultant LnK field captures more of the measured properties at the field site than that generated by fBm. In contrast, the LnK generated by fBm is more likely to form both high-K channels and low-K barriers. The fBm therefore predicts more extreme behaviours in flow and transport, including the preferential flow, low-concentration blocks and solute retention. The overall groundwater renewal period and solute travel time for the fLm simulation are slightly shorter than those for fBm. The impacts of the fLm and fBm models on the statistics of the resultant LnK fields and the dynamics of groundwater flow and solute transport revealed by this study shed light on the selection and evaluation of the fractional probability distribution models in capturing the K fields for low-K media.     

Lessons Learned from 25 Years of Research at the MADE Site

Field studies at well‐instrumented research sites have provided extensive data sets and important insights essential for development and testing of transport theories and mathematical models. This paper provides an overview of over 25 years of research and lessons learned at one of such field research sites on the Columbus Air Force Base in Mississippi, commonly known as the Macrodispersion Experiment (MADE) site. Since the mid‐1980s, field data from the MADE site have been used extensively by researchers around the world to explore complex contaminant transport phenomena in highly heterogeneous porous media. Results from field investigations and modeling analyses suggested that connected networks of small‐scale preferential flow paths and relative flow barriers exert dominant control on solute transport processes. The classical advection‐dispersion model was shown to inadequately represent plume‐scale transport, while the dual‐domain mass transfer model was found to reproduce the primary observed plume characteristics. The MADE site has served as a valuable natural observatory for contaminant transport studies where new observations have led to better understanding and improved models have sprung out analysis of new data.     

各向异性多孔介质数字岩心模型重构通用算法

多孔介质材料内部的微观孔隙结构直接决定着其宏观的各类物理特性, 岩石等天然多孔介质材料内部分布有复杂的孔隙结构, 且呈现出较明显的各向异性特征, 准确的描述和捕捉该特征对于岩石高精度数字岩心模型的构建以及后续的孔隙尺度模拟至关重要.基于先前所提出的重构加速框架, 选用单点概率函数、两点概率函数和线性路径函数等统计函数来对孔隙结构特征进行定量化描述, 通过对孔隙结构在不同方向上的分布特征进行独立的提取表征, 提出了一种能够适用于岩石等各类天然各向异性多孔介质数字岩心模型构建的改进模拟退火算法.以具有显著各向异性特征的人工合成二维切片图像及真实碳酸盐岩和砂岩三维CT扫描图像作为测试参考样本, 分别采用传统方法和本研究所提出方法构建了其相应的数字岩心模型, 并对两者的连通性及不同方向上孔隙结构的分布情况和渗透特性进行了全面的综合对比分析, 结果表明: 相比于传统方法, 本研究所提出方法在确保重构岩心模型具有良好连通性的前提下, 能够对不同方向上孔隙结构分布的差异性进行更加准确的捕获和表征, 且重构岩心模型各个方向上的渗透性也与参考真实岩心表现出了较高程度的一致性, 实现了高精度还原重建参考图像中孔隙结构各向异性特征的目标, 从而验证了本研究所提出的各向异性多孔介质数字岩心模型重构方法的准确性和有效性.