Analytical modeling and probabilistic evaluation of gas production from a hydraulically fractured shale reservoir using a quad-linear flow model

Journal of Petroleum Science and Engineering 184 (2020) 106516

Contents lists available at ScienceDirect

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

Analytical modeling and probabilistic evaluation of gas production from a hydraulically fractured shale reservoir using a quad-linear flow model Wenxi Ren, Hon Chung Lau * Department of Civil and Environmental Engineering, National University of Singapore, Singapore

A R T I C L E I N F O

A B S T R A C T

Keywords: Shale gas Secondary fractures Analytical model Global sensitivity analysis Probabilistic assessment

During hydraulic fracturing, secondary fractures are created around hydraulic fractures and serve as important channels connecting the shale matrix and the hydraulic fractures. We developed an analytical model to simulate shale gas production from a multi-stage fractured horizontal well. The analytical model is based on the trilinear flow model, but it is free from the assumption that the length of secondary fractures is equal to the half-distance between hydraulic fractures. This analytical model captures the fundamental characteristics of a hydraulically fractured shale reservoir by incorporating transient gas flow in the matrix, secondary fractures, and hydraulic fractures. Because this model describes four linear flow systems (unstimulated reservoir volume, matrix in the stimulated reservoir volume, secondary fracture, and hydraulic fracture), it is called the quad-linear flow model. Production data from two shale gas wells were used to validate the model. Moreover, using the variance-based method and quasi-Monte Carlo simulation, we developed a stochastic framework to perform a probabilistic assessment of shale gas production. We also conducted a global sensitivity study with nine uncertainty param­ eters to identify important parameters and their interactions that affect well productivity. The results indicate that the existence of dense secondary fractures (the distance between secondary fractures is small) results in higher initial production rate and ultimate recovery, but a sharper decline in production rate. Increasing the length of the secondary fractures can more effectively enhance gas production in the presence of dense secondary fractures. Moreover, the benefits of dense secondary fractures are mostly limited to a high porosity shale reservoir.

1. Introduction Hydraulic fracturing is a common practice to enhance the produc­ tivity of shale gas wells. In addition to initiating hydraulic fractures, this practice also activates natural fractures, leading to a complex fracture system including hydraulic fractures and secondary fractures (Cheng et al., 2015; Dahi-Taleghani and Olson, 2011; Gale et al., 2014; Gas­ parrini et al., 2014). In this paper, we use the term secondary fractures to mean induced micro-scale hydraulic fractures and reactivated natural fractures. Microseismic events, tracer data, and production history all indicate the existence of secondary fractures near hydraulic fractures (Sharma and Manchanda, 2015). To date, many researchers have attempted to model such fracture complexity. Brown et al. (2011) introduced an analytical trilinear flow solution to model the production performance of a multi-stage fractured horizontal well. They used a dual porosity model to simulate the fluid transfer from shale matrix to the natural fractures and assumed that the dual porosity region is

represented by alternative slab-shaped matrix blocks and parallel hori­ zontal fractures (natural fractures). The trilinear flow model is limited to the cases in which the stimulated reservoir region completely covers the region between adjacent hydraulic fractures. However, the model pro­ vides a practical method to investigate the effect of the fracture network on well productivity. Subsequently, the multi-linear flow approach was widely used to investigate pressure and rate-transient responses of shale gas horizontal wells. Dehghanpour and Shirdel (2011) developed a trilinear flow model by proposing that each matrix block consists of sub-matrices which are surrounded by microfractures. Stalgorova and Mattar (2012) extended the work of Brown et al. (2011) and developed a modified trilinear flow model in which secondary fractures partially cover the region between adjacent hydraulic fractures. Moreover, the secondary fractures around each hydraulic fracture is represented by an enhanced permeability region and the remaining reservoir is unstimu­ lated. Later, Stalgorova and Mattar (2013) developed a more flexible linear flow model which considers five linear flow regions instead of three. Further, Zeng et al. (2017, 2018) used a spherical particle model

* Corresponding author. E-mail address: [email protected] (H.C. Lau). https://doi.org/10.1016/j.petrol.2019.106516 Received 4 April 2019; Received in revised form 24 August 2019; Accepted 20 September 2019 Available online 21 September 2019 0920-4105/© 2019 Elsevier B.V. All rights reserved.

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

Nomenclature Symbols pini pw T

ρg

v2 t ϕ2 k2 Ct2 Ct1 Ctf CtF Cg Ct2a Ct1a

ρr μ

na n0 pL p2 s v1 ϕ1 k1 vf ϕf kf vF

ϕF kF wF xF nF num H hf hft hm x1 x2 q Q S1 ST N Sij

initial reservoir pressure, MPa wellbore pressure, MPa reservoir temperature, K gas density, mol⋅m 3 gas velocity in region 2, m⋅s 1 time, s matrix porosity (region 2) matrix permeability (region 2), nd total compressibility of the unstimulated reservoir, MPa 1 total compressibility of the matrix in the stimulated reservoir, MPa 1 total compressibility of secondary fracture, MPa 1 total compressibility of hydraulic fracture, MPa 1 gas compressibility, MPa 1 modified total compressibility of the unstimulated reservoir, MPa 1 modified total compressibility of the stimulated reservoir, MPa 1 rock density, kg⋅m 3 gas viscosity, Pa⋅s adsorption amount, mol⋅kg 1 maximum adsorption capacity, mol⋅kg 1 Langmuir pressure, MPa reservoir pressure (region 2), MPa Laplace variable gas velocity in region 1, m⋅s 1 matrix porosity (region 1) matrix permeability (region 1), nd gas velocity in secondary fracture, m⋅s 1 secondary fracture porosity secondary fracture permeability, md gas velocity in hydraulic fracture, m⋅s 1

ω’ μ1, μ2 σ 1, σ2 μ3 σ3

hydraulic fracture porosity hydraulic fracture permeability, md width of hydraulic fracture, m half-length of hydraulic fracture, m number of hydraulic fractures number of secondary fractures reservoir thickness, m width of secondary fracture, m total width of secondary fractures, m distance between secondary fractures, m length of secondary fracture, m half distance between hydraulic fractures, m gas production rate, 104 m3 d 1 cumulative gas production, 107 m3 first-order sensitivity index total sensitivity index sample size second-order sensitivity index weighted coefficient of the dual-Gaussian distribution mean of the dual-Gaussian distribution, nd standard deviation of the dual-Gaussian distribution, nd mean of the Gaussian distribution standard deviation of the Gaussian distribution

Abbreviations HF hydraulic fracture SF secondary fracture SRV stimulated reservoir volume URV unstimulated reservoir volume MC Monte Carlo QMC quasi-Monte Carlo VBSA variance-based sensitivity analysis DOE Design of Experiment SAFE Sensitivity Analysis For Everybody PDF probability distribution function

to represent the stimulated reservoir region and developed a seven linear flow model which considers fluid supply from upper and lower reservoirs. Heidari Sureshjani and Clarkson (2015) and Al-Rbeawi (2017) extended the trilinear flow model of Brown et al. (2011) by considering the asymmetrical placement of a hydraulic fracture in its drainage area. Wang et al. (2015) and Fan and Ettehadtavakkol (2017) introduced the fractal theory to characterize the heterogeneity of the stimulated reservoir region and developed fractal linear-flow models. Apart from the linear flow approach, the semi-analytical model based on the Green function and Newman product method has been used to model gas production from multiple fractured horizontal wells (He et al., 2017, 2019). Moreover, discrete fracture modeling (Jiang and Younis, 2016; Wang, 2018; Xu et al., 2017; Yu et al., 2018) and nodal analysis approach (Chen et al., 2016, 2018; Yang et al., 2016; Yu et al., 2017) with the aid of microseismic data have been used to simulate gas flow in fractured shale gas reservoirs. However, microseismic data only inform the modeler of the location of failures during hydraulic fracturing. It does not give direct information on the connected stimulated volume during production. Above all, the linear flow approach is a practical method for production evaluation, but neglects the nonlinearity of desorption and gas properties (compressibility and viscosity). The semi-analytical model, discrete fracture modeling, and nodal analysis approach are more sophisticated than the linear flow approach, but they are less attractive due to their computation cost. Since linear flow is the dominant flow regime in fractured shale gas wells (Brohi et al., 2011; Ozkan et al., 2011) and probabilistic evaluation of shale gas production involves intense computation, this work uses the linear flow approach.

