Dynamic pore network modelling of real gas transport in shale nanopore structure

Journal of Petroleum Science and Engineering 184 (2020) 106506

Contents lists available at ScienceDirect

Journal of Petroleum Science and Engineering journal homepage: www.elsevier.com/locate/petrol

Dynamic pore network modelling of real gas transport in shale nanopore structure

T

Wenhui Songa, Yao Juna,∗, Dongying Wanga, Yang Lib,∗∗, Hai Suna, Yongfei Yanga a b

School of Petroleum Engineering, China University of Petroleum, No. 66, Changjiang West Road, Huangdao District, Qingdao, 266580, China Department of Oilfield Exploration & Development, Sinopec, Beijing, 100029, China

A R T I C LE I N FO

A B S T R A C T

Keywords: Shale gas Transport mechanisms Dynamic pore network model Gas permeability

Gas transport in shale nanopores is controlled by the complex transport mechanisms and pore structure characteristics. So far much work has been done on the single pore and regular structure based pore network but little is known on real gas transport behavior on a realistic shale pore structure model. In this work, a dynamic pore network model is proposed to describe single component methane transport in nano-scale porous media. Gas transport behavior takes transport mechanisms of slip flow, transition flow, surface diffusion and ad/desorption into account. Real gas effect under high pressure and temperature is considered when calculating gas properties. A three dimensional pore network model is built from a three dimensional image that is reconstructed using multi-point statistics from a small organic-rich area on a two dimensional shale SEM image. This pore network model is used to analyze dynamic gas transport and pressure drop transmission process. Our simulated results reveal that there exists a time period that gas permeability on the pore network is influenced by the pressure drop transmission process and pore structure that flows on. When pressure drop reaches outlet, gas transport turns steady state. At steady state period gas permeability becomes constant and is influenced by the effective stress, pressure, temperature and shale rock property.

1. Introduction Gas transport in nano-scale porous media has drawn more and more attention because of recent success in shale gas development (Middleton et al., 2017; Yuan et al., 2015; Zou et al., 2016). Gas transport in shale nanopores is controlled by the complex transport mechanisms and pore structure characteristics. So far much work has been done on the single pore and regular structure based pore network but little is known on real gas transport behavior on a realistic shale pore structure model. The interaction between gas molecules and pore wall can not be neglected and the conventional Darcy law becomes invalid to describe shale gas transport because of the nanoscale shale pore size (Ho and Tai, 1998; Song et al., 2016; Yao et al., 2013). Knudsen number is generally given by the molecular mean free path divided by the pore radius. The corresponding gas transport regime in Fig. 1 differs at various Knudsen numbers. The organic matter is widely distributed in the shale formation and serves as the adsorbed gas accumulation place (Clarkson et al., 2013; Li et al., 2019; Saif et al., 2017; Sun et al., 2019; Wang et al., 2016a). Large amount of gas molecules are gathered inside the nanoporous organic matter with heterogeneous

∗

pore structure. The adsorbed gas molecules transport in the forms of surface diffusion along the pore surface caused by the molecular concentration gradients (Kim et al., 2017; Sun et al., 2015; Wang et al., 2016b). Attempts have been made to study the detailed shale gas transport behavior and assess the gas transport ability. A single pore gas transport model proposed by (Beskok and Karniadakis, 1999) covers flow regimes from continuum fluid flow to free molecular flow. The second order slip permeability solution for single pore gas transport under high Knudsen number condition were studied by (Ziarani and Aguilera, 2012) and (Wang et al., 2017b). The gas transport model in shale nanopores established by (Civan, 2010; Civan et al., 2013, 2011) considers the pore structure characteristics and slippage effects. Later surface diffusion of adsorbed gas and slippage of free gas are both considered (Xiong et al., 2012). On the other hand, a series of single pore gas transport models were developed by superimposing slip flow and Knudsen diffusion (Darabi et al., 2012; Javadpour, 2009; Javadpour et al., 2007; Singh, 2013; Singh and Javadpour, 2016). The gas transport models for shale porous media developed by (Song et al., 2016) accounts for different pore type based on the equivalent hydraulic radius with a given

Corresponding author. Corresponding author. E-mail addresses: [email protected] (J. Yao), [email protected] (Y. Li).

∗∗

https://doi.org/10.1016/j.petrol.2019.106506 Received 1 April 2019; Received in revised form 8 June 2019; Accepted 18 September 2019 Available online 18 September 2019 0920-4105/ © 2019 Elsevier B.V. All rights reserved.

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

Fig. 1. Gas flow regimes at various Knudsen numbers.

adopted in a wide range of previous studies (Cao et al., 2016; Sheng et al., 2018; Wang et al., 2018; Wu et al., 2016).

porosity and tortuosity. The shale pore structure heterogeneity is neglected in these single pore models and reduces the accuracy of predicting gas transport ability. Fractal theory has been applied to model gas transport by taking the critical shale pore structure parameter into account (Kou et al., 2009; Song et al., 2019). The ideal gas slippage fractal model was initially established using the simplified Beskok et al. model (Zheng et al., 2013). The gas slippage fractal model for real gas transport is derived based on the Extended Navier-Stokes Equation (ENSE) (Geng et al., 2016). The gas slippage model specifically for the nanoscale shale porous media is developed by considering both the second-order slip for free gas and surface diffusion for adsorbed gas (Song et al., 2018b). However, the pore size spatial variation and porethroat connectivity are neglected in the fractal theory. In order to assess the influence of the spatial pore structure on gas transport ability, pore network model has been applied to analyze gas transport in porous media, but current research on gas transport through pore network is mostly based on conceptual pore network model and only free gas with constant property is considered (Huang et al., 2016; Ma et al., 2014; Mehmani et al., 2013; Song et al., 2017; Zhang et al., 2015). The recent development of high resolution imaging technology and digital rock physics (DRP) approach have made it feasible to visualize the three dimensional pore structure and analyze the detailed pore structure characteristic (Bai et al., 2013; Curtis et al., 2010; Loucks et al., 2009; Ostadi et al., 2010). Numerical simulations can be conducted on the imaging sample to predict gas transport property and mechanical properties (Blunt et al., 2013; Jiang et al., 2007; Knackstedt et al., 2006; Oren et al., 1998). The objective of this work is to build a dynamic shale pore network model to study the dynamic single component methane transport process on the realistic pore structure and account for the influences of gas transport mechanisms, adsorption/desorption, real gas effect and gas critical property change. To the best knowledge of the authors, this is the first time dynamic gas transport process is studied in detail based on the pore network extracted from the realistic pore structure.

