# Transient Pressure and Rate Decline Analysis for Horizontal Well in Stress-Sensitive Composite Reservoir

Authors: Yu Huang; Xiaoping Li; Xiaohua Tan
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

The existence of stress sensitivity effect made the percolation mechanism of low-permeability reservoirs significantly complex. Further, numerous reservoirs have composite stratum properties in actual development procedure. This paper based on the concept of permeability stress sensitivity presents an unsteady flow model for horizontal well taking both stress sensitivity and composite reservoir into account. Analytic solutions for the transient pressure and the rate decline behaviors are obtained by comprehensive utilization of regular perturbation method, Laplace transformation, orthogonal transformation, and Stehfest numerical inversion. The example analysis verifies that the proposed model is reliable and practical. Likewise, there is a discussion of the influence of permeability modulus and other relevant parameters on the transient pressure and the rate decline for horizontal well in stress-sensitive composite reservoir. The work of this paper improved the previous researches and provided a more accurate basis for transient flow analysis and formation evaluation of this typical reservoir.


---

## Body

## 1. Introduction

The development of low-permeability gas reservoirs has received considerable attention because of the decrease of conventional oil and gas resources. However, the existence of “stress sensitivity” in low-permeability reservoirs made its percolation mechanism significantly complex. Furthermore, the formation acidizing or the characteristics of the reservoir itself cause the formation parameters near wellbore to be different from the area far from the wellbore. This kind of reservoir is called composite gas reservoir. However, researches taking both stress sensitivity and composite reservoir into account have been rarely seen.Ren and Guo (2014) described the relationship between permeability and pressure as a power function and built a numerical solution model for transient pressure analysis with consideration of stress sensitivity [1]. Li et al. (2014) presented an eigenvalue method of well test analysis for horizontal gas well with consideration of different inner and outer boundaries [2]. Zhang (2015) built an unstable seepage model of horizontal well with consideration of gas slippage effect and analyzed the impact of slippage effect on transient pressure performance of unstable seepage [3]. Guo et al. (2015) proposed a semianalytical solution method to analyze the fluid flow in the fractured horizontal well and discussed the influence of various parameters on the transient pressure behavior [4]. Zhao et al. (2016) assumed that permeability is an exponential function of pressure and built a nonlinear model to analyze the behavior of horizontal well in permeability stress-sensitive reservoir [5]. Daolun et al. (2016) presented an unsteady flow model with consideration of starting pressure gradient for low-permeability reservoir by using the numerical method [6]. Li-Na et al. (2016) built a method with consideration of starting pressure for horizontal well transient rate decline analysis by employing an orthogonal transformation method [7]. Ren et al. (2017) offered a semianalytical solution to the pressure distribution for fractured horizontal well of stress-sensitive oil reservoirs [8]. Luo et al. (2017) proposed a new semianalytical method to simulate the flow behavior of vertical fractured wells and analyzed the effect of permeability modulus on pressure and pressure derivative curve [9].From all the researches above, we discovered that although some scholars have studied the percolation mechanism in low-permeability gas reservoir regarding the stress sensitivity, these researches did not consider the situation of formation heterogeneity.Wang et al. (2010) through sink-source superposition and integral transformation obtained a new model with consideration of dual permeability in composite reservoir and discussed the influence of different outer boundaries on the transient pressure [10]. Imad et al. (2011) analyzed the pressure dynamic behaviors of fractured horizontal well in a dual-porosity inner zone and single-porosity outer zone composite reservoir [11]. Nie et al. (2011) presented a model for transient well test in a multiple-zone composite reservoir with consideration of negative skin effect [12]. Wang et al. (2013) deduced a nonlinear fluid flow model of multiple-zone composite reservoir with consideration of the quadratic pressure gradient term and plotted the pressure dynamic type curves of every zone [13]. Zhang et al. (2014) analyzed the continuous point source functions of composite reservoir under three kinds of outer boundaries and drew the type curves of flow behavior [14].Although the mathematical models for composite reservoirs are well developed, there is almost no report regarding the percolation analysis for the horizontal wells considering both formation heterogeneity and stress sensitivity.This paper proposed an unsteady flow model for horizontal well with consideration of stress sensitivity in composite reservoir and acquired the analytical solution of the model by utilizing fluid mechanics in porous medium, mathematical physical methods, and fathoming strategies for complex science. We discussed pressure transient behavior and watched its distinctive flow stages and analyzed the rate decline characteristics of horizontal well in composite gas reservoir. The proposed model is more in accordance with the actual conditions; it can provide a more accurate guidance for the formation evaluation and the development of this type of reservoirs.


## 2. Methodology

### 2.1. Physical Model

A physical model of a composite gas reservoir is shown in Figure1. To make the proposed model more specific and easier to understand, the following assumptions and descriptions are made:Figure 1
Physical models for horizontal well of composite reservoir.① The top and bottom boundaries of the reservoir are closed, the horizontal boundary is infinite, and the thickness ish. The inner radius is rm, and the initial formation pressure is pi.② The vertical and horizontal permeability areKv and Kh, respectively. The horizontal length is 2L, the horizontal section is parallel to the vertical boundaries of the gas reservoir, and its location is unlimited and represented by zw.③ A slightly compressible gas has a constant viscosity and compressibility factor. The single-phase gas flow obeys a Darcy law, and the rock permeability has stress sensitivity and it occurs mainly in the near wellbore zone.④ There is no additional pressure drop in the inner and outer zone interfaces. The impact of gravity and the capillary force are ignored, but the wellbore storage and the skin factor are taken into consideration.


### 2.2. Mathematical Model

