1 Introduction

Gas transport in tight rocks is of relevance in various geotechnical contexts such as carbon dioxide sequestration, nuclear waste disposal and the exploitation of shale gas and coalbed methane (Civan 2010a). Associated transport processes are typically controlled by different mechanisms. In dry porous systems, a major feature is the superposition of fluid-dynamic (i.e., slip flow) and rock mechanical stress (pore compressibility) effects (Fink et al. 2017b). Different flow regimes (Knudsen number) may prevail along a nanoporous flow path and structural pore deformation additionally affects the rheology of the gas phase. Therefore, accurate permeability data and a profound understanding of the influencing processes are essential for predicting fluid transport in low-permeable rocks.

Anomalies in the flow of gases in technically manufactured media, such as narrow capillaries or channels, have been observed and studied in earlier works (Maxwell 1879; Kundt and Warburg 1875; Warburg 1876; Christiansen 1890; Knudsen 1909), mainly in the context of vacuum physics and rarefied gases. Maxwell (1879) was the first to describe deviations of gas flow behavior near a solid wall. Kundt and Warburg (1875) and Knudsen (1909) already use the term “slippage” and “slip flow coefficient,” which was later popularized by Klinkenberg (1941), who successfully expanded the slip flow concept to gases at higher pressures and natural rocks with particular relevance to natural gas exploration and production. A practical consequence of slip flow in constrained matrix systems of tight rocks is the increase in transport rates.

In modern research, there is a clear shift toward smaller diameters in low- (e.g., micro-/nano-fluidics) and high-pressure (e.g., oil recovery) applications. Several studies (Arkilic et al. 1997; Maurer et al. 2003; Roy et al. 2003; Ewart et al. 2007; Gruener & Huber 2008; Graur et al. 2009; Yamaguchi et al. 2011) performed fluid flow experiments through unstressed tubes and channels in the µm to nm range with low pressures (< 0.1 MPa). They reported deviations from the conventional first-order slip model toward higher Knudsen numbers (0.1 < \(\mathrm{Kn}\hspace{0.17em}\)< 10, transitional flow regime) that were then implemented in higher-order slip models accordingly. These low-pressure studies on artificial materials were then applied and extended to high-pressure shale gas models (Javadpour 2009; Civan 2010b) without proper verification on tight rocks.

However, reliable studies about fluid dynamics on low-permeable rocks at elevated pressures are rare as this requires the separation of interfering effects, such as poro-elastic deformation, and high-quality data (Fink et al. 2017b). Recent discussions (Moghaddam 2018; Fink et al. 2018) have shown that interpretations of transport mechanisms are still controversial. Further issues that have to be considered are uncertainties in the characterization of the material properties of tight rocks that remains a key issue due to rock–fluid interactions (e.g., gas sorption) and the complexity, heterogeneity and anisotropy of pore structures in the µm- to nm-scale.

Another fluid transport phenomenon, the “helium anomaly,” has been documented independently in several studies (Sinha et al. 2013; Ghanizadeh et al. 2014a, b; Gensterblum et al. 2014; Fink et al. 2017a, b). Here, the intrinsic (slip flow-corrected) permeability coefficients measured with helium on different types of rocks (shales and coal) were found to be consistently and significantly higher by a factor of 1.1 to 63 than those measured with other gases (Ar, N2, CH4, C2H6, CO2). This effect cannot be readily explained by general concepts and no satisfactory and unambiguous explanation has been given so far.

The objective of this research was to investigate, without interference of other (in particular rock-mechanic) processes, the fluid-dynamic effects associated with gas transport processes in tight natural rocks. For this purpose, rather than using natural rock samples, we selected well-characterized, artificial materials as rock analogues. A similar approach was pursued by Welch et al. (2017) who compared the results of image analysis and numerical simulations of nanoporous ceramic materials with experimental results.

This study is divided into two parts based on the type of analogue material used in single-phase gas transport experiments: (a) µm-sized capillary tubes and (b) a synthetic ceramic material with nm-sized pores. The ceramic filter material, similar to the one used by Welch et al. (2017), is a homogeneous porous medium with a very narrow pore size distribution. Due to its high rigidity (low poro-elasticity), fluid-dynamic effects, in particular gas flow behavior in the transitional flow regime, can be studied without interference of poro-mechanical effects. The capillary tube experiments consider potential changes in the rheology of gases with decreasing diameters and can be considered as “viscosimeter” test in a well-defined confined space. By comparison with literature data, it was attempted to identify fluid flow effects that might help explain phenomena, such as the “helium anomaly”, frequently observed in fluid flow tests on tight rocks.

2 Theoretical Background

2.1 Flow Regimes and Knudsen Number (Fluid Dynamics)

Whereas porosity reflects the storage capacity of a porous medium, permeability describes its ability of transporting fluids through the interconnected pore network (Tiab and Donaldson 2015). Fluid transport through porous media generally involves advective and diffusive transport mechanisms, the former having a higher transport capacity. Based on fluid type, P–T conditions, pore or fracture size and flux, fluid-dynamic processes control the gas flow in porous media in different flow regimes. In the order of decreasing characteristic flow path size(s) within porous media, the flow regimes change from turbulent flow to viscous flow, slip flow, transitional flow and Knudsen diffusion. While turbulent flow occurs at high pressures and flow rates/velocities, viscous flow prevails in micro-fractures and macropores, and slip flow and Knudsen diffusion become dominant in meso- to micropores. When gas permeates through a heterogeneous, porous medium containing fractures, micro- (< 2 nm), meso- (2–50 nm) and macropores (> 50 nm; IUPAC classification), it may be subject to several different flow regimes along its flow path (Ghanizadeh et al. 2014a).

Conventionally, fluid flow regimes from viscous flow to Knudsen diffusion are classified by means of the Knudsen number \({\text{Kn}}\) (–), which is defined as the mean free path length \(\lambda\) (m) of the gas divided by the characteristic length (i.e., pore diameter) \(d\) (m) of the porous medium.

$${\text{Kn}} = \frac{\lambda }{d}.$$
(1)

Depending on the classification scheme, the transition from viscous flow to slip flow is assumed at \({\text{Kn}}\) = 0.01 (Karniadakis et al. 2006) or \({\text{Kn}}\) = 0.001 (Roy et al. 2003; Javadpour et al. 2007; Wang et al. 2008). In this study, we adopt the latter value (\({\text{Kn}}\) = 0.001). It should be kept in mind, however, that the transition is not abrupt but gradual. The mean free path length \(\lambda\) (m) describes the average distance a gas molecule travels until it collides or interacts with another molecule. According to the kinetic theory of gases, the mean free path length is inversely proportional to the gas pressure and can be calculated as (Cussler and Cussler 1997):

$$\lambda = \frac{{k_{{\text{B}}} T}}{{\sqrt 2 \pi d_{m}^{2} P}},$$
(2)

where \(k_{{\text{B}}}\) (J K−1) is the Boltzmann constant (1.381 · 10–23), \(T\) (K) is the temperature, \(d_{m}\) (m) is the kinetic diameter and \(P\) (Pa) is the pressure. The mean free path model involves strong simplifications since gas molecules are not “hard spheres” and their diameters are not well defined. Furthermore, attraction and repulsion forces act as a function of intermolecular distance during the collision. These effects are reflected in the dynamic viscosity \(\mu\) (Pa s) and gas density \(\rho\) (kg m−3), which are measurable quantities that can be used to directly express the mean free path length as follows (Bird 1983):

$$\lambda = \frac{\pi \mu }{{\rho \sqrt {2 \pi \frac{{k_{{\text{B}}} N_{{\text{A}}} }}{M}T} }},$$
(3)

where \(N_{{\text{A}}}\) (mol−1) is the Avogadro constant and \(M\) (kg mol−1) is the molar mass. At isothermal conditions, mean free path lengths increase with decreasing pressure (gas density) and, in consequence, Knudsen numbers increase while they decrease with increasing characteristic length.

2.1.1 Viscous, Continuum (Darcy) and Turbulent (Non-Darcy) Flow (Kn < 0.001)

Laminar, pressure gradient-driven volume flow (Darcy flow) is the most important transport mechanism in the matrix system of a porous medium (Gensterblum et al. 2015). For very small Knudsen numbers (< 0.001), the fundamental physical equation for the description of flow of viscous substances is the Navier–Stokes equation (Javadpour et al. 2007). Analytical solutions of the Navier–Stokes equation for flow through capillaries and rocks can be derived using simplifying assumptions such as laminar flow and certain geometrical shapes (Tiab and Donaldson 2015). Fluid flow through a cylindrical capillary tube can then be expressed by Hagen–Poiseuille’s law for incompressible fluids (e.g., water):

$$Q_{{\text{v}}} = - \frac{{\pi r^{4} }}{8 \mu }\frac{{{\text{d}}P}}{{{\text{d}}x}}$$
(4)

