# A Comprehensive Model for Estimating Stimulated Reservoir Volume Based on Flowback Data in Shale Gas Reservoirs

Authors: Qi Chen; Shaojun Wang; Dan Zhu; Guoxuan Ren; Yuan Zhang; Jinghong Hu
Journal: Geofluids
 (2020)

---

## Abstract

Stimulated reservoir volume (SRV) which is generated by horizontal drilling with multistage hydraulic fracturing governs the production in the shale gas reservoirs. Although microseismic data has been used to estimate the SRV, it is high-priced and sometimes overestimated. Additionally, the effect of stress sensitivity on SRV is not considered in abnormal overpressure areas. Thus, the objective of this work is to characterize subsurface fracture networks with stress sensitivity of permeability through the shale gas well production data of the early flowback stage. The flowback regions are first identified with the flowback data of two shale gas wells in South China. Then, we measured the permeability stress sensitivity of the core after fracturing, coupled to the dynamic relative permeability (DRP) calculation to obtain an accurate and simple DRP curve. After that, a comprehensive model is built considering dynamic two-phase relative permeability function and stress sensitivity. Finally, we compared the calculated results with the microseismic data. The results show that the proposed model could reasonably predict the SRV using the flowback data after fracturing. Additionally, compared with the microseismic data, the stress sensitivity should be included, especially in the abnormal overpressure block. It is believed that this mathematical model is accurate and useful. The work provides an efficient approach to estimate stimulated reservoir volume in the shale gas reservoirs.


---

## Body

## 1. Introduction

The development of shale gas has gained increasing attention with the decline of the production from conventional reservoirs [1, 2]. Due to the ultralow permeability of shales (from 10−23to10−17m2) [3], a combination of horizontal drilling and multistage hydraulic fracturing has been widely used to SRV, which increases effective contact significantly to improve gas production from an individual well. Therefore, characterization of the stimulated reservoir volume appears to be of vital importance to manage and predict shale reservoir performance [4, 5].While microseismic data has been extensively used to validate the stimulated reservoir volume, the SRV obtained directly from the microseismic data is generally overestimated. For example, SRV would include large unstimulated void regions by the convexity assumption [6–9]. Also, implementation of microseismic techniques in a field sometimes is not cost-effective, which limits a wide application in fields. Therefore, using early production data after hydraulic fracturing together with early flowback data of hydraulic fluids is fast becoming a key instrument in evaluating SRV [10–13], namely, transient analysis of flowback data.Two-phase flow appears to be a classic and representative flow behavior during the flowback process after hydraulic fracturing [2, 14]. For example, a large number of production data obtained from wells in the Horn River shale show that two-phase flow was observed for every single well during the flowback process after hydraulic fracturing [15, 16]. Clarkson [17] and Williams-Kovacs and Clarkson [18] also discovered two-phase flow within 48 hours of well opening and production in the Barnett and Marcellus shales. Similarly, the two-phase flow stage occurs in the Silurian Longmaxi Formation shale gas wells in South China after the well opening [19]. And an even greater advance in the research of production data patterns is that after processing the production data from shale gas wells in the Horn River Basin, Adefidipe et al. [20, 21] and Xu et al. [22, 23] divided the production data into two stages as early gas production (EGP) and late gas production (LGP) for the “V” shape of the relationship between the gas-water ratio (GWR) and the cumulative gas production. They believed that the EGP stage was special and worth studying because the LGP stage already had a matrix gas supply.To take advantage of the two-phase flow characteristics during the early gas production (EGP) stage, research studies have endeavored to estimate the SRV by coupling rock physics with various types of reservoir models. For example, Ezulike et al. [24–26] integrated an explicitly determined dynamic relative permeability (DRP) function in linear dual-porosity models to explore the two-phase flowback data. The function extends the single-phase water flow equation to gas-water two-phase flow which was calculated using cumulative production data of water and hydrocarbon compounds. Clarkson et al. [17, 27, 28] initially used an improved method for coalbed methane (CBM) flow to analyze multiphase flows. They assumed that the flow of gas and water through the generated fracture network was similar to the simultaneous flow of gas and water during long-term production from fractured coal reservoirs. Later, they proposed a large number of two-phase flow models based on this assumption to quantitatively calculate multiflow data, including tight oil and tight gas, to help analyze the flow characteristics of shale gas wells which are analogous to analyzing two-phase (gas and water) CBM production. And in the follow-up research work [29–32], they verified the analysis results from the above methods by integrating various techniques, including a modified material balance equation (MBE) applied to before-breakthrough (BBT) (before gas breakthrough) single-phase rate-transient analysis (RTA) and after-breakthrough (ABT) (after gas breakthrough) multiphase rate-transient analysis (RTA).While the aforementioned models provide insights in calculating SRV through multiflow data quantitatively, challenges of applying these models to the field remain because the Langmuir volume needs to be increased in the calculation process, which will lead to some extent uncertainty and inaccuracy. To overcome the complexity of applying these models in fields, Clarkson et al. [17, 28] and Fu et al. [33] established a linear matching curve relation to describe the straight line relationship between rate-normalized pressure (RNP) and material balance time (MBT). This is because fracture permeability and fracture volume can be obtained on the basis of the characteristics of the fitting curve. Also, Alkouh et al. [34] modified a log-log plot of RNP and RNP derivative against MBT for both the oil and gas cases and provided a straightforward approach to calculate the key fracture parameters covering fracture permeability and fracture half-length. However, the influence of fracture compressibility during the flowback process is not taken into consideration in the models, which likely leads to estimation error of the total compressibility.However, due to some defects in the hypothesis of DRP function, a nonlinear phenomenon occurs when the model is applied to some field cases. And the phenomenon of the nonzero intercept leads researchers to question the physical properties of these problems. It may be due to the inaccuracies caused by stress sensitivity in some areas. In particular, there is abnormal overpressure in shale gas reservoirs in southern China [35], so the influence of stress sensitivity, especially for the two-phase relative permeability, has to be taken into account [36]. Moreover, all the researchers did not make a production data verification and comparison with the microseismic data.We thus aimed to develop a stress-dependent two-phase relative permeability, which can be used to extend the existing single-phase model. Also, we proposed a new multiphase model to estimate SRV and fracture characteristics. Moreover, to verify our model, we compared the results with the microseismic data.In this work, we first processed the production data of shale gas wells and obtain the gas-water ratio characteristics. Then, we calculated the simple and practical two-phase dynamic relative permeability (DRP) based on the production data collected in the field belonging to shale gas wells in southern Sichuan. Subsequently, the stress-strain curves of permeability are obtained by using the laboratory experiments. Finally, the material balance equation and diffusion equation are used to process the production data to derive the model, thus calculating fracture parameters. We compared the calculation results with the microseismic data to verify the rationality and accuracy of the calculation results.


## 2. Data Preparation

### 2.1. Field Data from a Typical Shale Gas Reservoir in Southern China