(1)
The relationship between reservoir permeability and pressure can be expressed as(1)K=Kie-γpi-p (2)
Motion equation:The equation of motion without consideration of permeability stress sensitivity is as follows:(2)vr=-3.6Khμ∂p∂rvz=-3.6Kvμ∂p∂zThe equation of motion with consideration of permeability stress sensitivity is as follows:(3)vr=-3.6Khμe-γpi-p∂p∂rvz=-3.6Kvμe-γpi-p∂p∂z(3)
State equation:(4)φ=φ01+Cφp-p0ρ=ρ01+Cρp-p0ρ=MpRTZ(4)
Continuity equation:(5)∇·ρv⃑+∂ρφ∂t=0We only consider the effect of stress sensitivity in the inner zone since the stress sensitivity mainly occurs in the near wellbore area. Combining (3)–(5) the differential equation of percolation that accounts for the stress sensitivity of inner zone can be expressed as follows:Inner zone0≤r≤rm:(6)Kh1r∂∂rrp1μ1Ze-γpi-p1∂p1∂r+Kv1∂∂zp1μ1Ze-γpi-p1∂p1∂z=φ1μ1Ct13.6p1μ1Z∂p1∂tCombining (2) and (4)-(5) the differential equation of percolation that accounts for the stress sensitivity of outer zone can be expressed as follows:Outer zonerm≤r<∞: (7)Kh2r∂∂rrp2μ2Z∂p2∂r+Kv2∂∂zp2μ2Z∂p2∂z=φ2μ2Ct23.6p2μ2Z∂p2∂tDefine the pseudopressure:(8)mp=∫0p2pμZdpDefine the permeability modulus with pseudopressure: (9)K=Kie-γmmi-mIntroducing the definition of dimensionless variables (Table 1), the dimensionless equations of the differential equation of percolation of composite gas reservoir with consideration of stress sensitivity effect can be expressed as follows:Inner zone0≤rD≤rmD: (10)∂2m1D∂rD2+1rD∂m1D∂rD-γmD∂m1D∂rD2+LD2∂2m1D∂zD2-γmD∂m1D∂zD2=hDLD2eγmDm1D∂m1D∂tDOuter zonermD≤rD≤∞: (11)1rD∂∂rDrD∂m2D∂rD+LD2∂2m2D∂zD2=M12ω12hDLD2∂m2D∂tDInitial condition:(12)m1DtD=0=m2DtD=0=0Inner boundary condition:(13)limεD→0⁡limrD→0⁡∫zwD-εD/2zwD+εD/2rDe-γmDm1D∂m1D∂rDdzwD=-12zD-zwD≤εD2Interface condition:(14)m1DrD=rmD=m2DrD=rmD(15)∂m1D∂rDrD=rmD=1M12∂m2D∂rDrD=rmDInfinite horizontal boundary condition:(16)m2DrD→∞=0Closed top and bottom boundaries:(17)∂m1D∂zDzD=0,1=0,∂m2D∂zDzD=0,1=0Equations (10)–(17) are the mathematical model of composite gas reservoir with consideration of stress sensitivity.Table 1
Definitions of the dimensionless variables.
Dimensionless variable Equation Dimensionless variable Equation Time t D = 3.6 K h 1 t φ 1 μ 1 C t 1 r w 2 Storage ratio ω 12 = φ 1 φ 2 Radial distance r D = r L Mobility ratio M 12 = λ 1 λ 2 = K 1 K 2 Vertical distance z D = z h Storage coefficient C D = 0.159 C φ 1 C t 1 h L 2 Horizontal length L D = L h K v K h Permeability modulus γ m D = T q s c 78.489 K h 1 h γ m Reservoir thickness h D = h r w K h K v Pseudo pressure of inner zone m 1 D = 78.489 K h 1 h T q s c m i - m 1 Horizontal section position z w D = z w h Pseudo pressure of outer zone m 2 D = 78.489 K h 1 h T q s c m i - m 2


### 2.3. Solutions of the Mathematical Model

