## Abstract

A statistical technique for the pore-scale analyses of heterogeneity and representative elemental volume (REV) in unconventional shale rocks is hereby presented. First, core samples were obtained from shale formations. The images were scanned using microcomputed tomography (micro-CT) machine at 6.7 μm resolution with voxels of 990 × 990 × 1000. These were then processed, digitised, thresholded, segmented and features captured using numerical algorithms. This allows the segmentation of each sample into four distinct morphological entities consisting of pores, organic matter, shale grains and minerals.

#### Author

##### James O. Adeleye, Lateef T. Akanji

Petroleum Engineering Division, School of Engineering, University of Aberdeen, Aberdeen, UK.

Received: 2 February 2017 / Accepted: 23 July 2017

© The Author(s) 2017.

In order to analyse the degree of heterogeneity, Eagle Ford parallel sample was further cropped into 96 subsamples. Descriptive statistical approach was then used to evaluate the existence of heterogeneity within the subsamples. Furthermore, the Eagle Ford parallel and perpendicular samples were analysed for volumetric entities representative of the petrophysical variable, porosity, using corner point cropping technique. The results of porosity REV for Eagle ford parallel and perpendicular indicated sample representation at 300 μm voxel edge. Both pore volume distribution and descriptive statistical analyses suggested that a wide variation of heterogeneity exists at this scale of investigation. Furthermore, this experiment allows for adequate extraction of necessary information and structural parameters for pore-scale modelling and simulation. Additional studies focusing on re-evaluation at higher resolution are recommended.

## Introduction

Pore geometry, tortuosity, grains size and shape are properties that are important to describe and characterise fluid flow in shale rock. However, getting a single material point measurement at which this shale rock property is determined depends on the ability to accurately extract information from its structure.

Determination of the nature and extent of heterogeneities at pore scale can enhance fluid characterisation in porous media. In shale rock, variability is extreme, because shales are detrital sediments formed as a result of alteration of mud or clay deposits. They occur as fine-grained rich in illite and fragments with particle sizes generally less than 0.062 mm (Leith 2016). Hence, its petrophysical properties (lithology, porosity, permeability water saturation) become difficult to evaluate.

Another major challenge of shale petrophysical evaluation is its heterogeneity as this makes quantification of porosity and other properties difficult. These parameters depend on knowing the details of how shape or size grains (Milner et al. 2010), pore throats (Curtis et al. 2012) and tortuosity (Bai et al. 2013; Katsube et al. 1991) are distributed. Several studies have reported the complex microstructure of the shale rocks (Curtis et al. 2012), heterogeneity of shale rocks (Chen et al. 2013) as well as microcracks (Bai et al. 2013) without necessarily quantifying the degree of heterogeneity. Extensive research into the effects of particle shape and polydispersity on flow through porous media using synthetic models has shed more light on how pore geometry and complexity of different porous media affect porosity and permeability (Mota et al. 2005; Garcia et al. 2009). However, little has been performed to fully characterise same in unconventional shale geometries.

Microcomputed tomography (micro-CT) is one of the most powerful non-destructive imaging techniques that can be used to characterise the microstructure of porous materials at the microscopic scale (Rozenbaum and du Roscoat 2014). Over the years, this technique has contributed immensely to the understanding and characterisation of fluid flow through pore geometries of shale rocks (Watson and Mudra 1994; Boruah and Ganapathi 2015). Ma et al. (2016) recently reported that they were able to scan sample size of 100 µm and obtained spatial pixel size of 130 nm. This, however, brings about the question of the ability of the acquired images to adequately represent the microstructural characteristics of the data sets. In order to simulate fluid flow processes, an average of seven grains diameter is required for a fully developed representative elemental volume (REV) (Garcia et al. 2009).

Replacing a heterogeneous property of the rock with an equivalent homogeneous one through a continuum description could best be done through REV method, which represents appreciable property to capture the heterogeneity (Bear 1972). This analysis is carried out by plotting different sample volumes with their corresponding measured property. As shown in Fig. 1, the measured property (n) varies intensely with small changes in the sample length (L) and begins to reduce until measured property is relatively insensitive to small changes in volume or location.

**Fig. 1** Graphical representation of how representative elementary volume (REV) is determined for a specific property (Bear 1972).

Then, representative amount of porosity can be confidently determined (Bear 1972). According to Sahimi (1995), since the continuum approach breaks down if the correlation lengths in the system approach the size of the system, it is important to evaluate the effect of various length scales at which the system may be considered homogeneous.

In this paper, we developed a workflow for the characterisation of heterogeneity of unconventional shale rock samples. The degree of heterogeneity was evaluated using 1 microcomputed tomography (micro-CT) image from Eagle Ford parallel. The sample was divided into 96 for descriptive statistical analysis. Meanwhile, REV for porosity was evaluated using a set of two microcomputed tomography (micro-CT) images from Eagle Ford parallel and perpendicular. Finally, corner point reference techniques where each of the sample was cropped into 7 were used to achieve REV for porosity. The images were processed, digitised, thresholded, segmented and each feature captured using marching cube algorithm. This allows the segmentation of each sample into four distinct entities, consisting of pores, organic matter, shale grains and minerals.

## Image acquisition and processing

The 3D image of the samples used for the analysis was obtained using an industrial micro-CT device phoenix v|tome|x s. It has a 180-kV, micro- and nano-focus X-ray tube and digital detector array (1000 × 1000 pixels). The flow chart for this study is shown in Fig. 2. Full sample (Fig. 2a) was placed in the sample holder, and 6.7 mm size was scanned by rotating it 360°, and the signal-to-noise ratio used for the acquisition of these data was 37.8%. The resulting projections were converted into a 3D image stack using the PhoenixImaging resolution 3D reconstruction software developed from filtered backprojection Feldkamp algorithm (Feldkamp et al. 1984). The detailed procedure for the sample imaging can be found in Singhal et al. (2013).