Currently, a comprehensive analysis of the effect of secondary frac­ tures on production performance is scarce. Following Brown et al. (2011), we developed a quad-linear flow model in which the length of secondary fractures does not have to be half the distance between hy­ draulic fractures as assumed by Brown et al. (2011). In addition, based on variance-based analysis and quasi-Monte Carlo simulation, we developed an integrated probabilistic framework to investigate how production performance depends on reservoir and fracture parameters. 2. Model development The proposed model follows the trilinear flow model proposed by Brown et al. (2011). An outline of the proposed model is presented in the following section (for a detailed derivation please see the supplementary material). The linear flow approach originates from the analytical so­ lution of Cinco-Ley and Meng (1988) which is well accepted in the in­ dustry (Brohi et al., 2011; Brown et al., 2011; Fan and Ettehadtavakkol, 2016; Ozkan et al., 2011; Stalgorova and Mattar, 2013; Zeng et al., 2018). This analytical approach assumes a 1D flow in each flow region. The flux from the nearby region into this region is treated as a source term, regardless of flow direction. These assumptions, though approxi­ mate, are necessary to render the problem soluble analytically. If these assumptions are not made, the solution can only to be solved numeri­ cally which is not our approach in this paper. To simulate a shale gas reservoir depleted by a multi-stage fractured horizontal well, we assume that equal-length hydraulic fractures are uniformly spaced along the horizontal wellbore and penetrate the formation completely, as shown 2

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

in Fig. 1(a). Flow contribution from the region beyond the tip of hy­ draulic fracture is assumed to be negligible. During hydraulic fracturing, induced microfractures propagate in different directions and interact with reactivated natural fractures resulting from the tensile or shear failure of pre-existing natural fractures, thus creating a second fracture network around each hydraulic fracture, as shown in Fig. 1(a). This region is idealized as a dual porous medium consisting of matrix blocks and vertical secondary fractures. The secondary fractures are idealized as planar fractures and equally spaced by a distance, hm, along the hy­ draulic fracture. The dual-porosity assumption is only a first approxi­ mation, but a practical approach to model a fracture network in tight unconventional reservoirs (Brown et al., 2011; Ozcan et al., 2014). Ac­ cording to the state-of-the-art techniques of monitoring and modeling the fracture propagation during fracturing treatments, real fracture networks are more complex than planar shape (Li et al., 2018; Weng et al., 2014; Wu et al., 2012). But the complex fracture networks are usually not known a priori. Moreover, simplifications are needed to establish an analytical model for describing gas transport in complex fracture networks. For these reasons, we used the simple fracture network in this study. Some researchers (Albinali and Ozkan, 2016; Liu et al., 2018; Raghavan and Chen, 2013, 2019) suggested to use anom­ alous diffusion to characterize gas transport in fractured rocks. This method seeks to describe flow in a complex fracture system without knowing the intrinsic details of the system (Albinali and Ozkan, 2016; Liu et al., 2018). We reserve the investigation of anomalous diffusion for future research. We call the dual porous medium the stimulated reservoir volume (SRV), as shown in Fig. 1(b). Correspondingly, the unstimulated reser­ voir volume is abbreviated as URV. Moreover, gas flow is considered as single-phase, isothermal flow. In our model, the following assumptions of gas transport are made. Gas flows from the URV along the x*-direction into the matrix of the SRV, but no into the tip of the secondary fracture. In other words, the tip of the secondary fractures is sealed. Gas then flows from the matrix of the SRV into the secondary fracture in the y*direction. This gas then enters into the hydraulic fracture. There is no flow of gas from the matrix of the SRV directly into the hydraulic frac­ tures. These assumptions are consistent with those made by other au­ thors. (Brown et al., 2011; Chen et al., 2016). The direction of gas transport is shown in Fig. 1(b). Because we assumed 1D linear flow in the

four regions (URV, matrix in the SRV, secondary fracture, and hydraulic fracture), our model is called the quad-linear flow model. We couple the linear flow solutions for the SRV (region 1), URV (region 2), and hydraulic fracture region by imposing pressure and flow continuity boundary conditions. Note that we do not differentiate gas flow in organic matter and inorganic matrix, which is due to the fact that kerogen is highly dispersed in shale matrix and not expected to be continuous at the macroscopic scale (Chen, 2016; Olorode et al., 2017). Moreover, stress sensitivity and multiphase flow will be considered in future research. The diffusion-type equation describing gas flow in region 2 is � ∂ � ∂ ∂ ρv ¼ ϕ ρ þ ½ð1 ∂x g 2 ∂t 2 g ∂t v2 ¼

k2

(1)

rp2

(2)

p2 pL þ p2

(3)

μ

na2 ¼ n0

ϕ2 Þρr na2 �

The Langmuir equation is used to model gas adsorption/desorption on shale. Note that Darcy’s law rather than the slip-flow theory is used to describe gas transport in the shale matrix. This is because petroleum engineers and scientists have not reached a consensus on this issue. Supporters of the slip-flow theory believe that gas molecules slip on nanopore surfaces, which enhances gas transport (Cai et al., 2018; Darabi et al., 2012; Javadpour, 2009; Ren et al., 2016; Wu et al., 2015). However, some researchers indicate that the slip-flow theory is based on the rarefied gas dynamics (dilute gas), which is not suitable to describe gas flow in shale because shale gas is highly compressed in typical reservoir conditions (Chen and Shen, 2018a, b; Patzek et al., 2013). Thus, we reserve this issue for further research and use the Darcy’s law, as was done by Brown et al. (2011). Equation (1) can be rewritten in dimensionless term and converted to the Laplace domain for simplicity. The definition of dimensionless variables is given in the Appendix.

∂2 m2D s ¼ m ∂xD2 η2D 2D

(4)

Fig. 1. Plan schematics of (a) a multi-stage fractured horizontal well and (b) the quad-linear flow model. The inset in Fig. 1(b) shows a local-coordinate system describing gas flow from shale matrix to secondary fracture. Arrows indicate flow directions. Gas flows from the URV along the x-direction to the matrix in the SRV, as shown by the blue arrows. Gas flows from the matrix in the SRV along the y*-direction to the secondary fractures, as shown by black arrows. All gas entering the hydraulic fracture comes from the secondary fractures in which the flow direction is along the x-direction, as shown by the green arrows. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.) 3

W. Ren and H.C. Lau

η2D ¼

Journal of Petroleum Science and Engineering 184 (2020) 106516

k2 ϕf Ctf kf ϕ2 Ct2a

�

(5)

where m refers to pseudo pressure (Al-Hussainy et al., 1966), a hori­ zontal bar over a variable denotes the variable in the Laplace domain, the subscript D refers to dimensionless variables, and s is the Laplace variable. Cti (i ¼ 1, 2, f, or F) is the total compressibility. Assuming that the total compressibility is dominated by gas compressibility (Nobakht and Clarkson, 2012), Cti � Cg. Moreover, when considering gas desorp­ tion effect, Ct2 is modified as Ct2a ¼ Ct2 þ

