Transient Pressure Analysis of Inclined Well in Continuous Triple-Porosity Reservoirs With Dual-Permeability Behavior

Inclined wells has recently been adopted to develop fractured-vuggy carbonate hydrocarbon reservoirs (FVCHRs) composed by matrix, fracture and vug system. Therefore, it is significant for us to describe pressure transient of inclined well for FVCHRs. In this paper, it is assumed that vug and fracture system connect with wellbore, and inter-porosity flow from vug to fracture system appears. Therefore, the pressure transient responses model of inclined well with triple-porosity dual-permeability behavior was built. The model is solved by employing Laplace integral transform and finite cosine Fourier transform. Real-domain solution of the model is obtained by Stehfest inversion algorithm. On the basis model of the published paper, the solution of the simplified model of this paper was validated with horizontal and vertical well of FVCHRs with triple-porosity dual-permeability behavior, and results reach a good agreement. Type curve, according to pressure derivative curve characteristic, can be divided into eight flow regimes, which includes wellbore storage, skin reflect, early vertical radial flow, top and bottom boundary reflection, linear flow, inter-porosity flow from vug to fracture, inter-porosity from matrix to vug and fracture and pseudo-radial flow regimes. The influence of some vital parameters (inclination angle, inter-porosity flow coefficient, well length and permeability ratio etc.) on dimensionless pressure and its derivative curves are discussed in detail. The presented model can be used to understand pressure transient response characteristic of inclined wells in FVCHRs with dual-permeability behavior.


