Abstract

A dual fractal reservoir transient flow model was created by embedding a fracture system simulated by a tree-shaped fractal network into a matrix system simulated by fractal porous media. The dimensionless bottom hole pressure model was created using the Laplace transform and Stehfest numerical inversion methods. According to the model's solution, the bilogarithmic type curves of the dual fractal reservoirs are illustrated, and the influence of different fractal factors on pressure transient responses is discussed. This semianalytical model provides a practical and reliable method for empirical applications.

1. Introduction

Numerous researchers have simulated the nonuniform distribution of fractures using fractal networks and have studied fluid flow behavior in fractured reservoirs. By assuming that the fracture network is fractal, Camacho-Velázquez et al. [1] studied the production decline behavior in a naturally fractured reservoir. Zhang and Tong [2] introduced a stress-sensitive coefficient and built a transient pressure analysis model for fractal reservoirs which considers stress-sensitive effects.

Jafari and Babadagli [3] illustrated the 3D permeability distribution of a reservoir using outcrop, well log, and well test data which served as the basis for applying fractal networks to a reservoir. Zhang et al. [4] solved a nonlinear flow model for a stress-sensitive dual media fractal reservoir using a finite element method.

Previous scholars often embedded fractal structures into matrix networks by using straight or intersecting lines, but this did not correctly simulate well bottom radial flow. This type of radial flow has not been sufficiently studied in the underground seepage and oil development fields.

By referring to a plant lamina’s bifurcation structure, Wechsatol et al. [5] used a tree-shaped fractal structure that connected center points to different circles. The fractal network could simulate the radial flow tending toward the well bottom. Based on their research, the construction method and optimization rules for tree-shaped fractal structures [6] were formulated.

Xu and Yu [7] presented a tree-shaped fractal flow model that considered the dynamic behavior of branching tubes in a tree-shaped fractal network. Based on this model, our model’s transport properties and mass transfer capabilities [8, 9] were analyzed.

Given the capillary pressure effect, the starting pressure gradient influence, and pore fractal characteristics, Yun et al. [10] developed a fractal model that describes Bingham fluid flow in porous media. Based on Yun et al. [10], Wang et al. [11] proposed a tree-shaped fractal model that considered the influence of the starting pressure gradient on Bingham fluid seepage in a porous medium.

In this paper, fractures are simulated using a tree-shaped fractal network, as it accurately simulates radial flow tending to the well bottom, and the matrix system is simulated using fractal porous media. A transient flow model of dual fractal reservoirs is then presented by embedding the fracture network into a matrix system. Factors influencing the dynamic characteristics of transient pressure responses in dual fractal reservoirs are analyzed. This semianalytical model provides a practical and reliable method for empirical applications.

2. Physical Model

Figure 1 shows a well located in the reservoir center, , with a thickness, , and a well radius, . The physical model assumptions are as follows.(1)The fractal porous media is divided into M annular sections in a tree-shaped fractal network.(2)Each fracture section’s properties are different, but the fluid properties are identical. Matrix permeability is much lower than that of the fracture. (3)Rock and single-phase fluid are slightly compressible causing isothermal flow to be considered. Capillary pressure and gravity effects are neglected.(4)Fluid flows to the wellbore only through the fracture system. Fluid flow in the matrix and the fracture system of each section satisfies the linear flow rule.

3. Mathematical Model

3.1. Matrix System

According to fractal geometry theory, the fractal scaling law can be used to describe the cumulative size distribution of pores in matrix system [12]: where is pore diameter, is maximum pore diameter, and is the fractal dimension of the pore space.

Tortuosity is often used to describe flow path tortuosity, as flow in porous media is tortuous. The matrix flow path tortuosity is defined as [13]: where is the actual length of the tortuous flow path and is the straight length along the macroscopic pressure gradient.

The matrix system porosity and permeability are defined as [14]: where is the unit cell area.

3.2. Fracture System

Fracture system parameters are directly generated from a tree-shaped fractal network. Sets of branch structures form the tree-shaped fractal network. During network generation, the branches on each level must end up on the same circle with all circles having the same center of origin, . tubes which start at make up the tree-shaped fractal network. The tube’s initial length and diameter are and , respectively. The double branches (), whose angles are and total network branch levels are , are applied in this network. Furthermore, two scale factors are used in this fractal network, length ratio, , and diameter ratio, . The branch tube is assumed to be smooth, and tube wall thickness is ignored.

For the th level fracture, length is given by

For the th level fracture, diameter is

Distance from the well to each section’s boundary is defined as the radial distance, which is expressed by [8]:

Xu et al. [8] proposed the th level permeability expression of fracture system, which is expressed by