ϕ2 ÞZ ρr n0 pL RT

ð1

The modified total compressibility (Ct2a) and gas viscosity (μ) are evaluated at pm, which is equal to (pini þ pw)/2, as was done by Sure­ shjani et al. (2016) and Fuentes-Cruz and Valko (2015). The boundary conditions are � ∂m2D �� ¼0 (7) ∂xD � xD ¼ x2D � � m2D ��

xD ¼ x1D

� � � ¼ m1D � * � xD ¼ x1D

� � ∂ ∂ ∂ ρ v þ F21 ¼ ϕ ρ þ ½ð1 ∂y* g 1 ∂t 1 g ∂t

ϕ1 Þρr na1 �

� ∂ � ∂ ρv ¼ ϕρ ∂x g f ∂t f g

ω¼

ϕ1 Ct1a hm ϕf Ctf hf

η2D

η2D

Ct1a ¼ Ct1 þ

ð1

ϕ1 ÞZ ρr n0 pL RT 2

ϕ1 p1 ðpL þ p1 Þ

ρg k1 ∂p1 �� hf μ ∂y* � y* ¼ 1=2hm

ff ðsÞ ¼ 1 þ

(10)

(21)

λ 3s

(22)

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ! 3ωs 3ωs αtanh α λ λ

with the boundary conditions � ∂mfD �� ¼0 ∂xD � xD ¼ x1D � � mfD ��

xD ¼ 1=2wFD

� � ¼ mFD ��

xD ¼ 1=2wFD

(23)

(24) (25)

The solution to Eq. (22) with the boundary conditions given by Eqs. (24) and (25) is hpffiffiffiffiffiffiffiffiffiffi i � cosh ff ðsÞs ðxD x1D Þ � � � � �� (26) mfD ¼ mFD � pffiffiffiffiffiffiffiffiffiffi 1 xD ¼ 1=2wFD cosh ff ðsÞs 2wFD x1D

(13)

is

The diffusion equation describing gas flow in the hydraulic fracture � ∂ � ∂ ρv ¼ ϕ ρ þ FfF ∂y g F ∂t F g

(14)

k1 hm x2F λ ¼ 12 kf hf h2m

2

� ∂2 mfD � þ ff ðsÞs mfD ¼ 0 ∂x2D

where the subscript * refers to the local-coordinate system describing gas flow from the matrix in the SRV to the secondary fracture, as shown in Fig. 1(b). F21 characterizes the flow contribution from region 2 (URV region). Equation (10) can be converted to the Laplace domain and rewritten in terms of dimensionless variables to yield � � ∂2 m1D 3ω þ α (12) s m1D ¼ 0 *2 λ ∂yD � x2D Þ

(20)

where F1f characterizes the flow contribution from the matrix in region 1 (the matrix feeds the secondary fractures). The dimensionless form of Eq. (20) in the Laplace domain is

(11)

�rffiffiffiffiffiffi rffiffiffiffiffiffi s s tanh ðx1D

F1f

�

F1f ¼

�

ρg k2 ∂p2 �� x1 μ ∂x � x ¼ x1

1 h2m k2 4 k1 x1 xF

(18)

(8)

This solution is consistent with the solution of Stalgorova and Mattar (2012) given in Eq. (B-8). Region 2 continues to supply gas to region 1, which exerts a signif­ icant effect on the pressure distribution of region 1. Region 1 consists of secondary fractures and the matrix. The diffusion equation describing gas flow in the matrix is given by

α¼

� � � � � � m1D � * ¼ mfD �� hm � yD ¼ 1 � yD ¼ 2xF

When α is equal to 0 (neglecting the flow contribution from the re­ gion 2). Equation (19) reduces to the solution of Serra et al. (1983) given in Eq. (A-6). The diffusion-type equation describing gas flow in the secondary fracture is

The solution to Eq. (4) with the boundary conditions described by Eqs. (7) and (8) is hpffiffiffiffiffiffiffiffiffiffiffi i � � ðx cosh s= η x Þ D 2D 2D � h ffiffiffiffiffiffiffiffiffiffi i m2D ¼ m1D � * (9) � xD ¼ x1D cosh pffis= η2D ðx1D x2D Þ

F21 ¼

(17)

The modified total compressibility (Ct1a) is evaluated at pm, which is equal to (pini þ pw)/2. The solution to Eq. (12) with the boundary conditions given by Eqs. (17) and (18) is hpffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi i � � cosh 3ωs=λ α y*D � hpffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi i (19) m1D ¼ mfD �� hm � yD ¼ 2x cosh 3ωs=λ α F

(6)

ϕ2 p2 ðpL þ p2 Þ2

∂m1D �� ¼0 � ∂y*D � y*D ¼ 0

(15)

FfF ¼

(16)

kfb ¼

� 2kfb ρg ∂pf � � wF μ ∂x � x ¼ 1=2wF hft kf xF

(27) (28) (29)

where FfF characterizes gas transfer between secondary fracture and hydraulic fracture (the secondary fractures feed gas into the hydraulic fracture). Note that hft ¼ num � hf where num is the number of the sec­ ondary fractures (num ¼ xF/hm) and hf is the width of the secondary

with the boundary conditions

4

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

fracture. Using Laplace transforms, we solved Eq. (27) under a constant rate condition. The wellbore-pressure solution is � pffiffi � � π ¼ 0 pffiffi coth β (30) mFD �� yD ¼ 0 sCFD β 0

CFD ¼ 0

kf ¼ kf β¼2

kF wF 0 kf xF

(31)

hft H

(32)

� � pffiffiffiffiffiffiffiffiffiffi kfb xF pffiffiffiffiffiffiffiffiffiffi ff ðsÞstanh ff ðsÞs x1D k F wF

1 wFD 2

�� þ

s

ηFD

Table 1 Model parameters for calculating the dimensionless cumulative production. Value

Parameters

Value

pini T

30 MPa 350 K 2500 kg m 3 0.1 MPa 50 m 0.07 mol kg 3 16 MPa 570 nd 0.06 0.06 20 m

k2 x2 kf ϕf hf kF ϕF xF wF hm

570 nd 50 m 1 md 0.33 3 � 10 4 m 100 md 0.3 150 m 3 � 10 3 m 20 m

ρrock pw H n0 pL k1 ϕ1 ϕ2 x1

(33)

where CFD is the dimensionless hydraulic fracture conductivity. Equa­ tion (30) is similar to the solution of Brown et al. (2011) given in Eq. (41). Following Van Everdingen and Hurst (1949), the solution to Eq. (27) under a constant wellbore pressure is pffiffi 0 pffiffi � βCFD 1 � qD ¼ ¼ tanh β (34) � π s s2 mFD �� yD ¼ 0 0

Moreover, the cumulative gas production is pffiffi 0 pffiffi � βCFD qD QD ¼ nF ¼ nF tanh β s π s2

Parameters