with \(Q_{{\text{v}}}\) (m3 s−1) the (liquid) volumetric fluid flow rate, \(r\) (m) the radius of the capillary tube, \(\mu\) (Pa s) the viscosity of the permeating fluid and \(\frac{{{\text{d}}P}}{{{\text{d}}x}}\) (Pa m−1) denoting the pressure gradient.

Fluid flow through porous media, such as rocks, can be described phenomenologically by Darcy’s law for incompressible fluids, where the Darcy velocity \(u\) (m s−1) is proportional to the pressure gradient \(\frac{{{\text{d}}P}}{{{\text{d}}x}}\) (Ho and Webb 2006):

$$u = \frac{{Q_{{\text{v}}} }}{A} = - \frac{k}{\mu }\frac{{{\text{d}}P}}{{{\text{d}}x}} \leftrightarrow Q_{v} = - \frac{k A}{\mu }\frac{{{\text{d}}P}}{{{\text{d}}x}}.$$
(5)

Here, \(A\) (m2) is the cross-sectional area of the sample plug and \(k\) (m2) the permeability coefficient. Frequently, k is reported in units of Darcy (D), where 1 D = 0.987 × 10–12 m2.

On a laboratory scale, turbulent or inertial flow may occur at high absolute and/or differential pressures that lead to interferingly high flow rates, causing inertial effects and non-Darcy flow (Rushing et al. 2004; Ziarani and Aguilera 2011).

The onset of non-Darcy, turbulent flow is estimated using the critical Reynolds number \({\text{Re}}\), which is defined as:

$${\text{Re}} = \frac{{Q_{{\text{v}}} }}{A}\frac{d}{\mu }.$$
(6)

While Darcy’s law is valid at low Reynolds numbers (\({\text{Re}}\) < 1–2300; Javadi et al. 2014), it fails to describe flow properties when turbulent effects emerge. Turbulent flow can be described by Forchheimer’s equation, which is a second-order extension of Darcy’s law (Forchheimer 1901):

$$\frac{{{\text{d}}P}}{{{\text{d}}x}} = - \frac{\mu }{k}\frac{{Q_{{\text{v}}} }}{A} + \beta \rho \left( {\frac{{Q_{{\text{v}}} }}{A}} \right)^{2},$$
(7)

where \(\beta\) (–) is the inertial flow or non-Darcy coefficient. The Forchheimer equation expresses a nonlinear relationship between fluid flow rate and pressure gradient (Ghanizadeh et al. 2017). As inertial effects diminish (i.e., \(\beta\) or \(Q_{{\text{v}}}\) approaching zero), the second-order term becomes negligible (Ho and Webb 2006). The occurrence of turbulent flow can be assessed by means of the Reynolds number \({\text{Re}}\) and graphical representations.

2.1.2 Slip Flow (0.001 < Kn < 0.1)

Slip flow is a nonlinear, non-Darcy effect that occurs in small pores (e.g., in shales) at low gas pressures. As gas pressure decreases, the mean free path of gas molecules increases and may approach the characteristic diameter of the transport pores. At these conditions, collisions between molecules and pore walls (slippage) become more important than molecule–molecule collisions (Klinkenberg 1941; Ghanizadeh et al. 2014a) and the flow velocity of the individual molecules increases. Thus, the pressure-dependent apparent gas permeability coefficient \(k_{{{\text{gas}}}}\) increases with decreasing gas pressure. The contribution of gas slippage increases with decreasing pore size, i.e., in the meso- to microporous matrix system. Klinkenberg (1941) experimentally demonstrated the linear relationship between \(k_{{{\text{gas}}}}\) and the reciprocal mean gas pore pressure \(P_{{\text{p}}}\) (Pa) for dry rocks:

$$k_{{{\text{gas}}}} = k_{\infty } \left( {1 + \frac{b}{{P_{{\text{p}}} }}} \right).$$
(8)

Here, \(k_{\infty }\) (m2) is the Klinkenberg-corrected (“intrinsic”) gas permeability coefficient and is interpreted as the equivalent liquid or absolute permeability \(k\), and \(b\) (Pa) is the (first order) gas slippage factor. \(P_{{\text{p}}}\) is the arithmetic mean of the upstream and the downstream pressure:

$$P_{{\text{p}}} = \frac{{P_{{\text{u}}} + P_{{\text{d}}} }}{2}.$$
(9)

The coefficient \(k_{\infty }\) is the y-axis intercept of the regression line of measured \(k_{{{\text{gas}}}}\) plotted versus the reciprocal of \(P_{{\text{p}}}\) (Klinkenberg plot). Assuming flow through cylindrical capillaries, the gas slippage factor \(b\) can be related to the mean free path, the mean pore pressure and mean transport pore throat diameter \(d\) (m) by:

$$\frac{b}{{P_{{\text{p}}} }} = \frac{8 C \lambda }{d}$$
(10)

with \(C\) ≈ 0.9 being the dimensionless Adzumi constant (Adzumi 1937). The relationship indicates that \(b\) depends on the pressure-independent gas properties \(\lambda \cdot P_{{\text{p}}}\) (where \(\lambda \sim \frac{1}{{P_{{\text{p}}} }}\)) and the characteristic diameter of the transport pores (\(d\)).

2.1.3 Transitional Flow (0.1 < Kn < 10)

Transitional flow is described as an intermediate flow regime where both slip flow and diffusional transport can occur. Conventional fluid-dynamic equations (e.g., continuum equation) fail when gas molecules travel more independently through the pores (Knudsen flow) rather than interacting with each other (Ziarani and Aguilera 2011).

In order to account for deviations from the first-order slip flow concept, several higher-order slip flow extensions have been proposed in the past based on Maxwell’s theory (Maxwell 1879; Beskok and Karniadakis 1999). Similar slip flow modifications can be applied to apparent gas permeability in porous media (Civan 2010b). The second-order slip extension of the Klinkenberg equation (Eq. 8) used to describe fluid flow in the transitional flow regime can be written as (Tang et al. 2005):

$$k_{{{\text{gas}}}} = k_{\infty } \left( {1 + \frac{{b_{1} }}{{P_{{\text{p}}} }} + \frac{{b_{2} }}{{P_{{\text{p}}}^{2} }}} \right)$$
(11)

with \(b_{2}\) (MPa2) being the second-order gas slippage factor.

2.1.4 Knudsen Diffusion or Free Molecular Flow (Kn > 10)

At high Knudsen numbers, fluid transport in porous media is mainly controlled by Knudsen diffusion through voids. At  \(\mathrm{Kn}\) > 10, molecular interaction virtually no longer occurs since the mean free path is significantly greater than the diameter of the flow channel. As a result, individual molecules move autonomously and reflect diffusely (non-oriented) upon collision with the pore wall. Knudsen diffusion can co-occur with molecular and surface diffusion (Ziarani and Aguilera 2011).

3 Materials and Methods

The study on fluid flow in fused silica capillaries (Sect. 3.1) comprises fluid flow experiments (Sect. 3.7.1) and scanning electron microscopy combined with broad ion beam polishing (BIB-SEM, Sect. 3.6). For the study using the synthetic nanoporous ceramic material (Sect. 3.2), fluid flow experiments (Sect. 3.7.2) were performed and additionally the pore structure of the ceramic was thoroughly characterized using low-pressure nitrogen sorption (LPNS, Sect. 3.4), mercury intrusion porosimetry (MIP, Sect. 3.5) and BIB-SEM (Sect. 3.6).

3.1 Fused Silica Capillaries

The first set of fluid flow experiments was performed on fused silica capillaries with nominal inner diameters of 10, 5 and 2 µm (CS-Chromatographie GmbH, Germany). The outer diameter of all capillary columns is 150 µm. The outside of the capillary is coated with a polyimide polymer to add resistivity and flexibility, whereas the inside is uncoated fused silica.

3.2 Nanoporous Ceramic Membrane

The second set of fluid flow experiments was conducted on a synthetic sintered ceramic plug with a nominal mean pore size of 50 nm (Cobra Technologies B.V., Netherlands). The crystal white, cylindrical plug has a dry mass of 23.7597 g, a diameter of 38.1 mm and a length of 7.5 mm. It consists of pure aluminum oxide (> 99% α-Al2O3, corundum) and is very hard, rigid and brittle. Prior to each experiment, the sample material was dried in a vacuum oven at 105 °C for several days until constant weight was reached.

3.3 Helium Pycnometry—Porosity and Specific Pore Volume Determination on the Unconfined Sample

Porosity \(\Phi\) (%) and specific pore volume \(V_{{{\text{sp}}}}\) (m3 kg−1) were determined from the skeletal volumes, measured by helium pycnometry (expansion), and the bulk volumes calculated from the dimensions of the plug. A detailed description of this procedure can be found in Webb (2001). The specific pore volume \(V_{{{\text{sp}}}}\) is the pore volume \(V_{{\text{p}}}\) (m3) normalized to the sample mass \(m\) (kg). Bulk (\(\rho_{{{\text{bulk}}}}\)) density, grain (\(\rho_{{{\text{grain}}}}\)) density (kg m−3), porosity \({\Phi }\) and specific pore volume \(V_{{{\text{sp}}}}\) can be interconverted and are related as follows:

