Linking multi-scale 3D microstructure to potential enhanced natural gas recovery and subsurface CO₂ storage for Bowland shale, UK

Abstract

Injection of CO₂ into shale reservoirs to enhance gas recovery and simultaneously sequester greenhouse gases is a potential contributor towards the carbon-neutral target. It offers a low-carbon, low-cost, low-waste and large-scale solution during the energy transition period. A precondition to efficient gas storage and flow is a sound understanding of how the shale's micro-scale impacts on these phenomena. However, the heterogeneous and complex nature of shales limits the understanding of microstructure and pore systems, making feasibility analysis challenging. This study qualitatively and quantitatively investigates the Bowland shale microstructure in 3D at five length scales: artificial fractures at 10–100 μm scale, matrix fabric at 1–10 μm-scale, individual mineral grains and organic matter particles at 100 nm–1 μm scale, macropores and micro-cracks at 10–100 nm scale and organic matter and mineral pores at 1–10 nm-scale. For each feature, the volume fraction variations along the bedding normal orientation, the fractal dimensions and the degrees of anisotropy were analysed at all corresponding scales for a multi-scale heterogeneity analysis. The results are combined with other bulk laboratory measurements, including supercritical CO₂ and CH₄ adsorption at reservoir conditions, pressure-dependent permeability and nitrogen adsorption pore size distribution, to perform a comprehensive analysis on the storage space and flow pathways. A cross-scale pore size distribution, ranging from 2 nm to 3 μm, was calculated with quantified microstructure. The cumulative porosity is calculated to be 8%. The cumulative surface area is 17.6 m² g⁻¹. A model of CH₄ and CO₂ flow pathways and storage with quantified microstructure is presented and discussed. The feasibility of simultaneously enhanced gas recovery and subsurface CO₂ storage in Bowland shale, the largest shale gas potential formation in the UK, was assessed based using multi-scale microstructure analysis. The potential is estimated to store 19.0–21.2 Gt CO₂ as free molecules, together with 18.3–28.5 Gt CO₂ adsorbed onto pore surfaces, implying a theoretical maximum of 47.5–49.5 Gt carbon storage in the current estimate of 38 trillion cubic metres (∼1300 trillion cubic feet) of Bowland shale. Simple estimates suggest 6.0–15.8 Gt CO₂ may be stored in practice.

---

Broader context

During the energy transition from fossil fuel to green energy, relatively clean gas extraction and subsurface carbon sequestration are two of the major pathways towards the carbon-neutral target. Injection of CO₂ into shale reservoirs has the potential to integrate both pathways, to enhance gas recovery and sequester greenhouse gases simultaneously. It offers a low-carbon, low-cost, low-waste and large-scale solution. Shales have heterogeneous and complex microstructure, making the feasibility investigation and scientific understanding of CO₂ injection in shales highly challenging. This paper imaged the key features, including artificial fractures, matrix fabric, individual mineral grains, organic matter particles, and pores, at five different scales in 3D and quantitatively analysed their contribution in the natural gas extraction and CO₂ injection pathways. A model of CH₄ and CO₂ flow pathways and storage in shales across multiple scales is built based on the enhanced understanding of 3D quantified microstructure. Meanwhile, the feasibility of simultaneously enhanced gas recovery and subsurface CO₂ storage in Bowland shale, the largest shale gas potential formation in the UK, was assessed. Both theoretical maximum storage and the technical and economical storage were estimated, providing crucial guidance for energy development in energy transition periods.

1 Introduction

Climate change has resulted in many countries aiming to reach net-zero greenhouse gas emissions by 2050. In the energy transition from fossil-based fuel to zero-carbon, several strategies need to be considered. Natural gas produced from shales is considered to be a relatively ‘clean’ energy source, which could increase the energy supply and reduce CO₂ emissions significantly, when compared with oil.¹ Additionally, it is generally accepted that to meet net-zero, carbon capture, utilisation and storage is needed.

The concept of CO₂ injection into shale reservoirs has the potential to achieve these two goals of enhanced natural gas extraction and carbon sequestration simultaneously e.g.ref. 2 and 3. It could potentially protect fresh water, enhance gas recovery, and sequester CO₂ into the subsurface.³ The process of preferentially adsorbing CO₂ while desorbing CH₄ has been suggested in shale reservoirs previously.^4,5 Globally there has been interest from government, industry and the academic community, but high costs and outstanding research questions restrict wide-scale implementation of this technology e.g.6. A core issue for commercialisation is the limited understanding of pore networks, gas transport pathways and gas behaviours in shales.⁷ The heterogeneous nature, fine-grain size and small pores within shales makes characterization of microstructure, pore networks and gas transport challenging.^{ref. 8} Shale microstructure classically consists of pores and fractures distributed within a matrix in which non porous and porous grains are embedded at different scales. All can contribute differently to transport properties and storage properties of CO₂ and CH₄. Furthermore, recent experimental results by Song et al., (2019) suggest that using CO₂ as a fracturing fluid rather than water may have a significant impact driven by a lower fracturing pressure, while leading to a denser network of largely shear-mode fractures.

A range of techniques have been applied to quantify microstructure and pore networks in shales, for example, gas adsorption, mercury intrusion, small angle neutron or X-ray scattering e.g.ref. 9–11. They can provide data on pore size, porosity, adsorption, and permeability, but the pore structure, spatial distributions and connectivity cannot be observed or quantified directly. The complexity of shale makes it necessary to additionally quantify the microstructure and their pore system in 3D by means of advanced imaging techniques. High-resolution SEM imaging techniques, such as broad ion beam-scanning electron microscopy (BIB-SEM), can characterise the microstructure over a range of scales in 2D.¹² X-Ray tomography has been applied to fractures,¹³ organic matter,¹⁴ large mineral grains⁸ and clay minerals¹⁵ at micro-scale in 3D. Ga⁺ focused ion beam-scanning electron microscopy (FIB-SEM)¹⁶ and Xe⁺ plasma FIB (PFIB, or P-FIB)¹⁷ can be applied to resolve macro- and meso-pores at nm-scale.¹⁶ Nanopores within organic matter or between clay minerals can be imaged by Transmission Electronic Microscopy (TEM) tomography¹⁷ and Helium Ion Microscopy (HIM).¹⁸ Owing to the complex microstructure and high heterogeneous nature of shales, the assessment of potential of CO₂ storage and enhanced gas recovery in shales requires a multi-scale, multi-approach investigation of the open pore spaces for free molecules, the pore surfaces for adsorbed molecules and the flow pathways between them. Particularly, a detailed quantification is needed of the nm-scale pores and the μm-scale phases associated with these pores (e.g. minerals and organic matter) alongside with the mm-scale sample fabric and stimulated fractures. To this end, in this study we combine imaging across multiple scales with other laboratory measurements, for example, porosity, permeability and adsorption, to characterize shale microstructure and pore system in 3D and to link these with the potential for enhanced gas recovery and subsurface CO₂ storage.

As the largest shale gas resource in the UK, the Bowland shale is estimated to contain a gas resource of 1300 trillion cubic feet (at surface temperature and pressure).¹⁹ The gas resource is comparable to those of major producing shales in North America, such as Barnett, Woodford and Fayetteville Shales.²⁰ However, previous studies mostly focused on the shale gas potential, while the carbon sequestration potential has not been considered. The aims of the research presented here are: (1) to characterise the shale microstructure across a range of scales (nm- to mm-) and to assess its heterogeneity; (2) to quantify the pore structure and network to investigate the possible gas storage and pathways for methane and CO₂ flow; and (3) to link the microstructure and reactive transport to further assess the potential for simultaneous enhanced gas recovery and subsurface CO₂ storage.

2 Materials and methods

2.1 Sample selection

Recent studies of the Carboniferous Bowland Shale Formation (herein referred to as the Bowland Shale) e.g.ref. 20 and 21, have demonstrated that the formation contains organic-rich horizons. Samples extracted from these horizons contain greater than 1% TOC, are quartz-rich (>50 wt%) and clay-poor (<20 wt%). Based on detailed core observation and microfacies analysis, three typical samples were selected in this study. These consist of (1) one typical organic-rich Bowland Shale (BOR; TOC 6.1 wt%) sample, refer to B6 in ref. 22 from the upper Bowland Shale Formation, which was selected for cross-scale 3D microstructure characterisation and adsorption experiments (Table 1). This sample was taken from a depth of 2344.7 m in the Preese Hall-1 borehole which was the UK's first shale gas exploration well (Clarke et al., 2018). (2) One typical organic-lean (BOL; TOC 1.0 wt%) sample refer to B8 in ref. 22 at a depth of 2495.3 in Preese Hall-1, which was used in adsorption experiments for comparison. (3) An adjacent sample (BF) to BOR with similar compositions in the same core (2354.5 m) was selected for fracturing characterisation owing to the limited volume of available BOR material.

Table 1 Samples compositions based on TOC and XRD analysis (wt%)

Sample name	TOC	Quartz	Calcite	Ankerite	Pyrite	Muscovite	Kaolinite	Albite
BOR	6	67	6	2	1	9	7	2
BOL	1	51	3	11	2	5	18	9
BF	1	69	7	5	1	6	8	3

---

2.2 Imaging techniques

The selected typical organic-rich sample was firstly examined in 2D for microstructure and pore network characterization and then in 3D for cross-scale quantification.

2.2.1 2D imaging (SEM and HIM). High resolution scanning electron microscopy (SEM) and Helion Ion Microscopy (HIM) were used for multi-scale 2D characterization of the microstructure and pore system. HIM observations were made at the Molecular Foundry, Lawrence Berkeley National Laboratory, USA. The pixel size for SEM ranges from 1 nm to 500 nm, and that for HIM image is 0.2 nm.

