Next Article in Journal
Predictive Modeling for Microchannel Flow Boiling Heat Transfer under the Dual Effect of Gravity and Surface Modification
Next Article in Special Issue
Induced Casing Deformation in Hydraulically Fractured Shale Gas Wells: Risk Assessment, Early Warning, and Mitigation
Previous Article in Journal
Optimizing Microwave-Assisted Extraction from Levisticum officinale WDJ Koch Roots Using Pareto Optimal Solutions
Previous Article in Special Issue
Preparation and Application of CO2-Resistant Latex in Shale Reservoir Cementing
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

The Inversion Method of Shale Gas Effective Fracture Network Volume Based on Flow Back Data—A Case Study of Southern Sichuan Basin Shale

1
Shale Gas Research Institute, CNPC Southwest Oil and Gas Field Company, Chengdu 610051, China
2
Sichuan Key Laboratory of Shale Gas Evaluation and Exploitation, Chengdu 610051, China
3
National Key Laboratory of Oil and Gas Reservoir Geology and Exploitation, Southwest Petroleum University, Chengdu 610500, China
*
Author to whom correspondence should be addressed.
Processes 2024, 12(5), 1027; https://doi.org/10.3390/pr12051027
Submission received: 14 April 2024 / Revised: 15 May 2024 / Accepted: 16 May 2024 / Published: 18 May 2024

Abstract

:
Fracture network fracturing is pivotal for achieving the economical and efficient development of shale gas, with the connectivity among fracture networks playing a crucial role in reservoir stimulation effectiveness. However, flow back data that reflect fracture network connectivity information are often ignored, resulting in an inaccurate prediction of the effective fracture network volume (EFNV). The accurate calculation of the EFNV has become a key and difficult issue in the field of shale fracturing. For this reason, the accurate shale gas effective fracture network volume inversion method needs to be improved. Based on the flow back characteristics of fracturing fluids, a tree-shaped fractal fracture flow back mathematical model for inversion of EFNV was established and combined with fractal theory. A genetic algorithm workflow suitable for EFNV inversion of shale gas was constructed based on the flow back data after fracturing, and the fracture wells in southern Sichuan were used as an example to carry out the EFNV inversion. The reliability of the inversion model was verified by testing production, cumulative gas production, and microseismic results. The field application showed that the inversion method proposed in this paper can obtain tree-shaped fractal fracture network structure parameters, fracture system original pressure, matrix gas breakthrough pressure, fracture compressibility coefficient, reverse imbibition index, equivalent main fracture half length, and effective initial fracture volume (EIFV). The calculated results of the model belong to the same order of magnitude as those of the HD model and Alkouh model, and the model has stronger applicability. This research has important theoretical guiding significance and field application value for improving the accuracy of the EFNV calculation.

1. Introduction

The economical and effective development of shale gas plays an important role in the world’s energy structure. The shale gas revolution in the United States has achieved great success and thus the world energy pattern has been reshaped. In 2010, the first shale gas fracturing well was implemented in China. After 10 years of development, shale gas production has achieved a new leap [1]. The practice of shale gas fracturing has found that the information on the characteristics of complex fracture network hidden in the flow back data is very important. In the past, fracturing fluid flow back data were mostly ignored. In recent years, the shale gas effective fracture network volume and fracture characteristics information included in the fracturing fluid flow back data to evaluate hydraulic fracturing effects have attracted attention [2]. Due to the special two-phase seepage characteristics of shale gas, the fracturing evaluation method based on flow back is still in its infancy, and there are very few studies on the effective fracture network volume prediction of shale.
The fracture network that develops post-fracturing serves as the primary pathway for fluid flow back and gas production. It is of great significance to the evaluation of fracturing effect and productivity prediction. Effective fracture network volume information can be obtained by accurately interpreting flow back data [3]. It is of great significance to the evaluation of fracturing effect and productivity prediction. At present, microseismic monitoring is often used to analyze the effect of reservoir stimulation. The volume of reservoir stimulation determined by microseismic analysis rarely provides the actual recoverable reservoir volume of gas shale reservoirs. Reservoir stimulation volumes are overestimated or underestimated due to the occurrence of non-contributing microseismic events and the natural fractures [4]. A relatively new concept used to quantify the actual reservoir stimulated volume is the effective reservoir stimulated volume (ESRV) [5,6]. The ESRV is critical for estimating ultimate recovery, stimulated volume, optimal fracture length, and spacing. Some scholars have carefully identified high-frequency flow back data (flow rate and pressure) as a possible method to extract the properties of hydraulic fractures. A new mathematical model was used to approximate the effective volume of fractures and elucidate the mechanisms of early gas production [7] compared to previous flow back models [8,9,10,11]. The water and natural gas production data were analyzed in the Horn River shale, and the estimation of the effective fracture volume was based on the assumption of a two-phase tank model for the fracture system. In addition, due to the high uncertainty introduced by complex flow back models in the evaluation of fracture pore volume, a two-phase tank model was built for reducing parameter uncertainty in order to estimate the fracture pore volume accurately [12]. Comparing the estimated effective fracture pore volume with fracture design parameters (soaking time and proppant concentration) revealed that fracture closure drives single-phase flow back during further analysis of the effective fracture volume. Therefore, the effective fracture pore volume largely depends on the compressibility of the fracture [13]. The evaluation of fracture compressibility is essential for calculating the effective fracture volume, as well as evaluating the fracture volume change. It was found that fracture compressibility depends on how fracture porosity and pore size vary with effective stress. The fracture compressibility is mainly affected by porosity changes in propped fractures [14]. The harmonic decline (HD) model was used to estimate the water phase flow back and the initial effective fracture volume to study the change in the effective fracture volume. According to the findings, pressure depletion and fracture closure during early-time water flow back can lead to a loss of up to 30% of the effective fracture volume. When gas from the matrix enters the fracture network and offers adequate pressure support during late-time flow back, the rate of fracture volume loss reduces [15]. A gas–water two-phase flow model was established for a fractured horizontal well considering the influence of water saturation on gas–water permeability. Sensitive studies were performed to investigate the impact of various parameters on the productivity of a horizontal well [16,17]. The flow back period can be used to quantitatively identify the characteristics of the formed fracture network and calculate the effective fracture network volume; however, such data analysis is not commonly reported in the literature. Therefore, it is necessary to develop an appropriate model to match the data and acquire valuable information on the volume of the fracture network. In this paper, a tree-shaped fractal fracture network two-phase flow back model was first established based on the flow back data, and then the effective fracture network volume was inverted. And the reliability of the model was verified by cumulative production and microseismic data. The shale gas wells were taken as an example in the southern Sichuan Basin, to obtain the effective fracture network volume, which is useful for evaluating post-fracture effects. The flow back models are summarized in Table 1.

2. Model Development

2.1. The Model of Two-Phase Flow Back

