Abstract
In this paper, a new rate decline analysis model of horizontal wells with variable conductivity and uneven distribution of multiple fractures is proposed. By Laplace transformation, point source integration, and superposition principle, solutions of multiple infinite conductivity fractures in closed reservoirs are obtained. By coupling Fredholm integral equation of variable conductivity, linear equations of variable conductivity fractures in Laplace space are obtained. Gauss-Newton iteration, Duhamel convolution, and Stehfest numerical inversion method are used to obtain the bottom hole production solution. The accuracy of the results is verified by comparing with Eclipse software simulation. Then, the influence of some important reservoir and fracture parameters on the production is analyzed. The calculative results show that the smaller the fracture spacing is, the earlier the fracture begins to decline, the more the production will decrease; the change of different fracture length with the total fracture length unchanged has almost no effect on the production; the angle between fracture and -axis has an important effect on the production; the smaller the angle between fracture and -axis is, the stronger the interference between fractures is, the higher the production; the initial fracture conductivity affects the early production behavior, and the higher the initial fracture conductivity, the higher the production; the larger the fracture declines index, the lower the production, but the decreasing range gradually decreases with the increase of the decline index; the larger the reservoir drainage radius, the later the energy depletion stage, the higher the production. At last, a good fitting effect is obtained by fitting an example of oil field. The model proposed in this paper enriches the model base of rate decline analysis of fractured horizontal wells and lays a theoretical foundation for efficient development and practice of tight reservoirs.
1. Introduction
Multistage fracturing horizontal well is an important measure for the stimulation of low permeability reservoirs and especially for tight reservoirs. Multiple fractures generated by hydraulic fracturing can improve the near wellbore seepage mode and increase oil and gas seepage channels. Based on the theory of transient seepage flow, it is of much significance to guiding production practice to establish a productivity decline analysis model and method for multistage fractured horizontal wells and finally obtain reservoir parameters and well-controlled reserves using the typical curve fitting technology [1].
The transient well test analysis of hydraulic fractures reached the peak in the 1970s and 1980s and laid a theoretical foundation for the subsequent development of models for multistage fractured horizontal wells. Gringarten [2] firstly proposed practical typical curves of testing of vertically fractured wells with infinite conductivity and applied the curves. Cinco-Ley et al. [3] proposed a semianalytical model for the pressure of vertical fractures with finite conductivity in real space, which laid the theoretical foundation for test interpretation of vertically fractured wells. Cinco-Ley et al. [4] developed a dual-medium fractured well model and established typical curves. Some scholars have extended Cinco-Ley method to the pressure and production instability model for multistage fractured horizontal wells. Raghavan and Joshi [5] gave a steady-state productivity analysis method for multifracture systems with infinite conductivity based on the steady-state seepage theory. They conducted detailed theoretical derivation using the methods such as sink-source integration and pressure drop superposition. Wan and Aziz [6] proposed a new model for fractured horizontal wells. In their model, the included angle between an infinite conductivity fracture and a horizontal well can be arbitrary. Al-Kobaisi et al. [7] proposed a hybrid model for fractured horizontal wells with finite conductivity, which is very effective for early productivity analysis. Demarchos et al. [8] proposed a series of production optimization charts for fractured horizontal wells based on the idea of Economides, which are of much significance to the design of the construction of fractured horizontal wells. Brown [9] considered the “volume fracturing” characteristics of hydraulic fractures in tight reservoirs and represented the area between fractures as a double-porosity three-flow model. Luo et al. [10] proposed a fracture single-wing model and used it to solve the problem of noncoplanar fracture pressure instability analysis of multistage fractured horizontal wells. Based on predecessors’ achievements, some scholars later did some extended research [11–16]. Yuan et al. [17] proposed a workflow to analyze hydraulic fracture effect on hydraulic fractured horizontal well production in composite formation system. Zhao and Du analyzed the performance of horizontal wells in composite tight gas reservoirs considering stress sensitivity.
Great progress has been made in the pressure analysis of multistage fractured horizontal wells. Most scholars assume that the conductivity is constant, and fractures have equal length and are uniformly distributed. This will bring large errors to the production data analysis of fractured wells with variable conductivity in general. Domestic and foreign scholars have conducted a lot of proppant conductivity experiments to study the variable conductivity effect of fractures [18–21], where there are 3 spatial variable conductivity forms such as linear, exponential, and logarithmic ones. However, at present, there are relatively few rate decline models that consider both variable conductivity and uneven distribution of multiple fractures in horizontal wells. Sun et al. [22] established the single-fracture rate decline analysis method considering the variable conductivity of fractures, but it fails to consider uneven distribution of multiple fractures. The purpose of this paper is to establish a new rate decline model taking account of both the variable conductivity of fractures and uneven distribution of multiple fractures. In this paper, a point source seepage model for closed reservoirs has been firstly established; moreover, an infinite conductivity model with uneven distribution of multiple fractures has been obtained by means of Laplace transformation, point source integration, and superposition principle. Then, the Fredholm integral equation has been established, which takes account of the change of conductivity. By coupling the infinite conductivity multifracture model with the variable conductivity Fredhom integral equation, the system of linear equations of fracture flow rate and horizontal well pressure in Laplace space has been obtained. The real space rate solution has been obtained based on the convolution principle and Stehfest numerical inversion algorithm. Finally, the type curves and flow characteristics have been analyzed using the obtained solution, thereby establishing a theoretical foundation for the rate decline analysis of tight reservoirs.
2. Productivity Model for Horizontal Wells with Variable Conductivity Fractures Distributed Nonuniformly
2.1. Physical Model
Figures 1(a) and 1(b) are the 3D and 2D schematic diagrams of a model for a horizontal well with variable conductivity fractures distributed nonuniformly. After multistage fracturing of the horizontal well in the tight reservoir, multiple fractures are formed in the well (Figure 1(a)). This model is different from the conventional equal conductivity model, and fractures gradually pinch out from the center to the two ends.
(a)
(b)
The hypothesis of fracture physical model is very reasonable and conforms to the real fracture spatial distribution. The fractures can be unevenly distributed (Figure 1(b)); fluids mainly flow from the reservoir to the hydraulic fractures, and finally, the fluids from the hydraulic fractures converge to the horizontal wellbore (Figures 1(a) and 1(b)). The physical model is simplified as follows: (1)All rock properties such as reservoir permeability, reservoir porosity, fracture permeability, and fracture porosity are all constant. The compressibility factor and viscosity of fluids and reservoir thickness are also a constant(2)The hydraulic fracture formed by staged fracturing of the horizontal well is assumed to run through the whole reservoir; in addition, the length of each fracture is , and the fracture width varies on - plane(3)The reservoir belongs to the radial coordinate. Assume that fluids can spread to the boundaries after flowing for a certain period. Therefore, it can be assumed that the horizontal boundaries are circular and closed, and the upper and lower boundaries are impermeable(4)Isothermal Darcy seepage occurs in the reservoir; moreover, consider that the seepage in the fractures is of finite conductivity(5)The total production of all hydraulic fractures is a constant (6)The pressure is distributed uniformly in the initial state and equals to the initial pressure (7)All hydraulic fractures can be unevenly distributed, and their spacing, length, and angle can all be different
2.2. Point Source Model with Multiple Unevenly Distributed Fractures
According to the physical model assumptions, a tight reservoir subjected to multistage fracturing stimulation of horizontal wells is known; by combining the continuity equation, velocity equation, and state equation, the obtained point source seepage control system of equations for the whole reservoir is as follows (for the following symbols, see their annotations):
The initial condition of the reservoir is:
The closed external boundary condition of the reservoir is:
The production of the point source is changed, so the internal boundary condition of the point source is below:
The dimensionless definition expressions of the physical quantities are as follows:
After nondimensionalization and Laplace transform of equations (1) to (4), the governing equations and boundary condition equations can be rewritten. The governing equations can be rewritten as follows:
The internal boundary condition of the point source can be rewritten below:
The external boundary condition can be rewritten below:
By jointly solving equations (5)–(7), the point source equation can be generally solved as follows:
For obtaining the solution to multiple unevenly distributed fractures at any angle, assume that the position of any point of the reservoir and the point source is, respectively, (, ) and (, ), and the fracture reference length is . The maximum length of the fractures is defined as the reference length (a fracture is from the intersection of the horizontal well with the fracture to fracture pinch out); the length of each fracture is , and the distance from the point source to any point is . Then, the following geometrical relationship is established:
Substitute equation (10) into equation (9) to obtain:
where where is the Laplace space-time variable. Assuming that the number of fractures in the multistage fractured horizontal well is , the pressure drop of each section of fractures can be obtained as follows by means of point source integration, superposition principle, and coordinate transformation:
where is the included angle between the th fracture and the -axis; is the position of the fracture source; (, ) are the coordinates of the starting endpoint of the fracture. If each fracture is divided into sections, the pressure drop formula for the th section of the th fracture can be written as follows:
2.3. Hydraulic Fracture Flow Model
The flow in fractures can be assumed to be one-dimensional linear flow. According to the idea of Cinco-Ley [4], the dimensionless governing equation of flow in the single wing of a one-dimensional fracture can be written as follows:
The initial condition of a hydraulic fracture can be written as:
The inner boundary condition can be written as:
Integrate equation (15) from 0 to while using boundary equations (17) and (18) to obtain:
Substitute equation into equation (19) to obtain:
At the wellbore:
Substitute equation (21) into equation (20) to obtain:
Carry out the Laplace transformation to obtain:
Assume that a multistage fractured horizontal well consists of multiple single-wing fractures. If the th fracture wing is divided into sections, the pressure of the th section of fracture can be obtained from the discrete equation (24):
where , , and are the pressure and flow rate of the th section of the th fracture wing.
The specific expression of variable conductivity of a vertical fracture must use the laboratory experiment method or the fracturing experience formula of an actual reservoir. In this paper, the expression uses an exponential empirical formula referring to reference [19], which is as follows:
Figure 2 shows the variation of the fracture conductivity with the along-fracture distance at different exponent coefficients. As the along-fracture distance increases, the fracture conductivity decreases gradually; in addition, as the decline exponent increases, the decline amplitude gradually increases, and the fracture conductivity decreases sharply.
2.4. Coupling Solution to the Reservoir Model and Fracture Model
Assume that fracture wings are formed during horizontal well fracturing. If the th fracture wing is divided into sections, the total number of fracture sections is:
Based on the superposition principle, the pressure on the fracture surface shall be equal to that on the hydraulic fracture section, then:
Now, equations can be listed, but the Laplace space horizontal well pressure in equation (25) is also unknown, and there are m +1 unknowns. Assuming that the horizontal well is produced at a fixed flow rate, the sum of the flow rates of all fracture sections shall be equal to the fixed flow rate of the horizontal well; then:
By solving equations (27) and (28), flow values and the bottom hole pressure in the horizontal well in the Laplace space can be obtained. To get the production solution, it is needed to use the Duhamel convolution principle; then:
The dimensionless bottom hole pressure can be obtained through Stehfest numerical inversion.
3. Validation of the Presented Model
In this section, Eclipse commercial numerical simulation software is compared with our model to validate it. Some important parameters are as follows: reservoir area , reservoir permeability , fracture half-length , fracture permeability , number of fractures formed by fracturing , fracture spacing 100 m, and horizontal well pressure difference 2 MPa. The fractures are evenly distributed and are of equal length, and the initial fracture conductivity is . In order to ensure that the dimensionless drainage area is equal to the circle, , and the conductivity decline exponent is . The simulation validation result is shown in Figure 3. As shown in the figure, the result from the commercial software is in good agreement with the presented model, thus verifying the reliability of the model result.
4. Analysis and Application of Production Influence Factors
4.1. Influence of Fracture Spacing on Production
In the presented model, the dimensionless fracture length is , and the included angle between different fractures and the -axis of the reservoir is 90°, 90°, 90°, 90°, and 90°, respectively. The coordinate positions of different fractures are shown in Figure 4. The reservoir radius is ; the number of fractures is ; the initial fracture conductivity is ; the fracture decline exponent is , and the dimensionless fracture length is , 3, and 5, respectively. Figure 4 shows the pressure distribution caused under different fracture spacings at time . Due to the constant pressure production of the bottom hole of the horizontal well, the dimensionless pressure is lower at the position more distant from the pressure.
(a)
(b)
(c)
As shown in the figure, the smaller the fracture spacing, the darker the horizontal wellbore pressure color, the more intense the pressure change, and the more severe the pressure interference. Figure 5 shows the variation of production with time at different fracture spacings. As shown in the figure, all the production curves gradually decline with time. The three production curves of fractures at different fracture spacings coincide because of no interference between fractures when the initial is less than 0.1. When the dimensionless time , the production curves begin to show different decline trends. The smaller the fracture spacing is, the earlier the fracture starts to decline, and the more seriously the production decreases. After the pressure wave reaches the boundary, the larger the fracture spacing, the faster the production decline. This is because after the pressure wave reaches the boundary, the fracture at the far end is close to the boundary, so pressure attenuation is more serious.
4.2. Influence of Fracture Length on Production
In the presented model, the dimensionless fracture spacing is ; the coordinates of the fracture center point are (0,0), (2,0), (4,0), (6,0), and (8,0), respectively, and the included angle between different fractures and the -axis of the reservoir is 90°, 90°, 90°, 90°, and 90°, respectively. In addition, the reservoir radius is , and the number of fractures is . Actual fractures are defined as half of fractures. Therefore, in the simulation process, and the total length of these fractures is set to 5. Figure 6 shows the fracture distribution coordinates and the length of each fracture. The initial fracture conductivity is ; the fracture decline exponent is ; the number of fractures is , and the dimensionless total length of fractures is set to 5. Figure 6 shows the pressure field distribution under different fracture lengths at time . The figure shows the distribution of the equipotential lines of influence at the constant total fracture length. The more uneven the fracture length, the more irregular the shape of the equipotential lines. Figure 7 shows the influence of different fracture lengths on production at the constant total length. As shown in the figure, when the total fracture length is constant, different fracture lengths have almost no influence on the change of production.
(a)
(b)
(c)
4.3. Influence of the Dip Angle θ of the Fracture and the -Axis on Production
In the presented model, the dimensionless fracture spacing is , and the coordinates of the fracture center point are (0,0), (2,0), (4,0), (6,0), and (8,0), respectively. There are three cases of the included angle between different fractures and the -axis of the reservoir, such as 90°, 60°, and 30°. In addition, the reservoir radius is , and the number of fractures is . Actual fractures are defined as half of fractures. Therefore, in the simulation process. The position and angle of fractures in different cases are shown in Figure 8. Figure 8 shows the influence of different fracture dip angles on dimensionless production. As shown in the figure, as the fracture dip angle gradually decreases, the equipotential lines around the fracture become increasingly irregular, indicating that the direction of pressure wave propagation has changed. As the fracture dip angle becomes small, the contact area between the fractures and the reservoir in the direction gradually increases, and the interference between the fractures becomes small. Figure 9 shows the influence of the include angle between the fracture and the -axis on production. With the passage of time, the production gradually decreases. When the dimensionless time is more than 1000, pressure disturbances propagate to the boundary, and production rapidly declines. Throughout the flow period, as the fracture angle becomes small, the fracture production increases. This is because after the fracture angle becomes small, the interference between the fractures becomes small, thus greatly increasing the production. In the actual fracturing production process, the interference between fractures can be reduced to increase the production by forming low-angle fractures to the greatest extent possible.
(a)
(b)
(c)
4.4. Influence of Initial Fracture Conductivity on Production
In the presented model, the dimensionless fracture length is , and the included angle between different fractures and the -axis of the reservoir is 90°, 90°, 90°, 90°, and 90°, respectively. The reservoir radius is ; the number of fractures is ; the fracture decline exponent is , and the initial fracture conductivity is , 50, and 100, respectively. Figure 10 shows the influence of initial fracture conductivity on dimensionless production. As shown in the figure, the initial fracture conductivity mainly affects the early flow stage. When the dimensionless time is less than 0.01 in the early stage, the initial fracture conductivity affects the production to a very large extent. The higher the initial fracture conductivity, the higher the production. However, all production curves coincide in the midlate stage, and the initial fracture conductivity has almost no influence on the production.
4.5. Influence of Fracture Decline Exponent on Production
In the presented model, the dimensionless fracture length is , and the included angle between different fractures and the -axis of the reservoir is 90°, 90°, 90°, 90°, and 90°, respectively. The reservoir radius is ; the number of fractures is ; the initial fracture conductivity is , and the different fracture decline exponents are , 2, and 4. Figure 11 shows the influence of different fracture decline exponents on dimensionless production. As shown in the figure, the fracture decline exponent has an influence on the production in the whole flow period of production decline. The larger the fracture decline exponent, the lower the dimensionless production, but its decline amplitude gradually decreases.
4.6. Influence of Boundary Distance on Production
In the presented model, the dimensionless fracture length is , and the included angle between different fractures and the -axis of the reservoir is 90°, 90°, 90°, 90°, and 90°, respectively. The number of fractures is ; the initial fracture conductivity is ; the conductivity decline exponent is , and the dimensionless drainage radius is , 40, and 60. The magnitude of dimensionless drainage radius reflects that of reservoir drainage area, so the influence of dimensionless drainage radius on productivity reflects that of reservoir drainage area on productivity. According to the production curves in Figure 12, a larger reservoir drainage radius leads to later entering of energy depletion stage and higher production. As the drainage radius decreases, some flows in characteristic sections cannot appear and even will affect the linear flow period.
4.7. Example Analysis
Through fitting of the site-measured data of pressure and pressure derivatives based on the model in the paper, the parameters of a fractured reservoir after multistage fracturing of a horizontal well can be obtained as follows: matrix permeability , reservoir thickness , oil viscosity , compressibility coefficient , fracture decline exponent , fracture conductivity , initial formation pressure , bottom hole pressure , fracture half-length , fracture spacing , fracture number 5, and included angle between different fractures and the -axis of the reservoir is 60°, 60°, 60°, 60°, and 60°, respectively. The matched parameters are matrix permeability , fracture decline exponent , fracture conductivity , fracture half-length , fracture spacing, and included angle between different fractures and the -axis of the reservoir. Figure 13 shows the result of the comparison of actual production data with the presented model. According to the result, they tally with each other well, and the accuracy of the model presented in this paper has been verified.
5. Conclusions
(1)This paper has presented a new model for production calculation of horizontal wells with multiple unevenly distributed variable conductivity fractures in a closed reservoir by means of Laplace transformation, point source integration, and superposition principle(2)The accuracy of the proposed semianalytical solution has been verified through comparison with Eclipse commercial numerical simulation software(3)The influence of some reservoir and fracture sensitivity parameters on typical curves has been analyzed in detail, mainly including fracture spacing, nonuniform fracture length, the included angle between fractures, and the -axis, initial fracture conductivity, fracture decline exponent, and fracture drainage radius(4)Finally, the accuracy of the presented model has been verified through the analysis of an example, thereby providing theoretical guidance to the productivity decline analysis of multistage fractured horizontal wells in fractured reservoirs
Abbreviations
| : | Rectangular coordinates, m |
| : | Point source coordinates, m |
| : | Radial distance, m |
| : | Hydraulic fracture number |
| : | Reservoir pressure, MPa |
| : | Dimensionless reservoir pressure and hydraulic fracture pressure |
| : | Total productivity of horizontal wells, m3/d |
| , : | Reservoir permeability and hydraulic fracture permeability, respectively, m2 |
| : | Reservoir thickness, m |
| : | Drainage radius, m |
| : | Reference length, m |
| : | Fracture half-length, m |
| : | Fracture spacing, m |
| : | Initial reservoir pressure, Pa |
| : | Hydraulic fracture pressure, Pa |
| μo: | In-place oil viscosity, Pa·s |
| : | Volume factor of in-place oil, m3/m3 |
| : | Dimensionless initial fracture conductivity |
| : | Dimensionless fracture conductivity |
| : | Fracture conductivity index |
| qwf: | Bottom hole production of oil well, m3/d |
| : | Time, d |
| : | Laplace space time variable |
| : | Dimensionless time |
| : | Dimensionless hydraulic fracture surface flow rate. |
Data Availability
No data were used to support this study.
Conflicts of Interest
No conflicts to declare.