Eq. (10), the differential equation of percolation with consideration of stress sensitivity, is a nonlinear equation and cannot be solved directly. Therefore, we use Pedrosa substitution and regular perturbation method (Guo 2009) [15] to eliminate the nonlinearity of the equation.Pedrosa substitution:(18)m1DrD,zD,tD=-1γmDln⁡1-γmDξ1DrD,zD,tDwhereξ1DrD,zD,tD is the perturbed deformation function.Hence, we can achieve the following expressions:(19)e-m1DγmD=1-γmDξ1D∂m1D∂rD=11-γmDξ1D∂ξ1D∂rD∂2m1D∂rD2=11-γmDξ1D∂2ξ1D∂rD2+γmD1-γmDξ1D2∂ξ1D∂rD2∂m1D∂zD=11-γmDξ1D∂ξ1D∂zD∂2m1D∂zD2=11-γmDξ1D∂2ξ1D∂zD2+γmD1-γmDξ1D2∂ξ1D∂zD2∂m1D∂tD=11-γmDξ1D∂ξ1D∂tDAccording to regular perturbation method, the following items of dimensionless permeability modulus can be expanded in power series:(20)ξ1D=ξ1D0+γmDξ1D1+γmD2ξ1D2+⋯11-γmDξ1D=1+γmDξ1D+γmD2ξ1D2+⋯-1γmDln1-γmDξ1D=ξ1D+12γmDξ1D2+⋯Since the dimensionless permeability modulus is usually very small (γmD≪1), the zeroth order perturbation solution can meet the requirements of the engineering precision. Combining (19)–(20), the mathematical model (10)–(17) can be modified as(21)1rD∂∂rDrD∂ξ1D0∂rD+LD2∂2ξ1D0∂zD2=hDLD2∂ξ1D0∂tD0≤rD≤rmD1rD∂∂rDrD∂m2D∂rD+LD2∂2m2D∂zD2=M12ω12hDLD2∂m2D∂tDrmD≤rD≤∞ξ1D0tD=0=m2DtD=0=0limεD→0limrD→0⁡∫zwD-εD/2zwD+εD/2rD∂ξ1D0∂rDdzwD=-12zD-zwD≤εD2ξ1D0rD=rmD=m2DrD=rmD∂ξ1D0∂rDrD=rmD=1M12∂m2D∂rDrD=rmDm2DrD→∞=0∂ξ1D0∂zDzD=0,1=0,∂m2D∂zDzD=0,1=0Combining the initial condition, the mathematical model in the Laplace domain (Van Everdingen and Hurst 1949) [16] can be written as (22)1rD∂∂rDrD∂ξ-D10∂rD+LD2∂2ξ-D10∂zD2=hDLD2sξ-D100≤rD≤rmD1rD∂∂rDrD∂m-2D∂rD+LD2∂2m-2D∂zD2=M12ω12hDLD2sm-2DrmD≤rD≤∞limεD→0limrD→0⁡∫zwD-εD/2zwD+εD/2rD∂ξ-1D0∂rDdzwD=-12szD-zwD≤εD2ξ-1D0rD=rmD=m-2DrD=rmD∂ξ-1D0∂rDrD=rmD=1M12∂m-2D∂rDrD=rmDm-2DrD→∞=0∂ξ-1D0∂zDzD=0,1=0,∂m-2D∂zDzD=0,1=0Combining the closed top and bottom boundaries and the achievements of Wu (2008) [17] and Hu (2013) [18], the integral transformation kernel of Z direction is cos⁡nπzD.Defineξ--1D0 as the orthogonal transformation function of ξ-1D0, and define m--2D as the orthogonal transformation function of m-2D, then we can obtain (23) and (24)(23)ξ--1D0=∫01ξ-1D0rD,zDcos⁡nπzDdzD(24)m--2D=∫01m-2DrD,zDcos⁡nπzDdzDHence, the orthogonal transformation of the mathematical model (22) can be modified as (25)1rD∂∂rDrD∂ξ--1D0∂rD-hDLD2s+n2π2LD2ξ--1D0=00≤rD≤rmD1rD∂∂rDrD∂m--2D∂rD-hDLD2M12ω12s+n2π2LD2m--2D=0rmD≤rD≤∞limrD→0⁡rD∂ξ--1D0∂rD=-cos⁡nπzwD2sξ--1D0rD=rmD=m--2DrD=rmDξ--1D0∂rDrD=rmD=1M12∂m--2D∂rDrD=rmDm--2DrD→∞=0Letting λ1n=hDLD2s+n2π2LD2, λ2n=hDLD2M12/ω12s+n2π2LD2, the differential equations can be written as(26)∂2ξ--1D0∂rD2+1rD∂ξ--1D0∂rD-λ1n2ξ--1D0=00≤rD≤rmD(27)∂2m--2D∂rD2+1rD∂m--2D∂rD-λ2n2m--2D=0rmD≤rD≤∞Equations (26) and (27) are both Bessel functions, where the general solution of them can be expressed respectively as follows: (28)ξ--1D0=AnK0λ1nrD+BnI0λ1nrD(29)m--2D=CnK0λ2nrD+DnI0λ2nrDwhere I0 is the first zero order modified Bessel function, K0 is the second zero order modified Bessel function, and An, Bn, Cn, Dn are constants.Substituting (28) into inner boundary condition of (25),(30)limrD→0⁡rD∂ξ--1D0∂rD=limrD→0⁡rDλ1n-AnK1λ1nrD+BnI1λ1nrD=-cos⁡nπzwD2sAccording to the properties of Bessel function (Liu 2002; Duan 2018) [19, 20], ifz→0, then I1(z) = 0, K1(z) = 1/z, and we can obtain the value of An: (31)An=cos⁡nπzwD2sSubstituting (28) and (29) into interface condition of (25):(32)AnK0λ1nrmD+BnI0λ1nrmD=CnK0λ2nrmD+DnI0λ2nrmD(33)-Anλ1nK1λ1nrmD+Bnλ1nI1λ1nrmD=-CnM12λ2nK1λ2nrmD+DnM12λ2nI1λ2nrmDCombining infinite horizontal boundary condition of (25), (29), and I0∞=∞, we can obtain the value of Dn: (34)Dn=0Combining (31)–(34), we can obtain the value of Bn and Cn: (35)Bn=λ1n/λ2nK1λ1nrmDK0λ2nrmD-1/M12K0λ1nrmDK1λ2nrmDλ1n/λ2nI1λ1nrmDK0λ2nrmD+1/M12I0λ1nrmDK1λ2nrmD·An(36)Cn=AnK0λ1nrmD+BnI0λ1nrmDK0λ2nrmDSubstituting An, Bn into (28) and taking the inverse orthogonal transformation, we can achieve the pressure solution in the Laplace domain: (37)ξ-1D0=A0K0λ10rD+B0I0λ10rD+2∑n=1∞AnK0λ1nrD+BnI0λ1nrD·cos⁡nπzDwhere zD=zwD+rwDLDAccording to the source function idea method, the horizontal line source well can be achieved by integrating the horizontal section of the wellbore. Therefore, we can acquire the bottom pseudopressure solution for the horizontal well of composite gas reservoir:(38)
ξ
-
1
w
D
N
=
A
0
∫
-
1
1
K
0
λ
10
x
D
-
α
2
d
α
+
B
0
∫
-
1
1
I
0
λ
10
x
D
-
α
2
d
α
+
2
∑
n
=
1
∞
A
n
∫
-
1
1
K
0
λ
1
n
x
D
-
α
2
d
α
+
B
n
∫
-
1
1
I
0
λ
1
n
x
D
-
α
2
d
α
·
cos
⁡
n
π
z
DAccording to Duhamel’s theorem (Van Everdingen and Hurst) [21], take into consideration the wellbore storage coefficient and the skin factor, and the pseudopressure can be expressed as follows: (39)ξ-wD=sξ-1wDN+Ss+CDs2sξ-1wDN+SThrough numerical inversion of Stehfest (1970) [22], the bottom dimensionless pseudopressure solution for horizontal well with consideration of stress sensitivity, storage coefficient, and the skin factor in composite gas reservoir can be obtained:(40)m1D=-1γmDln⁡1-γmDξwDThe dimensionless rate equation in the Laplace space can be determined as follows: (41)q-Ds=1s2m-1DsCombining (38)–(41), we can acquire the dimensionless rate solution for horizontal well under the influence of the stress sensitivity, storage coefficient, and the skin factor for composite gas reservoirs.


## 2.1. Physical Model

