Fault Slip Risk Assessment and Treating Parameter Optimization for Casing Deformation Prevention: A Case Study in the Sichuan Basin

College of Petroleum Engineering, China University of Petroleum (Beijing), Beijing 102249, China Drilling Technology Institute, CNPC Engineering Technology R&D Company Limited, Beijing 102206, China Key Laboratory of Geotechnical & Underground Engineering of Ministry of Education, Tongji University, Shanghai 200092, China Department of Geotechnical Engineering, College of Civil Engineering, Tongji University, Shanghai 200092, China Shale Gas Research Institute, PetroChina Southwest Oil & Gasfield Company, Chengdu 610051, China Chongqing Gas Field, PetroChina Southwest Oil & Gasfield Company, Chongqing 400707, China


Introduction
The Changning-Weiyuan national shale gas demonstration area is located in the Sichuan Basin, southwest China. Since 2009, horizontal drilling and hydraulic fracturing have been used to develop shale gas in this block. However, serious casing deformation problems occurred during hydraulic fracturing operations. As of December 2018, among 113 wells that have completed fracturing in the Changning block, casing deformation occurred in 39 wells, accounting for 34%. Simi-larly, as of October 2018, among 74 wells that have completed fracturing in the Weiyuan block, casing deformation occurred in 36 wells, accounting for 48%. As a result, casing deformation has become a key engineering problem that needs to be resolved urgently for shale gas development in the Changning-Weiyuan block. Not only does casing deformation decrease the number of viable fracturing stages but also reduces the overall life cycle of the well, which in turn seriously restricts the efficient development of shale gas resources.
In the Changning-Weiyuan block, cluster perforation and bridge plug techniques are used for multistage fracturing. The drifting gauge can pass the well smoothly before fracturing; however, the bridge plug is often unable to pass through after fracturing. Therefore, it can be inferred that casing deformation occurred during hydraulic fracturing operations. Fracturing is the engineering factor for causing casing deformation. The complex geological conditions in the Sichuan Basin have resulted in extensive faults and fractures, which are prone to cause casing shear deformation due to fault slip [1][2][3][4][5][6][7][8]. The main mechanism is believed to be that the reduction of effective stress caused by the hydraulic fracturing fluid leads to fault slip. If the wellbore crosses faults, casing deformation occurs due to fault slip.
Many scholars have explored how to mitigate casing deformation from different aspects. Xi et al. [9] established the relationships among microseismic moment magnitude, slip distance, and the reduction of casing inner diameter. Their results show that keeping the designed horizontal segment of the well trajectory either away from the fracturedeveloped area or parallel to natural fractures can decrease the fault slip distance. Yin et al. [10] developed a 3D nonlinear finite element model to simulate the mechanical behavior of casing crossing slip deformation. Simulation results demonstrate that decreasing the crossing angle between faults and casing, employing low Young's modulus cement or without cementing are the effective countermeasures to prevent casing failure during hydraulic fracturing. Liu et al. [11] put forward a semianalytical model to calculate the fault slip distance and suggested that the wells drilled through the fault should be stimulated below a certain net pressure to avoid slippage while stimulating nearby stages. Zhang et al. [12] used the discrete element method to simulate the fault slip induced by hydraulic fracturing. The simulation results show that lowing pumping rate and fluid viscosity can effectively decrease the shear displacement of the fault plane and help mitigate casing deformation. These studies have documented a lot of valuable suggestions to mitigate casing deformation, including optimizing the well trajectory, reducing the elastic modulus of cement sheath, no cementing, and reducing treating parameters. However, firstly, some measures are limited to theoretical analyses and lack verification by field practices. The feasibility and rationality need to be further demonstrated. Secondly, the proposed measures generally are incapable of predicting casing deformation, which makes the process from prediction to control incomplete. Consequently, limitations still exist by applying these measures in the field.
In this study, we first statistically analyze casing deformation in the Ning209 block and investigate the mitigation effect of reducing treating parameters on casing deformation. Then, the geological structure and in-situ stress characteristics are analyzed to explore the root cause of casing deformation concentration. Finally, based on the fault likelihood and 3D in-situ stress model, a geomechanical probability model of fault slip is established to identify the areas with high slip risk to explain the distribution of casing deformation and evaluate the impact of reducing pumping rate on the faults with different slip risks.

