An Analytical Method for Parameter Interpretation of Fracture Networks in Shale Gas Reservoirs considering Uneven Support of Fractures

In shale gas reservoirs, the production data analysis method is widely used to invert reservoir and fracture parameter, and productivity prediction. Compared with numerical models and semianalytical models, which have high computational cost, the analytical model is mostly used in the production data analysis method to characterize the complex fracture network formed after fracturing. However, most of the current calculation models ignore the uneven support of fractures, and most of them use a single supported fracture model to describe the flow characteristics, which magnifies the role of supported fracture to a certain extent. Therefore, in this study, firstly, the fractures are divided into supported fractures and unsupported fractures. According to the near-well supported fractures and far-well unsupported fractures, the SRV zone is divided into outer SRV and inner SRV. The four areas are characterized by different seepage models, and the analytical solutions of the models are obtained by Laplace transform and inverse transform. Secondly, the material balance pseudotime is introduced to process the production data under the conditions of variable production and variable pressure. The double logarithmic curves of normalized production rate, rate integration, the derivative of the integration, and material balance pseudotime are established, and the parameters are interpreted by fitting the theoretical curve to the measured data. Then, the accuracy of the method is verified by comparison the parameter interpretation results with well test results, and the influence of parameters such as the half-length and permeability of supported and unsupported fractures on gas production is analyzed. Finally, the proposed method is applied to four field cases in southwest China. This paper mainly establishes an analytical method for parameter interpretation after hydraulic fracturing based on the production data analysis method considering the uneven support of fractures, which is of great significance for understanding the mechanism of fracturing stimulation, optimization of fracturing parameters, and gas production forecast.


Introduction
China is rich in shale gas resources and has a broad prospect of gas reserves and production. Shale formations are so tight that fracturing of horizontal wells is the key technique for their efficient development. And the parameter interpretation of stimulated reservoir and hydraulic fractures is the key to hydraulic fracturing evaluation [1].
Parameter inversion methods can be divided into pressure dynamic inversion (well testing analysis method) and production dynamic inversion (production data analysis method). The shale formations are so tight that the well test requires a long time to shut in the well and measure the pressure, which seriously affects production [2][3][4][5]. The production data analysis method is similar to the well testing analysis method. It is based on seepage theory, introduces the material balance time and normalized rate to process production and pressure data, and inverts the reservoir and fracture parameters via theoretical curve and measured curve fitting [6][7][8][9][10][11][12]. The production data analysis method does not require closing the well. Through comprehensive analysis of the daily production data, fracture network parameters equivalent to the well test data can be obtained, which cannot only reduce production loss but can also make full use of the production data accumulated during the development process.
After fracturing, parameters such as length, azimuth, and dip angle of each fracture are different in the stimulated reservoir volume (SRV) [13]; so, it is impossible to accurately describe each fracture. The first priority of parameter inversion is to describe the complex fracture network and establish a corresponding seepage model to analyze the production or pressure of fractured horizontal wells [14,15]. According to the solution method, seepage models can be classified into three types of numerical solution model, semianalytical solution model, and analytical solution model. The stimulated reservoir area is divided into many small grids by the numerical model, and then the finite volume, finite element, and other numerical methods are used to solve the model. Discrete fracture model (DFM) [16][17][18][19] and embedded discrete fracture model (EDFM) [20][21][22][23] are the most widely used numerical models. The DFM deals with the fractures by the local mesh refinement approach. The large-scale fractures are extracted individually and dealt with discrete media when multiscale fractures are simulated in the EDFM. In most semianalytical models of fractured horizontal wells, source function models are used to treat fractures as linear sources, and the point source function is used to solve the model. Then, the superposition principle is applied to deal with multistage fractures [24][25][26][27][28][29]. However, when a large number of fractures are generated, the numerical models and semianalytical models require a large amount of computation and have low computational efficiency. At present, most parameter inversion methods based on production data analysis mostly adopt analytical models. In most analytical models, the SRV are divided into different areas to characterize the heterogeneity of the fractured reservoir and solved through coupling of the production and pressure. Bello and Wattenbarger [30] used the fracture matrix dual model to approximate the productivity of multistage fractured horizontal shale gas wells. Based on the dual media model, Alahmadi and Wattenbarger [31] further considered the secondary fractures in the reservoir and established a three-hole model. Brown et al. [32] considered the contribution of the unreconstructed area outside of the reformed area to the production and established a classic trilinear flow model. Su et al. [33] regarded the flow perpendicular to the hydraulic fracture as compound flow with internal double holes and homogeneous external areas, and they extended the three linear flow model to the four linear flow model. Stalgorova and Louis [34] further considered the unreconstructed area between the fractures, divided the model into five areas, and established the five linear flow model. Zeng et al. [35] extended the five-linear flow model and established the seven linear flow model. These analytical models take into account the heterogeneity of the reservoir after hydraulic fracturing, but ignore the uneven support of fractures. The length and conductivity of the fractures obtained based on these seepage models are the equivalent parameters of supported fractures and unsupported fractures. In fact, affected by factors such as ground stress, rock strength, and the degree of development of natural fractures, various fractures including proppant-filled fractures and unfilled fractures are formed after hydraulic fracturing [36,37]. Previous studies on the conductivity of propped fractures by means of experiments and numerical simulations have shown that the conductivity of propped fractures plays a major role in the productivity of gas wells [38,39]. However, an increasing number of studies have shown that unpropped fractures have conductivity even after closure, which also plays an important role in productivity [40]. And the more complex the fracture network is, the smaller the average proppant concentration in the fracture, and the majority of fractures are unpropped.
In this study, the complex fracture network after fracturing is equivalent to the SRV area, and the SRV area is divided into the outer SRV and inner SRV based on the uneven support of the fractures. First, a seepage model considering the shale gas adsorption and gas physical character nonlinearity is established. The analytical solution of the model is obtained by Laplace transform and Stehfest numerical inversion. Then, a corresponding interpretation method of key parameters after fracturing based on production data analysis method is established. Finally, the proposed method was applied to four field examples which is provided to demonstrate the application of the proposed analytical method for parameter interpretation, and six parameters were obtained, including the half-lengths and permeabilities of the supported and unsupported fractures and the permeabilities of the outer and inner SRV, which can quantitatively distinguish between supported fractures and unsupported fractures.