$$\Phi = 1 - \frac{{\rho_{{{\text{bulk}}}} }}{{\rho_{{{\text{grain}}}} }} = 1 - \frac{1}{{\rho_{{{\text{grain}}}} V_{{{\text{sp}}}} + 1}}.$$
(12)

3.4 Low-Pressure Nitrogen Sorption (LPNS)

Low-pressure nitrogen sorption (LPNS) measurements were taken on a Micromeritics Gemini VII 2390t instrument (Micromeritics Instrument Corporation, USA). An aliquot of the sample was crushed and sieved to grain fractions varying from 63 to 300 µm. Around 0.5 g of sample material was placed in a measuring cell and degassed using a VacPrep 061 (Micromeritics Instrument Corporation) at 378 K for at least 10 h. Subsequently, adsorption and desorption were measured in a liquid N2 bath at 101 relative pressure steps between 0.001 and 0.995 p/p0 during which the amount of adsorbed nitrogen was recorded. Measured isotherms were analyzed and interpreted by means of different theories and models. The total pore volume was calculated by the Gurvich rule at a p/p0 of 0.995 (Gurvich 1915), specific surface area was derived by the Brunauer–Emmett–Teller (BET) (Brunauer et al. 1938) method and the pore size distribution by the Barrett–Joyner–Halenda (BJH) (Barrett et al. 1951) theory. The theoretical background and applicability of the different methods were recently discussed by Bertier et al. (2016) and Seemann et al. (2017).

3.5 Mercury Intrusion Porosimetry (MIP)

The MIP measurement was performed using a Micromeritics Autopore IV 9500 series instrument (Micromeritics Instrument Corporation). An aliquot of the sample was crushed into small fragments with lengths ranging from 1 to 5 mm to fit into the penetrometer. Mercury was injected continuously up to the final pressure of 413 MPa (60,000 psi). The amount of injected Hg was recorded as a function of injection pressure and the equivalent pore (throat) radius [\(r_{{{\text{eq}}}}\) (m)] distribution was determined using the capillary pressure equation (Washburn 1921) for cylindrical capillaries:

$$r_{\text{eq}} = - \frac{2 \gamma_{\text{Hg/air}} \cos \theta }{P_{\text{cap}}} .$$
(13)

Here, \(\gamma\) (0.485 N/m) is the Hg/air interfacial tension, \(\theta\) (140°) is the contact angle between the fluid (Hg) and solid phase, \(P_{\text{cap}}\) (Pa) is the mercury (capillary) pressure. The smallest equivalent pore radius that can theoretically be accessed (at a Hg pressure of 413 MPa) is approximately 1.8 nm.

3.6 BIB-SEM

3.6.1 Fused Silica Capillary

Twenty capillaries, aligned side by side, were glued onto a cover glass (2 × 8 × 10 mm) using epoxy resin. The capillaries were impregnated with low-viscosity epoxy resin to prevent contamination with sputter dust during subsequent ion-milling. The top surface, perpendicular to the long axes of the capillaries, was ion-milled, using a Leica TIC3X ion-mill at 5 kV acceleration voltage for 2.5 h producing a smooth and dust-free cross section of around 1 mm2. Subsequently, the cross section was sputter-coated with a 3 nm thick layer of Tungsten and then imaged in the SEM using secondary electrons. The sizes of seven capillaries in the center of the cross section were measured with a pixel resolution of 6 nm.

3.6.2 Ceramic Disk

The ceramic sub-sample was ion polished using a JEOL (SM-09010) cross section polisher at 6 kV acceleration voltage for 8 h producing a planar and damage-free cross section of around 0.5 mm2. After tungsten coating (7.5 nm thick), twenty random locations of the BIB-cross section were imaged in the SEM at a magnification of 40,000 (pixel size of 7.3 nm). The images were segmented by thresholding to quantify the porosity and to obtain the pore space geometries.

3.7 Fluid Flow Experiments

For the non-steady-state and steady-state fluid flow tests in capillaries (Sect. 3.7.1) and the ceramic disk (Sect. 3.7.2), different experimental setups (Fig. 1) were used.

Fig. 1
figure 1

Schematic sketches of a capillary non-steady-state fluid flow setup for the determination of gas viscosity, b steady-state permeability setup with counteracting piston pumps and c non-steady-state permeability setup with additional reservoir volumes (Vu, Vd). All setups are temperature controlled (thermocouples type K, NiCr–Ni) and equipped with high-precision pressure transducers (KELLER AG für Druckmesstechnik; model PAA-33X) (e.g., Pu, Pd) and valves (V1–V4). For fluid flow experiments performed on the ceramic plug in the isostatic cell (b, c), the confining pressure Pc (stressed condition) can be adjusted

3.7.1 Experimental Viscosity Tests (FS Capillaries)

Single-phase gas flow experiments through µm-diameter fused silica capillary tubes were conducted at unstressed conditions with different gases (He, Ar, N2, H2, CH4, C2H6 and CO2) using the non-steady-state pressure pulse decay (PPD) technique in a closed system.

Each capillary was placed between an upstream and a downstream reservoir, both of known volumes (Vu and Vd reservoir volumes between valve V1 and top or bottom end of the capillary, respectively, see Fig. 1a). Each reservoir was equipped with a pressure sensor. Due to low gas flow rates through the capillaries, the system was carefully leak-tested and the reservoir volumes were kept very small (< 1 ml) for increased sensitivity.

Prior to each measurement, the system was set to a desired downstream pressure (Pd) and equilibrated. Hereafter, the upstream pressure (Pu) was increased. The differential pressure (initial differential pressure < 0.7 MPa) results in a pressure gradient across the capillary tube that declines over time until pressure equilibrium is reached between the two compartments. The pressures in the upstream and downstream compartments were monitored continuously. Mean gas pressures ranged from approx. 0.2–3.5 MPa. The pressure transients and temperature data were then used to calculate the viscosity coefficients (Sect. 3.7.3).

3.7.2 Permeability Under Controlled Effective Stress (Ceramic Disk)

Single-phase gas permeability measurements under confined stress were conducted using the isostatic flow cell setup shown in Fig. 1b, c. The cylindrical sample plug was placed between two stainless steel pistons, each with a capillary tube inlet for gas. Radial and circular grooves on the base of the pistons acted as flow distributors. The sample surfaces were sealed by a lead foil and a double layer of rubber to prevent fluid bypass and to separate the plugs from the surrounding confining fluid (water).

Gas permeability tests were performed with different gases (He, Ar, N2, CH4, H2 and CO2) and gas mixtures (Ar–N2, ratio of 1:1). Measurements were taken successively at different effective stresses, with \(\sigma\) being defined according to Terzaghi’s principle (1923) as the difference between the isostatic confining pressure \({P}_{\mathrm{c}}\) and pore pressure \({P}_{\mathrm{p}}\). To achieve the desired effective stress levels (\(\sigma \hspace{0.17em}\) = 10, 20, 30 and 40 MPa), \({P}_{\mathrm{c}}\) was adjusted between 10.2 and 50.5 MPa and \({P}_{\mathrm{p}}\) from 0.2 to 30.5 MPa (Table 1). Confining and pore pressures were adjusted by piston pumps (260D syringe pump, Teledyne ISCO). The measuring sequence for a given gas and stress level was started at low pore pressures, which were then successively increased to higher values. Throughout all measurements, the gas pressure difference across the samples was kept small (approx. 0.1 MPa) to avoid turbulent flow. Measurements for a given gas type were either conducted as steady-state experiments or non-steady-state pressure pulse decay (PPD) experiments, depending on various experimental considerations (e.g., safety issues, suitability of technical components).

Table 1 Protocol of permeability measurements on the ceramic plug at 35 °C

Steady-state permeability tests were performed with counteracting piston pumps connected to the upstream and downstream sides of the sample, respectively (Fig. 1b). The control systems of the gas-filled pumps ensured that the desired pressure difference was kept constant. After applying the differential pressure, steady-state conditions adjusted within a few minutes and the resulting volumetric gas flow rate \({Q}_{\mathrm{v}}\) and the corresponding pressures (Ppump,u and Ppump,d) were recorded continuously by both pumps. The true \({P}_{\mathrm{u}}\) and \({P}_{\mathrm{d}}\) (up- and downstream pressure, respectively) were recorded by independent pressure sensors connected directly to the pistons.

Non-steady-state permeability experiments were conducted by means of the PPD technique in a closed system with calibrated reservoirs on the upstream and downstream sides (Vu, Vd) of the sample (Fig. 1c). The experimental approach is the same as described in Sect. 3.7.1 for the capillaries. Unlike the capillary experiments, where flow rates were low and reservoir volumes were kept small (< 1 ml), flow rates through the ceramic filter disk were very high and, thus, additional reservoirs of approximately 300 ml were installed to prolong the pressure decay and improve data quality.

