Use of Gas Adsorption and Inversion Methods for Shale Pore Structure Characterization

Abstract

The analysis of porosity and pore structure of shale rocks has received special attention in the last decades as unconventional reservoir hydrocarbons have become a larger parcel of the oil and gas market. A variety of techniques are available to provide a satisfactory description of these porous media. Some techniques are based on saturating the porous rock with a fluid to probe the pore structure. In this sense, gases have played an important role in porosity and pore structure characterization, particularly for the analysis of pore size and shapes and storage or intake capacity. In this review, we discuss the use of various gases, with emphasis on N₂ and CO₂, for characterization of shale pore architecture. We describe the state of the art on the related inversion methods for processing the corresponding isotherms and the procedure to obtain surface area and pore-size distribution. The state of the art is based on the collation of publications in the last 10 years. Limitations of the gas adsorption technique and the associated inversion methods as well as the most suitable scenario for its application are presented in this review. Finally, we discuss the future of gas adsorption for shale characterization, which we believe will rely on hybridization with other techniques to overcome some of the limitations.

1. Introduction

Shale rocks have received special attention since they arose as an important source of hydrocarbons as well as for their potential for decades ahead [1,2,3]. Although characterizing shale is not new, the most significant efforts were originally devoted to the understanding of their characteristics in conventional basins, where they play a role in wellbore stability or as sealing rocks [4,5,6]. For shale reservoirs, the need for a deeper understanding of the rock properties and their economic viability have become the most important research targets. Different factors can determine if a shale reservoir is economically viable. Among them, the quality and quantity of the organic matter (OM), mineralogy, brittleness, and petrophysical properties such as porosity, including pore structure, and permeability are the most relevant [7,8]. Therefore, laboratory analysis methods to estimate and evaluate these properties are in demand. Several techniques have been borrowed from those used for characterization of conventional reservoirs. However, the complexity of shales requires adaptation of these techniques, the incorporation of analysis methods from other fields as well as the implementation of new technologies [7,9]. It is worth noticing that even though all properties/characteristics of shale relate to one another, there is no established protocol that would allow their combined analysis.

Organic matter plays an important role in the evaluation of a shale reservoir, as it relates to the kerogen and the hydrocarbons in the rock. Generally, thermal maturity, related to vitrinite reflectance, the type of kerogen and the mass percentage also known as Total Organic Carbon (TOC), are key features that can determine the quality of a reservoir [2,7,8,10]. Productive reservoirs have been observed to have kerogen types I, II and II, with a TOC (%) greater than 2% and a vitrinite reflectance greater than 1.2. The analysis of these properties is not new, since these rocks were originally studied as source rocks [11,12,13,14]. To analyze the organic component of the reservoir, several techniques are used, e.g., the type of kerogen commonly classified through Rock-Eval pyrolysis analysis and the hydrogen index/oxygen index, combustion infrared detection used for TOC determination and microscopy used for the determination of vitrinite reflectance [9,10,15,16].

In quite a few cases, the inorganic component of the rock in shale reservoirs frequently differs from that in conventional reservoir. However, its relative high content of clay minerals has been the reason to classify them as “mudstones” [17]. The rock mineral analysis often relies on X-Ray Diffraction (XRD) and X-Ray Fluorescence (XRF) [7,9,18,19,20,21,22,23]. These two techniques can provide a detailed description of the minerals in the rock, though they are typically semi-quantitative techniques regarding their ability to determine mass percentages, depending on the processing method used. There are also physical or mechanical properties related to the inorganic component of the rock, a list that includes strength, elasticity, and extensional and compressive strength among others [24,25].

Although many properties of shale reservoirs are relevant and cross-correlated, the growing literature on shales in the last decade shows that petrophysical properties remain of special interest and play a determinant role in understanding unconventional systems [26,27,28,29,30,31,32,33,34,35]. It is fair to say that petrophysics is responsible for differentiating a shale from a conventional reservoir, e.g., permeabilities are often in the nano or millidarcies regime and the porosity often has a significant representation in the micro and mesopore regions. Most studies in the literature also provide information regarding the dependence of the petrophysical properties on shale characteristics as inorganic and organic phases. The term petrophysical properties can include different rock characteristics, but its analysis mainly focuses on permeability and porosity, and several methods are available to study them. The analysis of shale permeability remains, to an extent, like the analysis of conventional rocks such as sandstones and carbonates. Some of the most used techniques include Mercury Intrusion Porosimetry (MIP) and the use of different gases and differential pressure, while most of the techniques to estimate permeability rely on the use of Darcy’s law. Although permeability has been widely studied, it still represents a challenge as shown by several works in the literature [36,37,38,39,40,41,42].

Porosity and pore architecture [43] are as important for a shale reservoir as permeability is. The study of shales pore structure has received special attention due to the complexity observed in this type of rock matrix. Even though shale reservoirs can contain fractures, microfractures, and large pores, a significant fraction of the porosity in shale is in the nanometer range and often gets expressed through heterogeneities throughout the reservoir. The nanopores in shales have been classified in different ways. According to the IUPAC [44], pores in the nanometer regime can be divided into micropores (<2 nm), mesopores (between 2 and 50 nm), and macropores (>50 nm). Due to the significant fraction of nanoporosity in shale, the search for techniques that allow the analysis of the pore structure at this scale has become prevalent. According to Anovitz et al., 2015 [45], there are several methods available, often classified as radiation methods such as optical microscopy, SEM, AFM, TEM, FID-SEM, and SANS/USANS SAX. In contrast, the so-called saturation, or fluid-based methods, such as Mercury Intrusion Porosimetry and Helium Pycnometry, are also available. The saturation methods list includes the gas adsorption category, which implies the use of a variety of gases, e.g., N₂ and CO₂, CH₄, Ar.

