Analysis of Fracture Roughness Control on Permeability Using SfM and Fluid Flow Simulations: Implications for Carbonate Reservoir Characterization

Zambrano, Miller; Pitts, Alan D.; Salama, Ali; Volatili, Tiziano; Giorgioni, Maurizio; Tondi, Emanuele

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

Geofluids

On this page

Abstract Introduction Results Discussion Conclusions Data Availability Conflicts of Interest Acknowledgments References Copyright Related Articles

Special Issue

Geofluids and Energy for the XXI Century

View this Special Issue

Research Article | Open Access

Volume 2019 | Article ID 4132386 | https://doi.org/10.1155/2019/4132386

Analysis of Fracture Roughness Control on Permeability Using SfM and Fluid Flow Simulations: Implications for Carbonate Reservoir Characterization

Miller Zambrano,^1,2,3Alan D. Pitts,^1,3Ali Salama,¹Tiziano Volatili,^1,2,3Maurizio Giorgioni,⁴and Emanuele Tondi^1,2,3

Academic Editor: Jaewon Jang

Received30 Nov 2018

Accepted14 Feb 2019

Published23 Apr 2019

Abstract

Fluid flow through a single fracture is traditionally described by the cubic law, which is derived from the Navier-Stokes equation for the flow of an incompressible fluid between two smooth-parallel plates. Thus, the permeability of a single fracture depends only on the so-called hydraulic aperture which differs from the mechanical aperture (separation between the two fracture wall surfaces). This difference is mainly related to the roughness of the fracture walls, which has been evaluated in previous works by including a friction factor in the permeability equation or directly deriving the hydraulic aperture. However, these methodologies may lack adequate precision to provide valid results. This work presents a complete protocol for fracture surface mapping, roughness evaluation, fracture modeling, fluid flow simulation, and permeability estimation of individual fracture (open or sheared joint/pressure solution seam). The methodology includes laboratory-based high-resolution structure from motion (SfM) photogrammetry of fracture surfaces, power spectral density (PSD) surface evaluation, synthetic fracture modeling, and fluid flow simulation using the Lattice-Boltzmann method. This work evaluates the respective controls on permeability exerted by the fracture displacement (perpendicular and parallel to the fracture walls), surface roughness, and surface pair mismatch. The results may contribute to defining a more accurate equation of hydraulic aperture and permeability of single fractures, which represents a pillar for the modeling and upscaling of the hydraulic properties of a geofluid reservoir.

1. Introduction

Fractures exert a vital contribution on determining the migration and storage for geofluids, such as groundwater, and hydrocarbons. Thus, the analysis and modeling of fractures are imperative for characterizing reservoirs and simulating their behavior during the production stage. Fluid flow through fractures is traditionally described by the cubic law, derived from the Navier-Stokes equation for the flow of an incompressible fluid between two smooth-parallel plates [1]. Thus, the permeability (intrinsic permeability) of a single fracture may be represented by the equation where corresponds to the hydraulic aperture. Since fractures are normally rough, the hydraulic aperture differs from the mechanical aperture (separation between the two fracture wall surfaces). The hydraulic aperture is in fact an equivalent value which can be derived from field analysis like tracer tests and hydraulic tests [2] and laboratory fluid flow experiments (e.g., [3]).

Several authors have studied the effect of roughness of the walls on fracture permeability working with various materials, such as glass [4], rocks [5], and concrete [6]. These authors included a correction term, a so-called friction factor (), on the permeability of rough fractures: where the friction factor is defined by the common base equation where is the difference between the highest peak and the lowest valley of the physical roughness, and both terms and are constants with slightly different values depending on the author. Thus, the friction factor depends only on the relative hydraulic roughness ignoring the frequency (or wavelength) of the asperities.

Another widely used methodology derives the hydraulic aperture from the mechanical aperture, , and the joint roughness coefficient as proposed by Barton et al. [7]: where is derived by comparing the fracture profile obtained with the Barton Comb with the standard tables provided by Barton and Choubey [8]. This methodology is perhaps the simplest and cheapest way to obtain fracture surface roughness values and has been widely used in outcrop studies [9–13]. The disadvantages of this method are related to the moderate resolution (about 1 mm) and the inaccuracy of equation (4) at relatively wide apertures (with respect to the JRC value). For instance, considering a fracture with 100 μm mechanical aperture and a JRC equal to 2, equation (4) gives a hydraulic aperture equal to 1767 μm.

Considering the previous arguments, the main objective of this work is to find empirical equations that describe the effect of fracture roughness on permeability at different apertures. In order to reach this goal, some problems which should be overcome are (i) to develop a protocol for mapping the fracture surface, (ii) to evaluate the fracture roughness as a function of the wavelength of the asperities, and (iii) to validate the relationships using a significant number of samples, roughness values, and aperture scenarios.

Various approaches have been reported in the literature for mapping the surface of fractures and faults in the field or laboratory involving the use of Lidar [14, 15], laboratory profilometers [16–18], and SfM photogrammetry (e.g., [12, 19]). Corradetti et al. [19] applied SfM photogrammetry for mapping fracture surfaces obtaining 3D reconstructions with point-cloud densities of equal quality to Lidar-derived data. Among these methods, SfM photogrammetry shows great promise as it is inexpensive (photo-camera and processing software) and extremely flexible regarding the scales and conditions (applicable in the field and laboratory). For example, SfM photogrammetry has been successfully used as an analytical tool to gather geologic data from outcrop studies [20–23], as well as at smaller scales on fault surfaces [19] and fossilized human footprints [24].

Evaluation of fracture roughness is achieved by implementing the power spectral density (PSD), which provides a more objective description based on the frequency distribution of the asperities in the Fourier domain. This approach has been successfully applied by previous authors for describing the roughness of fractures (e.g., [25]) and fault surfaces [14, 15, 18, 19]. To increase the statistical significance of the results, approximately 2000 fractures were modeled using the software SYNFRAC [25–27] creating computer-generated synthetic fractures, using input parameters of the fractal dimension, fracture roughness (an output of the PSD analysis of real fractures), and the standard deviation of the asperities height.