2.2.2 X-Ray computed tomography (XCT). To characterise the heterogeneous microstructure in 3D, multi-scale X-ray computed tomography (XCT) was used in this study.
(1) Synchrotron source macro-CT. BF sample was imaged in I12 Beamline (EE17606) at Diamond Light Source (DLS, UK) to characterize the fractures created by high-pressure fluids during injection experiments and under confinement [details of the experiment can be found in ref. 23]. At 51 MPa confining pressure, the fracture opens at 53 MPa injection pressure.
(2) Synchrotron source micro-XCT. BOR sample was cut into a 2 × 2 × 10 mm³ stick and scanned in I13 Beamline (MT15506-1) at DLS using PCO Edge 4× camera and Pink beam. The energy was 20 keV and exposure time is 0.04 s. 4000 projections were taken for each scan with voxel size 1.6 μm. Four tomograms were taken at different vertical locations with 10% overlap to access the larger range of microstructure variation.
(3) Lab source nano-XCT. BOR sample was cut into a 60 μm diameter cylinder using a laser cutter and scanned in a Zeiss Ultra 810 nano-CT scanner in Manchester X-ray Imaging Facility (MXIF). Two scans using phase contrast mode and adsorption mode were performed for the same volume. These two images were registered together for segmentation. 8 kV monochromatic beam was used for the scan with 160 s exposition time. The voxel sizes of both images were 64 nm.

2.2.3 FIB-SEM. To characterize the nanoscale features in shales, advanced Focused Ion Beam Scanning Electron Microscopy (FIB-SEM) is used to perform high-magnification 3D imaging. Specifically, Xe⁺ plasma FIB (PFIB) and Ga⁺ FIB are both used in this study at diffferent scales.¹⁷
(1) Plasma FIB-SEM. The BOR sample scanned in nano-CT was used for plasma FIB. Pore and organic matter cannot be distinguished in the nano-CT images, but are identifiable in the PFIB data. Secondary electron images were taken at 1 kV and the organic matter and pores were quantified after image processing. The SEM images were taken at 3449 × 3735 voxels with 33 × 42 nm² resolution. 815 slices were collected during rock milling with a milling distance 50 nm.
(2) Ga⁺ FIB-SEM. Two sites were selected to be imaged with Ga⁺ FIB-SEM, including one for mineral matrix and one for a single organic matter particle. The organic matter particle site was slightly smaller than the mineral matrix site (Table 2), and both of the data set were collected at 2 kV.

Table 2 Image techniques used for 3D microstructure of BOR in this study

Imaging technique	Synchrotron macro-CT	Synchrotron micro-CT	Nano-CT	Plasma-FIB	Ga FIB
Corresponding scale	10–100 μm	1–10 μm	100 nm–1 μm	10–100 nm	1–10 nm
Target feature	Stimulated fractures	Non-fracture texture	Fine-grain matrix	Macropores	Pores associated with organic matter particles	Pores associated with mineral matrix
Sample volume	∼20 mm diameter and 20 mm length cylinder	∼2 × 2 × 10 mm^3	∼60 μm diameter and 60 mm length cylinder	∼60 × 60 × 40 μm^3	∼10 × 8 × 4 μm^3	∼10 × 10 × 7 μm^3
Voxel size	16 × 16 × 16 μm^3	1.6 × 1.6 × 1.6 μm^3	64× 64× 64 nm^3	33 × 42 × 50 nm^3	10 × 10 × 18 nm^3	7 × 7 × 14 nm^3
Facility modal	I12 beamline	I13–2 beamline	Zeiss ultra 810	FEI Helios PFIB	FEI Nova 600	FEI Nova 600
Facility location	Diamond Light Source	Diamond Light Source	MXIF, UoM	Henry Royce Institute, UoM	Electron Microscopy Centre, UoM	Electron Microscopy Centre, UoM

---

2.2.4 Image processing. X-Ray tomography images required specific image processing to align the images into the same locations. Four micro-CT tomograms were registered into a whole image using the overlap volumes of the top and bottom of adjacent tomograms in Avizo (Thermo Fisher Scientific, United States). Two nano-CT tomograms of adsorption and phase-contrast scans were registered into the same location for an enhanced quantification of minerals and organic matter. Specifically, granular minerals and clay minerals were segmented in the phase-contrast image while the heavy minerals and organic matter were segmented at the adsorption image. PFIB and both Ga⁺ images were aligned (rotation) and sheared to correct the drift during the imaging.

Subsequently, all images were filtered and segmented using varied filtering (e.g. non-local means, FFT, median filters) and segmentation methods (e.g. thresholding, tophat, watershed) in different series of images to present the identified features at the corresponding scale, as a standard workflow.²⁴ Features below 27 voxels (3 × 3 × 3 for XYZ) for XCT data and 9 voxels (3 × 3 for XY plane) for FIB data were not considered. The size, volume, surface area, volume fractures of each feature were measured in voxels using Avizo (Standard and Fire versions, FEI, Hillsboro, United States). The sizes of pores or particles are measured in sphere equivalent diameter in this paper. This is defined as the diameter of a sphere with the equivalent volume.

2.2.5 Heterogeneity analysis. The heterogeneity of the microstructure and pore system in the 3D images was quantified with respect to vertical variation, fractal dimensions and degree of anisotropy.
(1) Vertical variation. Horizontal bedding was recognized in the Bowland shale samples (e.g. 1–2 mm), and therefore vertical variations of different phases were measured for heterogeneity analysis of microstructure. In each image dataset, the volume fractions of identified phases at this scale were quantified in all horizontal slices and were then plotted with the corresponding vertical positions. The maximum, minimum, median values and standard deviations for each phase at each scale were also identified. The variations of phase distribution at different scales were then quantified.
(2) Fractal dimensions. The fractal dimension was determined by a box-counting algorithm to analyse the spatial arrangement of porous phases or pores within the 3D binary images using ‘BoneJ’ plug-in in ImageJ.²⁵ Certain sizes of boxes with the identified feature scan the whole image and the number of the boxes were counted. The box counting starts from the size equalling to 80% of the side length of the original image and ends at the 27 voxels (3 × 3 × 3). Features smaller than 27 voxels are not considered owing to the small sizes. The scaling factor is set to be 1.2 for each of the image to get enough steps for the analysis, which means the box size was divided by 1.2 after each iteration.
(3) Degree of anisotropy. The degree of anisotropy is calculated to quantify the directionality of the shale components and pores using the ‘BoneJ’ plugin in ImageJ.²⁶ It consists of four steps to determine the degree of anisotropy including (1) find mean intercept length (MIL) of vectors in random directions; (2) plot MIL vectors from background to foreground into a point cloud; (3) solve the equation of an ellipsoid that best fits the point cloud and (4) calculate the degree of anisotropy from the radii of the ellipsoid (longest/shortest). The degree of anisotropy ranges from 0–1 and the closer to 1 means higher anisotropy of the identified feature.

2.3 Laboratory bulk property measurement

Bulk laboratory measurements were conducted to provide the essential information for the assessment of gas recovery and CO₂ sequestration along with the image results. These include (i) high-pressure adsorption measurements (CO₂ and CH₄); (ii) low-pressure N₂ and CO₂ adsorption measurements for pore size distribution analysis; and (iii) porosity and permeability.

(1) High-pressure adsorption. CO₂ and CH₄ adsorption measurements were made BOR and BOL at 7 MPa and 80 °C (reservoir conditions). The measurements have been carried out on powdered samples (2 g) in a Rubotherm Magnetic Suspension Balance following the procedure described in ref. 27. Samples were dried at 120 °C under vacuum for at least 12 hours in situ prior to the measurements. The obtained excess amount adsorbed, n^ex, was converted to absolute amount adsorbed, n^a, using the following equation: n^a = n^ex/(1 − ρ/ρ^a), where ρ is the bulk density of the given gas at the selected pressure and temperature, while ρ^a is the adsorbed-phase density, presumed to be constant at 21.1 mol L⁻¹ for CO₂ and 26.3 mol L⁻¹ for CH₄.

(2) Pore-size distribution. Cryogenic adsorption of N₂ at 77 K and CO₂ at 273 K was conducted on about 1 g of both samples using a Micromeritics 3Flex Surface Characterisation Analyzer. Samples were first degassed at 120 °C under vacuum ex situ for at least 12 hours and for a further 3 hours in situ at the same temperature. Gases were probed in the pressure range 2 × 10⁻⁶–1 × 10⁻¹ MPa. The distribution of pore sizes above 2 nm (mesopores and macropores) was obtained by the application of the BJH theory to the adsorption branch of the N₂ isotherms²⁸ while the micropore (pores less than 2 nm in size) volume was obtained by the application of the linearised Dubinin–Radushkevich (DR) equation to the CO₂ adsorption isotherm.²⁹ The Harkins–Jura equation was used to describe the thickness of the adsorbed layer.³⁰ A liquid density of 21.1 mol L⁻¹ was assumed for CO₂.

(3) Porosity and permeability measurement. Porosity was measured by helium porosimeter and gas permeability was measured by pore pressure oscillation³¹ using Argon gas, and in two orientations, parallel and normal to bedding. Six permeability measurements were recorded in each orientation to determine a mean value and its uncertainty with at pore pressure 5 MPa and confining pressure 10.5 MPa.

2.4 Porosity and surface area upscaling from 1–100 nm to 10–100 μm scale

The porosity, pore size distribution and surface area measured at 1–100 nm from images and lab measurements were upscaled to pore associated phases (organic matter and minerals) at μm-scale; and, then to the sample fabric at 10–100 μm scale using an image-based upscaling method that has been shown to be applicable to shales.²⁴