INTRODUCTION
Since the 1960s, numerous scholars mainly studied the flow characteristic of the fluid in naturally fractured reservoirs by the method of the pressure transient analysis (Warren and Root, 1963;Pruess and Narasimhan, 1985;Kuchuk and Biryukov, 2014;Kui et al., 2017;Tan et al., 2018;Yin, 2018;Song et al., 2019;Pang et al., 2020;Tontiwachwuthikul et al., 2020;Wang et al., 2020). Since a large number of the fractured-vuggy carbonate hydrocarbon reservoirs (FVCHRs) were found around the world, most of scholars have focused on studies of FVCHRs to understand flow characteristic of the fluid in the FVCHRs (Abadassah and Ershaghi, 1986;Liu et al., 2003;Sun et al., 2018). Compared with naturally fractured reservoirs, the complex structures of the FVCHRs posed significant research challenges towards scholars and scientists about pressure transient response of the FVCHRs. Therefore, the mathematic model descripting flow process of the fluid in FVCHRs is established to evaluate the properties and pressure transient characteristics of this type reservoirs.
During the last 20 years, pressure transient analysis is used widely to obtain some vital parameters of reservoirs. With the help of these parameters, reservoir engineers can evaluate oil and gas reservoirs effectively and develop a reasonable development plan for the FVCHRs. The source function method plays a vital role in the process of well test model solution. Gringarten and Ramey (1973) first put forward the real-domain point source analytical solution under the conditions of the different outer boundary. On the basis of the real-domain point source solution, Ozkan and Raghavan (1991; applied Green's function theory to expand the point source to the Laplace domain and the corresponding calculation method is given. Most of the world's giant oil ans gas produce from naturally fractured and vuggy carbonate reserviors with complex pore system. Vug porosity is common pore system in FVCHRs, and the size of the vug have an imoprtant influence on pressure and production performance. The FVCHRs can be subdivided into triple-continuum medium model and discontinuous medium according to the size of the vug. For the research of the pressure and production of triple-continuum medium model, Abadassah and Ershaghi (1986) first presented triple-porosity single-permeability well-test mathematical model for vertical well, and semi-log and log-log curves were drawn to analysis flow characteristic of the triple-continuum medium. Camacho-Velázquez et al. (2005) established and solved mathematical model of triple-porosity dual-permeability of the vertical well, which assumed that fluid of the natural fracture and vuggy porosity flow into wellbore directly, and the characteristic of pressure transient and decline curve are analyzed. Wu et al. (2004;2011) solved a triple-continuum pressure transient model of the multi-phase fluid by employing analytical and numerical method. Corbett et al. (2010) studied on fractured carbonate rocks by adopting the numerical well test, and summarized that numerical well testing has its limitations. Guo et al. (2012) established a triple-porosity single-permeability flow mathematical model of the horizontal well in a triple-continuum FVCHRs under different boundary condition, and standard pressure and its derivative log-log curves were drawn. Nie et al. (2012) presented pressure response mathematical model of the horizontal well in coal seam. Pressure and its derivative response curves are drawn in log-log plot, and influence of some vital parameters on type curves were analyzed. Jia et al. (2013) studied on flow issues of a porous-vuggy carbonate reservoir without fracture system. It is assumed that vug system is spherical and pressure distribute symmetrically. Corresponding mathematical model was established and solved by employing Laplace transform. Wang et al. (2018) developed a pressure transient mathematical model of acid fracturing wells in FVCHRs. Therefore, a composite system is considered, the inner region is made up with finite number of artificial fractures and wormholes, and the outer region is made up triple-porosity medium. Wang et al. (2014) derived well test mathematical model of hydraulic fractured well in FVCHRs by point sink integral method. This model assumed that fluid of matrix and natural fractures system flow into wormholes directly. The existence of large-scale cavities will have a noticeable influence on the transient pressure behavior. For the research of the pressure and production of discontinuous medium model, Gao et al. (2016) and Xing et al. (2018) presented a flow mathematical model, it is assumed that well is drilled through large-scale cavities and fluid flow into a filled cavity. Wellbore storage and skin factor were considered in their model. Du et al. (2019) presented an analytical well test model of fractured vuggy carbonate reservoirs. Pressure difference between the wellbore and large vug was considered in their model.
Because the reservoir geological condition becomes more and more complex, the technology of inclined wells has recently been adopted to solve the problem of onshore oil and gas reservoirs extraction. To estimate pressure transient behavior of inclined well, many scholars established pressure response mathematical model of inclined well. Cinco-Ley et al. (1975;2007) presented an analytical solution of inclined well in real domain and shown the existence of the pseudo-radial flow period at large time values. Yildiz and Ozkan Ozkan et al., 1999;Ozkan and Raghavan, 2000) studied transient pressure behavior and rate distributions of slanted well and perforated wells. Corresponding calculation method of slanted and perforated wells was shown to improve the computational efficiency and increase the accuracy. Feng and Liu (2015) established a well test mathematical model of inclined well with an impermeable fault by employing point source theory. Dong et al. (2018) presented a skin factor model of partially penetrated directionally-drilled wells in anisotropic reservoirs to study the stimulation effect of slanted.
This paper summarizes the published literature about well test model of the FVCHRs and slanted well, absorbs their highlights and studies on pressure transient behavior of the inclined well in FVCHRs. Therefore, the objective of this paper is to establish and solve well test mathematical model of inclined well in FVCHRs with triple-porosity dual-permeability behavior by employing Laplace transform and finite cosine Fourier transform. The fluid of natural fractures and vug system can flow into wellbore directly. Inter-porosity flows from matrix system to both natural fracture and vug system, and inter-porosity flow from vug system to a natural fracture system were also considered in this paper. Pressure and its derivative type curves under condition of constant production by adopting the Stehfest numerical inversion are drawn. Influence of some vital parameters on type curves are discussed in detail. fractures and vug system flow into the inclined wellbore directly. There is an inclined well in top and bottom closed boundary in FVCHRs. The length of the well is denoted by h w , inclination angle of the well is denoted by θ w , and height of reservoirs is denoted by h, the coordinates of point source in the reservoir is M(r, 0, z), the coordinates of any point on the well is M ′ (r ′ , 0, z ′ ) ( Figure 2). The end view of an inclined well opened practically in FVCHRs is shown in Figure 3. The fluid of the fracture and vug systems flow into the inclined wellbore directly, Inter-porosity flows from matrix system to both natural fracture and vug systems, and inter-porosity flow from vug system to natural fracture system were considered. This process is named tripleporosity dual-permeable flow behavior. In order to better establish mathematical model of an inclined well in FVCHRs, the other assumptions are as follows: (1) The production of the inclined well is denoted by q sc in FVCHRs of closed the top and bottom boundaries.
(2) The rock and fluid are slightly compressible, and the fluid flowing in the fractures and reservoirs is isothermal seepage.
(3) The inclined well is assumed to be penetrated fully or partially, and formation thickness and well length meet the relationships: The fluid in the natural fracture and vug system can flow into the inclined well and inter-porosity flow from vug system to a natural fracture system were considered. The fluid in the matrix system flow into inclined wellbore through vug and fracture system (see Figure 3A). (5) The permeability of the matrix system is k m , the permeability of the natural fracture is k f , the permeability of the vug is k f , and the permeability anisotropy of the horizontal and vertical directions is considered in FVCHRs. (6) The original formation pressure of the matrix, fracture and vug system is equal to initial pressure p e .