**Fig. 2** a Eagle Ford (EGF) shale rock sample obtained parallel to the bedding plane. b Micro-CT image of the shale sample at resolution of 6.7 μm showing 3D volume raw data with artefacts. c The 3D volume after artefacts have been removed. d Graphical annotation design of six (6) REV image windows of 50, 100, 200, 300, 400, 500 and 600 μm pixel length. e Graphical annotation design of ninety-six (96) 1 mm × 1 mm × 1 mm subsamples for REV and heterogeneity analyses, respectively. f Sample description used for heterogeneity analysis. Note: The length of each subsample is denoted by lx where x∈{1,2,3,4,⋯,n} and total length L=∑^{n}_{x=1}lx.

In order to accurately estimate petrophysical properties of the samples, proper segmentation of solid(s) and void phases takes a priority. Segmentation involves the segregation of the grey-level voxels of the 3D image into distinct phases. The presence of artefacts such as streaks, brightness, non-uniformity, and phase-contrast fringes at edges and/or noise would reduce the accuracy of segmentation of these 3D images (Ketcham and Carlson 2001) and subsequently lead to misidentification of shale components (Fig. 2b). Due to complex intrinsic shale properties, advanced image processing such as artefact removal and multiband thresholding is necessary (Iassonov et al. 2009). In addition to acquiring high-quality images, the images were cropped to remove edges artefacts present (Fig. 2b). This will also reduce the processing power requirements of the computing machine.

## Polygon surface representation: marching cubes algorithm

The marching cube algorithm is applied to reconstruct a surface CT scan volumetric data sets. It generates inter-slice connectivity, through linear interpolation of the scan line to calculate triangle vertices. The algorithm works in two steps: initially the surface that matches to the user-specified value is located, and subsequently, the normal to the surface of each vertex of each triangle is calculated.

The algorithm basically operates by determining how the surface intersects one cube and then proceeds (marches) to the subsequent cube. Full details of the implementation of the marching cube algorithm can be found in Lorensen and Cline (1987). In summary, the marching cubes algorithm for creating a polygonal surface representation of an isosurface is given as follows:

- Four slices are read into the memory;
- Two slices are scanned to form a cube from neighbouring slices, and four others from the next slice;
- By contrasting between the density values of the eight vertices while keeping the surface constant, the cube index can be calculated;
- From the index value obtained, using a pre-calculated table, a list of edges is determined;
- Applying linear interpolation to find the surface edge intersection using the densities at each vertex edge;
- The unit normal at each cube vertex is calculated using central difference approach and interpolating the normal to each triangle vertex by using the equations:

where D(i,j,k) is the density at pixel (i, j) in slice k and Δx,Δy,Δz are the lengths of the cube edges. When the gradient is divided by its length, it produces the unit normal at the vertex required for rendering.

7. Output the vertex normal and triangle vertices.

### Representative elementary volume (REV) approach

In pore-scale modelling, smaller subsamples can be treated as REV in order to be used in the evaluation of Darcy fluxes. REV is assumed to be obtained when the computed variable (in this case porosity) plotted against increasing sample size (Fig. 2d) does not change significantly at a plateau value (e.g. Figure 1). This is particularly useful when determining what sample to be used in modelling and simulating fluid flow through them.

In order to evaluate deterministic REV, each sample was further incrementally cropped into seven subsamples using corner point reference technique (see Fig. 2d). From the origin (at x, y, z = 0), we extract seven (7) three-dimensional subsample image windows of 50, 100, 200, 300, 400, 500 and 600 μm pixel lengths. The length of each subsample is denoted by l_{x} where x∈{1,2,3,4,…,n} and totallength,L=∑^{n}_{x=1}l_{x} (Fig. 2d). Similar REV data to the conceptual diagram shown in Fig. 1 were obtained for each sample, for porosity variable. REV was determined using the approaches of Costanza-Robinson et al. (2011) and was taken as the minimum window length scale (l) at which the absolute value of the relative gradient error (ε_{g}) in the porosity (∅) remained below 0.003:

where l is the window increment identity/number and δ is the magnitude of l_{x}. In physical terms, the relative gradient error REV criterion requires changes in the measured variable over a given length-scale increment to be relatively small proportional to the increment size.

The choice of the ε criterion (e.g. <0.002) depends on region II analysis given by (Bear 1972). Costanza-Robinson (2011) states that the semi-quantitative ε approach used here is advantageous because it makes REV estimation to be automated and reproducible across numerous images. This allows quantitative relationships between REV and other variables to be readily evaluated; however, it is not applicable to shale rocks which are highly heterogeneous.

### Image processing and heterogeneity analysis technique

In order to investigate the degree of sample heterogeneity, each image sample was subdivided into 1 mm × 1 mm × 1 mm as shown in Fig. 2e. Hence, x, y, z represent row (R), layer (L), and column number (N), respectively. Pore volume was measured for each of the subcropped samples after processing (Fig. 3).

**Fig. 3** Illustration of calculation of absolute value of the relative gradient error.

Each of the subcropped images was morphologically processed. A morphological filtering algorithm opens up holes and gaps in a mask to get rid of small and potentially spurious features. This is to enhance component segmentation (foreground from background as the case may be, Fig. 4a, b).

**Fig. 4** Micro-CT image of Eagle Ford shale rock sample obtained parallel to the bedding plane at resolution of 6.7 µm and physical extent is 1 mm × 1 mm. a 2D view of typical subsample before morphological processing. b 2D view of typical subsample after morphological processing.