3.7.3 Evaluation of Fluid Flow Experiments

Due to their high compressibility, gases expand as they move under the effect of the pressure gradient. Thus, the volumetric flow rate \(Q_{{\text{v}}}\) (m3 s−1) increases along the flow path. In order to take this effect into account, the mass flow rate \(Q_{{\text{m}}}\) (kg s−1) is used for evaluating gas flow data, as it is constant due to mass balance constraints.

For evaluating non-steady-state gas transport through a single capillary, the following “pseudo-steady state” formula, based on the Hagen–Poiseuille equation (Eq. 4), was used to calculate the experimental viscosity:

$$\mu_{{{\text{exp}}}} = \frac{{\pi r^{4} \left( {P_{{\text{u}}}^{2} - P_{{\text{d}}}^{2} } \right) M}}{{16 Q_{{\text{m}}} \left( {P_{{\text{p}}} ,T} \right)L z\left( {P_{{\text{p}}} ,T} \right) R T}}.$$
(14)

Here, \(r\) (m) represents the inner capillary radius, \(P_{{\text{u}}}\) and \(P_{{\text{d}}}\) (Pa) are the pressures at the upstream and downstream ends of the capillary, \(M\) (kg mol−1) is the molar mass, \(Q_{{\text{m}}}\) (kg s−1) is the mass flow rate, \(L\) (m) is the length of the capillary tube, \(z\) (–) is the compressibility factor, \(R\) (J mol−1 K−1) is the universal gas constant (8.314) and \(T\) (K) denotes the temperature. The derivation of Eq. 14 is given in “Appendix 1.” Pressure and temperature recordings over defined time intervals were used to calculate \(\rho_{{{\text{gas}}}} \left( {P_{{\text{p}}} ,T} \right)\) (kg m−3). The mass flow rate \(Q_{{\text{m}}}\) was determined from the changes of \(\rho_{{{\text{gas}}}}\) over time and the known compartment volumes \(V_{{\text{u}}}\) and \(V_{{\text{d}}}\) [m3]:

$$Q_{{\text{m}}} = \frac{{{\text{d}}\rho_{{{\text{gas}}}} \left( {P_{{\text{p}}} ,T} \right)}}{{{\text{d}}t}}V_{{\text{u,d}}}.$$
(15)

\(\rho_{{{\text{gas}}}} \left( {P_{{\text{p}}} ,T} \right)\) and \(z\left( {P_{{\text{p}}} ,T} \right)\) were calculated using the GERG-2008 Wide-Range Equation of State (EoS) for Natural Gases and Other Mixtures (Kunz and Wagner 2012).

The steady-state gas permeability coefficients \(k_{{{\text{gas}}}}\) were calculated using Darcy’s law based on the mass flow rate \(Q_{{\text{m}}}\):

$$k_{{{\text{gas}}}} = Q_{{\text{m}}} \left( {P_{{\text{p}}} ,T} \right)\frac{{z\left( {P_{{\text{p}}} ,T} \right)R T}}{M}\frac{{2 \mu \left( {P_{{\text{p}}} ,T} \right) L }}{{A \left( {P_{{\text{u}}}^{2} - P_{{\text{d}}}^{2} } \right)}}.$$
(16)

The derivation of Eq. 16 is given in “Appendix 2.” The dynamic viscosity \(\mu \left( {P_{{\text{p}}} ,T} \right)\) for each individual mean pressure step was obtained from the NIST chemistry webbook (Lemmon et al. 2021). The piston pumps continuously recorded the flow rate \(Q_{{\text{v}}}\) (m3 s−1) over time. \(Q_{{\text{m}}}\) was calculated based on \(\rho_{{{\text{gas}}}} \left( {P_{{\text{p}}} ,T} \right)\) (kg m−3) and \(Q_{{\text{v}}}\).

$$Q_{{\text{m}}} = Q_{{\text{v}}} \rho_{{{\text{gas}}}} \left( {P_{{\text{p}}} ,T} \right).$$
(17)

The gas permeability coefficients \(k_{{{\text{gas}}}}\), from the non-steady-state experiments, were calculated according to the following equation (cf. Cui et al. 2009; Ghanizadeh et al. 2014a):

$$k_{{{\text{gas}}}} = - \frac{{s_{1} \mu L c_{{\text{g}}} }}{{f_{1} A\left( {\frac{1}{{V_{{\text{u}}} }} + \frac{1}{{V_{{\text{d}}} }}} \right)}}.$$
(18)

Here, \(V_{{\text{u}}}\) and \(V_{{\text{d}}}\) are the calibrated upstream and downstream volumes, respectively. Parameter \(c_{{\text{g}}}\) (Pa−1) is the isothermal gas compressibility which was computed numerically from the gas density change with pressure \(\left( {\frac{{{\text{d}}\rho_{{{\text{gas}}}} }}{{{\text{d}}P}}} \right)_{T}\). Parameter \(s_{1}\) is the slope of the late-time solution of the mass transfer of the amount of substance \(n\) between the upstream and downstream side (\(\ln \left( {n_{{\text{u}}} \left( t \right) - n_{{\text{d}}} \left( t \right)} \right)\) vs. time). In order to evaluate permeability on this basis, the amount of substance \(n\) was determined from the recorded \(P_{{\text{u}}} \left( T \right)\) and \(P_{{\text{d}}} \left( T \right)\) by means of the EoS. Parameter \(f_{1}\) is the mass flow correction factor that depends on the storativity of the sample’s pore volume \(V_{{\text{p}}}\) with respect to the reservoir volumes (\(V_{{\text{u}}} + V_{{\text{d}}}\)) (Cui et al. 2009). If \(V_{{\text{u}}} + V_{{\text{d}}}\) is much larger than \(V_{{\text{p}}}\), as it is the case in our experiments, the contribution of \(f_{1}\) can be neglected.

4 Results and Discussion

4.1 Capillary Tube Experiments: Experimental Versus Literature Viscosity Data

In order to ensure laminar flow of gas in the capillary tube experiments throughout the experimental pressure range, a linearity analysis was performed by plotting the Darcy pressure drop versus the mass flow rate/flux. This is exemplary shown for argon and nitrogen for the capillaries of 10 and 2 µm in diameter (Fig. 2). Maximum mean pore and differential pressure applied were 3.5 and 0.7 MPa, respectively. The trends for both capillary sizes are linear and pass through the origin, hence, indicating laminar viscous flow and no influence of inertial effects. Previous studies (e.g., Rushing et al. 2004; Ghanizadeh et al. 2017) have shown that turbulent flow would lead to a concave upward deviation since the mass flow rate/flux decreases as a result of obstructively high pressures. A more detailed explanation can be found in Rushing et al. (2004) and Ghanizadeh et al. (2017).

Fig. 2
figure 2

Darcy pressure drop versus mass flux displaying linear trends for both capillary sizes indicating linear Darcy flow behavior. Maximum differential pressures for the 10 µm (a) and 2 µm (b) capillaries were 0.8 and 1.2 MPa, respectively

Based on the fluid flow experiments, viscosity was initially calculated using the “nominal” inner capillary diameters provided by the manufacturer (10, 5 and 2 µm). This resulted in significant discrepancies between the calculated and published viscosity data (NIST chemistry webbook; Lemmon et al. 2021), particularly for the 2 µm capillary. As these differences could also be attributed to inaccurate diameter data, the fluid flow data in the laminar flow regime were matched with the literature viscosity data of all gases by adjusting the capillary diameter. The adjusted “true” inner diameters of 10.02, 4.95 and 1.59 µm (Fig. 3) resulted in better agreement of experimental and literature data for all gases. Evidently, the difference between “nominal” and “true” diameter is relatively large (~ 20%) for the smallest capillary.

Fig. 3
figure 3

Viscosity was measured by flow through capillary tubes with different nominal, inner diameters: a 2 µm, b 5 µm and c 10 µm. Left: Experimental viscosity values (solid symbols) plotted against mean pore pressure to show deviations from the corresponding viscosity values from the NIST chemistry webbook (dashed lines). Right: experimental viscosity plotted against the Knudsen number to indicate flow regime-dependent changes. Note: dnominal is the nominal inner capillary diameter reported by the manufacturer, dfit is the adjusted inner diameter for consistent results and L depicts the capillary length

The adjusted capillary diameters were subsequently determined independently using scanning electron microscopy (BIB-SEM). Microphotographs of the capillary with a nominal inner diameter of 2 µm are shown in Fig. 4. The average inner diameter of the capillary cross section determined by this method was 1.56 µm (number = 7, standard error = 0.01 µm, median = 1.55 µm, mode = 1.58 µm, standard deviation = 0.03 µm). This is in excellent agreement with the adjusted diameter of 1.59 µm deduced from capillary flow experiments and supports the validity and accuracy of this approach.

Fig. 4
figure 4

BIB-SEM cross section of the capillary with a “nominal” diameter of 2 µm. a Overview showing a single silica capillary including outer polymer coating and inner diameter hole. b Detail showing the “true” measured capillary inner diameter of 1.56 µm

