Pressure Transient Analysis for a Horizontal Well in Heterogeneous Carbonate Reservoirs Using a Linear Composite Model

State Key Laboratory of Oil and Gas Reservoir Geology and Exploitation, Southwest Petroleum University, Chengdu 610500, Sichuan, China Department of Petroleum Engineering, College of Engineering and Applied Science, University of Wyoming, 1000 East University Avenue, Laramie, WY 82071-2000, USA Engineering Technology Research Institute, Petro China Southwest Oil & Gasfield Company, Chengdu 610031, Sichuan, China


Introduction
Carbonate reservoirs play a very important role in the world's gas and oil production. Most of the world's giant fields produce from naturally fractured and vuggy carbonate reservoirs that have complex pore structure. Drilling horizontal wells would be challenging because of unfavorable geomechanical properties [1,2]. e majority of carbonate reservoirs are very heterogeneous consisting of different pore classes as well as being fractured. One of the most prominent problems of carbonate reservoirs is the strong heterogeneity and multiporosity network characteristic, which makes the fluid flow system complicated and causes relatively low oil recovery. Many mature theories and practices of clastic rocks in reservoir development cannot be applied to carbonate reservoirs. It is generally considered that carbonate reservoirs have their own characteristics, that is, well-developed connected fracture systems, most of which also contain tiny pore systems connected with fracture systems. Unlike the conventional single-porosity sandstone reservoirs, these kinds of reservoirs need to consider the permeability and porosity of the rock mass as well as the exchange of quality between the fractures and matrix in addition to considering the permeability and porosity of the fracture network as the main seepage channel. Different from the common conceptual dual-porosity coal reservoirs where gas is mostly stored in the coal matrix and Darcy fluid flow occurs in the natural fracture system [3,4] and shale gas reservoirs where gas is trapped in poorly connected micropores by absorption on and in particles of organic materials and clay minerals in the matrix of the host rock [5], fractures have high permeability and low porosity, while matrix has low permeability and high porosity. e low flow in the matrix mainly comes from the fluid exchange between matrix and fracture network. e presence of fractures and vugs significantly influences the dynamic of fluid flow in carbonate reservoirs and the transient pressure analysis of wells drilled in these reservoirs.
Barenblatt et al. [6] pointed out that it is not promising to study the nonsteady-state seepage in fractured rocks by assuming a system of fractures, and therefore, they proposed the concept of dual porosity for fluid flow in naturally fractured reservoirs for the first time. ey mentioned that the interconnected fractures and matrix rock in fractured reservoirs as two continuous media are spatially superimposed on each other. e permeability of matrix rock is far less than that of the fracture system, but the porosity of matrix rock is far higher than that of the fracture system. Based on this hypothesis, Warren and Root [7] developed an idealized model for studying the characteristic behavior of a porous medium, which contains regions that contribute significantly to the pore volume of the system but contribute negligibly to the flow capacity. Motivated by the early model by Warren and Root, the unsteady flow in dual-porosity media was mathematically modeled and solved analytically [8,9]. en, with the development of modern well testing methods, the effect of matrix-block-size distribution on pressure transient response for naturally fractured reservoirs was discussed [10][11][12][13]. Considering both pseudosteady state and unsteady-state fluid transfer from matrix blocks into fracture network, Raghavan and Ohaeri [14] delineated the characteristics of a well flowing at a constant pressure in a naturally fractured reservoir.
However, most of the fractured carbonate reservoirs in the world have more complex pore structure. Vugs contribute to both effective porosity and permeability. Matrix and vugs do not have the same interaction with the fracture network [15]. Abdassah and Ershaghi [16] first proposed a triple-porosity/single-permeability model, which considers an unsteady-state interporosity flow between the fractures and two types of matrix blocks and a primary flow only through the fracture network. Jalali and Ershaghi [17] analyzed the unsteady pressure behavior of such reservoirs. New dual-fracture models [18] considering two types of flow units were proposed to differentiate between the microfractures and the macrofractures. en, a triple porosity model [15] was built to address high secondary porosity and mainly vuggy porosity in naturally fractured reservoirs.
Horizontal and multifractured horizontal wells can improve the flow conductivity and enhance well production and recovery factor by increasing the exposure of wells to formation [19,20]. Traditional pressure transient models assume that the entire length of a horizontal or multilateral well remains in the same formation with uniform properties. In reality, most horizontal and multilateral wells extend through the reservoirs consisting of layers, sections, or pockets with different properties [21]. Bowman and Crawford [22] presented a sufficient number of tables to calculate the pressure transient in a linear semi-infinite water-driven reservoir, where the properties in the oil and water zones were different. Bixel et al. [23] proposed a detailed treatment of the pressure transient behavior for a well located near a linear discontinuity. Strltsova and Mckinley [24] discussed the buildup pattern for linear discontinuities and the effect of flow time on quantitative assessment of reservoir parameters. Ambastha et al. [25] employed analytical solutions to obtain pressure-transient behavior for a line-source well in a composite reservoir. Based on Ambastha's model, an analytical solution for a three-zone laterally composite reservoir model [26] was found in the Laplace-Fourier space, and Kuchuk and Tarek [27] extended this model to multizone case, explaining the changes in porosity, permeability, and compressibility. Based on the material balance equation and considering the nonlinearity caused by gas desorption and diffusion, and the complex flow of MFHW, Wei et al. [28,29] developed a numerical method for multifractured horizontal wells in shale gas reservoirs. Medeiros et al. [21] presented a semianalytical model for the pressure-transient analysis of horizontal wells in linear composite reservoirs and provided the detailed solutions for a two-zone, laterally composite reservoir. If the reservoir is divided into more than two sections, their solution methodology will not work. e reason is that the coupling of middle sections of multizone reservoir is much more complicated than two-zone reservoir. Dejam et al. [30][31][32] used semianalytical solutions for pressure transient analysis of a hydraulically fractured vertical well and multifractured horizontal well in a bounded dual-porosity reservoir. As carbonate reservoirs usually have strong anisotropy, a horizontal well completed in these reservoirs may extend through multiple zones, including homogeneous, dual-porosity, and triple-porosity formations. Traditional well test models are inapplicable here. erefore, it is necessary to build new models to satisfy the need for well testing in a horizontal well in heterogeneous carbonate reservoirs.
is study is developed based on Medeiros et al.'s model and works for a horizontal well extending through a multizone reservoir. Considering different combinations of natural fractures, rock matrix, and vugs in carbonate reservoirs, various physical and mathematical models are established and solved by using the point source function, Green's function, and coupling of multiple reservoir sections. e corresponding type curves are obtained by computer programming, and sensitivity analysis are carried out.

