Open Access
Issue
Oil Gas Sci. Technol. – Rev. IFP Energies nouvelles
Volume 76, 2021
Article Number 32
Number of page(s) 13
DOI https://doi.org/10.2516/ogst/2021013
Published online 11 May 2021

© C.D. Tsakiroglou et al., published by IFP Energies nouvelles, 2021

Licence Creative CommonsThis is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1 Introduction

Tight and shale reservoir systems pose a tremendous potential resource for future development, and the determination of their microscopic and macroscopic properties has long attracted the attention of scientific society. In particular, shale gas reservoirs possess many so-called “unconventional” features, the traditional approaches are unable to determine the reservoir permeability, and both the micro-scale analysis of flow regimes and macroscopic estimation of extremely low permeability pose emerging challenges (Hill and Nelson, 2000; Javadpour et al., 2007).

Depending on the pore size range, and prevailing conditions, the flow of natural gas through conventional porous rocks, tight formations and shale systems, is controlled by a variety of flow regimes, including Knudsen, transition, slip and viscous Darcy flow (Ziarani and Aguilera, 2012). Having known and understood the rock properties, including pore geometry, flow/transport parameters, the determination of hydrocarbon storage and recovery from shales becomes tractable. The microstructural characterization of shale samples could be achieved by combining Mercury Intrusion Porosimetry (MIP), with low‐field Nuclear Magnetic Resonance (NMR), N2 adsorption/desorption, and high resolution Focused Ion Beam–Scanning Electron Microscopy (FIB–SEM) images (Al Hinai and Rezaee, 2015).

The gas permeability of shales can be correlated with the capillary pressure data from mercury injection measurements, by using theoretical and empirical equations (Al Hinai and Rezaee, 2015). On the other hand, the pore network approach is a practical way to explain the macro flow behavior of porous media from a microscopic point of view, so that the shale matrix is simulated by 3-dimensional pore network models that include typical bimodal pore size distribution, anisotropy and low connectivity (Zhang et al., 2017). Moreover, the digital shale rock could be developed by 3-Dimensional (3D) images obtained from micro/nano Computed Tomography (CT) and FIB-SEM images reconstructed from 2D SEM images of the pore structure, and using the Lattice Boltzmann method to compute the intrinsic permeability, porosity and tortuosity; these parameters are used to calculate the apparent permeability by accounting for different gas transport mechanisms (Sun et al., 2017).

Three-dimensional nano-porous structures of shales reconstructed from SEM images of shale samples enabled the computational analysis of structural characteristics and calculation of tortuosity, effective Knudsen diffusivity and intrinsic permeability by Lattice-Boltzman method. It was found that Knudsen diffusion is always of key importance on shale gas transport mechanisms in the reconstructed shales, while no Darcy (viscous) flow regime was observed (Chen et al., 2015). The FIB-SEM technique was employed to observe and characterize the morphology of the nanopores in shales, Lattice Boltzmann Method (LBM) was used to simulate high-Knudsen gas flows, and the nano-scale pore morphology was found to play an important role on the gas transport properties of shale matrix (Zheng et al., 2017). A method was developed to identify organic and inorganic pores from SEM images of shale samples collected from the wells of a productive formation, and it was found that: (i) most of inorganic pores are mesopores; (ii) the sizes for organic pores tend to be lognormal distributions with peak in the range 10–30 nm; (iii) the gas storage space obtained by SEM is some orders of magnitude larger than that obtained by gas adsorption method; (iv) the permeability is strongly anisotropic, and a significant non-Darcy effect is evident (Cao et al., 2019). The image-based fractal characteristics of shale pore structure were examined from 100 SEM images, at multiscale resolutions ranging from 15 to 420 nm, their effect on the accurate prediction of gas permeability was analyzed, and it was found that the box counting method enables the accurate computation of fractal dimension, compared to Sierpinski carpets analytical solution, and the gas permeabilities approach to the exact values when the resolution is better than 50 nm (Song et al., 2019). Four major techniques, including N2 adsorption, SEM, MIP, and NMR, were employed to characterize the Pore Size Distributions (PSDs) of a shale oil reservoir, and it was found that there is a wide range of pore sizes ranging from micropores to macropores (>1 μm), interconnected by small-scale throats, while the comprehensive characterization of shale oil reservoir pore structures requires the combination of more than one techniques, like NMR and MIP (Zhang et al., 2019). A sophisticated algorithm of pore-and-throat network reconstruction from 3D SEM images with resolution of 4 nm was developed by accounting for complex geometry and extremely small sizes of the pore-bodies and throats, and validated for a Marcellus shale core sample (Zheng et al., 2019).