Gas adsorption includes the execution of static and dynamic experiments. The dynamic experiments are often used for the analysis of flow properties (permeability). On the other hand, the static experiments are related to the determination of storage capacity and pore architecture characterizations. Significant efforts related to gas adsorption has been made to understand the behavior of gases such as, CO₂ and CH₄, in the complex shale porous media [46,47]. The analysis of these of these gases is related to the fact that CO₂ is an environmental solution for formation stimulation and CH₄ is the most produced gas from unconventional reservoirs. However, the use this technique for pore architecture characterization reflects some problems and less progress. Therefore, this is the main scope for this review. Gas adsorption has been widely used for shale pore structure characterization over the last decade. Some work on the use of this technique can be traced back to the late 1960s and 1970s [48,49,50]. However, this method was not originally designed for the analysis of shale but for the analysis of mesoporous materials widely used in the chemical engineering field (e.g., zeolites and carbons). The advantage of the use of gas adsorption relies on its relatively low cost and its ability to provide useful information regarding pore architecture and a description of the pore structure (pore-size distribution). However, the fidelity of this method to describe a porous medium relies on the selected inversion methods. Whereas classical methods include the use of BET (Brunauer–Emmett–Teller) for surface area and BJH (Barrett–Joyner–Halenda) for pore-size distribution methods or FFH (Frenkel–Halsey–Hill) for fractal dimension analysis, more advanced interpretation techniques involve the use of DFT (Density Functional Theory), NLDFT (Non-Local Density Functional Theory), Molecular Dynamics and Montecarlo simulations. Although gas adsorption and the associated inversion methods have been widely used and accepted, the available literature does not reflect efforts to improve the performance on the technique or an agreement regarding the selection of an adsorbate or an inversion method, specifically for shale pore architecture characterization. We collated a significant fraction of the available literature (200+ published works) in the past 10 years to draw the main directions on shale characterization. The search of the available literature was performed mainly through Web of Science database, using keywords such as: shale, unconventional, gas adsorption, physisorption, pore structure, pore structure characterization, pore size distribution, nitrogen (N₂), carbon dioxide (CO₂), argon (Ar), krypton (Kr). The obtained list of references was sorted out to fit the scope of this review. Some additional works from other databases were considered due to the importance for the review.

The main purpose is to guide the selection of gas adsorption toward pore structure characterization of shale systems and the available inversion methods and the objective of this review is threefold. First, we summarize the state of the art and the theory behind gas adsorption in the context of shale rock. This serves to highlight advantages of these techniques in contrast with other characterization methods used for this type of rock. We also describe drawbacks that derive either from the technique itself or the interpretation methods. Second, we provide a potential workflow to assist users of the technique to select combinations of gases that offer the best potential to draw information on the rock system, in our opinion. Thirdly and finally, we speculate on potential improvements that could potentially resolve some of the issues associated with the use of these techniques, including hybridization with other techniques. In this sense, the article is organized the following way. We start by providing a brief description of the adsorption mechanism on a solid surface and summarize the theory behind physisorption. Then, we include a description of the most used gases for shale characterization, followed by the most used inversion method to describe the pore architecture. Finally, we include a discussion regarding the advantages and drawbacks of gas adsorption with respect to other techniques and posit our ideas of where the future of the technique could be.

2. Adsorption Mechanisms

The determination of the pore volume or the pore-size distribution through gas-adsorption isotherm methods relies on mechanisms driving the adsorption of a gas. In the case of shales, this method is dominated by physisorption or physical adsorption of the commonly used gases, e.g., CO₂, N₂, Ar [51,52,53,54]. The adsorption process consists of the accumulation of gas molecules on a surface (usually a solid) through weak van der Waals and electrostatic interactions, giving gas molecules the ability to form a single layer or multiple ones on the surface. The density of one or multiple layers will commonly differ from the bulk phase (Figure 1). The physisorption process can be connected to condensation in some cases [55,56,57,58,59].

Results from a physisorption experiment are commonly represented in terms of the amount of gas adsorbed by one gram of solid at a certain pressure at constant temperature, reason why it is consequently denominated an isotherm. Several techniques allow the measurement of adsorbed gas, but the most common methods are either the volumetric (or manometric) or gravimetric methods [51,58,60]. The volumetric method is based on the measurement of the gas pressure in a calibrated, constant volume at a known temperature. Meanwhile, the gravimetric approach relies on different devices to estimate changes in the mass of the samples at different pressures. Simplified schematics of these two devices are shown in Figure 2. Rouquerol et al., 2014, provide a detailed description and the math behind each method. However, additional works in the available literature discuss the techniques to build the isotherms [61]. Although gas adsorption has been used for decades and a standard for the analysis of certain material exists [62], in the world of shales there is no agreement on the use of one between the two methods (volumetric and gravimetric).

