Geofluids

Volume 2019, Article ID 6810467, 13 pages

https://doi.org/10.1155/2019/6810467

## Impact of Synthetic Porous Medium Geometric Properties on Solute Transport Using Direct 3D Pore-Scale Simulations

^{1}National Research Council of Italy, Water Research Institute, Area della Ricerca di Roma 1-Montelibretti, Strada Provinciale 35d, km 0.7 Montelibretti Roma, Italy^{2}Institute of Geochemistry and Petrology, ETH Zurich, Clausiusstrasse 25, CH-8092 Zurich, Switzerland^{3}Department of Earth, Environmental and Planetary Sciences, Brown University, Providence, 02912 RI, USA^{4}Department of Computational Hydrosystems Helmholtz, Centre for Environmental Research, UFZ Permoserstr. 15 04318 Leipzig, Germany

Correspondence should be addressed to Nicolas Guyennon; ti.rnc.asri@nonneyug

Received 17 July 2018; Revised 30 November 2018; Accepted 22 January 2019; Published 3 April 2019

Academic Editor: Andri Stefansson

Copyright © 2019 Paolo Roberto Di Palma et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

Transport processes in porous media have been traditionally studied through the parameterization of macroscale properties, by means of volume-averaging or upscaling methods over a representative elementary volume. The possibility of upscaling results from pore-scale simulations, to obtain volume-averaging properties useful for practical purpose, can enhance the understanding of transport effects that manifest at larger scales. Several studies have been carried out to investigate the impact of the geometric properties of porous media on transport processes for solute species. However, the range of pore-scale geometric properties, which can be investigated, is usually limited to the number of samples acquired from microcomputed tomography images of real porous media. The present study takes advantage of synthetic porous medium generation to propose a systematic analysis of the relationships between geometric features of the porous media and transport processes through direct simulations of fluid flow and advection-diffusion of a non-reactive solute. Numerical simulations are performed with the lattice Boltzmann method on synthetic media generated with a geostatistically based approach. Our findings suggest that the advective transport is primarily affected by the specific surface area and the mean curvature of the porous medium, while the effective diffusion coefficient scales as the inverse of the tortuosity squared. Finally, the possibility of estimating the hydrodynamic dispersion coefficient knowing only the geometric properties of porous media and the applied pressure gradient has been tested, within the range of tested porous media, against advection-diffusion simulations at low Reynolds (<10^{-1}) and Peclet numbers ranging from 10^{1} to 10^{-2}.

#### 1. Introduction

A quantitative relation between the structure of porous media and its effect on solute transport is fundamental to our understanding of groundwater contamination dynamics and the development of pollutant remediation strategies [1]. Transport processes in porous media have been traditionally studied through the parametrization of macroscale properties (e.g., hydraulic conductivity or hydrodynamic dispersion coefficients) by means of volume-averaging or upscaling methods over a representative elementary volume.

However, these upscaling approaches are not valid when considering length scales on the order of a single pore or the complexity of processes (i.e., physical, chemical, and biological) that take place at the pore scale [2]. Developing averaging approaches from the pore-scale approach can provide further information on the macroscopic medium properties and then enhance the understanding of transport effects which exhibit at a greater scale (e.g., field or laboratory scale) [3–6].

The development of digital image processing led to a large expansion of the quantitative analyses of rock and soil structure. Various tomographic techniques are nowadays available for a three-dimensional non-invasive visualization of structural properties with a spatial resolution down to microns [7–11]. Even with the availability of these techniques, the cost of imaging limits the number of samples that can be analyzed [12]. As an alternative, a statistical approach to generate media with given imposed structural features can be a cost-effective means for producing large numbers of realizations that represent those features of natural samples deemed relevant. The two approaches that are commonly used to generate synthetic porous media can be divided into statistical models (e.g., multipoint statistics and autocorrelation functions) and object-based methods, where geometrical objects are randomly placed in a domain to simulate diagenesis processes [8].

The well-established use of pore-scale approach consists of estimating averaged properties, such as permeability through direct single-phase fluid flow simulations [6, 13–18] and relative permeability for multiphase fluid flows [4, 19–22]. Recently, some authors focused on the relationship between velocity distribution and pore structure [23–28] as well as transport processes [2, 24, 29–36].