kins =

⎛1 − φ2 r22 = ⎝ 2 φ1 r1 ⎛1 − ⎝

3 σeff 2 m

( ) ⎞⎠ ( ) ⎞⎠ p1

3 σeff 1 m p1

(4)

Porosity is directly proportional to the third power of pore radius and Eq. (4) can be reorganized as (Tiab and Donaldson, 2015): 5

r2 r15

⎛1 − = ⎝ ⎛1 − ⎝

3 σeff 2 m

( ) ⎞⎠ ( ) ⎞⎠ p1

3 σeff 1 m p1

(5)

Eq. (5) can be rearranged into the following expression:

⎛1− r = r0 ⎜ ⎜1 − ⎝

m

( ) ( ) σeff p1

σeff 0 p1

⎞ ⎟ m ⎟ ⎠

3 5

(6)

Previous study shows that the shale nanopore shape exhibits different pore morphology (Afsharpoor and Javadpour, 2016; Song et al., 2018c). The circle pore and square pore shown in Fig. 2 are the most commonly observed pore shape during shale pore structure characterization (Yang et al., 2017). Therefore we consider the gas flow conductance in these two different pore shape. The absorbed gas molecules reduce the pore space available to the free gas molecules transport. The surface coverage of the absorbed gas molecules on the pore wall can be defined in Eq. (7) based on the Langmuir monolayer adsorption.

We first investigate the influence of stress on the single shale nanopore size. Gangi's model (Gangi, 1978; Kwon et al., 2004) is applied to correlate the gas permeability with effective stress: m 3

(1)

σeff = pcon − αs p

(3)

Where ϕ is the porosity, τ is the tortuosity. In this work, we use the change of effective hydraulic radius at different stress conditions to study gas permeability change at different confining stress. When the effective stress changes from σeff1 to σeff2, intrinsic permeability changes from kins1 to kins2. The limitation of the effective hydraulic radius model is that it only works for isotropic confining stress condition and the contribution of stress condition to the calculated permeability may not be accurate if the confining stress is not loaded isotropically. More accurate simulation of the pore structure change at different confining stress will be considered in future work. Assuming the tortuosity is constant value (Wu et al., 2016), the relationship between shale nanopore radius and effective stress can be derived according to Eq. (1) and Eq. (3):

2. Gas flow conductance in single shale nanopore

σeff km = k 0 ⎛⎜1 − ⎜⎛ ⎟⎞ ⎟⎞ ⎝ p1 ⎠ ⎠ ⎝

φr 2 8τ

(2)

Where k0 is the permeability at zero effective stress (μm ), σeff is the effective stress on the shale nanopores (MPa), m is associated with the surface roughness of pores, p1 is the effective stress when the shale nanopores are closed completely (MPa), pcon is the confining pressure (MPa), αs is the effective stress coefficient, p is the pore pressure (MPa). The relationship between the nanopore radius and intrinsic permeability can be expressed in Eq. (3) (Civan et al., 2013), which has been 2

θ=

p/Z pL + p / Z

(7)

Where Z is the gas compressibility factor, PL is the Langmuir pressure (MPa). The critical pressure Pc and critical temperature Tc change with the pore size when the pore size is less than 10 nm. Islam's equations (Islam et al., 2015) are applied to correlate the critical temperature and 2

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

(21)

Y = 2.447 − 0.2224X

The unified Hagen–Poiseuille-type equation established by Beskok et al. is applied to describe free gas transport in circle pores. The free gas volumetric flux for circle pore can be expressed as (Beskok and Karniadakis, 1999):

q = f (Kncir )

4 πreff ΔP 8μ l

(22)

Where l is the pore length, △P is the pressure drop, f(Kn) is the flow condition function and can be given by (Beskok and Karniadakis, 1999):

4Kncir ⎞ ⎟ f (Kncir ) = (1 + α cir Kncir ) ⎛⎜1 + 1 − βKncir ⎠ ⎝

Fig. 2. Illustration of different pore geometry in this study.

The value of slip coefficient β (−1) results in an accurate model of the velocity distribution for a wide range of Knudsen number based on the Direct simulation Monte Carlo (DSMC) and linearized Boltzmann solutions (Karniadakis et al., 2006). The influence of slip coefficient on gas transport will be discussed in Section 4. αcir is the dimensionless rarefaction coefficient and can be written as:

pressure in nano-pores with pore size in Eq. (8) and Eq. (9).

Tc =

8 ⎡ σ σ a − 2σ 3εN 2 ⎛2.6275 − 0.6743 ⎞ ⎤ 27bR ⎣ r⎝ r ⎠⎦

(8)

Pc =

8 ⎡ σ σ a − 2σ 3εN 2 ⎛2.6275 − 0.6743 ⎞ ⎤ 27b2 ⎣ r⎝ r ⎠⎦

(9)

α cir =

Where a is the vdw energy parameter (Pa dm6/mol2), b is vdW energy parameter (m3/mol), σ is the Lennard-Jones size parameter (m), ε is the Lennard-Jones energy parameter. The gas compressibility factor under the high pressure condition is calculated by (Mahmoud, 2014):

Ppr =

p Pc

(10)

Tpr =

T Tc

(11)

(12) The effective capillary radius in consideration of adsorbed gas layer can be expressed as:

gfree _ cir =

Kncir =

λ reff

(14)

Knsqu =

λ 2reff

(15)

Q = C (AR)

K=

pM ZT

(9.379 + 0.01607M ) T1.5 (209.2 + 19.26M + T )

X = 3.448 +

986.4 + 0.01009M T