The correlation coefficients (i.e. volume and surface area ratios) of pore and pore associated phases were quantified at 1–100 nm. Those were then assigned into the correlation coefficients of pore associated phases and sample matrix at μm-scale, and finally to the whole sample at 10–100 μm scale. Subsequently, pore size distribution was normalised as number of pores per μm³, and as surface areas per μm³ for the comparisons between scales.

2.5 CO₂ storage calculation

2.5.1 Derivation of a relationship between pore and reservoir scales. The volumes of storage for CO₂ in shales are calculated through eqn (3) and (4), which are modified based on the US-DOE-NETL method [eqn (1) and (2)];³²

G_CO2 = A_tE_aH_gE_h[ρ_CO2ΦE_Φ + ρ_sCO2(1 − Φ)E_s] (1)

= V_e[ρ_CO2ΦE_Φ + ρ_sCO2(1 − Φ)E_s] (2)

where G_CO2 is the mass of CO₂ that is stored in the available volume; A_t and H_g are the area and thickness of gas shale; ρ_CO2 is the density of liquid or supercritical CO₂ at downhole conditions; Φ is the porosity for free gas; 1 − Φ is the solid fraction; and ρ_sCO2 is the mass of CO₂ adsorbed per unit volume of solid rock (corresponding to 1 − Φ). E_A and E_h are the field-scale efficiency factors (ratios of efficient and total values) for the area and thickness. E_Φ and E_s are pore-scale efficiency factors for pore volume and sorbed volume, which are used to correct values measured in the laboratory, including effects arising from inaccessible porosity, variable mineralogy and the presence of moisture.

Furthermore, E_s can be calculated by the equation below:

E_s = E_s(d) × E_s(m) (3)

where E_s(d) refers to the adsorbed CO₂ at dry conditions and accounts for the variability in sample mineralogy. E_s(m) is the moisture equilibrated adsorption efficiency, and refers to the ratio of the adsorbed CO₂ at moisture-equilibrated and dry conditions. It additionally accounts for a reduction in gas uptake due to the presence of moisture.

Therefore eqn (1) and (2) can be further re-written as

G_CO2 = VE_v[ρ_CO2ΦE_Φ + ρ_sCO2(1 − Φ)E_s(d) × E_s(m)] (4)

where V is the efficient volume for shale gas, E_v is the fraction of shale volume used for CO₂ storage, which also refers to the practical efficiency compared with the theoretical maximum value.

2.5.2 Parameters used in the equations. Chen et al. 2019³³ report that pores below 2 nm aperture are predominantly filled by adsorbed gas and the majority of free gas resides in pores >5 nm aperture (which make up 77.7% of total pore volume). In this study, Φ is therefore computed as the cumulative pore space found in pores larger than 2 nm and 5 nm and these two thresholds are used to estimate to find the range of the stored gas volume. V, ρ_CO2, and E_v were determined as the average values reported in the previous studies of Bowland shale.^19,20

E _v, the practical efficiency, is a reservoir scale efficiency factor not measured in this study and therefore set to be 1 for the theoretical maximum value estimation. Different values of E_v were used for practical estimation in the next step. ρ_sCO2 was estimated from the adsorption experiments by the gas uptake measured on a dry shale. E_s(d) was estimated from measurements carried out on BOR and BOL samples. It is expected that the primary formation water reduces the adsorption capacity, and this effect is included in E_s(m). The efficiency factor E_s was considered at moisture equilibrated status in the subsurface. The adsorption was measured on dry samples in the laboratory for E_s(d). In reality, the exiting of primary formation water in the surface reduces the adsorption ability, which affects E_s(m). The water saturation in Bowland shale is relatively low with an average value 32.4% across all formations.²⁰ Therefore, the final adsorbed values of Bowland shale in the subsurface are assumed based on moisture equilibrated status according to other studies with a similar water content, for example, 22% less CO₂ adsorbed in a 32.10% moisture sample than the dried samples e.g.ref. 34. Therefore E_s(m) in eqn (3) is set to be 0.78 (= 1–0.22) in this study.

Then, two simple calculations were performed for practical volume estimations. The first calculation was based on an estimation that technically accessible CO₂ storage was around 32% in some commercial shale plays³⁵ and therefore E_v was set to be 0.32. The second estimation was performed using P10 (10 percent probability that an estimated value is less than the P10) and P90 (90 percent probability that the estimated value is less than the P90 value) of E_v for free gas (0.15; 0.36) and adsorbed gas (0.11; 0.24), which were extracted from the modelling results of different shale reservoirs³⁶ for practical storage potentials estimation of stored CO₂ volume.

3 Results

Artificial fractures, sample fabric, mineral matrix and organic matter, macro-pores, and nano-pores were imaged and quantified at different scales. Specifically, artificial fractures were quantified at 10–100 μm scale; sample fabric, including the distribution of granular minerals, clay minerals, organic matter particles and heavy minerals, were characterised at low-resolution micro-scale (1–10 μm scale); the grain sizes and higher-resolution distribution of these porous phases were further characterised and quantified at high-resolution micron-scale (100 nm–1 μm scale); macropores (>50 nm) within these phases were quantified at low-resolution nano-scale (10–100 nm scale) and mesopores (2–50 nm) were characterised at high-resolution nano-scale (1–10 nm scale). Additionally, micropores (<2 nm) were measured by sorption experiments.

3.1 Identification of microstructures

The major compositions are identified as granular minerals (GM), clay minerals (CM), organic matter (OM) and heavy minerals (HM). The distribution of these compositions and pores were recognized as sample microstructure. The overall microstructure of the sample is relatively similar horizontally with clusters of large grains and matrix (Fig. 1a). Granular minerals and heavy minerals scatter in the sample, while clay minerals and organic matter distribute partially along the beddings in the matrix (Fig. 1b).

Fig. 1 SEM images of the Bowland shale. (a) Overview of the sample microstructure; (b) higher resolution overview of the sample microstructure; identification of (c) granular minerals and clay minerals; (d) heavy minerals and organic matter particles; (e) minerals pores and micro-cracks; (f) organic matter pores; (g) inter-mineral pores (imaged with HIM); (h) organic matter pores (imaged with HIM). GM – granular mineral, CM – clay minerals, OM – organic matter, HM – heavy minerals. Qz – quartz, Alb – albite, Cal – calcite, Ank – ankerite, Py – pyrite, Mus – muscovite, Kao – kaolinite. HIM – Helion ion microscopy.

3.1.1 Bowland shale components:.
(1) Granular minerals. Quartz, calcite, ankerite and albite are observed as granular minerals in these Bowland Shale samples (Fig. 1a–d). Small pores exist in the crystal gaps, with spherical or irregular polygon shapes, which are unlikely to have an interconnected pore network (Fig. 1e).
(2) Clay minerals. Kaolinite and muscovite dominate in clay minerals in these samples, partly orientated to the bedding plane (Fig. 1a–d). Inter-mineral pores are present between clay mineral grains or between clay mineral and other mineral grains with elongated wedge shapes (Fig. 1e). The latter of the two are generally slightly larger, especially when fine-grain granular minerals inlay the clay mineral matrix.
(3) Organic matter. Organic matter particles are present scattered in the matrix (Fig. 1a–d). Some are porous while others are non-porous. Many of the porous organic matter particles are large and spherical in shape. Some are lamellar, interbedded with chlorite or other clay minerals. Non-porous organic matter is dominated by elongate shape or irregular polygons. Some inter particle pores are present in the interface of non-porous organic matter and adjacent minerals. They are commonly present as curved flakes or possessing a crack-like geometry, aligned along the edge of non-porous organic matter particles (Fig. 1f).
(4) Heavy minerals. The major heavy mineral detected in the Bowland Shale is pyrite, present in the forms of framboids and single crystals (Fig. 1a–d). Framboids host some polygon shaped pores while single crystals do not host pores.

3.1.2 Microcracks and pores. Microcracks and pores were observed in and between minerals and organic matter (Fig. 1e and f). Microcracks in the matrix follows the mineral boundaries primarily, while those in organic matter have different distributions. Some of cracks present at the interface of organic matter and clay minerals, and some cracks open in the middle or organic matter particles. Pores are observed as mineral pores and organic matter pores. The minerals pores are generally irregular in shape with wedge-shape corners while the organic matter pores are predominately spherical. Pores and microcracks are not distinguished in this study owing to the ambiguous definitions and comparable roles they are playing for carbon storage and sequestration.

3.2 Multi-scale 3D quantifications

The spatial distribution of the four major components was quantified at 1–10 μm scale (1.6 μm voxel size) and the individual grains of these components were characterized at 100 nm–1 μm scale (0.1 μm voxel size). Macropores were imaged at 10–100 nm scale, and then mesopores in minerals and organic matter then quantified respectively at 1–10 nm scale. To assess the heterogeneity, the volume fractions of each phase were plotted with the vertical distance from the top of the images at LR-μm, HR-μm and 10–100 nm scales. This was not present either for the 10–100 μm scale as the matrix is not able to be identified or for 1–10 nm scale owing to the small sample volume below the representativity elementary volume (REV).

3.2.1 Fractures at 10–100 μm scale. Stimulated fractures opened at 53 MPa injection pressure from a borehole oriented normal to bedding, and a confining pressure of 51 MPa. Two vertical (cross-bedding; along the injection direction) fractures and one horizontal (along bedding; vertical to the injection direction) were identified in the image (Fig. 2). The fractures are relatively planar across the sample section. This image was captured immediately after the fracture opening and the width and volume of the fractures may vary over time, therefore the quantification of fracture apertures and volumes is not presented here.

Fig. 2 Artificial fractures within BF sample at 53 MPa and injection pressure and 51 MPa confining pressure.