The present research is aimed at analyzing the relationships between geometric features of different porous media and transport processes through direct simulations of fluid flow and advection-diffusion of non-reactive solute using the lattice Boltzmann method (LBM). The LBM is selected for its ability to easily model different phenomena and its suitability to be parallelized and to exploit the computing power of HPC systems [37]. The question we address here is which geometric properties of the porous media control fluid flow and solute transport? In order to perform a sensitivity analysis on a wide class of porous microstructures, multiple synthetic porous media are generated using a geostatistical approach. We also consider two available sandpack samples (https://www.imperial.ac.uk/earth-science/research/research-groups/perm/research/pore-scale-modelling/micro-ct-images-and-networks/, last access 16/01/2019) to compare results for synthetic and natural porous media. The choice of sandpack images is motivated because alluvial aquifers are the most commonly contaminated systems [1]. Each porous medium is described by the Minkowski functionals: porosity, specific surface area, mean curvature, and Euler characteristic [38].

#### 2. Materials and Methods

To identify the dominant pore-scale parameters that control the transport of non-reactive solutes, LBM simulations (Section 2.2) are run over synthetic porous media and micro-CT images of sandpacks (Section 2.1). A specific simulation setup is designed to first consider separately the impact of advective and diffusive transport (Section 2.3). Our objective is to reduce the number of variables to consider for the advection-diffusion simulations, which are subsequently carried out on a selection of porous media. In fact, a full sensitivity analysis between the correlation lengths () set for the generation model, the Minkowski functionals for porous media characterizations (porosity (), specific surface area (SSA), mean curvature (MC), and Euler number ()), and tortuosity () over different transport regimes would require a computational effort that goes beyond our current capacity. Finally, the metrics used to quantify the simulation results over the selected parameters are described in Section 2.4.

Acronyms and abbreviations adopted throughout the paper are summarized in Supplementary Material S3.

##### 2.1. 3D Porous Medium

The procedure to generate the porous media is described in details in Di Palma et al. [23]. Here, we only summarize the main steps as follows: (i) the 3D domain and resolution of the porous medium are set; (ii) the geostatistical model is assigned, setting the type of semivariogram (exponential and Gaussian) and the correlation length . Pore-scale variability (i.e., for multiple length scale) is neglected; (iii) a 3D Gaussian random field (zero mean and variance of 1) following the assigned spatial autocorrelation function (semivariogram) is generated with the *TFracGen* code [39]; (iv) a threshold is imposed to the normal cumulated distribution function to obtain an indicator field of pores and grains with an assigned target porosity *φ*; and (v) the spatial autocorrelation of the generated indicator field is estimated by fitting the Matérn covariance model [40, 41] as a function of the smoothness parameter and the associated practical range. This step is necessary as the transformation of a Gaussian field into an indicator field does not preserve the characteristics of the spatial autocorrelation [12, 42]. The resulting smoothness parameters of the Matèrn model for the porous media generated by exponential and Gaussian variograms are and , respectively. These values are consistent with the findings of Di Palma et al. [23]. However, for the sake of simplicity, we refer for the rest of the paper to EXX and GXX to porous media generated by exponential and Gaussian functions, respectively, where XX refers to the correlation length of each medium.

In this work, each porous medium is generated at a spatial resolution of 10 *μ*m for a domain size of 3 mm^{3}, resulting in a number of voxels equal to 300^{3}. A greater calculation domain would have involved a huge computationally cost, especially for advection-diffusion simulations. Porous media are generated in a range between 3% and 11% of the ratio correlation length over domain size to guarantee a suitable balance between the spatial resolution impact on short correlation lengths and the representativeness of porous medium generation for great correlation lengths compared to domain size [23]. It is worth stressing that the synthetic porous media are isotropic (i.e., show only one independent length scale) and statistically homogeneous. We are aware that this modelling assumption constitutes a limitation of our approach, e.g., for sedimentary porous media that often show pronounced anisotropy between the horizontal and vertical directions. However, we think that this is justified by the reduction of the numbers of descriptors required to describe each medium. Such limitation does not significantly affect the interpretability of our results.

In order to compare the analysis of the synthetic porous media with real samples, we performed additionally flow and solute transport simulations on two different sandpacks (i.e., F42C and LV60A) with the same spatial resolution of 10 *μ*m and 300^{3} voxels, whose micro-CT images are available from the pore-scale modelling Consortium of Imperial College (https://figshare.com/articles/F42C_sandstone/1189261, https://figshare.com/articles/LV60A_sandpack/1153795, last access 15/01/2019). Both the X-ray scanned sandpack media have been verified to be isotropic and homogenous [17].

To describe the spatial characteristics of the generated porous media, we adopted the four Minkowski functionals and the tortuosity. The Minkowski functionals provide a detailed description of the geometric characteristics of a porous medium and are widely used in mathematical morphology [12, 38]. We briefly summarize the general definitions for each of them but refer interested readers to Vogel et al. [38] for a complete treatment of the Minkowski functional.

The first functional is the porosity: defined as the volume of the voids over the total volume of the medium ().

The second functional is the specific surface area, as the interstitial surface area of the voids and solids: where is the surface area of the grains.

The third functional is the mean curvature (MC) representing a measure of the average pore size of the sample over the total volume:

Here, is related to the mean width and is defined as the mean value of the distance between a pair of parallel support planes [12, 43].

Finally, the fourth functional is the total curvature, also known as the Euler number that measures the connectivity of the pore space. This is more easily defined by Vogel [38, 44, 45] as a topological parameter: where are the isolated objects, are interconnected regions or loops, and are the completely enclosed cavities. Broadly speaking, a negative corresponds to high pore space connectivity [45].

We calculate the tortuosity following the definition of hydraulic tortuosity in Clennell [46] as , being the average length of the fluid paths and is the straight-line length in the direction of flow [46]. is obtained from steady-state flow simulations by counting the number of streamlines passing through the porous medium. Streamlines are calculated using MATLAB® function “stream3,” which is based on an explicit Euler integration of the ordinary differential equation for a given set of starting points , where , , and represent the components of the velocity.

In our 3D media, the porosity and correlation length are imposed in the generation process, whereas the Minkowski functionals and tortuosity are calculated after the porous medium realization and fluid flow simulations, respectively.

The imposed geometric properties of the generated porous media and associated nomenclature (ID) are reported in Table 1, together with the corresponding estimated geometric properties for the scanned samples.