In the present work, with the aid of percolation theory, analytic models of the N2 adsorption–desorption isotherms and Hg intrusion curve are formulated as functions of the pore and throat radius distributions and pore network connectivity (Tsakiroglou et al., 2004). Then, a stepwise methodology is adopted to estimate the pore network parameters with inverse modeling of N2 sorption and Hg intrusion curves and applied to datasets of several shale samples. The critical path analysis of percolation theory is used to calculate the liquid permeability and Knudsen diffusion coefficient of porous media explicitly from the aforementioned pore structure parameters and compare the predictions with experimental measurements.

2 Materials and methodology

2.1 Shales and transport properties

The Perth Basin is a 100 000 km2 area covering the Western Australian margin between Augusta and Geraldton. Samples were collected at various depths from the same well of Carynginia shale gas formation (Al Hinai and Rezaee, 2015). The N2 adsorption/desorption isotherms were measured with Micromeritics TriStar II 3020 (pore sizes probed: 1.3–200 nm), whereas the Hg intrusion curve was measured with Micrometrics Autopore IV Porosimeter (pressure range: 1–60 000 psi; pore sizes probed: 400 μm–3.5 nm). The liquid permeability and all pertinent Klinkenberg parameters were measured with a modified steady-state method, where the sample was embedded in a resin disc, the gas flow rate was measured at the outlet, the average gas pressure was in the range of 1–7 bar, and the net confining pressure was 70 bar (Al Hinai et al., 2013). The properties of samples, measured by the various techniques, are shown in Table 1.

Table 1

Properties of shale samples.

2.2 Pore-scale models of drainage and imbibition

The actual pore space of sedimentary rocks, including shales (Fig. 1a) is complex, and it is commonly characterized in terms of pore-and-throat networks (Fig. 1b). In the present work, the porosity of shales is regarded as a stochastic network of pores (sites)-and-throats (bonds), the sizes of which are sampled from bimodal site f s (r) and bond f b (r) radius distributions, both composed of log-normal component functions, namely:

(1a)

(1b)where <r s1>, σ s1, <r s2>, σ s2 are the mean radii and standard deviations of the two component site distributions with contribution fractions to the overall one c s and 1 − c s , respectively, and <r b1>, σ b1, <r b2>, σ b2 are the mean radii and standard deviations of the two component bond distributions with contribution fraction to the overall one c b and 1 − c b , respectively.

thumbnail Fig. 1

(a) SEM image of AC1. (b) Representation of the pore structure as a pore-and-throat network.

The throats (bonds) are considered as narrow constrictions without volume, while the pores are assigned the entire pore volume according to the relation:

(2)where β s is a volume shape factor (β s  ≥ 0). The parameter value β s  = 2.0, implies a network of prismatic pores with their volume proportional to r 2. On the other hand, β s  = 1.0, implies a network of slit-shaped pores with their volume proportional to r 1. So any value in the range 1 < β s  < 2 corresponds to a mixture of pore shapes between the foregoing limits.

The representation of the pore shapes in reservoir rocks with angular geometries contributes to their realistic characterization in terms of pore-and-throat networks (Tsakiroglou and Payatakes, 1993). On the other hand, angular geometries allow the rational computation of the oil and water flow rates, and wettability effects on simulated relative permeability curves (Firincioglu et al., 1999; Tsakiroglou, 2014; Valvante et al., 2005; Zhou et al., 1997). The cross-section of each pore is modeled by a high-order regular polygon with the number of sides (angles), area, and perimeter equal to n s , A, and P, respectively (Fig. 2a). The pore shape factor, G, is given (Tsakiroglou et al., 2004) by,

(3)and the pore radius, r, coincides to the radius of the inscribed circle (Fig. 2a), given by,

(4)

thumbnail Fig. 2

(a) Polygonal pore; (b) adsorbed layer along with condensed N2 along a triangular pore; (c) capillary equilibrium at a pore corner.

2.3 Models of N2 adsorption/desorption in pore-and-throat networks