The literature shows that volumetric methods are the most common for low-pressure adsorption and therefore used for porosity and pore structure analysis. Proof of that is that the most used instruments are from Micromeritics and Quantachrome, which rely on volumetric expansion for determining the volume of gas adsorbed. On the other hand, gravimetric methods are more commonly used for high-pressure adsorption experiments, which for shales commonly relate to storage capacity. However, some publications report the use of both methods for shale characterization [63,64]. It is important to consider that both techniques rely on a corrected version of the gas law to determine changes in gas volume related to changes in pressure.

The result of the physisorption experiment is the isotherm, which provides preliminary information on the adsorbent, e.g., surface area and porosity. The shape of the isotherm depends on the interaction of the gas molecules and the solid surface as well as on their physical characteristics. The International Union of Pure and Applied Chemistry (IUPAC) [65] provides a classification of different types of isotherms as well as an interpretation to relate them with characteristics of the porous system. The classification shows six different types of isotherms. Type I isotherms can be represented by the Langmuir equation, meaning that the solid reaches maximum monolayer adsorption. This type of isotherms is often related to microporous solids with a small surface area. Type I isotherms often exhibit a linear relationship at low pressure and reach a plateau as pressure increases [56]. Type II represents systems where more than one layer of the absorbed gas is formed on the solid surface. Type III has the peculiarity of exhibiting a concave shape over the entire range of relative pressure that indicates a weak adsorbent-adsorbate interaction. Type IV is commonly observed on mesoporous material, being Type II-like, showing an inflection point as well, which suggests the formation of multi-layers. The main characteristic of the Type IV isotherm is the presence of a hysteresis loop that is indicative of capillary condensation. Type V is similar to Type III, but it also exhibits a hysteresis loop. Type VI shows as stepwise multilayer adsorption and is the less common of all. Several works in the literature provide a more detail description of each type of isotherm in addition to the information provided by the IUPAC [51,58,66]. In addition to the classification of the isotherms, the IUPAC provides an additional classification for isotherms with the presence of a hysteresis loop. The classification goes from type H1 to type H4, where H1 is related to porous material formed by agglomerates of approximately uniform spheres. H2 is commonly related to a system with not well-defined pore structures, but in the past, it has been related to ink-bottle shaped pores. Type H3 and H4 have been shown to relate to solid formed by aggregates with plate-like shapes that derive into slit-like pores.

In the case of shales, which is the focus of this work, different gases have been used for characterization purposes, e.g., CH₄, CO₂, N₂, Ar. CH₄ and CO₂ are often used for the determination of storage capacity, through high pressure adsorption experiments. These experiments usually do not include desorption analysis and do not imply any condensation process. The most used gas for the acquisition of total pore volume and pore-size distribution is N₂. From the classification presented by Rouquerol, 2014 [58], it is possible to say that the most observed isotherm corresponds to a type IV (due to the presence of the hysteresis loop), with a hysteresis loop of the types H3 or H4 (Figure 3).

As was previously described, isotherms in Figure 3 correspond to pore structures formed by aggregates with plate-like shapes giving some slit-like pores. Shales have pores of different types and shapes, but the significant amount of clays, crystals and other plate-like minerals provide a significant amount of slit pores [26]. There are additional classifications of types of hysteresis loops that relate the final isotherm with a specific pore shape [67,68]. In the case of shales, some works in the literature have used the classification reported by Boer 1958 [67], where five different types of hysteresis loop are shown (Figure 4). Type A loop is represented mainly by pore shaped as cylinders, Type B by slit pores, Types C and D by wedge pores, and Type E bottle neck [69,70].

From this last classification, the most common isotherms observed for shale systems are of type A and B. However, the complexity of shale pore structure complicates the proposal of a specific isotherm that can be related to a specific geometry, and it is possible to say that the porosity of these systems is a combination of multiple pore shapes. Despite the complication involved with the complexities of shales, similarities to existing isotherms allow researchers to select an adequate inversion method to obtain information as pore volume and pore-size distribution. Inversion methods will be covered later in this review.

3. Adsorptive Capabilities and Selection

As described in the previous section, different gases can be used for shale characterization. However, choosing an appropriate gas is critical to obtain the correct information.

Table 1. Gases used for shale characterization.

Name	Formula	Kinetic Diameter (nm)	Common Use in Shales
Methane	CH4	0.380	Storage capacity, Permeability
Helium	He	0.260	Porosity, void volume
Krypton	Kr	0.360	Permeability, Porosity, void volume
Argon	Ar	0.340	Permeability, Porosity, Pore size distribution
Nitrogen	N2	0.364	Permeability, Porosity, Pore size distribution, Surface area
Carbon Dioxide	CO2	0.330	Storage capacity, Permeability, Porosity, Pore size distribution

shows a list of gases that have been used for material characterization as well as some examples of how those gases have been used for shale characterization purposes.

Methane is one of the most important gases used for characterization and is closely related to shale, since over 60% of natural gas is produced from unconventional reservoirs [71], where CH₄ is the main component. Most measurements are addressed to understand the behavior of this gas in the reservoir. Storage capacity and permeability measurements are the main properties of interest, meaning that performed experiments that involve the use of methane are related to the analysis of the behavior of this gas in porous media [34,68,69,70,71,72,73,74,75,76,77,78,79,80,81,82,83,84]. Although CH₄ is not directly used for porosity or pore structure characterization, there is evidence that it can provide information regarding pore architecture alteration [43,85]. Furthermore, in Medina-Rodriguez et al. [43], we demonstrated that if complementary techniques are combined, namely, TD-NMR and gas adsorption, it is possible to analyze changes in the surface properties of the rock wettability.