Methodology
2.1. Physical Model. Due to the limitations of fracturing technology, it is difficult to evenly distribute the proppant in the fractures. Current studies believe that the near-well zone is the main fracture, which is filled with proppant, and most of the fractures in the far-well zone are not filled with proppant [41]. For shale gas wells in the early stage of production, the unreconstructed area outside of the stimulated area has a small contribution to gas production. In order to reduce the computation cost of parameter inversion and the multiplicity of solutions, the model in this study does not consider the unreconstructed region. According to the distribution of proppants in the hydraulic fracture, we divide the fractures into supported fracture, which is filled with proppant and unsupported fracture, which is unfilled with proppant. And the SRV area is divided into the outer SRV and inner SRV based on the near-well supported fractures and the far-well unsupported fractures (Figure 1), and we assume that the flow in each area is linear. The other assumptions are as follows. (1) The reservoir is closed vertically, and the horizontal well is located in the center of the reservoir. (2) The gas flows in a single phase and in an unsteady state within the stimulated reservoir and fractures.

Diffusivity Equation in SRV and Fractures
2.2.1. Outer SRV. In the outer SRV, the continuity equation is given below: In Equation (1), q a is desorption rate of adsorbed gas per unit volume per unit time, v om is the gas velocity, ρ is gas density, and they are defined as Volume factor is defined as The adsorption model obeys the Langmuir isotherm: Pseudopressure is defined as There is not only free gas but also much absorbed gas in the shale reservoirs. c tom is the modified compressibility defined by Gerami et al. [42] considering the desorption

Dimensionless variables Definition
Dimensionless rate at constant pressure 3 Geofluids effect, and it is defined as follows: In Equation (7), c om is the rock compressibility, c g is the gas compressibility, and c d is the desorption compressibility. Their definition is as follows: Substituting Equations (2)-(11) into Equation (1) and adopt the engineering unit system, to arrive at The fluid in the outer SRV flows into the unsupported fracture, and the boundary conditions are   1E-04  The boundary conditions are 2.2.3. Unsupported Fracture. Shale gas flow linearly from the outer SRV into the unsupported fracture and the continuity equation of the unsupported fracture are given below: Velocity is defined as Equation (14) is the mass of the gas flowing into the unsupported fracture from per unit volume of the outer SRV per unit time, and its definition is Substituting Equations (4), (7), (15), and (16) into Equation (14), the diffusivity equation becomes The gas flow from unsupported fracture into supported fracture and the unsupported and supported fractures are coupled through the rate and pressure. Therefore, the boundary conditions are as follows: Similarly, the diffusivity equation in the supported fracture can be obtained as follows: The horizontal well is produced at constant pressure, and the boundary conditions are as follows: When the shale gas fractured horizontal well is produced at constant pressure, the productivity is 2.3. Dimensionless Mathematical Models. The dimensionless mathematical models of the outer and inner SRV and unsupported and supported fractures are shown in Equations (24)- (27), and the dimensionless parameters are shown in Table 1.