The graual saturation of pores with liquid nitrogen takes place through two parallel mechanisms (Fig. 2b): (1) gas adsorption in multiple layers across the pore-walls, and (2) capillary condensation of liquid at the corners of adsorbed layer. The thickness of the adsorbed layer, t c , is governed by the generalized Frenkel–Halsey–Hill (FHH) equation (Gregg and Sing, 1982),

(5)where x is the relative vapor pressure, σ is the diameter of N2 molecule (σ = 3.54 A), the parameter b depends on the properties of adsorbent and adsorptive, and s is an index with values ranging from 2 to 3. The radius of curvature of the condensate/vapor interfaces, r c , is given by Kelvin equation (Gregg and Sing, 1982)

(6)

To get a visual view of the calculated values and fluid saturation, we first assume that β s  = 2.0 (prismatic pores) and then we generalize the so-produced equations to any β s value. The liquid N2-saturation in a pore includes the adsorbed and condensed phases and is given by,

(7)where F LG is the fraction of the pore cross-section area occupied by liquid N2, and depending on the liquid/gas/solid contact angle, θ LG (Fig. 2c; we assume that θ LG = 0°, for liquid N2, unless otherwise stated). If equation (7) is generalized for any β s value, then it is written as,

(8)

With regard to the fluid saturation history, the N2 capillary condensation resembles the imbibition of a wetting fluid (liquid N2), and the critical radius of curvature for condensation, r c  = r con, is given by,

(9a)

(9b)

Equations (9a) and (9b), in conjunction with equations (5) and (6), yield the critical radius r = r s , given by,

(10a)

(10b)

The saturation of a pore network with liquid N2 can be expressed as a function of the relative vapor pressure by the relation,

(11)

The capillary evaporation of liquid N2 from a pore resembles the drainage of the wetting fluid (liquid N2), and hence the critical radius of curvature for evaporation, r evp , is given by,

(12a)

(12b)which in combination with equations (5) and (6) lead to,

(13a)

(13b)where r b is the critical radius of throats which are “allowable” to N2 evaporation (r ≥ r b ). N2 empties a pore if its radius is less than the critical size (r b or r s ) and is accessible either to the bulk vapor phase (primary desorption-drainage) or to isolated vapor pockets (secondary desorption-drainage). During the N2 desorption process, the fraction of bonds (throats) q b and sites (pores) q s which are “allowable” to the vapor phase are given by,

(14)respectively. In percolation theory, the accessibility function depends on pore network, topology, size, and statistical spatial pore size correlations, and is defined as the fraction of allowable to a fluid pores that are accessible to this fluid (Tsakiroglou et al., 2004, 2009). In the specific case of N2 desorption, the accessibility function relates the fraction of pores which are accessible to vapor N2 with the fraction of throats that are allowable to vapor N2. Finally, the liquid N2 saturation during desorption is given by (Tsakiroglou et al., 2004),

(15)

2.4 Model of Hg intrusion in pore-and-throat network

In general, the value of pore volume calculated from N2 sorption does not coincide to the corresponding value calculated from Hg intrusion curve (Tab. 1). Although Hg is a non-wetting fluid for the majority of materials, by “convention” the Hg contact angle, θ Hg, is measured with respect to the wetting fluid (low pressure air) and is considered always lower than 90°. In order to estimate a unique pore volume value, V p (cm3/g), the following method is adopted. For a given Hg intrusion contact angle, the maximum intrusion capillary pressure is converted to an equivalent relative vapor pressure. This value is interpolated in N2 desorption curve to calculate the corresponding liquid N2 volume. The sum of the two volumes yields the total pore volume. The procedure is repeated for various values of θ Hg, and the produced datasets (V p vs. θ Hg) are fitted to a non-linear parametric function of the form,

(16)so that the values of A 1, A 2, θ 0, w are estimated.

The Hg intrusion in a pore-and-throat network is analogous to capillary evaporation for q si  = 0 (primary desorption) and is controlled by throat sizes. The mercury saturation in a pore network is given (Tsakiroglou et al., 2004, 2009) by,

(17)where,

(18)

(19a)

(19b)

2.5 Accessibility functions

The accessibility functions can be described satisfactorily by the following functional form,