Helium is probably the most used gas in the world for shale characterization. The ability of the He molecule to reach small pores and to interact insignificantly with the solid walls makes it the most appropriate gas for porosity analysis (pycnometry). Another reason why it is the most used gas relies on how it is used almost exclusively on every adsorption experiment to measure void volumes. For this reason, we are confident to say that Helium is present in all works considered for this review. Despite the importance of He for porosity analysis, it is not capable of providing a description of the pore structure.

Nitrogen has been accepted as a standard for surface area and pore structure analysis of mesoporous materials and the adsorption experiments are performed at 77 K. N₂ also has shown to not be the most adequate adsorbate for the analysis of microporous structures. The quadrupolar moment of the nitrogen molecule is responsible for the problem associated with the adsorption of N₂. If nitrogen is used to analyze polar materials in pores, the pore-filling pressure can be shifted to low relative pressure and will not correlate with micropore sizes, e.g., carbons and metal organic frameworks [86,87,88]. For shales, nitrogen adsorption has shown to be significantly impacted by the presence of organic matter and clays, where organic matter has shown to constrain accessibility of the gas molecules and the amount and type of clays will have a direct impact on the cation exchange capacity of the surface [89]. The constraints provided by the N₂ quadrupolar moment can be overcome by using argon as an adsorbate. Ar adsorption is not affected by any specific interaction with any type of solid surface. This characteristic, in addition to the slightly smaller kinetic diameter, allows argon molecules to access the micropores better and provides a superior description. Although there have been attempts to use Argon for shale pore characterization [53], it is not the popular adsorbate. This is probably due to the temperature needed to perform the experiment (87 K). However, argon is present in shale characterization in other ways. It has been used for ion milling/polishing associated with imaging techniques as well as for permeability measurements [90,91,92,93,94,95]. Finally, if the purpose is to acquire a better description of the micropores, the most appropriate adsorbate is carbon dioxide. The smaller kinetic diameter and the use of higher pressures and temperature (273 K) for adsorption makes it capable to access micropores easily. For this reason, CO₂ has been chosen as a standard adsorbate for microporous carbon characterization [87,96] and widely used for shales.

Considering the aforementioned statements, Figure 5 provides percentages of use of these gases in the world of shales. In order of importance, nitrogen shows the largest use at 64%. Carbon dioxide occupies the second place at 33%. However, it is worth noticing that CO₂ has been increasingly used particularly for storage purposes. Additionally, literature shows that N₂ and CO₂ are often used together to obtain full coverage of accessible micro and mesopores. Figure 5 does not separate the individual use of any gas.

4. Inversion Methods

Once an isotherm has been collected, the next step is to select the inversion technique that will support the interpretation of the porous medium characteristics. Not all the available techniques widely used have been exclusively developed for shale pore structure description but have proven to provide an acceptable description of unconventional porous media. For a better organization of this section, we decided to divide the available inversion methods into two groups. The first group is identified as Classical Methods, while the second group is labelled as Advance Techniques. In the first group, we briefly describe the traditional BET (Brunauer–Emmett–Teller), BJH (Barrett–Joyner–Halenda) and Dubinin methods (D-R, D-A). In the second group, we include Molecular Dynamics Simulations, Monte Carlo Simulations, DFT (Density Functional Theory), NLDFT (Non-Local Density Functional Theory), and QSDFT (Quenched Solid Density Functional Theory). Finally, we provide information on other methods that despite their relationship with the aforementioned two groups, their application for shale characterization is not as significant.

Figure 6 shows the percentage of application to shale of each inversion method category, corresponding to the percentage of articles in the open literature. It turns out that Classical Methods dominate notoriously over the Advanced Methods.

In addition, this section aims at providing the most important equations as well as highlighting references to guide the reader in case additional details are required. A table in the Supplementary Materials contains the references, selected for preparation of Figure 5 and Figure 6.

4.1. Classical Methods

○

Surface Area (BET)

The use of the BET theory for the determination of surface area is perhaps the most widely accepted and used method for shale characterization. At least 90% of published works in this field have reported its use. The method was originally developed by Bruanuer, Emmett, and Teller. It has been widely used to describe the multilayer adsorption of gas molecules on a surface. The method allows the determination of the surface area on the analyzed solid through the transformation of the isotherm to a linearized BET plot that allows the determination of the monolayer coverage. The traditional BET equation follows:

$$\frac{1}{v\left[\left(\frac{p_o}{p}\right)-1\right]}=\frac{C-1}{v_{m}C}\left(\frac{p}{p_o}\right)+\frac{1}{v_{m}C}$$

(1)

where $v$ is the adsorbed amount, p_o is the saturation pressure of the adsorbate, $p$ is the equilibrium pressure, $v_m$ is the monolater adsorbed gas, and $C$ is the BET constant that can be expressed as $C=\exp \left(\frac{E_1-E_L}{RT}\right)$; $E_1$ is the heat of the adsorption for the first layer, and $E_L$ is the heat of vaporization. Through the linear plot shown in Figure 7, we obtain the slope and the intercept.

Once linear-fit parameters are obtained, one can determine the total the volume of the monolayer and subsequently the total surface area, which if divided by the total mass of the samples will provide the BET surface area of the samples through Equation (2).