The th section tortuosity of fracture systems can be obtained by the following expression:

Substituting (6) and (9) into (8), the th section permeability in the fracture system is

For the th section, total system volume can be calculated by where is the reservoir thickness.

For the th section, pore volume of fracture systems can be calculated by

For the th section, total system volume, , is expressed as , where is total volume of the matrix system. is related to the pore volume of matrix system, , and expressed as , where is the matrix system porosity. Thus, can be calculated as

For the kth section, porosity of a fracture system can be obtained by dividing (12) by (11):

For the kth section, porosity of a matrix system can be obtained by the following expression:

The permeability of a fracture system, , and the porosity of a fracture system, , do not change with the radial distance, r, in traditional double porosity (fracture and matrix system) reservoir transient flow models [15]. In order to compare dual fractal reservoir transient flow models with double porosity reservoir transient flow models, we have to clarify how to keep and independent of r in dual fractal reservoir transient flow models.

Under the condition of double branches (), the permeability, , and porosity, , of every section in a fracture system are equal; that is,

Substituting (10) and (14) into (16), we can obtain: , , and   .

When is smaller than 0.707, permeability of a fracture system increases with the radius of a dual fractal reservoir. When is larger than 0.707, permeability of a fracture system decreases with the radius of a dual fractal reservoir.

3.3. Dual Fractal Reservoir

According to the physical model, the flow mathematical model of a dual fractal reservoir can be described as follows.

Governing differential equations in a dual fractal reservoir, we have the following.

For fracture system [15],

For matrix system [15],

Initial condition:

Interface connecting conditions of each zone, pressure continuity [16]:

Interface connecting conditions of each zone, rate continuity [16]:

Well production condition [17]:

External boundary condition (infinite):

External boundary condition (constant pressure):

External boundary condition (closed):

4. Mathematical Model Solution

To simplify the mathematical model and its solution, dimensionless parameters are defined as follows [18, 19].

The dimensionless pressure of the fracture system of the kth section:

The dimensionless pressure of the matrix system of the kth section:

The dimensionless effective radius:

The dimensionless effective interface radius:

The dimensionless effective time:

The dimensionless effective wellbore storage coefficient:

The fluid capacitance coefficient of the kth section:

The interporosity flow coefficient of the kth section:

Substituting (26)–(33) into (17)–(25), dimensionless mathematical models are obtained [20].

For a fracture system, the governing differential equation in a dual fractal reservoir is

For matrix system,

Initial condition:

Interface connecting conditions of each zone, pressure continuity:

Interface connecting conditions of each zone, rate continuity:

Well production condition:

External boundary condition (infinite):

External boundary condition (constant pressure):

External boundary condition (closed):

The flow mathematical model in Laplace space is obtained by taking the Laplace transformation of (34)–(42) based on /. The flow mathematical model is as follows [20].

For fracture system,

For matrix system,

Initial condition:

Interface connecting conditions of each zone, pressure continuity:

Interface connecting conditions of each zone, rate continuity:

Well production condition:

External boundary condition (infinite):

External boundary condition (constant pressure),

External boundary condition (closed).

In (43) and (44), the general solution of is calculated by

The derivative of in (52) is calculated by

Substitute (52) and (53) into well production condition equation (48);

Substitute (52) and (53) into interface connecting condition equation (46);

Substitute (52) and (53) into interface connecting condition equation (47);

Substitute (52) and (53) into external boundary condition equations (49)–(51);

, , and () can be obtained by solving the simultaneous equations (54)–(59). In (54)–(59), the tree-shaped fractal network parameters can be directly used to express the parameters of the permeability ratio, , the dimensionless effective interface radius, , and the function, .

The dimensionless effective interface radius expression can be derived by substituting (7) into (29):

The permeability ratio can be calculated using (10):

The expression of the function is as follows: where

The interporosity flow coefficient, , can be obtained by substituting (4) and (10) into (33), and the fluid capacitance coefficient, , can be obtained by substituting (3), (14), and (15) into (32).

5. Analysis of Type Curve Characteristics

Dimensionless bottom hole pressure in Laplace space, , is obtained by solving the linear equations (54)–(57) using the Stehfest numerical inversion method. The bilogarithmic type curves of the dual fractal reservoirs can then be illustrated.