In this work, we selected the flowback data from two production wells located in southern China: the shale gas formations of the Lower Silurianare, a typical marine shale gas reservoir. A large number of micron-nanoscale pores are developed in the mineral grains and organic matter [37]. These pores are the main channels for shale gas storage and migration, which play an important role in improving shale gas storage performance and providing good storage space for marine shale gas accumulation and enrichment [38]. The depth of the shale gas reservoir is over 2000 m at subsurface with a formation pressure coefficient up to 2.0. The thickness of the entire shale reservoir is about 100 m with porosity ranging from 2% to 4% and matrix permeability of 0.001 mD. The production at the early stage was about 10~20×104m3/d with a stable production around at 5~8×104m3/d.The east-west anticline structure belt with few faults is arranged in the left echelon row as the main geological structure of this area [39]. This drilling platform is located at the end of the anticline structure. Although two compresso-shear faults exist around the well pad, the overall upper lateral sealing performance leads to an effective accumulation and preservation of shale gas [37, 38]. Figure 1 shows the layout of the horizontal well group. Due to the low quartz content in the high-quality shale segment at the bottom of the Longmaxi Formation, natural fractures are not developed throughout the reservoir [40]. In the whole horizontal well group, Well 1, Well 4, and Well 5 are located in one formation, while Well 2, Well 3, Well 6, and Well 7 belong to another formation.Figure 1
Well pad schematic. Seven MFHWs drilled in the southern Sichuan Basin.In this work, we used data from Well 3 and Well 4 to test our model for the following two reasons. Firstly, the two wells were drilled through the shale gas reservoir in different directions and are located in two different formations. Secondly, microseismic monitoring is conducted in the two wells together with a comprehensive drilling, completion data, and well test data. Figure2 shows gas and water production rates of the two wells. The two wells both experienced gas breakthrough immediately after well opening and production, meaning that there was no single-phase flow. The difference is that Well 3 experienced a precipitous drop in production at one point (80 hrs), while Well 4 experienced a more gentle fluctuation, possibly due to different formation conditions. We are sorry that there is no further research on the causes of this phenomenon in this paper on account of insufficient high-frequency production data or complete adjacent test well data.Figure 2
Diagnostic plots for two wells in shale gas formations belonging to southern Sichuan. (a) Production rate plot for Well 3. (b) Production rate plot for Well 4.
(a)(b)


### 2.2. Producing Gas-Water Ratios (GWR)

In order to explore the trend of GWR of production data, this section simply processes and analyzes the production data to observe whether there is an immediate gas breakthrough after the production of shale gas wells [17, 27], which means that SRV can be calculated by establishing relevant models based on production data in this area.This is largely because the effective fracture network system is saturated with both the gas and water phases after two shut-ins. The gas source here is assumed to be from three aspects: (1) the originally existing initial gas in the active natural fracture, (2) the gas displaced by fracturing fluid under the influence of the strong countercurrent water imbibitions into the shale matrix due to the huge pressure difference during the first shut-in period, and (3) the gas accumulation that resulted from spontaneous imbibitions of fracturing fluid during the second shut-in period [41, 42].Similarly, before establishing the model, we also need to use the gas-water yield ratio to gain the V-shaped trend. Thus, we processed the initial two-phase production data of Well 3 and Well 4 and finally obtained the trend of the gas-water ratio (GWR). Figure3 shows a V-shaped trend GWR over time as a feature of early two-phase flowback data. Although the inflection points in the GWR curves of the two wells occur at different times due to different formation conditions and production systems, the same V-shaped trend still appears obviously. Ghanbari et al. [15] and Abbasi [16] also found the V-shaped GWR in the wells drilled in the Horn River shales in Canada.Figure 3
Diagnostic plots for two wells in shale gas formations belonging to southern Sichuan. (a) GWR plot for Well 3. (b) GWR plot for Well 4.
(a)(b)The gas-water ratio decreases and then increases, and we can study it from the most basic theory. Assume that the gas-water two-phase flow satisfies Darcy’s law:(1)Q=−kAμdPdL,where Q is the fluid flow rate under pressure P, m3/s; k is the fluid permeability, μm2; A is the fluid flow area, m2; μ is the fluid viscosity, mP·s; P is the pressure, MPa; and L is the fluid flow length, m.After ignoring the capillary force in the fracture network system, the gas-water ratio becomes(2)qgqw=μwμgPkrgSgkrwSw∂Pg∂Pw≈μwμgPkrgSgkrwSw,where qg and qw are the flow rates of the gas phase and water phase, respectively, m3/d; μw and μgP are the viscosity of the gas phase and water phase under P pressure, respectively, mP·s; krgSg and krwSw are the relative permeability of the gas phase and water phase under their respective saturation, μm2; Sg and Sw are the saturation of the gas phase and water phase, dimensionless; and Pg and Pw are the pressure of the gas phase and water phase, respectively, MPa.After the well opening for production, the gas viscosity decreases with the decrease of pressure, while the water viscosity remains relatively unchanged. According to the change of the gas-water ratio over time, the ratio of gas-water relative permeability decreases. Therefore, on the basis of the relationship between permeability and saturation, it can be inferred that the saturation ratio has a corresponding variation trend, meaning that the initial gas saturation (Sgi) is not zero, which also proves the existence of initial free gas in the fracture before well opening. When the wellbore storage effect disappears, the gas flows from the matrix to the fracture network, resulting in the replenishment of gas in the fracture system. Accordingly, the gas saturation increases, so the gas-water relative permeability increases and the gas-water ratio decreases. This is why the slope of the GWR curve drops first and then rises.In general, the negative slope on the GWR diagnostic plot is called early gas production (EGP) and the rise of the GWR diagnostic plot is called late gas production (LGP). It is generally believed that the EGP region is the stage of wellbore storage effects and both the gas and water productions come from the effective fracture connected with the horizontal well. Also, the increase of water relative permeability is greater than the reduction of gas viscosity. The LGP phase is the result of the matrix gas transfer to the fracture network after wellbore effects becoming negligible [10, 11, 17, 25–27]. According to the above analysis, we have figured out that in the EGP stage, the ground flow of flowback is only from the effective fracture system, while after the LGP stage beginning, the matrix gas initiates participation in the flow. Therefore, we will establish an analysis model for the EGP phase to calculate the SRV of shale gas wells after multistage fracturing.


## 2.1. Field Data from a Typical Shale Gas Reservoir in Southern China

In this work, we selected the flowback data from two production wells located in southern China: the shale gas formations of the Lower Silurianare, a typical marine shale gas reservoir. A large number of micron-nanoscale pores are developed in the mineral grains and organic matter [37]. These pores are the main channels for shale gas storage and migration, which play an important role in improving shale gas storage performance and providing good storage space for marine shale gas accumulation and enrichment [38]. The depth of the shale gas reservoir is over 2000 m at subsurface with a formation pressure coefficient up to 2.0. The thickness of the entire shale reservoir is about 100 m with porosity ranging from 2% to 4% and matrix permeability of 0.001 mD. The production at the early stage was about 10~20×104m3/d with a stable production around at 5~8×104m3/d.The east-west anticline structure belt with few faults is arranged in the left echelon row as the main geological structure of this area [39]. This drilling platform is located at the end of the anticline structure. Although two compresso-shear faults exist around the well pad, the overall upper lateral sealing performance leads to an effective accumulation and preservation of shale gas [37, 38]. Figure 1 shows the layout of the horizontal well group. Due to the low quartz content in the high-quality shale segment at the bottom of the Longmaxi Formation, natural fractures are not developed throughout the reservoir [40]. In the whole horizontal well group, Well 1, Well 4, and Well 5 are located in one formation, while Well 2, Well 3, Well 6, and Well 7 belong to another formation.Figure 1
Well pad schematic. Seven MFHWs drilled in the southern Sichuan Basin.In this work, we used data from Well 3 and Well 4 to test our model for the following two reasons. Firstly, the two wells were drilled through the shale gas reservoir in different directions and are located in two different formations. Secondly, microseismic monitoring is conducted in the two wells together with a comprehensive drilling, completion data, and well test data. Figure2 shows gas and water production rates of the two wells. The two wells both experienced gas breakthrough immediately after well opening and production, meaning that there was no single-phase flow. The difference is that Well 3 experienced a precipitous drop in production at one point (80 hrs), while Well 4 experienced a more gentle fluctuation, possibly due to different formation conditions. We are sorry that there is no further research on the causes of this phenomenon in this paper on account of insufficient high-frequency production data or complete adjacent test well data.Figure 2
Diagnostic plots for two wells in shale gas formations belonging to southern Sichuan. (a) Production rate plot for Well 3. (b) Production rate plot for Well 4.
(a)(b)