Mathematical Model Linearization. As in Equation
By substituting Equation (24) into Equation (10), the diffusivity equation becomes The dimensionless pseudotime is defined as Therefore, the linearized dimensionless outer SRV math-ematical model becomes The linearized dimensionless inner SRV model is 3. Inversion Method for the Fracture Network Parameters 3.1. Derivation of the Mathematical Model. In order to facilitate the derivation, the model was transformed using the Laplace transform, and the solution was obtained in the Laplace domain [47,48]. The solution of the dimensionless rate in the Laplace domain at constant bottom hole flow pressure is α om , α Im , α oF , andα IF are defined as follows: Furthermore, the Stehfest numerical inversion method [49] was used to obtain the solution in the time domain.

Parameter
Interpretation. The general procedure of the production data analysis method is as follows: firstly, a 6 Geofluids seepage model in the SRV area is established according to actual situation of reservoir, and then the production data is fitted based on the established seepage model to obtain parameters of stimulated reservoir and fractures. During the production of shale gas wells, the gas production rate and bottom hole pressure are constantly changing. However, the model proposed above was established under the condition of a constant flowing bottom hole pressure; so, the actual rate curve and the theoretical curve cannot be directly fitted to interpret the parameters [50][51][52]. Blasingame et al. [53] introduced the material balance pseudotime to establish the equivalent relationship between the constant and the variable production rate and found that the curve for the constant pressure case and the curve for the constant rate case overlap after the material balance time conversion. Therefore, the material balance pseudotime can be used to convert the variable rate and the variable pressure production curves into constant pressure curves for parameter interpretation. We substituted the material balance pseudotime for the actual time, substituted the normalized rate for the actual rate, and established data processing methods for the variable rate or pressure cases.
The material pseudotime is defined as follows:    Figure 6: Effects of the permeability of supported fractures on gas production rate (a) and cumulative production (b).

Geofluids
In Equation (15), Z * is the modified gas compressibility factor considering the influence of desorption, which is defined as follows: The pseudopressure is used to normalize the rate: In order to reduce the influence of the noise of the production data and to reduce the nonuniqueness of the     Figure 7: Effects of the permeability of outer SRV on gas production rate (a) and cumulative production (b). 8 Geofluids solution, the dimensionless rate integration and the derivative of the integration are applied, which are defined as follows: Establish the chart of the relationship between normalized production rate, rate integration, and the derivative of the integration and material balance pseudotime. The parameters are interpreted by curve fitting. The theoretical curve is drawn according to the proposed mathematical model. The production data are processed, and then the measured curve is drawn.

Model Validation. A shale gas well YY-1 in southwestern
China was taken as an example for verification. The well test analysis method assumes that the fractures are all supported fractures. In order to compare the parameters interpreted by the proposed method with the well test method, the uneven support of fractures is ignored, and the SRV model proposed in Figure 1 is simplified (Figure 2). YY-1 horizontal well section is 1502 m long, with 63 clusters of 23 fracturing stages. According to the gas testing, the initial formation pressure is 75 MPa. The thickness of high-quality shale is 38 m. Water saturation is 56%, matrix porosity is 5.15%, and pore compressibility is 3 × 10 − 4 MPa -1 . Langmuir pressure is 6.5 MPa, and Langmuir volume is 3 m 3 /t obtained by the adsorption isotherm experiment. The proposed parameter interpretation method is used to analyze the gas production rate and flowing bottom hole pressure of well YY-1 (Figure 3), and three parameters were obtained through  Figure 10: Effects of the half-length of supported fractures on gas production rate (a) and cumulative production (b).
q Di ⨯100, q Did /100 t D Figure 11: Effects of the x IF/x F on normalized gas production. The fitting curve is shown in Figure 4, in which the abscissa is the pseudotime of the material balance, and the ordinate is the normalized rate, rate integral, and integral derivative. By comparing the inversion parameters with the values obtained from well testing (Table 2), it can be seen that the half-length and permeability of the fracture and the permeability of the stimulated reservoir obtained through inversion are in good agreement with the values obtained from well testing, and the simplified model fits the production data well. Therefore, the parameter   10 Geofluids interpretation method proposed in this study can be applied to the quantitative characterization of parameters after hydraulic fracturing in shale gas reservoirs.

