TY - JOUR AU - Guo,, Xiao AB - Abstract Flow in extremely low permeability shale gas reservoirs undergoes a transition from a Darcy regime to other flow regimes including slip flow, transition flow and free molecular flow, due to the significant effect of molecular collisions with pore walls on gas transport. These various flow regimes in shales and their effect on actual gas production are not clearly understood, and multi-fractured horizontal wells are crucial for the economic production of unconventional resources. Therefore, a numerical model which is valid for the entire Knudsen range (continuum flow, slip flow, transition flow and free molecular flow) in shale gas reservoirs has been developed, with the effect of gas–water flow and the simulation of hydraulic fracturing cracks taken into consideration as well. Furthermore, the effect of different flow regimes on both the dynamic permeability of shales and then the actual gas production rate are analyzed thoroughly, and we study the influencing factors such as gas–water flow and fracturing treatment. The simulation results of the developed model are validated with other numerical models and field data. It is concluded that different flow regimes, considering slippage and diffusion in nanopores, have a greater impact on gas production when shale permeability is lower than 1 × 10−4 mD. The presented numerical model can provide a theoretical basis for evaluating and predicting the productivity of shale gas wells. shale gas, multi-scale, gas–water flow, dynamic permeability, multi-stage fractured horizontal well Nomenclature Latin b slippage coefficient Bg gas volume factor Bw water volume factor D depth of formation, m DK Knudsen diffusion coefficient, m2 s–1 g gravitational acceleration, m s–2 h thickness of grid block, m i, j, k coordinates of grid block K apparent permeability k permeability measured in laboratory KeH equivalent permeability of hydraulic fracture KH permeability of hydraulic fractures considering stress-dependence effect KH0 initial permeability of hydraulic fractures Kn Knudsen number krg relative permeability of gas phase krw relative permeability of water phase K0 absolute permeability M molecular mass, kg mol–1 p pressure, Pa pg pressure of gas phase, Pa pgi initial pressure of gas phase, Pa pL Langmuir’s pressure, Pa pw pressure of water phase, Pa pwf pressure of wellbore, Pa qg gas volume flux per unit volume of shale and per unit time qvg produced gas volume flux per unit volume and per unit time in well-block qvw produced water volume flux per unit volume and per unit time in well-block qw water volume flux per unit volume of shale and per unit time R universal gas constant, 8.314 J K−1 mol−1 r pore radius, m re equivalent radius of well-block, m rw wellbore radius, m S skin factor sg gas saturation sw water saturation T temperature at formation condition, K v gas flow rate, m s–1 Va volume of adsorbed gas (standard condition) under formation pressure, m3 kg–1 VL Langmuir’s volume at standard condition, m3 kg–1 wH width of hydraulic fracture, m x, y, z distance coordinates, m ΔxH length of the refined grid along the well, m Greek δ(i,j,k) Dirac delta function α rarefication coefficient αs stress-sensitivity coefficient for hydraulic fracture, Pa−1 ξ permeability correction factor λ¯ mean free path of gas molecule, m μ viscosity, Pa · s μg gas viscosity, Pa · s μw water viscosity, Pa · s ρbi bulk density of shale at initial reservoir pressure, kg m–3 ρg gas density, kg m–3 ρw water density, kg m–3 τ external boundary of reservoir ϕ shale porosity 1. Introduction Unconventional gas resources play an important role in the world energy supply due to the decline in conventional gas production. Core analysis shows that the permeability of shale is in the range 1 × 10−9–1 × 10−3μm2, and gas flow in extremely low permeability shales undergoes a transition from the Darcy regime to other regimes owing to the significant effect of collisions between molecules and pore walls on gas transport (Deng et al2014). The Darcy formula cannot describe all the flow regimes in shale gas reservoirs (Javadpour et al2007), so the dynamics of gas flow in shales with wide pore-size distribution has become an important research topic during the current decade in the oil and gas industry. Gas transfer in nanopores is a complex process under the effects of pore diffusion, slip flow and gas desorption (Brown et al1946, Roy et al2003, Javadpour 2009, Li et al2014, Wu et al2016a and 2016b). Numerical models, which consider Knudsen diffusion and slip flow in nanopores for shale gas reservoirs, have been analyzed by many researchers (Ozkan et al2010, Shabro et al2011, Swami and Settari 2012, Huang et al2015a, 2015b). However, pore structure in shale is complicated because of its wide pore-size distribution, and Knudsen diffusion and the slippage effect in nanopores cannot describe all the flow regimes in multi-scale pores. So, a Knudsen number is used to classify the flow regimes in pores (Roy et al2003, Ziarani and Aguilera 2012), and Beskok and Karniadakis (1999) developed a rigorous equation, which is applicable for the entire Knudsen range, to account for all flow regimes including continuum flow, slip flow, transition flow and molecular flow. Ziarani and Aguilera (2012) used the Beskok–Karniadakis equation to correct the permeability of shales and validated it with data from the Mesaverde formation. Wang et al (2015) presented the apparent permeability for gas transport in the nanopores of shale reservoirs based on the Beskok–Karniadakis equation. Yuan et al (2016) used the Beskok–Karniadakis equation to build an analytical model of apparent gas permeability for tight porous media. Besides this, the economic production of unconventional resources also relies heavily on advanced completion technology such as horizontal wells with multi-stage hydraulic fracture simulations. Regarding the research on the rate transient analysis of multi-fractured horizontal wells, an analytical trilinear flow model was developed to study the transient behavior of fractured wells (Ozkan et al2010, Brown et al2011). A semi-analytical model based on Green’s function and the source/sink method has also been presented to facilitate the transient pressure/rate for a fractured well (Guo et al2012, Zhao et al2013, Wang 2014). However, many simplified assumptions are needed to solve the analytical or semi-analytical solution for shale gas reservoirs with complex storage and seepage mechanisms, which are not fully built on the real physical model compared with the numerical models. Furthermore, the performance of multi-fractured horizontal wells has been modeled by many researchers (Wu et al2009, Cipolla et al2010, Clarkson et al2012, Li et al2015 and 2016b), employing the method of numerical simulation or using commercial software. However, little work has been done in previous literature (Wu et al2009, Cipolla et al2010, Ozkan et al2010, Brown et al2011, Guo et al2012, Clarkson et al2012, Zhao et al2013, Wang 2014, Li et al2015, Wang et al2015, Yuan et al2016) to simultaneously incorporate the complex storage and water–gas flow mechanisms in shales and simulate multi-stage fractured horizontal wells using numerical models in shale gas reservoirs. In this paper, we provide a way of thinking or a method that researchers can use to combine different flow regimes in microscale pores, gas–water two-phase flow and a multi-stage fractured horizontal well in a macroscale shale gas reservoir together, and the simulation results are validated with other numerical models and field data. The presented numerical models can provide a theoretical basis for evaluating and predicting the productivity of shale gas wells. 2. Multi-scale seepage nonlinear model in shale gas reservoirs 2.1. Knudsen number and multi-scale flow regimes in shale gas reservoirs Shales consist of large amounts of nanoscale pores, a certain number of microscale pores and micro-fractures (figure 1). In microscale pores or microfractures, the Darcy flow is the dominant flow transport mechanism; however, when the pore radius is as small as a few nanometers, gas flow in nanoscale pores must consider the effect of collision between gas molecules and solid walls as well as the slippage phenomenon (Roy et al2003, Javadpour et al2007). So, a rigorous approach is needed to characterize the permeability of shales when nanoscale pores, microscale pores and micro-fractures co-exist. Figure 1. View largeDownload slide Microscopic pore structure of shales. Reproduced from Huang et al (2017). CC BY 4.0. Figure 1. View largeDownload slide Microscopic pore structure of shales. Reproduced from Huang et al (2017). CC BY 4.0. The Knudsen number Kn is defined as the ratio of the molecular mean free path λ and pore radius r, which is a widely recognized dimensionless parameter to determine the degree of appropriateness of the continuum model Kn=λ¯r. 1 Gas flow regimes can be classified into four categories (Roy et al2003) based on the Knudsen number (table 1): (1) continuum flow; (2) slip flow; (3) transition flow; (4) free-molecule flow. In the continuum flow regime, the no-slip boundary condition is valid and gas flow is linear. As the Knudsen number increases, the rarefaction effects become more pronounced and the continuum assumption breaks down eventually. So for flow regimes other than the continuum flow, the traditional equation of Darcy’s law is not applicable anymore. Table 1. Classification of gas flow regimes based on the Knudsen number. Knudsen number Kn ≤ 0.001 0.001 < Kn ≤ 0.1 0.1 < Kn ≤ 10 Kn > 10 Flow regime Continuum flow Slip flow Transition flow Free-molecule flow Knudsen number Kn ≤ 0.001 0.001 < Kn ≤ 0.1 0.1 < Kn ≤ 10 Kn > 10 Flow regime Continuum flow Slip flow Transition flow Free-molecule flow View Large Table 1. Classification of gas flow regimes based on the Knudsen number. Knudsen number Kn ≤ 0.001 0.001 < Kn ≤ 0.1 0.1 < Kn ≤ 10 Kn > 10 Flow regime Continuum flow Slip flow Transition flow Free-molecule flow Knudsen number Kn ≤ 0.001 0.001 < Kn ≤ 0.1 0.1 < Kn ≤ 10 Kn > 10 Flow regime Continuum flow Slip flow Transition flow Free-molecule flow View Large To analyze the flow regimes in shale gas reservoirs, figure 2 presents the Knudsen number for different pore sizes between 1 nm–50 μm and different pressures ranging from atmospheric to 100 MPa. The Knudsen number increases when the pressure drops and the pore throat diameter decreases. As illustrated in figure 2, we find that the gas flow in microscale pores (pore radius: 12–800 μm) and the micro-fracture is mainly a continuum flow, which can be presented by the Darcy formula. However, when the pore radius decreases from the microscale to the nanoscale, gas flow undergoes a transition from the Darcy regime to slip flow and transition flow. Only when the pressure is below 1 MPa and the pore radius is about 1 nm will the gas flow in a shale gas reservoir have free-molecule flow. Thus, we can conclude that gas flow in shale gas reservoirs is a multi-scale flow process mainly including continuum flow, slip flow and transition flow. Figure 2. View largeDownload slide Knudsen number at different pressures and pore radii. Reproduced from Huang et al (2017). CC BY 4.0. Figure 2. View largeDownload slide Knudsen number at different pressures and pore radii. Reproduced from Huang et al (2017). CC BY 4.0. 2.2. Dynamic permeability model valid for different flow regimes Beskok–Karniadakis (1999) developed a rigorous equation, valid for the entire Knudsen range, which accounts for different flow regimes including the continuum flow, slip flow, transition flow and free-molecule flow. The relationship between the flow rate and pressure gradient can be described as (Civan 2010, Ziarani and Aguilera 2012): v=−K0μ(1+αKn)1+4Kn1−bKndpdx. 2 So the apparent permeability and permeability correction factor can be written as: K=K0ξ 3 ξ=(1+αKn)1+4Kn1−bKn. 4 We can see from equations (3) and (4) that the apparent permeability is nearly equal to the value of absolute permeability when Kn approaches 0. When the Knudsen number gets bigger, the flow in the microtube is no longer Darcy flow, which means the apparent permeability should be corrected. The rarefication coefficient is defined by Beskok and Karniadakis (1999): α(Kn)=12815π2tan−1(4Kn0.4). 5 The gas phase molecular mean free path is defined as (Guggenheim 1960): λ¯=πRT2Mμp. 6 The Knudsen diffusion coefficient mainly changes with pore diameter, and a detailed derivation process is given by Javadpour et al (2007), which is written as: DK=2r38RTπM0.5. 7 The Knudsen number can be represented by the Knudsen diffusion coefficient according to equations (6) and (7): Kn=3π8r2μpDK. 8 So equation (3) can be rewritten as: K=K01+α3πμ8r2pDK1+12πμDK8r2p−3πbμDK. 9 According to the Hagen–Poiseuille model, the relationship between permeability and the pore radius in porous media can be written as: K0=ϕr28. 10 Putting equation (10) into equation (9), we can get the dynamic apparent permeability of shales during the gas production period: K=K0ξ=K01+α3πμϕ64K0pDK×1+12πμϕDK64K0p−3πbμϕDK 11 where the dynamic permeability is relevant to the absolute permeability, porosity, pressure, slippage factor and Knudsen diffusion coefficient. The pore radius can also be calculated through the correlations listed in table 2 (Ziarani and Aguilera 2012): Table 2. Various correlations for calculating pore radius. Model Correlation Assumption Heid et al (1950) r=8.886×10−2kϕ0.5 The core is homogeneous Kolodzie (1980) r=5.395k0.588(100ϕ)0.864 Obtained through mercury injection capillary pressure test Aguilera (2002) r=2.665k100ϕ0.45 At 35% mercury saturation Model Correlation Assumption Heid et al (1950) r=8.886×10−2kϕ0.5 The core is homogeneous Kolodzie (1980) r=5.395k0.588(100ϕ)0.864 Obtained through mercury injection capillary pressure test Aguilera (2002) r=2.665k100ϕ0.45 At 35% mercury saturation Units: r (μm), k (md) and ϕ (fraction). View Large Table 2. Various correlations for calculating pore radius. Model Correlation Assumption Heid et al (1950) r=8.886×10−2kϕ0.5 The core is homogeneous Kolodzie (1980) r=5.395k0.588(100ϕ)0.864 Obtained through mercury injection capillary pressure test Aguilera (2002) r=2.665k100ϕ0.45 At 35% mercury saturation Model Correlation Assumption Heid et al (1950) r=8.886×10−2kϕ0.5 The core is homogeneous Kolodzie (1980) r=5.395k0.588(100ϕ)0.864 Obtained through mercury injection capillary pressure test Aguilera (2002) r=2.665k100ϕ0.45 At 35% mercury saturation Units: r (μm), k (md) and ϕ (fraction). View Large 3. Simulation of shale gas reservoir considering gas–water flow mechanism 3.1. Mathematical model of gas and water transport in multi-scale shale gas reservoirs In order to deduce and develop the seepage equations for shale gas reservoirs, some assumptions have been made, which include the isothermal reservoir, non-Darcy flow in shales and Darcy flow in hydraulic fractures. The adsorbed gas only desorbs from the pore walls within the shales. Based on the above assumptions, the material balance equations for shale gas reservoirs are derived as follows, according to the law of conservation of matter. 3.1.1. Gas flow in multi-scale shales The mass balance equation for gas flow in multi-scale shales per unit time per bulk volume of shale can be derived, considering the adsorbed gas on pore walls and the multi-scale flow regimes in shales: ∇⋅K⃑krgBgμg(∇pg−ρgg∇D)−qg=∂∂tϕsgBg+ρbiVLpgpL+pg 12 where the adsorbed gas can be calculated by using the Langmuir isotherm equation: Va=VLpgpL+pg. 13 Putting equation (11) into equation (12), we can get the governing equation of the gas phase: ∇⋅K0⃑ξkrgBgμg(∇pg−ρgg∇D)−qg=∂∂tϕsgBg+ρbiVLpgpL+pg 14 where ξ in equation (14) can be calculated by equation (11), which is used to correct the apparent permeability of shales caused by the Knudsen diffusion and slippage effect in nanoscale and microscale pores: ξ=1+α3πμϕ64K0pDK1+12πμϕDK64K0p−3πbμϕDK. 3.1.2. Water flow in multi-scale shales The mass balance equation for water flow in multi-scale shales per unit time per bulk volume of shale is: ∇⋅K⃑krwBwμw(∇pw−ρwg∇D)−qw=∂∂tϕswBw. 15 We can get the governing equation of the water phase: ∇⋅K0⃑krwBwμw(∇pw−ρwg∇D)−qw=∂∂tϕswBw. 16 3.1.3. Gas and water flow in hydraulic fracture Experimental results indicate that the effective conductivity of both propped and unpropped fractures decreases with a decrease in effective stress (Wang et al2016). According to the experimental results of Zhang et al (2013), the relationship between the effective permeability and gas pressure can be written as: KH=KH0⋅e−αs(pgi−pg). 17 By applying the theory of equivalent fracture conductivity (Wang et al2016), the permeability of the hydraulic fracture considering the stress-dependence effect can be calculated by: KeH=KH⋅wHΔxH=KHi⋅wHΔxH⋅e−αs(pgi−pg). 18 Gas flow in the hydraulic fracture can be described by the Darcy equation, so the following equation is achieved: −∇pl=μKeHkrlvl(l=waterorgas). 19 Outer boundary condition: We often consider the closed boundary for reservoirs without edges or bottom water, so the outer boundary condition can be written as: ∂pg∂nτ=0. 20 Inner boundary condition: If the well is producing at a constant rate, we can employ the Dirac function to deal with the well point if the well radius is very small compared with the well spacing. The production rate can be represented by: qg(i,j,k,t)=qg(t)δ(i,j,k) 21 where δ(i, j, k) = 1 when a production well exists in grid block (i, j, k). If the well is producing at a constant wellbore pressure, the gas flow in this well-block is assumed to be at a pseudo-steady state and the production rate can be determined by the Peaceman model (Peaceman 1983): qvg=2πhK0krgBgμg(lnre/rw+S)(pgi,j,k−pwf)δ(i,j,k) 22 qvw=2πhK0krwBwμw(lnre/rw+S)(pwi,j,k−pwf)δ(i,j,k). 23 For anisotropic formation, the permeability of shale in equations (22) and (23) can be calculated by equation (24), and the equivalent radius of the well-block can be described by equation (25) (the direction of permeability refers to table 3): K0=klkm 24 re=0.28[(kl/km)1/2Δm2+(km/kl)1/2Δl2]1/2(kl/km)1/4+(km/kl)1/4. 25 Table 3. Coordinate transformation. Horizontal well Coordinates Horizontal section parallel to X axis Horizontal section parallel to Y axis Vertical well L Y X X M Z Z Y N X Y Z Horizontal well Coordinates Horizontal section parallel to X axis Horizontal section parallel to Y axis Vertical well L Y X X M Z Z Y N X Y Z View Large Table 3. Coordinate transformation. Horizontal well Coordinates Horizontal section parallel to X axis Horizontal section parallel to Y axis Vertical well L Y X X M Z Z Y N X Y Z Horizontal well Coordinates Horizontal section parallel to X axis Horizontal section parallel to Y axis Vertical well L Y X X M Z Z Y N X Y Z View Large 3.2. Validation of the numerical model A multi-stage fractured horizontal well in a three-dimensional shale gas reservoir is made to characterize and simulate the gas flow from nano-micro-scale pores in shales to hydraulic fractures. Figure 3 shows a simplified flow chart of our three-dimensional, multi-phase, multi-mechanism simulator. The dimension of the reservoir is 1500 × 750 × 50 m and the grid size is 114 × 21 × 1 in the x, y and z direction. Figure 3. View largeDownload slide Three-dimensional multi-phase multi-mechanism simulator flow chart. Figure 3. View largeDownload slide Three-dimensional multi-phase multi-mechanism simulator flow chart. To approve the accuracy of the developed numerical model, the simulation results are validated with the simulation results of the model solved by the fully implicit method (figure 4). The results show that the difference in production rates between the two models is less than 5%. The basic parameters for the numerical simulation are listed in table 4. Then, the simulation results are validated with actual field data, and the parameters of the shale gas reservoir and multi-stage fractured horizontal well refer to the literature (Bello and Wattenbarger 2010). The simulated gas production rate is proved to match well with the field data, which is shown in figure 5. Figure 4. View largeDownload slide The matching results of the reservoir simulators solved by different methods. Figure 4. View largeDownload slide The matching results of the reservoir simulators solved by different methods. Table 4. Basic parameters for the numerical simulation. Parameter Value Unit Formation depth 2180 m Formation thickness 50 m Initial pressure 29.65 MPa Initial temperature 80 °C Shale porosity 0.08 Shale permeability 1.5 × 10−4 mD Initial density of rock 2500 kg m−3 Rock compressibility 1.45 × 10−4 MPa−1 Langmuir pressure 4.55 MPa Langmuir volume 4.199 cm3 g−1 Initial gas saturation 0.6 Gas specific gravity 0.65 Horizontal well length 1050 m Hydraulic fracture spacing 130 m Half-length of hydraulic fracture 125 m Wellbore pressure 16 MPa Wellbore radius 0.08 m Parameter Value Unit Formation depth 2180 m Formation thickness 50 m Initial pressure 29.65 MPa Initial temperature 80 °C Shale porosity 0.08 Shale permeability 1.5 × 10−4 mD Initial density of rock 2500 kg m−3 Rock compressibility 1.45 × 10−4 MPa−1 Langmuir pressure 4.55 MPa Langmuir volume 4.199 cm3 g−1 Initial gas saturation 0.6 Gas specific gravity 0.65 Horizontal well length 1050 m Hydraulic fracture spacing 130 m Half-length of hydraulic fracture 125 m Wellbore pressure 16 MPa Wellbore radius 0.08 m View Large Table 4. Basic parameters for the numerical simulation. Parameter Value Unit Formation depth 2180 m Formation thickness 50 m Initial pressure 29.65 MPa Initial temperature 80 °C Shale porosity 0.08 Shale permeability 1.5 × 10−4 mD Initial density of rock 2500 kg m−3 Rock compressibility 1.45 × 10−4 MPa−1 Langmuir pressure 4.55 MPa Langmuir volume 4.199 cm3 g−1 Initial gas saturation 0.6 Gas specific gravity 0.65 Horizontal well length 1050 m Hydraulic fracture spacing 130 m Half-length of hydraulic fracture 125 m Wellbore pressure 16 MPa Wellbore radius 0.08 m Parameter Value Unit Formation depth 2180 m Formation thickness 50 m Initial pressure 29.65 MPa Initial temperature 80 °C Shale porosity 0.08 Shale permeability 1.5 × 10−4 mD Initial density of rock 2500 kg m−3 Rock compressibility 1.45 × 10−4 MPa−1 Langmuir pressure 4.55 MPa Langmuir volume 4.199 cm3 g−1 Initial gas saturation 0.6 Gas specific gravity 0.65 Horizontal well length 1050 m Hydraulic fracture spacing 130 m Half-length of hydraulic fracture 125 m Wellbore pressure 16 MPa Wellbore radius 0.08 m View Large Figure 5. View largeDownload slide The matching result of the gas production rate for shale gas reservoirs. Figure 5. View largeDownload slide The matching result of the gas production rate for shale gas reservoirs. 4. Results and discussion According to the validated simulator for shale gas reservoirs, the influence of submicron effects in multi-scale shales on both the dynamic permeability and then the actual gas production rate are analyzed thoroughly, and the effects of gas–water flow and the main factors of the multi-stage fractured horizontal well on production rate and cumulative production are studied as well. 4.1. Influence of submicron effects in multi-scale shales on dynamic permeability According to the gas phase pressure of each grid block in a shale gas reservoir (figure 6(a)), we can obtain the gas molecular mean free path through equation (6). The pore radius can be calculated through equation (10) or the correlations listed in table 2, as the absolute permeability and porosity of shales can be measured in a laboratory. So, we can calculate the Knudsen number of every grid block in a shale gas reservoir (figure 6(b)) and classify the gas flow regimes, then determine the permeability correction factor of each grid block (figure 6(c)) according to equation (11). The gas flow in a shale gas reservoir is mainly the transition flow and slip flow in figures 6–8, but the pressure of formation will decrease with the continuous production of shale gas, and the Knudsen number will increase, leading to a higher permeability correction factor. So, the gas flow regime is different in formation at different production periods as the Knudsen number is changing all the time, which leads to the dynamic permeability of the shale gas reservoir. Figure 6. View largeDownload slide A distribution map of the parameters in a shale gas reservoir when K0 equals 1 × 10−4 mD after 600 days (nx and ny represent the numbers of grid blocks in the x direction and y direction respectively). Figure 6. View largeDownload slide A distribution map of the parameters in a shale gas reservoir when K0 equals 1 × 10−4 mD after 600 days (nx and ny represent the numbers of grid blocks in the x direction and y direction respectively). Figure 7. View largeDownload slide A distribution map of parameters in a shale gas reservoir when K0 equals 1 × 10−3 mD after 600 days (nx and ny represent the numbers of grid blocks in the x direction and y direction respectively). Figure 7. View largeDownload slide A distribution map of parameters in a shale gas reservoir when K0 equals 1 × 10−3 mD after 600 days (nx and ny represent the numbers of grid blocks in the x direction and y direction respectively). Figure 8. View largeDownload slide A distribution map of parameters in a shale gas reservoir when K0 equals 1 × 10−2 mD after 600 days (nx and ny represent the numbers of grid blocks in the x direction and y direction respectively). Figure 8. View largeDownload slide A distribution map of parameters in a shale gas reservoir when K0 equals 1 × 10−2 mD after 600 days (nx and ny represent the numbers of grid blocks in the x direction and y direction respectively). Figure 6(a) shows the distribution of gas phase pressure when shale permeability equals 1 × 10−4 mD after 600 days, and the pressure of the hydraulic fracture drops faster compared with a shale gas reservoir with higher permeability (figures 7(a) and 8(a)). The lower the permeability of the shale, the lower the pressure in the hydraulic fractures becomes. This makes the gas stored in the formation far from the fractures harder to exploit, leading to a rapid pressure drop nearby. For shale gas reservoirs with lower permeability, the pressure drop funnel cannot spread over a large range, and the pressure of the hydraulic fractures drops faster. This is because it is easier for gas to flow in hydraulic fractures compared with nano-micro-scale pores in shales. 4.2. Effects of dynamic permeability considering submicron effects on gas production Figure 6(b) presents the Knudsen number of each grid block in shale gas reservoirs after 600 days of production. We can conclude that the closer one gets to the hydraulic fractures, the bigger the Knudsen number becomes. Because the Knudsen number is negatively correlated to pressure, so the Knudsen number increases with the decrease of pressure. When the permeability of shales becomes lower, the Knudsen number gets bigger as the pressure drop of the fractures is bigger (figures 6(b), 7(b) and 8(b)). This leads to a rapid increase in the Knudsen number the closer one gets to the hydraulic fractures. The distribution of the permeability correction factor is similar to the Knudsen number (figures 6(c), 7(c) and 8(c)), and the permeability correction factor increases with the Knudsen number, because it is positively related to the Knudsen number according to equation (4). Figure 9 shows the effect of different flow mechanisms on gas production rates. The gas production rate, which considers gas desorption, Knudsen diffusion and the slippage effect, matches well with the field data. If the gas desorption, Knudsen diffusion and slippage effects existing in the nanoscale pores are neglected, the gas production rate of multi-fractured horizontal wells is underestimated. This is the reason why shale gas production modeled with conventional simulators/models is much lower than the actually observed field data (Swami and Settari 2012). Figure 9. View largeDownload slide The effect of different flow mechanisms on gas production rate and cumulative production. Figure 9. View largeDownload slide The effect of different flow mechanisms on gas production rate and cumulative production. Figure 10 presents the gas production rate and cumulative production of multi-fractured horizontal wells under the different permeability of shales. It is shown that the gas production rate and cumulative production both increase with the permeability of shales. However, the increasing rate of production and cumulative production, which consider slippage and Knudsen diffusion in nanoscale pores within shales, decreases as the permeability of shales increases, compared with the production rate and cumulative production without considering slippage or diffusion. This means that the effects of different flow regimes on production rate and cumulative production, considering slippage and diffusion in nanopores, become significant as shale permeability gets lower. Figure 10. View largeDownload slide The effects of slippage and diffusion in nanoscale pores within a shale matrix for different permeability of shales. (The solid lines in this figure represent gas production rate (figure 10(a)) and cumulative production (figure 10(b)), considering the effects of slippage and diffusion. The dashed lines in this figure represent the gas production rate (figure 10(a)) and cumulative production (figure 10(b)) without considering the effects of slippage or diffusion.) Figure 10. View largeDownload slide The effects of slippage and diffusion in nanoscale pores within a shale matrix for different permeability of shales. (The solid lines in this figure represent gas production rate (figure 10(a)) and cumulative production (figure 10(b)), considering the effects of slippage and diffusion. The dashed lines in this figure represent the gas production rate (figure 10(a)) and cumulative production (figure 10(b)) without considering the effects of slippage or diffusion.) Figure 11 illustrates the ratios of the production rate (cumulative production) calculated by the model considering submicron effects in the nanopores (solid lines in figures 10(a) and (b)) and the model without considering submicron effects (dashed lines in figures 10(a) and (b)) over a different production period. The lower the absolute permeability of shales, the bigger the production ratio of the two models. The difference between the ratio of the production rate and the ratio of cumulative production thus increases with the decrease of shale permeability. This is because the influence of submicron effects is bigger in shales with lower permeability. When the shale permeability decreases to 1 × 10−4 mD, the ratio of the production rate increases first and then declines, but is still bigger than the ratios when shale permeability is higher. Gas stored in hydraulic fractures is easier to extract, so pressure declines faster with the increasing difference between the permeability of the hydraulic fractures and the shale matrix. With gas produced continuously, gas stored in the nanopores of the shale matrix starts to be extracted, and adsorbed gas desorbs from the pore walls, which can supplement the gas to pores, so the decrease in pressure becomes slower. The permeability correction factor caused by submicron effects is negatively related to the pressure and pore radius according to equations (1), (4) and (6), so the production rate ratio (when shale permeability equals 1 × 10−4 mD) increases first and then declines, because the pressure declines quickly at first and then the decrease of pressure becomes slower. Figure 11. View largeDownload slide The ratio of production rates (cumulative production) calculated by models considering and without considering submicron effects in nanopores based on figure 10. Figure 11. View largeDownload slide The ratio of production rates (cumulative production) calculated by models considering and without considering submicron effects in nanopores based on figure 10. 4.3. Effect of gas–water flow on gas production Because shale gas wells produce no or little water, so we mainly analyze the effect of initial water saturation, sw, on the production rate and the cumulative production of multi-stage fractured horizontal wells in shale gas reservoirs, as shown in figure 12. Greater initial water saturation means that less free gas is stored in reservoirs, so the production rate declines a lot at the beginning of exploitation as initial water saturation increases. Besides this, greater water saturation causes the lower relative permeability of the gas phase, which is not favorable for gas flow in reservoirs either. Figure 12. View largeDownload slide The effect of initial water saturation on the gas production rate and cumulative production (declining curves represent gas production rate; progressively increasing curves represent cumulative production). Figure 12. View largeDownload slide The effect of initial water saturation on the gas production rate and cumulative production (declining curves represent gas production rate; progressively increasing curves represent cumulative production). Figure 13 shows the relative permeability curves of shales proposed by Li et al (2016a), and figure 14 illustrates the effect of relative permeability on the gas production rate and cumulative production. The production rates are the same at the initial production period of the fractured horizontal well, so the production rate declines faster as the relative permeability decreases (nanopore radius decreases). This is because the initial gas and water stored in reservoirs has not changed, but the gas production rate will reduce as the relative permeability of the gas phase decreases. Figure 13. View largeDownload slide Relative permeability curves with a different nanopore radius R. Reprinted from Li et al (2016a), copyright 2016, with permission from Elsevier. Figure 13. View largeDownload slide Relative permeability curves with a different nanopore radius R. Reprinted from Li et al (2016a), copyright 2016, with permission from Elsevier. Figure 14. View largeDownload slide The effect of the relative permeability on the gas production rate and cumulative production (the declining curves represent the gas production rate; the progressively increasing curves represent cumulative production). Figure 14. View largeDownload slide The effect of the relative permeability on the gas production rate and cumulative production (the declining curves represent the gas production rate; the progressively increasing curves represent cumulative production). 4.4. Effect of main factors of multi-stage fractured horizontal wells on gas production The effect of hydraulic fracture conductivity on the production rate and cumulative production is shown in figure 15. The gas production rate and cumulative production increase with fracture conductivity, and the production rate declines quickly at the beginning of production as the fracture conductivity decreases. The effect of fracture length on the production rate and cumulative production is presented in figure 16, which shows that the gas production rate and cumulative production increase with fracture length. The production rate is the same at the beginning of production, but then declines faster as the fracture length decreases. This is because the fracture conductivity affects the gas flow in the near-wellbore area, but the fracture length starts to affect production only when the pressure spreads to both sides of the hydraulic fractures (the horizontal well goes through the middle of them). Figure 15. View largeDownload slide The effect of hydraulic fracture conductivity on the production rate and cumulative production (the declining curves represent gas production rate; the progressively increasing curves represent cumulative production). Figure 15. View largeDownload slide The effect of hydraulic fracture conductivity on the production rate and cumulative production (the declining curves represent gas production rate; the progressively increasing curves represent cumulative production). Figure 16. View largeDownload slide The effect of hydraulic fracture length, Lf, on the production rate and cumulative production (declining curves represent gas production rate; progressively increasing curves represent cumulative production curves). Figure 16. View largeDownload slide The effect of hydraulic fracture length, Lf, on the production rate and cumulative production (declining curves represent gas production rate; progressively increasing curves represent cumulative production curves). 5. Conclusions This paper presents a numerical model which is valid for the entire Knudsen range (continuum flow, slip flow, transition flow and free molecular flow) in shale gas reservoirs. The effect of gas–water flow and simulation of the hydraulic fracturing cracks were taken into consideration as well. The simulation results were validated with other numerical models and actual field data, and the influencing parameters were analyzed thoroughly. The following conclusions can be drawn: (1) The flow regimes in shale gas reservoirs were analyzed by calculating the Knudsen number under different pressure and different pore throat diameters. Gas flow in shale gas reservoirs with nano-micro-scale pores is mainly continuum flow, slip flow and transition flow. So, for flow regimes other than continuum flow, the traditional equation of Darcy’s law is not applicable anymore. (2) A new numerical model considering multi-scale flow in nano-micro-scale pores within shales was developed based on the Beskok–Karniadakis equation, which is valid for all flow regimes. The simulation results show that the numerical model considering multi-scale flow and gas desorption matches well with the on-site production data, so gas desorption and submicron effects in nanopores cannot be neglected in the simulation. (3) The production rate of multi-stage fractured horizontal wells increases with shale permeability, but the effects of different flow regimes on the production rate and cumulative production, considering slippage and diffusion in nanopores, become significant as the shale permeability gets lower. (4) The initial water saturation has a greater impact on gas production than gas–water relative permeability. Greater initial water saturation means that less free gas is stored in reservoirs, so the production rate declines a lot at the beginning of exploitation as the initial water saturation increases. Greater water saturation causes the lower relative permeability of the gas phase, which is not favorable for gas flow in reservoirs either. Acknowledgments The authors would like to thank the Hubei Cooperative Innovation Center of Unconventional Oil and Gas for supporting our work. This project was partly supported by the National Natural Science Foundation of China (grant nos. 51404037 and 51704032). The authors would like to acknowledge the financial support of Hubei Provincial Department of Education (grant no. Q20171303), the Science and Technology Department of Hubei Province (grant no. 2017CFB146) and the College Students’ Innovation and Entrepreneurship Training Project (grant no. 2017090). References Aguilera R . , 2002 Incorporating capillary pressure, pore throat aperture radii, height above free water table, and Winland r35 values on Pickett plots , AAPG Bull. , vol. 86 (pg. 605 - 624 ) Bello R O , Wattenbarger R A . , 2010 Multi-stage hydraulically fractured shale gas rate transient analysis SPE North Africa Technical Conf. and Exhibition Cairo, Egypt 14–17 February pg. SPE 126754 https://doi.org/10.2118/126754-MS Beskok A , Karniadakis G . , 1999 A model for flows in channels pipes, and ducts at micro and nano scales , Microscale Thermophys. Eng. , vol. 3 (pg. 43 - 77 ) https://doi.org/10.1080/108939599199864 Google Scholar Crossref Search ADS Brown G P , DiNardo A , Cheng G K , Sherwood T K . , 1946 The flow of gases in pipes at low pressure , J. Appl. Phys. , vol. 17 (pg. 802 - 813 ) https://doi.org/10.1063/1.1707647 Google Scholar Crossref Search ADS Brown M , Ozkan E , Raghavan R , Kazemi H . , 2011 Practical solutions for pressure-transient responses of fractured horizontal wells in unconventional shale reservoirs , SPE Reservoir Eval. Eng. , vol. 14 (pg. 663 - 676 ) https://doi.org/10.2118/125043-PA Google Scholar Crossref Search ADS Cipolla C L , Lolon E P , Erdle J C , Rubin B . , 2010 Reservoir modeling in shale-gas reservoirs , SPE Reservoir Eval. Eng. , vol. 13 (pg. 638 - 653 ) https://doi.org/10.2118/125530-PA Google Scholar Crossref Search ADS Civan F . , 2010 Effective correlation of apparent gas permeability in tight porous media , Transport Porous Med. , vol. 82 (pg. 375 - 384 ) https://doi.org/10.1007/s11242-009-9432-z Google Scholar Crossref Search ADS Clarkson C R , Nobakht M , Kaviani D , Ertekin T . , 2012 Production analysis of tight-gas and shale-gas reservoirs using the dynamic-slippage concept , SPE J. , vol. 17 (pg. 230 - 242 ) https://doi.org/10.2118/144317-PA Google Scholar Crossref Search ADS Deng J , Zhu W , Ma Q . , 2014 A new seepage model for shale gas reservoir and productivity analysis of fractured well , Fuel , vol. 124 (pg. 232 - 240 ) https://doi.org/10.1016/j.fuel.2014.02.001 Google Scholar Crossref Search ADS Guggenheim E A . , 1960 , Elements of the Kinetic Theory of Gases Oxford Pergamon Guo J , Zhang L , Wang H , Feng G . , 2012 Pressure transient analysis for multi-stage fractured horizontal wells in shale gas reservoirs , Transport Porous Med. , vol. 93 (pg. 635 - 653 ) https://doi.org/10.1007/s11242-012-9973-4 Google Scholar Crossref Search ADS Heid J G , McMahon J J , Nielsen R F , Yuster S T . , 1950 Study of the permeability of rocks to homogeneous fluids Drilling and Production Practice (New York, 1 January) API-50-230 Huang T , Guo X , Chen F . , 2015a Modeling transient flow behavior of a multiscale triple porosity model for shale gas reservoirs , J. Nat. Gas Sci. Eng. , vol. 23 (pg. 33 - 46 ) https://doi.org/10.1016/j.jngse.2015.01.022 Google Scholar Crossref Search ADS Huang T , Guo X , Chen F . , 2015b Modeling transient pressure behavior of a fractured well for shale gas reservoirs based on the properties of nanopores , J. Nat. Gas Sci. Eng. , vol. 23 (pg. 387 - 398 ) https://doi.org/10.1016/j.jngse.2015.02.020 Google Scholar Crossref Search ADS Huang T , Tao Z , Li E , Lyu Q , Guo X . , 2017 Effect of permeability anisotropy on the production of multi-scale shale gas reservoirs , Energies , vol. 10 pg. 1549 https://doi.org/10.3390/en10101549 Google Scholar Crossref Search ADS Javadpour F . , 2009 Nanopores and apparent permeability of gas flow in mudrocks (shales and siltstone) , J. Can. Petroleum Technol. , vol. 48 (pg. 16 - 21 ) https://doi.org/10.2118/09-08-16-DA Google Scholar Crossref Search ADS Javadpour F , Fisher D , Unsworth M . , 2007 Nanoscale gas flow in shale gas sediments , J. Can. Petroleum Technol. , vol. 46 (pg. 55 - 61 ) https://doi.org/10.2118/07-10-06 Kolodzie S Jr . , 1980 Analysis of pore throat size and use of the Waxman-Smits equation to determine OOIP in Spindle Field, Colorado SPE Annual Technical Conf. and Exhibition Dallas, Texas 21–24 September pg. SPE 9382 https://doi.org/10.2118/9382-MS Li D , Xu C , Wang Y J , Lu D . , 2014 Effect of Knudsen diffusion and Langmuir adsorption on pressure transient response in shale gas reservoir , J. Petrol. Sci. Eng. , vol. 124 (pg. 146 - 154 ) https://doi.org/10.1016/j.petrol.2014.10.012 Google Scholar Crossref Search ADS Li D , Zhang L , Lu D . , 2015 Effect of distinguishing apparent permeability on flowing gas composition, composition change and composition derivative in tight-and shale-gas reservoir , J. Petrol. Sci. Eng. , vol. 128 (pg. 107 - 114 ) https://doi.org/10.1016/j.petrol.2015.02.002 Google Scholar Crossref Search ADS Li T , Song H , Wang J , Wang Y , Killough J . , 2016a An analytical method for modeling and analysis gas-water relative permeability in nanoscale pores with interfacial effects , Int. J. Coal Geol. , vol. 159 (pg. 71 - 81 ) https://doi.org/10.1016/j.coal.2016.03.018 Google Scholar Crossref Search ADS Li D , Zhang L , Wang Y J , Lu D , Du J . , 2016b Effect of adsorption and permeability correction on transient pressures in organic rich gas reservoirs: vertical and hydraulically fractured horizontal wells , J. Nat. Gas Sci. Eng. , vol. 31 (pg. 214 - 225 ) https://doi.org/10.1016/j.jngse.2016.02.033 Google Scholar Crossref Search ADS Ozkan E , Raghavan R S , Apaydin O G . , 2010 Modeling of fluid transfer from shale matrix to fracture network SPE Annual Technical Conf. and Exhibition Florence, Italy 19–22 September pg. SPE 134830 https://doi.org/10.2118/134830-MS Peaceman D W . , 1983 Interpretation of well-block pressures in numerical reservoir simulation with nonsquare grid blocks and anisotropic permeability , SPE J. (pg. 531 - 543 ) https://doi.org/10.2118/10528-PA Roy S , Raju R , Chuang H F , Cruden B A , Meyyappan M . , 2003 Modeling gas flow through microchannels and nanopores , J. Appl. Phys. , vol. 93 (pg. 4870 - 4879 ) https://doi.org/10.1063/1.1559936 Google Scholar Crossref Search ADS Shabro V , Torres-Verdin C , Javadpour F . , 2011 Numerical simulation of shale-gas production: from pore-scale modeling of slip-flow, Knudsen diffusion, and Langmuir desorption to reservoir modeling of compressible fluid SPE North American Unconventional Gas Conf. and Exhibition Woodlands, Texas, USA 14–16 June pg. SPE 144355 https://doi.org/10.2118/144355-MS Swami V , Settari A . , 2012 A pore scale gas flow model for shale gas reservoir Americas Unconventional Resources Conf. Pittsburgh, Pennsylvania, USA 5–7 June pg. SPE 155756 https://doi.org/10.2118/155756-MS Wang H T . , 2014 Performance of multiple fractured horizontal wells in shale gas reservoirs with consideration of multiple mechanisms , J. Hydrol. , vol. 510 (pg. 299 - 312 ) https://doi.org/10.1016/j.jhydrol.2013.12.019 Google Scholar Crossref Search ADS Wang J , Liu H , Wang L , Zhang H , Luo H , Gao Y . , 2015 Apparent permeability for gas transport in nanopores of organic shale reservoirs including multiple effects , Int. J. Coal Geol. , vol. 152 (pg. 50 - 62 ) https://doi.org/10.1016/j.coal.2015.10.004 Google Scholar Crossref Search ADS Wang J , Luo H , Liu H , Cao F , Li Z , Sepehrnoori K . , 2016 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. , vol. 22 (pg. 244 - 264 ) https://doi.org/10.2118/174996-PA Google Scholar Crossref Search ADS Wu K , Chen Z , Li X . , 2016a A model for multiple transport mechanisms through nanopores of shale gas reservoirs with real gas effect—adsorption-mechanic coupling , Int. J. Heat Mass Transf. , vol. 93 (pg. 408 - 426 ) https://doi.org/10.1016/j.ijheatmasstransfer.2015.10.003 Google Scholar Crossref Search ADS Wu K , Li X , Guo C , Wang C , Chen Z . , 2016b A unified model for gas transfer in nanopores of shale-gas reservoirs: coupling pore diffusion and surface diffusion , SPE J. , vol. 21 (pg. 1583 - 1611 ) https://doi.org/10.2118/2014-1921039-PA Google Scholar Crossref Search ADS Wu Y , Moridis G , Bai B , Zhang K . , 2009 A multi-continuum model for gas production in tight fractured reservoirs SPE Hydraulic Fracturing Technology Conf. Woodlands, Texas, USA 19–21 January pg. SPE 118944 https://doi.org/10.2118/118944-MS Yuan Y , Doonechaly N G , Rahman S . , 2016 An analytical model of apparent gas permeability for tight porous media , Transport Porous Med. , vol. 111 (pg. 193 - 214 ) https://doi.org/10.1007/s11242-015-0589-3 Google Scholar Crossref Search ADS Zhang J , Kamenov A , Zhu D , Hill A D . , 2013 Laboratory measurement of hydraulic fracture conductivities in the Barnett Shale SPE Hydraulic Fracturing Technology Conf. Texas, USA 4–6 February pg. SPE 163839 https://doi.org/10.2118/163839-MS Zhao Y , Zhang L , Zhao J , Luo J , Zhang B . , 2013 ‘Triple porosity’ modeling of transient well test and rate decline analysis for multi-fractured horizontal well in shale gas reservoir , J. Petrol. Sci. Eng. , vol. 110 (pg. 253 - 262 ) https://doi.org/10.1016/j.petrol.2013.09.006 Google Scholar Crossref Search ADS Ziarani A S , Aguilera R . , 2012 Knudsen’s permeability correction for tight porous media , Transport Porous Med. , vol. 91 (pg. 239 - 260 ) https://doi.org/10.1007/s11242-011-9842-6 Google Scholar Crossref Search ADS © 2018 Sinopec Geophysical Research Institute TI - A nonlinear seepage model of gas and water transport in multi-scale shale gas reservoirs based on dynamic permeability JO - Journal of Geophysics and Engineering DO - 10.1088/1742-2140/aaae8d DA - 2018-08-01 UR - https://www.deepdyve.com/lp/oxford-university-press/a-nonlinear-seepage-model-of-gas-and-water-transport-in-multi-scale-m2aB5G0l0r SP - 1255 VL - 15 IS - 4 DP - DeepDyve ER -