Solution Methodology
First, in order to solve the problem under study, where a horizontal well traverses two or more regions with distinct properties, Green's function formulation of the pressure transient solution for a locally homogeneous reservoir section with arbitrary flux conditions at the boundaries [21] is applied. en, the coupling of different sections of the reservoir is explained.

Source Functions for the Inner and Outer Boundaries.
Green's function for a rectangular parallelepiped satisfies the adjoint diffusion equation with its flux vanishing at the boundaries of the domain.
is Green's function corresponds to the point source solution in the Laplace transform domain for a rectangular parallelepiped one in Figure 1 [33]. Where where D is the domain consisting of porous medium, eD is the outer boundary of the domain, l is the characteristic reference length in the system, K is the uniform permeability, K z is the principal permeability in the z direction, and μ is the viscosity of a slightly compressible fluid of constant compressibility.
In equation (1) and equations (5)-(7), the new parameter u was introduced to extend the solution into the dual-porosity or triple-porosity idealization of naturally fractured reservoirs [33,34] and is given by u � s, for homogeneous reservoirs, sf(s), for naturally fractured reservoirs.
Appropriate expressions for f(s) are related to the type of formation [34,35].
is paper mainly considers the conventional dual-porosity and the triple-porosity models, where the seepage type is triporate-uniphase parallel flow.
e corresponding expressions for f(s) are given, respectively, as follows: where w is the dimensionless storativity and λ is the dimensionless transfer coefficient. e subscript N f denotes the natural fracture system, v refers to the vug, m represents the matrix, and s stands for the Laplace variable. Based on Medeiros and Ozkan's model, the inner boundary, B w , is a line equal to the length of the horizontal well, L h . en, the source function in the Laplace domain, S α , given by Medeiros and Ozkan can be expressed as (11) where x wi , y wi , and z wi are the coordinates of the midpoint and L hi is the length of the i th segment of the horizontal well. Substituting G given by equation (1) into equation (11) and evaluating the integral yield, ch ε j,n y D1i + ch ε j,n y D2i where e outer boundary consists of the planar surfaces of the rectangular parallelepiped representing the domain. In this study, the source function for a segment on the boundary perpendicular to the x-axis (a planar source in the y-z plane) is given by S ej � y ej +y pj /2 where y ej and z ej are the coordinates of the midpoint and y pj and z pj are the lengths of the sides of the planar segment.
Substituting equation (1) into equation (14) and evaluating the internal result in S ej M D , s � Cz pj y ej +y pj /2 cos kπx D x eD cos kπx wD x eD y ej +y pj /2 cos kπx D x eD cos kπx wD x eD y ej +y pj /2 ch ε k,n y D1j + ch ε k,n y D2j ε k,n sh ε k,n y eD dy ′ .