4.2.
Effect of the Uneven Support of Fractures. The halflength and conductivity of hydraulic fractures and permeability of the SRV area are the key parameters for evaluation of fracturing result, which also have a significant impact on shale gas productivity forecast. In this section, we analyzed the influence of six parameters, including the half-lengths and permeabilities of the supported and unsupported fractures and the permeabilities of the outer and inner SRV, on the average gas production rate in the first year and the cumulative gas production in twenty years.
The variable values are shown in Table 3. A single factor analysis method is used to analyze the influence of each parameter on the production. The default value of halflengths and permeabilities of the supported and unsupported fractures and the permeabilities of the outer and inner SRV are 100 m, 30 m, 10 mD, 1000 mD, 2 × 10 −5 mD, and 1 × 10 −4 mD, respectively.
The effects of six parameters on gas production are shown in Figures 5~10. Fracture conductivity or permeability is a very important parameter in fracturing. When the permeability of unsupported fracture varies in a wide range, it has a great influence on productivity. As it increases, the gas production increases, but the increasing amplitude slows down ( Figure 5). As permeability of supported fracture increases, gas production increases with small gradient. This is because a given minimum supported fracture permeability is sufficient for fluid transport (Figure 6). An important purpose of fracturing stimulation is to increase the extent and permeability of the SRV. It can be seen from Figures 7 and  8, with the increase of the permeability of the outer SRV or inner SRV, the gas production increased, but the increase rate tends to be reduced. The half-length of fracture determines the range size of the SRV area, which has a great influence on productivity. As shown in Figures 9 and 10, along with the increase of the half-length of supported or unsupported fracture, the gas production rate and cumulative production increase. In particular, gas production increases almost linearly as the half-length of the supported fracture increases.
In order to further study the effect of the half-length of supported and unsupported fractures on gas production, the half-length of fracture, which is the sum of the halflength of supported and unsupported fractures, is intro-duced, and the logarithmic curve of normalized production rate, rate integration, and the derivative of the integration and material balance pseudotime is drawn. Assuming that the half-length of the fracture is 100 m, the half-length of supported fracture is set at 20%, 40%, 60%, 80% and 100% of the half-length of the fracture, respectively, in Figure 11. Assuming that the half-length of the supported fracture accounts for 50% of the half-length of the fracture, the half-length of the fracture is set at 80 m, 100 m, 120 m, 130 m, and 150 m, respectively, in Figure 12. As can be seen from the logarithmic curve, the fracture linear flow, matrix fracture bilinear flow, matrix linear flow, transition flow, and system boundary control flow appear from left to right. The value of x IF /x F mainly affects the matrix fracture bilinear flow stage, the matrix linear flow stage, and the boundary control flow stage. It mainly affects the occurrence time of the boundary control flow stage and has a greater impact on the gas production. When x IF /x F is constant, the value of x F affects the matrix linear flow stage and the boundary control flow stage, mainly affecting the appearance time of the matrix linear flow stage, and also has a great influence on the production.

Field
Examples. The stimulated reservoir and fracture parameter inversion method proposed in this study, which considers the uneven support of fractures, are used to fit the production data of four fractured horizontal wells in Southwest China. In order to illustrate the difference in fitting results of different models to the same production data, this section compares the two models: considering and  11 Geofluids ignoring the uneven support of fractures. The gas production and bottom hole flow pressure data of the four wells are shown in Figure 13, the production fitting results are shown in Figures 14~17, and the parameter inversion results are shown in Tables 4 and 5.
It should be noted that when the uneven support of fractures is ignored, the integral derivative data of production is too noisy; so, only the normalized production and integral derivative curves of production are fitted. As we can see from the fitting curve ( Figures 14~17), when the production time is more than 100 days, except for the well YY-2 without considering the uneven support of fractures, the two models can obtain better fitting results in the matrix linear flow and boundary control flow stages. The fitting effect of production data in the early stage of fracture linear flow is not good, which is mainly related to the flowback of fracturing fluid in the fracture at the initial stage of production. However, the model considering the uneven support of fracture has a better fitting effect. Comparing the parameter inversion results in Tables 4 and 5, it can be seen that the inverted fracture parameters are the equivalent parameters of supported fractures and unsupported fractures when the uneven support of fractures is not considered. Although different models can fit the same production data, the fitting effects are different, and the key parameters of the inversion are also different. In order to accurately invert reservoir and fracture parameters and further understand the mechanism of fracturing stimulation measures, it is necessary to consider the influence of the uneven support of fractures in the production analysis of shale gas fractured horizontal wells.  Figure 15: Type curve fitting of YY-3 production data considering (a) and ignoring (b) the uneven support of fractures.  Figure 16: Type curve fitting of YY-4 production data considering (a) and ignoring (b) the uneven support of fractures.