(20)where is the initial fraction of sites occupied by the non-wetting fluid; for Hg intrusion (primary drainage, i = 0) and for nitrogen desorption (secondary drainage initiating before the completion of adsorption, i = 1). Evidently, two different sets of parameters (a 0, b 0) and (a 1, b 1) are required to quantify the accessibility functions Y s0(q b ) and Y s1(q b ) that specify the Hg intrusion and N2 desorption in pore networks, respectively. A “physical” parameter characterizing the network accessibility function is the “percolation threshold”, q b  = q bci , which corresponds to the inflection point of the curve and is obtained from the numerical solution of the following equation,

(21)

2.6 Critical Path Analysis (CPA)

The transport properties of porous materials, characterized in the terms of a parametric pore-and-throat network, with pore sizes spanning some orders of magnitude, can be computed either explicitly with simulation of the corresponding process (e.g. single phase flow) in a reconstructed network (Jones et al., 2018), or implicitly with the use of approximate methods, such as the CPA of percolation theory (Tsakiroglou and Ioannidis, 2008; Tsakiroglou et al., 2009). The main idea hidden behind the CPA is that for a pore system with a broad pore size distribution most of the flow is transferred by a network spanning cluster consisting of pores with conductance greater than a critical value, g c (Katz and Thompson, 1986). This critical conductance g c depends on pore space characteristics (e.g. pore shape, pore size distribution, connectivity, dimensionality, etc), and CPA enables the calculation of the corresponding critical pore dimensions. The transport properties of the critical path may be governed by universal scaling laws allowing the approximate computation of accessibility function at the vicinity of percolation threshold (Tsakiroglou et al., 2009). In mathematical terms, the viscous flow or Knudsen flow problem reduces to the calculation of the critical pore dimensions that maximize the total conductance of the critical path (Hunt, 2001; Tsakiroglou and Ioannidis, 2008; Tsakiroglou et al., 2009). In general, the dependence of the local viscous flow conductance on pore size may be different from that of the local Knudsen diffusion conductance, so that the critical pore dimension may depend on the transport mechanism (Tsakiroglou et al., 2009). The gas molar flux in a porous medium for the transition region between Knudsen diffusion and viscous flow is given by the sum of the individual contributions, namely,

(22)where k L is the liquid permeability, and D kn is the Knudsen diffusion coefficient of the critical path. According to the CPA (Tsakiroglou et al., 2009), these properties are approximated by the relationships,

(23)

(24)where φ is the porosity, t is the universal scaling exponent of conductivity (t = 2.0 for 3-D systems), and mean pore radius, <r sq > can be expressed as function of throat radius, r, by,

(25)

The values of k L and D kn are maximized for r = r max,f and r = r max,d , respectively, and are considered approximately equal to the permeability and Knudsen diffusion coefficient of the investigated porous medium.

For low permeable porous media, the measured gas permeability is pressure-dependent and is commonly expressed by the relationship,

(26)where b is the Klinkenberg factor, and P m is the mean gas pressure (Al Hinai et al., 2013). For the steady-state gas flow through a porous medium, the integrated form of Darcy equation, including the Klinkenberg effect, takes the form,

(27)where u 2 , P 2 are the flow velocity and pressure at the outlet, L is the sample length, μ g is the gas viscosity, and ΔP is the pressure drop across the sample. By combining equation (22) with equation (27), it is obtained,

(28)

2.7 Demonstration of methodology

For the non-linear parameter estimation, numerical codes in the software platform of Athena Visual Studio (Stewart and Caracotsios, 2008) were developed, and a step-wise procedure was adopted (Fig. 3). First, the parameters n s , β s , and f s (r; <r s1>, σ s1, <r s2>, σ s2, c s ) were estimated with inverse modeling of N2 adsorption isotherm. Then, the foregoing parameters were fixed, and f b (r; <r b1>, σ b1, <r b2>, σ b2, c b ) along with Y si (q b ) were estimated with inverse modeling of N2 desorption isotherm. Then, Y s0 (q b ) was estimated from inverse modeling of the Hg intrusion curve. Finally, the throat size distribution, secondary and primary pore accessibility functions, and Hg contact angle were re-estimated, with the simultaneous inverse modeling of N2 adsorption/desorption isotherms and Hg intrusion curve.

thumbnail Fig. 3

Multistep methodology of pore structure characterization.

3 Results and discussion