constant pressure (Beygi and Rashidi, 2011; Fuentes-Cruz and Valko, 2015; Sureshjani et al., 2016). This analytical solution is an approximate one. But in practical cases, the wellbore pressure is not likely to be at ambient atmosphere pressure (0.1 MPa), which means that the realistic change of the total compressibility and gas viscosity is insignificant as assumed in the simulation case. Thus, our model is acceptable for simulating shale gas production. Secondly, our analytical solution was history matched against pub­ lished production data from two wells. One was a shale gas well in Marcellus (Meyer et al., 2010). This was a horizontal well with a seven stage hydraulic fracturing. However, stage two showed a negligible contribution to total production. Thus, only six hydraulic fractures (each stage is treated as a hydraulic fracture) were considered in history matching. The production data covering a period of 200 days are shown in Fig. 3(a). Matrix porosity and permeability, and fracture parameters were adjusted to achieve a match. The model parameters are listed in Table 2. Fig. 3(a) shows that the model results agree well with the field data. The fitted value of k1 is 575 nd, which is consistent with the result of Meyer et al. (2010) (580 nd). The fitted value of xF is 100 m, which is close to the result of Meyer et al. (2010) (97.54 m). Moreover, the values of the other fitted parameters are within reasonable ranges. This verifies the practical utility of the analytical model. The other well was from China (Xu et al., 2015). The specific location of this well is not reported in the literature. This horizontal well was stimulated with three stages hydraulic fracturing and had produced for more than 500 days. Table 2 lists the fracture and reservoir parameters. Fig. 3(b) shows a comparison of model results with the production data. The fitted values of matrix porosity and permeability, and fracture pa­ rameters are within reasonable ranges. Note that although a good match is obtained, the match is non-unique, as is true for many history matching exercises in reservoir engineering. Therefore, the engineer should exercise caution in the history matching process, such as using various combinations of parameters for the history matching. Engi­ neering judgement and comparison with analog data may provide quality check on the history matching process.

(35)

where nF is the number of hydraulic fractures. Equations (34) and (35) are transformed into real-time by using the Gaver-Stefhest algorithm (Davies, 2001). 3. Validation Firstly, we checked whether the present model satisfies the material balance requirements. We plotted the dimensionless cumulative pro­ duction (defined by the ratio of cumulative production to gas-initiallyin-place) as a function of time, as shown in Fig. 2. The model parame­ ters are listed in Table 1. Fig. 1 shows that the dimensionless cumulative production for the case without desorption converges to 0.96, and the dimensionless cumulative production for the case with desorption con­ verges to 0.93. This is due to the fact that the linear flow approach ne­ glects the nonlinearity of desorption and gas properties, namely, the total compressibility and gas viscosity are usually evaluated at a

Variance-based sensitivity analysis We investigated the sensitivity of shale gas production using the variance-based sensitivity analysis (VBSA). Unlike the conventional sensitivity analyses (local methods), VBSA is a global sensitivity analysis method, which provides first-order and higher-order (interactions) sensitivity indices (Gatelli et al., 2009). Moreover, unlike Design of Experiment (DOE) in which the effects of input parameters are estimated over levels, VBSA uses the whole input space by using stochastic methods, such as Monte Carlo (MC) method (Saltelli et al., 2010). More details on VBSA can be found in Saltelli et al. (2008). We developed a Matlab toolbox for VBSA based on the VBSA module of SAFE (Sensitivity Analysis For Everybody) toolbox (Pianosi et al., 2015), which is freely available for non-commercial users. We made three improvements to the current VBSA module. First, quasi-Monte Carlo (QMC) sampling method

Fig. 2. Dimensionless cumulative production with time. The dimensionless cumulative production is the ratio of cumulative production to gas-initially-inplace. The model parameters are listed in Table 1. 5

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

Table 2 Model parameters for history matching. Parameters

Marcellus shale a

pini T pw H kF nF wF xF ϕF k1 k2 ϕ1 ϕ2 kf hf hm x1 x2 ϕf n0 pL a b c

32.6 352.6a 3.69a 49.38a 300c 6a 3 � 10 100c 0.33c 575c 575c 0.042a 0.042a 5.85c 3 � 10 20c 16c 50c 0.21c 0.17c 8.85c

China shale b

33.6 360b 16b 49b 650c 3b 3 � 10 200b 0.33c 150c 150c 0.03c 0.03c 1.4c 3 � 10 10c 43c 50a 0.21c 0.08b 8b

3

4c

Unit MPa K MPa m md

3c

m m nd nd

4c

md m m m m mol⋅kg MPa

1

Value is taken from Meyer et al. (2010). Value is taken from Xu et al. (2015). Fitting parameter.

interaction between inputs xi and xj on model output. The Ishigami function (Ishigami and Homma, 1990) and the Sobol function (Saltelli and Sobol, 1995) were used to test the toolbox. The two functions exhibit strong non-linearity and non-monotonicity and are widely used as examples for sensitivity analysis methods. The Sobol function is fS ðXÞ ¼

8 Y j4Xi i¼1

(37)

where a1 ¼ 0, a2 ¼ 1, a3 ¼ 4.5, a4 ¼ 9, a5 ¼ 99, a6 ¼ 99, a7 ¼ 99, a8 ¼ 99, and Xi are input parameters and independent random variables uni­ formly distributed in [0, 1]. For this function the first-order sensitivity indices can be calculated analytically (Gatelli et al., 2009). The theo­ retical values are shown in Table 3. A base sample size of 25 � 104 is used, which takes about 31.36 s on a desktop computer with an Intel Xeon E5-1650 processor at 3.20 GHz and 16 GB of RAM. The calculation results are consistent with the theoretical values, as shown in Table 3. The Ishigami function is

Fig. 3. History matching results with the analytical model. The field data are from (a) a shale gas well in Marcellus and (b) a shale gas well in China. The parameter set is given in Table 2.

(38)

fI ðxÞ ¼ sinðx1 Þ þ asin2 ðx2 Þ þ bx43 sinðx1 Þ

where a ¼ 7, b ¼ 0.1, x1, x2, and x3 are input parameters and indepen­ dent random variables uniformly distributed in [ π, π]. A base sample size of 20 � 104 is used, which takes about 1.50 s on the desktop com­ puter. The calculation results are in good agreement with the theoretical values, as shown in Table 4.

was used to generate sample matrices. QMC computations are more efficient than MC (Sobol and Levitan, 1999). We used the well-known Sobol sequence to generate quasi-random samples, as was done by Sobol and Levitan (1999). Second, first-order sensitivity index Si and total effect index STi were estimated using the Saltelli’s scheme (Saltelli et al., 2010), which is more efficient than the Sobol scheme (Saltelli and Sobol, 1993) and the Homma scheme (Homma and Saltelli, 1996) used in SAFE (Saltelli et al., 2010). Si and STi measure the first-order and total effects (first order and higher order effects) of input xi on model output, respectively. STi is given by (Saltelli et al., 2008) X X STi ¼ Si þ Sij þ Sijl þ … (36) i6¼j

2j þ ai 1 þ ai

Table 3 Comparison of calculated sensitivity indices with exact values.

i6¼j6¼l

where Sij and Sijl denote the second- and third-order sensitivity indices, P P Sijl respectively. On the right-hand side of Eq. (36), the term ( Sij þ þ …) represents all higher-order effects due to interactions. If Si ¼ STi, there is no interaction among inputs. Third, we developed a code to calculate the second-order sensitivity index Sij (current SAFE only pro­ vides Si and STi) and integrate the code in SAFE. Sij measures the effect of

Sensitivity index

Theoretical value

Calculated

S1 S2 S3 S4 S5 S6 S7 S8

0.72a 0.18a 0.02a 0.01a 0a 0a 0a 0a

0.72 0.18 0.02 0.01 0 0 0 0

a

6

Value is taken from Gatelli et al. (2009).

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

Table 4 Comparison of calculated sensitivity indices with exact values. Sensitivity index

Theoretical value

S1 S2 S3 S12 S23 S13 a

Calculated

a

0.31 0.44a 0a 0a 0a 0.24a

0.31 0.44 0 0 0 0.24

Value is taken from Homma and Saltelli (1996).