In order to establish the flow back model, the assumptions were as follows:
(1)
The effective crack network of 1/2 single cluster is equivalent to the tree-shaped fractal fracture network shown in Figure 1.
(2)
Consider the reverse imbibition effect and the redistribution of original free gas in activated natural fractures.
(3)
Ignore the capillary pressure in the effective fracture network system and ignore the influence of gravity.
(4)
The matrix is only considered as the gas source. Material exchange occurs between the matrix system and the effective fracture network system through the cross-flow equation.
(5)
The effective fracture network system is an elastic porous medium, and assuming that its compressibility is much greater than that of the shale matrix, the tree-shaped fractal fracture permeability and effective fracture volume are pressure-dependent variables.
As shown in Figure 1, the flow rate of rectangular fractures at level k of 1/2 single-cluster tree-shaped fractal fracture network is described by Equation (1) according to the Hagen–Poiseuille equation:
Q k = W fk 3 h fk 12 μ Δ P k l k
where ΔPk is the pressure difference of the k-th level fracture, Pa; μ is fluid viscosity, Pa s; lk, Wfk, hfk are the length, width, and height of the k-th level branch fracture, m.
The length, width, and height of cracks at level k are obtained by Equation (2).
l k = l 0 R L k W fk = W f 0 R W k h fk = h f 0 R h k
where l0, Wf0, and hf0 are the length, width, and height of the initial tree-shaped fractal fractures, m; RL, RW, and Rh are ratio of length, width, and height of fractures, dimensionless.
From Equation (1), it can be known that the viscous resistance of fluid flow in a single fracture is
R k = 12 μ l k β W fk 3 h fk = 12 μ R L k l 0 β R W 3 k R h k W f 0 3 h f 0
where β is correction coefficient, which can be corrected according to the conductivity experiment in practical applications.
According to the principle of fluid pressure drop in parallel and in series [18], the parallel model was used to convert each branch into a single channel, and then the series model was applied to convert all single channels into a channel. In this way, the total viscous resistance of the network can be calculated, and the total flow resistance of the tree-shaped fractal fracture network can be expressed as Equation (4):
R = k = 0 m R k N k = R 0 1 R L n R W 3 R h m + 1 1 R L n R W 3 R h
where
R 0 = 12 μ l 0 β W f 0 3 h f 0
N k = n k
where n is the total number of branches fractal fracture, n = 2.
Then the flow rate is described by Equation (7):
Q = Δ P R = β W f 0 3 h f 0 1 R L n R W 3 R h Δ P 12 μ l 0 1 R L n R W 3 R h m + 1
where ΔP is the total pressure difference of fractal fracture network, MPa.
Δ P = P f P wf
As the fracturing fluid is produced, the fracture pressure decreases, and the fracture will be compressed under the closure stress. At this time, assuming that its height and length remain unchanged, and the width of Wfk at the k-th level has a width of Wfkc after compression, then the volume change is shown in Equation (9):
V fk V fkc = C f V fk Δ P f
where Cf is fracture compression coefficient; Vfk is the k-level single fracture volume under original fracture pressure; Vfkc is the k-level single fracture volume under current fracture pressure; and ΔPf is fracture system pressure drop, ΔPf = PfiPf.
V fk = W fk h fk l k
V fkc = W fkc h fk l k
where Wfkc is the width of the k-th level single fracture under the current fracture pressure.
Equations (10) and (11) are substituted into Equation (9), and then Equation (12) can be obtained.
W fkc = 1 C f Δ P f W fk
The flow rate of 1/2 single-cluster tree-shaped fractal fracture network is expressed in Equation (13) by considering the fracture closure effect.
Q = Δ P R = β W f 0 c 3 h f 0 1 R L n R W 3 R h Δ P 12 μ l 0 1 R L n R W 3 R h m + 1
where
W f 0 c = 1 C f Δ P f W f 0
The two-phase flow depends on the relationship between gas–water relative permeability and water saturation in the fractal fracture system. In the flow model of a single-phase tree-shaped fracture network, the relative permeability is taken into account. Equation (15) illustrates the flow rate calculation model for a 1/2 single-cluster tree-shaped fracture network in gas/water two-phase flow.
Q i = β K ri S w W f 0 3 h f 0 1 R L n R W 3 R h P f P wf 12 μ i B i l 0 1 R L n R W 3 R h m + 1
where Pf is the average pressure in the fracture system, which changes as fluid is produced, MPa; Pwf is the wellbore bottom flow pressure in a horizontal well, MPa; I is water or gas; Bi is volume coefficient, m3/m3; and Kri(Sw) is the relative permeability of gas/water, dimensionless.
In this paper, the linear relative permeability model was adopted:
K rw = S w
K rg = 1 S w
where Sw is the water saturation in the fracture.
Gas production superposition model is as follows:
Q g = j = 1 2 N f Q g j
Water production superposition model is as follows:
Q w = j = 1 2 N f Q w j
where Nf is the total clusters in horizontal well, which satisfies the following relationship:
N f = n f n CL
where nf is the number of operation stages, dimensionless; and nCL is the clusters number, dimensionless.

2.2. Calculation of Fracture Network Volume

The volume of the fractal fracture network is
V fi = 2 N f k = 0 m N k V k = 2 N f k = 0 m n k W fk h fk l k = 2 N f V 0 1 n R W R h R L m + 1 1 n R W R h R L
Recently, the initial effective fracture volume Vfi has been used internationally as an important parameter for preliminary evaluation of shale gas fracture network fracturing [19,20]. In the above formula,
V 0 = W f 0 h f 0 l 0
where Vfi is initial volume of effective fracture, m3; and V0 is initial fracture volume, m3.
It can be seen from Figure 1 that the longitudinal equivalent fracture half-length is
x f = l 0 1 + R L 1 R L m cos θ 1 R L
where θ represents the bifurcation fractures, which varies from 0 to π/2.
As shown in Figure 1, the lateral extension width is given by
w f = l 0 R L sin θ 1 R L m 1 R L
The effective fracture network volume can be obtained by Equation (25):
EFNV = N f w f x f h f 0 V b _ overlap
where xf denotes the longitudinal expansion degree and wf denotes the lateral expansion degree. It is an important parameter commonly used in quantitative evaluation of shale gas fracture network fracturing [21,22], where
w f = l 0 R L sin θ 1 R L m 1 R L
V b _ overlap = 2 S overlap h f 0 N f 1 w f > L w N f 0 w f L w N f
where Vb_overlap signifies the volume of overlap area within the EFNV.
In the matrix, the gas is entered the fracture network through the reverse imbibition of the fracturing fluid and the redistribution of the original free gas in the activated natural fractures during the process of fracturing and soaking. The reverse imbibition index (Iimb) is adopted to represent the combined effect. The primary origin of initial free gas in the effective fracture system is the reverse imbibition process [23]. The shale reverse imbibition index Iimb represents the proportion of free gas in the fracture network to the fracture volume, 0 ≤ Iimb ≤ 1. Obviously, there is the following relationship, the gas saturation of the original fracture is
S gi = I imb
The volume of free gas within the subsurface fractures is
V gfi = I imb V fi
Then the volume occupied by water within the fracture is
V wi = 1 I imb V fi
Water saturation of fractures is as follows:
S wi = V wi V fi = 1 I imb
where Sgi represents the original gas saturation within fractures, dimensionless; Vgfi represents the volume of free gas within the subsurface fractures, m3; Vwi is water volume in fractures, m3; and Swi is water saturation in fractures under initial conditions, dimensionless.

2.3. Material Balance Equation

The change of fracture volume is very complicated in the process of flow back. The pressure of the fracture system drops from the original fracture pressure Pfi to the current fracture pressure Pf, and the fracture pressure drop is ΔPf = PfiPf when a certain amount of fracturing fluid (Wp) is discharged from the fracture. The reduction of the fracture volume, the expansion of the free gas volume in the fracture, and the entry of matrix gas into the fracture will all reduce the volume of the fracturing fluid.
(1)
The reduction of fracture volume is calculated by
Δ V f = V fi C f Δ P f
  • (a)
    In fracture network, the expansion of free gas can be expressed as
Δ V gf = V gf V gfi = I imb V fi B gf B gfi 1
where Vgf is the gas volume in the current fracture, m3; Vgfi is the gas volume in the original fracture, m3; Bgf is the gas volume factor in the fracture, m3/m3; Bgfi is the original gas volume factor in the fracture, m3/m3; Pfi is original fracture pressure, MPa; and Pf is current fracture pressure, MPa.
  • (b)
    Intrusion of matrix gas