3.1 Inverse modeling of Hg intrusion curve and N2 sorption isotherms

It is well-known that for low porosity samples, surface irregularities (e.g. nooks, crannies, crevices) and inter-granular pores that are filled with Hg over the low pressure range are interpreted as macro-pores. It is of key importance to exclude such characteristics from the deconvolution of Hg intrusion and N2 isotherms, and pay attention mostly on mesopores and micropores. In this context, the Hg intrusion curves of all samples were corrected by removing the part of the curve extending over low pressures (Fig. 4).

thumbnail Fig. 4

Corrected Hg intrusion curves for (a) AC1, (b) AC5, (c) AC8.

For the three shale samples (Tab. 1), the Hg intrusion data were combined with N2 desorption isotherms to calculate the total pore volume as a function of the Hg contact angle, and the results were fitted with equation (16) to estimate all pertinent parameters (Tab. 2, Fig. 5). The V p (θ Hg) relationship was incorporated into the set of equations modeling Hg intrusion and N2 sorption processes.

thumbnail Fig. 5

Correlation of the total pore volume with Hg contact angle by matching Hg intrusion with N2 desorption datasets (Tab. 2).

Table 2

Estimated parameters of equation (16).

The application of the methodology of pore structure characterization from N2 adsorption/desorption isotherms and Hg intrusion curves was applied to the datasets of the three shale samples (Tab. 1), and the results are shown in Table 3.

Table 3

Estimated parameter values for shale samples AC1, AC5, AC8.

The measured N2 isotherms are predicted satisfactorily by simulated ones for all cases (Figs. 6a, 7a, 8a and 8b), indicating that the estimated statistics of pore and throat sizes over the low size range (<0.5 μm) are consistent with the estimated accessibility function. On the other hand, regarding the Hg intrusion curves, there is a discrepancy between the measured and predicted curve over the low and very high pressure ranges (Figs. 6b, 7b and 8c). It seems that, in spite of the large number of “unknown” parameters and the iterative method of inverse modeling, we are unable to track precisely the N2 adsorption/desorption curves, and simultaneously the Hg intrusion curve. The initial guesses of TSD, PSD, and ACF are obtained from N2 adsorption (PSD) and desorption (TSD and ACF) curves, and unavoidably the finally estimated PSD and TSD might have been biased toward the low radii sizes. On the other hand, large throat radii that are reflected in the low pressure range of Hg intrusion curve (Figs. 6b, 7b and 8c) are overlooked since they are not “shown” at all in the N2 desorption curve. Other factors that prevent the simultaneous prediction of all datasets and favor the shift of PSD and TSD toward small sizes are: (i) the higher contribution fraction of N2 sorption to the estimation process (2 datasets against 1 dataset); (ii) the uncertainty of the “fixed values” of parameters included in N2 adsorption/desorption equations (e.g. liquid N2 contact angle); (iii) the lack of sufficient information concerning the large pore sizes that are not detected at all by N2 adsorption isotherm, and are shadowed in Hg intrusion curve; such large pore sizes might be correlated with the Hg retraction curve (Tsakiroglou et al. 2009) which was overlooked in the present approach.

thumbnail Fig. 6

Measured vs. predicted (a) N2 sorption isotherms, (b) Hg intrusion curve for sample AC1.

thumbnail Fig. 7

Measured vs. predicted (a) N2 sorption isotherms, (b) Hg intrusion curve for sample AC8.

thumbnail Fig. 8

Measured vs. predicted (a) N2 sorption isotherms (simulated structure S1), (b) N2 sorption isotherms (simulated structure S2), (c) Hg intrusion curve (simulated structures S1 & S2) for sample AC5.

In the case of AC5, the significant discrepancy observed between observed and simulated Hg intrusion curve (Fig. 8c) led us to repeat the last step of parameter estimation [f b (r), Y si (q b ), Y s0(q b ), θ Hg] by modifying the initial guesses. Finally, the estimated TSD widens toward large radii, and ACFs shift to higher qb values, characterized by higher percolation thresholds. Not having imposed any limitation on the ratio of pore to throat sizes, the overlapping between PSD and TSD may be discernible (Figs. 9a9c), with the mean value of estimated TSD exceeding that of estimated PSD (Tab. 3, AC5, S1/S2). The significant overlapping of TSD and PSD (Fig. 9) leads to comparable values of q bc0 and q bc1 (for AC1 and AC8), while in the case of AC5, with the throat sizes being greater than the pore sizes, q bc0 is much less than q bc1, which is opposite to conventional predictions of percolation-type approaches. The physical interpretation of such a quantitative description is that the critical capillary sizes for Hg penetrating in throats are comparable to those of Hg penetration in interconnecting pores. Probably, a simpler model like a capillary channel network, described by only one pore size distribution and one accessibility function, might be more convenient than the present model for the representation of the pore space.