5. Results and discussion VBSA with QMC was performed on nine uncertain parameters to determine the key-input parameters that influence cumulative gas pro­ duction (model output). The input parameters were taken from field practice and available experimental data (Table 5). The matrix perme­ ability data were from Javadpour et al. (2007), who measured the permeability of 152 shale samples. The measurements are shown in Fig. 4 and the distribution can be approximated by a dual-Gaussian distribution. Although the dual-Gaussian distribution underestimates the proportion of high permeability samples (>300 nd), Javadpour et al. (2007) indicated that 90% of the measurements are less than 150 nd. Thus, this deviation does not limit our model. The other model param­ eters are listed in Table 6. The present study is not aimed at optimizing completion designs but evaluating the effect of secondary fractures on shale gas production. Thus, we fixed the region drained by secondary fractures (SRV þ URV), i.e., that the half-length of hydraulic fracture xF and the half distance between hydraulic fractures x2 were set to be constant. The input parameters listed in Table 5 are independent of each other and sampled by using QMC. The generated samples were fed into the analytical model, which calculated cumulative gas production at different times. When all simulations were completed, the output was used to calculate sensitivity indices. Note that the sample size (N) is a key parameter for VBSA. Increasing the sample size results in more precise estimates but increases the computation cost. A convergence check gave N to be 1 � 106, as shown in Fig. 5, S1 and S2 in the sup­ plementary material. One million model runs were performed based on the parameters listed in Tables 5 and 6. Fig. 6 shows the cumulative gas production as a function of time. It can be seen that the uncertainty spans four orders of magnitude. Cumulative gas production at 1 year ranges from 1.54 � 105 m3 to 27.96 � 107 m3 with a P50 of 6.68 � 107 m3. The maximum 30-year cumulative gas production is 74.61 � 107 m3, which is about 1000 times greater than the minimum value (6.36 � 105 m3). Results in Fig. 6 were then used to calculate the P10/P90 ratio. Since

Fig. 4. Probability versus shale permeability. The measured permeability dis­ tribution is approximated with a dual-Gaussian distribution, with fitting pa­ rameters given in Table 5.

Fig. 5. Estimates of the first-order and total sensitivity indices of hm and x1 as a function of sample size. The model output is cumulative gas production at 1 year. The parameter set is given in Tables 5 and 6. The first-order and total sensitivity indices of other parameters show the similar tendency and thus are not shown in this figure. Table 6 Other model parameters for production evaluation and sensitivity analysis.

Table 5 Uncertain parameters for production evaluation and sensitivity analysis. Parameter

PDF

Range

kf kF hf hm x1 k1a ϕ1b n0c pLc

Uniform Uniform Uniform Uniform Uniform Dual-Gaussian Gaussian Uniform Uniform

0.01–50 md 100–1 � 104 md 1–500 μm 2–75 m 5–45 m 0–400 nd >0 0.14–0.23 mol kg 14–20 MPa

Parameters

Value

Parameters

Value

pini T

33.6 MPa 360 K 2500 kg m 16 MPa 49.38 m

nF wF xF ϕf x2

25 3 � 10 150 m 0.32 50 m

ρrock pw H

3

3

m

P90 and P10 are low and high estimates, respectively. The P10/P90 ratio characterizes the difference in production performance between the best and worst cases. The P10/P90 ratio is useful in quantifying the uncer­ tainty in production performance. A large ratio is indicative of a large difference in production performance between the best and worst cases, implying higher uncertainty (for a small ratio it is vice versa). Fig. S3 shows a decrease of the P10/P90 ratio with time, which is expected because our simulation geometry is a closed reservoir. Fig. S3 also shows

1

a We used a dual-Gaussian distribution (μ1 ¼ 58.46 nd, σ 1 ¼ 23.77 nd, μ2 ¼ 94.44 nd, σ2 ¼ 85.11 nd, ω’ ¼ 0.54) to approximate the measured perme­

ability distribution. b The range and PDF (μ3 ¼ 0.057, σ 4 ¼ 0.021) of porosity were taken from Kale et al. (2010). c Langmuir parameters were taken from Yang et al. (2015).

7

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

Fig. 6. Cumulative gas production from a multi-stage fractured horizontal well over 30 years. The parameter set is given in Tables 5 and 6. P10, P50, and P90 for cumulative gas production were calculated based on the results of the one million simulations.

Fig. 8. First-order and total sensitivity indices of cumulative gas production at (a) 1 year, (b) 10 years, and (c) 30 years.

The first-order sensitivity index measures the influence of an individual parameter on model output. The total sensitivity index, including the influence of parameter interactions, quantifies the total contribution of the individual parameter to model output. A parameter with a small first-order sensitivity index but a large total sensitivity index, influences model output primarily through interaction with other parameters. Fig. 8(a) shows that at 1 year of production the most significant parameter is the distance between secondary fractures (hm). The second most significant parameter is the permeability of secondary fracture (kf), followed by the width of secondary fracture (hf), the length of secondary fracture (x1), matrix porosity (ϕ1), and permeability (k1). Hydraulic fracture permeability (kF), the maximum adsorption capacity (n0), and the Langmuir pressure (pL) are insignificant parameters. Moreover, the first-order and total sensitivity indices for all the parameters are not equal, which means that the parameters interact with each other. Spe­ cifically, the first-order and total sensitivity indices for hm are 0.51 and 0.65, respectively, which means that the sum of higher-order sensitivity indices is 0.14 based on Eq. (36). Fig. 8(b) shows that at 10 years of production the most significant parameter is still hm. ϕ1 becomes the second most significant parameter, followed by kf, x1, hf, and k1. The influences of kF, n0 and pL on cumulative gas production are still insig­ nificant. At 30 years of production, ϕ1 becomes the most significant parameter, and hm becomes the second most significant parameter, as shown in Fig. 8(c). This is because secondary fractures serve as channels connecting shale matrix and hydraulic fractures. Free gas from hydraulic fractures and secondary fractures are the dominant contributions to gas production at initial production stage. With the depletion of free gas in hydraulic fractures and secondary fractures, the matrix starts to supply gas into secondary fractures, which is limited by the low matrix permeability. Thus, at the late production stage, the effect of secondary fractures decreases. Fig. 8 also shows that the most significant fracture parameter affecting production is hm. To clarify this further, we calculated cumu­ lative gas production by only changing hm. Fig. 9(a) shows that the early time drainage is significantly affected by hm. When the distance between secondary fractures is relatively narrow (hm ¼ 3 m), the gas recovery is much faster compared to that when hm ¼ 75 m (relatively sparse sec­ ondary fractures). For example, if hm ¼ 3 m, 69.85% of the 30-year cu­ mulative gas production is produced in the first five years. In contrast, when hm ¼ 75 m, only 28.09% of the 30-year cumulative gas production is produced in the first five years. The 30-year cumulative gas produc­ tion is generally regarded as an approximation of estimated ultimate recovery (EUR). Thus, the existence of dense secondary fractures (small hm) can result in higher initial gas production rate and more effective

that the P10/P90 ratio decreases sharply at an early stage of production and starts to flatten after 5 years. Specifically, P10/P90 ratio is 12.54 at 1 year of production, 6.47 at 5 years of production, and 3 at 30 years of production, which indicates that shale gas production is governed by different factors at different production stages. This issue will be detailed later. The 30-year cumulative gas production distributions for the cases with and without secondary factures are shown in Fig. 7. For the case without secondary fractures, all gas entering the hydraulic fracture comes from the region 1. The region 2 feeds gas into the region 1. The two distributions can be approximated by a Gaussian distribution. The 30-year cumulative gas production distribution for the case with sec­ ondary fractures shows a shift toward right with an increasing vari­ ability. This rightward shift shows that secondary fractures can enhance gas production significantly. The mean value for the case with secondary fractures is about 25.70 � 107 m3, which is around 1.11 times that for the case without secondary fractures. Next, the simulation results were used to calculate the sensitivity indices. The first-order and total sensitivity indices are shown in Fig. 8.

Fig. 7. 30-year cumulative gas production distributions for the cases with and without secondary factures. HF refers to the case without secondary fractures, and HF þ SF refers to the case with secondary fractures. 8

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