wh3 ΔP 6Kn ⎞ ⎛ ⎞ (1 + αduct Kn) ⎜⎛1 + ⎟ 12μ ⎝ l ⎠ 1 − βKn ⎠ ⎝

192(AR) ⎡ C (AR) = ⎢1 − π5 ⎣

(26)

∞

∑ i = 1,3,5...

tan(iπ /2(AR)) ⎤ ⎥ i5 ⎦

(27)

The duct width equals to duct height for square pore case and the dimensionless rarefaction coefficient αduct can be expressed as (Karniadakis et al., 2006):

(16)

μ = (1 × 10−4) K exp(XρY )

(25)

Where w is the duct width (m), h is the duct height (m), αduct is the dimensionless rarefaction coefficient, AR is the aspect ratio, C (AR) is the correction factor and can be given as (White and Corfield, 2006):

Where M is the gas molecular weight (g/mol), R is the ideal gas constant, μ is the gas viscosity (Pa·s) and can be given by Lee et al. (1966) method, which has been adopted for predicting gas viscosity in confined pore space (Bui et al., 2016; Kim et al., 2016; Landry et al., 2016; Wang et al., 2017a):

ρ = 1.4935 × 10−3

8μl

The volumetric gas flux in a duct can be expressed as (Karniadakis et al., 2006):

Where λ is the mean free path and can be expressed as:

πZRT μ 2M p

(24)

4 πreff f (Kncir )

(13)

Where r is the inscribed circle radius on the pore cross section (m). Therefore Knudsen number in circle pore Kncir and square pore Knsqu can be respectively given as:

λ=

128 tan−1 [4.0Kncir 0.4] 15π 2

The pressure drop on the pore network is very small based on the typical pressure gradient at subsurface conditions. Given the definition of Knudsen number in Eq. (14) and Eq. (15), the change of Knudsen number with the pressure drop along the conduits can be neglected and Knudsen number is mainly dependent on the pore radius. As the gas hydraulic conductance depends on the pore radius and Knudsen number change, therefore the length effect on the calculation of hydraulic conductance is usually neglected and the gas flow in each pore and throat on the pore network is generally fully developed, as is also commonly adopted in previous studies (Mehmani et al., 2013; Mehmani and Prodanović, 2014; Tian and Daigle, 2018). According to Eq. (22), the free gas flow conductance in circle pore can be written as:

2 2 Z = (0.702e−2.5Tpr ) Ppr − 5.524e−2.5Tpr Ppr + 0.044Tpr − 0.164Tpr + 1.15

reff = r − dm θ

(23)

αduct = 1.7042

2 0.5 tan−1 [8Knsqu ] π

(28)

According to Eq. (26), free gas flow conductance in square pore can be written as:

(17)

gfree _ squ = C (AR)

(18)

6Knsqu ⎞ wh3 (1 + αduct Knsqu ) ⎜⎛1 + ⎟ 1 12μ βKnsqu ⎠ − ⎝

(29)

It should be noted that the value of slip coefficient in Eq. (23) and Eq. (29) is valid within the full range of flow regimes based on the DSMC and Boltzmann solutions (Karniadakis et al., 2006). The surface diffusion process of adsorbed gas molecules can be modelled as the general diffusion process. The molar flow rate per unit area of the

(19) (20) 3

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

Fig. 3. Two dimensional polished shale SEM image used in this study, unsegmented image on the left and segmented image on the right).

VAsqu =

4M dθ 2 2 dp (r − reff ) Ds Ca max dx ρ dp

(34)

According to methane adsorption experimental data and molecular simulation data (Song et al., 2018a), the surface diffusion coefficient at zero gas coverage can be expressed in Eq. (35) based on the Hwang and Kammermeyer's model (Guo et al., 2008; Hwang and Kammermeyer, 1966): 0.8

ΔH ⎞ Ds0 = 3.1321 × 10−7T 0.5 exp ⎛− ⎝ RT ⎠ ⎜

⎟

(35)

Where △H is the isosteric adsorption heat at the zero surface gas coverage (J/mol). Ds0 is obtained under the low pressure condition. In order to describe the gas surface diffusion in shale nanopores under the high pressure condition, the influence of gas coverage on surface diffusion is considered in Eq. (36) by kinetic method (Chen and Yang, 1991): κ

Ds = Ds0