## 2.2. Producing Gas-Water Ratios (GWR)

In order to explore the trend of GWR of production data, this section simply processes and analyzes the production data to observe whether there is an immediate gas breakthrough after the production of shale gas wells [17, 27], which means that SRV can be calculated by establishing relevant models based on production data in this area.This is largely because the effective fracture network system is saturated with both the gas and water phases after two shut-ins. The gas source here is assumed to be from three aspects: (1) the originally existing initial gas in the active natural fracture, (2) the gas displaced by fracturing fluid under the influence of the strong countercurrent water imbibitions into the shale matrix due to the huge pressure difference during the first shut-in period, and (3) the gas accumulation that resulted from spontaneous imbibitions of fracturing fluid during the second shut-in period [41, 42].Similarly, before establishing the model, we also need to use the gas-water yield ratio to gain the V-shaped trend. Thus, we processed the initial two-phase production data of Well 3 and Well 4 and finally obtained the trend of the gas-water ratio (GWR). Figure3 shows a V-shaped trend GWR over time as a feature of early two-phase flowback data. Although the inflection points in the GWR curves of the two wells occur at different times due to different formation conditions and production systems, the same V-shaped trend still appears obviously. Ghanbari et al. [15] and Abbasi [16] also found the V-shaped GWR in the wells drilled in the Horn River shales in Canada.Figure 3
Diagnostic plots for two wells in shale gas formations belonging to southern Sichuan. (a) GWR plot for Well 3. (b) GWR plot for Well 4.
(a)(b)The gas-water ratio decreases and then increases, and we can study it from the most basic theory. Assume that the gas-water two-phase flow satisfies Darcy’s law:(1)Q=−kAμdPdL,where Q is the fluid flow rate under pressure P, m3/s; k is the fluid permeability, μm2; A is the fluid flow area, m2; μ is the fluid viscosity, mP·s; P is the pressure, MPa; and L is the fluid flow length, m.After ignoring the capillary force in the fracture network system, the gas-water ratio becomes(2)qgqw=μwμgPkrgSgkrwSw∂Pg∂Pw≈μwμgPkrgSgkrwSw,where qg and qw are the flow rates of the gas phase and water phase, respectively, m3/d; μw and μgP are the viscosity of the gas phase and water phase under P pressure, respectively, mP·s; krgSg and krwSw are the relative permeability of the gas phase and water phase under their respective saturation, μm2; Sg and Sw are the saturation of the gas phase and water phase, dimensionless; and Pg and Pw are the pressure of the gas phase and water phase, respectively, MPa.After the well opening for production, the gas viscosity decreases with the decrease of pressure, while the water viscosity remains relatively unchanged. According to the change of the gas-water ratio over time, the ratio of gas-water relative permeability decreases. Therefore, on the basis of the relationship between permeability and saturation, it can be inferred that the saturation ratio has a corresponding variation trend, meaning that the initial gas saturation (Sgi) is not zero, which also proves the existence of initial free gas in the fracture before well opening. When the wellbore storage effect disappears, the gas flows from the matrix to the fracture network, resulting in the replenishment of gas in the fracture system. Accordingly, the gas saturation increases, so the gas-water relative permeability increases and the gas-water ratio decreases. This is why the slope of the GWR curve drops first and then rises.In general, the negative slope on the GWR diagnostic plot is called early gas production (EGP) and the rise of the GWR diagnostic plot is called late gas production (LGP). It is generally believed that the EGP region is the stage of wellbore storage effects and both the gas and water productions come from the effective fracture connected with the horizontal well. Also, the increase of water relative permeability is greater than the reduction of gas viscosity. The LGP phase is the result of the matrix gas transfer to the fracture network after wellbore effects becoming negligible [10, 11, 17, 25–27]. According to the above analysis, we have figured out that in the EGP stage, the ground flow of flowback is only from the effective fracture system, while after the LGP stage beginning, the matrix gas initiates participation in the flow. Therefore, we will establish an analysis model for the EGP phase to calculate the SRV of shale gas wells after multistage fracturing.


## 3. Methods

### 3.1. Material Balance Equation (MBE) and Diffusion Equation for the Fracture System