A physical model of a composite gas reservoir is shown in Figure1. To make the proposed model more specific and easier to understand, the following assumptions and descriptions are made:Figure 1
Physical models for horizontal well of composite reservoir.① The top and bottom boundaries of the reservoir are closed, the horizontal boundary is infinite, and the thickness ish. The inner radius is rm, and the initial formation pressure is pi.② The vertical and horizontal permeability areKv and Kh, respectively. The horizontal length is 2L, the horizontal section is parallel to the vertical boundaries of the gas reservoir, and its location is unlimited and represented by zw.③ A slightly compressible gas has a constant viscosity and compressibility factor. The single-phase gas flow obeys a Darcy law, and the rock permeability has stress sensitivity and it occurs mainly in the near wellbore zone.④ There is no additional pressure drop in the inner and outer zone interfaces. The impact of gravity and the capillary force are ignored, but the wellbore storage and the skin factor are taken into consideration.


## 2.2. Mathematical Model

(1)
The relationship between reservoir permeability and pressure can be expressed as(1)K=Kie-γpi-p (2)
Motion equation:The equation of motion without consideration of permeability stress sensitivity is as follows:(2)vr=-3.6Khμ∂p∂rvz=-3.6Kvμ∂p∂zThe equation of motion with consideration of permeability stress sensitivity is as follows:(3)vr=-3.6Khμe-γpi-p∂p∂rvz=-3.6Kvμe-γpi-p∂p∂z(3)
State equation:(4)φ=φ01+Cφp-p0ρ=ρ01+Cρp-p0ρ=MpRTZ(4)
Continuity equation:(5)∇·ρv⃑+∂ρφ∂t=0We only consider the effect of stress sensitivity in the inner zone since the stress sensitivity mainly occurs in the near wellbore area. Combining (3)–(5) the differential equation of percolation that accounts for the stress sensitivity of inner zone can be expressed as follows:Inner zone0≤r≤rm:(6)Kh1r∂∂rrp1μ1Ze-γpi-p1∂p1∂r+Kv1∂∂zp1μ1Ze-γpi-p1∂p1∂z=φ1μ1Ct13.6p1μ1Z∂p1∂tCombining (2) and (4)-(5) the differential equation of percolation that accounts for the stress sensitivity of outer zone can be expressed as follows:Outer zonerm≤r<∞: (7)Kh2r∂∂rrp2μ2Z∂p2∂r+Kv2∂∂zp2μ2Z∂p2∂z=φ2μ2Ct23.6p2μ2Z∂p2∂tDefine the pseudopressure:(8)mp=∫0p2pμZdpDefine the permeability modulus with pseudopressure: (9)K=Kie-γmmi-mIntroducing the definition of dimensionless variables (Table 1), the dimensionless equations of the differential equation of percolation of composite gas reservoir with consideration of stress sensitivity effect can be expressed as follows:Inner zone0≤rD≤rmD: (10)∂2m1D∂rD2+1rD∂m1D∂rD-γmD∂m1D∂rD2+LD2∂2m1D∂zD2-γmD∂m1D∂zD2=hDLD2eγmDm1D∂m1D∂tDOuter zonermD≤rD≤∞: (11)1rD∂∂rDrD∂m2D∂rD+LD2∂2m2D∂zD2=M12ω12hDLD2∂m2D∂tDInitial condition:(12)m1DtD=0=m2DtD=0=0Inner boundary condition:(13)limεD→0⁡limrD→0⁡∫zwD-εD/2zwD+εD/2rDe-γmDm1D∂m1D∂rDdzwD=-12zD-zwD≤εD2Interface condition:(14)m1DrD=rmD=m2DrD=rmD(15)∂m1D∂rDrD=rmD=1M12∂m2D∂rDrD=rmDInfinite horizontal boundary condition:(16)m2DrD→∞=0Closed top and bottom boundaries:(17)∂m1D∂zDzD=0,1=0,∂m2D∂zDzD=0,1=0Equations (10)–(17) are the mathematical model of composite gas reservoir with consideration of stress sensitivity.Table 1
Definitions of the dimensionless variables.
Dimensionless variable Equation Dimensionless variable Equation Time t D = 3.6 K h 1 t φ 1 μ 1 C t 1 r w 2 Storage ratio ω 12 = φ 1 φ 2 Radial distance r D = r L Mobility ratio M 12 = λ 1 λ 2 = K 1 K 2 Vertical distance z D = z h Storage coefficient C D = 0.159 C φ 1 C t 1 h L 2 Horizontal length L D = L h K v K h Permeability modulus γ m D = T q s c 78.489 K h 1 h γ m Reservoir thickness h D = h r w K h K v Pseudo pressure of inner zone m 1 D = 78.489 K h 1 h T q s c m i - m 1 Horizontal section position z w D = z w h Pseudo pressure of outer zone m 2 D = 78.489 K h 1 h T q s c m i - m 2


## 2.3. Solutions of the Mathematical Model