Fig. 4. Reconstructed three dimensional shale binary model (pore phase is marked by blue color, matrix phase is marked by red color, physical size: 4.8 μm × 4.8 μm × 4.8 μm, voxel size: 400 × 400 × 400). (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

concentration gradient within the adsorbed monolayer is expressed as (Cunningham and Williams, 1980):

dCa dx

(30)

Where Ds is the surface diffusion coefficient (m /s), Ca is adsorbed gas concentration in the gas adsorption layer (mol/m3) and can be calculated by: (31)

Where Camax is the maximum absorbed gas concentration inside the gas adsorption layer (mol/m3) and can be expressed as (Wasaki and Akkutlu, 2014):

Ca max =

κb κm

(32)

Where VL is the Langmuir volume (m /kg), ρsc,gas is the gas density at standard condition (kg/m3), ρgrain is the rock density (kg/m3), εks is the total organic grain volume per total grain volume, dimensionless. Combining Eq. (30)-Eq. (32), adsorbed gas volumetric flux for the circle pore and square pore are given below:

M dθ 2 dp ) Ds Ca max π (r 2 − reff dx ρ dp

(37)

(38)

gads _ cir =

1M dθ 2 Ds Ca max π (r 2 − reff ) l ρ dp

(39)

gads _ squ =

4M dθ 2 2 Ds Ca max (r − reff ) l ρ dp

(40)

3

VAcir =

(36)

The total gas flow conductance can be derived in Eq. (41) and Eq. (42) based on the superposition of free gas flow conductance and adsorbed gas flow conductance:

VL ρsc, gas ρgrain εks M

κ 2

When κm is the rate constant for forward migration in surface diffusion, κb is the rate constant for blockage in surface diffusion. When κm is larger than κb, surface diffusion occurs. When κm is smaller than κb, gas molecules are blocked and surface diffusion stops. According to Eq. (33) and Eq. (34), adsorbed gas flow conductance in square pore and circle pore can be expressed respectively as:

2

Ca = Ca max θ

(1 − θ + θ )

2

H (1 − κ ) = 0, κ ≥ 1; 1,0 ≤ κ ≤ 1 κ=

Ja = Ds

κ

(1 − θ) + 2 θ (2 − θ) + {H (1 − κ )}(1 − κ ) 2 θ 2

gcir =

4 πreff f (Kncir )

8μl

gsqu = C (AR) +

(33) 4

+

1M dθ 2 Ds Ca max π (r 2 − reff ) l ρ dp

(41)

6Knsqu ⎞ wh3 (1 + αduct Knsqu ) ⎜⎛1 + ⎟ 1 12μ βKnsqu ⎠ − ⎝

4M dθ 2 2 (r − reff ) Ds Ca max l ρ dp

(42)

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

Fig. 8. Coordination number distribution. Fig. 5. Extracted shale pore network (throats are represented by gray color, pores are represented by yellow color, physical size: 4.8 μm × 4.8 μm × 4.8 μm). (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

Fig. 9. Calculation unit of gas flow conductance (single throat with two connected pores). Table 1 Parameters for modelling shale gas transport on pore network. Reservoir properties Temperature T (K) Pore pressure p (MPa) Pressure gradient dp/dx (MPa/m) Total organic grain volume per total grain volume εks Rock density ρgrain (kg/m3) Gangi's model parameter Surface roughness of pores m Effective stress when the shale nanopores are closed p1 (MPa) Effective stress coefficient αs Confining pressure Pcon (MPa) Gas properties Langmuir pressure pL (MPa)

Fig. 6. Shape factor distribution.

Langmuir volume VL (m3/Kg) Isosteric adsorption heat at zero gas coverage △H (J/mol) Ratio of the rate constant for blockage to the rate constant for forward migration κ Molecular weight M (kg/mol) Phase change parameters (Islam et al., 2015) vdW energy parameter a (m6·Pa/mol2) vdW energy parameter b (m3/mol) Lennard-Jones size parameter σ (m) Lennard-Jones energy parameter ε

Fig. 7. Pore and throat size distribution.

5

400 30 0.1 0.01 (Wasaki and Akkutlu, 2015) 2.66 × 103 (Wasaki and Akkutlu, 2015) 0.5 179.263682 (Wasaki and Akkutlu, 2015) 0.5 (Wasaki and Akkutlu, 2015) 40 13.789514 (Wasaki and Akkutlu, 2015) 0.0031 (Wasaki and Akkutlu, 2015) 16000 0.5 (Wu et al., 2016) 0.016 0.22998 4.28 × 10−5 3.73 × 10−9 2.0434 × 10−21

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

Fig. 10. Pressure drop distributions on the shale pore network during gas transport process (a) t = 5 × 10−9 s; (b) t = 10−7 s; (c) t = 2 × 10−7 s; (d) t = 4.8626 × 10−7 s.

Fig. 11. Decreasing gas permeability during dynamic gas transport process.

Fig. 12. Decreasing gas permeability and constant Darcy permeability on the shale pore network model versus pore pressure.

3. Gas transport on the pore network model phase by using median filter and then Otsu segmentation (Jassim and Altaani, 2013; Walls and Sinclair, 2011). Then a three dimensional shale binary model is reconstructed in Fig. 4 by applying the MultiplePoint Statistics (MPS) method (Okabe and Blunt, 2005) from the segmented shale SEM image in Fig. 3. The shale pore network model in Fig. 5 is extracted from 3D shale binary model according to the

3.1. Construction of a three dimensional shale model and pore network model A high resolution SEM image obtained from organic-rich Eagle Ford shale sample in Fig. 3 is first segmented into pore phase and matrix 6

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

1 1 1 1 = + + gi, j gi gt gj

(43)

The equivalent pore radius in consideration of shape factor for each individual element is expressed as (Ma et al., 2014; Oren et al., 1998; Valvatne and Blunt, 2004):

2Ap

r=

Pd

= 2 Ap G

(44) 2

Where Ap is the pore cross section area (m ), Pd is the pore perimeter (m). For single pore i at time step k to k+1, the gas mass change inside the pore i equals to the mass flux from pore i to all connected pores j: k k k k ⎛ pi + pj ⎞ ⎛ pi + pj ⎞ ⎛ k ρ p − pjk ⎞⎟ Δt ⎟⎜ i ⎟ ⎜ 2 2 ⎝ ⎠ ⎠ ⎠ ⎝ ⎝

n

mik − mik + 1 =

∑ g⎜ j=1

(45)

According to real gas equation of state, the real gas density can be given as:

Fig. 13. The slight increase of gas permeability on the shale pore network model versus formation temperature.

ρ=

maximal ball fitting method (Blunt et al., 2013; Dong and Blunt, 2009). The shape factor G is defined as the pore cross section area divided by the perimeter squared in each pore and throat spatial location on the shale pore network and its distribution is shown in Fig. 6. The pore shape are recognized as the circle pore and square pore according to the shape factor value. The fraction of pores and throats that are circular is about 17.5%. In this work, we follow the commonly used pore network model assumption that a circular cross section pore in 3D have cylindrical shape (Afsharpoor and Javadpour, 2016; Blunt, 2001; Yang et al., 2019). The calculated pore size distribution and throat size distribution are shown in Fig. 7. Coordination number is defined as the number of pore throats connected to a single pore and reflects the spatial connectivity of the pore space. The calculated coordination number distribution is shown in Fig. 8.

pM Z (p) RT

(46)

Therefore, gas mass in single pore at time step k and k+1 can be written as:

mik =

pik M Z (pik ) RT

mik + 1 =

Vi

(47)

pik + 1 M Z (pik + 1 ) RT

Vi

(48)

The volume of pore i can be expressed as:

Vi =

4 3 πri 3

(49)

Combining Eqs. (33)–(37) yields:

pik + 1

3.2. Pressure distributions on the shale pore network model

Z (pik + 1 )

For each unit comprising of a single throat with two connected pores in Fig. 9, the gas flow conductance gi,j is calculated by Eq. (43) according to the three individual element conductances (Valvatne, 2004).

=

pik Z (pik )

−

3 4πri3

n

k k k k ⎛ pi + pj ⎞ pi + pj ⎛ k k⎞ ⎜pi − pj ⎟ Δt ⎟ pik + pjk 2 ⎝ ⎠ ⎛ ⎞ ⎠ 2Z ⎜ ⎝ ⎟ 2 ⎝ ⎠

∑ g⎜ i=1

(50) A new definition to represent pressure in consideration of gas compressibility factor is regulated below:

Fig. 14. The impact of surface chemical property on gas permeability (a) maximum absorbed gas concentration (b) Ratio of the rate constant for blockage to the rate constant for forward migration. 7

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

Fig. 15. The impact of rock physical property on gas permeability (a) rock density (b) surface roughness. Ni

∑ Qij = 0 j=1

(54)

Qij = gij (pi − pj )

(55)

Because the gas viscosity changes with the pore size dependent critical properties (Pc,Tc), apparent gas permeability is calculated based on the product of gas flux and viscosity for every inlet pore voxel to compute the sum of the total flux at the inlet face in Eq. (56). Gas flux in each inlet pore is obtained by Eq. (55). The pressure value in Eq. (55) is obtained from the dynamic pore network calculation results at each time step or steady state pore network calculation results. N

kapp =

p Z

AΔp

(56)

The pressure drop on the pore network is 0.48Pa according to the model parameters in Table 1. Gas flows into the pore network from the inlet and gradually reaches the outlet during the pressure drop transmission process. The dynamic pressure distributions before gas reaching outlet and steady state pressure distribution when gas reaches outlet at 4.8626 × 10−7 s are shown in Fig. 10. Note that the simulation stops when gas reaches outlet and pressure distribution becomes stable.

Fig. 16. Decreasing trend of average pore radius versus effective stress at surface roughness m = 0.5.

m (p) =

L ∑i =inlet gi Δpi μinlet 1

(51)

Finally according to Eq. (51), Eq. (45) can be rearranged as:

m (pik + 1 ) = m (pik ) −

3 4πri3

n

4. Results and discussion

k k k k ⎛ pi + pj ⎞ ⎛ pi + pj ⎞ ⎛ k m p − pjk ⎞⎟ Δt ⎟⎜ i ⎟ ⎜ 2 2 ⎠ ⎠⎝ ⎠ ⎝ ⎝

∑ g⎜ i=1

Fig. 11 shows gas permeability notably decreases during the dynamic gas transport process. The initial permeability is 1.8 × 10−3 μm2, which is three order higher than the final steady state permeability 1.946 × 10−6 μm2. When the pressure drop initially imposes on the pore network, the pressure drop in each inlet pore △Pi equals to the total pressure drop △P on the pore network. Therefore, the initial gas permeability can be written in Eq. (57) and equals to 1.8 × 10−3 μm2. During the pressure drop transmission, △Pi significantly decreases which can be clearly seen in Fig. 10(a) and (b). The value of △Pi/△P is far less than one and gas permeability in Eq. (56) decreases with the increase of pressure drop transmission area. When the pressure drop reaches outlet, △Pi and gas permeability turns constant value. As is shown in Fig. 12, steady state gas permeability through pore network decreases with the increase of pore pressure. The simulated decreasing trend is in accordance with experimental data (Soeder, 1988). The Darcy permeability calculated based on the pore network model is compared with the apparent permeability obtained

(52) The range of time step is determined based on Courant number. Courant number reflects the relation between time step and space step and generally should be smaller than 1 in order to stabilize calculation:

vΔt = Δl

g ⎜⎛ ⎝

pik + pjk

⎞ ⎛p k − pk ⎞⎟ Δt i j ⎠⎝ ⎠ ≤1 2 πri (Li + Lt + Lj ) 2

⎟⎜

(53)

When the pressure drop reaches outlet, the pressure distribution on the pore network doesn't change with time and gas transport turns steady state. The steady state gas transport is defined by a system of the mass conservation equations (Eq. (54) on the whole network and momentum equations (Eq. (55)) for every network segment comprising a throat and two connected pores. 8

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

Fig. 17. The contribution of surface diffusion and free gas flow on total gas permeability at different pore pressure and effective stress (a) surface diffusion (b) free gas flow.

permeability. The increase of rock density increases the maximum absorbed gas concentration in Eq. (32) and the gas permeability slightly increases in Fig. 15(a). Gas permeability notably decreases with the increase of effective stress in Fig. 15(b). This can be attributed to the fact that the larger value of effective stress causes the further shrinkage of shale nanopore radius and the average pore radius decreases from 17.5 nm to 12.5 nm in Fig. 16. The small surface pore space caused by the surface roughness can behave as the effective gas transport path especially at high effective stress conditions with shrinked pore size. Therefore gas permeability decreases more slowly with the increase of effective stress at larger surface roughness value in Fig. 15(b). Fig. 17 shows the contribution of surface diffusion and free gas flow on total gas permeability at different pore pressure and effective stress. Fig. 17 (a) indicates that the contribution of surface diffusion on total gas permeability decreases from about 3.85% at 5 MPa to 0.007% at 50 MPa in low effective stress (σeff = 25 MPa) and free gas flow dominates the total gas permeability in Fig. 17(b). The decrease of pore size at high effective stress (σeff = 65 MPa) reduces the relative area for free gas flow. Therefore the contribution of surface diffusion on total gas permeability notably increases with the increase of effective stress and can reach 7.5% at 5 MPa. As free gas flow dominates the total gas permeability, we further analyze the influence of slip coefficient on gas permeability in Fig. 18. The gas permeability slightly increases with the increase of the slip coefficient when the pore pressure is less than 20 MPa and almost stays constant at 30 MPa because of the small Knudsen number in Eq. (14) and Eq. (15) at high pore pressure.

Fig. 18. The increasing trend of gas permeability versus slip coefficient at different pore pressure.

from our model. It is can be seen in Fig. 12 that the gap between kDarcy and kapp is larger in low pressure than that in high pressure. This is because gas flow conductance in Eq. (41) and Eq. (42) is larger in low pore pressure than that in high pore pressure. The gas mean free path in Eq. (16) increases with the increase of formation temperature. Therefore Knudsen number and gas permeability in Fig. 13 slightly increases as a consequence.

5. Conclusion

N

kinitial =

L ∑i =inlet gi μinlet 1 A

The dynamic single component methane transport process on the realistic shale pore structure was studied in detail by a dynamic pore network model. Gas transport mechanisms considered slip flow, transition flow, surface diffusion and ad/desorption. Gas critical property change induced by confined pore space and real gas effect were considered when calculating gas properties. The analysis results reveal that gas permeability is influenced by the dynamic pressure transmission process. Gas permeability notably decreases during dynamic gas transport process before pressure drop reaching outlet. When the pressure drop reaches outlet, gas permeability is influenced by the effective stress, pressure, temperature and shale rock property. Our model can be applied to better understand dynamic shale gas transport process during field production and laboratory core experiment.

(57)

Fig. 14 shows gas permeability change at different maximum absorbed gas concentration and ratio of the rate constant for blockage to the rate constant for forward migration. The increase of maximum absorbed gas concentration enhances the surface diffusion capacity by increasing the surface diffusion conductance in Eq. (41), Eq. (42). The gas permeability consequently increases with the increase of maximum absorbed gas concentration in Fig. 14(a). The surface diffusion coefficient decreases with the increase of the ratio of the rate constant for blockage to the rate constant for forward migration according to Eq. (36)-Eq. (38) and causes the slight decrease of gas permeability in Fig. 14(b). Fig. 15 shows the impact of the rock physical property on gas 9

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

Javadpour, F., Fisher, D., Unsworth, M., 2007. Nanoscale gas flow in shale gas sediments. J. Can. Pet. Technol. 46. Jiang, Z., Wu, K., Couples, G., Van Dijke, M.I.J., Sorbie, K.S., Ma, J., 2007. Efficient extraction of networks from three‐dimensional porous media. Water Resour. Res. 43. Karniadakis, G.E., Beskok, A., Aluru, N., 2006. Microflows and Nanoflows: Fundamentals and Simulation. Springer Science & Business Media. Kim, C., Jang, H., Lee, Y., Lee, J., 2016. Diffusion characteristics of nanoscale gas flow in shale matrix from Haenam basin, Korea. Environ. Earth Sci. 75, 350. Kim, T.H., Cho, J., Lee, K.S., 2017. Evaluation of CO2 injection in shale gas reservoirs with multi-component transport and geomechanical effects. Appl. Energy 190, 1195–1206. Knackstedt, M., Arns, C., Ghous, A., Sakellariou, A., Senden, T., Sheppard, A., Sok, R., Averdunk, H., Pinczewski, W.V., Padhy, G.S., 2006. 3D Imaging and Flow Characterization of the Pore Space of Carbonate Core Samples. Kou, J., Wu, F., Lu, H., Xu, Y., Song, F., 2009. The effective thermal conductivity of porous media based on statistical self-similarity. Phys. Lett. A 374, 62–65. Kwon, O., Kronenberg, A.K., Gangi, A.F., Johnson, B., Herbert, B.E., 2004. Permeability of illite‐bearing shale: 1. Anisotropy and effects of clay content and loading. J. Geophys. Res. Solid Earth 109. Landry, C.J., Prodanović, M., Eichhubl, P., 2016. Direct simulation of supercritical gas flow in complex nanoporous media and prediction of apparent permeability. Int. J. Coal Geol. 159, 120–134. Lee, A., Gonzalez, M., Eakin, B., 1966. The viscosity of natural gases. J. Pet. Technol. 18, 997–1000. Li, Z., Yao, J., Ren, Z., Sun, H., Zhang, L., Yang, Y., Fan, J., Kou, J., 2019. Accumulation behaviors of methane in the aqueous environment with organic matters. Fuel 236, 836–842. Loucks, R.G., Reed, R.M., Ruppel, S.C., Jarvie, D.M., 2009. Morphology, genesis, and distribution of nanometer-scale pores in siliceous mudstones of the Mississippian Barnett Shale. J. Sediment. Res. 79, 848–861. Ma, J., Sanchez, J.P., Wu, K., Couples, G.D., Jiang, Z., 2014. A pore network model for simulating non-ideal gas flow in micro-and nano-porous materials. Fuel 116, 498–508. Mahmoud, M., 2014. Development of a new correlation of gas compressibility factor (Zfactor) for high pressure gas reservoirs. J. Energy Resour. Technol. 136, 12903. Mehmani, A., Prodanović, M., 2014. The application of sorption hysteresis in nano-petrophysics using multiscale multiphysics network models. Int. J. Coal Geol. 128, 96–108. Mehmani, A., Prodanović, M., Javadpour, F., 2013. Multiscale, multiphysics network modeling of shale matrix gas flows. Transp. Porous Media 99, 377–390. Middleton, R.S., Gupta, R., Hyman, J.D., Viswanathan, H.S., 2017. The shale gas revolution: barriers, sustainability, and emerging opportunities. Appl. Energy 199, 88–95. Okabe, H., Blunt, M.J., 2005. Pore space reconstruction using multiple-point statistics. J. Pet. Sci. Eng. 46, 121–137. Oren, P.-E., Bakke, S., Arntzen, O.J., 1998. Extending predictive capabilities to network models. SPE J. 3, 324–336. Ostadi, H., Rama, P., Liu, Y., Chen, R., Zhang, X.X., Jiang, K., 2010. 3D reconstruction of a gas diffusion layer and a microporous layer. J. Membr. Sci. 351, 69–74. Saif, T., Lin, Q., Butcher, A.R., Bijeljic, B., Blunt, M.J., 2017. Multi-scale multi-dimensional microstructure imaging of oil shale pyrolysis using X-ray micro-tomography, automated ultra-high resolution SEM, MAPS Mineralogy and FIB-SEM. Appl. Energy 202, 628–647. Sheng, G., Javadpour, F., Su, Y., 2018. Effect of microscale compressibility on apparent porosity and permeability in shale gas reservoirs. Int. J. Heat Mass Transf. 120, 56–65. Singh, H., 2013. Nonempirical Apparent Permeability of Shale. Unconventional Resources Technology Conference (URTEC). Singh, H., Javadpour, F., 2016. Langmuir slip-Langmuir sorption permeability model of shale. Fuel 164, 28–37. Soeder, D.J., 1988. Porosity and permeability of eastern Devonian gas shale. SPE Form. Eval. 3, 116–124. Song, W., Wang, D., Yao, J., Li, Y., Sun, H., Yang, Y., Zhang, L., 2019. Multiscale imagebased fractal characteristic of shale pore structure with implication to accurate prediction of gas permeability. Fuel 241, 522–532. Song, W., Yao, B., Yao, J., Li, Y., Sun, H., Yang, Y., Zhang, L., 2018a. Methane surface diffusion capacity in carbon-based capillary with application to organic-rich shale gas reservoir. Chem. Eng. J. 352, 644–654. Song, W., Yao, J., Li, Y., Sun, H., Yang, Y., 2018b. Fractal models for gas slippage factor in porous media considering second-order slip and surface adsorption. Int. J. Heat Mass Transf. 118, 948–960. Song, W., Yao, J., Li, Y., Sun, H., Zhang, L., Yang, Y., Zhao, J., Sui, H., 2016. Apparent gas permeability in an organic-rich shale reservoir. Fuel 181, 973–984. Song, W., Yao, J., Ma, J., Couples, G., Li, Y., 2017. Assessing relative contributions of transport mechanisms and real gas properties to gas flow in nanoscale organic pores in shales by pore network modelling. Int. J. Heat Mass Transf. 113, 524–537. Song, W., Yao, J., Ma, J., Li, A., Li, Y., Sun, H., Zhang, L., 2018c. Grand canonical Monte Carlo simulations of pore structure influence on methane adsorption in micro-porous carbons with applications to coal and shale systems. Fuel 215, 196–203. Sun, F., Yao, Y., Li, G., Zhang, S., Xu, Z., Shi, Y., Li, X., 2019. A slip-flow model for oil transport in organic nanopores. J. Pet. Sci. Eng. 172, 139–148. Sun, H., Yao, J., Fan, D., Wang, C., Sun, Z., 2015. Gas transport mode criteria in ultratight porous media. Int. J. Heat Mass Transf. 83, 192–199. Tiab, D., Donaldson, E.C., 2015. Petrophysics: Theory and Practice of Measuring Reservoir Rock and Fluid Transport Properties. Gulf professional publishing. Tian, X., Daigle, H., 2018. Permeability prediction from a pore-scale network model

Declarations of interest None. Acknowledgments This project was supported by the Major Projects of the National Science and Technology (2016ZX05061), the Fundamental Research Funds for the Central Universities (18CX06007A, 18CX06008A, 17CX05003), National Natural Science Foundation of China (No. 51504276, No. 51674280), Shandong Provincial Natural Science Foundation, China (ZR2019JQ21), Key Research and Development Plan of Shandong Province (2018GSF116009). Appendix A. Supplementary data Supplementary data to this article can be found online at https:// doi.org/10.1016/j.petrol.2019.106506. References Afsharpoor, A., Javadpour, F., 2016. Liquid slip flow in a network of shale noncircular nanopores. Fuel 180, 580–590. Bai, B., Elgmati, M., Zhang, H., Wei, M., 2013. Rock characterization of Fayetteville shale gas plays. Fuel 105, 645–652. Beskok, A., Karniadakis, G.E., 1999. Report: a model for flows in channels, pipes, and ducts at micro and nano scales. Microscale Thermophys. Eng. 3, 43–77. Blunt, M.J., 2001. Flow in porous media—pore-network models and multiphase flow. Curr. Opin. Colloid Interface Sci. 6, 197–207. Blunt, M.J., Bijeljic, B., Dong, H., Gharbi, O., Iglauer, S., Mostaghimi, P., Paluszny, A., Pentland, C., 2013. Pore-scale imaging and modelling. Adv. Water Resour. 51, 197–216. Bui, B.T., Liu, H.-H., Chen, J., Tutuncu, A.N., 2016. Effect of capillary condensation on gas transport in shale: a pore-scale model study. SPE J. 21, 601–612. Cao, P., Liu, J., Leong, Y.-K., 2016. General gas permeability model for porous media: bridging the gaps between conventional and unconventional natural gas reservoirs. Energy Fuels 30, 5492–5505. Chen, Y.D., Yang, R.T., 1991. Concentration dependence of surface diffusion and zeolitic diffusion. AIChE J. 37, 1579–1582. Civan, F., 2010. Effective correlation of apparent gas permeability in tight porous media. Transp. Porous Media 82, 375–384. Civan, F., Devegowda, D., Sigal, R.F., 2013. Critical evaluation and improvement of methods for determination of matrix permeability of shale. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. Civan, F., Rai, C.S., Sondergeld, C.H., 2011. Shale-gas permeability and diffusivity inferred by improved formulation of relevant retention and transport mechanisms. Transp. Porous Media 86, 925–944. Clarkson, C.R., Solano, N., Bustin, R.M., Bustin, A.M.M., Chalmers, G.R.L., He, L., Melnichenko, Y.B., Radliński, A.P., Blach, T.P., 2013. Pore structure characterization of North American shale gas reservoirs using USANS/SANS, gas adsorption, and mercury intrusion. Fuel 103, 606–616. Cunningham, R.E., Williams, R.J.J., 1980. Diffusion in Gases and Porous Media. Springer. Curtis, M.E., Ambrose, R.J., Sondergeld, C.H., 2010. Structural characterization of gas shales on the micro-and nano-scales. In: Canadian Unconventional Resources and International Petroleum Conference. Society of Petroleum Engineers. Darabi, H., Ettehad, A., Javadpour, F., Sepehrnoori, K., 2012. Gas flow in ultra-tight shale strata. J. Fluid Mech. 710, 641–658. Dong, H., Blunt, M.J., 2009. Pore-network extraction from micro-computerized-tomography images. Phys. Rev. E 80, 36307. Gangi, A.F., 1978. Variation of whole and fractured porous rock permeability with confining pressure. In: International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts. Elsevier, pp. 249–257. Geng, L., Li, G., Tian, S., Sheng, M., Ren, W., Zitha, P., 2016. A fractal model for real gas transport in porous shale. AIChE J. 63, 1430–1440. Guo, L., Peng, X., Wu, Z., 2008. Dynamical characteristics of methane adsorption on monolith nanometer activated carbon. Huagong Xuebao/J. Chem. Ind. Eng. 59, 2726–2732. Ho, C.-M., Tai, Y.-C., 1998. Micro-electro-mechanical-systems (MEMS) and fluid flows. Annu. Rev. Fluid Mech. 30, 579–612. Huang, X., Bandilla, K.W., Celia, M.A., 2016. Multi-physics pore-network modeling of two-phase shale matrix flows. Transp. Porous Media 111, 123–141. Hwang, S., Kammermeyer, K., 1966. Surface diffusion in microporous media. Can. J. Chem. Eng. 44, 82–89. Islam, A.W., Patzek, T.W., Sun, A.Y., 2015. Thermodynamics phase changes of nanopore fluids. J. Nat. Gas Sci. Eng. 25, 134–139. Jassim, F.A., Altaani, F.H., 2013. Hybridization of OTSU Method and Median Filter for Color Image Segmentation. arXiv Prepr. arXiv1305.1052. Javadpour, F., 2009. Nanopores and apparent permeability of gas flow in mudrocks (shales and siltstone). J. Can. Pet. Technol. 48, 16–21.

10

Journal of Petroleum Science and Engineering 184 (2020) 106506

W. Song, et al.

Wu, K., Chen, Z., Li, X., Guo, C., Wei, M., 2016. A model for multiple transport mechanisms through nanopores of shale gas reservoirs with real gas effect–adsorptionmechanic coupling. Int. J. Heat Mass Transf. 93, 408–426. Xiong, X., Devegowda, D., Villazon, M., German, G., Sigal, R.F., Civan, F., 2012. A fullycoupled free and adsorptive phase transport model for shale gas reservoirs including non-Darcy flow effects. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. Yang, R., Hao, F., He, S., He, C., Guo, X., Yi, J., Hu, H., Zhang, S., Hu, Q., 2017. Experimental investigations on the geometry and connectivity of pore space in organic-rich Wufeng and Longmaxi shales. Mar. Pet. Geol. 84, 225–242. Yang, Y., Wang, K., Zhang, L., Sun, H., Zhang, K., Ma, J., 2019. Pore-scale simulation of shale oil flow based on pore network model. Fuel 251, 683–692. Yao, J., Sun, H., Huang, Z., Zhang, L., Zeng, Q., Sui, H., Fan, D., 2013. Key mechanical problems in the development of shale gas reservoirs. Sci. Sin. Phys. Mech. Astron. 43, 1527. Yuan, J., Luo, D., Feng, L., 2015. A review of the technical and economic evaluation techniques for shale gas development. Appl. Energy 148, 49–65. Zhang, P., Hu, L., Meegoda, J.N., Gao, S., 2015. Micro/nano-pore network analysis of gas flow in shale matrix. Sci. Rep. 5. Zheng, Q., Yu, B., Duan, Y., Fang, Q., 2013. A fractal model for gas slippage factor in porous media in the slip flow regime. Chem. Eng. Sci. 87, 209–215. Ziarani, A.S., Aguilera, R., 2012. Knudsen's permeability correction for tight porous media. Transp. Porous Media 91, 239–260. Zou, Y., Yang, C., Wu, D., Yan, C., Zeng, M., Lan, Y., Dai, Z., 2016. Probabilistic assessment of shale gas production and water demand at Xiuwu Basin in China. Appl. Energy 180, 185–195.

constrained by low-pressure nitrogen sorption isotherms. J. Pet. Sci. Eng. 162, 554–566. Valvatne, P.H., 2004. Predictive Pore-Scale Modelling of Multiphase Flow. Imperial College London. Valvatne, P.H., Blunt, M.J., 2004. Predictive pore‐scale modeling of two‐phase flow in mixed wet media. Water Resour. Res. 40. Walls, J.D., Sinclair, S.W., 2011. Eagle Ford shale reservoir properties from digital rock physics. First Break 29, 97–101. Wang, J., Li, Y., Yang, Z., Gong, H., Dong, M., 2016a. Measurement of dynamic adsorption–diffusion process of methane in shale. Fuel 172, 37–48. Wang, J., Luo, H., Liu, H., Cao, F., Li, Z., Sepehrnoori, K., 2017a. An integrative model to simulate gas transport and production coupled with gas adsorption, non-Darcy flow, surface diffusion, and stress dependence in organic-shale reservoirs. SPE J. 22, 244–264. Wang, J., Yang, Z., Dong, M., Gong, H., Sang, Q., Li, Y., 2016b. Experimental and numerical investigation of dynamic gas adsorption/desorption–diffusion process in shale. Energy Fuels 30, 10080–10091. Wang, S., Lukyanov, A.A., Wang, L., Wu, Y.-S., Pomerantz, A., Xu, W., Kleinberg, R., 2017b. A non-empirical gas slippage model for low to moderate Knudsen numbers. Phys. Fluids 29, 12004. Wang, S., Shi, J., Wang, K., Sun, Z., Miao, Y., Hou, C., 2018. Apparent permeability model for gas transport in shale reservoirs with nano-scale porous media. J. Nat. Gas Sci. Eng. 55, 508–519. Wasaki, A., Akkutlu, I.Y., 2015. Permeability of organic-rich shale. SPE J. 20, 1–384. Wasaki, A., Akkutlu, I.Y., 2014. Permeability of organic-rich shale. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. White, F.M., Corfield, I., 2006. Viscous Fluid Flow. McGraw-Hill, New York.

11