Fig. 9. Effect of hm on (a) cumulative gas production and (b) production rate. The reservoir parameters are listed in Table 6. The other model parameters are kf ¼ 10 md, kF ¼ 1000 md, hf ¼ 50 μm, ϕ1 ¼ ϕ2 ¼ 4%, x1 ¼ 30 m, k1 ¼ k2 ¼ 100 nd, n0 ¼ 0.15 mol kg 1, pL ¼ 18 MPa.

drainage. In other words, decreasing hm not only accelerates gas re­ covery, but also increases EUR. However, the existence of dense sec­ ondary fractures results in a sharper decline in production rate, as shown

in Fig. 9(b). After 15 years, the production rate for hm ¼ 3 m is 5.29 � 103 m3 d 1, which is even lower than that for hm ¼ 75 m. Above all, the geometry of secondary fractures is important for shale

Fig. 10. Second-order sensitivity indices of cumulative gas production at (a) 1 year, (b) 10 years, and (c) 30 years. 9

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

When hm ¼ 50 m cumulative gas production at 1 year is 1.14 � 107 m3 for x1 ¼ 40 m, which is 1.13 times more than that for x1 ¼ 10 m. Thus, increasing x1 can more effectively enhance gas production in the pres­ ence of dense secondary fractures. At 10 years of production, the most significant interaction is between hm and ϕ1, as shown in Fig. 10(b). Fig. 11(b) shows that for ϕ1 ¼ 0.09, the cumulative gas production at 10 years decreases with increasing hm; while for ϕ1 ¼ 0.03, the cumulative gas production at 10 years decreases slightly with increasing hm. This is due to the fact that increasing hm can more effectively accelerate gas production from a high porosity reser­ voir than from a low porosity reservoir. From a practical viewpoint, when designing hydraulic fracturing, it is important to identify the lo­ cations of high porosity shale along horizontal wellbore by well logging. The typical industry practice is to create equally-spaced hydraulic fractures along the horizontal wellbore. An improvement may be the use of gas content and rock mechanic properties, such as porosity and brittleness index, to select the sweet spots for multi-stage fracturing (Stegent et al., 2013). Fig. 10(c) shows that at 30 years of production the most significant interaction is between hm and ϕ1, which is the same as that at 10 years. Fig. 11(c) shows that the benefits of dense secondary fractures (small hm) are mostly limited to a high porosity shale reservoir, which also reflects the fact that reservoir parameters gradually control shale gas production with increasing time.

gas production. However, current methods cannot accurately describe the fracture geometry. The commonly used microseismic analysis cannot characterize fracture connectivity (Kumar and Sharma, 2018). Tracer-test method and flow-back analysis can provide the information of the secondary fracture network connected to hydraulic fractures (Zolfaghari et al., 2016) but pay little attention to the distribution of secondary fractures. Previous researchers have developed different models to interpret flowback data (Xu et al., 2016). Magnetic nano­ particles are considered to be potential nano-sensors for detecting fracture networks in unconventional reservoirs (An et al., 2017; Yekeen et al., 2019). The magnetic nanoparticles are used as contrast agents. The information of fracture networks may be obtained from the differ­ ence in the magnetic susceptibility before and after injecting magnetic nanoparticles. kf is the second most significant fracture parameter affecting pro­ duction. Dense secondary fractures can effectively increase the SRV. However, the aperture of these secondary fractures (micron-sized) is much smaller than that of the hydraulic fractures (millimeter-sized). Regular proppants used in hydraulic fracturing are 100 mesh (149 μm) or larger and are too large to enter and support the secondary fractures. Research should be directed to come up with micro-sized proppants that can enter and support the secondary fractures (Lau et al., 2019). Our results also show that the initial shale gas production is mainly controlled by fracture parameters (hm, x1, kf, hf, and kF), while at the late stage gas production is mainly controlled by reservoir parameters (ϕ1, n0, and pL). This is because fracture parameters affect well productivity but not affect the gas-initially-in-place, which is only dependent on reservoir parameters. Good completion quality can accelerate gas pro­ duction, which also means faster decline of production rate. Here, completion quality refers to fracture complexity, which is related to the fracture parameters hm, x1, kf, hf, and kF. Thus, the influence of fracture parameters on cumulative gas production decreases with increasing time. But the influence of gas desorption on cumulative gas production increases over time. This is because adsorbed gas gradually desorbs with depletion of reservoir pressure. Fig. 10 shows the interactions among the nine parameters with a color legend. Different colors are used to indicate different values of Sij. Fig. 10(a) shows that at 1 year of production, the interaction between hm and x1 is more important than others. Fig. 11(a) shows that for small hm (<10 m), the cumulative gas production is sensitive to x1, while for large hm (>10 m) the cumulative gas production is insensitive to x1. Specifically, when hm ¼ 3 m, the cumulative gas production at 1 year is 10.21 � 107 m3 for x1 ¼ 40 m, which is 1.37 times that for x1 ¼ 10 m.

6. Conclusions We developed an analytical quad-linear flow model to simulate shale gas production. We also developed a global sensitivity analysis toolbox based on the variance based method and quasi-Monte Carlo simulation. By integrating the toolbox and the analytical model, we performed a probabilistic evaluation and global sensitivity analysis of shale gas production. One million model runs were run spanning the probability distribution of model parameters. Simulation results indicate that: (1) The creation of secondary fractures consisting of secondary hy­ draulic fractures and reactivated natural fractures around hy­ draulic fractures can increase the SRV significantly, thus leading to both higher initial production rate and ultimate recovery. (2) Initial production rate and EUR can be improved with dense secondary fractures (i.e. secondary fractures close to each other), especially for a higher porosity (e.g. ϕ1 > 5%) shale reservoir.

Fig. 11. (a) Interaction between hm and x1 on cumulative gas production at 1 year. The reservoir parameters are listed in Table 6. The other model parameters are kf ¼ 25 md, kF ¼ 1000 md, hf ¼ 50 μm, hm ¼ 10 m, ϕ1 ¼ ϕ2 ¼ 5%, k1 ¼ k2 ¼ 30 nd, n0 ¼ 0.18 mol kg 1, pL ¼ 15 MPa. (b) Interaction between hm and ϕ1 on cumulative gas production at 10 years, and (c) interaction between hm and ϕ1 on cumulative gas production at 30 years. The reservoir parameters are listed in Table 6. The other model parameters are kf ¼ 25 md, kF ¼ 1000 md, hf ¼ 50 μm, hm ¼ 10 m, x1 ¼ 20 m, k1 ¼ 30 nd, n0 ¼ 0.18 mol kg 1, pL ¼ 15 MPa. 10

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

(3) The initial gas production is mainly controlled by fracture pa­ rameters (hm, x1, kf and kF) while the gas production at late time is controlled by reservoir porosity, the Langmuir pressure, and the maximum adsorption capacity. (4) The developed probabilistic framework combing quasi-Monte Carlo simulation and the quad-linear flow model can be used to evaluate shale gas production and determine key parameters that affect well productivity. This, in turn, provides valuable infor­ mation to improve stimulation design.

Acknowledgments We specially thank Prof. Andrea Saltelli (University of Bergen) for helpful discussions on VBSA. We thank Prof. Salam Al-Rbeawi (METU Northern Cyprus Campus) and Dr. Jie Zeng (University of Western Australia) for many discussions on the linear-flow model.

Appendix xD ¼

x xF

(A-34)

x1D ¼

x1 xF

(A-35)

x2D ¼

x2 xF

(A-36)

yD ¼

y xF

(A-37)

y*D ¼

y* 1=2hm

(A-38)

kF ϕf Ctf kf ϕF CtF