thumbnail Fig. 9

The pore- and throat-radius distributions obtained from inverse modeling for (a) AC1, (b) AC8, (c) AC5.

3.2 Prediction of transport properties

The parameter values of simulated throat-and-pore networks (Tab. 3, Fig. 9) were used to compute k L and D kn as functions of throat radius, r, through the approximate equations (23) and (24), respectively (Figs. A1A3). The maximum k L,max or D kn,max occur at the minimum throat radius, over which the interconnected critical path of throats and pores transfers most of the viscous (Figs. A1a, A2a and A3a) or Knudsen (Figs. A1b, A2b and A3b) flow, respectively. The so-predicted permeability and Knudsen diffusion coefficient are compared with corresponding experimental values (Tab. 4). It seems that D kn is predicted satisfactorily for almost all cases (Tab. 4), while the permeability is almost predictable only for sample AC5, while it fails for AC1 and AC8 (Tab. 4).

Table 4

Comparison of predicted vs. measured flow properties of shales.

The Knudsen number is commonly used to classify flow regimes in fine pores where deviation from continuum flow is important (Ziarani and Aguilera, 2012). It is defined as the ratio of molecular mean free path, λ (nm) to the mean pore radius, r s (nm), and is given by,

(29)

If K n  < 0.01, then viscous flow dominates, and the conventional fluid dynamics equations and Darcy equation are applicable. If 0.01 < K n  < 0.1, then slip flow regime occurs and Darcy equation with Klinkenberg or Knudsen modification is applicable. If 0.1 < K n  < 10, then the transition flow regime prevails, where both slip (continuum) and diffusion flows can occur, and the traditional fluid dynamics equations start to fail. The higher the Knudsen number, the higher the chance of failure. Even though conventional equations could be applied (i.e., Darcy with Knudsen correction), the validity of such a formulation is questionable. It is safer to use Knudsen diffusion equations, especially for flow at higher Knudsen numbers (close to 10). If K n  > 10, then the Knudsen (free molecular) flow regime prevails, where the continuum fluid flow equations fully break down, and diffusion-based formulations must be applied. Looking at the post-calculated Knudsen number, K n (Tab. 4), the transition (AC5) and Knudsen (AC1, AC8) flow regimes (AC1, AC8) prevail. Therefore, Knudsen flow is dominant over all cases (AC1, AC5, AC8), but only in the one case (AC5), viscous flow becomes evident. Based on the aforementioned interpretation, failing to predict the value of kL for AC1 and AC8, and succeeding to predict satisfactorily D kn for all cases seem reasonable.

4 Conclusion

Analytical models of Hg intrusion and N2 adsorption-desorption in porous materials are used for the inverse modeling of experimental datasets of shale samples and estimation a full set of parameters for the pore space, regarded as a stochastic pore-and-throat network. Using the critical path analysis of percolation theory, analytical relationships are developed for the explicit computation of the Klinkenberg parameters (liquid permability and Knudsen diffusion coefficient) of the gas permeability from pore space parameters. The methodology is evaluated with demonstration to the datasets of three shale samples:

  • The estimated N2 adsorption/desorption curves are predicted satisfactorily, while the predicted Hg intrusion curve deviates from the experimental one, over low and high pressure range.

  • The initial guesses of parameters are based on N2 sorption isotherms, so that the estimated throat and pore size distributions are biased over the low sizes.

  • The overlapping between the pore and throat size distributions is significant for the investigated shales.

  • Application of the methodology to datasets of shales indicates that, based on the estimated pore space parameters, the liquid permeability is predictable at moderate K n values, over which the transition regime is dominant, while the Knudsen diffusion coefficient is also predictable at high K n values, over which Knudsen flow is dominant.

Acknowledgments