3.2.2 Sample fabric at 1–10 μm scale. The four phases in the selected volume of the sample, (organic matter, granular minerals, clay minerals, and heavy minerals), were imaged using synchrotron-based X-ray tomography at 1–10 μm scale (Fig. 3). Micro-cracks or pores are not distinguishable from the above four phases at this scale. The quantification uncertainty of volume percentages is ±4 vol% at this scale granular minerals occupy 71 vol% while clay minerals account for 26 vol%. Organic matter and heavy minerals only contribute 2 and 2 vol% in volume fractions respectively. Granular minerals show a decreasing trend upwards which ranges from 85% to 48% with a median value 72%, while the clay mineral increase from 12% to 45% with a median value 25%. Organic matter slightly increases from 1% to 3%, correspondingly with clay minerals, with a median value 14%. Heavy minerals distribute relatively evenly, and the volume fraction ranges from 1% to 3%, with a median value 2%. In this study, note here that the organic matter content is able to be higher than 1% in local intervals of 500 μm to 1.5 mm length.

Fig. 3 3D microstructure of BOR sample (3 × 1 × 1 mm³) at 1–10 μm scale. (a) Slice view, (b) organic matter, (c) granular minerals, (d) clay minerals, (e) heavy minerals, (f) mixed phases, (g) volume fractions of different compositions, (h) X axis is the zoom-in view of (g) for the organic matter and heavy minerals curves.

3.2.3 Porous phases at 100 nm–1 μm scale. The four major components, organic matter, granular minerals, clay minerals, and heavy minerals, were imaged at high-resolution micro-scale using nano-CT (Fig. 4) for higher resolution quantification at 100 nm–1 μm scale. The quantification uncertainty of volume percentages is ±2 vol% at this scale. Granular minerals, clay minerals, organic matter and heavy minerals occupy 71, 17, 10 and 2 vol% respectively. Pores and micro-cracks are not distinguished at this scale, but are separated in the same sample at the next scale.

Fig. 4 3D microstructure of BOR sample at 100 nm–1 μm scale. (a) Slice view, (b) granular minerals, (c) organic matter, (d) clay minerals, (e) heavy minerals, (f) mixed phases, (g) volume fractions of different compositions.

The vertical variation for granular minerals ranges from 46% to 52% with a median value 46% and for clay minerals ranges from 1% to 2% with a median value 1%. Organic matter varies from 43% to 8% with a median value 7%. Heavy minerals max at 15% and min at 0% with a median 4%.

3.2.4 Pore size, structure, and networks at nm-scale.
3.2.4.1 Macro-pores at 10–100 nm scale. The sample scanned in PFIB was part of the one scanned in nano-CT, and macropores (>100 nm) are distinguishable from the matrix (Fig. 6a), which occupies 0.64% (Fig. 5b). The quantification uncertainty of volume percentages is ±0.5 vol% at this scale. 16 056 matter particles and 5832 pores were identified in the PFIB images at 10–100 nm scale. The vertical variation of the porosity ranges from 0.0 to 3.3% at this scale with the step of Z axis 64 nm (Fig. 5c). A good correlation between organic matter and pore volume fraction was found with a correlation coefficient³⁷ value 0.84. The largest pore sizes were quantified to be 2.8 μm with bimodal distribution peaking at 800 nm and 2 μm (Fig. 8b).

Fig. 5 3D microstructure of BOR sample at 10–100 nm scale. (a) Vertical slice view, (b) organic matter and pores distribution, blue-organic matter, red-pores, (c) volume fractions of pores by horizontal slice.

3.2.4.2 Meso-pores at 1–10 nm scale. Based on the observations under high-resolution SEM, organic matter and minerals are two typical phases hosting pores (Fig. 1). Imaging and quantification of the pore system was performed in an organic matter dominated site and a mineral dominated site (Fig. 6a and d). Due to the image resolution limitation (equivalent pixel size ∼12 nm), only pores above 30 nm aperture are quantified and described in this section. 5574 mineral pores and 3353 organic matter pores were identified in the two FIB images respectively at HR-nanoscale. The quantification uncertainty of volume percentages is ±0.3 vol% at this scale. The normalized cumulative porosity and surface area of the whole range of pores is shown in Fig. 8.

Fig. 6 Mesopores observation and quantification in FIB images of BOR sample. (a) Slice view of mineral dominated site, (b) segmented mineral pores, (c) pore network model of mineral pores, (d) slice view of organic matter dominated site, (e) segmented organic matter pores, (f) pore network model of organic matter pores, (g) number of organic matter pores and mineral pores per μm³, (h) incremental porosity of organic matter and minerals by size, (i) cumulative porosity of organic matter and pores by size, (j) incremental and cumulative surface area of mineral pores and organic matter pores. Colors in (b and e) refer to the unconnected pores, and colors in (c and f) refer to the size.

(1) Mineral pores. Pores inside minerals (intra-mineral pores) or between minerals (inter-mineral pores) are both observed in these images (Fig. 6a). Intra-mineral pores are commonly found in granular minerals such as quartz, calcite, ankerite and albite, with wedge shapes. Small pores exist in the crystal gaps, with polygon shapes. They only occupy 8.5% of the total mineral pores. Inter-mineral pores are observed between many mineral types, and account for 91.5% of total mineral pores. Pores between granular minerals display irregular shapes, sometimes featuring sharp corners, while pores between clay mineral grains or between clay minerals and other mineral grains with wedge shapes (Fig. 6b).

The pore network shows a low connectivity at the resolution of the study and a relatively homogenous texture (Fig. 6c). The mineral pores occupy 1.4% of the volume within the mineral phase, while accounting for 65.8% of total pores.

The mineral pore range from 30 nm to 900 nm (equivalent diameter). The major peak in number of pores appears at 30 nm (Fig. 6g), and in terms of incremental porosity, pores between 200–700 nm equivalent diameter make the majority contribution (Fig. 6h). Cumulative mineral porosity from FIB images is 2.3% (Fig. 6i) and cumulative surface area is 520 μm² μm⁻³ (≈1.325 × 10³ μm² g⁻¹) (Fig. 6j).

(2) Organic matter pores. Organic matter pores are mainly spherical, with variable sizes (Fig. 6d and e). The distribution of OM pores is more heterogeneous compared with the mineral pores, and the size ranges from 30 nm to 3 μm in this image, with a major peak at 40 nm and 200 nm. They form a high-density map of pores although not a globally connected network was observed in the FIB image (Fig. 6f). Organic matter pores occupy 9.0% within the organic matter phase.

Although there is a higher density and larger total volume of organic matter pores than mineral pores in the single image dataset, they tend to be lower in total number of pores per μm³ when the whole sample is taken into account (Fig. 6g) owing to the lower volume fraction of organic matter overall. Number of pores show a bimodal distribution, with peaks at 30 nm and 200 nm (Fig. 6g). The incremental porosity peaks at 2000 nm, which is greater in sizes compared with the mineral pores (Fig. 6h). Cumulative organic matter pores in FIB image is 1.2% (Fig. 6i) and the cumulative surface area is 21 μm⁻¹ (≈54.6 × 10³ μm² g⁻¹) (Fig. 6j).

3.2.5 Heterogeneity across scales. The fractal dimensions and the degree of anisotropy of each feature were obtained at multiple scales for heterogeneity studies. Generally, fractal dimensions decrease with higher resolution and the degree of anisotropy increases for each phase due to the sample's limited size (Table 3). At 1–10 μm scale, the fractal dimensions of four components range in a small span of only less than 0.1, from 2.546 (OM) to 2.642 (GM). The degrees of anisotropy for GM and CM are relatively higher with values of 0.601 and 0.635 and those for OM and HM are only 0.435 and 0.272 respectively. When resolution improves to 100 nm–1 μm scale, the individual grains become visible. Compared with 1–10 μm scale, fractal dimensions at this scale decrease to the range of 2.350–2.523. The degrees of anisotropy for HM, GM and OM decrease to 0.171, 0.181 and 0.342 respectively while that for CM increase to 0.758.

Table 3 Statistical analysis of the vertical variation of the volume percentages of each imaged phase of BOR sample at all scales (GM – granular mineral, CM – clay mineral, HM – heavy mineral, and OM – organic matter)

Scale	Parameter of volume percentages	GM	CM	HM	OM	Macropores	Mesopores (OM)	Mesopores (mineral)
1–10 μm	Max	0.853	0.460	0.030	0.036
	Min	0.488	0.121	0.011	0.006
	Median	0.723	0.246	0.018	0.014
	Standard deviation	0.096	0.089	0.003	0.008
	Fractal dimension (FD)	2.642	2.622	2.558	2.546
	R ^2 for FD	0.989	0.991	0.990	0.990
	Degree of anisotropy	0.601	0.635	0.272	0.435
100 nm–1 μm	Max	0.518	0.020	0.149	0.535
	Min	0.074	0.001	0.000	0.044
	Median	0.460	0.016	0.048	0.067
	Standard deviation	0.043	0.003	0.042	0.035
	Fractal dimension (FD)	2.523	2.350	2.409	2.510
	R ^2 for FD	0.993	0.991	0.992	0.994
	Degree of anisotropy	0.181	0.758	0.171	0.342
10–100 nm	Max				0.039	0.033
	Min				0.000	0.000
	Median				0.004	0.004
	Standard deviation				0.008	0.006
	Fractal dimension (FD)				2.458	2.709
	R ^2 for FD					0.997
	Degree of anisotropy					0.429
1–10 nm	Max						0.045	0.043
	Min						0.079	0.000
	Median						0.119	0.024
	Standard deviation						0.022	0.008
	Fractal dimension (FD)						2.316	2.195
	R ^2 for FD						0.996	0.992
	Degree of anisotropy						0.736	0.504

---