(A-39)

kf t ϕf Ctf μx2F

(A-41)

ηFD ¼ tD ¼

mID ¼

qD ¼

πkf hft Tsc qsc psc T

1 qsc psc T ¼ mFD πkf hft Tsc ðmFi

QD ¼

(A-42)

ΔmI I ¼ m; f or F

(A-43)

mF Þ

Qsc psc T

πϕf Ctf hft Tsc x2F ðmFi

wFD ¼

(A-44)

mF Þμ

wF xF

(A-45)

Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi.org/10.1016/j.petrol.2019.106516.

References

Brown, M., Ozkan, E., Raghavan, R., et al., 2011. Practical solutions for pressuretransient responses of fractured horizontal wells in unconventional shale reservoirs. SPE Reserv. Eval. Eng. 14 (6), 663–676. Cai, J., Lin, D., Singh, H., et al., 2018. Shale gas transport model in 3D fractal porous media with variable pore sizes. Mar. Pet. Geol. 98, 437–447. Chen, C., 2016. Multiscale imaging, modeling, and principal component analysis of gas transport in shale reservoirs. Fuel 182, 761–770. Chen, K.P., Shen, D., 2018. Drainage flow of a viscous compressible fluid from a small capillary with a sealed end. J. Fluid Mech. 839, 621–643. Chen, K.P., Shen, D., 2018. Mechanism of fluid production from the nanopores of shale. Mech. Res. Commun. 88, 34–39. Chen, Z., Liao, X., Yu, W., et al., 2018. Pressure-transient behaviors of wells in fractured reservoirs with natural-and hydraulic-fracture networks. SPE J. 24 (1), 375–394. Chen, Z., Liao, X., Zhao, X., et al., 2016. A semianalytical approach for obtaining type curves of multiple-fractured horizontal wells with secondary-fracture networks. SPE J. 21 (2), 538–549. Cheng, W., Jin, Y., Chen, M., 2015. Reactivation mechanism of natural fractures by hydraulic fracturing in naturally fractured shale reservoirs. J. Nat. Gas Sci. Eng. 23, 431–439.

Al-Hussainy, R., Ramey Jr., H., Crawford, P., 1966. The flow of real gases through porous media. J. Pet. Technol. 18 (5), 624–636. Al-Rbeawi, S., 2017. How much stimulated reservoir volume and induced matrix permeability could enhance unconventional reservoir performance. J. Nat. Gas Sci. Eng. 46, 764–781. Albinali, A., Ozkan, E., 2016. Anomalous diffusion approach and field application for fractured nano-porous reservoirs. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. An, C., Yan, B., Alfi, M., et al., 2017. Estimating spatial distribution of natural fractures by changing NMR T2 relaxation with magnetic nanoparticles. J. Pet. Sci. Eng. 157, 273–287. Beygi, M.E., Rashidi, F., 2011. Analytical solutions to gas flow problems in low permeability porous media. Transp. Porous Media 87 (2), 421–436. Brohi, I.G., Pooladi-Darvish, M., Aguilera, R., 2011. Modeling fractured horizontal wells as dual porosity composite reservoirs-application to tight gas, shale gas and tight oil cases. In: SPE Western North American Region Meeting. Society of Petroleum Engineers.

11

W. Ren and H.C. Lau

Journal of Petroleum Science and Engineering 184 (2020) 106516

Cinco-Ley, H., Meng, H.-Z., 1988. Pressure transient analysis of wells with finite conductivity vertical fractures in double porosity reservoirs. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. Dahi-Taleghani, A., Olson, J.E., 2011. Numerical modeling of multistranded-hydraulicfracture propagation: accounting for the interaction between induced and natural fractures. SPE J. 16 (3), 575–581. Darabi, H., Ettehad, A., Javadpour, F., et al., 2012. Gas flow in ultra-tight shale strata. J. Fluid Mech. 710, 641–658. Davies, B., 2001. Integral Transforms and Their Applications, third ed. Springer. Dehghanpour, H., Shirdel, M., 2011. A triple porosity model for shale gas reservoirs. In: Canadian Unconventional Resources Conference. Society of Petroleum Engineers. Fan, D., Ettehadtavakkol, A., 2016. Transient shale gas flow model. J. Nat. Gas Sci. Eng. 33, 1353–1363. Fan, D., Ettehadtavakkol, A., 2017. Analytical model of gas transport in heterogeneous hydraulically-fractured organic-rich shale media. Fuel 207, 625–640. Fuentes-Cruz, G., Valko, P.P., 2015. Revisiting the dual-porosity/dual-permeability modeling of unconventional reservoirs: the induced-interporosity flow field. SPE J. 20 (1), 124–141. Gale, J.F., Laubach, S.E., Olson, J.E., et al., 2014. Natural fractures in shale: a review and new observations. AAPG Bull. 98 (11), 2165–2216. Gasparrini, M., Sassi, W., Gale, J.F., 2014. Natural sealed fractures in mudrocks: a case study tied to burial history from the Barnett Shale, Fort Worth Basin, Texas, USA. Mar. Pet. Geol. 55, 122–141. Gatelli, D., Kucherenko, S., Ratto, M., et al., 2009. Calculating first-order sensitivity measures: a benchmark of some recent methodologies. Reliab. Eng. Syst. Saf. 94 (7), 1212–1219. He, Y., Cheng, S., Li, S., et al., 2017. A semianalytical methodology to diagnose the locations of underperforming hydraulic fractures through pressure-transient analysis in tight gas reservoir. SPE J. 22 (3), 924–939. He, Y., Cheng, S., Qin, J., et al., 2019. Interference testing model of multiply fractured horizontal well with multiple injection wells. J. Pet. Sci. Eng. 176, 1106–1120. Heidari Sureshjani, M., Clarkson, C.R., 2015. An analytical model for analyzing and forecasting production from multifractured horizontal wells with complex branchedfracture geometry. SPE Reserv. Eval. Eng. 18 (3), 356–374. Homma, T., Saltelli, A., 1996. Importance measures in global sensitivity analysis of nonlinear models. Reliab. Eng. Syst. Saf. 52 (1), 1–17. Ishigami, T., Homma, T., 1990. An Importance Quantification Technique in Uncertainty Analysis for Computer Models, First International Symposium on Uncertainty Modeling and Analysis. IEEE. Javadpour, F., 2009. Nanopores and apparent permeability of gas flow in mudrocks (shales and siltstone). J. Can. Pet. Technol. 48 (8), 16–21. Javadpour, F., Fisher, D., Unsworth, M., 2007. Nanoscale gas flow in shale gas sediments. J. Can. Pet. Technol. 46 (10), 55–61. Jiang, J., Younis, R.M., 2016. Hybrid coupled discrete-fracture/matrix and multicontinuum models for unconventional-reservoir simulation. SPE J. 21 (3), 1009–1027. Kale, S.V., Rai, C.S., Sondergeld, C.H., 2010. Petrophysical Characterization of Barnett Shale, SPE Unconventional Gas Conference. Society of Petroleum Engineers. Kumar, A., Sharma, M.M., 2018. Diagnosing fracture-wellbore connectivity using chemical tracer flowback data. In: Unconventional Resources Technology Conference. Society of Petroleum Engineers. Lau, H.C., Radhamani, A.V., Ramakrishna, S., 2019. Maximizing production from shale reservoir by using micro-sized proppants. In: International Petroleum Technology Conference. Society of Petroleum Engineers. Li, J., Yu, W., Wu, K., 2018. Analyzing the impact of fracture complexity on well performance and wettability alteration in Eagle Ford shale. In: Unconventional Resources Technology Conference. Society of Petroleum Engineers. Liu, S., Li, H., Valk� o, P.P., 2018. A Markov-Chain-based method to characterize anomalous diffusion phenomenon in unconventional reservoir. In: SPE Canada Unconventional Resources Conference. Society of Petroleum Engineers. Meyer, B.R., Bazan, L.W., Jacot, R.H., et al., 2010. Optimization of multiple transverse hydraulic fractures in horizontal wellbores. In: SPE Unconventional Gas Conference. Society of Petroleum Engineers. Nobakht, M., Clarkson, C.R., 2012. A new analytical method for analyzing linear flow in tight/shale gas reservoirs: constant-flowing-pressure boundary condition. SPE Reserv. Eval. Eng. 15 (3), 370–384. Olorode, O.M., Akkutlu, I.Y., Efendiev, Y., 2017. Compositional reservoir-flow simulation for organic-rich gas shale. SPE J. 22 (6), 1963–1983. Ozcan, O., Sarak, H., Ozkan, E., et al., 2014. A trilinear flow model for a fractured horizontal well in a fractal unconventional reservoir. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. Ozkan, E., Brown, M.L., Raghavan, R., et al., 2011. Comparison of fractured-horizontalwell performance in tight sand and shale reservoirs. SPE Reserv. Eval. Eng. 14 (2), 248–259. Patzek, T.W., Male, F., Marder, M., 2013. Gas production in the Barnett Shale obeys a simple scaling theory. Proc. Natl. Acad. Sci. U. S. A. 110 (49), 19731–19736. Pianosi, F., Sarrazin, F., Wagener, T., 2015. A Matlab toolbox for global sensitivity analysis. Environ. Model. Softw 70, 80–85. Raghavan, R., Chen, C.-C., 2013. Fractional diffusion in rocks produced by horizontal wells with multiple, transverse hydraulic fractures of finite conductivity. J. Pet. Sci. Eng. 109, 133–143.