MATHEMATICAL MODELING
Based on the scheme of an inclined well as shown in Figure 2, the mathematical model of the inclined well is established in FVCHRs. In order to help reader to understand coordinate   system as shown in Figure 2, three different views of an inclined well in top and bottom closed boundary are given in Figure 4. According to the state equation, the motion equation and the continuity equation, the differential equation of matrix, vug and natural fracture systems is established in FVCHRs (Camacho-Velázquez et al., 2005). Inclined well could be regarded as a line source in 3D space. Therefore, to obtain solutions for inclined well, we must establish a 3D point source model. Then we can get a surface source solution through point source integral method.
(1) The two-dimensional seepage differential equation of natural fracture system can be described as: (2) The two -dimensional seepage differential equation of vug system can be described as: (3) The seepage differential equation of matrix system can be described as: (4) Initial condition can be described as: (5) Inner condition boundary can be described as: (6) Outer condition boundary can be described as: where: p V is pressure of vug system, MPa; p f is pressure of natural fracture system, MPa; p m is pressure of matrix system, MPa; p e is initial reservoirs pressure, MPa; k fh is horizontal permeability of natural fracture system, mD; k fv is vertical permeability of natural fracture system, mD; k mh is horizontal permeability of matrix system, mD; k Vh is horizontal permeability of vuggy system, mD; C tf is total compressibility factor of natural fracture system, MPa −1 ; C mf is total compressibility factor of matrix system, MPa −1 ; C Vf is total compressibility factor of vuggy system, MPa −1 ; φ f is porosity of natural fracture system, dimensionless; φ m is porosity of matrix system, dimensionless; φ V is porosity of vuggy system, dimensionless; t is production time, s; μ is fluid viscosity, cp;q is strength of continuous point source of reservoirs, m 2 /s; α is shape factor; x, y, z is Descartes coordinates, m; ε, δ is infinitesimal distance, m; q sc is production rate, m 3 /s; r w is wellbore radial, m. For convenience of derivation, the relevant dimensionless variables are defined in Table 1: Laplace transformation with respect to t D is then adopted to invert the above seepage model into Laplace domain (Van Everdingen and Hurst, 1949).
where: p jD is dimensionless pressure in Laplace domain; According to the dimensionless definition, the dimensionless forms of the differential equation and boundary conditions can be given by: where: k D is permeability ratio between vug and fracture, dimensionless; κ is permeability ratio, dimensionless; λ mf is inter-porosity flow coefficient from matrix to fracture, dimensionless; λ mv is inter-porosity flow coefficient from matrix to vuggy, dimensionless; λ vf is inter-porosity flow coefficient from vuggy to fracture, dimensionless; ω f is fracture storativity ratio, dimensionless; ω m is matrix storativity ratio, dimensionless; ω v is vuggy storativity ratio, dimensionless; s is Laplace variables; subscripts D is dimensionless variables; superscriptsis variables in the Laplace domain; By employing the finite cosine Fourier transform with respect to z D in Eq. 10, and expressing dimensionless pressures of natural fracture and vug systems by dimensionless matrix system pressure in Laplace domain after the finite cosine Fourier transform, Eq. 10 can be described as: where: , (n 0, 1, 2, ..., n); It is useful to express the Eq. 12 in polar coordinates, let the polar coordinates of the points (x D , y D ) and (x ′ D , y ′ D ) be given by (r D , θ) and (r ′ D , θ) respectively according to Raghavan and Ozkan (2017). We can do integration to the point source solution along inclined well in polar coordinates system as shown in Figure 4C. Therefore, the solution of the inclined well is obtain as following.
p wDincline 1 sh wD sin ψ ′ where: K 0 (x) is the second kind modified Bessel function; Gringarten et al. (1974) pointed out that the wellbore pressure of horizontal well and infinite conductivity hydraulic fractured well can be calculate by an equivalentpressure point x D 0.732, there is difference between horizontal well and inclined well. Fortunately, Cinco-Ley et al. (1975) given the equivalent-pressure point of inclined well in the polar coordinates. As following: In order to enhance calculation speed of the inclined well, Ozkan and Raghavan (2000) provides an efficient algorithm to compute transient pressure responses of inclined wells in Laplace domain. The principle of convolution can be used to obtain solution with the effects of wellbore storage and skin factor, which could be given by where: p wD is dimensionless wellbore pressure with wellbore storage and skin; p wDincline is dimensionless wellbore pressure without wellbore storage and skin; C D is dimensionless wellbore storage; S is skin factor. The solution of inclined well in real domain for FVCHRs with triple-porosity dual-permeability behavior is obtained by employing Stehfest (1970) numerical inversion algorithm. Dimensionless pressure and its derivative type curves are drawn analyzed.