In a condition of closed top and bottom boundary, the transient flow process, which has six flow regimes, can be clearly shown (Figure 2). The full and dashed lines represent pressure and pressure derivative curves, respectively. Regime 1 is the pure wellbore storage regime. Pressure and its derivative curves appear as upward straight lines with a slope of 1. Regime 2 is the transition flow regime. The shape of the derivative curve looks like an “arch,” which is influenced by the wellbore storage coefficient and skin factor. Regime 3 is the fracture system inter-porosity flow regime. The pressure derivative curve is V shaped, which depicts the response of inter-porosity flow between the fractures that are heterogeneously distributed. This inter-porosity flow regime is caused by spread of the pressure wave through the fracture system. Regime 4 is the fracture system radial flow regime. Slope of the pressure derivative curve is zero. In this scenario, the pressure wave spreads through the whole fracture system and begins to spread to the matrix system. Regime 5 is the inter-porosity flow regime of matrix system to fracture system. The pressure derivative curve is also V shaped. However, it is influenced by spread of the pressure wave through the matrix system. Regime 6 is the total system radial flow regime. The pressure derivative curve converges to a horizontal line, which depicts the response of the pressure dynamic balance state in the whole system.

Figure 3 shows the type curve characteristics affected by pore tortuosity, . As increases, the two V shaped curves deepen and shift to the right, which indicates a longer lasting inter-porosity flow regime and a delayed occurring time of the radial flow regime in the fracture system. It also indicates a delayed inter-porosity flow regime from matrix to fracture. An increase in the initial branch number, N, has a similar influence on flow regimes 3, 4, and 5 as an increase in , but the V shaped curves associated with the pressure derivative become shallower and shift to the right (Figure 4).

Figure 5 exhibits the type curve characteristics affected by pore fractal dimension, . Permeability of matrix system, , and porosity of matrix system, , increase with an increase in . A larger leads to greater flow capacity in the matrix system with an earlier transition to regime 5. A larger leads to a larger supplying capacity in the matrix system, and regime 5 occurs earlier and lasts longer. As increases, regime 5 occurs earlier and lasts longer, which is depicted as a deeper and wider second V shaped pressure derivative curve.

Figure 6 shows the type curve characteristics affected by branch angle, θ. Permeability of fracture system, , decreases with an increase in θ. A smaller leads to lower flow capacity in the fracture system with an earlier transition to regime 5. Porosity of fracture system, , increases with an increase in θ. A larger leads to greater supplying capacity in the fracture system, and regime 5 occurs later with a shorter duration. When θ increases, it has the opposite effect on flow regime 5 as an increase in .

Figure 7 exhibits the type curve characteristics affected by the length ratio, . decreases with an increase in α. A smaller leads to a lower supplying capacity in the fracture system and an earlier transition to a longer lasting regime 5. When a large increases, regime 5 occurs earlier, which manifests in a deeper and wider second V shaped type curve.

Figure 8 exhibits the type curve characteristics affected by diameter ratio, . and increase with an increase in β. A larger leads to larger flow capacity in the fracture system and a later transition to regime 5. A larger leads to greater supplying capacity in the fracture system with a later and shorter regime 5. When β is smaller than 0.707, increases with an increase in r, and when β is greater than 0.707, decreases with r. β affects all regimes except for pure wellbore storage and transition flow regime, which are not affected by . A larger β leads to a lower location of the dimensionless pressure curve, and regime 5 occurs later, resulting in a shallower and narrower second V shaped type curve. When β equals the critical value of 0.707 (16), the horizontal line representing regime 6 equals 0.5.

Figures 9 and 10 exhibit the type curve characteristics affected by total branch level, M, when diameter ratio, , is 0.65 and 0.75, respectively. Radius, , increases with an increase in M. Additionally, r enhances the type curve characteristics affected by . When is less than 0.707, a large M leads to a higher dimensionless pressure curve, and, vice versa, when is greater than 0.707, a large M leads to a lower dimensionless pressure curve.

6. Conclusions

The transient flow model for pressure responses in dual fractal reservoirs is established and solved, type curves are illustrated, and dual fractal flow behavior characteristics are analyzed. The following conclusions were obtained.(1)Fracture and matrix systems can be simulated using a tree-shaped fractal network and fractal porous media, respectively.(2)Six flow regimes for pressure type curves can be established. Type curves are dominated by fracture and matrix fractal parameters. These various parameters affect type curves differently.(3)Type curves are dominated by external boundary conditions, fractal parameters the fluid capacitance coefficient and the inter-porosity flow factor.(4)Semianalytical dual fractal modeling is suitable for various naturally fractured oil or gas reservoirs and provides a practical method to solve empirical cases.

Acknowledgments

The authors are grateful for financial support from the National Science Fund for Distinguished Young Scholars of China (Grant no. 51125019), the National Key Basic Research and Development Program of China (Grant no. 2011CB201005), and the Science and Technology Innovation Fund of Southwest Petroleum University (Grant no. GIFSB0701).