Using the adjusted capillary diameters, the experimentally determined viscosity of all gases agrees well with the literature data, especially at pressures above approximately 0.1 MPa. Deviations primarily occur at lower pressures and for smaller capillary diameters, i.e., at higher Knudsen numbers. The Knudsen numbers were calculated according to Eqs. 1 and 3 using the adjusted capillary diameters.

Theoretically, flow regimes covered by the set of experimental conditions range from the Darcy (viscous) to the slip flow regime. The transition between both flow regimes is indicated in Fig. 3 for the 10 and 5 µm capillaries. The experiments performed on the 2 µm capillary cover only the nominal “slip flow regime.” The viscosities measured with the 2 µm capillary show significant scatter, in particular for helium and argon. This indicates a decrease in accuracy of the measurements at very low mass transfer rates.

4.2 Flow Experiments with Synthetic Porous Ceramic Disk

4.2.1 Comparison of Pore Structure Information from He-Pycnometry, MIP, LPNS and BIB-SEM

In an attempt to relate the fluid flow behavior of gases with the pore structure, the synthetic porous ceramic material was thoroughly characterized using He-pycnometry, MIP, LPNS and BIB-SEM. The results of these measurements are summarized in Table 2. The bulk densities determined by caliper measurements and from immersion into mercury during MIP were 2.85 and 2.81 g cm−3, respectively. Grain densities from helium pycnometry and MIP were 3.99 g cm−3 and 3.89 g cm−3, respectively, and, thus, are very close to the literature density values of Al2O3 (3.95 g cm−3) and corundum (4.02 g cm−3). The specific surface area determined by LPNS and evaluated using the BET method was 4.93 m2 g−1.

Table 2 Summary of helium pycnometry, mercury porosimetry (MIP), low-pressure nitrogen sorption (LPNS) and BIB-SEM results.

Porosity was also calculated from the specific pore volume (“Gurvich volume”) of the LPNS measurement and the grain density from helium pycnometry. Porosity values determined by helium pycnometry, MIP and LPNS are very similar, with 28.64%, 27.66% and 28.93%, respectively (Table 2). The variation of 1% (absolute) is within the range of inherent measuring inaccuracy.

The BIB-SEM images clearly show that the sample is highly porous with evenly distributed and well-connected pores of similar size (Fig. 5a). These pores are irregularly shaped with no preferred orientation and have neither slit nor circular pore geometries as commonly assumed for pore models. Instead, the pore space consists of a complex network of void spaces with connecting pathways. After segmentation of pores (Fig. 5b), the two-dimensional (2D) porosity and pore size distributions were computed by thresholding from a series of images scanned at a magnification of 40,000×. The “visible” porosity from different image locations ranges from 13.01 to 18.08% with an average value of 15.38%. This is much lower than the true porosity values determined by the conventional methods helium pycnometry, MIP and LPNS (Table 2).

Fig. 5
figure 5

BIB-cross-sectioning and segmentation of the same image. a Typical 2D pore structure (40,000×) in the BIB-section showing “visible” pore space (dark gray to black) and Al2O3 matrix (light gray). Within big pores, actual 3D pore space can be observed as the grayscale changes due to the incident electron beam. b Segmentation to differentiate between matrix (blue) and pore space (white). Pore geometries and pore size distribution were obtained by fitting coextensive ellipses (black)

The substantial underestimation of the true porosity by BIB-SEM image analysis may have different causes: Firstly, within the larger pores, shades with changing gray scales can be observed (Fig. 5a). The intensity of these shades changes as a function of the incident electron beam of the SEM and therefore is an indicator for the depth of the pores and of actual 3D porosity. Missing these pore volume fractions during pore segmentation may contribute to an underestimation of total porosity, but it cannot explain a porosity discrepancy of almost 50%. A second reason for this discrepancy could be the process of sample preparation. BIB sputtering can lead to minor re-deposition of Al2O3 dust in the pores that reduces the pore volume. Due to the hardness of the material, the polishing process requires more time, which subsequently may lead to even more re-deposition. Since SEM cannot distinguish between the monomineralic α-Al2O3 solid and dust phase, actual pores that are constricted or clogged by re-deposited material would be misinterpreted as matrix. Due to the small size of the selected BIB-section area (< 100 µm2), it might also be that the section is not representative for the bulk sample.

Welch et al. (2017) observed similar discrepancies between gravimetrically (Archimedes’ principle) measured porosity (15.2%) of nanoporous ceramics and porosity determined from FIB-SEM image analysis (initially 6.6%) They report that by eliminating “high grayscale artifacts from the pore walls” the porosity from image analysis could be increased to 16.5% to closely match the gravimetric porosity. However, they do not describe the procedure of artifact elimination in detail and it appears that the “true” porosity has to be known beforehand in order to achieve a realistic best fit. Definitely, the assessment of porosity by micro-imaging techniques still requires further studies and critical scrutiny.

Pore size distributions and mean pore sizes can be derived by MIP, LPNS and BIB-SEM despite principal differences in the measuring techniques. The fluid invasion techniques MIP and LPNS are in very good agreement and show a sudden increase in cumulative pore volume in a very narrow pore diameter range from approximately 45–105 nm (Fig. 6a). This is very characteristic for the synthetic, homogeneous, porous material, with a narrow unimodal PSD, used in this study. The total pore volume of around 0.1 cm3 g−1 is mainly located in macropores (> 50 nm), whereas only a small portion lies within the mesoporous range (2–50 nm). Micropores (< 2 nm; according to the IUPAC classification) are completely absent. Figure 6b illustrates that the pore size distribution derived from both methods is so narrow that the pore-filling process is resolved by only four data in the MIP measurement and 11 data points in the LPNS. The most prominent equivalent pore diameters (mean value of the pore size distribution) are 69 and 68 nm for MIP and LPNS (Table 2), respectively.

Fig. 6
figure 6

Results of MIP, LNPS and BIB-SEM measurements. a Cumulative pore volumes plotted as a function of the logarithm of the equivalent pore diameter. b Pore size distribution indicating the most prominent pore diameter. The majority of pores are macro- and mesopores, micropores are absent. c Frequency distribution of quantitative BIB-SEM analysis and MIP data. Frequency of MIP was derived from a simplified capillary bundle model (see Sect. 4.2.5). Normalized pore count fraction is plotted as a function of equivalent circular pore diameter. The amount of same-sized pores (bins) was normalized to the total pore count

For easier comparison, porosity fractions derived from BIB-SEM image analysis have been converted into specific pore volumes using the grain density derived from helium pycnometry (3.995 g cm−3). Results in Fig. 6a show a gradual increase in cumulative pore volume in the pore diameter interval from approx. 400 down to 30 nm. In accordance with the aforementioned porosity difference, the total specific pore volume of 0.046 cm3 g−1 accounts for only about half of the volume obtained from the other methods (0.1 cm3 g−1). The pore size distribution (Fig. 6b) is broad and flat with a peak located at around 90 nm and a most prominent equivalent pore diameter of 109 nm.

However, presenting the BIB-SEM data as frequency distribution based on the total count of pores (Fig. 6c), the distribution of most abundant pores collectively shifts to smaller equivalent pore diameters, as small pores are naturally more frequent than large ones. The mode of distribution is represented by the peak at 60 nm equivalent pore diameter (11%), which is in the same range as MIP and LPNS. In comparison, the frequency distribution of MIP, estimated by the capillary bundle model in Sect. 4.2.5, is very narrow and its maximum at 70 nm is substantially higher (41%), as the total pore size range is smaller (50–90 nm).

It is evident that BIB-SEM results do not only deviate from MIP and LPNS results in terms of absolute porosity values but also in terms of pore size distributions and equivalent pore diameters. As stated above, the reasons for this discrepancy can lie in the manner of presenting the results. It can be seen from Fig. 5b that pore areas are complex and often have irregular shapes that were mathematically approximated by ellipses of equal surface area. The segmentation algorithm interprets a cluster of connected pores in the 2D image as one large pore and, thus, distorts the true pore size and collectively shifts the distribution to higher average pore sizes. However, BIB-SEM quantifies actual pore bodies in 2D opposed to 3D pore throats with MIP.

Comparison of the intrusion and extrusion curve of the MIP measurement (Fig. 6a) shows a hysteresis of 0.062 cm3 g−1, i.e., approximately 63% of the injected mercury remains trapped in the pore space of the sample. The volume of the extruded mercury can be assumed to represent the pore space responsible for fluid flow (“transport porosity”; “effective porosity”), which in this case amounts to a porosity of about 10%.

4.2.2 Non-Darcy Turbulent Flow and Stress-Sensitivity Analysis

Inertial effects (turbulent flow) and stress effects are commonly observed in gas permeability measurements on rocks. Stress-induced permeability changes as a result of pore compressibility have been reported extensively by various authors for different types of sedimentary rocks (Dong et al. 2010; Chalmers et al. 2012; Ghanizadeh et al. 2014a, b; Heller et al. 2014; Letham and Bustin 2016; Fink et al. 2017a, b).