Eq. (10), the differential equation of percolation with consideration of stress sensitivity, is a nonlinear equation and cannot be solved directly. Therefore, we use Pedrosa substitution and regular perturbation method (Guo 2009) [15] to eliminate the nonlinearity of the equation.Pedrosa substitution:(18)m1DrD,zD,tD=-1γmDln⁡1-γmDξ1DrD,zD,tDwhereξ1DrD,zD,tD is the perturbed deformation function.Hence, we can achieve the following expressions:(19)e-m1DγmD=1-γmDξ1D∂m1D∂rD=11-γmDξ1D∂ξ1D∂rD∂2m1D∂rD2=11-γmDξ1D∂2ξ1D∂rD2+γmD1-γmDξ1D2∂ξ1D∂rD2∂m1D∂zD=11-γmDξ1D∂ξ1D∂zD∂2m1D∂zD2=11-γmDξ1D∂2ξ1D∂zD2+γmD1-γmDξ1D2∂ξ1D∂zD2∂m1D∂tD=11-γmDξ1D∂ξ1D∂tDAccording to regular perturbation method, the following items of dimensionless permeability modulus can be expanded in power series:(20)ξ1D=ξ1D0+γmDξ1D1+γmD2ξ1D2+⋯11-γmDξ1D=1+γmDξ1D+γmD2ξ1D2+⋯-1γmDln1-γmDξ1D=ξ1D+12γmDξ1D2+⋯Since the dimensionless permeability modulus is usually very small (γmD≪1), the zeroth order perturbation solution can meet the requirements of the engineering precision. Combining (19)–(20), the mathematical model (10)–(17) can be modified as(21)1rD∂∂rDrD∂ξ1D0∂rD+LD2∂2ξ1D0∂zD2=hDLD2∂ξ1D0∂tD0≤rD≤rmD1rD∂∂rDrD∂m2D∂rD+LD2∂2m2D∂zD2=M12ω12hDLD2∂m2D∂tDrmD≤rD≤∞ξ1D0tD=0=m2DtD=0=0limεD→0limrD→0⁡∫zwD-εD/2zwD+εD/2rD∂ξ1D0∂rDdzwD=-12zD-zwD≤εD2ξ1D0rD=rmD=m2DrD=rmD∂ξ1D0∂rDrD=rmD=1M12∂m2D∂rDrD=rmDm2DrD→∞=0∂ξ1D0∂zDzD=0,1=0,∂m2D∂zDzD=0,1=0Combining the initial condition, the mathematical model in the Laplace domain (Van Everdingen and Hurst 1949) [16] can be written as (22)1rD∂∂rDrD∂ξ-D10∂rD+LD2∂2ξ-D10∂zD2=hDLD2sξ-D100≤rD≤rmD1rD∂∂rDrD∂m-2D∂rD+LD2∂2m-2D∂zD2=M12ω12hDLD2sm-2DrmD≤rD≤∞limεD→0limrD→0⁡∫zwD-εD/2zwD+εD/2rD∂ξ-1D0∂rDdzwD=-12szD-zwD≤εD2ξ-1D0rD=rmD=m-2DrD=rmD∂ξ-1D0∂rDrD=rmD=1M12∂m-2D∂rDrD=rmDm-2DrD→∞=0∂ξ-1D0∂zDzD=0,1=0,∂m-2D∂zDzD=0,1=0Combining the closed top and bottom boundaries and the achievements of Wu (2008) [17] and Hu (2013) [18], the integral transformation kernel of Z direction is cos⁡nπzD.Defineξ--1D0 as the orthogonal transformation function of ξ-1D0, and define m--2D as the orthogonal transformation function of m-2D, then we can obtain (23) and (24)(23)ξ--1D0=∫01ξ-1D0rD,zDcos⁡nπzDdzD(24)m--2D=∫01m-2DrD,zDcos⁡nπzDdzDHence, the orthogonal transformation of the mathematical model (22) can be modified as (25)1rD∂∂rDrD∂ξ--1D0∂rD-hDLD2s+n2π2LD2ξ--1D0=00≤rD≤rmD1rD∂∂rDrD∂m--2D∂rD-hDLD2M12ω12s+n2π2LD2m--2D=0rmD≤rD≤∞limrD→0⁡rD∂ξ--1D0∂rD=-cos⁡nπzwD2sξ--1D0rD=rmD=m--2DrD=rmDξ--1D0∂rDrD=rmD=1M12∂m--2D∂rDrD=rmDm--2DrD→∞=0Letting λ1n=hDLD2s+n2π2LD2, λ2n=hDLD2M12/ω12s+n2π2LD2, the differential equations can be written as(26)∂2ξ--1D0∂rD2+1rD∂ξ--1D0∂rD-λ1n2ξ--1D0=00≤rD≤rmD(27)∂2m--2D∂rD2+1rD∂m--2D∂rD-λ2n2m--2D=0rmD≤rD≤∞Equations (26) and (27) are both Bessel functions, where the general solution of them can be expressed respectively as follows: (28)ξ--1D0=AnK0λ1nrD+BnI0λ1nrD(29)m--2D=CnK0λ2nrD+DnI0λ2nrDwhere I0 is the first zero order modified Bessel function, K0 is the second zero order modified Bessel function, and An, Bn, Cn, Dn are constants.Substituting (28) into inner boundary condition of (25),(30)limrD→0⁡rD∂ξ--1D0∂rD=limrD→0⁡rDλ1n-AnK1λ1nrD+BnI1λ1nrD=-cos⁡nπzwD2sAccording to the properties of Bessel function (Liu 2002; Duan 2018) [19, 20], ifz→0, then I1(z) = 0, K1(z) = 1/z, and we can obtain the value of An: (31)An=cos⁡nπzwD2sSubstituting (28) and (29) into interface condition of (25):(32)AnK0λ1nrmD+BnI0λ1nrmD=CnK0λ2nrmD+DnI0λ2nrmD(33)-Anλ1nK1λ1nrmD+Bnλ1nI1λ1nrmD=-CnM12λ2nK1λ2nrmD+DnM12λ2nI1λ2nrmDCombining infinite horizontal boundary condition of (25), (29), and I0∞=∞, we can obtain the value of Dn: (34)Dn=0Combining (31)–(34), we can obtain the value of Bn and Cn: (35)Bn=λ1n/λ2nK1λ1nrmDK0λ2nrmD-1/M12K0λ1nrmDK1λ2nrmDλ1n/λ2nI1λ1nrmDK0λ2nrmD+1/M12I0λ1nrmDK1λ2nrmD·An(36)Cn=AnK0λ1nrmD+BnI0λ1nrmDK0λ2nrmDSubstituting An, Bn into (28) and taking the inverse orthogonal transformation, we can achieve the pressure solution in the Laplace domain: (37)ξ-1D0=A0K0λ10rD+B0I0λ10rD+2∑n=1∞AnK0λ1nrD+BnI0λ1nrD·cos⁡nπzDwhere zD=zwD+rwDLDAccording to the source function idea method, the horizontal line source well can be achieved by integrating the horizontal section of the wellbore. Therefore, we can acquire the bottom pseudopressure solution for the horizontal well of composite gas reservoir:(38)
ξ
-
1
w
D
N
=
A
0
∫
-
1
1
K
0
λ
10
x
D
-
α
2
d
α
+
B
0
∫
-
1
1
I
0
λ
10
x
D
-
α
2
d
α
+
2
∑
n
=
1
∞
A
n
∫
-
1
1
K
0
λ
1
n
x
D
-
α
2
d
α
+
B
n
∫
-
1
1
I
0
λ
1
n
x
D
-
α
2
d
α
·
cos
⁡
n
π
z
DAccording to Duhamel’s theorem (Van Everdingen and Hurst) [21], take into consideration the wellbore storage coefficient and the skin factor, and the pseudopressure can be expressed as follows: (39)ξ-wD=sξ-1wDN+Ss+CDs2sξ-1wDN+SThrough numerical inversion of Stehfest (1970) [22], the bottom dimensionless pseudopressure solution for horizontal well with consideration of stress sensitivity, storage coefficient, and the skin factor in composite gas reservoir can be obtained:(40)m1D=-1γmDln⁡1-γmDξwDThe dimensionless rate equation in the Laplace space can be determined as follows: (41)q-Ds=1s2m-1DsCombining (38)–(41), we can acquire the dimensionless rate solution for horizontal well under the influence of the stress sensitivity, storage coefficient, and the skin factor for composite gas reservoirs.