Dr. C.D. Tsakiroglou acknowledges support of this work by the project POLITEIA II (MIS 5002478) which is implemented under the “Action for the Strategic Development on the Research and Technological Sector”, funded by the Operational Programme “Competitiveness, Entrepreneurship and Innovation” (NSRF 2014–2020) and co-financed by Greece and the European Union (European Regional Development Fund).

Appendix

thumbnail Fig. A1

(a) Liquid permeability, and (b) Knudsen diffusion coefficient of the critical path as a function of throat radius for the simulated pore structures of sample AC1.

thumbnail Fig. A2

(a) Liquid permeability, and (b) Knudsen diffusion coefficient of the critical path as a function of throat radius for the simulated pore structures of sample AC8.

thumbnail Fig. A3

(a) Liquid permeability, and (b) Knudsen diffusion coefficient of the critical path as a function of throat radius for simulated pore structures (S1, S2) of sample AC5.

References

  • Al Hinai A., Rezaee R. (2015) Pore geometry in gas shale reservoirs, in: Fundamentals of gas shale reservoirs, Rezaee R. (ed), John Wiley & Sons Inc., pp. 89–116. [Google Scholar]
  • Al Hinai A., Rezaee R., Saeedi A., Lenormand R. (2013) Permeability prediction from mercury injection capillary pressure: An example from the Perth Basin, Western Australia, APPEA J. 53, 31–36. [Google Scholar]
  • Cao G., Lin M., Ji L., Jiang W., Yang M. (2019) Characterization of pore structures and gas transport characteristics of Longmaxi shale, Fuel 258, 116146. [Google Scholar]
  • Chen L., Zhang L., Kang Q., Viswanathan H.S., Yao J., Tao W. (2015) Nanoscale simulation of shale transport properties using the lattice Boltzmann method: Permeability and diffusivity, Sci. Rep. 5, 8089. [PubMed] [Google Scholar]
  • Firincioglu T., Blunt M.J., Zhou D. (1999) Three-phase and wettability effects in triangular capillaries, Colloids Surf. 155, 259–276. [Google Scholar]
  • Gregg S.J., Sing K.S.W. (1982) Adsorption, surface area and porosity, Academic Press, London. [Google Scholar]
  • Hill D.G., Nelson C.R. (2000) Gas productive fractured shales: An overview and update, Gas TIPS 6, 3, 4–13. [Google Scholar]
  • Hunt A.G. (2001) Applications of percolation theory to porous media with distributed local conductances, Adv. Water Resour. 24, 279–307. [Google Scholar]
  • Javadpour F., Fisher D., Unsworth M. (2007) Nanoscale gas flow in shale gas sediments, J. Can. Pet. Technol. 46, 10, 55–61. [Google Scholar]
  • Jones K.L., Laudone G.M., Matthews G.P. (2018) A multi-technique experimental and modeling study of the porous structure of IG-110 and IG-430 nuclear graphite, Carbon 128, 1–11. [Google Scholar]
  • Katz A.J., Thompson A.H. (1986) Quantitative prediction of permeability in porous rock, Phys. Rev. B 34, 8179–8181. [Google Scholar]
  • Song W., Wang D., Yao J., Li Y., Sun H., Yang Y., Zhang L. (2019) Multiscale image-based fractal characteristic of shale pore structure with implication to accurate prediction of gas permeability, Fuel 241, 522–532. [Google Scholar]
  • Stewart W.E., Caracotsios M. (2008) Computer-aided modeling of reactive systems, John Wiley & Sons, Hoboken, NJ. [Google Scholar]
  • Sun H., Yao J., Cao Y.-C., Fan D.-Y., Zhang L. (2017) Characterization of gas transport behaviors in shale gas and tight gas reservoirs by digital rock analysis, Int. J. Heat Mass Transfer 104, 227–239. [CrossRef] [Google Scholar]
  • Tsakiroglou C.D. (2014) Computation of the two-phase flow properties of intermediate-wet porous media: A pore network approach, Can. J. Chem. Eng. 92, 515–523. [Google Scholar]
  • Tsakiroglou C.D., Payatakes A.C. (1993) Pore-wall roughness as a fractal surface and theoretical simulation of mercury intrusion/retraction in porous media, J. Colloid Interface Sci. 159, 287–301. [Google Scholar]
  • Tsakiroglou C.D., Ioannidis Μ.Α. (2008) Dual porosity modeling of the pore structure and transport properties of a contaminated soil, Eur. J. Soil Sci. 59, 744–761. [Google Scholar]
  • Tsakiroglou C.D., Burganos V.N., Jacobsen J. (2004) Pore structure analysis by using nitrogen sorption and mercury intrusion data, AIChE J. 50, 489–510. [Google Scholar]
  • Tsakiroglou C.D., Ioannidis M.A., Amirtharaj E., Vizika O. (2009) A new approach for the characterization of the pore structure of dual porosity rocks, Chem. Eng. Sci. 64, 847–859. [Google Scholar]
  • Valvante P.H., Piri M., Lopez X., Blunt M.J. (2005) Predictive pore-scale modeling of single and multiphase flow, Transp. Porous Media 58, 23–41. [Google Scholar]
  • Zhang P., Hu L., Meegoda J.N. (2017) Pore-scale simulation and sensitivity analysis of apparent gas permeability in shale matrix, Materials 10, 104, 1–13. [CrossRef] [Google Scholar]
  • Zhang P., Lu S., Li J. (2019) Characterization of pore size distributions of shale oil reservoirs: A case study from Dongying sag, Bohai Bay basin, China, Mar. Petrol. Geol. 100, 297–308. [Google Scholar]
  • Zheng J., Wang Z., Gong W., Ju Y., Wang M. (2017) Characterization of nanopore morphology of shale and its effects on gas permeability, J. Nat. Gas Sci. Eng. 47, 83–90. [Google Scholar]
  • Zheng D., Wang W., Reza Z. (2019) Pore-network extraction algorithm for shale accounting for geometry-effect, J. Pet. Sci. Eng. 176, 74–84. [Google Scholar]
  • Zhou D., Blunt M., Orr F.M. (1997) Hydrocarbon drainage along corners of noncircular capillaries, J. Colloids Interface Sci. 187, 11–21. [Google Scholar]
  • Ziarani A.S., Aguilera R. (2012) Knudsen’s permeability correction for tight porous media, Transp. Porous Media 91, 239–260. [Google Scholar]