$$S_t=\frac{v_{m}Ns}{V} \to S_{BET}=\frac{S_t}{m}\left[\frac{m^2}{g}\right]$$

(2)

where $S_t$ is the total surface area, $v_m$ is monolayer absorbed gas volume, $N$ is the Avogadro’s number, $s$ is cross-sectional area of the gas molecule, $V$ is the molar volume of adsorbed gas, and m is the sample mass.

It has been reported that traditional BET does not accurately represent the total surface area of porous media containing a significant amount of micropores. Therefore, a modified BET has been proposed, leading to Equation (3).

$$v=\frac{\left(\frac{v_m{C}_p}{p_o}\right)}{\left(1-\frac{p}{p_o}\right)\ast \left[1+\frac{\left(C-1\right)p}{p_o}\right]}+v_{mic}$$

(3)

where a new parameter is included, $v_{mic}$, which denotates the volume adsorbed by micropores. The BET constant C is not affected by micropores. The unknown parameters can be drawn using a non-linear fitting procedure that allows the determination of the total surface area. It is important to state that despite the error introduced by the traditional BET equation, it has been shown to be the first selection for the determination of surface area of shales.

○

Pore Size Distribution (BJH)

After BET, BJH is the next most frequently used method for isotherm analysis. The method allows the determination of the pore-size distribution of a mesoporous solid [58,60]. The physical principle behind the method relies on Kelvin’s equation for the analysis of condensation in pores and is expressed as

$$\ln \left(\frac{P}{P_o}\right)=-\frac{2\gamma V}{rRT}cos\theta$$

(4)

where $P$ is the equilibrium pressure in the pore, $P_o$ is the equilibrium saturation pressure in a planar surface, $r$ is the pore radius, $\gamma$ is the surface tension, $V$ is the molar volume of the liquid phase, and $\theta$ is the angle formed by the meniscus, traditionally denominate as the contact angle.

In general terms, the relationship between the gas adsorption and the pore structure of the solid develop by Barrett, Joyner, and Halenda can be summarized as follows:

$$v\left(\frac{p_j}{p_o}\right)={\displaystyle \sum }_{i=1}^j\Delta V_i\left(r_i\le r_c\left(\frac{p_j}{p_o}\right)\right)+{\displaystyle \sum }_{i=j+1}^n\Delta S_i{t}_i\left(r_i>r_c\left(\frac{p_j}{p_o}\right)\right)$$

(5)

where $v$ is the volume of adsorbate at a specific relative pressure, r represents the radius, $V$ is the pore volume, $S$ is the surface area, and $t$ is the thickness of the adjacent adsorbed layer formed on the pore wall. The expression can be divided into two parts. The first relates to the volume that condensates in all pores smaller than a specific size or radius (“$r_c$”) that depends on the relative pressure. The second part corresponds to a volume related to the adsorbed film that is estimated by the sum of the relation between the surface area and the thickness of the film in the pore. Some of the terms, such as the surface area and volumes, can depend on the assumed geometry [97]. However, in most cases, the pore shape is assumed cylindrical. The surface area is often determined by BET, and the thickness of the film can be obtained from the available models and equations (

Table 2. Thickness equation commonly used for pore structure characterization.

Harkins and Jura [98,99,100]	$t={\left(\frac{13.99}{0.034-{\log }_{10}\left(\frac{P}{P_o}\right)}\right)}^{\frac{1}{2}}$
Halsey [101]	$t=3.54{\left(-\frac{5}{\ln \left(\frac{P}{P_o}\right)}\right)}^{\frac{1}{3}}$
Broekhoff–de Boer [102]	${\log }_{10}\left(\frac{P}{P_o}\right)=-\frac{16.11}{t^2}+0.1682e^{-0.1137}$
Kruk, Jaroniec, and Sayari [103,104]	$t={\left(\frac{60.65}{0.03071-{\log }_{10}\left(\frac{P}{P_o}\right)}\right)}^{0.3968}$
Carbon Black STSA [105]	$t=2.98+6.45\left(\frac{P}{P_o}\right)+0.88{\left(\frac{P}{P_o}\right)}^2$

). Regardless of whether the adsorption or desorption isotherms are used, the calculation remains similar and will always start from the point where all pores are assumed to be filled.

It is important to know that estimating the real thickness of the adsorbed layer may have several implications on the analysis of several points on the solid surface. However, to simplify the process, the use of an average thickness for the entire surface provides some efficiency in the determination of the pore structure of a solid. In some case, corrections can be applied to the listed equations to improve their performance. For shale analysis, the most used thickness equations are the Harkins and Jura equation and the Halsey equation. There are additional thickness equations available in the literature, but they are excluded from this work due to their minimal appearance in shale-related works. The thickness equations are used in the t-plot analysis of microporosity. However, other methods (D-A, D-R) are preferred since they provide a better description of the micropores.

○

Dubinin equations.

Considering the restrictions and limitations that traditional BET and BJH have for describing the microporosity of a system, Dubinin equations represent a good alternative for an adequate description and are mostly used for CO₂ adsorption. These group of equations is formed by the Dubinin–Radushkevich (D-R) and the Dubinin–Astakhov (D-A) [106]. The latter equation can be considered a generalization of the former. Both expressions are shown in

Table 3. Dubinin equations.