To check for any of these interfering poro-elastic stress effects, permeability tests with helium were conducted at four different stress levels prior to the actual flow tests with other gases. As displayed in Fig. 7, apparent permeability coefficients measured over a large (Terzaghi) effective stress \(\sigma_{{{\text{eff}}}}\) range (10–40 MPa) show very good agreement and virtually no stress dependence irrespective of the mean pore pressure. This clearly proves that the ceramic material is resistant against mechanically induced stress. Thus, the interference of poro-elastic effects, commonly observed on natural rock samples, can be excluded. All permeability phenomena observed in this study can therefore be attributed unambiguously to fluid-dynamic effects.

Fig. 7
figure 7

Results of stress-sensitivity analysis of permeability coefficients measured with helium on the synthetic ceramic disk. No stress sensitivity can be observed and apparent permeability coefficients at a given mean pore pressure are identical over a large range of Terzaghi effective stress levels (10–40 MPa)

The plot of the Darcy pressure drop versus mass flow rate/flux in Fig. 8a shows a linear trend and no indication of deviations from viscous flow conditions. In Fig. 8b, the measured apparent gas permeability coefficients are plotted as a function of differential pressures at the largest mean pressure (30.5 MPa) used in these tests. The apparent permeability coefficients scatter around a value of 27.72 µD and do not decrease with increasing differential pressure (increased flow rates), as would be expected for inertial or turbulent flow effects. The absence of inertial (turbulent) flow effects is furthermore supported by a low Reynolds number, ranging between 10–9 and 10–5 (Table 3) far off the laminar-turbulent transition as defined by the critical Reynolds number range (from 1 to 2300, Javadi et al. 2014).

Fig. 8
figure 8

Test for non-Darcy turbulent flow. Turbulent flow can be excluded because a the Darcy pressure drop versus mass flux shows no deviation from linearity and b apparent gas permeability (He) does not decrease with increasing differential pressure at the largest mean pore pressure (30.5 MPa)

Table 3 Ranges of Reynolds numbers under the experimental conditions for all gases tested

4.2.3 Apparent Gas Permeability Data

Overall, the apparent gas permeability coefficients of the ceramic filter disk range between 20 and 530 µD at an effective stress of 20 MPa. The plots of the apparent permeability coefficients versus the reciprocal mean gas pressure (“Klinkenberg plot”) are linear for all gases and have very similar y-axis intercepts (Fig. 9a). However, slopes of the Klinkenberg trends differ significantly and decrease in the order He > H2 > Ar ≈ N2 ≈ Ar–N2 > CH4 > CO2. This corresponds directly to the order of the mean free path lengths \(\lambda\) of the gases at a given pressure (Hirschfelder et al. 1964; Bird 1983).

Fig. 9
figure 9

a Klinkenberg plots of the apparent permeability coefficients of different gases at the same effective stress of 20 MPa. Dashed lines represent the range of pore pressures used for extrapolation to the Klinkenberg-corrected permeability (y-axis intercept). Blue shaded area represents the transition between slip and transitional flow regime (Kn = 0.1; see also Table 4). b Enlarged view of a showing the deviant CO2 permeability trend. c Flow regimes in terms of ranges of Knudsen numbers covered in this study

Klinkenberg-corrected permeability coefficients \(k_{\infty }\) range from 20 to 25 µD (Table 4). Evidently, the \(k_{\infty }\) value measured with helium (25.04 µD) is significantly higher than those measured with all other gases (19.98–23.29 µD). This “He-anomaly” will be discussed in Sect. 4.3.

Table 4 Summary of characteristic parameters derived from permeability measurements on a synthetic nanoporous ceramic disk with different gases at 20 MPa effective stress

The apparent permeability coefficients measured with CO2 deflect significantly from a linear trend at \(1/P_{{\text{p}}}\) values < 0.2 MPa−1 (Fig. 9b). For the evaluation of intrinsic permeability coefficients in this study, only the linear portion of the CO2 trend was taken into account. The \(k_{\infty }\) obtained from this extrapolation is considerably lower (19.98 µD) than those obtained with the other gases. All CO2-related phenomena will be discussed in Sect. 4.2.6.

The gas slippage factor \(b\) for each gas was determined using Klinkenberg’s conventional first-order slip equation (Eq. 8). Helium exhibits the highest b-factor with 3.47 MPa and CO2 the lowest b-factor with 0.86 MPa (Table 4) as expected from the sequence of mean free path length \(\lambda\) at a given pressure.

Due to the large range of pore pressures used in this set of flow experiments (0.2–30.5 MPa), the permeability data covers the transition between slip and transitional flow regime. The corresponding Knudsen numbers range between 0.005 and 1.7 (Fig. 9c and Table 4). They were calculated according to Eqs. 1 and 3 using the average of the equivalent pore diameters from the MIP and LPNS measurements (Table 2).

4.2.4 Second-Order Gas Slip Analysis

The blue shaded area in the Klinkenberg plot (Fig. 9a) indicates the transition from slip flow to transitional flow regime around \({\text{Kn}}\) = 0.1 for the different gas types (Table 4). The linear regression of the experimental data by the first-order Klinkenberg equation has a coefficient of determination (R2) of 1.00 (Fig. 9a), i.e., there is no evidence for a deviation from a linear trend. This was further confirmed by applying the second-order slip function, which essentially results in \(b_{1}\) coefficients nearly identical to the first-order \(b\) values (\(b\) and \(b_{1}\) in Table 4) and \(b_{2}\) coefficients of essentially zero (\(b_{2}\) in Table 4).

In conclusion, deviations from the conventional first-order slip model, as reported by several studies (e.g., Maurer et al. 2003; Ewart et al. 2007; Graur et al. 2009) on micro- to nano-tubes and channels at sub-atmospheric pressures (< 0.1 MPa), could not be observed within the scope of this study on a nanoporous ceramic material at high pressures (see \({\text{Kn}}\) range in Table 4). Thus, second-order flow regime effects from the slip to transitional flow regime are negligible and do not have to be considered to describe permeability behavior even at comparatively high \({\text{Kn}}\) up to 1.

4.2.5 Transport Porosity

In studies on the characterization of porous substances by MIP (Purcell 1949) or sorption and capillary condensation (e.g., Barrett et al. 1951), simple capillary bundle models have been used to develop concepts, classify experimental data and standardize methods. Thus, the standard evaluation procedure of MIP is based on the model of cylindrical capillaries with non-uniform radii that are successively filled with mercury. Such “equivalent” models also form the basis of fluid transport concepts such as the Kozeny–Carman equation (cf. Carman 1956). It is generally understood that they constitute over-simplifications and cannot account for phenomena like, for instance, hysteresis. Nevertheless, they are useful for comparison and classification purposes.

In this study, we used a simplified capillary bundle model to estimate the transport porosity of the nanoporous disk. Combining the macroscopic Darcy law (Eq. 5) and the microscopic Hagen–Poiseuille formulation (Eq. 4) of fluid transport through an arrangement of \(n_{i}\) capillaries of diameter \(d_{i}\) distributed over a cross section area \(A\), one obtains:

$$k_{i} = \frac{{n_{i} \left( {\frac{{d_{i} }}{2}} \right)^{4} \pi }}{8 A}.$$
(19)

Here, \(k_{i}\) is the proportion of permeability contributed to the total permeability of the sample by the \(n_{i}\) capillaries (“pores”) of diameter \(d_{i}\).

Evidently,

$$\phi_{i} = \frac{{n_{i} \left( {\frac{{d_{i} }}{2}} \right)^{2} \pi }}{A} = \frac{{n_{i} r_{i}^{2} \pi }}{A}$$
(20)

is the proportion of porosity contributed to the total porosity of the sample by this set of capillaries. Equation 21 can thus be rewritten as:

$$k_{i} = \frac{{ \left( {\frac{{d_{i} }}{2}} \right)^{2} \phi_{i} }}{8 }.$$
(21)

Due to the very narrow (equivalent) pore size distribution of our sample, a uniform transport pore diameter of 68.5 nm (mean value of equivalent diameters from MIP and LPNS) can be assumed. With the mean value of the intrinsic permeability coefficient \(k_{\infty }\) of 23 µD, Eq. 21 yields a value of 15.7% for the transport porosity.

Taking explicitly into account the pore size distribution from the MIP measurement, this computation can be refined. For this purpose, the experimentally measured MIP pore size distribution was discretized into 1 nm steps in the pore diameter range from 40 to 120 nm. For each set of the 81 sets of capillaries, the \(k_{i}\) value was calculated according to Eq. 21 and the \(\phi_{i}\) variables adjusted such that the sum of the \(k_{i}\) values matched the observed intrinsic permeability coefficient \(k_{\infty }\). The sum of the \(\phi_{i}\) values obtained by this procedure amounted to 13.5% and was thus slightly lower than the value obtained assuming a uniform pore diameter. This is due to the disproportionate contribution of larger pores to the overall permeability coefficient.