All Tables

Table 1

Properties of shale samples.

Table 2

Estimated parameters of equation (16).

Table 3

Estimated parameter values for shale samples AC1, AC5, AC8.

Table 4

Comparison of predicted vs. measured flow properties of shales.

All Figures

thumbnail Fig. 1

(a) SEM image of AC1. (b) Representation of the pore structure as a pore-and-throat network.

In the text
thumbnail Fig. 2

(a) Polygonal pore; (b) adsorbed layer along with condensed N2 along a triangular pore; (c) capillary equilibrium at a pore corner.

In the text
thumbnail Fig. 3

Multistep methodology of pore structure characterization.

In the text
thumbnail Fig. 4

Corrected Hg intrusion curves for (a) AC1, (b) AC5, (c) AC8.

In the text
thumbnail Fig. 5

Correlation of the total pore volume with Hg contact angle by matching Hg intrusion with N2 desorption datasets (Tab. 2).

In the text
thumbnail Fig. 6

Measured vs. predicted (a) N2 sorption isotherms, (b) Hg intrusion curve for sample AC1.

In the text
thumbnail Fig. 7

Measured vs. predicted (a) N2 sorption isotherms, (b) Hg intrusion curve for sample AC8.

In the text
thumbnail Fig. 8

Measured vs. predicted (a) N2 sorption isotherms (simulated structure S1), (b) N2 sorption isotherms (simulated structure S2), (c) Hg intrusion curve (simulated structures S1 & S2) for sample AC5.

In the text
thumbnail Fig. 9

The pore- and throat-radius distributions obtained from inverse modeling for (a) AC1, (b) AC8, (c) AC5.

In the text
thumbnail Fig. A1

(a) Liquid permeability, and (b) Knudsen diffusion coefficient of the critical path as a function of throat radius for the simulated pore structures of sample AC1.

In the text
thumbnail Fig. A2

(a) Liquid permeability, and (b) Knudsen diffusion coefficient of the critical path as a function of throat radius for the simulated pore structures of sample AC8.

In the text
thumbnail Fig. A3

(a) Liquid permeability, and (b) Knudsen diffusion coefficient of the critical path as a function of throat radius for simulated pore structures (S1, S2) of sample AC5.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.