Organic matter is the only feature which can be characterized crossing three scales owing to the relatively large particles sizes, the fractal dimensions of which slightly decrease from 2.546 at 1–10 μm scale to 2.510 at 100 nm–1 μm scale and then to 2.458 at 10–100 nm scale. The degree of anisotropy of organic matter range between 0.435 to 0.342 and the standard deviations are 0.008 and 0.035. Fractal dimension of total pores is 2.709 and it is slightly higher for OM pores (2.315) than mineral pores (2.195). The degree of anisotropy increases from 0.429 at 10–100 nm scale to 0.736 for OM pores and 0.504 for mineral pores.

As a consequence, particle and pore anisotropy (a) significantly and quantitatively changes across the scales from the μm to nm, and (b) it is heterogeneous for the different phases at the same scale.

3.3 Laboratory bulk property measurements

3.3.1 CO₂ and CH₄ adsorption. CO₂ and CH₄ adsorption of BOR and BOL was measured at reservoir conditions (80 °C and 7 MPa). The CO₂ is present as supercritical phase (^scCO₂) under this condition. The amount of ^scCO₂ adsorbed in BOR sample was 4.5 kg per t while the CH₄ is only 0.6 kg per t. In contrast, BOL was found to only adsorbed 2.9 kg per t ^scCO₂ and 0.3 kg per t CH₄ (Fig. 7A).

Fig. 7 Adsorption measurements of BOR and BOL samples (A) and permeability measurements of BOR sample (B; note the log scale in Y axis).

3.3.2 Permeability. The permeability result shows two orders magnitude heterogeneity between flow orientations. The mean value is 3.9 × 10⁻²¹ m² in the direction of normal to beddings and 3.3 × 10⁻¹⁹ m² in the direction of parallel to beddings (Fig. 7B).

3.3.3 Pore size distribution (PSD) extracted from adsorption data. The PSD measured by adsorption shows that the major pores distribute between 2–5 nm with a small hump around 25 nm, which corresponds to the peaks shown in the image results (Fig. 8a). When this is converted into incremental pore volume, the peaks are found at 2–3 nm and 10–50 nm (Fig. 8b).

Fig. 8 Pore size distribution by numbers of pores (a), incremental porosity (b), incremental surface area (c) and cumulative porosity (d).

3.4 Upscaled pore size, volume and surface area

The porosity and surface area measured in the images at 1–100 nm were upscaled to 10–100 μm scale using an image based upscaling method.²⁴ The number of pores shows an increasing trend when pore sizes decrease, and the peak appears at 2 nm which is the limit of the techniques used in this study (Fig. 8a). The incremental porosity has the major peak at 200–500 nm and the second peak at 20–40 nm (Fig. 8b). The incremental surface area shows a continuously increasing trend with decreased pore size (Fig. 8c). The cumulative porosity is calculated to be ∼7.8% and the cumulative surface area is 1.79 m² g⁻¹ (Fig. 8d). Among them, pore above 100 nm contribute 47.5% in total porosity (3.7% of sample volume) and less than 0.1% for cumulative surface area (1.94 × 10⁻⁶ m² g⁻¹); pores between 10–100 nm contribute 35.9% in total porosity (2.8% of the sample volume) and less than 0.1% for cumulative surface area for cumulative surface area (7.67 × 10⁻⁴ m² g⁻¹); pore below 10 nm contribute 16.6% in total porosity (1.3% of the sample volume) and over 99.9% for cumulative surface area (1.79 m² g⁻¹) (Fig. 8d).

3.5 CO₂ storage calculation

The potential CO₂ storage volumes in shales are calculated through eqn (4), first for theoretical maximum values, and then the practical values using P10 and P90 values for free gas and adsorbed gas simulated in other studies.³⁶

3.5.1 Theoretical maximum values for CO₂ storage. V, the effective CH₄ shale volume, for the upper Bowland unit is 9.31 × 10¹¹ m³ and for the lower Bowland unit is 3.45 × 10¹² m³.¹⁹E_v, the practical efficiency is set to be 1 for theoretical maximum value. This assumption implies that the potential storage volumes estimated here should be thought of as upper bounds. The mean porosity, Φ, is 4.0%.²⁰ρ_CO2 is the density of the liquid or ^scCO₂ which is 134.1 kg m⁻³ at reservoir conditions (80 °C and 7 MPa). The ρ_sCO2 values measured in the laboratory on BOR (TOC 6%; 4.5 kg per t) and BOL (TOC 1%, 2.9 kt per t) at 7 MPa and 80 °C (Fig. 8a) were combined to estimate ρ_sCO2E_s. Specifically, the TOC content of the Bowland Shale varies mainly between 2 and 2.5 wt%. We assume here a shale composition of 75% BOL (TOC 1%) and 25% BOR (TOC 6%), yielding a mean TOC 2.25%, as observed in the field. The mean density of Bowland shale is 2.6 t per m³.¹⁹E_s(m), the moisture equilibrated adsorption efficiency factor, is not measured in this study and set to be 0.78 based on the measurements of other shale samples e.g.ref. 34 with a similar moisture equilibration status as Bowland shale (average 32.4% moisture in Bowland shale).

Therefore, ρ_sCO2E_s = ρ_sCO2E_s(d)E_s(m) = (4.5 kg per t × 0.25 + 2.9 kg per t × 0.75) × 2.6 t per m³ × 0.78 = 6.7 kg m⁻³.

91.5% of total mineral pores and all organic matter pores are CO₂ and CH₄ accessible, so E_Φ, the pore connectivity efficiency, are calculated to be 0.95 (= 1 − 0.082 (non-assessable pores in mineral pores) × 0.658 (mineral porosity in total porosity)).

Total porosity of BOR sample is calculated to be 7.8% of total sample volume (Fig. 8d). Specifically, pores above 2 nm and 5 nm occupied 94% and 85% of the total porosity respectively, corresponding to the min and max values of the free gas porosity (Section 2.5). When converted to mean porosity of Bowland shale (4%), Φ (above 2 nm) is 3.8% (= 4% × 94%) and Φ (above 5 nm) is 3.4% (= 4% × 85%). These values were used to define the boundaries to calculate the total amount of stored CO₂ using the eqn (5) and (6) which were derived from eqn (4).

In each 1 m³ shale with 32.4% water saturation, 4.3–4.8 kg CO₂ can be stored as free molecules and ∼6.4 kg CO₂ can be adsorbed, assuming full efficiency of CO₂ shale and sorbed volume. This can be converted into 2.0–2.2 kg per t and 3.1 kg per t respectively as the average density for shale is 2.6 kg m⁻³. Therefore, 19.0–21.2 Gt CO₂ is estimated to be stored as free gas and 28.3–28.5 Gt CO₂ is estimated to be adsorbed on pore surface, consisting of 47.5–49.5 Gt carbon storage as a theoretical maximum in Bowland shale. It is noted that these values were calculated based on the assumption of all gas shales in Bowland shale are used for CO₂ (i.e. E_v = 1) owing to the lack of reservoir-scale measurements.

3.5.2 A simple practical estimation for CO₂ storage. The practical storage potential will be smaller than the theoretical maximum values shown above. Here we provide two simple estimations for practical stored volumes. Estimate (1) was calculated on basis of E_v = 0.32.³⁵ Estimate (2) was calculated using the simulated P10 (0.15 for E_Φ and 0.11 for E_s) and P90 (0.36 for E_Φ and 0.24 for E_s) values of E_v after 60 years’ injection based on the studies from Myshakin et al., 2018 (Table 4).

Table 4 Practical estimates for CO₂ storage based on the simulated P10 and P90 values from Myshakin et al., 2018

		Total storage (Gt)	Free gas (Gt)	Adsorbed gas (Gt)
Theoretical maximum		47.5–49.5	19.0–21.2	28.3–28.5
Practical estimation (1)	32%	15.2–15.8	6.1–6.8	9.1
Practical estimation (2)	P10	6.0–6.3	2.9–3.2	3.1
	P90	13.5–14.3	6.8–7.6	6.7

---

The practical estimation of total gas storage is 6.0–15.8 Gt with 2.9–7.6 Gt free gas and 3.1–9.1 Gt adsorbed gas. It is noted that the exact practical values will need to be refined when the reservoirs-scale efficiency factors are measured in future studies and the ethically and economically stored volumes will consequently be left for those studies.

4 Discussion

4.1 Heterogeneity at multiple scales

The heterogeneity analysis relies on the targeting features and the corresponding scales. The variation of formation depths and thickness at basin-scale,¹⁹ lithofacies and microfacies at core-scale,²¹ and porosity, permeability and mineral compositions of bulk samples²⁰ can be found in many other studies. Based on these understandings, the heterogeneity is further extended to sub-cm-scales in this study through the distribution of mineral clusters and matrix at 1–10 μm scale, individual grains of minerals and organic matter at 100 nm–1 μm scale, macropores and microcracks at 10–100 nm scale, and mesopores at 1–10 nm scale. These features are crucial to understand the pore-structure, flow pathways, and gas behaviors for accurate prediction of carbon sequestration and enhanced gas recovery.

For the same component, increased resolution of images leads to a lower fractal dimension owing to more details of complexity and roughness quantified. At 1–10 μm scale, each component is present as assemblies without clear boundaries of individual grains, therefore the fractal dimensions do not vary dramatically between components (∼2.5 to 2.6). Individual grains of organic matter and minerals can be distinguished at 100 nm–1 μm scale, and the fractal dimensions, representing for the surface roughness of these grains, decrease to ∼2.3–2.5. Pores can be only observed at LR- and HR-nanoscale, and the fractal dimension of macropores and fractures is relatively high when the voxel size is 50 nm and it decreases to 2.2–2.3 and when the voxel sizes reach 20 nm and the mesopores are taken into account.