Engineering Background
In this study, we take the Ning209 block in the Changning-Weiyuan shale gas development area as the research object. Stimulation has been conducted in 13 pads of the Ning209 block so far. According to the statistics of issues during running bridge plug, a total of 24 casing deformation positions are identified. The spatial distribution of casing deformation is shown in Figure 1. Clearly, casing deformation locations are mainly concentrated in the four pads of H2, H4, H6, and H29 of the upper region. In the lower region, except for the H10 pad, almost no casing deformation occurred on the other pads. It can be seen that the spatial distribution of casing deformation shows clear characteristics of concentration in the Ning209 block.
Different from the typical hydraulic fracturing conducted in other blocks, treating parameter optimization was carried out in the Ning209 block to mitigate casing deformation. We use the fracturing operation design of Ning209 block to obtain the optimized parameters. These data determined how much pumping rate or fluid volume would be used for each stage during field fracturing operations. From the fracturing operation design, all stages near the faults were taken as risky stages for casing deformation so that measures such as reducing the pumping rate by 2 m 3 /min (from 16 m 3 /min to 14 m 3 /min) or the fluid volume by 300 m 3 (from 1800 m 3 to 1500 m 3 ) would be taken.
The summary of treating parameter optimization in the Ning209 block is provided in Table 1. 13 pads had reduced pumping rate, and 3 pads had reduced both pumping rate and fluid volume. The statistical results show that there are a total of 190 faults which intersect the treating parameters optimization stages in the Ning209 block, and 24 casing deformation occurred in total, with a casing deformation ratio of 12%. To compare, the Ning201 block, where 256 faults intersect wells, was not stimulated with optimization, and 57 casing deformation occurred in total, with a ratio of 22%. To date, reducing the pumping rate as the main optimization measure has indeed played a significant role in mitigating casing deformation. In addition, we find that the effect of reducing pumping rate on casing deformation in the lower region is more profound than in the upper region, as evidenced from Figure 1

Observation of Faults
In order to understand the concentration of casing deformation, the correlation between the faults and casing deformation is first explored. The fault likelihood is an attribute used to describe the distribution of faults based on seismic data. This method transforms the seismic data into a volume of local-fault-extraction (LFE) estimates that represents the likelihood that a given point lies on a fault surface [13,14]. It can not only accurately reflect the development characteristics of large-scale faults but also has a strong ability to distinguish the small-scale faults. The faults interpreted by fault likelihood in the Ning209 block are shown in Figure 2. And we statistics the position of casing deformation and the strike of fault nearest to 2 Geofluids casing deformation as shown in Table 2. We can see that most of the deformed locations agree with fault locations very well. However, the fault intensity shows no noticeable difference between the upper and the lower region, which could not explain the concentration of casing deformation yet. Figure 3 shows the seismic profile of fault systems in the Ning209 block. The lines with arrow represent the faults. From top to bottom, the red line represents the top of the Maokou formation, the pink line represents the bottom of the Liangshan formation, and the green line represents the bottom of the Longmaxi shale formation, which is the main pay zone. We can see that the faults developed in this block are mainly high angle faults that cross the Longmaxi formation.
Earlier studies have shown that the in-situ stress is one of the dominant factors controlling fault activation [15][16][17][18][19]. It is thus necessary to analyze the distribution of in-situ stress to explore whether a relationship exists between the in-situ stress and the concentration of casing deformation.