## 3. Type Curves and Discussions

As per (40), we plot the dimensionless pressure and dimensionless pressure derivative type curves under the influence of different parameters of horizontal well in stress-sensitive composite gas reservoir with tD/CD as X-axis and m1D, m1′D·tD/CD as Y-axis (Figures 2 and 3). As per (41), dimensionless rate decline and rate decline derivative type curves could be acquainted with tD/CD as X-axis and qD, q′D·tD/CD as Y-axis (Figures 4 and 5).Figure 2
Influence of permeability modulus on pressure type curve.Figure 3
Influence of different parameters on pressure type curve. (a) Influence of permeability modulus on pressure type curve (CD=10-5, LD=10, S=0.1, M12=1.6, rmD=15). (b) Influence of horizontal length on pressure type curve (CD=10-5, γmD=0.06, S=0.1, M12=1.6, rmD=15). (c) Influence of skin factor on pressure type curve (CD=10-5, γmD=0.06, LD=10, M12=1.6, rmD=15). (d) Influence of mobility ratio on pressure type curve (CD=10-5, γmD=0.06, S=0.1, LD=10, rmD=15). (e) Influence of inner zone radius on pressure type curve (K1>K2) (CD=10-5, γmD=0.06, S=0.1, LD=10, M12=2.5). (f) Influence of inner zone radius on pressure type curve (K2>K1) (CD=10-5, γmD=0.06, S=0.1, LD=10, M12=0.5).
(a) (b) (c) (d) (e) (f)Figure 4
Influence of permeability modulus on rate decline type curve.Figure 5
Influence of different parameters on rate decline type curve. (a) Influence of permeability modulus on rate decline type curve (CD=10-5, LD=10, S=0.1, M12=1.6, rmD=15). (b) Influence of horizontal length on rate decline type curve (CD=10-5, γmD=0.06, S=0.1, M12=1.6, rmD=15). (c) Influence of skin factor on rate decline type curve (CD=10-5, γmD=0.06, LD=10, M12=1.6, rmD=15). (d) Influence of mobility ratio on rate decline type curve (CD=10-5, γmD=0.06, S=0.1, LD=10, rmD=15). (e) Influence of inner zone radius on rate decline type curve (K1>K2) (CD=10-5, γmD=0.06, S=0.1, LD=10, M12=2.5). (f) Influence of inner zone radius on rate decline type curve (K2>K1) (CD=10-5, γmD=0.06, S=0.1, LD=10, M12=0.5).
(a) (b) (c) (d) (e) (f)


### 3.1. Dimensionless Transient Pressure Analysis

Figure2 presents the influence of permeability modulus on dimensionless pressure (solid line) and dimensionless pressure derivative (dotted line) type curves. The type curve can be used to analyze the characteristics of the flow behavior and determine the property parameters of the formation. As demonstrated in Figure 2, there are seven flow stages of the dimensionless type curve:I.
Pure wellbore storage flow stage: in this stage, the curves of the pressure and pressure derivative coincide with a straight line with the slope equal to 1.II.
Transition flow stage: this stage is influenced by both wellbore storage coefficient and skin factor, pressure curve increases slowly after stage 1, and the pressure derivative curve rises to the maximum value and then gradually decreases.III.
Early radial flow stage: during this stage, the pressure derivative curve turns into a horizontal line, and the value is relative to the horizontal length. This stage shows the feature of radial flow which is perpendicular to the horizontal axis before pressure wave reaches to the top and bottom boundary.IV.
Mid-term linear flow stage: the pressure derivative curve is a straight line and the slope is 1/2. In this stage, the pressure wave has already reached the top and bottom boundary.V.
Pseudoradial flow stage of inner zone: at the point whenγmD=0, the pressure derivative curve (black dotted line) is a horizontal line with a value of 0.5.VI.
Transition flow stage from inner zone pseudoradial flow to outer zone pseudoradial flow.VII.
Pseudoradial flow stage of outer zone: at the point whenγmD=0, the pressure derivative curve (black dotted line) is a horizontal line, and its value is related to the mobility ratio M12. This stage demonstrates the characteristics of the pseudoradial flow of the whole system.Figure3 shows the effect of different parameters on pressure type curve. As demonstrated by Figure 3(a), the influence of permeability modulus on pressure curve is for the most part mirrored in stage V to stage VII. The higher the permeability modulus, the higher position of dimensionless pressure and derivative curve, and the greater pressure drop of a horizontal well.Figure3(b) shows the effect of horizontal length on pressure type curve, and it is mainly reflected in stage III. The longer the horizontal length, the shorter the duration of stage III. This is because the longer horizontal length means the slower the pressure drop speed of the gas reservoir. With the increase of horizontal length, the position of pressure type curve decreased, and the values of the horizontal line in pressure derivative curve (the dotted lines) are equal to “1/(4LD)”. However, the decrease of the pressure type curve is getting smaller with the increase of the length of the horizontal; therefore, the horizontal length is not as longer as better.As demonstrated by Figure3(c), for the most part, the influence of skin factor on pressure type curve is reflected in stage II. With the increase of skin factor, the position of the pressure type curve increased, and the “Hump” value of pressure derivative curve (the dotted lines) also increased. This is because the greater the skin factor, the more serious the pollution in the near area in wellbore, and the greater additional resistance of the gas flow.Figure3(d) shows the effect of mobility ratio on pressure type curve, and it is mainly reflected in stage VI to stage VII. When the permeability of the inner zone is superior to the outer zone (M12>1), the position of the second step is higher than the first step, and the greater the mobility ratio, the higher the position of the second step (the red and the blue dotted lines). Likewise, when the permeability of the outer zone is superior to the inner zone (M12<1), the position of the second step is lower than the first step (the black dotted line). This is the most obvious feature of the flow behavior in composite reservoir. Overall, with the increase of mobility ratio, the position of the pressure type curve also increased. The reason is that the greater mobility ratio means that the permeability of the outer zone is smaller as compared to the inner zone, and the worse is the flow ability of the fluid in stage VII.Figures3(e) and 3(f) show the effect of inner zone radius on pressure type curve; the influence of inner zone radius on pressure curve is for the most part mirrored in stage V to stage VII. At the point when the permeability of the inner zone is superior to the outer zone, the position of pressure type curve decreased with the increase of inner zone radius (Figure 3(e)). At the point when the permeability of the outer zone is superior to the inner zone, the position of pressure type curve increased with the increase of inner zone radius (Figure 3(f)). The greater the inner zone radius, the longer the duration of stage V (the dotted lines), the later the appearance of stage VII. This is because the greater the inner zone radius is, the longer the time it takes for the pressure wave to spread to the inner boundary. Therefore, it can be used to determine the location of the boundary between inner zone and outer zone.