The transport porosity values thus obtained are in the same order as the 10.2% transport porosity derived from the comparison of the intrusion and extrusion curves (hysteresis) of the MIP measurement, indicating that the simple models used here yield (qualitatively and semiquantitatively) consistent results.

It can be observed that the most prominent pore diameter of 68.5 nm (MIP, LPNS) does not reflect the maximum permeability contribution (Fig. 10). Instead, the maximum is located at a slightly greater pore size of 73 nm. This shows that the maximum in the pore size distribution may not represent the maximum of the flow contribution, which is shifted to larger pore sizes. Therefore, using MIP data as input to estimate the flow regime based on the Knudsen number is likely to be in error and result in too high Knudsen numbers and overestimation of the slip flow contribution (Moghaddam and Jamiolahmady 2016).

Fig. 10
figure 10

Permeability contribution of pore sizes based on the MIP and LPNS pore size distribution. The orange curve represents the incremental permeability contribution per nm, whereas the blue curve displays the cumulative sum of all capillary bundles with increasing size. The maximum cumulative permeability corresponds to the experimentally deduced Klinkenberg-corrected permeability k (23 µD). The dotted line indicates the mean of the most prominent equivalent pore diameter (68.5 nm) derived from LPNS and MIP and the dashed line indicates the capillary pore diameter (73 nm) with the highest incremental permeability contribution

4.2.6 Non-ideal Gas Effects of CO2 on Permeability Near the Critical Point

In contrast to the other gases that essentially show a linear Klinkenberg trend, tests with CO2 exhibit a permeability peak around the critical pressure of 7.38 MPa (1/Pp = 0.14 MPa−1) (Fig. 9b). This nonlinearity coincides with the non-ideal nature and dramatic change in thermodynamic properties of CO2 in this pressure range. CO2 is supercritical above its critical pressure and temperature of 7.38 MPa and 30.98 °C, respectively. Storage- and transport-related properties, such as fluid density, isothermal compressibility and viscosity, are not well defined close to the critical point (Span and Wagner 1996) and are very sensitive to pressure and temperature fluctuations.

However, isothermal compressibility \({c}_{\mathrm{g}}\) and dynamic viscosity \(\mu\) are both required for the computation of permeability coefficients from non-steady permeability tests and inaccuracies will have direct impact on the determination of CO2 permeability coefficients (Eq. 18). Viscosity shows a strong increase with pressure around the critical point. However, the viscosity change is poorly correlated with permeability and can, therefore, not be responsible for the observed permeability deviations (Fig. 11a). In contrast, the trends in apparent permeability and isothermal compressibility coefficients are very similar (Fig. 11b), indicating that compressibility dominates the nonlinear deviations of the permeability trend.

Fig. 11
figure 11

Zoom-in on the Klinkenberg plot of CO2 permeability at 35 °C versus a dynamic viscosity and b isothermal compressibility, respectively. Note: permeability and compressibility peaks are located near 1/Pp of 0.14 MPa−1 (equivalent to Pp of 7.38 MPa, critical pressure)

To demonstrate the influence of compressibility input data on the permeability trend, compressibility values used for permeability evaluation were adjusted to obtain a linear Klinkenberg trend. The basic assumption behind this approach is that changes in apparent permeability are solely governed by the first-order slip flow equation and that deviations are due to inaccuracies (uncertainties) in the compressibility values.

For this approach, only the low-pressure experimental permeability data (\({P}_{\mathrm{p}}\hspace{0.17em}\)< 1 MPa) were used to construct a Klinkenberg diagram by linear regression (Fig. 12a). Based on this regression, the compressibility coefficients required for permeability calculation (Eq. 18) were adjusted over the entire experimental pore pressure range to match the regression line. The relative deviation between the original (at 35 °C) and the adjusted compressibility coefficients has a pronounced peak around the critical region (Fig. 12b). Comparing this deviation with the influence of temperature variations on the compressibility coefficients in this temperature and pressure range shows that a temperature error of less than 1 °C is already sufficient to cause the observed offsets in the CO2 permeability data around the critical region.

Fig. 12
figure 12

a Klinkenberg plot displaying the original, experimental CO2 permeability data (at 35 °C) over the entire pressure range and linear regression line for 1/Pp > 1 MPa−1 (black dashed line). b Relative deviation of CO2 compressibility (normalized to the experimental temperature of 35 °C) as a function of fluid pressure. Isotherms of 34 and 36 °C (dashed lines) represent lower and upper bound, respectively. Adjusted compressibility values (black circles) divert quite significantly from the experimental compressibility (black line, 35 °C) around the critical region

However, two outliers at \({P}_{\mathrm{p}}\hspace{0.17em}\)> 12 MPa cannot be readily explained as the adjusted compressibility values are unrealistic. This is further supported by the Klinkenberg plot (Fig. 12a) where these two apparent permeability points reveal an uncharacteristic and steep downward trend. A possible explanation is discussed below.

Evidently, permeability measurements close to the critical point involve uncertainties related mainly to the strong fluctuations in isothermal compressibility and, to a lesser extent, viscosity. This problem is amplified by insufficient temperature control, i.e., in our experiments, the temperature was not directly measured within the capillary tubes and close to the sample. Therefore, rapid temperature fluctuations within the gas stream, by, e.g., the Joule–Thomson effect, could not be detected. A decrease in temperature can lead to phase changes and the formation of immobile liquid films that can affect permeability measurements. This might explain the unrealistic deviations of the two outliers at pressures above 12 MPa.

4.3 Helium Anomaly

Differences between the Klinkenberg-corrected permeability coefficients (\({k}_{\infty }\)) of rocks measured with helium as compared to other gases (“helium anomaly”) have been commonly observed in previous studies (Sinha et al. 2013; Ghanizadeh et al. 2014a, b; Gensterblum et al. 2014; Fink et al. 2017a; Shabani et al. 2020). Permeability coefficients measured with helium can be several times higher than those measured with other gases (Fig. 13). While for the majority of samples deviations are less than a factor 2, factors of up to 63 were observed in extreme cases, e.g., for the TOC-rich (45% TOC) Kimmeridge clay (Gaus 2020). For the essentially inert nanoporous ceramic material used in this study, the Klinkenberg-corrected permeability coefficient \({k}_{\infty }\) for helium was also found to be consistently higher (by a factor of 1.13) than those measured with other gases at the same conditions (Table 4). Generally, the difference between permeability coefficients measured with helium and other gases tends to increase with decreasing permeability. This holds true for shales with the exception of carbonate-rich samples from Iran (Shabani et al. 2020).

Fig. 13
figure 13

Log-linear plot of the intrinsic helium permeability normalized to other gases (argon, nitrogen or methane) vs. the intrinsic helium permeability coefficients \({k}_{\infty }\). The dashed line represents the case that the permeability of helium and the other gases are identical. The presented samples include shales (circles), carbonates (diamonds) and the synthetic ceramic disk (red X)

This effect is probably controlled by a complex interaction of multiple rock- and fluid-specific properties, such as rock composition, pore/matrix structure (e.g., pore size, fracture, cleats, cementation), mechanical stress behavior (e.g., rigid framework or creep effects), maturation, fluid type (e.g., sorbing gas, polarizability) and sorption/swelling effects (e.g., or organic matter and/or clay). Figure 13 clearly shows that elevated helium permeability values persist throughout various rock types and permeability ranges.

As opposed to natural rocks, sorption and swelling effects of the ceramic disk can be excluded due to its composition (α-Al2O3, corundum) that is considered non-reactive.

Flow regimes of gases depend on the pore size with regard to the mean free path lengths of gas molecules. Molecular sieving effects can occur in micropores (< 2 nm, IUPAC classification) due to selectivity and rejection of certain molecules based on their molecular size (George and Thomas 2001). Among the gases used in this study, helium and hydrogen have the smallest molecular diameters of 0.258 and 0.297 nm, respectively (Hirschfelder et al. 1964). Nevertheless, the discrepancy between helium and hydrogen permeability coefficients cannot be readily explained by accessibility or molecular sieving effects since the ceramic disk does not have micropores. As determined by MIP and LPNS, the smallest pores are around 40 nm and therefore more than a hundred times larger than the kinetic diameters of the molecules.

Another aspect of this study was the examination of dynamic gas viscosity data. In terms of the determination of permeability coefficients, viscosity is the only input parameter taken from literature databases such as NIST (Lemmon et al. 2021). The NIST database itself is based on models that employ an extensive collection of viscosity data. As viscosity can be determined in different ways (e.g., transpiration method, rotary viscometer), the objective was to compare and verify viscosity data for potential anomalies. As shown in Sect. 4.1, after fitting the capillary diameters, the overall consistency of literature and experimental viscosity data for all gases used in this study could be verified for pressures above approximately 0.1 MPa. Deviations below this threshold can readily be attributed to slip flow effects that increase with decreasing pressure. Gas viscosity measured under conditions where slip flow occurs is referred to as “apparent viscosity,” as it does not reflect the true viscosity. Despite the limited accuracy of our viscosity measurements in this pressure range, it can be stated that the helium viscosity values measured with the 2 µm capillary were consistently higher than the expected/predicted values and further investigations are required to resolve this problem.