The degree of anisotropy varies among different components. Clay mineral shows higher degree of anisotropy at both LR- and 100 nm–1 μm scales compared with other components, owing to the elongate geometry and the preferred orientation to the bedding planes. It is as high as 0.758 at high-resolution micro-scale, much higher value than 0.181 for granular minerals and 0.171 heavy minerals which have relatively regular crystal geometries. Organic matter pores show slightly higher degree of anisotropy than minerals pores, owing to the different mechanism of origin.

4.2 Potential flow pathways and storage of CO₂ and CH₄

A full understanding of the pore system is the basis of understanding flow pathways and storage volumes. The mechanism of CH₄ and CO₂ flow transport in shale reservoirs can take several forms (e.g. Darcy flow, surface diffusion and molecular diffusion), and each depends on the pore sizes and pore-network size.^38,39 CO₂ is present as supercritical state at the reservoir conditions in Bowland shale, which has low viscosity and high diffusivity like gas but density similar to liquid CO₂.⁴⁰ The variations of subsurface condition factors such as temperature, pressure and moisture status will impact the storage volume, flow pathways and fluid behaviours in the subsurface. The adsorption and permeability were both measured at reservoir pressures to simulate the subsurface environment, in order to provide better estimates. The pore systems of Bowland Shale samples are primarily associated with mineral matrix (60.2%) with wedge shapes or within organic matter (34.2%) with spherical shapes in terms of total pore numbers. Three size ranges of pores were detected in this study playing different roles for the fluid pathways. The micro-cracks and pores above 100 nm, shown in the PFIB images, contribute 47.5% in total pore volume and many of them exist at the interface of organic matter and minerals with relatively elongate shape (Fig. 5). These pores have extremely low connectivity themselves but they can be connected to smaller pores. They may act as flow porosity, in the form of Darcy flow, allowing the CH₄ and CO2 exchange between different compositions. The CH₄/CO₂ displacement efficiency in these pores changes relatively rapid.⁴¹ Therefore once the CO₂ reaches these pores, the CH₄ adsorbed on the pore surface can be extracted more easily. The injected CO₂ flow travels between mineral matrix and organic matter through these pores and the extracted CH₄ gas (both free and adsorbed gas) is transported from the organic matter to mineral matrix. The pores between 10 nm to 100 nm, majority of which were shown in the FIB images, contribute 35.9% for the pore volume and therefore provide the storage space for the free molecules of CH₄ and CO₂ and the pathways to the smaller pores (less than 10 nm size). Pores below 10 nm, identified by the nitrogen adsorption, SEM and HIM images, are associated with organic matter and clay minerals. They contribute over 99.9% of cumulative surface area and may provide the potential pathways in the alternative forms of molecular diffusion cross pores and adsorption–desorption at the surface of pores.^42–44 The adsorption capacity is also the highest in these pores. The capillary condensation phenomenon may be able to be characterised in these pores using relevant techniques. Additionally, 5.6% of pores (0.4% of the sample volume) are identified as intra-minerals (Fig. 1e), which are bounded in the mineral crystal, are likely to contribute little to the overall storage and flow transport for either CH₄ or CO₂ due to the extremely low connectivity relationship to other pores. In addition, the artificial fractures may connect many pores with different sizes and provide Darcy flow pathways for CH₄ and CO₂.

The CO₂ and CH₄ flow transport in shale reservoirs are complex. Based on the previous understanding e.g.ref. 2, 45 and 46–48 and the results in this study, a flow pathways and storage mechanism is proposed below (Fig. 9). At the beginning of the CO₂ injection, high concentration of CO₂ will start the transport from the induced fractures (20–80 μm) to the pre-existing microcracks (1–20 μm) and macropores (20 nm–1 μm), as a form of free molecules. Then CO₂ transports into the mesopores and micropores (<20 nm) in the forms of free molecules and surface diffusion, stored in the open pore space and adsorbed to the pore surface. As the concentration of CO₂ increases, CH₄ desorbs and CO₂ adsorbs, and CO₂ penetrates to smaller pores gradually. The surface diffusion takes more weight in the finer pores. When the transport reaches the equilibrium, the CO₂ may be stored in shales as the following forms: hydrodynamic trapping of buoyant CO₂ remains as a mobile fluid, adsorption trapping to the surface of the host rock, solution trapping of dissolving CO₂ in the formation water, mineral trapping due to mineral precipitation; and the residual or capillary trapping into an immobile fraction.

Fig. 9 Schematic diagram of the storage space and flow pathways for CO₂ and CH₄.

4.3 Feasibility of simultaneous CH₄ extraction and CO₂ storage in Bowland shale

The Bowland Shale meets the criteria for CO₂ storage in shale proposed by the US-DOE:³² (1) hydrocarbons can be commercially produced from the shale through stimulated fracture networks, which has been established by numerous studies e.g.ref. 19 and 20. The gas in place is estimated to be 822 (P90) – 1329 (P50) – 2281(P10) tcf.¹⁹ (2) The Bowland shale reservoir is found at a depth >1500 m, which can maintain CO₂ in a supercritical state at the reservoir's pressure and temperature conditions. These conditions were attained in this study for the adsorption experiments with CO₂ and CH₄. (3) A suitable seal system is in place to prevent the vertical migration of CO₂ to the surface. The extremely low vertical permeability of Bowland Shale measured in this study, 3.9 × 10⁻²¹ m² (∼3.9 nD, two orders of magnitudes smaller than the horizonal permeability), indicate that CO₂ will be retained in the subsurface once injected.

Sufficient pore volume for free molecules and pore surface areas for adsorption are two important parameters to assess the potential of CO₂ storage in shales. Based on the estimation in this study, ∼40% of the CO₂ will be stored as free molecules and ∼60% will be stored in as adsorbed molecules in the Bowland shale. The adsorption of these gases largely depends on the organic and clay content in shale, which at the nanoscale is reflected in distinct surface functionality.⁴⁹ Oxygen-containing groups including –COOH, –H and –OH bonds, have been shown to have a severe impact on the adsorption of CO₂ and CH₄ on organic matter and mineral surface at reservoir conditions.⁵⁰ Inorganic vibrational bands of quartz, illite and kaolinite (Si–O, Al–O–H and Si–O–Al bonds) as well as carbonates can exhibit higher CO₂ adsorption than hydrocarbons (C–H) on mineral surfaces.^51–53 CO₂ injection causes a reduction in the intensity of C–O and C–H bonds within organic matter in shales,⁵⁴ leading to an increase in hydrocarbon recovery from organic matter.⁵⁵ In practice, clay minerals have a significant impact on gas adsorption on shale,⁵⁶ although in the presence of water the organic content reveals a stronger correlation with gas adsorption than clay minerals.⁵⁷ In the Bowland shale, the clay mineral content is relatively low with an average value 20% and the TOC value normally varies between 1% and 6%. The selected organic-rich (silicate and carbonate 75%, kaolinite 7%, muscovite 10% and TOC 6%) and lean (silicate and carbonate 65%, kaolinite 18%, muscovite 5% and TOC 1%) samples are typical for the whole Bowland shale formation.^20,21 In our study, we observed that CO₂ adsorption is approximately 8 and 10 times larger than CH₄ adsorption in BOR and BOL, respectively. This is in line with other studies reporting 3–10 times greater adsorption of CO₂ than CH₄ in shales.^58,59 Importantly for this application, enhanced CH₄ gas recovery could result from CO₂ injection in Bowland shale. Furthermore, the low clay mineral content in Bowland Shale is expected to minimise the impact of clay swelling following water injection, meaning that the differences between water and CO₂ injection observed in the Longmaxi shale formation (Song et al., 2019) are likely to be significantly smaller here.

The cumulative porosity (including micro-cracks) is calculated to be 7.8% and the cumulative surface area is 1.79 m² g⁻¹ in the selected samples (Fig. 8c). We have estimated a theoretical maximum storage capacity in Bowland shale of 47.5–49.5 Gt, including 19.0–21.2 Gt CO₂ (38–45%) stored in the open pore space and 28.3–28.5 Gt CO₂ (55–62%) adsorbed onto the pore surfaces. Additional porosity can be produced by artificial fracturing (4.8% in the 20 mm diameter sample considered in this study), leading to an increase of storage space for free molecules. The enhancement of the storage capacity will therefore also depend on the implemented fracturing job.

4.4 Overall assessment of simultaneous CH₄ extraction and CO₂ storage in different shale plays

Compared to many successfully commercial gas plays in North America, Bowland shale is similar in composition (low clay mineral content), depth (>1000 m), and water saturation (<35%). While the porosity of Bowland shale is lower than Marcellus and Haynesville shale, it is similar to values reported for Ohio shale, New Albany shale and Barnett shale. This may lead to a reduced gas recovery using conventional extraction methods (e.g. hydraulic fracturing), but the overall CO₂ storage volume is still high (Table 5).

Table 5 CO₂ storage potentials of Bowland shale compared with other shale plays in the world

	Total estimated CO2 storage volume (Gt)	Free gas (Gt)	Adsorbed gas (Gt)	Estimated hydrocarbon GIP (TCM)	Thickness (m)	Area (km^2)	TOC (%)	Ro (%)	Dominated compositions	Clay content (%)	Porosity (%)	Water content (%)	Permeability (mD)	Ref.
Bowland shale	47.5–49.5	19.0–21.2	28.3–28.5	12.4	Up to 3000	>1000	2.5	1.7–2.4	Mixed clastic–carbonate	35	3.7–7.8	12.3–55.2	0.00003	19 and 20
Marcellus	171^3	72	99	410	12–270	27, 510	3.8	1.56	Mixed clastic–carbonate	5–40	1.95 to 7.56%	12–35	0.0003–0.95	3, 60, 64 and 65
	10–18^63