Raghavan, R., Chen, C.-C., 2019. Evaluation of fractured rocks through anomalous diffusion. In: SPE Western Regional Meeting. Society of Petroleum Engineers. Ren, W., Li, G., Tian, S., et al., 2016. An analytical model for real gas flow in shale nanopores with non-circular cross-section. AIChE J. 62 (8), 2893–2901. Saltelli, A., Annoni, P., Azzini, I., et al., 2010. Variance based sensitivity analysis of model output. Design and estimator for the total sensitivity index. Comput. Phys. Commun. 181 (2), 259–270. Saltelli, A., Ratto, M., Andres, T., et al., 2008. Global Sensitivity Analysis: the Primer. John Wiley & Sons. Saltelli, A., Sobol, I., 1993. Sensitivity analysis for non-linear mathematical models. Math Model. Comput. Exp. 1 (4), 407–414. Saltelli, A., Sobol, I.M., 1995. About the use of rank transformation in sensitivity analysis of model output. Reliab. Eng. Syst. Saf. 50 (3), 225–239. Serra, K., Reynolds, A., Raghavan, R., 1983. New pressure transient analysis methods for naturally fractured reservoirs (includes associated papers 12940 and 13014). J. Pet. Technol. 35 (12), 2,271-272,283. Sharma, M., Manchanda, R., 2015. The role of induced un-propped (IU) fractures in unconventional oil and gas wells. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. Sobol, I., Levitan, Y.L., 1999. On the use of variance reducing multipliers in Monte Carlo computations of a global sensitivity index. Comput. Phys. Commun. 117 (1), 52–61. Stalgorova, E., Mattar, L., 2012. Practical analytical model to simulate production of horizontal wells with branch fractures. In: SPE Canadian Unconventional Resources Conference. Society of Petroleum Engineers. Stalgorova, K., Mattar, L., 2013. Analytical model for unconventional multifractured composite systems. SPE Reserv. Eval. Eng. 16 (3), 246–256. Stegent, N.A., Wagner, A., Stringer, C., et al., 2013. Engineering approach to optimize development strategy in the oil segment of the eagle ford shale: a case study. SPE Prod. Oper. 28 (3), 226–234. Sureshjani, M.H., Behmanesh, H., Soroush, M., et al., 2016. A direct method for property estimation from analysis of infinite acting production in shale/tight gas reservoirs. J. Pet. Sci. Eng. 143, 26–34. Van Everdingen, A., Hurst, W., 1949. The application of the Laplace transformation to flow problems in reservoirs. J. Pet. Technol. 1 (12), 305–324. Wang, H., 2018. Discrete fracture networks modeling of shale gas production and revisit rate transient analysis in heterogeneous fractured reservoirs. J. Pet. Sci. Eng. 169, 796–812. Wang, W., Shahvali, M., Su, Y., 2015. A semi-analytical fractal model for production from tight oil reservoirs with hydraulically fractured horizontal wells. Fuel 158, 612–618. Weng, X., Kresse, O., Chuprakov, D., et al., 2014. Applying complex fracture model and integrated workflow in unconventional reservoirs. J. Pet. Sci. Eng. 124, 468–483. Wu, K., Li, X., Wang, C., et al., 2015. A model for gas transport in microfractures of shale and tight gas reservoirs. AIChE J. 61 (6), 2079–2088. Wu, R., Kresse, O., Weng, X., et al., 2012. Modeling of interaction of hydraulic fractures in complex fracture networks. In: SPE Hydraulic Fracturing Technology Conference. Society of Petroleum Engineers. Xu, J., Guo, C., Wei, M., et al., 2015. Production performance analysis for composite shale gas reservoir considering multiple transport mechanisms. J. Nat. Gas Sci. Eng. 26, 382–395. Xu, Y., Cavalcante Filho, J.S., Yu, W., et al., 2017. Discrete-fracture modeling of complex hydraulic-fracture geometries in reservoir simulators. SPE Reserv. Eval. Eng. 20 (2), 403–422. Xu, Y., Ezulike, O.D., Zolfaghari, A., et al., 2016. Complementary surveillance microseismic and flowback data analysis: an approach to evaluate complex fracture networks. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. Yang, F., Ning, Z., Zhang, R., et al., 2015. Investigations on the methane sorption capacity of marine shales from Sichuan Basin, China. Int. J. Coal Geol. 146, 104–117. Yang, R., Huang, Z., Yu, W., et al., 2016. A comprehensive model for real gas transport in shale formations with complex non-planar fracture networks. Sci. Rep. 6, 36673. Yekeen, N., Padmanabhan, E., Idris, A.K., et al., 2019. Nanoparticles applications for hydraulic fracturing of unconventional reservoirs: a comprehensive review of recent advances and prospects. J. Pet. Sci. Eng. 178, 41–73. Yu, W., Wu, K., Sepehrnoori, K., et al., 2017. A comprehensive model for simulation of gas transport in shale formation with complex hydraulic-fracture geometry. SPE Reserv. Eval. Eng. 20 (3), 547–561. Yu, W., Xu, Y., Weijermars, R., et al., 2018. A numerical model for simulating pressure response of well interference and well performance in tight oil reservoirs with complex-fracture geometries using the fast embedded-discrete-fracture-model method. SPE Reserv. Eval. Eng. 21 (2), 489–502. Zeng, J., Wang, X., Guo, J., et al., 2017. Composite linear flow model for multi-fractured horizontal wells in heterogeneous shale reservoir. J. Nat. Gas Sci. Eng. 38, 527–548. Zeng, J., Wang, X., Guo, J., et al., 2018. Composite linear flow model for multi-fractured horizontal wells in tight sand reservoirs with the threshold pressure gradient. J. Pet. Sci. Eng. 165, 890–912. Zolfaghari, A., Dehghanpour, H., Ghanbari, E., et al., 2016. Fracture characterization using flowback salt-concentration transient. SPE J. 21 (1), 233–244.

12