MODEL VERIFICATION
Compared with vertical wells and horizontal wells, wellbore pressure solution of the inclined well is more universal, so our work expands analytical solution of vertical well and horizontal well in FVCHRs with triple-porosity dual-permeability behavior. However, in order to verify the correctness of our model and solution, the simplified model of this paper is compared with analytical solution of horizontal well for FVCHRs with tripleporosity dual-permeability behavior  and vertical well for FVCHRs with triple-porosity dualpermeability behavior (Camacho-Velázquez et al., 2005). Guo et al. (2012) established the mathematical model of the horizontal well under condition of different outer boundary in FVCHRs with triple-porosity dual-permeability behavior. What is different with our model is that they assumed that fluid of the naturally fractured and the matrix system are directly flow into the wellbore, inter-porosity flow from matrix to fracture and inter-porosity flow from vug to matrix and fracture is considered. In order to make verification result reliable, if parameters λ Vf , λ Vm , λ mf of the Guo' model are replaced by parameters λ mf , λ mV , λ Vf of this paper correspondingly and corresponding value of these parameters is same. On the basis of replacement, we can prove that our model is same with Guo' model. If θ w is 90°in our paper, and h w 1, 000 of our paper, inclined well of this paper will become horizontal well, which is similar to the Guo et al.' model with horizontal well length is 500. Other common parameters are same and their values are: h 100, r w 0.1, k hv 1, k D 1, κ 0.7, ω f 0.001, ω V 0.1, λ mf 10 − 12 , λ mV 10 − 10 , λ Vf 10 − 8 . As shown in Figure 5A, result reaches a good agreement. Camacho-Velázquez et al. (2005) established the mathematical model of the vertical well in FVCHRs with triple-porosity dual-permeability behavior. It is assumed that fluid of the fracture and the vug system are directly flow into the wellbore. Inter-porosity flow from vug to fracture and interporosity flow from matrix to vug and fracture is considered, which is same with our assumption. If inclination angle of wellbore is less than 5°and well is penetrated fully, we can think that our model is equal to vertical well. Other common FIGURE 8 | Type curve affected by length of inclined well (C D 0.1, S 0.1, h 100, r w 0.1, θ w 60 o , k hv 1, k D 1, κ 0.7, ω f 0.01, ω V 0.1, λ mf 10 − 12 , λ mV 10 − 10 , λ Vf 10 − 8 ) FIGURE 9 | Type curve affected by inclination angle of well (C D 0.1, S 0.1, h 100, h w h/cos(ψ′), r w 0.1, k hv 1,k D 1, κ 0.7, ω f 0.01, ω v 0.1, λ mf 10 − 12 , λ mV 10 − 10 , λ Vf 10 − 8 ) Frontiers in Energy Research | www.frontiersin.org December 2020 | Volume 8 | Article 601082 7 parameters are same and their values are: h 100, r w 0.1, k hv 1, k D 1, κ 0.7, ω f 0.001, ω V 0.1, λ mf 10 − 12 , λ mV 10 − 10 , λ Vf 10 − 8 . As shown in Figure 5A, result reaches a good agreement.
By comparing with horizontal and vertical well in FVCHRs with triple-porosity dual-permeability behavior, the results reached a good agreement, which verifies the correctness of our model and solution.