The ground volume of intrusion of matrix gas channeling into the tree-shaped fractal fracture network is obtained using Equation (34) by considering the compressibility of matrix pores and the desorption effect of adsorbed gas in shale matrix.
G mf = V b ϕ m B gmi V b ϕ m 1 C m Δ P m B gm + V b 1 ϕ m V Ei V E
where Bgm and Bgmi are matrix gas volume coefficient under the condition of current matrix pressure and original matrix pressure, m3/m3; ϕm is matrix porosity of shale reservoir, dimensionless; Cm is rock compressibility coefficient of shale matrix, 1/MPa; VEi and VE are unit shale adsorbed gas volume under original matrix pressure and current matrix pressure, m3/m3; and ΔPm is matrix pressure drop, MPa.
Vb is the effective fracture network volume and the calculation method is shown in Equation (25). The average pressure drop of the matrix system and the desorption of matrix adsorbed gas are calculated as follows:
Δ P m = P mi P m
V Ei = V L P mi P L + P mi
V E = V L P m P L + P m
The underground volume of matrix shale gas intrusion can be obtained as Equation (38) by combining it with the gas volume coefficient under the fracture pressure.
V mf = G mf B gf
where Pmi is shale original matrix pressure, MPa; Pm is shale current matrix pressure, MPa; VL is Langmuir volume, sm3/m3; PL is Langmuir pressure, MPa; and Gmf is matrix shale surface volume of intrusion of gas channeling into the tree-shaped fractal fracture network, m3.
  • (c)
    Produced free gas (underground volume)
Δ V gp = G p B gf
The gas storage capacity of fractures is obtained by combining Equations (33), (38), and (39).
Δ V g = Δ V gf + V mf Δ V gp
where Gp is cumulative gas production, m3.
(2)
The remaining volume of fracturing fluid in the fracture system
The volume of fracturing fluid in fractures will all be reduced due to the reduction of fracture pore volume, the expansion of fracturing fluid in fractures, and the gas storage capacity in fractures. Therefore, when the fracture original pressure Pfi decreases to Pf, the fracture fracturing fluid volume is
V w = V wi Δ V f Δ V g
Then, Equation (42) can be obtained by combining Equations (32), (40), and (41):
V w = V wi V fi C f Δ P I imb V fi B gf B gfi 1 G mf B gf + G p B gf
The remaining fracturing fluid volume is converted to surface conditions:
W res = V w B wf = V wi V fi C f Δ P I imb V fi B gf B gfi 1 G mf B gf + G p B gf / B wf
The basic form of the fracturing fluid material balance equation in the fracture is as follows:
V wi B wi = W p + W res
where Bwi is the fracturing fluid volume coefficient under original pressure, dimensionless; Wp is current cumulative fracturing fluid volume of flow back m3; and Wres is the remaining fracturing fluid volume, m3.
The remaining fracturing fluid volume Wres is inputted into the above formula:
V wi B wi = W p + V wi V fi C f Δ P I imb V fi B gf B gfi 1 G mf B gf + G p B gf / B wf
Then,
W p B wf + G p B gf = V wi B wf B wi B wi + V fi C f Δ P f + I imb V fi B gf B gfi 1 + G mf B gf
where Bwf is the fracturing fluid volume coefficient under the current fracture pressure, m3/m3; and Bwfi is the fracturing fluid volume coefficient under the original fracture pressure, m3/m3.
The compressibility coefficients of gas and water are
C g = B gf B gfi B gfi Δ P f C w = B wf B wfi B wfi Δ P f
Then,
W p B wf + G p B gf = V fi Δ P 1 I imb C w f + C f + I imb C g + G mf B gf
The function h representing the fracture pressure and matrix pressure at the k + 1th time step is defined as follows:
h P f k + 1 , P m k + 1 = V fi P fi P f k + 1 1 I imb C w f + C f + I imb C g + G mf B gf P f k + 1 W p B wf P f k + 1 G p B gf P f k + 1
where
B gf P f k + 1 = p sc Z sc T sc Z f T i P f k + 1
where Cwf is fracturing fluid compression coefficient in fractures, dimensionless; Psc is ground pressure, MPa; Zsc is gas deviation factor in ground conditions, dimensionless; Tsc is ground temperature, K; Zf is gas deviation factor in fractures, dimensionless; and Ti is original formation temperature of reservoir, K.
Water saturation of fractures under current formation conditions is as follows:
S w = V w V w + V gfi + Δ V g
The equation can be simplified to
S w = V w V fi Δ V f
The equations for fracture and matrix pressure at time k + 1 can be derived by setting Equation (49) equal to zero.
h P f k + 1 , P m k + 1 = 0
The material balance equation needs to be established of the shale matrix system to solve this equation due to there being two unknowns, Pfk+1 and Pmk+1.
If the fracture pressure (Pf) falls below the matrix gas breakthrough pressure (PBT), it will result in channeling between the matrix and fractures. The pressure PBT is related to the porosity and permeability of the matrix. The channeling equation [24] is described by Equation (54).
q m = k m α mf P m P f μ g
If PBTPfk, then the channeling from the matrix to the fracture system from time k to k + 1 is 0:
Δ G mf = 0
If PBT > Pfk, the channeling from the matrix to the fracture system from time k to k + 1 is
Δ G mf = Δ t V b ϕ m q m = Δ t V b ϕ m k m α mf P m k P f k + 1 μ g
where ΔGmf is the channeling of gas from the shale matrix to the fracture system from time k to time k + 1, m3; Δt—the time from time k to time k + 1, s; Vb is shale reservoir stimulation volume, m3; ϕm is matrix porosity of reservoir, dimensionless; qm is gas channeling flow from matrix to fracture system per unit time, s−1; km is matrix permeability, um2; αmf is matrix to fracture channeling factor [25], m−2; Pmk is average pressure of shale matrix system at time k, MPa; Pfk+1 is average pressure of fracture system at time k + 1, MPa; and ug is shale gas viscosity, Pa·s.
The channeling flow diffusion from the matrix to the fracture system from time k to k + 1 can also be expressed as
Δ G mf = G mf k + 1 G mf k
where
G mf k = V b ϕ m B gmi V b ϕ m 1 C m P mi P m k B gm k + V b 1 ϕ m V Ei V E P m k
G mf k + 1 = V b ϕ m B gmi V b ϕ m 1 C m P mi P m k + 1 B gm k + 1 + V b 1 ϕ m V Ei V E P m k + 1
V E = V L P P + P L
where VL is Langmuir volume, sm3/m3; PL—Langmuir pressure, MPa. Bgm and Bgmi are matrix gas volume coefficient under current matrix pressure and original matrix pressure, m3/m3; Cm is rock compressibility coefficient of shale matrix, 1/MPa; and VEi and VE are the unit shale adsorbed gas volume under original matrix pressure and current matrix pressure, m3/m3.
Equations (58) and (59) are substituted into Equation (57), then Equation (61) can be obtained.
Δ G mf = G mf k + 1 G mf k = V b ϕ m 1 C m P m P m k B gm k V b ϕ m 1 C m P m P m k + 1 B gm k + 1 + V b 1 ϕ m V E P m k V E P m k + 1
The function g representing the fracture and matrix pressure at the k + 1 time step is derived from Equation (62) by combining Equations (56) and (61):
g P m k + 1 , P f k + 1 = ϕ m 1 C m P mi P m k B gm k ϕ m 1 C m P mi P m k + 1 B gm k + 1 + 1 ϕ m V E P m k V E P m k + 1 Δ t ϕ m k m α mf P m k P f k + 1 μ g
The equations for the fracture network and matrix pressure at time k + 1 can be derived by setting Equation (62) to zero:
g P f k + 1 , P m k + 1 = 0
There are two unknowns, Pfk+1 and Pmk+1, in above equation, and the fracture pressure Pfk+1 and matrix pressure Pmk+1 at time k + 1 can be obtained by solving Equations (53) and (63), and the fracture pressure Pfk and matrix pressure Pmk at time k are known.