The 3D Geomechanical Model
4.1. Model Construction. By integrating seismic, geological structure, logs, and core data, a 3D geomechanical model is built for the Ning209 block. The 3D geomechanical model includes 3D anisotropic mechanical properties, 3D pore pressure, and 3D in-situ stress field. The modeling process is shown in Figure 4, including the following steps: (i) Based on the laboratory experiment and measured data, a 1D geomechanical model is established by using the logging data of acoustic wave, density, and gamma ray (ii) On the basis of geological modeling, a 3D mechanical property model is developed by combining well and seismic data (iii) The 3D pore pressure prediction is conducted with seismic data and calibrated against pressure measurements, mud-logging data, and flowback data (iv) Taking strain as the boundary condition, the in-situ stress interpreted by well data is fitted to obtain the 3D in-situ stress field It has been widely reported that the mechanical properties of shale are anisotropic [20][21][22][23]. Recent studies show that the Longmaxi shale in the Sichuan Basin exhibits clear anisotropy [24][25][26][27]. To build a reliable geomechanical model, it is necessary to capture the anisotropy of the formation. This anisotropy is reflected in the differences of mechanical properties in different directions, including the differences of Young's modulus and Poisson's ratio in the transverse and longitudinal direction. The anisotropy of the formation is captured at wellbores and propagated to the 3D space guided by seismic data as shown in Figure 5. The vertical depth of the shale reservoir is about 3300 m. At this depth, the Young's modulus varies in the range of 30~70 GPa, and the Poisson's ratio varies in the range of 0.21~0.31.
Based on the 3D anisotropic mechanical properties model, finite element simulation is carried out to construct the 3D in-situ stress field by taking strain as the boundary condition. Strain can be determined through the poroelastic horizontal strain equation by applying two constraints simultaneously. First, ensure that the minimum horizontal stress matches the collected closure-pressure data with the given set of strain. Second, ensure that the borehole failure (both shear failure and tensile failure) predicted by the 1D  4 Geofluids geomechanical model is consistent with borehole failure observed from the caliper and image logs [28][29][30].
where σ h is the minimum horizontal stress (in MPa); E h and E v are Young's moduli in the horizontal and vertical directions, respectively (in GPa); υ h and υ v are Poisson's ratio in the horizontal and vertical directions, respectively; σ v is the overburden stress (in MPa); α is the Biot coefficient; P p is the pore pressure (in MPa); and ε h and ε H are the horizontal strain in minimum and maximum horizontal stress direction, respectively. The 3D in-situ stress field of the Ning209 block is shown in Figure 6. The vertical stress distribution is approximately in the 50~90 MPa range (Figure 6(a)), the maximum horizontal stress is approximately in the 50~110 MPa range ( Figure 6(b)), and the minimum horizontal stress is approximately in the 40~80 MPa range (Figure 6(c)). The stress regime over the region is mainly strike-slip faulting. The azimuth of the maximum horizontal stress is approximately 109°N.

The Relationship between Horizontal Stress Difference and Casing Deformation.
In the classical Mohr-Coulomb failure theory, a larger stress difference will enlarge the Mohr circle which describes the regional stress state and then cause the stress state of the fault to be closer to the critical stress state, which increases fault slip potential under hydraulic fracturing. The distribution of horizontal stress difference in the Ning209 block is displayed in Figure 7. It can be seen that the stress difference in the upper region is significantly greater than that in the lower region, which indicates that the faults in the upper region have higher potential to slip. Casing deformation can be caused by fault slip, so the greater the stress difference, the more severe the casing deformation. From Figure 7, most of the casing deformation is concentrated in the region with relatively large stress difference. The distribution characteristics of stress difference are consistent with the actual casing deformation, which explains why casing deformation is concentrated in the upper region to some extent.