Coupling of Multiple Reservoir Sections.
Medeiros et al. [21] developed a semianalytical model for the pressuretransient analysis of horizontal wells in linear composite reservoirs and provided the detailed solution for a two-zone reservoir. In their study, all section boundaries except the interface between sections are assumed impermeable. Physically, this case corresponds to a horizontal well that penetrates a closed and rectangular reservoir with two homogeneous sections of different properties. However, their model cannot provide solution if the horizontal well traverses more than two sections. e reason is that the coupling of middle sections of a multizone reservoir is much more complicated than a two-zone reservoir.
To demonstrate the coupling of multiple reservoir sections by using the solution given in equation (1), two rectangular reservoir sections in series in the x-direction are considered as the basis for the derivation. en, the models for three reservoir sections and even multiple reservoir sections models are developed. For the sake of simplicity, the discretization scheme shown in Figures 2 and 3, which includes three horizontal well segments and two interface segments between the sections, is used.
From equation (1), the pressure drop at the center of a well or interface segment i in section k is given by where M Di indicates the midpoint of segment i, q wj and q el are the fluxes on all segments of the inner and outer boundaries, respectively. Using equation (16) at the center of all well and interface segments shown in Figure 3 yields the following set of 14 linear equations with 28 unknowns: for i, k � 1, 2. To equate the number of equations with the number of unknowns, the condition of continuity of pressure and flux at the interface between sections 1, 2, and 3 is implemented. e continuity of pressure and flux requires is reduces the number of unknowns by twelve. If it is also assumed that there is infinite conductivity in the wellbore, then five additional equations are obtained: Finally, the mass balance requires that the sum of fluxes entering the wellbore is equal to the production rate at the sandface q: e linear system defined by equation (17)-(25) now has 16 equations with 16 unknowns.
For models addressing multiple reservoir sections, the discretization scheme shown in Figures 4 and 5, which includes n horizontal well segments and n − 1 interface segments between the sections, is used.
q w1~n q w2 Figure 5: Discretization scheme to demonstrate the coupling of n reservoir sections penetrated by a horizontal well. 6 Mathematical Problems in Engineering is reduces the number of unknowns by 4n. If it is also assumed that there is infinite conductivity in the wellbore, then 2n − 1 additional equations are obtained: Finally, mass balance requires that the sum of fluxes entering the wellbore is equal to the production rate at the sandface q: e linear system defined by equations (26)-(34) now has 2(3n − 1) equations with 2(3n − 1) unknowns.

Solution of the System of Linear Equations.
e system of linear equations, given by equations (16)
(40) e solution of this system of linear equations can be accomplished by using any matrix inverter. It is, however, in the Laplace transform domain and needs to be numerically inverted into the real-time domain. As it is common in the petroleum engineering literature, the Stehfest algorithm [36] is used for the numerical inversion of Laplace transforms in this study.

