Abstract

The nonlinear dual-porosity flow model, specifically considering the quadratic pressure gradient term, wellbore storage coefficient, well skin factor, and interporosity flow of matrix to natural fractures, was established for well production in a naturally fractured formation and then solved using a semianalytical method, including Laplace transform and a transformation of the pressure function. Analytical solution of the model in Laplace space was converted to numerical solution in real space using Stehfest numerical inversion. Nonlinear flow process for well production in a naturally fractured formation with different external boundaries was simulated and analyzed using standard pressure curves. Influence of the quadratic pressure gradient coefficient on pressure curves was studied qualitatively and quantitatively in conditions of a group of fixed model parameters. The research results show that the semianalytical modelling method is applicable in simulating the nonlinear dual-porosity flow behavior.

1. Introduction

The fluid flow in porous media was proved to be of nonlinearity [13]. The nonlinear flow models with the quadratic pressure gradient term for well production in underground formations were especially investigated in the past. Early in 1989, a steady-state and semisteady-state flow model with the quadratic pressure gradient term was presented [4]. In 1991, a transient flow model with the quadratic gradient term was established and solved [5]. In 1993, the pressure distribution of a radial flow model with the quadratic gradient term was studied and plotted [6, 7]. In 1994, a nonlinear flow model for a dual-porosity formation was established [8, 9]. In 1996, the analytic solution to the nonlinear diffusion equation with the quadratic pressure gradient by considering a constant compressibility of fluid was derived [10]. In 1998, the solutions between linear and nonlinear models were compared [11] and a special nonlinear model for variable-rate well-tests was especially researched [12]. In 2004, the exact solution of a flow model with the quadratic gradient term by considering wellbore storage effect was deduced [13]. In 2005, a fractal model with the quadratic gradient term was presented [14]. In 2008, a class of nonlinear type equations for porous flow was especially proved to be of existence using mathematical theory [2]. In 2009, the nonlinear equation of pressure diffusion was solved using a Hopf-Cole transformation [15]. In 2010, a spherical flow model for a partial perforation well in a formation with a larger thickness was established [16]. In 2013, the effects of the quadratic gradient term on the pressure curves and pressure derivative curves were analyzed qualitatively and quantitatively for homogenous and multiple-zone composite reservoirs [1719]. Despite copious literatures on the nonlinear flow subject with the quadratic pressure gradient term, most of them were focused on the single porosity medium and only two papers [8, 9] studied the dual-porosity media; however the nonlinear dual-porosity flow model of Bai suffered some limitations in face of the real world: well skin and wellbore storage effects that actually exist in every real well were not considered in his model; a standard set of log-log type curves for nonlinear flow process analysis was not developed; quantitative difference between nonlinear and linear models was not compared. Therefore, the main task of this paper is to address the three issues and clearly show the readers the nonlinear dual-porosity flow behavior in a different way.

2. Physical Model

Physical model assumptions are as follows.

For a single well production at constant rate in an underground dual-porosity formation saturated by a single-phase liquid (oil or water), the external boundary of formation may be infinite or closed or constant pressure.

Underground dual-porosity formation is constructed by natural fracture system and matrix system. It is supposed that fracture permeability is far larger than matrix permeability, so the flow pathway directly connected with wellbore is considered as fracture system.

Slightly compressible rock and liquid with a constant compressibility are considered.

The porous flow is isothermal and conforms to Darcy’s law and the gravity and capillary forces are ignored for the flow, and the effect of pressure on fluid viscosity is also ignored for the flow.

Skin effect is considered (near the wellbore where the formation could be damaged by drilling and completion operations there would be an additional pressure drop during production, with the “skin” being a reflection of additional pressure drop).

Wellbore storage effect is considered (in the beginning of opening well, the fluid stored in wellbore starts to flow, and the fluid in the formation does not flow).

At time , pressure is uniformly distributed in formation, equal to initial pressure ().

3. Mathematical Model

3.1. Establishment of Mathematical Model

The nonlinear governing equation of fluid flow in the radial cylindrical system is where is radial cylindrical coordinate; is matrix permeability, μm2; is fracture permeability, μm2; is matrix pressure, MPa; is fracture pressure, MPa; μ is viscosity, mPa·s; is matrix porosity, fraction; is fracture porosity, fraction; is liquid compressibility, MPa−1; is total compressibility of rock and liquid of matrix system, MPa−1; is total compressibility of rock and liquid of fracture system, MPa−1; is geometric shape factor of matrix block, m−2; is well production time, ; the subscripts “m” and “f” represent “matrix” and “fracture,” respectively.