In order to facilitate the establishment of the mathematical model for the EGP phase, the fracture network around the shale gas fractured well is simplified into the SRV region composed of the matrix system and fracture system. As shown in Figure4, the effective fracture system is made up of a fracture section and a matrix section. The fracture section consists of artificial hydraulic fractures, the secondary fractures generated by fracturing, and active natural fractures. The matrix section is the shale matrix connected to the fracture section. It is worth noting that the desorption and adsorption processes were ignored in the material balance equation (MBE) [17, 43].Figure 4
Schematic diagram of the horizontal well with multistage hydraulic fracturing for the development of material balance equation. Dashed arrows show fluid flow direction, which is sequentially from the matrix to fractures and fractures to the wellbore.This model simplifies complex, active natural and secondary fractures, as well as artificial hydraulic fractures, into a simple fracture system. The length of the artificial hydraulic fracture is used as the width of the entire stimulated reservoir volume, and the length of the horizontal wellbore is used as the length. In the whole stimulated reservoir area, the height of the major fracture, including the matrix part connected to the fracture system, participates in the flow. In the equivalent fracture system, the fracture is saturated with fracturing fluid (water phase) and natural gas (gas phase). It is assumed that no matrix gas is involved in the flow during the early gas production (EGP).It is assumed that the fracture system can be approximated as a homogeneous/closed/tank system. And the fluid flow from the fracture to the horizontal well is assumed to be linear. The mechanism driving the gas-water flow includes two aspects: (1) fracture closure and (2) expansion of the fluid (gas-water phase). Kuchuk et al. [44] also considered the assumption that fracture closure and expansion of water and free gas in the fracture system drive the accumulation of natural gas and water at the surface suggesting the occurrence of pseudo-steady-state flow in fractured reservoirs.In general, for the purpose of facilitating the establishment and solution of the model, we made the following assumptions: (1) capillary pressure in fracture systems is ignored; (2) the gas from the matrix is negligible in the EGP stage; (3) the fracture system is approximated as a homogeneous/closed system; (4) Darcy’s law applies to fluid flow; and (5) the effective fracture system is saturated with fracturing fluid (water phase) and natural gas (gas phase) initially.Effective compressibility terms are defined using simplified gas material balance equations (MBEs) [45–48]. In this work, the material balance equation (MBE) of the gas phase in the fracture network is simplified. The effects of gas and water (fracturing fluid) expansion and fracture closure are classified as the “effective compressibility” term given by
(3)C~t=1−GpGfiBgBgiSgiCg+1−WpWfiSwiCw+1Vfi∂Vf∂Pf,where Gp and Gfi are the cumulative gas production and initial volume of gas in the fracture, respectively, m3; Bg and Bgi are the gas formation volume factor and gas formation volume factor at initial conditions, respectively, m3/m3; Sgi and Swi are the initial gas phase saturation and initial water phase saturation, respectively, dimensionless; Cg and Cw are the gas compressibility and water compressibility, respectively, Pa-1; Wp and Wfi are the cumulative water production and initial volume of water in the fracture, respectively, m3; Vf and Vfi are the volume of effective fractures and volume of effective fractures at initial conditions, m3; and Pf is the fracture pressure, Pa.C~t is analogous to the total compressibility term proposed for conventional multiphase well testing. This is represented as a function of measurable flowback parameters, including cumulative gas production and cumulative water production. Each term in its expression outlines different driving mechanisms in the fracture system: (1) expansion of the gas phase, (2) expansion of liquid phase, and (3) fracture closure. Among them, 1/Vfi∂Vf/∂Pf is similar to the formation compressibility coefficient in traditional conventional material balance analysis [49]. When dealing with fractures, it refers to the inverse of the fracture stiffness determined by the elastic or strain energy required to keep the hydraulic fracture open [50].Ignoring the gas flow rate from the matrix system to the fracture, the gas phase material balance equation is(4)0−qgρg0=∂∂tVgtρgR,where ρg0 and ρgR are the density of gas at surface conditions and reservoir conditions, respectively, kg/m3; and Vgt is the volume of gas in the fracture system at any time, m3.According to the relationship between gas volumeVg, fracture volume Vf, and water volume Vw in the fracture system, the definition of gas compressibility is obtained as
(5)−qgBg=Vgcg∂Pf∂t−∂Vw∂t+∂Vf∂t,where Vf, Vg, and Vw are the fracture volume, gas volume, and water volume in the fracture system, respectively, m3.By substituting the gas productionGp and water production Wp, then using the chain rule, the total volume of fractures is derived:
(6)−1VfiqgBg+qwBw=Gfi−GpBgcgGfiBgi/Sgi∂Pf∂t+Wfi−WpBwcwWfiBwi/Swi∂Pf∂t+1Vfi∂Vf∂Pf∂Pf∂t,where Bw is the water formation volume factor at initial conditions, m3/m3.In the early gas production (EGP), it is assumed thatBw≈Bwi, the final material balance equation, can be obtained by defining a compressibility term and the total fluid flow rate by simplifying the above equation is
(7)∂Pf∂t=−qtC~tVfi,where qt=qgBg+qwB is total producing rate, m3/d.Given that single-phase, steady-state flow can be described using the continuity equation and Darcy’s law, the single-phase gas diffusion equation in the fracture system is given by Zhang and Winter [51]:
(8)∇PfμgZ∇Pf=ϕfKfPfZC~t∂Pf∂t,where Z is the gas compressibility factor, dimensionless; φf is the porosity for the fracture system, dimensionless; and ∇ is the gradient operator.Define the pseudopressure and pseudotime functions [46, 52] as
(9)ψPf=∫0Pf2PfμgZ∂Pf,(10)ta=∫0tkgtμgC~t∂t,where kgt is the relative permeability of the gas at a given time, μm2; and C~t is the total effective compressibility, Pa-1.Then, the governing equation of single-phase gas flow in the fracture system is(11)∂2ψPf∂y2=ϕfKf∂ψ∂ta.


### 3.2. Coupling Stress Sensitivity with Dynamic Relative Permeability Function

In Equation (10), kgt is the variable function of single-phase gas with time. By introducing the dynamic relative permeability (DRP) function krgt, the diffusion equation of single-phase gas is transformed into the diffusion equation of the two-phase system [53].The method in this paper is similar to that in Ezulike and Dehghanpour’s study [53], but due to the geological conditions of formations belonging to southern Sichuan being different from those of North America, the dynamic relative permeability (DRP) function of the gas phase is obtained by using the method of fitting field data. In the previous study, we have assumed that the fracture system is saturated with gas and water. The relationship between cumulative gas and water production and time can be reported through field data. Therefore, the gas-water saturation in the fracture system can be established as a function of gas-water production with time. Given the characteristics of abnormal overpressure in this area, the stress sensitivity in this area cannot be ignored [36, 54, 55]. Therefore, the stress-sensitive model was considered in the dynamic relative permeability function. The specific process is shown in Figure 5.Figure 5
Procedure to estimate the dynamic relative permeability (DRP) curve considering stress sensitivity.In Figure6, the relative permeability of the cores after fracturing in the research area is tested through laboratory experiments. The curve of pressure change with time during the flowback period is shown in Figure 7. Stress-strain experiments were performed to describe the curves between dimensionless permeability and effective stress by stabilizing the confining pressure and reducing the internal pressure. This is more in line with the actual formation flowback pressure changes. The early flowback stage of shale gas wells after well opening is significantly short; therefore, it can be assumed that the formation pressure has not changed much, but the inner pressure begins to decrease after the fluid flows out of the cores.Figure 6
Gas-water relative permeability curve of the cores from Well 3 and Well 4.Figure 7
Casing pressure and effective stress changed with time for two wells in shale gas formations belonging to southern Sichuan. (a) Well 3 and (b) Well 4.
(a)(b)Take the stress sensitivity experiments of cores belonging to the stimulated area of Well 3 as an example. Both the matrix and the fractures are most likely contained in these cores. We had required as many experiments as possible to obtain the stress-sensitive test data in the stimulated region after fracturing as accurately as possible. However, subject to the insufficient samples or the different distances of samples from the horizontal well, these data maybe cannot represent the whole area in fact but still have considerable reference and research value.The specific experimental procedures are as follows: (1) the initial confining pressure was set as the original formation pressure which is 38 MPa, and the internal pressure was 23 MPa; (2) the confining pressure was increased to 58 MPa slowly which is the formation pressure before well opening for flowback, and the internal pressure was increased to 43 MPa at the same time to keep the effective stress constant; and (3) the internal pressure was reduced to different pressure points to increase the effective stress, and the gas permeability of the sample was measured after each pressure point was stabilized.Dimensionless permeability is defined as(12)KK0=ae−bpi−pf,where K is the relative permeability of the gas at a given time, μm2; K0 is the initial gas relative permeability, μm2; and a are b are stress sensitivity coefficients by experiments, dimensionless.The experimental results are reported in Figure8. According to the calculation method in Figure 5, the results in Figure 9 were finally obtained. Then, according to the fitting curve, the function of gas-phase relative permeability changing with time is
(13)krgt=αt−β.Figure 8
Stress-strain experiment results of Well 3 and Well 4.Figure 9
Dynamic relative permeability function over time considering stress sensitivity for two wells in shale gas formations belonging to southern Sichuan. (a) DRP fitting curve for Well 3. (b) DRP fitting curve for Well 4.
(a)(b)Since the above equation is obtained based on flowback data and core experiment relative permeability,α and β are called the flowback data coefficient and flowback data intercept, respectively, which are determined by the data fitting curve.


### 3.3. Two-Phase Flowback Model for Early Gas Production (EGP)