Homogeneous-Dual-Porosity (Fracture-Matrix)-Homogeneous Model.
It is assumed that a horizontal well traverses a reservoir consisting of homogeneous, dual-porosity (fracture-matrix), and homogeneous sections with the same permeability. e physical model is illustrated in Figure 6.
In Figure 7, the solid curves represent the type curves of flow stage division for a horizontal well traversing a three-section reservoir including homogeneous, dualporosity (fracture-matrix), and homogeneous sections. e triangles represent the type curves for a dual-porosity infinite-conductivity model. e stages of pseudosteady interporosity flow from matrix into fracture, fracture pseudoradial flow, system linear flow, system pseudoradial flow, and pseudosteady flow occur in sequence. In the early stage of pseudosteady interporosity flow from matrix into fracture, a valley appears on the pressure derivative curve. Due to the effect of the homogeneous sections, the valley of the three-section model is shallower. en, the fracture pseudoradial flow happens with a constant pressure derivative value. In the stage of system linear flow, the derivative curve exhibits a straight line with a slope of 0.5, and thereafter, the system pseudoradial flow occurs. Finally, there is a pseudosteady flow stage with a slope of 1 on the pressure derivative curve.
It can be seen from Figure 8 that, in the early stage, due to the interporosity flow from matrix into fracture, the flow rate of the two segments of wellbore in the dual-porosity section is larger than the flow rate of the two segments of wellbore in the homogeneous sections. With the advent of pseudoradial flow regime, the rates for all six segments of the wellbore become equal. At late times, the flux distribution tends to be stable and the segment closer to the heel of the horizontal well gets higher flux. Based on this model, keeping the other parameters unchanged, when the permeability of two homogeneous sections decreases to 10 mD, the pressure drop and pressure derivative curves become different, as shown in Figure 9.

Mathematical Problems in Engineering
In Figure 9, the solid curves represent the original type curves of flow stage division for a horizontal well traversing a three-section reservoir including homogeneous, dual-porosity (fracture-matrix), and homogeneous sections. e triangles represent the type curves for the model with the decreased permeability for the two homogeneous sections. When the permeability of two homogeneous sections decreases from 100 mD to 10 mD, the values of pressure drop and pressure derivative increase. In addition, the stage of the interporosity flow from matrix into fracture becomes easy to recognize in the early period. e valley becomes deeper. Other flow regimes remain unchanged. In Figure 10, the flow rate of the dual-porosity section is larger than that of the two homogeneous sections because the permeability of section 2 is higher than that of sections 1 and 3. In the early stage, due to the influence of interporosity flow, the wellbore flow in the second section gradually decreases until the stage of the fracture pseudoradial flow appears. With the increase of time, the rates for all six segments of the wellbore reach a steady value. Based on the model, when the permeability of sections 1 and 2 keeps unchanged and the permeability of the homogeneous reservoir in section 3 increases from 10 mD to 40 mD and then to 100 mD, the wellbore pressure drops and pressure derivative values decreases as shown in Figure 11. When the permeability of homogeneous section 3 becomes closer to that of section 2 (dual-porosity), the dip of the interporosity flow stage becomes shallower. Figure 12 shows the flux distribution along the horizontal well when the permeability of section 3 is 40 mD. e  x e = 400m, y e = 400m, h = 40m,   x e = 400m, y e = 400m, h = 40m, Figure 9: Type curves of flow stage division for a horizontal well traversing a three-section reservoir including homogeneous, dualporosity (fracture-matrix), and homogeneous sections (solid curves). e triangles represent the type curves for the model with the decreased permeability for the two homogeneous sections. flow rate of the two segments of wellbore in section 2 with the highest permeability is greater than that of the two segments of wellbore in other two sections.