A key benefit of incorporating computer-generated synthetic fractures is the capability to work with a large amount of fracture data to perform direct fluid flow simulations, such as (i) the finite difference method (e.g., [28]), (ii) the finite element method (e.g., [29, 30]), (iii) the finite volume method (e.g., [31]), and the lattice-Boltzmann method (e.g., [32]). The lattice-Boltzmann method (LBM) describes the flow of many particles interacting with a medium and among themselves following the Navier-Stokes equation at the macroscopic scale (e.g., [33]). The LBM has been implemented to compute permeability using 3D images of rocks and soft sediments obtained by micro-CT imaging techniques [34–38] and from reconstructed models [39–41] generating results consistent with laboratory measurements [40, 42]. The simplest LBM is based on the Bhatnagar-Gross-Krook (BGK) collision operator, which consists of a single relaxation time approximation [43]. Despite its widespread use, the BGK-LBM brings some issues, for example, the computed permeability values may be viscosity-dependent [44]. An alternative approach involves the use of multiple relaxation times (MRT) methods, which solve the drawbacks of the BGK method and are characterized by more stability [45–47].

In this study, the permeability values of single isolated fractures (synthetic and natural) were calculated via LBM, using the PALABOS open source library [48]. The method and the code have been previously implemented by Zambrano et al. [38] for quantifying the permeability in deformed porous carbonates using X-ray microtomography synchrotron-based images [49]. Following these authors, rather than the BGK method (e.g., [34]), the MRT approach has been adopted in the present study to assure that permeability values are viscosity-independent.

The selected study area, the Roman Valley Quarry (Figure 1(a)), is an inactive quarry located at the northern termination of the Majella Mountain (Abruzzo region, Italy). The Majella Mountain is the orogenic expression of a thrust-related anticline, with internal deformation characterized by high-angle normal, strike-slip and oblique-slip normal faults, small folds, multiple sets of opening-mode fractures, pressure solution seams, and deformation bands [50–57].

(a)

(b)

The Roman Valley quarry has been heavily studied by previous authors focusing on the structural, sedimentological, and diagenetic properties [54, 56], and the fluid flow behavior of the fractured carbonates at the macroscale [13]. Here, the bitumen distribution suggests that the main hydrocarbon flow occurred through the damage zones of the principle NW-SE-oriented oblique slip faults [54] (Figure 1(b)). The distribution of major lithofacies at the Roman Valley Quarry is another element influencing the presence of bitumen which has been previously described by Rustichelli et al. [56] (see Table 1). Concerning the fractures, the most pervasive ones are represented by both pressure solution seams (often sheared with sliding/tearing mode displacement) and joints (opening mode fractures).

Considering the significance of these fractures, this work focused on investigating both cases of open mode and sheared fractures with a small sliding/tearing mode displacement, in the order of millimeters, allowing the assumtion of a negligible wall wearing. For this last case, the mismatch between opposite walls was also computed due to its importance as a mechanism for maintaining fracture openings even at reservoir depths.

2. Methodology

In this work, we present a multiphase integrated methodology for characterizing fracture surfaces and their effect on permeability. This approach combines fracture surface scanning using structure from motion photogrammetry, a statistical and spectral description of individual natural fracture surfaces, modeling of synthetic fractures, and computational calculation of permeability by fluid flow simulation.

2.1. Sample Collection

During the summer of 2018, a suite of oriented hand samples was collected from the study site comprising three (i.e., Au, B, and C) of the four major lithofacies of the Bolognano Formation present in the quarry (Figure 2, Table 1). The field sampling procedure involved manually removing blocks containing fracture surfaces that showed minimal signs of physical and chemical weathering. Sampling targeted these specific lithofacies based on accessibility, quality of well-developed fracture surfaces, and the fact that they have been well documented in previous studies focused on fracture distribution and mechanical properties [13, 54, 56]. After removal from the outcrop but prior to analysis, the rock samples were cleaned using a soft-bristled brush to remove debris and other obstructions but without abrading the surface.

(a)

(b)

(c)

(d)

2.2. Mapping Surface Topography

The workflow for mapping surface topography involves the following key stages: (1) fracture surface image set acquisition and (2) image alignment and three-dimensional digital rock model creation using SfM.

2.2.1. Image Set Acquisition

SfM photo scanning was performed at the University of Camerino photogrammetry laboratory (Figure 3(a)). We used a tripod-mounted Canon EOS 100D with the standard kit lens fixed at 55 mm. Images of samples were recorded inside a photo-lightbox to maintain full control of camera positions and ambient lighting and to reduce shadows and glare ensuring high image quality. To achieve the desired >70% inter-photo overlap, fracture samples were placed on a 360-degree rotating stage and manually rotated at 10-degree angular increments between photos. After completing a full orbit of the object, the camera was reset at a new vertical position, and the next orbit was conducted adhering to the same 10-degree increments but offset one position from the initial starting point. This procedure was repeated along three horizontal rotations from different vertical positions followed by 2-3 photos directly above the object oriented normal to the fracture surface.

(a)

(b)

(c)

(d)

2.2.2. Image Alignment and Three-Dimensional Model Creation