By substituting Equation (13) into Equation (10), the governing equation of single-phase gas flow in the fracture system is changed into the two-phase flow model. This means that when gas-water two-phase flow occurs in the fracture system, the new real pseudotime function is transformed into
(14)ta=∫0tkrgtμgC~t∂t.In order to establish the relationship between pseudopressure function and pseudotime function, we make the following transformation:(15)∂ψ∂ta=∂Pf∂t×∂ψ∂Pf×∂t∂ta.And∂Pf/∂t is derived by the final material balance equation (Equation (7)); in a similar way, we can gain ∂ψ/∂Pf and ∂t/∂ta by the defined pseudopressure (Equation (9)) and new pseudotime equation (Equation (14)), respectively. Finally, Equation (15) becomes
(16)∂ψ∂ta=−2VfiqtkrgtPfZ.Define the equivalent gas rate as [48]
(17)qg∗=1krgtqgBgi+qwBw.SubstitutingPf/Z in Equation (16) using the real gas law, then combine it with the equivalent gas rate (Equation (17)):
(18)∂ψ∂ta=−2VfiPiZiqg∗.By substituting∂ψ/∂ta (Equation (18)) in Equation (11), we can establish the two-phase diffusion equation:
(19)∂2ψPf∂y2=ϕfKf−2VfiPiZiqg∗.Then, the following inner and outer boundary conditions can be solved:y=Ye, ∂ψPf/∂y=0, and y=0, ψPf=ψPwf.Define the fracture storage coefficient as(20)Cst=VfiZi2Pi.Then, the following relationship is gained:(21)ψP¯f−ψPwfqg∗=ϕfKf1CstYe23.Substitute Equation (20) into Equation (18) and integrate to obtain
(22)ta=CstψPi−ψP¯fqg∗.By combining Equation (22) with Equation (15), we can modify the final two-phase flow model of the early gas production stage as follows:
(23)ψPi−ψPwfqg∗=1Cstta+ϕfKf1CstYe23.Theoretically, a plot of the rate-normalized pseudopressure (RNP) vs. the pseudotime should yield a straight line relationship. With the pseudotime function as an independent variable and the left side of the equation as a dependent variable, the slope and intercept can be obtained according to the fitting curve, and the relationship between the equivalent fracture porosity/half-length of the effective fracture system and effective fracture system permeability can be described as follows:(24)Cst=1slope=VfiZi2Pi,(25)YD=y‐axisinterceptslope=ϕfKfYe23.


### 3.4. Analysis Procedure

We propose the following analysis procedure:(1)
Obtain and process water and gas flowback data to explore a V-shaped gas-water ratio trend (see Figures2 and 3 and Equation (2))(2)
The early flowback period (EGP) is distinguished from the late flowback period (LGP) according to the V-shaped trend of the GWR curve (see Figure3)(3)
Conduct a simpler fracture network system model for the EGP (see Figure4)(4)
Calculate effective compressibility by Equation (3) to deduce the two-phase material balance equation (MBE) (Equation (7)) for the fracture system(5)
Define the pseudopressure and pseudotime functions (Equations (9) and (10)) to deduce the gas diffusion equation (Equation (11)) of the fracture system(6)
Plot the gas-water relative permeability curve (see Figure5) and the curves of dimensionless permeability with effective stress (see Figure 8) by experiments(7)
Calculate gas DRP coupling stress sensitivity (see Figure9) following the steps by Figure 6 to transform the gas diffusion equation (Equation (11)) into the two-phase diffusion equation (Equation (19))(8)
Modify the final two-phase flow mode (Equation (23)) by combining Equation (22) with Equation (15)(9)
Plot rate-normalized pressure change with pseudotime (see Figure10)(10)
Calculate SRV by Equation (24) and effective fracture system permeability by Equation (25)(11)
Verify the analytical model against microseismic data (see Figure11)Figure 10
Analysis of EGP data of two wells belonging to southern Sichuan: rate-normalized pressure change with pseudotime. (a) Well 3 and (b) Well 4.
(a)(b)Figure 11
(a) Renderings of microseismic detection technology for Well 3 and Well 4. (b) Comparison between the analytical model and the microseismic data.
(a)(b)The flowback data we need to obtain include production rates and pressure and cumulative production data profiles. Then, we got a V-shaped trend in the gas-water ratio curve by processing the flowback data. Thus, an analytical model was established for the EGP stage. Finally, the calculated results were compared with the microseismic data to validate the mathematical results.


## 3.1. Material Balance Equation (MBE) and Diffusion Equation for the Fracture System

In order to facilitate the establishment of the mathematical model for the EGP phase, the fracture network around the shale gas fractured well is simplified into the SRV region composed of the matrix system and fracture system. As shown in Figure4, the effective fracture system is made up of a fracture section and a matrix section. The fracture section consists of artificial hydraulic fractures, the secondary fractures generated by fracturing, and active natural fractures. The matrix section is the shale matrix connected to the fracture section. It is worth noting that the desorption and adsorption processes were ignored in the material balance equation (MBE) [17, 43].Figure 4
Schematic diagram of the horizontal well with multistage hydraulic fracturing for the development of material balance equation. Dashed arrows show fluid flow direction, which is sequentially from the matrix to fractures and fractures to the wellbore.This model simplifies complex, active natural and secondary fractures, as well as artificial hydraulic fractures, into a simple fracture system. The length of the artificial hydraulic fracture is used as the width of the entire stimulated reservoir volume, and the length of the horizontal wellbore is used as the length. In the whole stimulated reservoir area, the height of the major fracture, including the matrix part connected to the fracture system, participates in the flow. In the equivalent fracture system, the fracture is saturated with fracturing fluid (water phase) and natural gas (gas phase). It is assumed that no matrix gas is involved in the flow during the early gas production (EGP).It is assumed that the fracture system can be approximated as a homogeneous/closed/tank system. And the fluid flow from the fracture to the horizontal well is assumed to be linear. The mechanism driving the gas-water flow includes two aspects: (1) fracture closure and (2) expansion of the fluid (gas-water phase). Kuchuk et al. [44] also considered the assumption that fracture closure and expansion of water and free gas in the fracture system drive the accumulation of natural gas and water at the surface suggesting the occurrence of pseudo-steady-state flow in fractured reservoirs.In general, for the purpose of facilitating the establishment and solution of the model, we made the following assumptions: (1) capillary pressure in fracture systems is ignored; (2) the gas from the matrix is negligible in the EGP stage; (3) the fracture system is approximated as a homogeneous/closed system; (4) Darcy’s law applies to fluid flow; and (5) the effective fracture system is saturated with fracturing fluid (water phase) and natural gas (gas phase) initially.Effective compressibility terms are defined using simplified gas material balance equations (MBEs) [45–48]. In this work, the material balance equation (MBE) of the gas phase in the fracture network is simplified. The effects of gas and water (fracturing fluid) expansion and fracture closure are classified as the “effective compressibility” term given by
(3)C~t=1−GpGfiBgBgiSgiCg+1−WpWfiSwiCw+1Vfi∂Vf∂Pf,where Gp and Gfi are the cumulative gas production and initial volume of gas in the fracture, respectively, m3; Bg and Bgi are the gas formation volume factor and gas formation volume factor at initial conditions, respectively, m3/m3; Sgi and Swi are the initial gas phase saturation and initial water phase saturation, respectively, dimensionless; Cg and Cw are the gas compressibility and water compressibility, respectively, Pa-1; Wp and Wfi are the cumulative water production and initial volume of water in the fracture, respectively, m3; Vf and Vfi are the volume of effective fractures and volume of effective fractures at initial conditions, m3; and Pf is the fracture pressure, Pa.C~t is analogous to the total compressibility term proposed for conventional multiphase well testing. This is represented as a function of measurable flowback parameters, including cumulative gas production and cumulative water production. Each term in its expression outlines different driving mechanisms in the fracture system: (1) expansion of the gas phase, (2) expansion of liquid phase, and (3) fracture closure. Among them, 1/Vfi∂Vf/∂Pf is similar to the formation compressibility coefficient in traditional conventional material balance analysis [49]. When dealing with fractures, it refers to the inverse of the fracture stiffness determined by the elastic or strain energy required to keep the hydraulic fracture open [50].Ignoring the gas flow rate from the matrix system to the fracture, the gas phase material balance equation is(4)0−qgρg0=∂∂tVgtρgR,where ρg0 and ρgR are the density of gas at surface conditions and reservoir conditions, respectively, kg/m3; and Vgt is the volume of gas in the fracture system at any time, m3.According to the relationship between gas volumeVg, fracture volume Vf, and water volume Vw in the fracture system, the definition of gas compressibility is obtained as
(5)−qgBg=Vgcg∂Pf∂t−∂Vw∂t+∂Vf∂t,where Vf, Vg, and Vw are the fracture volume, gas volume, and water volume in the fracture system, respectively, m3.By substituting the gas productionGp and water production Wp, then using the chain rule, the total volume of fractures is derived:
(6)−1VfiqgBg+qwBw=Gfi−GpBgcgGfiBgi/Sgi∂Pf∂t+Wfi−WpBwcwWfiBwi/Swi∂Pf∂t+1Vfi∂Vf∂Pf∂Pf∂t,where Bw is the water formation volume factor at initial conditions, m3/m3.In the early gas production (EGP), it is assumed thatBw≈Bwi, the final material balance equation, can be obtained by defining a compressibility term and the total fluid flow rate by simplifying the above equation is
(7)∂Pf∂t=−qtC~tVfi,where qt=qgBg+qwB is total producing rate, m3/d.Given that single-phase, steady-state flow can be described using the continuity equation and Darcy’s law, the single-phase gas diffusion equation in the fracture system is given by Zhang and Winter [51]:
(8)∇PfμgZ∇Pf=ϕfKfPfZC~t∂Pf∂t,where Z is the gas compressibility factor, dimensionless; φf is the porosity for the fracture system, dimensionless; and ∇ is the gradient operator.Define the pseudopressure and pseudotime functions [46, 52] as
(9)ψPf=∫0Pf2PfμgZ∂Pf,(10)ta=∫0tkgtμgC~t∂t,where kgt is the relative permeability of the gas at a given time, μm2; and C~t is the total effective compressibility, Pa-1.Then, the governing equation of single-phase gas flow in the fracture system is(11)∂2ψPf∂y2=ϕfKf∂ψ∂ta.