Conclusion
In this study, a simple and general inversion method of reservoir and fracture parameters is established based on the production data analysis method considering the uneven support of fractures. We can draw some conclusions as follows: (1) Well test analysis and production data analysis are the two most commonly used method for fracture network parameter inversion. The production data analysis method can be better used for parameter interpretation because of its easy data collection and simple process in shale gas reservoirs (2) Six parameters can be interpreted by curve fitting, including the half-lengths and permeabilities of the supported and unsupported fractures and the permeabilities of the outer and inner SRV, and thus the method proposed can quantitatively distinguish between supported fractures and unsupported fractures (3) In the parameter inversion, the production data can be well fitted either with or without uneven support of fractures, but the model considering the uneven support of fracture has a better fitting effect. If the uneven support of fracture is not considered, the predicted output in the early stage of production will be much higher than the field data (4) The pseudotime and pseudopressure can deal with the nonlinearity of gas physical characteristics in the model. And a method for processing production data with variable production and pressure can be established by introducing material balance pseudo time and normalized production rate Nomenclature B g : Formation volume factor, m 3 /m 3 h: Effective formation thickness, m x oF : Half-length of the unsupported fracture, m x IF : Half-length of the supported fracture, m k om : Permeability of the outer SRV, md k Im : Permeability of the inner SRV, md k oF : Permeability of the unsupported fracture, md k IF : Permeability of the supported fracture, md C m : Rock compressibility, MPa -1 C g : Gas compressibility, MPa -1 C d : Desorption compressibility, MPa -1 C t : Total compressibility, MPa -1 C tom : Total compressibility of the outer SRV, MPa -1 C tIm : Total compressibility of the inner SRV, MPa -1 C toF : Total compressibility of the unsupported fracture, MPa -1 C tIF : Total compressibility of the supported fracture, MPa -1 Z: Gas compressibility factor, dimensionless p om : Pressure in the outer SRV, MPa p Im : Pressure in the inner SRV, MPa p sc : Pressure at standard condition, MPa p i : Initial pressure, MPa p L : Langmuir pressure, MPa q sc : Wellbore gas rate, m 3 /day q N : Normalized rate, m 3 /day/MPa 2 (mPa ⋅ s) t: Production time, day t a : Pseudoproduction time, day t ca : Material balance time, day T: Formation temperature, K T sc : Temperature at standard condition, K V: Adsorbed gas volume, m 3 V L : Langmuir volume, m 3 /m 3 ψ om : Pseudopressure in the outer SRV, MPa 2 /mPa ⋅ s  Figure 17: Type curve fitting of YY-5 production data considering (a) and ignoring (b) the uneven support of fractures.
13 Geofluids ψ oF : Pseudopressure in the unsupported fracture, MPa 2 /mPa ⋅ s ψ Im : Pseudopressure in the inner SRV, MPa 2 /mPa ⋅ s ψ IF : Pseudopressure in the supported fracture, MPa 2 /mPa ⋅ s ω oF : The width of the unsupported fracture, m ω IF : The width of the supported fracture, m ϕ om : Porosity of the outer SRV, dimensionless ϕ oF : Porosity of the unsupported fracture, dimensionless ϕ Im : Porosity of the inner SRV, dimensionless ϕ IF : Porosity of the supported fracture, dimensionless ρ: Gas density, kg/m 3 μ g : Gas viscosity, mPa•s.

Data Availability
The data used to support the findings of this study are included within the article.

Conflicts of Interest
The authors declare that they have no conflicts of interest.