The second power of the pressure gradient in (1) is called the quadratic pressure gradient term.

Initial condition is where is initial formation pressure, MPa.

Well production condition at constant rate production based on effective radius is where is fluid volume factor, dimensionless; is wellbore storage coefficient, m3/MPa; is wellbore pressure, MPa; is well rate at wellhead, m3/d; is effective wellbore radius, m; is formation thickness, m.

The effective wellbore radius is defined as [20] where is real wellbore radius, m; is skin factor, dimensionless: where is external boundary radius, .

The following dimensionless definitions are introduced to solve the mathematical model:dimensionless radial distance ;skin factor ; here is additional pressure drop near wellbore;dimensionless production time ;dimensionless fracture pressure ;dimensionless fracture pressure ;dimensionless wellbore storage coefficient ;dimensionless time for dual media reservoir ;interporosity flow factor of matrix to fracture systems ;fluid capacitance coefficient of fracture subsystem ;fluid capacitance coefficient of matrix subsystem , ;dimensionless coefficient of nonlinear term .

The dimensionless mathematical model is as follows.

The dimensionless governing equation of fluid flow in a radial cylindrical system is

In real cases, cannot be equal to 0 for the nonlinear flow model. If we take a limit of 0 to (), the nonlinear model will be reduced to the conventional linear model.

Initial conditions are

Well production condition at constant rate is

External boundary conditions are

3.2. Solution to Dimensionless Mathematical Model

Because (7) is a nonlinear equation with two unknown functions of and , the first- and second-order derivatives of to , the quadratic power of the first-order derivative, and the first-order derivative of to , it is hard to linearize the nonlinear equation. Therefore, we have to seek the analytical solution of the nonlinear model in the other way.

We introduce the Laplace transform on the basis of : where is the variable in real space; is the variable in Laplace space; is the dimensionless time in real space; is the time in Laplace space.

The dimensionless mathematical model in Laplace space is as follows.

The dimensionless governing equation of fluid flow in a radial cylindrical system is

By (14),

Substituting (15) into (13),

If we set and , the model is reduced to the single porosity medium (homogeneous) model, so (17) can be changed by

In other words, under the same conditions of well production and formation boundary, the difference of dual-porosity model with homogeneous model in Laplace space only exhibits the function expression of , which means the flow governing equation of homogeneous model in Laplace space is also (16). Therefore if we obtain the solution of homogeneous model in Laplace space firstly, we will easily obtain the solution of dual-porosity model in Laplace space by changing the function of .

Because (16) is a nonlinear equation with , the first- and second-order derivatives of to , and the quadratic power of the first-order derivative in Laplace space, it is hard to linearize the nonlinear equation. Therefore, we have to seek the analytical solution of homogeneous model back to real space firstly.

For homogeneous model, the flow governing equation in real space can be expressed by

Well production and formation external boundary conditions are the same as (9)–(11).

The following variable modification [11, 16] is introduced to linearize (19):

Then

Substitute (20)-(21) into (19), (9)–(11):

The homogeneous model in Laplace space via Laplace transform can be written by where is substitution variable of dimensionless fracture pressure in Laplace space and is the value at the wall of wellbore .

According to the previous derivation, we can obtain the following governing equation of dual-porosity model by changing from (18) to (17) for (23):

The general solution of (28) is

Substitute (29) into (28) and (24):

Substitute (29) into (25)–(27): where and are undetermined coefficients; is modified Bessel function of the first kind, zero order; is modified Bessel function of the first kind, first order; is modified Bessel function of the second kind, zero order; is modified Bessel function of the second kind, first order.

In (30)-(31), there are three unknown numbers and three equations, and solutions to the model in Laplace space can be easily obtained by using linear algebra, such as a Gauss-Jordan reduction.

In real space, and the derivative can be obtained using a Stehfest numerical inversion [21] to convert back to , and then dimensionless wellbore pressure and the derivative can be obtained by substituting into (20). The standard log-log type curves of well-test analysis of and () versus can then be obtained.