Homogeneous-Triple-Porosity (Fracture-Vug-Matrix)-Homogeneous Model.
It is assumed that a horizontal well traverses a reservoir consisting of homogeneous, tripleporosity (fracture-vug-matrix), and homogeneous sections with the same permeability. e physical model is illustrated in Figure 13.
In Figure 14, the solid curves represent the type curves of flow stage division for a horizontal well traversing a threesection reservoir including homogeneous, triple-porosity (fracture-vug-matrix), and homogeneous sections. e triangles represent the type curves for a triple-porosity infinite conductivity model. e stages of pseudosteady interporosity flows from vug and matrix into fracture, fracture pseudoradial flow, system linear flow, system pseudoradial flow, and pseudosteady flow occur sequentially. In the early stage of pseudosteady interporosity flow from vug and matrix into fracture, two valleys appear on the pressure derivative curve. Due to the effect of the homogeneous sections, the valleys of the three-section model are shallower. en, the fracture pseudoradial flow happens with a constant pressure derivative value. In the stage of system linear flow, the derivative curve exhibits a straight line with a slope of 0.5, and thereafter, the system pseudoradial flow occurs. Finally, there is a pseudosteady flow stage with a slope of 1 on the pressure derivative curve.
It can be observed from the flux distribution in Figure 15 that due to the two interporosity flows of the triple-porosity section in the early stage, the flow rate of the two segments of wellbore in the triple-porosity section is larger than the flow rate of the two segments of wellbore in the homogeneous sections. With the appearance of pseudoradial flow regime, the rates for all six segments of the wellbore become equal. At late times, the flux distribution tends to be stable, and the segment closer to the heel of the horizontal well gets higher flux. Based on this model, keeping the other parameters unchanged, when the permeability of two homogeneous sections decrease to 10 mD, the pressure drops and pressure derivative curves become different, as shown in Figure 16.
In Figure 16, the solid curves represent the original type curves of flow stage division for a horizontal well traversing a three-section reservoir including homogeneous, triple-porosity (fracture-vug-matrix), and homogeneous sections. e triangles represent the type curves for the model with the decreased permeability for the two homogeneous sections. When the permeability of homogeneous sections decreases from 100 mD to 10 mD, the values of pressure drop and pressure derivative increase and the stage of the interporosity flow from both vug and matrix into fracture becomes easier to recognize in the early period. e dips become deeper. Other flow regimes remain unchanged.
It can be seen from Figure 17 that the flow rate of the two segments of wellbore in section 2 is higher than that of the two segments of wellbore in other two sections.
at is because the permeability of the triple-porosity section is   higher than that of the first and third sections. In the early stage, due to the influence of interporosity flows, the wellbore flow in the second section gradually decreases until the stage of the fracture pseudoradial flow appears. With the increase of time, the rates for all six segments of the wellbore reach a steady value.

Homogeneous-Dual-Porosity (Fracture-Matrix)-Triple-
It is assumed that a horizontal well traverses a reservoir consisting of homogeneous, dual-porosity (fracture-matrix), and triple-porosity (fracture-vug-matrix) sections with the same permeability. e physical model is illustrated in Figure 18.
In Figure 19, the solid curves represent a horizontal well traversing a three-section reservoir including homogeneous, dual-porosity (fracture-matrix), and triple-porosity (fracture-vug-matrix) sections. e triangles represent the type curves for a horizontal well traversing a three-section reservoir including homogeneous, triple-porosity (fracturevug-matrix), and homogeneous sections. When a horizontal well traverses a reservoir consisting of homogeneous, dualporosity (fracture-matrix), and triple-porosity (fracturevug-matrix) sections with equal permeability, due to the influence of the first homogeneous section and the middle dual-porosity section, the interporosity flow from vug into fracture becomes weaker than the interporosity flow from matrix into fracture. e dip of the interporosity flow in dual-porosity section becomes easier to recognize; however, the general flow regimes remain unchanged.
In the flux distribution in Figure 20, due to the two interporosity flows of the triple-porosity section and the interporosity flow of the dual-porosity section in the early stage, the flow rate of the two segments of wellbore in the triple-porosity section and the flow rate of the two segments of wellbore in the dual-porosity section are larger than the flow rate of the two segments of wellbore in the homogeneous section. With the appearance of pseudoradial flow regime, the rates for all six segments of the wellbore become equal. At late time, the flux distribution tends to be stable and the segment closer to the heel of the horizontal well gets higher flux. Based on this model, when the permeability of homogeneous section decreases from 100 mD to 40 mD and then to 10 mD, the values of pressure drop and pressure derivative increase ( Figure 21). In addition, the stage of the interporosity flow from both vug and matrix into fracture system becomes harder to recognize at early times. e dip of the interporosity flow in dual-porosity section becomes shallower. Other flow regimes remain unchanged.
When the permeability of dual-porosity section decreases from 200 mD to 100 mD and then to 50 mD, the values of pressure drop and pressure derivative increase as it can be observed in Figure 22. In addition, the stage of the interporosity flow from both vug and matrix into fracture becomes harder to recognize at early times. However, the dip of the interporosity flow in dual-porosity section   Figure 19: Type curves of flow stage division for a horizontal well traversing a three-section reservoir including homogeneous, dualporosity (fracture-matrix), and triple-porosity (fracture-vug-matrix) sections (solid curves). e triangles represent the type curves for a horizontal well traversing a three-section reservoir including homogeneous, triple-porosity (fracture-vug-matrix), and homogeneous sections. becomes easier to identify. Other flow regimes remain unchanged. Figure 23 shows that when the permeability of tripleporosity section decreases from 200 mD to 100 mD and then to 50 mD, the values of pressure drop and pressure derivative increase. In addition, the stage of the interporosity flow from both vug and matrix into fracture becomes easier to recognize ate early times. Moreover, the dips of the interporosity flows in triple-porosity section become deeper. Other flow regimes remain unchanged.