Dubinin–Radushkevich	$\frac{V}{V_o}=\exp \left[-{\left(\frac{Rt}{\beta E}\ln \left(\frac{P_o}{P}\right)\right)}^2\right]$
Dubinin–Astakhov	$\frac{V}{V_o}=\exp \left[-{\left(\frac{Rt}{\beta E}\ln \left(\frac{P_o}{P}\right)\right)}^n\right]$

.

In these equations, $V$ is the volume of adsorbed gas, $V_o$ is the micropore volume, β is the gas affinity coefficient, T is the temperature, $R$ is the gas constant, and $E$ is the characteristic pressure. $P$ and $P_o$ are the equilibrium and saturation pressures, respectively. The reason why D-A is considered a generalized version of the D-R equation is because the exponent of the equation can change and be optimized to improve accuracy. Some works show that these expressions can be used to obtain a micropore size distribution [107]. The theory and principles of these methods have been extensively discussed in the literature [58,60,87,108]. Therefore, additional details will not be discussed in this work.

For shale micropore description, the following references should be consulted [109,110,111,112,113]. These techniques have been used for the determination of micropore volume and to obtain a pore-size distribution of the micropores [114].

○

Advantages and Disadvantages.

The main advantage provided by the classical inversion methods is their relative simplicity and the reduced computational work required, considering they are based on simpler thermodynamics, the law of gases and traditional concepts of phase transition and capillary condensation. However, the strength of these methods is related to their disadvantages. Although the methods are capable of describing different pore shapes, they can be highly affected by heterogeneities [115,116]. In addition, traditional methods do not provide an accurate description of the energetic interactions between the pore walls and the walls with the fluid [116].

4.2. Advanced Methods

Among the methods based on Statistical Thermodynamics, some stand out, such as Molecular Dynamic Simulation, Grand Canonical Monte Carlo Simulation, and Density Functional Theory in all its forms (DFT, NLDFT and QSDFT). These methods focus on determining the distribution of gas molecules in a pore or system that is at equilibrium. These methods are known for providing a better representation of a porous medium. It is important that most of the related works that imply the use of these methodologies are based on the Leonard–Jones 12–16 potential equation.

$$U_{ff}=4{\epsilon }_{ff}\left[{\left(\frac{{\sigma }_{ff}}{r}\right)}^{12}-{\left(\frac{{\sigma }_{ff}}{r}\right)}^6\right]$$

(6)

In this equation, ${\sigma }_{ff}$ and ${\epsilon }_{ff}$ represent the collision diameter and the well depth of the interaction potential, respectively. The collision can be defined as the distance where the energy of interaction is zero and the well depth represent how strong is the attraction between the particles. Therefore, the expression can be divided is two parts where the first can be defined as the repulsion term and the second as the attraction term [117].

Molecular Dynamic simulations consist of tracking each gas molecule at different time intervals and calculating the velocity and the position of the gas molecules. The method considers the interaction not just between gas molecules but also with the simulated solid surface. The accumulation of gas molecules near the solid surface makes it possible to determine the amount of gas adsorbed. The calculations behind it are relatively simple as long as the force field and the interaction models are available for the system of interest. However, the computational work required to model a system containing a large number of gas molecules as well as long times makes the application of this method time-consuming and less efficient. Despite its disadvantages, the literature shows multiple examples of the application of Molecular Dynamics for shale pore structure characterization, e.g., [71]. The method has been used mainly to understand the behavior and the diffusion of different gasses in porous media that could impact the adsorption capacity or related flow properties [47,118,119,120,121,122].

The Monte Carlo method requires less computationally intensive work compared to MD simulation. It consists of assessing a random initial configuration, in this case to gas molecules or particles, and then calculating the equilibrium distribution after a random movement of a selected particle. The equilibrium of the system is analyzed through the calculation of the total energy of the system after the random movement. The process is repeated on each particle until the energy of the system is reduced to a minimum. In case a particle displacement results in an increment of the total energy, the probability of that event is calculated to penalize it, which serves as a baseline to avoid it during the next simulation step. Similarly to MD simulations, the literature contains evidence of the use of Monte Carlo Simulation for shale pore structure characterization [123]. However, in a similar way, the preferred use is to understand the adsorption mechanism of gases like methane and carbon dioxide in shales [81,119,124,125,126,127,128].

The most used advance inversion method for pore structure characterization of shales is DFT and its derived methods (NLDFT and QSDFT). Its popularity has increased in the last few years, due to the relatively low computational requirement and good accuracy in comparison with MD and MC methods. In relation to the use of traditional method as BJH, Figure 8 shows how the use of BJH and DFT has changed in the last 10 years. It is important to clarify that not every work presented in the literature has been accounted for. However, we expect that the selected work to be representative of the general trend. As can be noticed, the use of DFT has increased in popularity, without being able to completely replace the use of traditional BJH.

Perhaps the most widely used DFT model for pore structure characterization is that developed by Tarzona and Tarzona et al. back in 1985 and 1987 [129,130,131]. However, in shale related works, the most common references include the works of Lastoskie et al., 1993 [132], Oliver and Occelli 2001 [133], Do and Do 2003 [117], Ravikovitch and Neimark 2006 [134], Landers et al. 2013 [135], and Zhang 2013 [136]. All these works describe the theory and calculations behind the process used for the development of a DFT model.

The basis of a DFT model is the calculation of the Grand Potential energy that is described as a functional of the density distribution of a single particle. To provide a better description, it should be said that a pore of a specific geometry (cylindrical, slit, etc.) is introduced in a homogeneous fluid. Once immersed, the fluid interacts with the pore wall until the system can reach equilibrium at a given chemical potential. Equation (7) describes the energy of the system.