The Relationship between Pore Pressure and Casing
Deformation. Pore pressure is another key factor affecting fault slip. Based on the effective stress principle, the effective normal stress of a fault plane would be at a lower level under a higher pore pressure, which leads to the weakening of fault shear strength. The pore pressure distribution in the Ning209 block is shown in Figure 8, which is similar to the characteristics of stress difference distribution. The abnormal high pore pressure is concentrated in the upper region. Consequently, under the combined effects of high pore pressure and high stress difference, the faults in the upper region have higher slip   The characteristics of in-situ stress distribution in the Ning209 block show that the high stress difference and high pore pressure distribution might be the main factors causing the concentration of casing deformation. In the next section, we use geomechanics probability analysis to evaluate the slip risk of faults and further investigate the distribution characteristics of casing deformation.

Quantitative Risk Assessment Method.
In the context of Coulomb faulting theory, whether or not a fault will slip in response to fluid injection depends on several factors, including the magnitude and orientation of the stress field and the orientation of the fault, pore pressure, and material parame-ters such as the coefficient of friction. Because there is uncertainty in each of these parameters, Quantitative Risk Assessment (QRA) can be used to evaluate the slip risk of faults in response to a given pore pressure increase. QRA is a Monte Carlo method used to evaluate the probability of an uncertain outcome when the input parameters are also uncertain. In QRA, uncertainties in input parameters are populated by iteratively running the model using inputs randomly sampled from probability density functions. The output is a cumulative distribution function describing the cumulative probability of fault slip as a function of pore pressure applied to the fault [31][32][33]. The risk assessment process of fault slip is demonstrated in Figure 9, including the following steps:    9 Geofluids probability is the frequency of the fault stress state exceeds the failure line in the probability model. Essentially, it is the relationship between the stress state of the fault and the failure line on Mohr's diagram. In the probability model, each input parameter has a central value within an error range; thus, there is also a central value describing the stress state of the fault in simulation results, which determines the location of faults with different slip risks on Mohr's diagram. If the slip probability of the fault is 50%, it means that the cases of slip and nonslip in the simulation results account for 50% each; then, the fault on Mohr's circle corresponding to the exact position meets the critical stress state, as shown by the yellow dot in Figure 10. Therefore, it is not difficult to see that the degree of slip risk is controlled by the distance between the central value of stress state and the failure line. Taking the intersection point of the Mohr's circle and the failure line as the center, it can be divided into three regions: high risk (slip probability > 70%), medium risk (70% > slip probability > 30%), and low risk (30% > slip probability), as illustrated in Figure 10.