### 3.2. Dimensionless Rate Decline Type Curve Analysis

As demonstrated in Figure4, the rate decline (solid line) and rate decline derivative (dotted line) type curves of composite reservoir are also composed of seven flow stages, and the flow stages are similar to the pressure type curve. The rate decline type curve can be used to estimate reserve evaluation and predict gas well production.Figure5 shows the influence of different parameters on rate decline type curve. As demonstrated by Figure 5(a), the influence of permeability modulus on rate decline curve is for the most part mirrored in stage V to stage VII (the solid lines). Overall, the position of the rate decline curve gradually falls down, with the increase of permeability modulus, demonstrating that the production rate is lower.As demonstrated in Figure5(b), the effect of horizontal length on rate decline type curve is mainly reflected in Stage III. The production rate is higher because of the increment in the horizontal length. The reason is that the longer horizontal length means the larger gas drainage area. However, it is not the longest horizontal length at the higher production rate, because the longer the horizontal length is, the greater it will bring friction resistance.As demonstrated in Figure5(c), for the most part, the influence of skin factor on rate decline type is reflected in stage II. The position of rate decline curve falls down, and the production rate is lower because of the increment of the skin factor. The reason is that the greater skin factor means more serious pollution in the vicinity of the wellbore. After stage II, the effect of skin factor on rate decline will gradually decrease and finally disappear.Figure5(d) shows the influence of mobility ratio on rate decline type curve, and it is mainly reflected in stage VI to stage VII. The position of the rate decline curve falls down and the production rate is lower because of the increase of the mobility ratio. The reason is that the greater the mobility ratio, the worse the flow ability of the fluid in stage VII, and the lowest production rate in the late period.Figures5(e) and 5(f) show the effect of inner zone radius on rate decline type curve; the influence of inner zone radius on rate decline curve is for the most part mirrored in stage V to stage VII. When the permeability of the inner zone is superior to the outer zone, the production rate increases continuously with the expansion of the inner zone radius (Figure 5(e)). When the permeability of the outer zone is superior to the inner zone, the production rate diminishes continuously with the expansion of the inner zone radius (Figure 5(f)).


## 3.1. Dimensionless Transient Pressure Analysis

Figure2 presents the influence of permeability modulus on dimensionless pressure (solid line) and dimensionless pressure derivative (dotted line) type curves. The type curve can be used to analyze the characteristics of the flow behavior and determine the property parameters of the formation. As demonstrated in Figure 2, there are seven flow stages of the dimensionless type curve:I.
Pure wellbore storage flow stage: in this stage, the curves of the pressure and pressure derivative coincide with a straight line with the slope equal to 1.II.
Transition flow stage: this stage is influenced by both wellbore storage coefficient and skin factor, pressure curve increases slowly after stage 1, and the pressure derivative curve rises to the maximum value and then gradually decreases.III.
Early radial flow stage: during this stage, the pressure derivative curve turns into a horizontal line, and the value is relative to the horizontal length. This stage shows the feature of radial flow which is perpendicular to the horizontal axis before pressure wave reaches to the top and bottom boundary.IV.
Mid-term linear flow stage: the pressure derivative curve is a straight line and the slope is 1/2. In this stage, the pressure wave has already reached the top and bottom boundary.V.
Pseudoradial flow stage of inner zone: at the point whenγmD=0, the pressure derivative curve (black dotted line) is a horizontal line with a value of 0.5.VI.
Transition flow stage from inner zone pseudoradial flow to outer zone pseudoradial flow.VII.
Pseudoradial flow stage of outer zone: at the point whenγmD=0, the pressure derivative curve (black dotted line) is a horizontal line, and its value is related to the mobility ratio M12. This stage demonstrates the characteristics of the pseudoradial flow of the whole system.Figure3 shows the effect of different parameters on pressure type curve. As demonstrated by Figure 3(a), the influence of permeability modulus on pressure curve is for the most part mirrored in stage V to stage VII. The higher the permeability modulus, the higher position of dimensionless pressure and derivative curve, and the greater pressure drop of a horizontal well.Figure3(b) shows the effect of horizontal length on pressure type curve, and it is mainly reflected in stage III. The longer the horizontal length, the shorter the duration of stage III. This is because the longer horizontal length means the slower the pressure drop speed of the gas reservoir. With the increase of horizontal length, the position of pressure type curve decreased, and the values of the horizontal line in pressure derivative curve (the dotted lines) are equal to “1/(4LD)”. However, the decrease of the pressure type curve is getting smaller with the increase of the length of the horizontal; therefore, the horizontal length is not as longer as better.As demonstrated by Figure3(c), for the most part, the influence of skin factor on pressure type curve is reflected in stage II. With the increase of skin factor, the position of the pressure type curve increased, and the “Hump” value of pressure derivative curve (the dotted lines) also increased. This is because the greater the skin factor, the more serious the pollution in the near area in wellbore, and the greater additional resistance of the gas flow.Figure3(d) shows the effect of mobility ratio on pressure type curve, and it is mainly reflected in stage VI to stage VII. When the permeability of the inner zone is superior to the outer zone (M12>1), the position of the second step is higher than the first step, and the greater the mobility ratio, the higher the position of the second step (the red and the blue dotted lines). Likewise, when the permeability of the outer zone is superior to the inner zone (M12<1), the position of the second step is lower than the first step (the black dotted line). This is the most obvious feature of the flow behavior in composite reservoir. Overall, with the increase of mobility ratio, the position of the pressure type curve also increased. The reason is that the greater mobility ratio means that the permeability of the outer zone is smaller as compared to the inner zone, and the worse is the flow ability of the fluid in stage VII.Figures3(e) and 3(f) show the effect of inner zone radius on pressure type curve; the influence of inner zone radius on pressure curve is for the most part mirrored in stage V to stage VII. At the point when the permeability of the inner zone is superior to the outer zone, the position of pressure type curve decreased with the increase of inner zone radius (Figure 3(e)). At the point when the permeability of the outer zone is superior to the inner zone, the position of pressure type curve increased with the increase of inner zone radius (Figure 3(f)). The greater the inner zone radius, the longer the duration of stage V (the dotted lines), the later the appearance of stage VII. This is because the greater the inner zone radius is, the longer the time it takes for the pressure wave to spread to the inner boundary. Therefore, it can be used to determine the location of the boundary between inner zone and outer zone.