Utica	340	—	—	103	100–250	149, 995	1.1	1.17–2.12	Carbonate	20	2.6–3.9	—	<0.1	3 and 66
Antrim	72	—	—	21	95	95, 203	0.25–2.4	0.4–0.6	Quartz	25–50	9	—	<0.1	3, 60, 67 and 68
Ohio shale and New Albany Shale	21	—	—	11	200	118, 736	1.2–13	0.35–1.41	Mixed clastic–carbonate	16–42	4.7	—	—	3, 47, 67 and 69
Barnett	2.1–3.1	—	—	43.38	300	10, 554	3.2	2.25	Mixed	24.2	4–5	25%	0.01	60 and 70–72
Devonian black shales	27.7	—	—	7.4	175	22, 468	0.69–14	0.78–1.59	Clastic	39	4.3	—	0.0005	58 and 67
Haynesville	—	—	—	74.41	250	9256	4.2	2.37	Silicate	20–20	8–9%	15–20%	0.001–0.000001	17, 60, 67 and 73
Longmaxi shale	>1	—	—	40.02	50–600	1300	0.4–18.4	2.1–3.0	Mixed clastic–carbonate	33.87	1.2–10.8	1–8	0.25–1.737	74 and 75

---

The theoretical CO₂ storage value of Bowland shale (47.5–49.5 Gt CO₂) is comparable to many commercial shale plays in the United States. For Marcellus shale, the CO₂ storage capacity was estimated to be 171 Gt (99 Gt of adsorbed gas and 72 Gt of free gas), while 340 Gt, 72 Gt and 21 Gt CO₂ were estimated fro Utica shale, Antrim Shale and Devonian Ohio shale, respectively.³ A total of approximately 28 Gt of CO₂ was estimated to be potentially stored in Ohio shale and New Albany Shale.⁵⁸ The average CO₂ storage capacity of Bowland shale is 5.4 CO₂ kg per t, which represent the lower limit of the range reported for many shale plays (5 to 10 kg per t).⁶⁰ Additionally, the average recovery rate of shale gas is 25% with a CO₂ enhanced rate 7% which makes the final CO₂ recovery gas 32%. Meanwhile, the technically accessible CO₂ storage is around 32% in Marcellus shale.³⁵ The estimated total gas in-place estimate (P50) in Bowland shale is 1329 tcf (∼37.6 tcm). Based on technical values estimated above, 438.8 tcf (∼12.4 tcm) methane may be recovered by CO₂ injection and 6.5 (P10)–15.8 (P90) Gt CO₂ may be stored practically. This is equivalent to over 10 years' CO₂ emission in the UK according to the statistical figure (0.364 Gt = 80% × 0.455 Gt) in 2019 published by UK government.⁷⁶

Some inherent limitations to the extraction operation include (i) the low viscosity of ^scCO₂ (1/4 to 1/5 of water) leading to a poor sand-carrying ability;⁴⁰ and (ii) the low matrix permeability (3.9 × 10⁻²¹ m² normal to beddings and 3.3 × 10⁻¹⁹ m² parallel to beddings) which introduces severe mass transfer limitations.⁶¹ Technical developments such as cyclic injection of carbon dioxide, and huff and puff operations, may be considered to further increase the efficiency in field operations. Compared with water, CO₂ injection require a lower pressure for fracturing and generates more complex fracturing networks,⁶² while variations between induced microseismicity during water or CO₂ injection requires careful studies.

5 Conclusions

(1) This paper qualitatively and quantitatively investigates the 3D microstructure and pore system of the Bowland shale and characterise the possible flow pathways and storage capacity of CH₄ and CO₂ in shales. It demonstrates that the multi-scale 3D characterization is a powerful tool to assess the feasibility of CO₂ storage in shale and enhanced gas recovery and this can be used for further technical development to enhance the energy efficiency.

(2) The representative features at corresponding scales were identified, including fractures induced by fluid injection, sample fabric, individual grains, macropores, mesopores and micropores at 10–100 μm scale, 1–10 μm, 100 nm–1 μm, 10–100 nm, and 1–10 nm respectively. The particle and pore anisotropy is heterogeneous for the different phases and pore sizes at the same scale, and significantly changes across the scales from the μm to nm.

(3) The cumulative porosity is calculated to be 7.8% while the cumulative surface area is 1.79 m² g⁻¹. Pores associated with mineral matrix within organic matter occupy 60.2% and 34.2% of the total porosity respectively with another 5.6% intra-mineral pores which are not connected to the whole pore network. A pore size distribution spanning the range of 2 nm to 3 μm shows pore above 100 nm, between 10 to 100 nm and below 10 nm contribute 3.7%, 2.8% and 1.3% of the sample volume respectively, and pore below 10 nm contribute over 99.9% of the cumulative surface area.

(4) At the study conditions, a theoretical maximum of 49.5–49.5 Gt carbon storage can be stored in Bowland shale, consisting of ∼40% as free molecules in the open pore space and ∼60% as adsorbed molecules to the pore surfaces. The technical and economical storage of CO₂ in Bowland shales may be less (e.g. 6.0–15.8 Gt based on simple estimations) and the exact numbers can be estimated in future studies when reservoirs-scale factors are measured. The practical value might be equivalent to over 10 years* CO₂ emission in the UK.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This project was funded in part by the UK-NERC (NE/M001458/1 and NE/R013527/1) and the European Union's Horizon 2020 716 Research and Innovation Programme under the ShaleXenvironmenT project, (grant no. 640979). The authors thank the Diamond Light Source Beamline I13-2 for beamtimes (EE13824, MT9866), the Molecular Foundry (#5837), and the Research Complex at Harwell and Henry Mosely for laboratory space. National Research Facility for Lab X-ray CT (NXCT) through EPSRC grant EP/T02593X/1 is acknowledged for imaging processing and quantification support.

References