## 3.2. Coupling Stress Sensitivity with Dynamic Relative Permeability Function

In Equation (10), kgt is the variable function of single-phase gas with time. By introducing the dynamic relative permeability (DRP) function krgt, the diffusion equation of single-phase gas is transformed into the diffusion equation of the two-phase system [53].The method in this paper is similar to that in Ezulike and Dehghanpour’s study [53], but due to the geological conditions of formations belonging to southern Sichuan being different from those of North America, the dynamic relative permeability (DRP) function of the gas phase is obtained by using the method of fitting field data. In the previous study, we have assumed that the fracture system is saturated with gas and water. The relationship between cumulative gas and water production and time can be reported through field data. Therefore, the gas-water saturation in the fracture system can be established as a function of gas-water production with time. Given the characteristics of abnormal overpressure in this area, the stress sensitivity in this area cannot be ignored [36, 54, 55]. Therefore, the stress-sensitive model was considered in the dynamic relative permeability function. The specific process is shown in Figure 5.Figure 5
Procedure to estimate the dynamic relative permeability (DRP) curve considering stress sensitivity.In Figure6, the relative permeability of the cores after fracturing in the research area is tested through laboratory experiments. The curve of pressure change with time during the flowback period is shown in Figure 7. Stress-strain experiments were performed to describe the curves between dimensionless permeability and effective stress by stabilizing the confining pressure and reducing the internal pressure. This is more in line with the actual formation flowback pressure changes. The early flowback stage of shale gas wells after well opening is significantly short; therefore, it can be assumed that the formation pressure has not changed much, but the inner pressure begins to decrease after the fluid flows out of the cores.Figure 6
Gas-water relative permeability curve of the cores from Well 3 and Well 4.Figure 7
Casing pressure and effective stress changed with time for two wells in shale gas formations belonging to southern Sichuan. (a) Well 3 and (b) Well 4.
(a)(b)Take the stress sensitivity experiments of cores belonging to the stimulated area of Well 3 as an example. Both the matrix and the fractures are most likely contained in these cores. We had required as many experiments as possible to obtain the stress-sensitive test data in the stimulated region after fracturing as accurately as possible. However, subject to the insufficient samples or the different distances of samples from the horizontal well, these data maybe cannot represent the whole area in fact but still have considerable reference and research value.The specific experimental procedures are as follows: (1) the initial confining pressure was set as the original formation pressure which is 38 MPa, and the internal pressure was 23 MPa; (2) the confining pressure was increased to 58 MPa slowly which is the formation pressure before well opening for flowback, and the internal pressure was increased to 43 MPa at the same time to keep the effective stress constant; and (3) the internal pressure was reduced to different pressure points to increase the effective stress, and the gas permeability of the sample was measured after each pressure point was stabilized.Dimensionless permeability is defined as(12)KK0=ae−bpi−pf,where K is the relative permeability of the gas at a given time, μm2; K0 is the initial gas relative permeability, μm2; and a are b are stress sensitivity coefficients by experiments, dimensionless.The experimental results are reported in Figure8. According to the calculation method in Figure 5, the results in Figure 9 were finally obtained. Then, according to the fitting curve, the function of gas-phase relative permeability changing with time is
(13)krgt=αt−β.Figure 8
Stress-strain experiment results of Well 3 and Well 4.Figure 9
Dynamic relative permeability function over time considering stress sensitivity for two wells in shale gas formations belonging to southern Sichuan. (a) DRP fitting curve for Well 3. (b) DRP fitting curve for Well 4.
(a)(b)Since the above equation is obtained based on flowback data and core experiment relative permeability,α and β are called the flowback data coefficient and flowback data intercept, respectively, which are determined by the data fitting curve.


## 3.3. Two-Phase Flowback Model for Early Gas Production (EGP)