Field Application
To show the advantage of the proposed model with a horizontal well traversing multisection formation with different properties, a field case in Tarim Basin is analyzed. e input parameters of the well are as follows. e length of the well is 590 m. Before the pressure buildup test, this well has produced about 150 days with an average daily production of 80 m 3 /d. en, the pressure buildup test has been carried out for about 180 h. e log-log fitting curves are shown in Figure 24. We select the model with a horizontal well traversing a threesection reservoir including homogeneous, dual-porosity, and homogeneous sections to fit the measured data. e permeability of the dual-porosity section is greater than that of homogeneous sections. We set the main parameters as x 1 � 200 m, x 2 � 160 m, and x 3 � 300 m and the corresponding permeability of each section as K 1 � 15 mD, K 2 � 45 mD, and K 3 � 15 mD to fit the data with the model. It can be seen that the agreement between the measured data and the fitting curve is excellent, except in the wellbore storage period. Because the manometers may be installed downward the middle of the formation, the wellbore storage period lasts for a short time and is not represented on the curve.
e main fitting results are as follows: K N f � 45 mD, K m � 15 mD, and λ � 2.66 × 10 − 7 . From the fitting results and curves, it is shown that there is a dual-porosity section along the horizontal well. According to the microseismic monitoring results, it can be verified that there is a dual-porosity section in the formation consisting of fractures and pores.

Conclusions
A well test model for pressure transient analysis of horizontal wells extending through a carbonate reservoir  consisting of natural fractures, rock matrix, and vugs with different properties was presented in this study. Considering the pseudosteady interporosity flows from rock matrix and vugs into fractures, a multizone triple-porosity model was established and solved by using the point source function, Green's function, and coupling of multiple reservoir sections. e corresponding type curves were developed, and sensitivity analysis was carried out. is approach satisfies the need for pressure transient analysis for a horizontal well that traverses two or more regions with distinct properties in heterogeneous carbonate reservoirs. e following conclusions can be drawn: (1) A horizontal well traversing a three-section reservoir including homogeneous and dual-porosity (fracturematrix) sections identifies the stages of pseudosteady interporosity flow from matrix into fracture, fracture pseudoradial flow, system linear flow, system pseudoradial flow, and pseudosteady flow occur in sequence. e greater the difference of permeability between the dual-porosity section and the two homogeneous sections, the deeper the valley on the pressure derivative curve. (2) A horizontal well traversing a three-section reservoir including homogeneous and triple-porosity (fracture-vug-matrix) sections recognizes the stages of pseudosteady interporosity flow from vug and matrix into fracture, fracture pseudoradial flow, system linear flow, system pseudoradial flow, and pseudosteady flow occur sequentially. e greater the permeability of the triple-porosity section, the deeper the dips of the interporosity flows in this section on the pressure derivative curve.
(3) It is shown that, in the early stage, the flow rate in the dual-porosity/triple-porosity section is larger than the flow rate in the homogeneous section. With the advent of the pseudoradial flow regime, the rates for all segments of the wellbore become equal. At late times, the flux distribution tends to be stable and the segment closer to the heel of the horizontal well gets higher flux.
ere are still some issues to be further explored. Acid fracturing techniques and hydraulic fracturing techniques have been widely practiced to enhance the well productivity and ultimate recovery for low permeability carbonate reservoirs and shale gas reservoirs. We can extend our model from the horizontal well to the multifractured horizontal well to study the effect of different fracture networks. Based on the pressure transient analysis of the multifractured horizontal well, we can further study the production decline analysis to improve the well productivity and ultimate recovery.