Fracture surface models were aligned and processed using Agisoft PhotoScan Pro (http://www.agisoft.com). For each fracture surface, approximately 63 photos were used as input images to create the digital point cloud model. We follow the procedure described by Carrivick et al. [58] and Zimmer et al. [24] for camera settings and photo procedure, and Pitts et al. [22] for image alignment and point cloud generation.

Agisoft-generated coded targets were placed inside the scene to aid in the imagery processing; these coded targets are automatically recognized by the software and help build connection points between the image sets (Figure 3(b)). Additionally, a 5-centimeter scale was placed on the sample surface to calibrate the spatial reference.

As a measure to define the error of the model, we follow the methodology established by Corradetti et al. [19]. This calls for modeling a piece of graph paper under the same condition as fracture imaging (same light, relative distance, and number of photos). Under the assumption that the graph paper is perfectly planar, the standard deviation of the height of the scanned asperities is considered as the error of the model [19]. In our case, this value is approximately 20 μm, whereas the point density (points/area) is 34 points/mm.

2.3. Fracture Surface Processing and Analysis

2.3.1. Extraction and Processing of 3D Fracture Surface Model

Trimmed 3D point clouds were exported from Agisoft as “.xyz” text files. Then, a rectangular subregion of each fracture surface of interest was extracted from the point cloud and processed to remove undesirable trend and eventual noise (Figure 4). This technique, previously documented by Corradetti et al. [19], consists of (i) removing the artificial trend of the surface, (ii) surface interpolation, and (ii) sampling in a regular grid.

(a)

(b)

2.3.2. Fracture Roughness Assessment

A complete description of the fracture roughness is given by the specification of two functions: the probability density function (depending on the media and standard deviation) for heights and the PSD [59] (Figure 5).

(a)

(b)

The Fourier power spectrum, , is defined as the square of the modulus of the Fourier transform [60]. Considering a cross section of the rock, this profile can be represented as a summation of sinusoidal components, each with its own amplitude, wavelength, and phase. The squared amplitude of each sinusoid component is referred to as its “power” (Figure 5(b)). The power spectrum regulated in an appropriate manner is referred to as the PSD, , and it provides a valuable definition of surface roughness. The PSD as a function of in a bi-logarithmic scale graph of a self-affine function exhibits an apparent linear slope, which is defined from the following power law equation: where the exponent of the power law is related with the fractal dimension, [59], as follows:

From a physical aspect, the fractal dimension shows the proportion of high-frequency to low-frequency sinusoid components (roughness). High values are related to greater surface roughness. By stacking and normalizing the power spectra, it is possible to reduce the noise associated with a single profile and create a single spectrum which represents the entire rough surface in each direction (e.g., [15, 19]). The MATLAB script used to perform the procedures described above is partially modified from Corradetti et al. [19].

2.4. Single Fracture Modeling

Since a limited number of natural fracture surfaces were available, additional synthetic fracture surfaces were used to strengthen the statistical significance of the results. Following the procedure described by Ogilvie et al. [25], more than 2000 computer-generated synthetic fractures were created using the software SYNFRAC [25–27]. SYNFRAC is based on a mathematical model of a rough surface reported by Brown [59]. The software can model open fractures by introducing mismatch values with the spatial and spectral roughness parameters. For the scope of the present study, the mismatch values were not considered for surface modeling, but were measured after the fracture modeling (see Section 2.5.1).

The individual fracture surfaces (natural and synthetic) were used to model dilation associated with (Figure 6): (i) opening mode displacement (joint and/or opened pressure solution seam), various ranges of aperture, and (ii) sliding/tearing mode displacement (sheared joint and/or sheared pressure solution seam). In both cases, it is assumed that both walls are identical. In the second case, because the shear process is minimum, and the displacement is in the order of mm, it is assumed that no physical wearing of the fracture surfaces has occurred. To illustrate different scenarios, a wide range of displacements (opening and sliding/tearing mode) were considered (the PYTHON code for fracture modeling is available at https://github.com/superrostom/synthetic-fracture).

(a)

(b)

2.5. Lattice-Boltzmann Method and Permeability Computation

Lattice-Boltzmann simulations were performed using the open-source computational fluid dynamics software PALABOS [48] following the methodology described by Zambrano et al. [38].

This procedure consists of imposing a single-phase fluid flow through a 3D porous media maintaining a fixed pressure gradient between the inlet and outlet opposing faces of the model; the rest of the faces were padded (Figure 7). A bounce-back boundary condition was assigned to the fracture surfaces. An MRT collisional operator [45, 46], with a D3Q19 lattice, is used instead of the popular BGK [43] as in Degruyter et al. [34]. Moreover, the geometry of the media is obtained by the SfM photogrammetry outputs and modeling differently than Degruyter et al. [34] and Zambrano et al. [38] who used X-ray micro CT images.

(a)

(b)

The simulation ended once the imposed steady-state condition was reached (standard deviation of the average after 1000 steps). Then, the permeability component parallel to the imposed flow was calculated applying Darcy’s law, where is the pressure gradient, is the fluid kinematic viscosity, and is the average fluid velocity per unit of area. The permeability was calculated, using the same procedure, in two orthogonal directions: along strike and along dip ( and , respectively). In the case of sliding/tearing mode displacement, the -direction corresponds to the slip direction. All the variables are handled in lattice units prior to permeability calculation. For convenience, permeability values were converted to millidarcy which is the most commonly used permeability unit in the oil industry. All the provided values of permeability correspond to a volume of ().

2.5.1. Mismatch Evaluation

The mismatch between the opposite fracture walls is of extreme importance since this factor may keep fractures open even at reservoir depths. Since the mismatch was not imposed during fracture modeling, it was measured after the generation of the synthetic fractures. The mismatch was evaluated only for the sliding/tearing mode fractures, whereas it was unnecessary in the case of opening mode fractures since the aperture is constant. For the evaluation of the mismatch value, the methodology of the power spectral density ratio (PSDr), introduced by Ogilvie et al. [25], was followed. The methodology consists of obtaining a relationship between the PSD of the aperture and both surfaces of the fracture, as follows:

The results of this calculation can be represented in a graph where the parameters associated with the mismatch and the degree of mismatch between the surfaces at different wavelengths can be obtained (Figure 8).

Following the definition of Ogilvie et al. [25], these parameters are the following: (i)Minimum mismatch length (ML_min): wavelength at which the fractures start to match, indicated by the wavelength where the PSD ratio values fall below its maximum value (PSDr_max)(ii)Maximum mismatch length (ML_max): wavelength at which the fracture opposing surfaces reach the maximum matching, thus the minimum value of PSD ratio (PSDr_min)

The calculation of these parameters was made using a MATLAB code. In this case, (ML_min) is considered as the only reliable indicator of the mismatch since (ML_max) often falls outside the scale of the study (Figure 8).

3. Results

The results of this work consist of an analysis of surface topography performed on fracture samples from three lithofacies (Au, B, and C), and the computed permeability in function of the fracture properties, including fractal dimension, opening and sliding/tearing displacement, and minimum mismatch length.

3.1. Fracture Surface Properties

In Table 2, the values correspond to the fractal dimension (), the average height of asperities, and their standard deviation. The traditional roughness measurement, JRC, was added to compare both techniques and the results with previous works in the same outcrop. For the fracture description, we followed Agosta et al. [54]. The set PS2a corresponds to pressure solution seams (generated during background deformation) often sheared, with normal or left lateral kinematics, often impregnated by tar. The set PS2b corresponds to pressure solution seams (generated during background deformation) which are generally un-sheared.

3.2. Permeability Results

The results of the work are presented in Figures 9–14. In each graph, the different surface roughness values (expressed in fractal dimension) are illustrated. When pertinent, error bars are added to show the variability of the results.

(a)

(b)

(a)

(b)

(a)

(b)

(a)

(b)

(a)

(b)

In the case of opening mode displacement, the results indicate that the permeability increased proportionally to the mechanical aperture following a positive power-law relationship (Figure 9(a)). Similarly, the hydraulic aperture (derived from equation (1)) is related with the mechanical aperture by means of a positive power law with exponents varying from 1.6 to 1.8 for fractal dimension () values between 1.6 and 2.5, respectively (Figure 9(b)).

With respect to the case of sliding/tearing displacement, the results indicate that the permeability component parallel to the displacement () increases proportionally to the sliding/tearing displacement () following a positive power-law relationship (Figure 10(a)). The permeability component perpendicular to the displacement () is also related by a power law with the sliding/tearing displacement (). The anisotropy ratio, , is generally higher for low values of fractal dimension (Figure 10(b)). The anisotropy ratio value tends to decrease as a function of the sliding/tearing displacement () following a negative power-law relationship. The highest value of anisotropy was near 2.6 for fractal dimension () equal to 2.0 and sliding/tearing displacement () equal to 0.5 mm.

The fracture roughness, expressed in terms of fractal dimension (), showed a different influence in the permeability on single fractures depending on their kinematic: opening or sliding/tearing mode. For the opening mode displacement case, the permeability is inversely proportional to the fractal dimension (Figure 11(a)). This relationship follows a negative power law with the slope depending on the opening mode displacement values. A higher roughness (fractal dimension) implies a decrease in permeability. For the simulated scenarios, rough fractures () showed permeability values between 45% and 65% lower than smooth fractures (). For fractures with the sliding/tearing displacement, the permeability is proportional to the fractal dimension following a positive power law (Figure 11(b)). In this case, an increment of the fractal dimension (roughness) from 2.0 to 2.5 represents an enhancement of the permeability of approximately one order of magnitude.

It is expected that this positive relationship between displacement and permeability should stabilize at a certain point, as the permeability and porosity cannot increase indefinitely. This behavior is observed when the porosity is evaluated at higher displacement values (Figure 11). In this case, thanks to the simplified calculation of porosity in comparison with permeability, a large volume of data was considered (more than 2000 fractures). The porosity is proportional to the sliding/tearing displacement following a nonlinear relationship with variable slope. The most significant change of slope occurs near 10 mm and after approximately 30 mm of sliding/tearing displacement, where the porosity seems to stabilize at values between 3% and 4%. The fractures with higher fractal dimension (roughness) tend to have higher porosity, which agrees with the permeability results. As previously mentioned, these results do not consider the possible physical wearing of the surfaces due to shearing particularly at high sliding/tearing displacement.

The fracture permeability is related to the porosity following a power law (Figure 13), which differs from the theoretical relationship based on the smooth parallel-plate equation; here, we also assumed smooth parallel plates for the porosity. The power-law relationship seems to be unaffected by the roughness (fractal dimension). Similar behaviors were obtained for both permeability components and .

The minimum mismatch length, ML_min, was evaluated as a function of the displacement (Figure 14(a)). Results indicate that within the evaluated range (a maximum displacement of 5 mm), the minimum mismatch is linearly proportional to the displacement. The permeability is proportional to the minimum mismatch length following a power-law relationship (depending on the fractal dimension). In other words, higher values of fractal dimension imply a higher permeability for similar mismatch values.

4. Discussion

The present work evaluates the effect of fracture surface features such as roughness, aperture, and mismatch on permeability using fracture surface scanning by SfM photogrammetry, numerical modeling, and lattice-Boltzmann fluid flow simulation.

4.1. SfM Photogrammetry Surface Scanning

The results of this study demonstrate the versatility of the SfM procedure as an analytical tool which can be applied at a wide range of scales including millimeter-scale features such as fracture surfaces. The controlled conditions in the photogrammetry laboratory allowed a highly detailed scan and extraction of the micro surface topography of samples sized producing a point cloud with a density of 34 points/mm and an estimated error of 20 μm. This method produces more realistic and applicable results than the traditional Barton Comb, with results comparable to those reported by Candela et al. [14], Renard et al. [18], and Corradetti et al. [19] using both Lidar or laser profilometers. However, the SfM methodology is several orders of magnitude more cost-effective and is readily accessible. A future implementation of this study could include developing a workflow suitable for in situ field studies. However, more variables need to be controlled and results yielding lower accuracy are expected.

4.2. Fracture Roughness Characterization

This methodology proved to be highly efficient in expressing the fracture roughness allowing a more accurate and representative measure with respect to the relative hydraulic roughness [4–6] or the JRC [7, 8]. However, more data is necessary for evaluating control of the lithofacies in the fracture roughness, as it is observed for other properties such as distribution and spacing (e.g., [61]).

Another important aspect of these results is the reproducibility of synthetic fractures with similar characteristics. This step was important to increase the data volumes leading to a greater statistical significance of the results and validity of the inferred relationships. The lattice-Boltzmann procedure also played a key role in this study as it allows the estimation of permeability values for controlled scenarios with different imposed properties (i.e., roughness, opening mode, and sliding/tearing displacement). This permits evaluation of the relationship between permeability, porosity, mismatch, and other imposed properties. The computed permeability may present some inaccuracy in low-resolution models as previously reported by Zambrano et al. [38]. Therefore, permeability results may be considered not as the real values but only as approximations.

4.3. Permeability in Function of Fracture Properties

Two situations were considered to explain the presence of open fractures: (i) dilation due to opening mode displacement (joint or opened pressure solution seam) and (ii) dilation due to mismatch caused by shearing and sliding/tearing displacement.

In the first case, the results followed the expectation and confirmed previous interpretations: (i) permeability tends to increase with opening following a nonlinear relationship, (ii) a higher fractal dimension (greater roughness) correlates to lower permeability, and (iii) the effect of roughness is less significant at greater opening values. It is expected that higher roughness (higher frequencies of asperities) may expose a wider area in contact with the migrating fluid, diminishing its velocity due to friction. Evidently, at higher opening values, this effect should be less evident because it is the specific area (area/volume) which has a control on the permeability, as has been previously reported by Zambrano et al. [38] for porous media.

The second case (sliding/tearing displacement) creates an aperture due to the mismatch between the opposite walls of the fracture. We found significant differences between these two cases concerning the effect of the roughness of the permeability. In fact, the effect of roughness on permeability is inverse. Given the same displacement, fractures with higher roughness values permit the creation of larger voids and therefore enhance the fluid flow. So, the effect of friction exerted by the roughness on fluid flow has a secondary role in the case of mismatch due to sliding/tearing displacement. The continuous dilatancy of the fracture due to sliding/tearing displacement should cease at a certain value depending on the asperity frequencies present in the fracture. Nevertheless, it is difficult to verify this behavior for real fractures, where a 20 mm sliding/tearing displacement likely leads to fracture wall wearing and the generation of cataclastic material [54], eventually reducing permeability. However, our model has greater applicability to small displacements where the damage of the fracture walls is negligible.

The permeability anisotropy in fractures with sliding/tearing displacement is significant and dependent on the roughness (fractal dimension). For low displacement (0.5 mm), the anisotropy can reach values up to 2.6 for fractures with a fractal dimension of 2.0. For the same displacement, fractures with high roughness () showed lower values of permeability anisotropy near 1.5. For higher sliding/tearing displacement, the permeability anisotropy decreases approaching the value 1.

The mismatch itself has a positive control on the permeability. The importance of this result is that the mismatch could also be produced by diagenetic processes (e.g., cementation, dissolution) and shearing. Zambrano et al. [38], using X-ray microtomographic images in shear compactive bands hosted in porous carbonates, showed the existence of complex channelized porous networks along the shear surface within these structures.

4.4. Consequences to Reservoir Modeling

The results agree with the macroscale observations of previous authors in the study area, where both opened/sheared pressure solution seams and fault-related joints present the greatest values of aperture and the most important bitumen impregnation [13, 54, 56].

The relationship between permeability and porosity for rough fractures clearly deviates from the ideal smooth parallel plate case (for the studied scenarios). Fracture permeability is lower for the same porosity range (<0.2%) in comparison to the theoretical values. Instead, the power-law slope is higher, indicating a more important control of porosity as it was expected. The equation itself may be useful to estimate the permeability of fractures if the fracture porosity is known.

After their formation, both closing and opening mode fractures are often subjected to a shear process, and even with a small imperceptible sliding/tearing displacement, they cannot be modeled as simple opening mode fracture. At reservoir depth, preexisting fractures (joints and pressure solution seams) that are favorably oriented to be sheared (accordingly to the orientations of the stress field which affected the area) may be characterized by a mismatch between the fracture walls enhancing the fracture opening. Therefore, the findings of this work may have a significant impact on fracture modeling workflows for subresolution faults (e.g., [62]). In fractures with small shear displacement, contrary to the open mode case, the roughness may influence positively the permeability. Sliding/tearing displacement also may enhance the fracture permeability eventually decelerating when the cataclasis process starts. If the shear displacement is small (), the permeability anisotropy is significant enough to be considered in the fracture modeling workflow.

5. Conclusions

We presented a new multifaceted approach to characterize surface fracture roughness by SfM photogrammetry, numerical modeling, and computational fluid dynamics simulation. This methodology provides a better quantification of surface parameters that are not possible to obtain using former surface roughness measurement and analysis tools.

In addition, this study illustrates the crucial relationships between permeability and other fracture properties, such as roughness, porosity, opening mode-sliding/tearing displacement, and mismatch. The obtained relationships pointed out the following statements: (i)In joints (opening mode fractures) and/or opened pressure solution seams, the roughness tends to reduce the permeability. Thus, the permeability is inversely proportional to the fractal dimension(ii)In sheared joints and/or pressure solution seams (assuming an insignificant surface wearing), the sliding/tearing mode displacement may cause mismatch and therefore enhance the porosity and permeability. The validity of this behavior may depend on the point that displacement starts to produce cataclastic material. Small shear displacements and mismatch may be extremely important to guarantee storage and migration of geofluid at depth thanks to asperities-supported aperture. Permeability anisotropy is very significant for small shear displacements, characterized by higher values of permeability component perpendicular to the shear displacement(iii)Porosity exerts a more important control on permeability in rough fractures (higher power-law slope). The empiric relationship may result in greater utility for estimating the fracture permeability if the fracture porosity is known

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This research was supported by the FAR Project 2014 “Characterization and modeling of natural reservoirs of geofluids in fractured carbonate rocks,” funded by the University of Camerino, Principal investigator Emanuele Tondi, and the Reservoir Characterization Project (http://www.rechproject.com). We also acknowledge and thank Amerigo Corradetti for generously sharing his code for power spectral density (PSD) analysis of surfaces and thoughtful conversations that helped move this work forward.

References

D. Snow, A parallel Plate Model of Fractured Permeable Media, [Ph.D. Thesis], University of California, Berkeley, 1965.
Y. W. Tsang, “Usage of “equivalent apertures” for rock fractures as derived from hydraulic and tracer tests,” Water Resources Research, vol. 28, no. 5, pp. 1451–1455, 1992.
View at: Publisher Site | Google Scholar
P. A. Witherspoon, J. S. Y. Wang, K. Iwai, and J. E. Gale, “Validity of cubic law for fluid flow in a deformable rock fracture,” Water Resources Research, vol. 16, no. 6, pp. 1016–1024, 1980.
View at: Publisher Site | Google Scholar
G. M. Lomize, Flow in Fractured Rocks, Gosenergoizdat, Moscow, 1951.
C. Louis, “A study of ground water flow in jointed rock and its influence on the stability of rock masses,” Rock Mechanics Research Report, Imperial College, London, 1969, Rock Mechanics Research Report 10.
View at: Google Scholar
E. F. de Quadros, Determinação das características do fluxo de água em fraturas de rochas, Department of Civil Construction Engineering, Polytechnic School, University of Sao Paulo, 1982.
N. Barton, S. Bandis, and K. Bakhtar, “Strength, deformation and conductivity coupling of rock joints,” International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts, vol. 22, no. 3, pp. 121–140, 1985.
View at: Publisher Site | Google Scholar
N. Barton and V. Choubey, “The shear strength of rock joints in theory and practice,” Rock Mechanics, vol. 10, no. 1-2, pp. 1–54, 1977.
View at: Publisher Site | Google Scholar
M. Zambrano, E. Tondi, I. Korneva et al., “Fracture properties analysis and discrete fracture network modelling of faulted tight limestones, Murge Plateau, Italy,” Italian Journal of Geosciences, vol. 135, no. 1, pp. 55–67, 2016.
View at: Publisher Site | Google Scholar
E. Panza, F. Agosta, M. Zambrano et al., “Structural architecture and discrete fracture network modelling of layered fractured carbonates (Altamura Fm., Italy),” Italian Journal of Geosciences, vol. 134, no. 3, pp. 409–422, 2015.
View at: Publisher Site | Google Scholar
E. Panza, F. Agosta, A. Rustichelli et al., “Fracture stratigraphy and fluid flow properties of shallow-water, tight carbonates: the case study of the Murge Plateau (southern Italy),” Marine and Petroleum Geology, vol. 73, pp. 350–370, 2016.
View at: Publisher Site | Google Scholar
D. H. Kim, I. Gratchev, and A. Balasubramaniam, “Determination of joint roughness coefficient (JRC) for slope stability analysis: a case study from the Gold Coast area, Australia,” Landslides, vol. 10, no. 5, pp. 657–664, 2013.
View at: Publisher Site | Google Scholar
T. Volatili, M. Zambrano, A. Cilona et al., “From fracture analysis to flow simulations in fractured carbonates: the case study of the Roman Valley Quarry (Majella Mountain, Italy),” Marine and Petroleum Geology, vol. 100, pp. 95–110, 2019.
View at: Publisher Site | Google Scholar
T. Candela, F. Renard, M. Bouchon et al., “Characterization of fault roughness at various scales: implications of three-dimensional high resolution topography measurements,” Pure and Applied Geophysics, vol. 166, no. 10-11, pp. 1817–1851, 2009.
View at: Publisher Site | Google Scholar
T. Candela, F. Renard, Y. Klinger, K. Mair, J. Schmittbuhl, and E. E. Brodsky, “Roughness of fault surfaces over nine decades of length scales,” Journal of Geophysical Research, vol. 117, no. B8, article B08409, 2012.
View at: Publisher Site | Google Scholar
F. Renard, C. Voisin, D. Marsan, and J. Schmittbuhl, “High resolution 3D laser scanner measurements of a strike-slip fault quantify its morphological anisotropy at all scales,” Geophysical Research Letters, vol. 33, no. 4, article L04305, 2006.
View at: Publisher Site | Google Scholar
F. Renard, K. Mair, and O. Gundersen, “Surface roughness evolution on experimentally simulated faults,” Journal of Structural Geology, vol. 45, pp. 101–112, 2012.
View at: Publisher Site | Google Scholar
F. Renard, T. Candela, and E. Bouchaud, “Constant dimensionality of fault roughness from the scale of micro-fractures to the scale of continents,” Geophysical Research Letters, vol. 40, no. 1, pp. 83–87, 2013.
View at: Publisher Site | Google Scholar
A. Corradetti, K. McCaffrey, N. De Paola, and S. Tavani, “Evaluating roughness scaling properties of natural active fault surfaces by means of multi-view photogrammetry,” Tectonophysics, vol. 717, pp. 599–606, 2017.
View at: Publisher Site | Google Scholar
S. P. Bemis, S. Micklethwaite, D. Turner et al., “Ground-based and UAV-based photogrammetry: a multi-scale, high-resolution mapping tool for structural geology and paleoseismology,” Journal of Structural Geology, vol. 69, pp. 163–178, 2014.
View at: Publisher Site | Google Scholar
K. J. W. McCaffrey, R. R. Jones, R. E. Holdsworth et al., “Unlocking the spatial dimension: digital technologies and the future of geoscience fieldwork,” Journal of the Geological Society of London, vol. 162, no. 6, pp. 927–938, 2005.
View at: Publisher Site | Google Scholar
A. D. Pitts, C. I. Casciano, M. Patacci, S. G. Longhitano, C. Di Celma, and W. D. McCaffrey, “Integrating traditional field methods with emerging digital techniques for enhanced outcrop analysis of deep water channel-fill deposits,” Marine and Petroleum Geology, vol. 87, pp. 2–13, 2017.
View at: Publisher Site | Google Scholar
P. R. Nesbit, P. R. Durkin, C. H. Hugenholtz, S. M. Hubbard, and M. Kucharczyk, “3-D stratigraphic mapping using a digital outcrop model derived from UAV images and structure-from-motion photogrammetry,” Geosphere, vol. 14, no. 6, pp. 2469–2486, 2018.
View at: Publisher Site | Google Scholar
B. Zimmer, C. Liutkus-Pierce, S. T. Marshall, K. G. Hatala, A. Metallo, and V. Rossi, “Using differential structure-from-motion photogrammetry to quantify erosion at the Engare Sero footprint site, Tanzania,” Quaternary Science Reviews, vol. 198, pp. 226–241, 2018.
View at: Publisher Site | Google Scholar
S. R. Ogilvie, E. Isakov, and P. W. Glover, “Fluid flow through rough fractures in rocks. II: a new matching model for rough rock fractures,” Earth and Planetary Science Letters, vol. 241, no. 3-4, pp. 454–465, 2006.
View at: Publisher Site | Google Scholar
E. Isakov, S. R. Ogilvie, C. W. Taylor, and P. W. Glover, “Fluid flow through rough fractures in rocks I: high resolution aperture determinations,” Earth and Planetary Science Letters, vol. 191, no. 3-4, pp. 267–282, 2001.
View at: Publisher Site | Google Scholar
S. R. Ogilvie, E. Isakov, C. W. Taylor, and P. W. J. Glover, “Characterization of rough-walled fractures in crystalline rocks,” Geological Society, London, Special Publications, vol. 214, no. 1, pp. 125–141, 2003.
View at: Publisher Site | Google Scholar
M. J. Blunt, B. Bijeljic, H. Dong et al., “Pore-scale imaging and modelling,” Advances in Water Resources, vol. 51, pp. 197–216, 2013.
View at: Publisher Site | Google Scholar
M. B. Cardenas, “Three-dimensional vortices in single pores and their effects on transport,” Geophysical Research Letters, vol. 35, no. 18, article L18402, 2008.
View at: Publisher Site | Google Scholar
M. B. Cardenas, “Direct simulation of pore level Fickian dispersion scale for transport through dense cubic packed spheres with vortices,” Geochemistry, Geophysics, Geosystems, vol. 10, no. 12, article Q12014, 2009.
View at: Publisher Site | Google Scholar
B. Bijeljic, A. Raeini, P. Mostaghimi, and M. J. Blunt, “Predictions of non-Fickian solute transport in different classes of porous media using direct simulation on pore-scale images,” Physical Review E, vol. 87, no. 1, article 013011, 2013.
View at: Publisher Site | Google Scholar
M. C. Sukop, H. Huang, C. L. Lin, M. D. Deo, K. Oh, and J. D. Miller, “Distribution of multiphase fluids in porous media: comparison between lattice Boltzmann modeling and micro-x-ray tomography,” Physical Review E, vol. 77, no. 2, article 026710, 2008.
View at: Publisher Site | Google Scholar
A. J. C. Ladd, “Numerical simulations of particulate suspensions via a discretized Boltzmann equation: part 2. Numerical results,” Journal of Fluid Mechanics, vol. 271, pp. 311–339, 1994.
View at: Publisher Site | Google Scholar
W. Degruyter, A. Burgisser, O. Bachmann, and O. Malaspinas, “Synchrotron X-ray microtomography and lattice Boltzmann simulations of gas flow through volcanic pumices,” Geosphere, vol. 6, no. 5, pp. 470–481, 2010.
View at: Publisher Site | Google Scholar
F. Khan, F. Enzmann, M. Kersten, A. Wiegmann, and K. Steiner, “3D simulation of the permeability tensor in a soil aggregate on basis of nanotomographic imaging and LBE solver,” Journal of Soils and Sediments, vol. 12, no. 1, pp. 86–96, 2012.
View at: Publisher Site | Google Scholar
H. Andrä, N. Combaret, J. Dvorkin et al., “Digital rock physics benchmarks—part II: computing effective properties,” Computational Geosciences, vol. 50, pp. 33–43, 2013.
View at: Publisher Site | Google Scholar
S. M. Shah, F. Gray, J. P. Crawshaw, and E. S. Boek, “Micro-computed tomography pore-scale study of flow in porous media: effect of voxel resolution,” Advances in Water Resources, vol. 95, pp. 276–287, 2015.
View at: Publisher Site | Google Scholar
M. Zambrano, E. Tondi, L. Mancini et al., “Fluid flow simulation and permeability computation in deformed porous carbonate grainstones,” Advances in Water Resources, vol. 115, pp. 95–111, 2018.
View at: Publisher Site | Google Scholar
G. Jin, T. W. Patzek, and D. B. Silin, “Direct prediction of the absolute permeability of unconsolidated consolidated reservoir rock,” in Proceedings of SPE Annual Technical Conference and Exhibition, pp. 1–15, Houston, TX, USA, September 2004.
View at: Publisher Site | Google Scholar
Y. Keehm, T. Mukerji, and A. Nur, “Permeability prediction from thin sections: 3D reconstruction and lattice-Boltzmann flow simulation,” Geophysical Research Letters, vol. 31, no. 4, article L04606, 2004.
View at: Publisher Site | Google Scholar
K. Wu, M. I. J. Dijke, G. D. Couples et al., “3D stochastic modelling of heterogeneous porous media – applications to reservoir rocks,” Transport in Porous Media, vol. 65, no. 3, pp. 443–467, 2006.
View at: Publisher Site | Google Scholar
Y. Keehm, Computational Rock Physics: Transport Properties in Porous Media and Applications, [Ph.D. Thesis], Stanford University, Stanford, CA, USA, 2003.
P. L. Bhatnagar, E. P. Gross, and M. Krook, “A model for collision processes in gases. I. Small amplitude processes in charged and neutral one-component systems,” Physics Review, vol. 94, no. 3, pp. 511–525, 1954.
View at: Publisher Site | Google Scholar
A. Narvaez, T. Zauner, F. Raischel, R. Hilfer, and J. Harting, “Quantitative analysis of numerical estimates for the permeability of porous media from lattice-Boltzmann simulations,” Journal of Statistical Mechanics: Theory and Experiment, vol. 2010, no. 11, article P11026, 2010.
View at: Publisher Site | Google Scholar
D. d’Humières, “Generalized lattice-Boltzmann equations,” in Rarefied Gas Dynamics: Theory and Simulations, Progress in Astronautics and Aeronautics, B. D. Shizgal and D. P. Weave, Eds., vol. 159, pp. 450–458, AIAA, Washington, D. C, USA, 1992.
View at: Publisher Site | Google Scholar
D. d'Humières, I. Ginzburg, M. Krafczyk, P. Lallemand, and L.-S. Luo, “Multiple–relaxation–time lattice Boltzmann models in three dimensions,” Philosophical Transactions of the Royal Society A, vol. 360, no. 1792, pp. 437–451, 2002.
View at: Publisher Site | Google Scholar
C. Pan, L. S. Luo, and C. T. Miller, “An evaluation of lattice Boltzmann schemes for porous medium flow simulation,” Computers & Fluids, vol. 35, no. 8-9, pp. 898–909, 2006.
View at: Publisher Site | Google Scholar
J. Latt, Palabos, Parallel Lattice Boltzmann Solver, FlowKit, Lausanne, Switzerland, 2009, http://www.palabos.org/.
M. Zambrano, E. Tondi, L. Mancini et al., “3D pore-network quantitative analysis in deformed carbonate grainstones,” Marine and Petroleum Geology, vol. 82, pp. 251–264, 2017.
View at: Publisher Site | Google Scholar
F. Ghisetti and L. Vezzani, “Normal faulting, extension and uplift in the outer thrust belt of the central Apennines (Italy): role of the Caramanico fault,” Basin Research, vol. 14, no. 2, pp. 225–236, 2002.
View at: Publisher Site | Google Scholar
V. Scisciani, E. Tavarnelli, and F. Calamita, “The interaction of extensional and contractional deformations in the outer zones of the Central Apennines, Italy,” Journal of Structural Geology, vol. 24, no. 10, pp. 1647–1658, 2002.
View at: Publisher Site | Google Scholar
E. Tondi, M. Antonellini, A. Aydin, L. Marchegiani, and G. Cello, “The role of deformation bands, stylolites and sheared stylolites in fault development in carbonate grainstones of Majella Mountain, Italy,” Journal of Structural Geology, vol. 28, no. 3, pp. 376–391, 2006.
View at: Publisher Site | Google Scholar
A. P. Lavenu, J. Lamarche, R. Salardon, A. Gallois, L. Marié, and B. D. Gauthier, “Relating background fractures to diagenesis and rock physical properties in a platform–slope transect. Example of the Maiella Mountain (central Italy),” Marine and Petroleum Geology, vol. 51, pp. 2–19, 2014.
View at: Publisher Site | Google Scholar
F. Agosta, M. Alessandroni, E. Tondi, and A. Aydin, “Oblique normal faulting along the northern edge of the Majella Anticline, central Italy: inferences on hydrocarbon migration and accumulation,” Journal of Structural Geology, vol. 31, no. 7, pp. 674–690, 2009.
View at: Publisher Site | Google Scholar
A. Aydin, M. Antonellini, E. Tondi, and F. Agosta, “Deformation along the leading edge of the Maiella thrust sheet in central Italy,” Journal of Structural Geology, vol. 32, no. 9, pp. 1291–1304, 2010.
View at: Publisher Site | Google Scholar
A. Rustichelli, E. Tondi, F. Agosta, C. Di Celma, and M. Giorgioni, “Sedimentologic and diagenetic controls on pore-network characteristics of Oligocene–Miocene ramp carbonates (Majella Mountain, central Italy),” AAPG Bulletin, vol. 97, no. 3, pp. 487–524, 2013.
View at: Publisher Site | Google Scholar
A. Iadanza, G. Sampalmieri, and P. Cipollari, “Deep-seated hydrocarbons in the seep “brecciated limestones” of the Maiella area (Adriatic foreland basin): evaporitic sealing and oil re-mobilization effects linked to the drawdown of the Messinian salinity crisis,” Marine and Petroleum Geology, vol. 66, pp. 177–191, 2015.
View at: Publisher Site | Google Scholar
J. L. Carrivick, M. W. Smith, and D. J. Quincey, Structure from Motion in the Geosciences, John Wiley & Sons, 2016.
View at: Publisher Site
S. R. Brown, “Simple mathematical model of a rough fracture,” Journal of Geophysical Research: Solid Earth, vol. 100, no. B4, pp. 5941–5952, 1995.
View at: Publisher Site | Google Scholar
A. Schuster, “On the investigation of hidden periodicities with application to a supposed 26 day period of meteorological phenomena,” Journal of Geophysical Research, vol. 3, no. 1, p. 13, 1898.
View at: Publisher Site | Google Scholar
A. Rustichelli, F. Agosta, E. Tondi, and V. Spina, “Spacing and distribution of bed-perpendicular joints throughout layered, shallow-marine carbonates (Granada Basin, southern Spain),” Tectonophysics, vol. 582, pp. 188–204, 2013.
View at: Publisher Site | Google Scholar
M. Antonellini, A. Cilona, E. Tondi, M. Zambrano, and F. Agosta, “Fluid flow numerical experiments of faulted porous carbonates, northwest Sicily (Italy),” Marine and Petroleum Geology, vol. 55, pp. 186–201, 2014.
View at: Publisher Site | Google Scholar
J. Ahrens, B. Geveci, and C. Law, ParaView: an End-User Tool for Large Data Visualization, Visualization Handbook, Elsevier, 2005.

Copyright

Copyright © 2019 Miller Zambrano 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.

PDF Download Citation

Download other formats

Order printed copies

Views

2789

Downloads

2197

Citations