2.4. Model Solution

(1)
The tree-shaped fracture network structure parameters are preset, including l0, Wf0, hf0, RL, Rw, Rh, θ, m, n, and Cf, and initial fracture network and matrix average pressure (Pf(k = 1) = Pfi, Pm(k = 1) = Pmi) and the bottom hole flow pressure Pwf, then QW(k) and Qg(k) is calculated by Equation (15), k = 1, 2, …, Num;
(2)
The cumulative production of step k is solved, the cumulative flow back Wp = sum(Qw(k)), and the cumulative gas production Gp = sum(Qg(k));
(3)
If Pf(k) > =PBT, then there is no matrix gas channeling into the effective fracture network, Gmf = 0, and then Pm(k + 1) = Pm(k). Based on the dichotomy method, the fracture pressure Pf(k + 1) is calculated by combining the fracture material balance Equation (49), and the network average pressure Pf(k + 1) ranges from [0, Pfi]. When Pf(k) < PBT, the matrix gas channeling into fractures, the fracture pressure Pf(k + 1) range [0, Pfi], under the condition of known Pf(k + 1), the matrix pressure Pm(k + 1) range [0, Pmi]. The matrix pressure Pm(k + 1) is calculated in combination with the matrix material balance Equation (62), and then Gmf is calculated by using Equation (34), which is then brought into the Equation (49), and Pf(k + 1) is solved by the dichotomy method.
(4)
The Pf(k + 1) and Pm(k + 1) are assigned to Pf(k) and Pm(k), and steps (1), (2), and (3) are repeated until k = Num.

2.5. Flow Back Model Validation

The rate transient analysis (RTA) method was used to identify and calculate the flow back flow stage of early shale fracturing. Since a constant flow rate or bottom hole flow pressure was rarely used in the flow back, a small oil nozzle is usually used to control the discharge, step by step amplification, and adjust a stable drainage system in southern Sichuan. Therefore, the rate normalized pressure (RNP) and its derivative (RNP’) were used to observe the flow back characteristics. The equations are as follows [26]:
RNP = P fi P wf q w
RNP ' = dRNP dln t MB
t MB = W P q w
where Pfi is initial fracture network pressure, MPa; Pwf is bottom hole flow pressure, MPa; qw is water production rate, m3/h; tMB is equivalent time, h; and Wp is cumulative water production, m3.
The established shale gas fracturing two-phase flow back model was used to calculate the flow back data of 1000 h. The model parameter settings are shown in Table 2. As shown in Figure 2, in the initial state, the effective fracture network volume is completely filled with water (Iimb = 0) and the observed flow stages are consistent with field data [27] and analytical models [26] and Williams-Kovacs [28]. For the case of Iimb = 0, the two phases in RNP’ are characterized by a slope equal to 1. If the free gas is contained in the effective fracture network volume in the initial state, the more free gas content, the smaller the slope in stage 2 (for example, the slope of 0.70 at Iimb = 0.1 is greater than the slope of 0.65 at Iimb = 0.3), and the slope of stage 3 is always 1. Stage 2 and stage 3 last for a long time, and the characteristics of these two stages can be fully reflected by the flow back model proposed in this paper, so the established model is reasonable.
Some assumptions were made on the real seepage state of the fracturing fluid in order to be used in the inversion of the effective fracture network volume, so that the model will deviate from the real characteristics of the fluid flow back stage. But because the permeability of the effective fracture network is very large, the pressure wave propagates very fast, and this flow stage is usually impossible to observe in traditional low-frequency flow back data (transient flow in the fracture system may last only 10 s [29], and the absence of this stage has little error in the inversion). Stage 2 and stage 3 are the main stages in the production process of flow back, and the model proposed in this paper can fully reflect the characteristics of these two stages. It shows that the flow back model proposed in this paper is reasonable.

3. EFNV Inversion

In this section, a mathematical model was developed to invert the volume of an effective fracture network in shale gas using a tree-shaped fractal fracture flow back approach. Then, based on the flow back data of shale gas well, a set of genetic algorithm workflows suitable for the effective fracture network volume inversion of shale gas is developed and combined with a genetic algorithm. An example was taken from shale fracturing wells in southern Sichuan, where the inversion of effective fracture network volume was conducted for 162 wells. A basic data set of fracture network volume was obtained.

3.1. Statistical Evaluation Methods

The following four statistical descriptors are introduced to evaluate the reliability of the fitting between the flow back data (Qwsc and Qgsc) and in the predicted data of the model (Qw and Qg), yi is the observed value, fi is the predicted value, and N is the number of samples.
(1)
Coefficient of determination (R2)
The closer the R2 value is to 1, the better the predictive performance of the model.
R 2 = 1 i = 1 N y i f i 2 i = 1 N y i y ¯ 2
y ¯ = 1 N i = 1 N y i
(2)
Root mean square error (RMSE)
RMSE is the standard deviation of the error between the observed values (Qwsc and Qgsc) and the predicted values (Qw and Qg). The smaller the RMSE, the better the prediction performance of the model.
R M S E = 1 N i = 1 N y i f i 2
(3)
Slope of regression line (K)
The K value is the slope of the fitted linear regression equation between the observed values (Qwsc and Qgsc) and predicted values (Qw and Qg).
K = i = 1 N y i f i i = 1 N f i 2
(4)
Willmott’s index of agreement (IA)
The range of IA is [0,1]: IA = 1 indicates that the forecast is perfect, and IA = 0 indicates that the forecast is invalid.
I A = 1 i = 1 N y i f i 2 i = 1 N f i y ¯ + y i y ¯ 2
According to statistical recommendations, a reliable result can be evaluated by R2 > 0.64, 0.85 < K < 1.15, or IA > 0.80 [30].

3.2. Fitness Function