Geomechanical Probability Analysis.
Before performing geomechanical probability analysis, the fault model and insitu stress model need to be determined. The in-situ stress model has been established as shown in Figure 6. The fault model is established based on the fault likelihood as shown in Figure 11.
The panels containing the red line in Figure 12 show the distribution of geomechanical input parameters of fault #1 that is marked in Figure 11. The results are generated by random sampling on uniformly distributed (red line) samples. The bottom right panel shows the distribution of pore pressure perturbation values required for slip on fault #1 in the generated probability model. It can be seen that the proportion of the perturbation values reaching the highest is in the range of 16-22 MPa. It is worth noting that this distribution describes the minimum perturbation condition required for fault slip. That is when the perturbation value gradually increases to [16][17][18][19][20][21][22] MPa, the fault begins to have a larger slip possibility. This possibility will increase with the increase of perturbation value.
The cumulative distribution function (CDF) curve can be obtained by integrating the distribution of pore pressure perturbation value as shown in Figure 13. This curve directly describes the slip tendency probability of fault #1 with the increase of pore pressure. Based on this curve and pore pressure perturbation value, the slip probability of fault #1 can be evaluated.
Before probability assessment, it is notable that compared with typical hydraulic fracturing, the Ning209 block has taken measures to reduce pumping rate, which will inevitably impact the pore pressure perturbation value induced by hydraulic fracturing. Thus, it is necessary to consider two cases. The first case is the typical hydraulic fracturing, which gives the distribution of fault slip risk and explains the concentration of casing deformation. The second is the case with reduced pumping rate, which shows the influence of reducing pumping rate on fault slip risk.
The range of pore pressure perturbation value after fracturing can be approximately determined by the difference between the instantaneous shut-in pressure (ISIP) and initial pore pressure. The ISIP values in this study were obtained by the treating pressure records per second in field fracturing data [17]. We used the standard procedure in the industry to obtain the ISIP values. The ISIP value is the treating pressure corresponding to the pumping rate reduced to 0 in the treating pressure records per second. Figure 14 shows the comparison of treating pressure corresponding to different pumping rate. We can see that the pressure difference between typical hydraulic fracturing and reducing pumping rate is approximately 4 MPa.
The pore pressure perturbation value induced by typical hydraulic fracturing is approximately 25 MPa, and the slip probability of fault #1 is about 71%, as shown in Figure 13. The pore pressure perturbation value induced by hydraulic fracturing with reduced pumping rate is approximately 21 MPa, and the slip probability of fault #1 is about 50%. It can be seen that fracturing with reduced pumping rate decreases the slip probability by 20% for fault #1 compared with fracturing with typical pumping rate.
The accuracy of the probabilistic analysis may be verified by field microseismic data. Well H2-3 intersects with fault #1, as shown in Figure 15(a), where microseismic monitoring has been carried out. The microseismic events of stages 5, 8, and 9 are shown in Figure 15(b). Due to casing deformation, stage 6 and stage 7 were abandoned. Figure 16 shows the treating pressure of these stages. We observe that the microseismic events from these three stages are all far away from the perforation locations. The multistage events overlap and converge into a microseismic event band, which conforms to the characteristics of fault activation [34][35][36][37][38]. The microseismic event band coincides with fault #1 spatially, and its extension direction is also consistent with the strike of fault #1. Therefore, it could be considered as the same fault. Microseismic data indicates fault #1 activation during hydraulic fracturing, which shows that the probability result is relatively reliable.

Results
The slip probabilities of faults under typical hydraulic fracturing (pore pressure perturbation value of 25 MPa) in the Ning209 block are shown in Figure 17. A total of 160 faults in the upper region and 43 faults were in high-risk state, accounting for 27%; 76 faults were in medium-risk state, accounting for 47%; 41 faults were in low-risk state, accounting for 26%. A total of 198 faults in the lower region and 101 faults were in medium-risk state, accounting for 51%; 97 faults were in low-risk state, accounting for 49%. It can be seen that under high stress difference and high pore pressure conditions, the slip risks of faults in the upper region are clearly higher than that in the lower region, which is the reason that casing deformation locations are mainly concentrated in the upper region. 10 Geofluids Specifically, casing deformation mainly occurred near the medium-risk and high-risk faults in the upper region; thus, it can be concluded that the distribution of casing deformation is mainly dominated by the medium-risk and high-risk faults.
The slip probabilities of faults under reduced pumping rate (pore pressure perturbation value of 21 MPa) in the Ning209 block are provided in Figure 18. For the upper region, 112 faults were in medium-risk state, accounting for 70%; 48 faults were in low-risk state, accounting for 30%. We see that most of the faults in the upper region are changed from high-risk state to medium-risk state when fracturing using reduced pumping rate. For the lower region, 52 faults were in medium-risk state, accounting for 26%; 146 faults were in low-risk state, accounting for 74%. The faults in the lower region are close to the low-risk state, and the slip probabilities of faults are relatively low, which is the reason for the reduction of casing deformation in the lower region. Hence, it can be concluded that reducing pumping rate has a significant effect on mitigating the reactivation of medium-risk faults.
The influence of reducing pumping rate on fault slip risk can be explained from the stress state of faults. Reducing pumping rate usually decreases the pore pressure

Geofluids
perturbation by a small extent. For high-risk fault, the stress state exceeds the critical stress after hydraulic fracturing ( Figure 10). With a small pressure difference, the fault might be still in medium-risk state (Figure 19(a)). For medium-risk fault, the stress state is near to the critical stress after hydraulic fracturing (Figure 10). With a small pressure difference, the fault moves to a position below the failure line and can turn into low-risk state (Figure 19(b)). Therefore, reducing pumping rate has a significant effect on medium-risk faults, while it has a little effect on high-risk faults. Since the faults are mainly medium risk in the lower region of the Ning209 block ( Figure 17), after using reduced pumping rate, casing deformation hardly occurred. Clearly, the field practice results in the Ning209 block have well proved this point.