$$\mathsf{\Omega }=F\left(\rho \left(\mathit{r}\right)\right)+\int \rho \left(\mathit{r}\right)u_{ex}\left(r\right)dr-\int \rho \left(\mathit{r}\right)\mu dr$$

(7)

where the terms can be separated from left to right and presented as the Helmholtz free energy ($F$), the contribution of the external force and the contribution of the chemical potential ($\mu$). In addition, $\rho$ is the fluid density and $\mathit{r}$ is the position vector inside the pore [117,135]. Equation (7) is the basis as well for the NLDFT theory, which in contrast to the local DFT allows to analyze the short range correlations near the wall by the incorporation of smoothing approximation of the density of the fluid as well as the weighing functions [58]. One of the most significant drawbacks is the assumption of a smooth and homogeneous surface of the pores. To overcome these inconveniences, more complex models have been developed (QSDFT), where the grand potential of the fluid and the solid are considered. The solid density included in the QSDFT can consider heterogeneities of the pore wall. The change with respect to the previous mathematical expression is the inclusion of the potentials and densities of the solid. A complete description of the equation can be found in Landers et al., 2013 [135].

After the description of all these methods, we should say that any of them can be applied to obtain the pore size distribution of a porous media directly. In other words, all these methods are addressed to elaborate theoretical adsorption curves of a solid with a specific pore geometry. These theoretical isotherms allow the elaboration of kernels that can be translated into the pore-size distribution through the application of the generalized adsorption isotherm equation

$$N_{\mathrm{ex}}\left(\frac{p}{p_o}\right)={\int }_{w_l}^{w_u}{N}_t\left(\frac{p}{p_o},w\right)f\left(w\right)dw$$

(8)

where $N_{\mathrm{ex}}$ corresponds to the experimental isotherm values, $w_u$ and $w_l$ are the maximum and minimum size of the pore width (w), $N_t$ correspond to the theoretical isotherm obtained from the DFT simulation. The method assumes that the isotherm from the experiment is made up from the combination of individual isotherms from single pores.

To summarize the process of the acquisition of the pore-size distribution through the application of any of the described advanced methods, Figure 9 shows a modified version of the diagram presented by Rouquerol et al., 2014.

○

Advantages and Disadvantages

The clear advantage of the advanced techniques is their ability to provide a better description of the fluid’s behavior in porous media. If we have an ordered, known, and homogeneous porous solid, it should be possible to obtain a representative pore-size distribution. The most significant drawback for shale application is related to the availability of kernels as well as the complexity provided by the shale itself. Shale is characterized by having pores formed by different minerals and of different shapes that make it almost impossible to have an accurate kernel. However, the available kernels have been shown to provide an acceptable description of the porosity of shales. In addition, the computational work required to obtain a new kernel could still be a challenge.

4.3. Other Methods

The literature shows that other methods can be developed to draw information on the porous media of shales. Among these techniques, we find that the Horvath–Kawazoe method [137] allows users to obtain the distribution of the micropores for slit shape pores [138,139,140]. The method was later modified to work on cylindrical pores as well as on spherical pores [141,142]. This method has been used for shale pore structure analysis. To obtain a description of the external surface area and to calculate the micropore volume of materials that are mainly formed by micropores, the method known as the deBoer t-plot [67,143] appears to be an acceptable option and has been used in shales [123,144]. Finally, one of the methods that appears to be gaining popularity is the application of the Frenkel-Halsey–Hill (FHH) theory for the study of the fractal dimensions of shales [113,145,146].

The analysis of the fractal dimension of a porous material is related to the characterization of the roughness of the pore surface, and it is expressed by an dimensionless number often between two and three, where two is the representation of a perfectly smooth surface and three correspond to rough and disordered surface [147,148]. In the case of the FHH model, the fractal dimension is described by the number of gas molecules (spheres) that cover the solid surface and form a monolayer [148]. The literature shows that, in addition to roughness information, the analysis of fractal dimensions of shales can provide information regarding pore connectivity [149]. However, this field is still to be explored, especially due its potential to describe complex porous structures commonly present in shales.

5. Gas Adsorption for Shale Characterization

We have discussed the background on the use of gas adsorption for shale pore structure characterization, but not the experimental procedure itself. There are several steps that are involved in the execution of an experiment. For example, the sample preparation that includes the definition of particle size and amount of sample, the degas conditions and the selection of the adsorbate. However, as we have previously stated in this review, there is no agreement on a specific methodology used for this purpose. For this reason, we have decided that is worth providing a diagram summarizing the steps for an adsorption experiment addressed to determine PSD. The information in the diagram, shown in Figure 10, was prepared from the compilation of the different works collected for the preparation of this review. It is worth mentioning that our purpose is to provide a simple guide and that by no means we are defining a standard methodology.