The difference of the nonlinear model with the conventional linear model is that the nonlinear model contains the quadratic gradient term; therefore the solutions of the linear and nonlinear models are different and the effect of the quadratic gradient term exhibits the difference of solutions in pressure transients and pressure derivative transients controlled by . In the following, qualitative and quantitative analyses will be implemented to compare the solutions between the linear and nonlinear models.

4. Analysis of Nonlinear Flow Characteristics

4.1. Analysis of Nonlinear Flow Processes

Type curves reflect properties of underground formations. Type curves graphically show the process and characteristics of fluid flow in formations. Figures 14 show the type curves of pressure transients for well production in an underground dual-porosity media formation.

Figure 1 shows the type curves of the nonlinear dual-porosity flow model with an infinite formation boundary. The pressure transients were simulated using a set of fixed parameters and a group of varying nonlinear term coefficients (). The type curves were obviously controlled by . Pressure transients were simulated by setting as a limit of 0, 0.01, and 0.05. The curves of “” are the curves associated with the conventional linear flow model (see the derivative curve ① in Figure 1). The nonlinearity of fluid flow influenced the pressure transients positively. Different formation has different value of . We can calculate the value of according to the definition of , and usually the value of varies from 0 to 0.05; therefore we simulated the type curves using 0.01 and 0.05 in the context.

Five main flow regimes can be observed.

(i) Regime I, pure wellbore storage regime: there are no differences in type curves between the linear and nonlinear models because fluid in formation does not start to flow and the influence of the nonlinear quadratic pressure gradient term is only produced for the flow in formation. Wellbore pressure transients are not affected by the nonlinearity of fluid flow in this regime.

(ii) Regime I, skin effect regime: there are little differences in type curves between the linear and nonlinear models. The curves of the nonlinear model gradually deviate from those of the linear model with time elapsing. The deviation of pressure derivative curve is more obvious than that of pressure curve for the same . A larger means a stronger nonlinearity on the type curves.

(iii) Regime III, early fracture radial regime: fluid in the fracture system of formation radially flows into the wellbore and the fluid in the matrix system of formation does not start to flow. The differences in type curves associated with the linear and nonlinear models are more obvious than those of regime II. A larger results in a lager curve deviation of the nonlinear model from the linear model.

(iv) Regime IV, interporosity flow regime of matrix system to fracture system: the pressure derivative curve is V-shaped, which is the reflection of the interporosity flow of matrix to fracture. The differences in type curves associated with the linear and nonlinear models are more obvious than those of regimes II and III. The influences the location of type curves more heavily.

(v) Regime V, whole radial flow stage of fracture and matrix systems: the interporosity flow regime had ended. The differences in type curves associated with the linear and nonlinear models are the most obvious among the five flow regimes. The differences increase with the increase of . For conventional linear dual-porosity model, the slope of the pressure derivative curve in both regime III and regime V is zero; however, the pressure derivative curves of nonlinear dual-porosity model in regimes III and V are inclined instead of horizontal (see the derivative curves ② and ③ in Figure 1). As time elapsed, the pressure derivative curves of nonlinear model gradually deviate from the pressure derivative curve of linear model.

For constant pressure and closed boundaries, the type curves are similar to the nonlinear homogenous model of Guo and Nie [17] and are omitted here.

Except for nonlinear term coefficient (), type curves are sensitive to the other model parameters.

Figures 2 and 3 reflect the shape characteristics of type curves affected by interporosity flow factor of fluid from matrix system to fracture system . The pressure transients were simulated using a set of fixed parameters and changing from to , respectively, for and . Because represents the starting time of interporosity flow of matrix system to fracture system, therefore the bigger the is, the earlier the time of interporosity is, and as the value of decreases, the V-shaped derivative curve moves right (see the derivative curves ① and ② in Figures 2 and 3).

Figures 4 and 5 reflect the shape characteristics of type curves affected by interporosity flow factor of fluid capacitance coefficients of fracture and matrix systems ( and ). The pressure transients were simulated using a set of fixed parameters and changing from 0.1 to 0.01 and from 0.9 to 0.99, respectively, for and . Fluid capacitance coefficient of fracture system is coupled with fluid capacitance coefficient of fracture system , and they represent the relative fluid storage capacitance for fracture and matrix systems, respectively. A greater is the response of relatively more reserves in fracture system. As decreases, the V-shaped derivative curve becomes deeper and wider (see the derivative curves ① and ② in Figures 4 and 5).