1. P. H. B. Zárante and J. R. Sodré, Evaluating carbon emissions reduction by use of natural gas as engine fuel, J. Nat. Gas Sci. Eng., 2009, 1(6), 216–220.
2. C. Oldenburg, K. Pruess and S. M. Benson, Process modeling of CO2 injection into natural gas reservoirs for carbon sequestration and enhanced gas recovery, Energy Fuels, 2001, 15(2), 293–298.
3. M. Godec, et al., Enhanced gas recovery and CO2 storage in gas shales: a summary review of its status and potential, Energy Procedia, 2014, 63, 5849–5857.
4. K. Louk, et al., Monitoring CO2 storage and enhanced gas recovery in unconventional shale reservoirs: Results from the Morgan County, Tennessee injection test, J. Nat. Gas Sci. Eng., 2017, 45, 11–25.
5. A. Busch, et al., Carbon dioxide storage potential of shales, Int. J. Greenhouse Gas Control, 2008, 2(3), 297–308.
6. Y. Shi, et al., Potential evaluation on CO2-EGR in tight and low-permeability reservoirs, Nat. Gas Ind. B, 2017, 4(4), 311–318.
7. Y. Jiang, et al., Effects of supercritical CO2 treatment time, pressure, and temperature on microstructure of shale, Energy, 2016, 97, 173–181.
8. L. M. Keller, et al., Characterization of multi-scale microstructural features in Opalinus Clay, Microporous Mesoporous Mater., 2013, 170, 83–94.
9. C. R. Clarkson, et al., Pore structure characterization of North American shale gas reservoirs using USANS/SANS, gas adsorption, and mercury intrusion, Fuel, 2013, 103, 606–616.
10. S. Lee, et al., Dehydration effect on the pore size, porosity, and fractal parameters of shale rocks: Ultrasmall-angle X-ray scattering study, Energy Fuels, 2014, 28(11), 6772–6779.
11. L. Wang, et al., Review of multi-scale and multi-physical simulation technologies for shale and tight gas reservoirs, J. Nat. Gas Sci. Eng., 2017, 37, 560–578.
12. J. Klaver, et al., BIB-SEM characterization of pore space morphology and distribution in postmature to overmature samples from the Haynesville and Bossier Shales, Mar. Pet. Geol., 2015, 59(0), 451–466.
13. X. Zhang, et al., Experimental study on fracture initiation and propagation in shale using supercritical carbon dioxide fracturing, Fuel, 2017, 190, 370–378.
14. L. Ma, et al., Multi-scale 3D characterisation of porosity and organic matter in shales with variable TOC content and thermal maturity: Examples from the Lublin and Baltic Basins, Poland and Lithuania, Int. J. Coal Geol., 2017, 180, 100–112.
15. L. Ma, et al., Novel 3D centimetre-to nano-scale quantification of an organic-rich mudstone: The Carboniferous Bowland Shale, Northern England, Mar. Pet. Geol., 2016, 72, 193–205.
16. M. E. Curtis, R. J. Ambrose and C. H. Sondergeld, Structural Characterization of Gas Shales on the Micro- and Nano-Scales, in Canadian Unconventional Resources and International Petroleum Conference, Society of Petroleum Engineers, 2010.
17. L. Ma, et al., Hierarchical integration of porosity in shales, Sci. Rep., 2018, 8(1), 11683.
18. H. E. King Jr, et al., Pore architecture and connectivity in gas shale, Energy Fuels, 2015, 29(3), 1375–1390.
19. I. Andrews, The Carboniferous Bowland Shale gas study: geology and resource estimation, 2013.
20. H. Clarke, et al., Shale gas resources of the Bowland Basin, NW England: a holistic study, Pet. Geosci., 2018, 24(3), 287–322.
21. S. M. Newport, et al., Sedimentology and microfacies of a mud-rich slope succession: in the Carboniferous Bowland Basin, NW England (UK), J. Geol. Soc., 2017, jgs2017.
22. A.-L. Fauchille, et al., An enhanced understanding of the Basinal Bowland shale in Lancashire (UK), through microtextural and mineralogical observations, Mar. Pet. Geol., 2017, 86, 1374–1390.
23. M. Chandler, et al., In-situ synchrotron X-ray tomography of fluid injection experiments, AGUFM, 2019, vol. 2019, p. MR41D-0073.
24. L. Ma, et al., A novel upscaling procedure for characterising heterogeneous shale porosity from nanometer-to millimetre-scale in 3D, Energy, 2019, 1285–1297.
25. N. Fazzalari and I. Parkinson, Fractal dimension and architecture of trabecular bone, J. Pathol., 1996, 178(1), 100–105.
26. A. Odgaard, Three-dimensional methods for quantification of cancellous bone architecture, Bone, 1997, 20(4), 315–328.
27. H. Ansari, et al., Supercritical adsorption in micro- and meso-porous carbons and its utilisation for textural characterisation, Microporous Mesoporous Mater., 2020, 308, 110537.
28. E. P. Barrett, L. G. Joyner and P. P. Halenda, The determination of pore volume and area distributions in porous substances. I. Computations from nitrogen isotherms, J. Am. Chem. Soc., 1951, 73(1), 373–380.
29. H. Marsh and B. Rand, The characterization of microporous carbons by means of the Dubinin-Radushkevich equation, J. Colloid Interface Sci., 1970, 33(1), 101–116.
30. K. Iyer and A. Kunju, Extension of Harkins—Jura adsorption isotherm to solute adsorption, Colloids Surf., 1992, 63(3–4), 235–240.
31. R. Mckernan, et al., Influence of Effective Pressure on Mudstone Matrix Permeability: Implications for Shale Gas Production, in SPE/EAGE European Unconventional Resources Conference and Exhibition, Society of Petroleum Engineers, 2014.
32. J. S. Levine, et al., U.S. DOE NETL methodology for estimating the prospective CO2 storage resource of shales at the national and regional scale, Int. J. Greenhouse Gas Control, 2016, 51, 81–94.
33. G. Chen, et al., Investigation of pore size effects on adsorption behavior of shale gas, Mar. Pet. Geol., 2019, 109, 1–8.
34. F. Tian, et al., Adsorption Behavior of CH4, C2H6, and CO2 on Moisture-Equilibrated Shale, Energy Fuels, 2020, 34(8), 9492–9497.
35. M. Godec, et al., Potential for enhanced gas recovery and CO2 storage in the Marcellus Shale in the Eastern United States, Int. J. Coal Geol., 2013, 118, 95–104.
36. E. M. Myshakin, et al., Numerical estimations of storage efficiency for the prospective CO2 storage resource of shales, Int. J. Greenhouse Gas Control, 2018, 76, 24–31.
37. R. Taylor, Interpretation of the correlation coefficient: a basic review, J. Diagn. Med. Sonogr., 1990, 6(1), 35–39.
38. F. Javadpour, Nanopores and Apparent Permeability of Gas Flow in Mudrocks (Shales and Siltstone), J. Can. Pet. Technol., 2009, 48(8), 16–21.
39. C. H. Sondergeld, et al., Petrophysical Considerations in Evaluating and Producing Shale Gas Resources, in SPE Unconventional Gas Conference, Society of Petroleum Engineers, 2010.
40. R. S. Middleton, et al., Shale gas and non-aqueous fracturing fluids: Opportunities and challenges for supercritical CO2, Appl. Energy, 2015, 147, 500–509.
41. J. Liu, et al., Influence of pore structure on shale gas recovery with CO2 sequestration: insight into molecular mechanisms, Energy Fuels, 2019, 34(2), 1240–1250.
42. L. Mi, H. Jiang and J. Li, The impact of diffusion type on multiscale discrete fracture model numerical simulation for shale gas, J. Nat. Gas Sci. Eng., 2014, 20(0), 74–81.
43. S. R. Etminan, et al., Measurement of gas storage processes in shale and of the molecular diffusion coefficient in kerogen, Int. J. Coal Geol., 2014, 123, 10–19.
44. W. Yuan, et al., Experimental study and modelling of methane adsorption and diffusion in shale, Fuel, 2014, 117(part A), 509–519.
45. D. Liu, et al., CO2 sequestration with enhanced shale gas recovery, Energy Sources, Part A, 2019, 1–11.
46. Y. Gensterblum, A. Busch and B. M. Krooss, Molecular concept and experimental evidence of competitive adsorption of H2O, CO2 and CH4 on organic material, Fuel, 2014, 115, 581–588.
47. F. Liu, et al., Assessing the feasibility of CO2 storage in the New Albany Shale (Devonian–Mississippian) with potential enhanced gas recovery using reservoir simulation, Int. J. Greenhouse Gas Control, 2013, 17, 111–126.
48. N. R. Backeberg, et al., Quantifying the anisotropy and tortuosity of permeable pathways in clay-rich mudstones using models based on X-ray tomography, Sci. Rep., 2017, 7(1), 14838.
49. Y. Liu and J. Wilcox, CO2 adsorption on carbon models of organic constituents of gas shale and coal, Environ. Sci. Technol., 2011, 45(2), 809–814.
50. P. J. Jarboe, et al., Extraction of hydrocarbons from high-maturity marcellus shale using supercritical carbon dioxide, Energy Fuels, 2015, 29(12), 7897–7909.
51. R. Heller and M. Zoback, Adsorption of methane and carbon dioxide on gas shale and pure mineral samples, J. Unconv. Oil Gas Resour., 2014, 8, 14–24.
52. J. Hwang and R. Pini, Supercritical CO2 and CH4 uptake by illite-smectite clay minerals, Environ. Sci. Technol., 2019, 53(19), 11588–11596.
53. X. Hu, et al., Characterization of CO2/CH4 competitive adsorption in various clay minerals in relation to shale gas recovery from molecular simulation, Energy Fuels, 2019, 33(9), 8202–8214.
54. Y. Pan, et al., Experimental investigation of the geochemical interactions between supercritical CO2 and shale: implications for CO2 storage in gas-bearing shale formations, Energy Fuels, 2018, 32(2), 1963–1978.
55. A. Fatah, et al., A Review on the Influence of CO2/Shale Interaction on Shale Properties: Implications of CCS in Shales, Energies, 2020, 13(12), 3200.
56. H. T. Schaef, et al., CO2 utilization and storage in shale gas reservoirs: experimental results and economic impacts, Energy Procedia, 2014, 63, 7844–7851.
57. I. Klewiah, et al., Review of experimental sorption studies of CO2 and CH4 in shales, J. Nat. Gas Sci. Eng., 2020, 73, 103045.
58. B. C. Nuttal, et al., Analysis of Devonian black shales in Kentucky for potential carbon dioxide sequestration and enhanced natural gas production, Greenhouse Gas Control Technologies 7, Elsevier, 2005, pp. 2225–2228.
59. S. M. Kang, et al., Carbon dioxide storage capacity of organic-rich shales, SPE J., 2011, 16(4), 842–855.
60. R. Khosrokhavar, S. Griffiths and K.-H. Wolf, Shale gas formations and their potential for carbon storage: opportunities and outlook, Environ. Process., 2014, 1(4), 595–611.
61. R. Xu, et al., Assessing the feasibility and CO2 storage capacity of CO2 enhanced shale gas recovery using Triple-Porosity reservoir model, Appl. Therm. Eng., 2017, 115, 1306–1314.
62. X. Song, et al., Fracturing with carbon dioxide: from microscopic mechanism to reservoir application, Joule, 2019, 3(8), 1913–1926.
63. Z. Tao and A. Clarens, Estimating the carbon sequestration capacity of shale formations using methane production rates, Environ. Sci. Technol., 2013, 47(19), 11318–11325.
64. L. Song, et al., Porosity and storage capacity of Middle Devonian shale: A function of thermal maturity, total organic carbon, and clay content, Fuel, 2019, 241, 1036–1044.
65. J. Goral, et al., Confinement effect on porosity and permeability of shales, Sci. Rep., 2020, 10(1), 1–11.
66. O. H. Ardakani, et al., Do all fractions of organic matter contribute equally in shale porosity? A case study from Upper Ordovician Utica Shale, southern Quebec, Canada, Mar. Pet. Geol., 2018, 92, 794–808.
67. U. S. Gas and S. O. Plays, Review of Emerging Resources, 2011.
68. A. P. Ruotsala, Mineralogy of Antrim Shale, Michigan, Dow Chemical Co., Midland, MI, USA, 1980.
69. M. Mastalerz, et al., Porosity of Devonian and Mississippian New Albany Shale across a maturation gradient: Insights from organic petrology, gas adsorption, and mercury intrusion, AAPG Bull., 2013, 97(10), 1621–1643.
70. S. L. Montgomery, et al., Mississippian Barnett Shale, Fort Worth basin, north-central Texas: Gas-shale play with multi–trillion cubic foot potential, AAPG Bull., 2005, 89(2), 155–175.
71. K. A. Bowker, Barnett shale gas production, Fort Worth Basin: Issues and discussion, AAPG Bull., 2007, 91(4), 523–533.
72. R. W. Edwards, et al., A model to estimate carbon dioxide injectivity and storage capacity for geological sequestration in shale gas wells, Environ. Sci. Technol., 2015, 49(15), 9222–9229.
73. P. J. Dowey and K. G. Taylor, Diagenetic mineral development within the Upper Jurassic Haynesville-Bossier Shale, USA, Sedimentology, 2020, 67(1), 47–77.
74. B. Zhang, et al., Review of Formation and Gas Characteristics in Shale Gas Reservoirs, Energies, 2020, 13(20), 5427.
75. S. Jiang, et al., Geology and shale gas resource potentials in the Sichuan Basin, China, Energy Explor. Exploit., 2016, 34(5), 689–710.
76. Waite, Final UK greenhouse gas emissions national statistics: 1990 to 2019, GOV.UK, retrieved 30 March 2021.

---