By substituting Equation (13) into Equation (10), the governing equation of single-phase gas flow in the fracture system is changed into the two-phase flow model. This means that when gas-water two-phase flow occurs in the fracture system, the new real pseudotime function is transformed into
(14)ta=∫0tkrgtμgC~t∂t.In order to establish the relationship between pseudopressure function and pseudotime function, we make the following transformation:(15)∂ψ∂ta=∂Pf∂t×∂ψ∂Pf×∂t∂ta.And∂Pf/∂t is derived by the final material balance equation (Equation (7)); in a similar way, we can gain ∂ψ/∂Pf and ∂t/∂ta by the defined pseudopressure (Equation (9)) and new pseudotime equation (Equation (14)), respectively. Finally, Equation (15) becomes
(16)∂ψ∂ta=−2VfiqtkrgtPfZ.Define the equivalent gas rate as [48]
(17)qg∗=1krgtqgBgi+qwBw.SubstitutingPf/Z in Equation (16) using the real gas law, then combine it with the equivalent gas rate (Equation (17)):
(18)∂ψ∂ta=−2VfiPiZiqg∗.By substituting∂ψ/∂ta (Equation (18)) in Equation (11), we can establish the two-phase diffusion equation:
(19)∂2ψPf∂y2=ϕfKf−2VfiPiZiqg∗.Then, the following inner and outer boundary conditions can be solved:y=Ye, ∂ψPf/∂y=0, and y=0, ψPf=ψPwf.Define the fracture storage coefficient as(20)Cst=VfiZi2Pi.Then, the following relationship is gained:(21)ψP¯f−ψPwfqg∗=ϕfKf1CstYe23.Substitute Equation (20) into Equation (18) and integrate to obtain
(22)ta=CstψPi−ψP¯fqg∗.By combining Equation (22) with Equation (15), we can modify the final two-phase flow model of the early gas production stage as follows:
(23)ψPi−ψPwfqg∗=1Cstta+ϕfKf1CstYe23.Theoretically, a plot of the rate-normalized pseudopressure (RNP) vs. the pseudotime should yield a straight line relationship. With the pseudotime function as an independent variable and the left side of the equation as a dependent variable, the slope and intercept can be obtained according to the fitting curve, and the relationship between the equivalent fracture porosity/half-length of the effective fracture system and effective fracture system permeability can be described as follows:(24)Cst=1slope=VfiZi2Pi,(25)YD=y‐axisinterceptslope=ϕfKfYe23.


## 3.4. Analysis Procedure

We propose the following analysis procedure:(1)
Obtain and process water and gas flowback data to explore a V-shaped gas-water ratio trend (see Figures2 and 3 and Equation (2))(2)
The early flowback period (EGP) is distinguished from the late flowback period (LGP) according to the V-shaped trend of the GWR curve (see Figure3)(3)
Conduct a simpler fracture network system model for the EGP (see Figure4)(4)
Calculate effective compressibility by Equation (3) to deduce the two-phase material balance equation (MBE) (Equation (7)) for the fracture system(5)
Define the pseudopressure and pseudotime functions (Equations (9) and (10)) to deduce the gas diffusion equation (Equation (11)) of the fracture system(6)
Plot the gas-water relative permeability curve (see Figure5) and the curves of dimensionless permeability with effective stress (see Figure 8) by experiments(7)
Calculate gas DRP coupling stress sensitivity (see Figure9) following the steps by Figure 6 to transform the gas diffusion equation (Equation (11)) into the two-phase diffusion equation (Equation (19))(8)
Modify the final two-phase flow mode (Equation (23)) by combining Equation (22) with Equation (15)(9)
Plot rate-normalized pressure change with pseudotime (see Figure10)(10)
Calculate SRV by Equation (24) and effective fracture system permeability by Equation (25)(11)
Verify the analytical model against microseismic data (see Figure11)Figure 10
Analysis of EGP data of two wells belonging to southern Sichuan: rate-normalized pressure change with pseudotime. (a) Well 3 and (b) Well 4.
(a)(b)Figure 11
(a) Renderings of microseismic detection technology for Well 3 and Well 4. (b) Comparison between the analytical model and the microseismic data.
(a)(b)The flowback data we need to obtain include production rates and pressure and cumulative production data profiles. Then, we got a V-shaped trend in the gas-water ratio curve by processing the flowback data. Thus, an analytical model was established for the EGP stage. Finally, the calculated results were compared with the microseismic data to validate the mathematical results.


## 4. Results and Discussions

### 4.1. Stimulated Reservoir Volume

We apply the analytical model presented above to analyze the flowback data of Well 3 and Well 4. However, there are several issues that need to be addressed and discussed:(1)
It is difficult to gain an appropriate initial gas saturation of the fracture system from actual field data. Unlike conventional numerical simulation, this parameter is unknown in actual field data. After the fracturing operation is completed, the effective fracture system connecting the wellbore is filled with fracturing fluid approximately. Thus, after the well is opened for flowback operation and before production, the volume of the recovered fracturing fluid under the ground is filled with gas renewedly. Therefore, in this paper, we take a reasonable value which is the recovery percentage of the total injection amount of fracturing fluid as the initial fracture system gas saturation in calculation. According to Xu’s paper [48], there is a little significant difference between the actual value of initial gas saturation and the recovery percentage. However, this was still not rigorous enough. Then, other possible values are substituted into the calculation by us to get the new SRV. The results showed that the impact of this value on SRV estimation was not obvious within the error range, and the difference was less than 1% compared with the microseismic dataDue to that, fracturing fluid may leak off into the existing inactive natural fractures and into the matrix during injection; the fluid flowback percentage decreases, but the initial gas saturation in the fracture system increases. Hence, when a significant portion of the injected fluid volume does not contribute to create fracture volume, this value of the initial gas saturation for calculation is actually low. One direction of future work is to consider using the results to iteratively optimize the initial gas saturation value.(2)
The value of the fracture closure term in the total compressibility cannot be accurately expressed. Fortunately, it is found that this value has no obvious influence on the results in the calculation process. Since this value refers to the inverse of fracture stiffness when dealing with fractures, it is considered to give a reasonable value in the subsequent work from the perspective of rock fracture mechanicsThe final model calculation results are shown in Figure10. Since microseismic detection technology was used to estimate the final hydraulic fracture stimulated effect in the two wells. We compared the calculated results with the microseismic data and finally found that the difference between the calculated results of the analytical model and the microseismic data was minor, as shown in Figure 11. This reason is that quantities of fractures that do not connect effectively with fracture systems are counted in the final result when microseismic detection techniques are used. Those “dead” fractures skew the microseismic results.


### 4.2. Effective Fracture System Permeability of Stimulated Reservoir

As shown in Table1, after obtaining the slope of rate-normalized pressure change with pseudotime (see Figure 10), SRV was gained by Equation (17). On the basis of the average length of hydraulic fractures provided by the microseismic data of the drilling platform where Well 3 and Well 4 are located and the average porosity of the cores after fracturing obtained through experiments, the permeability of the effective fracture system in the stimulated region is calculated by using Equation (18). The results are shown in Table 1, where the absolute deviation refers to the difference between the SRV calculated by the analytical model and estimated by the microseismic data, and the relative deviation refers to the ratio of the absolute deviation to the SRV estimated by microseismic data. Numerical errors in porosity of the effective fracture system are responsible for the error of permeability.Table 1
Different results of the analytical model and microseismic data.
Parameter nameWell 3Well 4SRV calculated by the analytical model3.4811×107m33.8478×107m3SRV estimated by the microseismic data3.5966×107m34.1305×107m3The absolute deviation1.1150×106m32.8270×106m3The relative deviation3.2%6.8%Effective fracture system permeability0.301 mD0.244 mDAverage half-length of hydraulic fractures320 m285 mEffective fracture system porosity15.5%13.8%If the SRV estimated by microseismic data is taken as a reference for comparison, the specific calculation process is as follows:(26)dA=VA−VM,(27)dR=VA−VMVM,where dA and dR are the absolute deviation and the relative deviation between the SRV calculated by the analytical model and estimated by the microseismic data, respectively, m3; and VA and VM are the SRV calculated by the analytical model and the SRV estimated by the microseismic data, m3.


### 4.3. Analytical Model without Considering Stress Sensitivity