TYPE CURVE AND SENSITIVITY OF PARAMETERS
Type Curve Figure 6 shows type curves of dimensionless pressure and its derivative of inclined well with dual-permeability behavior. According to the characteristic of the dimensionless pressure derivative curve, the type curve of inclined well can be divided into eight flow regimes for FVCHRs with triple-porosity dualpermeability behavior, which mainly includes well storage regime, skin effect regime, early vertical radial flow regime ( Figure 7A), top and bottom boundary reflection ( Figure 7B), linear flow regime ( Figure 7C), inter-porosity flow and pseudoradial flow regime ( Figure 7D).
Regime①: Wellbore storage regime. During the pure wellbore storage, pressure and its derivative curves coincide and is a unitslope straight line.
Regime②: Skin effect regime. Due to consideration of skin, pressure derivative curve looks like a "hump." Regime③: Early vertical radial flow regime. Pressure derivative curve is a constant value horizontal line in log-log plot ( Figure 7A).
Regime④: Top and bottom boundary reflection regime. Pressure derivative curve take on a small step in log-log plot. However, this stage may disappear because of change of the length of wellbore and inclination angle of wellbore ( Figure 7B).
Regime⑤: Linear flow regime. Pressure derivative curve is 1/2-slope straight line in log-log plot. Shorter well length lead to disappearance of 1/2-slope straight line ( Figure 7C).
Regime⑥: Inter-porosity flow from vug to fracture. Because vug permeability is better than matrix permeability, interporosity flow from vug to fracture appear firstly. The curve of pressure derivative is V-shaped.
Regime⑦: Inter-porosity flow from matrix to vug and fracture. The curve of pressure derivative is V-shaped, which is controlled by parameters ω m and ω v .
Regime⑧: Pseudo-radial flow regime. Inter-porosity flows have already finished and the curve of pressure derivative is a 0.5-value horizontal line in log-log plot ( Figure 7D). Figure 8 shows type curve affected by length of inclined well in FVCHRs with triple-porosity dual-permeability behavior. Because short inclined well length can decrease contract area between reservoirs and wellbore, the short inclined well length can increase pressure loss that fluid flows into wellbore. Therefore, short inclined well length can lead to high position of pressure curve in the whole flow regime in log-log plot and make top and bottom boundary reflection regime disappear and spherical flow characteristic more obvious in pressure derivative curve. Figure 9 shows type curve affected by inclination angle in FVCHRs with triple-porosity dual-permeability behavior. It is assumed that length of inclined well keep constant when inclination angle changes. The large inclination angle change can only lead to disappearance of the top and bottom boundary reflection regime of pressure derivative curve without consideration of anisotropy ( Figure 9A). However, if the inclined well is penetrated fully, inclination angle has a more effect on pressure and its derivative curves than inclined well penetrated partly. With the increasing of inclination angle length of inclined well become length, which leads to low position of pressure curve in whole regime and low position of pressure derivative curve during early vertical radial flow regime ( Figure 9B). Figure 10 shows type curve affected by vug storativity ratio in FVCHRs with triple-porosity dual-permeability behavior. According the flow model shown Figure 3A, fluid of the vug system do not flow into fracture, but also flow into wellbore directly. Therefore, as shown in Figure 10, small vug storativity ratio can make the second V-shaped of the pressure derivative curve become more and more deep while the first V-shaped become more and more shallow. Figure 11 shows type curve affected by fracture storativity ratio in FVCHRs with triple-porosity dual-permeability behavior. According the flow model shown Figure 3A, fluid of fracture system only flow into wellbore directly. Therefore, as is shown in Figure 11, fracture storativity ratio has an influence on the first V-shaped. The first V-shaped become more and more deeply with the decrease of fracture storativity ratio. The position of pressure curve become higher during vertical radial flow regime and top and bottom boundary reflection regime. Figure 12 shows type curve affected by inter-porosity flow coefficient λ Vf in FVCHRs with triple-porosity dual-permeability behavior. According the flow model shown Figure 3A, fluid of fracture system only flow into wellbore directly. Therefore, as shown in Figure 12, inter-porosity flow coefficient λ Vf λ Vf has an influence on the first V-shaped. As λ Vf increase, starting time of inter-porosity flow from vug system to fracture system become earlier, which makes the first V-shaped of pressure derivative move towards left in log-log plot. Figure 13 shows type curve affected by inter-porosity flow coefficient λ mV in FVCHRs with triple-porosity dual-permeability behavior. As shown in Figure 13, inter-porosity flow coefficient λ mv has an influence on the second V-shaped. As λ mV increase, starting time of inter-porosity flow from matrix system to vug system become earlier, which makes the second V-shaped of pressure derivative move towards left in log-log plot. Of course, with the flow coefficient λ mV keep increasing, first V-shaped may be covered and first V-shaped become more deeply. Figure 14 shows type curve affected by inter-porosity flow coefficient λ mf in FVCHRs with triple-porosity dual-permeability behavior. As shown in Figure 14, inter-porosity flow coefficient λ mf has an influence on the second V-shaped. As λ mf increase, starting time of inter-porosity flow from matrix system to fracture system become earlier, which makes the second V-shaped of pressure derivative move towards left in log-log plot. Of course, because some of fluid of fracture system come from supplement of vug system, but also others come from supplement of matrix system. Therefore, when inter-porosity flow coefficient λ mf < λ vf , the change of inter-porosity flow coefficient λ mf has an influence on pressure derivative curve slightly. Figure 15 shows type curve affected by permeability ratio κ in FVCHRs with triple-porosity dual-permeability behavior. The large permeability ratio indicates that permeability of fracture system is higher than vug system, and vug storativity capacity is lower than fracture storativity capacity. Therefore, the first V-shaped of pressure derivative curve becomes more and more deeply with the increase of permeability ratio. Figure 16 shows type curve affected by permeability ratio k hv between horizontal permeability and vertical permeability in FVCHRs with triple-porosity dual-permeability behavior. The large k hv indicates that permeability of horizontal permeability is higher than vertical permeability. Because inclination angle of wellbore is considered in vertical direction, difference between vertical and horizontal flow is shown obviously during early vertical radial flow regime. The greater the value of k hv , the greater the loss of wellbore flow pressure in the vertical direction.