The information presented in Figure 10 corresponds to a simplified procedure for shale pore-structure characterization applying gas adsorption. However, changes can occur in the procedure, depending on the purpose of the author. The selection of the gas depends on which fraction of the porosity needs to be analyzed. If the purpose is to analyze mesopores and part of the macropores, then nitrogen and argon are the best options. If the objective is the micropore region, CO₂ appears to be the most suitable adsorbate. A combination of N₂ and CO₂ has been the most frequently used approach to gather a full spectrum of the pore structure, covering micro and mesopores. Regarding the particle size, the literature shows that it can go from a few microns to the use of full-size core-plugs, and the particle size closely relates to the amount of mass used. The number of isotherm points or pressure steps should be optimized for the sample being analyzed. More points can translate into a more detailed pore-size distribution. However, this can also imply longer experimental times, especially for those samples with large particle sizes. Finally, the selection of the inversion method should be discretized by the researcher. We consider that if additional information about the samples is available, e.g., imaging, mineral composition, and a suitable kernel that matches the additional information of the rock is available, DFT model should be selected to obtain a better description of the system. On the other hand, if not additional detail or appropriate kernels are available, traditional methods, e.g., BJH, D-A, should be considered [150]. Although these methods are known for having a questionable accuracy, the results can still be acceptable, especially if the purpose is to analyze changes over the same sample.

6. Comparison with Other Techniques

Several works can be found in the literature in which authors have compared the performance of gas adsorption with respect to other techniques. For example, it has been compared with radiation techniques, e.g., SEM, TEM, SANS, USANS, and Mercury Intrusion Porosimetry (MIP). All these comparisons have proven that even using classical inversion methods can still yield an acceptable description of the Pore-Size Distribution (less than 300 nm). An example of that is the work presented by Clarkson et al., 2012 and 2013 [31,151], where the authors were able to show that the result of simply using BJH and D-A can produce acceptable results in comparison with the results obtained from SANS and USANS. With respect to MIP, the literature shows that the obtained pore-size distribution differs between the two techniques [152]. However, instead of being a suitable replacement, the literature suggests that both techniques can be complemented to obtain a description of the pore structure that can go from less than 1 nm (if CO₂ is used) to several microns.

7. Discussion of Gas Adsorption Limitations

Gas Adsorption may not be the best technique to obtain an accurate pore-size distribution of shales. However, its efficacy and low-cost place it as a good option for researchers. The drawbacks and limitations of techniques are well known and can be easily found in the literature. To provide a brief description of the limitations of the technique, we can mention: the lack of a complete understanding on capillary condensation and gas adsorption can be affected by the composition of the rock, where clays and the organic matter appear to induce the most significant alteration. Some of inversion methods have shown to overestimate or underestimate the total porosity of the shale with respect to other techniques. The isotherm can be highly affected by the particle size, if smaller particles are used, so a fraction of the macropores is lost. The BJH from desorption should not be used for characterization of pores smaller than four (4) nm, due to the presence of the tensile strength effect. No appropriate kernels are available to provide an accurate description of the shale matrix [89,150,152,155]. In addition to the previously mentioned limitations, we should add the fact that most of the shale pore-size distribution data have been obtained using the incorporated software on the instrument used, with Micromeritics and Quantachrome instruments being of popular use. Therefore, the user should carefully consider the theory behind the software and in case that DFT is used, the use should analyze the kernels provided in order to make an appropriate selection. The main problem related to the use of pre-existing models is related to the geometry (or shape) of the pores. All inversion methods, classic or advanced, assume a specific geometry that is used to calculate the dimensions of the pores. These geometries include the assumption of cylindrical, spherical, ink-bottle-shaped and slit-shaped pores. However, it is possible to say that any of these basic geometries is able to accurate describe the shale porosity that is known for its large heterogeneity. Shale porosity often include the presence of angular pores and others of different shapes across all scales, where this diversity has a significant impact on the understanding of the physical interaction between the pore walls and the gas and in some cases the interaction between the pore walls. Therefore, it is important to address efforts in creating models that consider more complex geometries, where advanced methods (MD, GCMC, DFT) appear to be the most suitable models to address this issue. For this purpose, the combination with complementary techniques, i.e., imaging, represent a valuable input by providing a more detailed description of the pore structure that could be used to feed the models.

8. Conclusions

The limitations and drawbacks of gas adsorption for pore architecture characterization are known and understood. However, its easy accessibility and low application cost make it worth to invest in its performance improvement. In the meantime, the available models satisfy the necessity of gathering a description of the shale pore architecture, if complemented by other techniques. With this in mind, we feel confident to express three main points regarding the use of gas adsorption for shale characterization.

• Gas adsorption can be considered the fastest and most inexpensive approach to obtain a description of shale pore structures.
• Gas adsorption is not a standalone technique and should be complemented by other techniques such as imaging.
• The most suitable scenario for gas adsorption is to compare changes over the same sample, after it has been subject to an alteration.

The future of gas adsorption relies not just on the development of more complex inversion models to describe the pore architecture but also on its combination with complementary techniques. We define this process as the hybridization of gas adsorption. In previous research, we demonstrated that combination of gas adsorption and NMR can be used draw information regarding surface alteration [43]. Others have shown that these two techniques can be combined to understand the distribution of gasses in the pores and to gather information related to diffusion [153,154,155]. In addition, we are also confident that the hybridization of these two techniques can provide insights to better understand the capillary condensation process in shale complex porous media. The hybridization potential is not limited to the combination with NMR, so other possibilities exist such as combinations with imaging techniques, e.g., tomography. Even within gas adsorption, it is possible to combine models to obtain more details about the pore architecture. For example, the understanding of the fractal dimensions of shales can help in the selection of an inversion method for pore-size distribution or could as well help in the development of new models that better consider the heterogeneities of a complex matrix.