If the stress sensitivity effect is not considered in the process of the DRP calculation, the deviation of the DRP curve will occur, as shown in Figure12. And the curve with a large error was substituted into the step shown in Section 3.4 to get the result in Figure 13. As a result, both the slope and the intercept of the RNP curve had changed. We obtained the new slope and intercept of the curve and calculated a new SRV according to Equations (17) and (18). Then, the calculation results and corresponding deviations as shown in Figure 14 were obtained.Figure 12
Dynamic relative permeability function over time without stress sensitivity for two wells belonging to southern Sichuan. (a) DRP for Well 3. (b) DRP for Well 4.
(a)(b)Figure 13
Analysis without considering stress sensitivity of EGP data of two wells in shale gas belonging to southern Sichuan: rate-normalized pressure change with pseudotime. (a) Well 3 and (b) Well 4.
(a)(b)Figure 14
Comparison between the analytical model with or without considering stress sensitivity and the microseismic data.It can be seen that the relative deviation calculated without considering the stress sensitivity effect is greater than that calculated with considering the stress sensitivity effect. Moreover, because the real SRV has been overestimated by the microseismic data, the calculation results of the model without considering the stress sensitivity effect are higher than those of the microseismic data, which further indicates that the calculation results without considering the stress sensitivity effect in the high-pressure area will have a large error, which cannot be ignored. We also calculated the permeability of the effective fracture system, and the effective permeability of Well 3 and Well 4 was 0.432 mD and 0.283 mD, respectively. It can be seen that the calculation results are 43.50% and 15.75% more than those considering the stress sensitivity effect.


## 4.1. Stimulated Reservoir Volume

We apply the analytical model presented above to analyze the flowback data of Well 3 and Well 4. However, there are several issues that need to be addressed and discussed:(1)
It is difficult to gain an appropriate initial gas saturation of the fracture system from actual field data. Unlike conventional numerical simulation, this parameter is unknown in actual field data. After the fracturing operation is completed, the effective fracture system connecting the wellbore is filled with fracturing fluid approximately. Thus, after the well is opened for flowback operation and before production, the volume of the recovered fracturing fluid under the ground is filled with gas renewedly. Therefore, in this paper, we take a reasonable value which is the recovery percentage of the total injection amount of fracturing fluid as the initial fracture system gas saturation in calculation. According to Xu’s paper [48], there is a little significant difference between the actual value of initial gas saturation and the recovery percentage. However, this was still not rigorous enough. Then, other possible values are substituted into the calculation by us to get the new SRV. The results showed that the impact of this value on SRV estimation was not obvious within the error range, and the difference was less than 1% compared with the microseismic dataDue to that, fracturing fluid may leak off into the existing inactive natural fractures and into the matrix during injection; the fluid flowback percentage decreases, but the initial gas saturation in the fracture system increases. Hence, when a significant portion of the injected fluid volume does not contribute to create fracture volume, this value of the initial gas saturation for calculation is actually low. One direction of future work is to consider using the results to iteratively optimize the initial gas saturation value.(2)
The value of the fracture closure term in the total compressibility cannot be accurately expressed. Fortunately, it is found that this value has no obvious influence on the results in the calculation process. Since this value refers to the inverse of fracture stiffness when dealing with fractures, it is considered to give a reasonable value in the subsequent work from the perspective of rock fracture mechanicsThe final model calculation results are shown in Figure10. Since microseismic detection technology was used to estimate the final hydraulic fracture stimulated effect in the two wells. We compared the calculated results with the microseismic data and finally found that the difference between the calculated results of the analytical model and the microseismic data was minor, as shown in Figure 11. This reason is that quantities of fractures that do not connect effectively with fracture systems are counted in the final result when microseismic detection techniques are used. Those “dead” fractures skew the microseismic results.


## 4.2. Effective Fracture System Permeability of Stimulated Reservoir

As shown in Table1, after obtaining the slope of rate-normalized pressure change with pseudotime (see Figure 10), SRV was gained by Equation (17). On the basis of the average length of hydraulic fractures provided by the microseismic data of the drilling platform where Well 3 and Well 4 are located and the average porosity of the cores after fracturing obtained through experiments, the permeability of the effective fracture system in the stimulated region is calculated by using Equation (18). The results are shown in Table 1, where the absolute deviation refers to the difference between the SRV calculated by the analytical model and estimated by the microseismic data, and the relative deviation refers to the ratio of the absolute deviation to the SRV estimated by microseismic data. Numerical errors in porosity of the effective fracture system are responsible for the error of permeability.Table 1
Different results of the analytical model and microseismic data.
Parameter nameWell 3Well 4SRV calculated by the analytical model3.4811×107m33.8478×107m3SRV estimated by the microseismic data3.5966×107m34.1305×107m3The absolute deviation1.1150×106m32.8270×106m3The relative deviation3.2%6.8%Effective fracture system permeability0.301 mD0.244 mDAverage half-length of hydraulic fractures320 m285 mEffective fracture system porosity15.5%13.8%If the SRV estimated by microseismic data is taken as a reference for comparison, the specific calculation process is as follows:(26)dA=VA−VM,(27)dR=VA−VMVM,where dA and dR are the absolute deviation and the relative deviation between the SRV calculated by the analytical model and estimated by the microseismic data, respectively, m3; and VA and VM are the SRV calculated by the analytical model and the SRV estimated by the microseismic data, m3.


## 4.3. Analytical Model without Considering Stress Sensitivity

If the stress sensitivity effect is not considered in the process of the DRP calculation, the deviation of the DRP curve will occur, as shown in Figure12. And the curve with a large error was substituted into the step shown in Section 3.4 to get the result in Figure 13. As a result, both the slope and the intercept of the RNP curve had changed. We obtained the new slope and intercept of the curve and calculated a new SRV according to Equations (17) and (18). Then, the calculation results and corresponding deviations as shown in Figure 14 were obtained.Figure 12
Dynamic relative permeability function over time without stress sensitivity for two wells belonging to southern Sichuan. (a) DRP for Well 3. (b) DRP for Well 4.
(a)(b)Figure 13
Analysis without considering stress sensitivity of EGP data of two wells in shale gas belonging to southern Sichuan: rate-normalized pressure change with pseudotime. (a) Well 3 and (b) Well 4.
(a)(b)Figure 14
Comparison between the analytical model with or without considering stress sensitivity and the microseismic data.It can be seen that the relative deviation calculated without considering the stress sensitivity effect is greater than that calculated with considering the stress sensitivity effect. Moreover, because the real SRV has been overestimated by the microseismic data, the calculation results of the model without considering the stress sensitivity effect are higher than those of the microseismic data, which further indicates that the calculation results without considering the stress sensitivity effect in the high-pressure area will have a large error, which cannot be ignored. We also calculated the permeability of the effective fracture system, and the effective permeability of Well 3 and Well 4 was 0.432 mD and 0.283 mD, respectively. It can be seen that the calculation results are 43.50% and 15.75% more than those considering the stress sensitivity effect.


## 5. Conclusions

In this paper, a comprehensive model is developed to predict SRV in the shale gas reservoirs using the flowback data after fracturing. The stress sensitivity is included as well. The results of SRV are compared with the microseismic data. The following conclusions can be drawn:(1)
The flowback data show that the gas-water ratio is V-shaped, i.e., the early descending stage and the late ascending stage in this field, which can be used to estimate the volume of the effective fracture system(2)
The stress sensitivity is a key factor affecting the permeability of the effective fracture system as well as the SRV in the shale gas reservoirs. Once it is neglected, the estimation will be overestimated(3)
Stress sensitivity is taken into account to forecast SRV in this typical block. Results show that the relative deviation of stimulated reservoir volume calculated by this proposed model and the microseismic data is less than 10%, indicating that this method could provide reasonable prediction



---

*Source: 8886988-2020-09-09.xml*