Discussion
By exploring the reason of casing deformation distribution after fracturing with reduced pumping rate in the Ning209   12 Geofluids block, our results show that the faults in the upper region have a higher slip probability than faults in the lower region, which resulted in the concentration of casing deformation in upper region. The high stress difference and high pore pressure are the main factors leading to a higher slip probability of faults. Therefore, much attention needs to be paid to the areas with high stress difference, high pore pressure, and high slip probability faults during hydraulic fracturing, and the casing deformation prevention measures should be focused on these areas. The casing deformation prevention measures currently applied in the field include well trajectory optimization and treating parameters optimization. The fault slip risk assessment can generally predict the medium-and highrisk faults well. Hence, optimizing the well trajectory to avoid the medium-and high-risk slip faults might be the most appropriate measure for preventing casing deformation. However, the trajectory of shale gas well usually needs to drill along the direction of the minimum horizontal in-situ stress in order to create transverse fractures to maximize fracture surface area. In addition, the faults and natural fractures are extremely well developed in some blocks due to the complex geological conditions of the Sichuan Basin. Hence, there might not be enough room  Figure 14: The comparison of treating pressure corresponding to different pumping rate.

Geofluids
for operation to avoid the medium-and high-risk faults. Consequently, the well trajectory optimization still has limitations in the field of application.
Alternatively, the treating parameter optimization has been the main practice to prevent casing deformation in this field. Our results show that medium-and high-risk faults are the main factors causing casing deformation; thus, they should be the main target when optimizing treating parameters. For medium-risk faults, the slip probability can be limited at a low level via reducing pumping rate, and then, casing deformation can be mitigated, which has been proved by our field practice results from the Ning209 block. For high-risk faults, to limit the slip probability to a lower level, a large pumping rate reduction may be needed. However, this could have a serious impact on production or may not even be able to effectively stimulate the reservoir. Therefore, reducing pumping rate alone may not be appropriate for high-risk faults.

Conclusions
In this study, we take the Ning209 block as an example to analyze the characteristics of casing deformation distribution and the effect of mitigation strategy using reduced treating   parameters on casing deformation. The fault slip risk assessment is used to identify the fractured stages near the medium-and high-risk faults before hydraulic fracturing. The treating parameter optimization is employed to reasonably reduce the treating parameters of these fractured stages to mitigate casing deformation, which bridges the gap between prediction and control for casing deformation mitigation. The main conclusions can be summarized as follows.
(1) The treating parameter optimization conducted in the Ning209 block shows that reducing the pumping rate has played a significant role in mitigating casing deformation. This directly confirms that reducing treating parameters is a practical measure to mitigate casing deformation (2) Fault slip risk results show that the faults in the upper region in the block have a higher slip probability than the faults in the lower region, which is consistent with the observation that casing deformation mainly occurred in the upper region. The high stress difference and high pore pressure are the main factors causing a high slip probability of faults (3) After reducing pumping rate, slip probability was seen to change from medium risk to low risk for faults in the lower region, which led to a significant reduction of casing deformation in the lower region. Reducing pumping rate can play a significant role in mitigating casing deformation caused by mediumrisk faults, which has been proved by the field practices In summary, this paper proposes a comprehensive method for preventing casing deformation by combining the fault slip risk assessment with the treating parameter optimization. It may provide an effective method for mitigating casing deformation issues in the field. Importantly, the current research only shows that reducing pumping rate has a positive effect on mitigating casing deformation caused by medium-risk faults. The conclusion in this study may only be qualitative, and the treating parameter optimization for high-risk faults needs to be further explored.

Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
The authors declare no competing interest.