## 3.2. Dimensionless Rate Decline Type Curve Analysis

As demonstrated in Figure4, the rate decline (solid line) and rate decline derivative (dotted line) type curves of composite reservoir are also composed of seven flow stages, and the flow stages are similar to the pressure type curve. The rate decline type curve can be used to estimate reserve evaluation and predict gas well production.Figure5 shows the influence of different parameters on rate decline type curve. As demonstrated by Figure 5(a), the influence of permeability modulus on rate decline curve is for the most part mirrored in stage V to stage VII (the solid lines). Overall, the position of the rate decline curve gradually falls down, with the increase of permeability modulus, demonstrating that the production rate is lower.As demonstrated in Figure5(b), the effect of horizontal length on rate decline type curve is mainly reflected in Stage III. The production rate is higher because of the increment in the horizontal length. The reason is that the longer horizontal length means the larger gas drainage area. However, it is not the longest horizontal length at the higher production rate, because the longer the horizontal length is, the greater it will bring friction resistance.As demonstrated in Figure5(c), for the most part, the influence of skin factor on rate decline type is reflected in stage II. The position of rate decline curve falls down, and the production rate is lower because of the increment of the skin factor. The reason is that the greater skin factor means more serious pollution in the vicinity of the wellbore. After stage II, the effect of skin factor on rate decline will gradually decrease and finally disappear.Figure5(d) shows the influence of mobility ratio on rate decline type curve, and it is mainly reflected in stage VI to stage VII. The position of the rate decline curve falls down and the production rate is lower because of the increase of the mobility ratio. The reason is that the greater the mobility ratio, the worse the flow ability of the fluid in stage VII, and the lowest production rate in the late period.Figures5(e) and 5(f) show the effect of inner zone radius on rate decline type curve; the influence of inner zone radius on rate decline curve is for the most part mirrored in stage V to stage VII. When the permeability of the inner zone is superior to the outer zone, the production rate increases continuously with the expansion of the inner zone radius (Figure 5(e)). When the permeability of the outer zone is superior to the inner zone, the production rate diminishes continuously with the expansion of the inner zone radius (Figure 5(f)).


## 4. Model Application with Measured Data

A production well X1 is located in a low-permeability composite gas reservoir. The relevant parameters are as follows: the reservoir thickness is 18.2 m, the average porosity is 1.51%, the well radius is 0.081 m, the horizontal length is 357 m, the reservoir temperature is 95.2°C, and the gas viscosity is 0.021 mPa·s. This well has steadily produced for 660 hours at a rate of 32.79×104m3/d before the buildup test. As demonstrated in Figure 6, the dimensionless pressure and pressure derivative curve cannot be fitted very well by using composite reservoir model without considering stress sensitivity. The proposed model in this paper can fit the actual well test data very well and can use for well test analysis of horizontal well in stress-sensitive composite reservoir. The fitted results are shown in Table 2.Table 2
Fitted results of X1.
Parameter Value Wellbore storage coefficientC (m3/MPa) 4.26 Skin factorS −1.57 Vertical permeability of inner zoneKv1(μm2) 0.83×10−3 Horizontal permeability of inner zoneKh1(μm2) 1.16×10−3 Storage ratioω12 3.72 Mobility ratioM12 1.5 Inner zone radiusrm (m) 282.4 Initial reservoir pressurepi (MPa) 23.8 Dimensionless permeability modulusγmD 0.11Figure 6
Log-log fitted curves of actual measured data.


## 5. Conclusions

This paper proposed a new analytical solution to analyze pressure transient behavior and transient rate decline of horizontal well in stress-sensitive composite gas reservoir. The example analysis verifies that the proposed model is reliable and practical. The accompanying conclusions can be summarized as follows:(1) For the most part, the permeability modulus influences stage V to stage VII of type curves. The position of pressure type curve increased with increasing permeability modulus, and the production rate decreased with increasing permeability modulus.(2) The horizontal length is mainly reflected in Stage III of type curves. The pressure drop speed decreased and the production rate increased with increasing horizontal length, but the horizontal length is not as longer as better.(3) The skin factor is mainly reflected in stage II of type curves. The mobility ratio is mainly reflected in stage VI to stage VII of type curves. The position of pressure type curve increased with increasing skin factor and mobility ratio. The production rate decreased with increasing skin factor and mobility ratio. The greater skin factor means more serious pollution in the near area in wellbore. The pressure derivative type curves under the influence of mobility ratio will appear “step,” and it is the most obvious feature of the flow behavior in composite reservoir.(4) The inner zone radius is mainly mirrored in stage V to stage VII of type curves. When the permeability of the inner zone is superior to the outer zone, the position of pressure type curve decreased and the production rate increased with increasing inner zone radius; otherwise, the result is reversed. This parameter can be used to determine the location of the boundary between inner zone and outer zone.



---

*Source: 8672578-2018-10-31.xml*