Figures 6 and 7 reflect the shape characteristics of type curves affected by skin factor of well . The pressure transients were simulated using a set of fixed parameters and changing from 3 to 0, respectively, for and . Greater leads to higher location of dimensionless pressure curve.

4.2. Quantitative Analysis of Nonlinear Influence

“DV” and “RDV” are employed to show the quantitative differences between type curves [17, 19]: where DV is the differential value between linear and nonlinear models and RDV is the relative differential value between linear and nonlinear models.

Tables 1 and 2 show the quantitative differences of nonlinear influence on type curves for “” and “,” respectively. Dimensionless pressure values and dimensionless pressure derivative values in Tables 1 and 2 were calculated by fixing a group of parameters , and the corresponding type curves (curves ② and ③) are shown in Figure 1. The tables show that dimensionless pressure and its derivative are different between linear and nonlinear models.

The following quantitative differences can be seen.

For “” (see Table 1): when “” (regime I in Figure 1), DV and RDV of pressure are 0, and DV and RDV of pressure derivative are also 0; when “” (regime II in Figure 1), DV of pressure is 0.0375 and RDV of pressure is 1.90%, and DV of pressure derivative is 0.0451 and RDV of pressure derivative is 4.77%; when “” (regime III in Figure 1), DV of pressure is 0.3326 and RDV of pressure is 6.73%, and DV of pressure derivative is 0.0681 and RDV of pressure derivative is 13.42%; when “” (regime IV in Figure 1), DV of pressure is 0.6709 and RDV of pressure is 9.50%, and DV of pressure derivative is 0.0664 and RDV of pressure derivative is 18.51%; when “” (regime V in Figure 1), DV of pressure is 1.0570 and RDV of pressure is 11.77%, and DV of pressure derivative is 0.1108 and RDV of pressure derivative is 22.18%.

For “” (see Table 2): when “” (regime I in Figure 1), DV and RDV of pressure are 0, and DV of pressure derivative is 0.0001 and RDV of pressure derivative is 0.23%; when “” (regime II in Figure 1), DV of pressure is 0.1795 and RDV of pressure is 9.08%, and DV of pressure derivative is 0.2069 and RDV of pressure derivative is 21.88%; when “” (regime III in Figure 1), DV of pressure is 1.3082 and RDV of pressure is 26.49%, and DV of pressure derivative is 0.2381 and RDV of pressure derivative is 46.91%; when “” (regime IV in Figure 1), DV of pressure is 2.4193 and RDV of pressure is 34.24%, and DV of pressure derivative is 0.2060 and RDV of pressure derivative is 57.41%; when “” (regime V in Figure 1), DV of pressure is 3.5745 and RDV of pressure is 39.80%, and DV of pressure derivative is 0.3181 and RDV of pressure derivative is 63.68%.

As shown in the tables, DV and RDV increase with an increase in elapsed time. The RDV of pressure derivative is larger than that of pressure at a fixed time, such as the fact that when “” for “” in Table 1, RDV of pressure derivative is 13.42%, which is greater than that of pressure (6.73%). It can be observed from the tables that DV and RDV increase with an increase in , such as the fact that when “” RDV of pressure for “” is 6.73% and RDV of pressure for “” is 26.49%.

In general, qualitative and quantitative analyses show that the nonlinearity caused by the quadratic pressure gradient term influences the pressure transients positively.

5. Conclusions

A nonlinear flow model for well production in an underground dual-porosity media formation was derived. The transient flow behavior caused by well production was modeled. The type curves of pressure transients simulated using different values of model parameters showed the nonlinear dual-porosity processes and reflected the differences between nonlinear and linear models. The nonlinear term in the flow equation is suggested to be retained.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgments

The authors would like to acknowledge NSFC (National Natural Science Foundation of China) for supporting this paper through a Project (no. 51304164)—Research on the Pressure Dynamics of Multiple-Acidized-Fractured Horizontal Wells in Fractured-Vuggy Carbonate Formations. The paper also was financially supported by 973 Program (National Basic Research Program) of China under Grant no. 2014CB239205 and a basic research project under Grant no. 2015JY0132 from Science and Technology Department of Sichuan Province. The authors would also like to thank the reviewers and editors whose critical comments were very helpful in preparing this paper.