In summary, in these tests with an inert nanoporous medium the “helium anomaly” could clearly be observed. Compared to most natural rocks, however, the effect is much less pronounced here, indicating that other influential factors, such as adsorption, molecular sieving (size exclusion) and secondary effects such as swelling, play an important role in those instances.

4.4 Implications of Fluid Flow Tests on Analogue Porous Media for Tight Rocks

The aim of this study was to improve the understanding of fluid transport mechanisms in tight rocks and give recommendations for future geotechnical applications by using artificial porous media as rock analogues.

Based on its petrophysical properties (mean pore diameter 68.5 nm; intrinsic permeability coefficient ~ 23 µD), the ceramic disk is a good proxy for tight rocks, such as tight sandstones and shales with comparable transport pore sizes. The gas slip factor, a measure of the average transport pore size, determined with helium is 3.47 MPa and thus well within the range of values for shales reported in other studies. Slip factors range mostly from 1.5 to 6 MPa with the highest published value of 8.6 MPa (Heller et al. 2014; Ghanizadeh et al. 2014a, b; Fink et al. 2017a, 2018; Letham and Bustin 2018; Gaus et al. 2019; Nolte et al. 2019; Shabani et al. 2020). Of course, the results obtained with this “ideal” system should only be applied to natural tight rock systems that have similar transport properties.

The permeability tests on this specimen have shown that the first-order Klinkenberg model is continuously valid from the slip to transitional flow regime; second-order slip effects are practically nonexistent and can be neglected for Knudsen numbers up to at least 1. Although second-order slip effects have been measured on micro- to nano-tubes and channels (Arkilic et al. 1997; Maurer et al. 2003; Roy et al. 2003; Ewart et al. 2007; Gruener and Huber 2008; Graur et al. 2009; Yamaguchi et al. 2011), these concepts cannot simply be applied to tight rocks and high gas pressures (Javadpour 2009; Civan 2010b). Attempts to interpret permeability results by means of second-order slip extensions (Moghaddam and Jamiolahmady 2016) were based on permeability trends measured with nitrogen on shale samples. These exhibited “concave downward deviations” toward lower pore pressures (min. 1.72 MPa) and higher Knudsen numbers (max. 0.56). However, as discussed in previous contributions (Moghaddam 2018; Fink et al. 2018), the deviations are quite likely due to stress rather than fluid-dynamic effects. This is further supported by our results where the permeability trend measured (with N2) on the synthetic nanoporous specimen remains linear for Knudsen numbers up to 0.69. The recurring issue with tight rocks, especially shales, is the proper and challenging separation of interfering effects (e.g., mechanical stress, sorption/swelling), which makes reliable fluid-dynamic studies on tight rocks, covering the slip to transitional flow regime, very difficult. Based on the results of this study, we argue that second-order slip models are neither required nor meaningful in the evaluation of experimental gas permeability data and in modeling of gas transport in nanoporous shales, as long as no clear evidence for additional slip flow parameters exist.

Furthermore, the determination of reliable Knudsen numbers and flow regimes in tight rocks is hardly possible, as they often exhibit multimodal pore size distributions so that the “true” and representative pore size, which is representative for the averaged, macroscopic fluid flow behavior, is unknown. Especially in gas/oil shales, small pores (2–5 nm) mostly dominate in quantity but comparatively larger pores (20–30 nm) are will dominate the flow contribution (Zendehboudi and Bahadori 2016).

However, in a geotechnical context, it can be instrumental to estimate fluid flow regimes in order to predict production rates for example. Figure 14 provides a simplistic estimation of expected flow regimes based on the average pore diameter and pore pressure.

Fig. 14
figure 14

Average pore diameter versus pore pressure for methane at experimental (35 °C) and reservoir (120 °C) temperatures to indicate the covered range of fluid flow regimes by means of the Knudsen number (Kn) isolines. Red line represents the experimental range of this study (68.5 nm; 0.2 to 30.5 MPa) at 35 °C

For a tight rock under reservoir conditions (120 °C) with an average flow-dominated pore diameter of 10 nm, pore pressures of only 8.5 MPa are required to reach a Knudsen number of 0.1, which designates the boundary between slip and transitional flow regime. Further, only 0.85 MPa of fluid pressure are required to reach a Knudsen number of 1 that, based on this study, does not induce flow regime effects or any changes in flow behavior.

In most tight reservoirs, rocks pore pressures will exceed a value of 8.5 MPa. Typical depths for tight rock reservoirs can range from 1 to 5 km (Zendehboudi and Bahadori 2016). Reservoir pressures can be deduced as a function of depth by assuming a hydrostatic pressure gradient 10 MPa/km. Consequently, expected pore pressures for that range of depth are 10–50 MPa. Therefore, it can be stated that fluid-dynamic or second-order effects can be neglected for most tight rock plays in pores > 10 nm.

Permeability tests with CO2 on the dry ceramic disk have emphasized the complexity of non-ideal gas behavior. These permeability effects described in this study were measured on an “ideal” sample where poro-elastic stress, moisture effects or sorption effects were absent. Even for such an idealized system, the interpretation of the CO2 permeability data was challenging due to thermodynamic effects close to the critical point.

However, geologic (“natural”) systems of interest (e.g., for CO2 sequestration) are always mechanically stressed (lithostatic pressure), mostly have significant compressibility, contain water and often operate near the critical pressure and temperature of CO2. Therefore, caution is advised when interpreting fluid flow data measured on natural rock systems that are controlled by complex interactions of multiple rock- and fluid-specific properties.

5 Conclusion

This systematic fluid transport study with different gases was performed on homogeneous, artificial porous media that serve as simplified tight rock analogues to compare pore structure characterization methods and investigate fluid-dynamic phenomena that occur in the laboratory and nature. With this approach, we aimed at reducing the amount of unknown parameters and the interference of different processes (e.g., mechanical deformation, swelling and sorption) that are usually superposed. The monomineralic ceramic disk (> 99% α-Al2O3, corundum) used in these tests fulfills these requirements as it has a homogeneous pore structure and is fully resistant toward mechanically induced stress (poro-elasticity).

Results from pore structure information obtained by different methods show that helium pycnometry, MIP and LPNS are in good agreement in terms of porosity (~ 28%) and the most prominent pore size (~ 68.5 nm). This contrasts with observations for natural tight rocks where differences between data derived from such measurements are regularly observed. Transport porosity is estimated to range between ~ 10 and 13%, which amounts to approximately one-third of the total porosity. One should keep in mind that this material was optimized for filtration and high fluid transport efficiency and transport porosities of natural rocks may be by orders of magnitude smaller. BIB-SEM analysis of total porosity showed significant deviations from the true porosity and most prominent pore size, but the frequency distribution compares reasonably well with MIP and LPNS. Further work is planned by means of liquid metal injection (LMI) followed by BIB-SEM.

The “helium anomaly” (intrinsic permeability coefficients of He > other gases) is significantly less pronounced for the artificial ceramic compared to most tight rocks indicating that it is described by a complex interaction of multiple rock- and fluid-specific properties (e.g., sorption effects). Overall, the intrinsic permeability coefficient is ~ 23 µD (2.3 · 10–17 m2).

The CO2 permeability trend exhibits a nonlinearity (Klinkenberg plot) around its critical pressure (7.38 MPa) which coincides with the drastic change in thermodynamic properties of supercritical CO2 (e.g., density, compressibility, viscosity). Within this transition, the thermodynamic properties are highly sensitive toward temperature and pressure fluctuations causing significant error that cannot be fully avoided. This needs to be considered during interpretation of CO2 flow data. Moreover, as fluid flow experiments with CO2 in “dry” condition are already challenging close to the critical point, we expect that the interpretation of data from “multi-phase” systems (e.g., for modeling of CO2 flow in the subsurface) is even more biased.

For modeling of gas flow in tight rocks (e.g., gas shales), second-order slip flow effects are regularly implemented within the transitional flow regime based on theoretical considerations. However, we obtained linear (first order) apparent permeability trends for various gases (He, Ar, N2, Ar–N2 mix, CH4, H2) over an extremely wide range of pore pressures (0.2–30.5 MPa). This experimentally derived linear and continuous “transition” between slip and transitional flow regime (Kn from 0.001 to 1) shows that second-order slip effects are not necessarily relevant to describe gas flow in these types of porous media. This is also of relevance for shale gas systems because similar Knudsen numbers in the range of 0.1–10 (transitional flow regime) can also be realized at elevated reservoir pressures if pores are smaller compared to our experiments (e.g., < 10 nm) (\({\mathrm{Kn}} \sim \frac{1}{{d\!\,\cdot\,\!P}_{\mathrm{p}}}\)).