The goal was to find the most suitable fracture network structure through flow back data inversion, that is, the flow back water and gas production predicted by the model should be the same or very close to the actual monitored flow back data. Therefore, the fitness function was established based on the idea of maximizing the determination coefficient R2 between the predicted and the observed value. The fitness function is as follows:
F x = f x x f e a s i b l e r e g i o n f x + p e n a l t y x x f e a s i b l e r e g i o n
where
f x = i = 1 N Q wsci Q wpredi 2 i = 1 N Q wsci Q ¯ wsci 2 + i = 1 N Q gsci Q gpredi 2 i = 1 N Q gsci Q ¯ gsci 2
p e n a l t y x = { 10 8 p f   or   p m f e a s i b l e   r e g i o n 0 p f   or   p m f e a s i b l e   r e g i o n
The decision variables were as follows:
x = x 1 , x 2 , , x 12
where
x 1 = l 0 x 5 = R h x 9 = p fi x 2 = W f 0 x 6 = R L x 10 = α mf x 3 = m x 7 = θ x 11 = C f x 4 = R W x 8 = p BT x 12 = I imb
The optimization goal was as follows:
minimize F x
The variable upper and lower bounds were as follows:
L B i x i U B i i = 1 , 2 , , 12
The constraints were as follows:
V fi T I V
where
V fi = x 1 x 2 h f 0 1 n x 4 x 5 x 6 x 3 + 1 1 n x 4 x 5 x 6
where Penalty(x) is the penalty function, which is taken as 108 in this paper; Qwsci is the flow back fluid volume actually measured at the i-th moment, m3/min; Qwpredi is the flow back fluid volume calculated by the flow back model at the i-th moment, m3/min; Qgsci is the gas production volume actually measured at the i-th moment, m3/min; Qgpredi is the gas production volume calculated by the flow back model at the i-th time moment, m3/min; Q ¯ wsci is the average flow back fluid volume actually measured, m3/min; Q ¯ gsci is the average gas production volume actually measured by the shale gas well, m3/min; l0 is the initial fracture length, m; Wf0 is the initial fracture width, m is the number of fractures; dimensionality; RW is the fracture width ratio, dimensionless; Rh is the fracture height ratio, dimensionless; RL is the fracture length ratio, dimensionless; θ is the fracture bifurcation angle, rad; pBT is the matrix gas breakthrough pressure, MPa; pfi is the initial fracture pressure, MPa; αmf is the channeling factor, m−2; Cf is fracture compression coefficient, 10−9 Pa−1; Iimb is the shale reverse imbibition index, dimensionless; LBi is the lower limit of fitting parameters for flow back data inversion; UBi is the upper limit of fitting parameters for flow back data inversion; Vfi is the initial effective fracture volume, m3; TIV is the total injected fracturing fluid volume, m3; hf0 is the initial height of tree-shaped fractal fractures, m; and n is the bifurcation number of the tree-shaped fractal fracture, dimensionless.

3.3. Genetic Algorithm Workflow

The genetic algorithm (GA) was established by Holland [31] and inspired by natural selection in biological evolution. Three basic operators were used to generate high-quality solutions for optimization and search problems: selection, crossover, and mutation. The GA is used to optimize a given objective function (fitness function). The inversion work was conducted using shale gas well flow back data to minimize errors between the actual flow back data (Qwsc and Qgsc) and the predicted data from the fracture network model (Qw and Qg). Therefore, the fitness function included the optimization of the dual objectives of gas production and water production.
The genetic algorithm was used to predict the fracture network. The workflow for inversion of effective fracture network volume is shown in Figure 3. The inversion simulation process is as follows:
Step 1: The inversion simulation starts. Genetic algorithm parameters are shown in Table 3, and the upper and lower limits of the flow back model parameters are shown in Table 4;
Step 2: The basic parameters of shale gas wells are inputted, including flow back data, fracturing operation parameters, initial reservoir parameters, etc.;
Step 3: The bottom hole flowing pressure is calculated as the input value through the Beggs–Brill model [32];
Step 4: An initial total group of genetic algorithm fitting parameters is formed;
Step 5: The predicted values (Qw and Qg) are calculated through the flow back model;
Step 6: The fitness function value is calculated;
Step 7: Judging whether all the best fitting parameter values are found, or whether the termination condition is met;
Step 8: If the judgment result is “yes”, then the fitting parameters are stored; if the judgment result is “no”, then genetic algorithm selection operation, crossover operation, and mutation operation are performed to generate the next generation population of fitting parameters, and from step 5 one starts recounting.

3.4. Model Validation

3.4.1. Analysis of Production

The shale gas wells in southern Sichuan were evaluated by the established EFNV inversion method, and the EFNV data set of the block was obtained. The average EFNV was 1497.8 × 104 m3, the standard deviation was 517.0 × 104 m3. The relationship between EFNV and test production, 60-day cumulative gas production, 90-day cumulative gas production, 120-day cumulative gas production, and 20-year technically recoverable reserves were analyzed.
The relationships between the EFNV obtained from the inversion and the test production, 60-day cumulative gas production, 90-day cumulative gas production, 120-day cumulative gas production, and 20-year technically recoverable reserves are shown in Figure 4, Figure 5, Figure 6, Figure 7 and Figure 8. It can be seen that the correlation decreases in the order (test production (R = 0.684) > 60-day cumulative gas production (R = 0.576) > 90-day cumulative gas production (R = 0.566) > 120-day cumulative gas production (R = 0.541) > 20-year technically recoverable reserves (R = 0.371)). The reason is that the EFNV obtained based on the flow back in this paper was mainly obtained from the inversion of the data within the first few hundred hours to two thousand hours of the flow back production, so its correlation with the test production is the highest.
The longer the production time is, the correlation with cumulative gas production gradually decreases. This is due to the limitations of the flow back data. However, the positive correlation between the effective fracture network volume and test production, 60-day cumulative gas production, 90-day cumulative gas production, 120-day cumulative gas production, and 20-year technically recoverable reserves (EUR) shows that if a larger effective fracture network volume is obtained, the probability of obtaining a larger gas production is greater. From the perspective of fracturing engineering, which can be used as an important parameter to evaluate fracturing effects and production prediction, this proves the reliability of the inversion results.

3.4.2. Microseismic Monitoring

Some wells were monitored by microseismic monitoring during the fracturing process in southern Sichuan. The relationship between the stimulated reservoir volume (SRV) of 20 microseismic monitoring wells and the EFNV obtained by flow back inversion was statistically analyzed. As shown in Figure 9, the correlation coefficient between them was R = 0.705, which belongs to a moderate positive correlation (0.5 < R < 0.8). The more fractured the reservoir detected by microseismic monitoring, the larger the SRV, and the greater the probability that the corresponding effective fracture network volume will be larger, which further indicates that the EFNV inversion method is reliable based on flow back, and the obtained EFNV can be used as an important parameter for evaluating fracturing effects and optimizing fracturing operation parameters.

4. Field Application

We gathered flow back data post-fracturing to input into the established tree-shaped fractal fracture flow back mathematical model for EFNV inversion. Then, we utilized the genetic algorithm workflow to process the data and predict fracture network parameters for optimization. We applied the inversion model to shale gas wells in the field, such as those in southern Sichuan Basin, to calculate the EFNV and key parameters for optimizing fracturing operations.
The fracturing fluid flow back analysis model and the EFNV inversion workflow were used for a field application. The flow back production data of a shale gas well H3 in southern Sichuan (Figure 10) was taken as an example. The statistics of the basic parameters are shown in Table 5.
The hydraulic fracturing fracture network characteristics inversion was performed on the flow back production data of Well H3 based on the proposed flow back model and the flow back data inversion genetic algorithm. The water transient history fitting of well H3 is shown in Figure 11. Among them, R2(Qw) = 0.883, R2(Qg) = 0.927, IA(Qw) = 0.969, IA(Qg) = 0.982, K(Qw) = 1.03, and K(Qg) = 0.95; according to the statistical recommendation standard, the inversion of fracture characteristics of Well H3 was reliable. The root mean square error RMSE of the predicted value and the measured value of Qw and Qg were 1.54 m3/h and 0.14 × 104 m4/h, respectively. The changes in the bottom hole flow pressure, fracture network pressure, shale matrix pressure, and wellhead nozzle size during the flow back process are shown in Figure 12. It can be seen that the pressure drop in the fracture network is larger than that in the shale matrix after the well is opened. After the first well shut-in, the gas in the matrix channeled into the fractures due to the pressure difference between the fracture network system and the matrix system. In the network, the fracture network pressure gradually rises to be equal to the matrix pressure. After the second well opening, the fracture network pressure begins to decrease again until the second well shut-in, and the fracture pressure gradually rises again. The network pressure did not rise back to equal the matrix pressure due to the short shut-in time of the second well. It can be clearly seen from Figure 11 that the fitting effect of the flow back water data is average in the early part. It is speculated that the possible reason is that the flow back data records on site are not accurate due to frequent changes in the size of the choke nozzle, and the size of the choke nozzle changed four times before the first shut-in and changed nine times before the second shut-in (Figure 12).
The fracture characteristics of Well H3 obtained through the genetic algorithm workflow inversion are shown in Table 6, including tree-shaped fractal fracture network structure parameters, fracture system original pressure, matrix gas breakthrough pressure, fracture compressibility coefficient, reverse imbibition index, effective fracture network volume, etc. Two types of hydraulic fracture networks will be formed after fracture network fracturing. Some fractures contribute to flow conduction and others are isolated due to fracture closure and for other reasons. The fracture network with a conductivity contribution is a highly conductive fracture network that allows matrix shale gas to enter and flow into the wellbore along it; these are also called effective fractures. A law of harmonic decline (HD) from the relationship between water production and cumulative water production can be observed by Lee [33], which can estimate the effective fracture volume (Table 7 and Figure 13). It can be seen from Figure 13 that when the HD model was applied to Well H3, the previous water flow data were overestimated, so the effective initial fracture volume would be overestimated.
Alkouh [2] observed a linear relationship between flow normalized pressure (RNP) and mass balance time (tMB) during the flow back process, and the effective initial fracture volume could be obtained (Table 7 and Figure 14). However, the fracture compressibility coefficient was ignored in the Alkouh model when the total compressibility coefficient was calculated, so the effective fracture volume would be overestimated.
As shown in Figure 15, the effective initial fracture volume calculated by the method proposed in this paper, in which the gas-water two-phase flow and fracture compressibility were taken into account, and the calculation results were of the same order of magnitude as those of the HD model and the Alkouh model, and were smaller than those of the two models. It shows that the inversion model of fracture network characteristics based on flow back is reasonable, and the applicability of the model is stronger. Shale gas fracturing generally has the problem of a low flow back rate [34], and the final estimated flow back rate of Well H3 was 20.8%. Most of the fracturing fluid is locked in the isolated ineffective fracture network with no conductivity, and it is difficult to flow back. This is the main reason for the low flow back rate.
The inversion effective fracture network volume fitting curves of some shale gas fracture network fracturing wells in southern Sichuan are shown in Figure 16. The model established can fit the actual field data very well in this paper, and it can be used to invert the EFNV efficiently after shale gas reservoir stimulation in southern Sichuan.

5. Conclusions

The inversion method of shale gas effective fracture network volume based on flow back data was established. Several effective fracture network volume inversion cases were simulated. The main conclusions are summarized below:
(1)
A two-phase flow back model of the tree-shaped fractal fracture network was established based on the flow back data, and the reliability of the model was verified by combining the RNP and RNP’. The model is reasonable, reliable, and has strong field adaptability.
(2)
Based on the principle of the genetic algorithm (GA), the GA workflow of EFNV inversion was established and combined with the two-phase flow back model, and then the EFNV inversion was carried out for shale gas wells in southern Sichuan. The results showed that the inversion of EFNV had a positive correlation with production and microseismic monitoring results, which confirmed the reliability of the inversion results.
(3)
The calculation showed that the results were in the same order of magnitude as those of the HD model and the Alkouh model, and were smaller than the calculation results of the two models, indicating that the fracture network inversion model proposed is reasonable and more applicable.
(4)
The fracture network feature inversion method can obtain tree-shaped fractal fracture network structure parameters and effective fracture network volume. It has important theoretical guiding significance and field application value for improving the calculation accuracy of EFNV and fracturing effect evaluation.
(5)
The model does not consider the effects of interference between wells and well channeling from adjacent wells, and accordingly the prediction of water production after well strings is relatively weak. It is recommended that the model can further consider the impact of simultaneous flow back of multiple wells and layers and the salinity of flow back water in the future.

Author Contributions

D.T.: conceptualization, methodology; J.W.: writing—review and editing; J.Z.: writing—original draft preparation; B.Z.: data collection; Y.S.: visualization, investigation; C.S.: formal analysis; L.R.: visualization; Y.H.: original draft preparation; Z.W.: software. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Data Availability Statement

The original contributions presented in the study are included in the article, further inquiries can be directed to the corresponding author.

Acknowledgments

Thanks to all the authors for their contribution to the paper and their efforts.

Conflicts of Interest

Authors Dengji Tang, Jianfa Wu, Bo Zeng, Yi Song, Cheng Shen, and Yongzhi Huang are employed by the CNPC Southwest Oil and Gas Field Company; The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

References

  1. Zhao, J.; Ren, L.; Jiang, T.; Hu, D.; Wu, L.; Wu, J.; Yin, C.; Li, Y.; Hu, Y.; Lin, R.; et al. Ten years of gas shale fracturing in China: Review and prospect. Nat. Gas Ind. B 2022, 9, 158–175. [Google Scholar] [CrossRef]
  2. Alkouh, A.; McKetta, S.; Wattenbarger, R.A. Estimation of effective-fracture volume using water-flowback and production data for shale-gas wells. J. Can. Pet. Technol. 2014, 53, 290–303. [Google Scholar] [CrossRef]
  3. Ren, L.; Wang, Z.; Zhao, J.; Lin, R.; Wu, J.; Song, Y.; Tang, D. Shale gas effective fracture network volume prediction and analysis based on flow back data: A case study of southern Sichuan Basin shale. Geoenergy Sci. Eng. 2023, 228, 211963. [Google Scholar] [CrossRef]
  4. Umar, I.A.; Negash, B.M.; Quainoo, A.K.; Ayoub, M.A. An outlook into recent advances on estimation of effective stimulated reservoir volume. J. Nat. Gas Sci. Eng. 2021, 88, 103822. [Google Scholar] [CrossRef]
  5. Qu, H.; Zhou, F.; Peng, Y.; Pan, Z. ESRV and production optimization for the naturally fractured keshen tight gas reservoir. In Proceedings of the International Field Exploration and Development Conference 2017; Springer: Singapore, 2019; pp. 1576–1595. [Google Scholar]
  6. Hussain, M.; Saad, B.; Negara, A.; Sun, S. Understanding the True Stimulated Reservoir Volume in Shale Reservoirs. In Proceedings of the SPE Kingdom of Saudi Arabia Annual Technical Symposium and Exhibition, Dammam, Saudi Arabia, 24–27 April 2017; OnePetro: Richardson, TX, USA, 2017. [Google Scholar]
  7. Adefidipe, O.A.; Xu, Y.; Dehghanpour, H.; Virues, C.J. Estimating effective fracture volume from early-time production data: A material balance approach. In Proceedings of the SPE/CSUR Unconventional Resources Conference–Canada, Calgary, AB, Canada, 30 September–2 October 2014; OnePetro: Richardson, TX, USA, 2014. [Google Scholar]
  8. Crafton, J.W. Oil and gas well evaluation using the reciprocal productivity index method. In Proceedings of the Oklahoma City Oil and Gas Symposium/Production and Operations Symposium, Oklahoma City, Oklahoma, 9–11 March 1997; SPE: San Antonio, TX, USA, 1997; p. 37409. [Google Scholar]
  9. Abbasi, M.A.; Ezulike, D.O.; Dehghanpour, H.; Hawkes, R.V. A comparative study of flowback rate and pressure transient behavior in multifractured horizontal wells completed in tight gas and oil reservoirs. J. Nat. Gas Sci. Eng. 2014, 17, 82–93. [Google Scholar] [CrossRef]
  10. Williams-kovacs, J.; Clarkson, C.R. Stochastic modeling of two-phase flowback of multi-fractured horizontal wells to estimate hydraulic fracture properties and forecast production. In Proceedings of the SPE Unconventional Resources Conference/Gas Technology Symposium, The Woodlands, TX, USA, 10–12 April 2013; SPE: San Antonio, TX, USA, 2013; p. 164550. [Google Scholar]
  11. Clarkson, C.; Qanbari, F.; Williams-Kovacs, J. Innovative use of rate-transient analysis methods to obtain hydraulic-fracture properties for low-permeability reservoirs exhibiting multiphase flow. Lead. Edge 2014, 33, 1108–1122. [Google Scholar] [CrossRef]
  12. Ezulike, D.O.; Dehghanpour, H.; Virues, C.J.; Hawkes, R.V.; Jones, R.S., Jr. Flowback fracture closure: A key factor for estimating effective pore volume. SPE Reserv. Eval. Eng. 2016, 19, 567–582. [Google Scholar] [CrossRef]
  13. Fu, Y.; Ezulike, D.O.; Dehghanpour, H.; Steven Jones, R. Estimating effective fracture pore-volume from early single-phase flowback data and relating it to fracture design parameters. In Proceedings of the SPE/CSUR Unconventional Resources Conference, Calgary, AB, Canada, 20–22 October 2015; OnePetro: Richardson, TX, USA, 2015. [Google Scholar]
  14. Xu, Y.; Ezulike, O.; Dehghanpour, H. Estimating compressibility of complex fracture networks in unconventional reservoirs. Int. J. Rock Mech. Min. Sci. 2020, 127, 104186. [Google Scholar] [CrossRef]
  15. Xu, Y.; Dehghanpour, H.; Ezulike, O.; Virues, C. Effectiveness and time variation of induced fracture volume: Lessons from water flowback analysis. Fuel 2017, 210, 844–858. [Google Scholar] [CrossRef]
  16. Jia, C.; Huang, T.; Yao, J.; Xing, H.; Zhang, H. Effect of isolated fracture on the carbonate acidizing process. Front. Earth Sci. 2021, 9, 698086. [Google Scholar] [CrossRef]
  17. Yao, J.; Ding, Y.; Sun, H.; Fan, D.; Wang, M.; Jia, C. Productivity Analysis of Fractured Horizontal Wells in Tight Gas Reservoirs Using a Gas–Water Two-Phase Flow Model with Consideration of a Threshold Pressure Gradient. Energy Fuels 2023, 37, 8190–8198. [Google Scholar] [CrossRef]
  18. Xu, P.; Yu, B.; Feng, Y.; Liu, Y. Analysis of permeability for the fractal-like tree network by parallel and series models. Phys. A Stat. Mech. Its Appl. 2006, 369, 884–894. [Google Scholar] [CrossRef]
  19. Fu, Y.; Dehghanpour, H.; Ezulike, D.O.; Jones, R.S., Jr. Estimating effective fracture pore volume from flowback data and evaluating its relationship to design parameters of multistage-fracture completion. SPE Prod. Oper. 2017, 32, 423–439. [Google Scholar] [CrossRef]
  20. Moussa, T.; Dehghanpour, H.; Fu, Y.; Ezulike, O. Dynamic fracture volume estimation using flowback data analysis and its correlation to completion-design parameters. In Proceedings of the SPE Hydraulic Fracturing Technology Conference and Exhibition, The Woodlands, TX, USA, 5–7 February 2019; SPE: San Antonio, TX, USA, 2019; p. 194322. [Google Scholar]
  21. Ren, L.; Lin, R.; Zhao, J.; Wu, L. Cluster spacing optimal design for staged fracturing in horizontal shale gas wells based on optimal SRV. Nat. Gas Ind. 2017, 37, 69–79. [Google Scholar] [CrossRef]
  22. Lin, R. Dynamic Simulation of Stimulated Reservoir Volume (SRV) during Hydraulic Fracturing in Horizontal Shale Gas Well; Southwest Petroleum University: Chengdu, China, 2018. [Google Scholar]
  23. Ezulike, O.D.; Dehghanpour, H. Modelling flowback as a transient two-phase depletion process. J. Nat. Gas Sci. Eng. 2014, 19, 258–278. [Google Scholar] [CrossRef]
  24. Wang, F.; Li, B.; Zhang, Y.; Zhang, S. Coupled thermo-hydro-mechanical-chemical modeling of water leak-off process during hydraulic fracturing in shale gas reservoirs. Energies 2017, 10, 1960. [Google Scholar] [CrossRef]
  25. Bian, X.; Zhang, S.; Zhang, J.; Wang, F. A new method to optimize the fracture geometry of a frac-packed well in unconsolidated sandstone heavy oil reservoirs. Sci. China Technol. Sci. 2012, 55, 1725–1731. [Google Scholar] [CrossRef]
  26. Song, B.; Ehlig-Economides, C.A. Rate-normalized pressure analysis for determination of shale gas well performance. In Proceedings of the SPE Unconventional Resources Conference/Gas Technology Symposium, The Woodlands, TX, USA, 14–16 June 2011; SPE: San Antonio, TX, USA, 2011; p. 144031. [Google Scholar]
  27. Abbasi, M.A. A Comparative Study of Flowback Rate and Pressure Transient Behaviour in Multifractured Horizontal Wells; University of Alberta: Edmonton, AB, Canada, 2013. [Google Scholar]
  28. Williams-Kovacs, J.D. Quantitative Analysis of Multi-Phase Flowback from Multi-Fractured Horizontal Wells. Doctoral Thesis, University of Calgary, Calgary, AB, Canada, 2017. [Google Scholar]
  29. Wylie, E.B.; Streeter, V.L. Fluid Transients; McGraw-Hill International Book Co.: New York, NY, USA, 1978. [Google Scholar]
  30. Roy, P.P.; Roy, K. On some aspects of variable selection for partial least squares regression models. QSAR Comb. Sci. 2008, 27, 302–313. [Google Scholar] [CrossRef]
  31. Holland, J.H. Adaptation in Natural and Artificial Systems: An Introductory Analysis with Applications to Biology, Control, and Artificial Intelligence; MIT Press: Cambridge, MA, USA, 1992. [Google Scholar]
  32. Beggs, D.H.; Brill, J.P. A study of two-phase flow in inclined pipes. J. Pet. Technol. 1973, 25, 607–617. [Google Scholar] [CrossRef]
  33. Lee, W.J.; Wattenbarger, R.A. Gas Reservoir Engineering; Henry, L., Ed.; Doherty Memorial Fund of AIME; Society of Petroleum Engineers: San Antonio, TX, USA, 1996. [Google Scholar]
  34. Singh, H. A critical review of water uptake by shales. J. Nat. Gas Sci. Eng. 2016, 34, 751–766. [Google Scholar] [CrossRef]
Figure 1. Tree-shaped fractal network of 1/2 single-cluster.
Figure 1. Tree-shaped fractal network of 1/2 single-cluster.
Processes 12 01027 g001
Figure 2. The relationship between RNP’, RNP, and tMB for flow back of shale gas wells.
Figure 2. The relationship between RNP’, RNP, and tMB for flow back of shale gas wells.
Processes 12 01027 g002
Figure 3. Genetic algorithm workflow for inversion of effective fracture network volume.
Figure 3. Genetic algorithm workflow for inversion of effective fracture network volume.
Processes 12 01027 g003
Figure 4. Relationship between effective fracture network volume and test production.
Figure 4. Relationship between effective fracture network volume and test production.
Processes 12 01027 g004
Figure 5. Relationship between effective fracture network volume and 60-day cumulative gas production.
Figure 5. Relationship between effective fracture network volume and 60-day cumulative gas production.
Processes 12 01027 g005
Figure 6. Relationship between effective fracture network volume and 90-day cumulative gas production.
Figure 6. Relationship between effective fracture network volume and 90-day cumulative gas production.
Processes 12 01027 g006
Figure 7. Relationship between effective fracture network volume and 120-day cumulative gas production.
Figure 7. Relationship between effective fracture network volume and 120-day cumulative gas production.
Processes 12 01027 g007
Figure 8. Relationship between effective fracture network volume and 20-year technically recoverable reserves.
Figure 8. Relationship between effective fracture network volume and 20-year technically recoverable reserves.
Processes 12 01027 g008
Figure 9. Relationship between inversion EFNV and microseismic monitoring SRV.
Figure 9. Relationship between inversion EFNV and microseismic monitoring SRV.
Processes 12 01027 g009
Figure 10. Flow back production data of well H3.
Figure 10. Flow back production data of well H3.
Processes 12 01027 g010
Figure 11. The transient performance of water and gas fitted of well H3.
Figure 11. The transient performance of water and gas fitted of well H3.
Processes 12 01027 g011
Figure 12. Transient performance of pressure and choke size of well H3.
Figure 12. Transient performance of pressure and choke size of well H3.
Processes 12 01027 g012
Figure 13. HD model.
Figure 13. HD model.
Processes 12 01027 g013
Figure 14. Alkouh model [2].
Figure 14. Alkouh model [2].
Processes 12 01027 g014
Figure 15. Comparison of initial effective fracture network volume estimated by different models.
Figure 15. Comparison of initial effective fracture network volume estimated by different models.
Processes 12 01027 g015
Figure 16. Effective fracture network volume inversion fitting curve.
Figure 16. Effective fracture network volume inversion fitting curve.
Processes 12 01027 g016aProcesses 12 01027 g016b
Table 1. Summary of current flow back models.
Table 1. Summary of current flow back models.
AuthorYearModelWork SummaryInputsOutputs
Williams-Kovacs and Clarkson [10]2013Two-phase flow back of multi-fractured horizontal wells stochastic mode.Using high-frequency fluid production and flowing pressures (hourly or greater) to characterize hydraulic fracture or reservoir parameters.Fractured stages, perforation clusters, et al.Bulk permeability, effective fracture half-length, and production.
Alkouh et al. [2]2014Water flow back and long-term production simulated mode.Procedures and examples are presented, including water flow back and water-production data, in the analysis of shale-gas wells using rate transient analysis.Initial pressure, fracture porosity, etc.Gas production and water production.
Fu et al. [13]2015A flowback model.How flow back data can be interpreted to estimate effective fracture pore-volume and its relationship to fracture design parameters are illustrates.Flow back rate, pressure date et al.Effective fracture pore volume.
Ezulike et al. [12]2016A two-phase tank model for reducing parameter uncertainty.The contributions of various drive mechanisms during flow back (fracture closure, gas expansion, and water depletion) are investigated.Well length, total injected volume et al.Fracture pore volume (PV), half-length, and permeability.
Xu et al. [15]2017An open-tank modelA comparative analysis to investigate the time variation of effective fracture volume during water flow back was conducted by this model.Initial water saturation, rock compressibility et al.Fracture–matrix interface area.
Hussain et al. [6]2017A fully 3-D simulation approach to estimate the SRV.The real-time changes in the reservoir’s geomechanics as a function of fluid pressures are considered.Fluid density, number of perforations et al.SRV.
Qu et al. [5]2019Predict production model.An integrated multiphysical model was developed to optimize the ESRV and predict the production.Matrix permeability, bottom hole temperature, etc.Production.
Xu et al. [14]2020Fracture compressibility modelThe effects of rock and proppant parameters on fracture compressibility are investigated.Fracture conductivity, etc.The compressibility of fracture networks.
Umar et al. [4]2021A simple analytical model to compare the pressure/rate transient behavior of the three flow back cases.The physics of flow back are understood by constructing basic diagnostic plots using two-phase flow back data from three multi-fractured horizontal wells.Early-time flow back pressure and rate data, etc.Single phase water rate and pressure.
Lan Ren et al. [3]2023Load recovery modeA genetic expression programming is established to calculate flow back.Geological engineering parameters.Load recovery equation.
Table 2. Model calculation parameters.
Table 2. Model calculation parameters.
ParameterUnitValueParameterUnitValue
Pmi106 Pa50TiK365.15
Pfi106 Pa60PBT106 Pa50
IimbDimensionless0, 0.1, 0.3Cf10−9 Pa−120
Cm10−9 Pa−10.1Cw10−9 Pa−10.4
μw10−3 Pa·s0.88μg10−3 Pa·s0.02
TpcK201Ppc106 Pa4.5
PL106 Pa2.95VLsm3/m30.3
Lwm1200hfm50
nfDimensionless20nCLDimensionless3
ϕmDimensionless0.04Km10−9 μm2200
BwDimensionless1.05αmfm−2200
Wf0m0.009l0m30
RWDimensionless0.7RhDimensionless0.7
RLDimensionless0.7θradπ/6
mDimensionless20as100,000
Table 3. Parameter settings of genetic algorithm for fitting flow back model parameters.
Table 3. Parameter settings of genetic algorithm for fitting flow back model parameters.
Genetic Algorithm Parameters Value
Initial parametersPopulation size500
Maximum algebra1200
Population initialization methodRandom method
Fitness functionEquation (72)
Genetic selectionRoulette
Termination conditionMaximum algebra
Genetic operatorCrossover probability85%
Mutation probability15%
Table 4. Upper and lower limits of fitting parameters for flow back data inversion.
Table 4. Upper and lower limits of fitting parameters for flow back data inversion.
ParameterUnitLowerUpper
l0m1300
Wf0m0.0000010.02
mDimensionless1500
RWDimensionless0.000011
RhDimensionless0.000011
RLDimensionless0.000011
θrad0π/2
pBTMPa0pmi
pfiMPa02pmi
αmfm−201000
Cf10−9 Pa−10.11000
IimbDimensionless01
Table 5. Basic parameters of shale gas well H3.
Table 5. Basic parameters of shale gas well H3.
ParameterUnitValueParameterUnitValue
Pmi106 Pa58.79TiK367.8
TIVm351,062.1dcasingm0.1143
Lwm1500hf0m43
Cm10−9 Pa−10.3Cw10−9 Pa−10.46
μw10−3 Pa·s0.28μg10−3 Pa·s0.042
TpcK190Ppc106 Pa4.61
PL106 Pa2.98VLsm3/m30.29
nfDimensionless28nCLDimensionless3
ϕmDimensionless0.054Km10−9 μm2380
BwDimensionless1.05nDimensionless2
Table 6. Inversion results of genetic algorithm workflow.
Table 6. Inversion results of genetic algorithm workflow.
ParameterUnitValueParameterUnitValue
l0m50.9Wf0m0.0154
mDimensionless24RLDimensionless0.618
RWDimensionless0.614RhDimensionless0.616
θrad1.07PBT106 Pa58.79
Pfi106 Pa58.81αmfm−2575
Cf10−9 Pa−1179IimbDimensionless0
EFNV104 m31172xfm90.85
Table 7. Results comparison of different models for estimating the initial EIFV.
Table 7. Results comparison of different models for estimating the initial EIFV.
ModelUnitValue
HD modelm321,192
Alkouh modelm315,895
Established modelm310,643
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Tang, D.; Wu, J.; Zhao, J.; Zeng, B.; Song, Y.; Shen, C.; Ren, L.; Huang, Y.; Wang, Z. The Inversion Method of Shale Gas Effective Fracture Network Volume Based on Flow Back Data—A Case Study of Southern Sichuan Basin Shale. Processes 2024, 12, 1027. https://doi.org/10.3390/pr12051027

AMA Style

Tang D, Wu J, Zhao J, Zeng B, Song Y, Shen C, Ren L, Huang Y, Wang Z. The Inversion Method of Shale Gas Effective Fracture Network Volume Based on Flow Back Data—A Case Study of Southern Sichuan Basin Shale. Processes. 2024; 12(5):1027. https://doi.org/10.3390/pr12051027

Chicago/Turabian Style

Tang, Dengji, Jianfa Wu, Jinzhou Zhao, Bo Zeng, Yi Song, Cheng Shen, Lan Ren, Yongzhi Huang, and Zhenhua Wang. 2024. "The Inversion Method of Shale Gas Effective Fracture Network Volume Based on Flow Back Data—A Case Study of Southern Sichuan Basin Shale" Processes 12, no. 5: 1027. https://doi.org/10.3390/pr12051027

APA Style

Tang, D., Wu, J., Zhao, J., Zeng, B., Song, Y., Shen, C., Ren, L., Huang, Y., & Wang, Z. (2024). The Inversion Method of Shale Gas Effective Fracture Network Volume Based on Flow Back Data—A Case Study of Southern Sichuan Basin Shale. Processes, 12(5), 1027. https://doi.org/10.3390/pr12051027

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop