Factors affecting casing equivalent stress in multi-cluster fracturing of horizontal shale gas wells: finite element study on Weirong Block, southern Sichuan Basin, China

Multi-cluster fracturing technology was often used in horizontal well reservoir reconstruction to achieve production increase, which also affected casing equivalent stress distribution. This paper focuses on multi-cluster fracturing and establishes a fracturing model in line with the reality. The three-dimensional finite element model of multi-cluster fracture-formation-cement sheath-casing was proposed, the influence of cluster spacing and fracturing cluster number on casing equivalent stress was studied. On this basis, a single segment 8-cluster three-dimensional finite element model was developed. The influence of rock elastic modulus, casing inner wall pressure, geostress change and elastic modulus of cement sheath on casing equivalent stress was simulated from two aspects of uniform and non-uniform extrusion of wellbore. Actual data was used and analyzed for the fracturing section of a well in Weirong Block, southern Sichuan Basin, China. The results showed that the casing equivalent stress decreased with the increase of fracture dip angle. The casing equivalent stress increased with the increase of cluster spacing; however, it decreased with the increase of rock elastic modulus. The casing equivalent stress increased with the increase of casing wall pressure. Also, the cracks extrude the casing evenly did not affect the change on casing equivalent stress. It was also found that, when casing was uniformly squeezed by multiple fractures, the difference of ground stress had little effect on casing equivalent stress, while non-uniform extrusion had greater effect on casing equivalent stress. Further, when there was no wellhead pumping pressure, the casing equivalent stress increased with the increase of the elastic modulus of the cement sheath, and decreased on the contrary. The elastic modulus of rock was lower than that of cement sheath, and the casing equivalent stress increased with the increase of the elastic modulus of cement sheath, and decreased on the contrary. The research results had certain guiding significance for the prevention and control of casing damage in fracturing section.


E s
The elastic modulus of casing, GPa The half length of hydraulic fracture, m p n The net fracture pressure, MPa p o The vector sum of induced stress field and geostress field, and the vertical stress is σ v + σ H . The minimum horizontal geostress is σ n + σ min , the maximum horizontal geostress is σ m + σ max , MPa p i The pressure of casing inner wall, MPa r 1 The casing inner diameter, mm r 2 Casing outer diameter, mm r 3 Cement sheath outer diameter, mm r 4 Surrounding rock outer diameter respectively, mm μ c Poisson's ratio of cement sheath μ f Poisson's ratio of cement rock Zhao Yang and Rui Sun have contributed equally to this work.
* Rui Sun 954003199@qq.com 1 1 3 μ s Poisson's ratio of casing ν Poisson's ratio of rock θ The angle between fracture and wellbore θ 1 The angle between this position and the crack tip θ 2 The included angle between this position and the crack surface at the crack initiation posi-tion θ 3 The angle between this position and the crack tip σ H The vertical in-situ stress, MPa σ m The maximum horizontal induced stress, MPa σ max The horizontal in-situ stress, MPa σ min The minimum horizontal in-situ stress, MPa σ n The minimum horizontal induced stress, MPa σ r The radial stress of casing, MPa σ θ The tangential stress stress of casing, MPa σ t The axial stress expression of casing, MPa σ v The vertical induced stress, MPa

Introdction
Shale gas is stored in deep buried shale with tight lithology. Horizontal well fracturing technology is the main method to develop shale gas at present Zhou et al. 2016;Li et al. 2017). When horizontal well staged fracturing is implemented, superimposed induced stress will be generated in multiple segments and clusters of fractures, which will change the magnitude and direction of the in-situ stress around the casing, and the stress field of the casing will also change accordingly (Wang et al. 2012;Yin et al. 2012). In the actual fracturing process, casing deformation is particularly prominent, affecting the normal production of shale gas wells . During oil and gas development stages, the casing is subjected to extremely harsh environmental conditions for a long time, which could lead to deformation. It has been observed that casing deformation falls into four categories, namely, buckling, shear, extrusion, and pitting (Yan et al. 2019). A considerable amount of literature has been published on casing deformation during hydraulic fracturing. These studies indicate that the potential causes of casing damage are asymmetry of fracture propagation, stress concentration, temperature changes, and poor cementing quality (Fan et al. 2017;Yan et al. 2017;Shen and Zhang 2019;Mu et al. 2021;Liu et al. 2022;Yang et al. 2023). Yu et al. (2014) pointed out that asymmetric hydraulic fracturing will lead to unbalanced formation rock deformation, resulting in casing asymmetry, extrusion load or stress deficiency around the casing (Yu et al. 2014). Liu et al. (2017) used the finite element method to simulate the casing equivalent stress and deformation under the condition of asymmetric fracturing . The research results showed that the asymmetric reconstruction around the wellbore will lead to the overall lateral bending load of the casing, which will cause serious shear deformation of the casing in the presence of faults. Based on summary and analysis of casing deformation in the Sichuan Basin over the years, Zhang et al. (2021) proposed quantitative risk assessment approach (QRA) and applied finite element method to research the effect of multi-scale nature fractures (Zhang et al. 2021). A multifracture dynamic propagation model was established by Meng et al. (2022) to investigate the stress regime transformation in the case of the mutual interference of multiple fractures (Meng et al. 2022). Chang et al. (2022) set the interlayer interfaces and rock mechanical parameters according to the vertical stratification data. The complex fracture propagation process was simulated under different cluster numbers and cluster spacings. The multi-fracture competitive propagation during multi-clustered fracturing and the three-dimensional occurrence of the fracture system after fracturing were analyzed (Chang et al. 2022). Yin et al. (2022) combined the fracture-induced stress calculation model and the casing-cement-formation system stress distribution model to analyze the mechanical properties of the system during hydrofracturing. We used the model to calculate the casing stress distribution under fracture-induced stress and analyzed various influencing factors that affect the maximum equivalent stress (MES) on the casing (Yin et al. 2022).
In recent years, different scholars have studied many changes in casing equivalent stress by developing induced stress calculation models and discrete element analysis methods. Most of the research objects are single cluster fractures, which can not truly reflect multi-cluster fracturing conditions, and two-dimensional discrete element analysis methods can not reflect the morphological characteristics of real fractures. Based on the finite element analysis method, a three-dimensional multi-cluster finite element model close to the fracture characteristics of the fracturing section can be developed to study the influence of different factors on casing equivalent stress. However, due to the difficulty of modeling and calculation, there are few relevant research results. In this paper, the finite element analysis method is used to simulate the casing equivalent stress distribution by developing a three-dimensional finite element model of multi-cluster fracture-formationcement sheath-casing. The influence of cluster spacing, number of fracturing clusters, rock elastic modulus, casing inner wall pressure, changes in geostress field and elastic modulus of cement sheath on the casing equivalent stress distribution is studied and analyzed through the actual fracturing section (A certain section in staged fracturing of actual field horizontal wells). The writing idea of this article is shown in Fig. 1.

Crack induced stress
The induced stress field will have an impact on the adjacent reservoir, and the increase of the number of segments will have a cumulative effect on the stress deformation of the casing at the target point.

3
Based on the homogeneous and isotropic two-dimensional plane strain model, the geometric model of vertical crack induced stress field is developed in Fig. 2.
Assume that the ideal crack shape is a long narrow crack with small width and large length, and there is a certain net pressure in the two-dimensional plane. As the size of the wellbore is significantly lesser than the fracture length, the fracture-induced stress in the vicinity of the wellbore can be ignored. In this study, the induced stress at the center of the wellbore is taken as the far-field induced stress. According to Sneddon formula, the induced stress in the plane generated by the net pressure is calculated as follows (Sneddon and Elliot 1946): where σ n is the minimum horizontal induced stress, MPa; σ m is the maximum horizontal induced stress, MPa; σ v is the vertical induced stress, MPa; p n is the net fracture pressure, MPa; a, b and c are the distance between the shaft cross section and the crack at the analysis point, m; ν is Poisson's ratio of rock of rock; H is the half length of hydraulic fracture, m; σ min is the minimum horizontal in-situ stress, MPa; θ 1 is the angle between this position and the crack tip, °; θ 2 is the included angle between this position and the crack surface at the crack initiation position, °; θ 3 is the angle between this position and the crack tip, °.

Casing equivalent stress
After the cement slurry is solidified, the casing, cement sheath and formation will be consolidated into a combined elastomer as shown in Fig. 3. The horizontal well casing is subject to the combined effect of induced stress, in-situ stress and casing inner wall pressure generated by hydraulic fractures: In Fig. 2, p o is the vector sum of induced stress field and geostress field, and the vertical stress is σ v + σ H , the minimum horizontal geostress is σ n + σ min , the maximum horizontal geostress is σ m + σ max , MPa;p i is the pressure of casing inner wall, MPa; r 1 is the casing inner diameter, mm; r 2 is the Casing outer diameter, mm; r 3 is the Cement sheath outer diameter, mm; r 4 is the surrounding rock outer diameter respectively, mm.
According to Lame formula, the stress distribution of casing is (Li and Xi 2020): (1) n = P n b H H 2 ac 3 2 sin 2 sin 3 2 1 + 3 where f 1 ~ f 8 are coefficients; σ r is radial stress, σ θ is tangential stress; σ t is axial stress expression of casing, MPa; μ s is the Poisson's ratio of casing; μ c is the Poisson's ratio of cement sheath; μ f is the Poisson's ratio of rock; E s is the elastic modulus of casing, GPa; E c is the elastic modulus of cement sheath; E f is the elastic modulus of rock; t 12 = r 2 /r 1 , t 23 = r 3 /r 2 , t 34 = r 4 /r 3 .

Model assumptions
After cementing, the cement slurry solidifies to form a solid cement sheath, which tightly couples the casing and the formation. In order to make the casing equivalent stress analysis more reliable and reasonable, it is necessary to simplify the mechanical model and make the following basic assumptions (Lin et al. 2022): (1) The casing cement sheath formation is in close contact with no gap (2) The casing, cement sheath and wellbore are all ideal cylindrical (3) Casing, cement sheath and formation are ideal elastic materials (4) The influence of temperature change on material properties is not considered (5) The contact in the model does not consider the relative displacement

Model parameter
Based on the logging data of a well in Changning Weiyuan National Shale Gas Demonstration Area in Sichuan Basin, the geostress parameters can be obtained, where σ H = 35 MPa, σ Min = 29 MPa, σ max = 48 MPa (Jiang et al. 2015).
The casing size commonly used in horizontal well section of fracturing construction in shale gas demonstration area is selected to establish the model. According to Saint Venant's principle (Mises 1945), more than 5 times of the well diameter size is taken as the formation size to eliminate the influence of boundary effect on casing equivalent stress in the model. The size parameters are as follows: the stratum model (length, width and height) is 2 m × 2 m × 2 m, outer diameter of cement sheath: 215.9 mm; Casing outer diameter: 139.7 mm; Inner diameter of casing: 121.36 mm.
The elastic modulus and Poisson's ratio of the formation material are determined by referring to the rock mechanical properties of the shale layer in a well 3000 m deep in Weiyuan National Shale Gas Demonstration Area. The rock properties in this area are characterized by high elastic modulus and low Poisson's ratio. The casing material shall be P110 steel grade with a yield strength of 758 MPa. The casing material parameters shall be set with reference to the mechanical properties of this grade of steel. The setting of mechanical properties of cement sheath refers to the cement properties commonly used in shale gas production and cementing. Mechanical parameters related to specific materials are shown in Table 1. Perforation parameter setting: perforation density is 2 holes/m, and perforation phase angle is 180°.

Mesh generation and load application
The mesh type: all finite element mesh models adopt hexahedral structured mesh for calculation and analysis, since casing is the main research object, the mesh size of cracks is coarsened to improve the calculation efficiency of the model. The approximate global mesh size of formation and fracture is 0.5, cement sheath is 0.1, and casing is 0.05. The mesh number: the total number of cells is 68563.The mesh generation results are shown in Fig. 4a.
Fully constrain the bottom surface of the model, and fix the displacement of the model surface in the X and Y directions of the model. Considering the fracturing fluid with high pressure inside the fracture, according to the actual construction pressure parameter of a shale gas well is 55 MPa, the bottom hole pressure is set, 70 MPa is applied to the inner wall of the casing and the inside of the fracture. The three-dimensional finite element model of fracture-formation-cement sheath-casing is developed as shown in Fig. 4b, where a is the formation, b is the hydraulic fracture, and c is the cement sheath-casing combination.

Model validation
Due to the heterogeneity of rocks and the existence of faults in deep shale reservoirs, the perforation phase angle needs to be corrected according to the distribution of oil and gas in the reservoir. There will be two situations: uniform fracture squeezing the wellbore or uneven fracture squeezing the wellbore. In order to study the distribution of casing equivalent stress during uniform and non-uniform extrusion of the wellbore, the fracture inclination Angle. The Angle θ (θ is the Angle between fracture and wellbore, °) between the fracture and the wellbore, was set as 15°, 30°, 45°, 60°, 75°, and 90° respectively (hydraulic fracture: uniform extrusion of the wellbore). At the same time, considering that the fracture width increases with the increase of perforation aperture in the actual fracturing construction process, set the fracture width as 40 mm and 90 mm to study the influence of the fracture width on the casing equivalent stress distribution. The simulation results are verified based on the mechanical model of fracture-formation-cement sheath-casing. Some It can be seen from Fig. 5, when the cracks compress the casing unevenly, the stress on the outer wall of the casing at the inclined side of the cracks increases. As the uneven compression degree of the cracks weakens, the stress on both sides of the outer wall of the casing at the perforation site begins to be evenly distributed. As the cracks change from uniform extrusion to non-uniform extrusion, the casing equivalent stress starts to concentrate from the stress concentration on the inner wall of the whole perforation part to the inner wall at one end of the perforation part (the inclined direction of the cracks).
It can be seen from Fig. 6 that the equivalent stress of casing decreases with the increase of fracture dip angle. This is because the degree of non-uniform fracture squeezing the wellbore is also decreasing. The casing equivalent stress starts to be evenly distributed, which is conducive to the protection of the casing. When the fracture dip angle is 15°, the casing equivalent stress exceeds the yield strength, and the casing undergoes plastic deformation. Therefore, 15° is the risk dip angle. When the fracture angle is the same, the increase of the crack width increases the stress concentration of the crack, so the influence of the crack width on the casing equivalent stress cannot be ignored.
The trend of theoretical calculation and model simulation results is the same, but there is some error between them. The numerical error shows a decreasing trend with the increase of fracture dip, but the overall value is close. The theoretical formula verifies that the fracture-formationcement sheath-casing model can be used to simulate the influence of induced stress field generated by fractures on casing equivalent stress.

Multi-cluster fracturing
Due to the large horizontal geostress difference in deep shale gas reservoirs, it is not conducive to the formation of complex mesh fractures and requires fracturing fluid injection with large displacement and large sand volume. Therefore, the conventional perforation method is no longer applicable to the fracturing and fracturing of deep shale gas reservoirs (Lin et al. 2022). At present, the method of large pore diameter and short cluster spacing is often used to improve the fracturing effect of deep shale gas wells .
In order to simulate the influence of single section and multiple clusters on casing equivalent stress distribution, with the help of the developed three-dimensional finite element combination model of fracture-formation-cement sheath-casing. The model size is extended from 2 to 10 m, that is, a single section "multiple cluster" model with a length of 10 m, a width of 2 m and a thickness of 2 m is developed, as shown in Fig. 7.

Cluster spacing
The fractures in the perforation section will produce induced stress fields, while the occurrence of multiple clusters of fractures will cause multiple fracture induced stress fields in the perforation section. Therefore, it is necessary to study the interaction between multiple induced stress fields, and the distance between clusters will directly affect the action of induced stress fields. Therefore, with four clusters in a single section as the background, the distance between clusters will vary from 900 to 2400 mm, taking into account the actual fracturing construction project. It is common that fractures extrude the wellbore unevenly, so the phase angles are set as 15° and 45°. When the fracture angle is 45°, the three-dimensional equivalent stress nephogram and twodimensional section equivalent stress nephogram under different cluster spacing are shown in Fig. 8, and the curve of casing equivalent stress variation with cluster spacing under different perforation phase angle is shown in Fig. 9.
It can be seen from Fig. 8 that when the cluster spacing is small (the cluster spacing is 900 mm), a stress dispersion area (minimum stress area) is formed at the center between clusters (along the seam length). When there is no crack in the inclined direction of the crack, the stress presents a convex distribution. The casing equivalent stress increases with the closer to the perforation section, and the casing equivalent stress concentration area is mainly located on the inner wall of the casing in the perforation section. When the cluster spacing is large (the cluster spacing is 2400 mm), the stress dispersion area between clusters turns to the center area along the seam height plane; however, that does not change the trend that the closer the cluster is to the perforation location, the greater the stress is. The stress concentration location of the casing expands from the inner wall of the perforation section casing to the outer wall of the seam height surface casing, and the stress concentration range expands with the increase of the cluster spacing.
It can be seen from Fig. 9 that the casing equivalent stress tends to increase with the increase of cluster spacing. However, when the cluster spacing exceeds 2100 mm, the casing equivalent stress tends to be stable, which means that the smaller the cluster spacing is, the more the induced stress field between clusters will offset; that is, the interference effect between cracks is obvious. With the increase of cluster spacing, the induced stress field changes from stress interference to partial superposition, which increases the casing equivalent stress; however, the cluster spacing increases to a certain distance. The influence between cluster induced stress fields tends to be minimized, and the change of cluster spacing has less influence on casing equivalent stress. When the crack dip angle is 15° and the cluster spacing is 1200 mm, the casing equivalent stress breaks through the yield strength, the casing undergoes plastic deformation, and the casing may be damaged. This is because when the cluster spacing is the same, the smaller the crack dip angle is, which enlarges the non-uniform compression degree of the crack and increases the stress concentration degree of the casing. Therefore, when considering the impact of cluster spacing on the casing damage, the non-uniform distribution of the cracks should also be considered in the actual construction. The occurrence of non-uniform extrusion of cracks will greatly increase the probability of casing damage (Chang et al. 2022).

Number of fracturing clusters
Under the condition of ensuring that each cluster of fractures in the section can be pressed open, increasing the number of perforating clusters can effectively improve the fracture control reserves, increase the fracturing volume, and increase the oil and gas production. At present, the fracturing construction mostly uses single segment and multiple clusters to increase production. Therefore, we need to study the impact of cluster number on the casing equivalent stress distribution. The single fracturing section generally contains 2 ~ 8 perforation clusters. At the same time, because the cluster number and cluster spacing are inseparable, the cluster spacing is still set at 900 mm ~ 2400 mm. The cluster number is simulated for two cases of single section with two clusters and single section with four clusters. The fracture is horizontal. When the cluster spacing is 2400 mm, the corresponding casing equivalent stress nephogram of different cluster numbers is shown in Fig. 10, and the curve of casing equivalent stress variation with cluster number under different cluster spacing is shown in Fig. 11.
It can be seen from Fig. 10 that when the fracture is horizontal, the stress concentration position of the casing is located on the inner and outer wall of the casing at the perforation position. The increase of the number of clusters only increases the casing equivalent stress value but does not change the casing equivalent stress distribution.
It can be seen from Fig. 11 that the increase of cluster number will not change the change of casing equivalent stress with cluster spacing. When cluster spacing is less than 1500 mm, the increase of cluster number will have little impact on casing equivalent stress. While when cluster spacing exceeds 1500 mm, the casing equivalent stress increases with the increase of cluster number, which indicates that the interference effect between clusters has reached a stable level when cluster spacing is small. The increase in the number of clusters does not affect the stress field near the wellbore, while the cluster spacing increases, the interference stress between clusters decreases. The increase in the number of clusters to a certain extent produces the additional stress effect, making the casing equivalent stress begin to increase with the increase in the number of clusters (Chang et al. 2022).

Factors affecting casing equivalent stress
With the help of the developed multi-cluster fracture-formation-cement sheath-casing 3D finite element combination model, the model size is changed from 10 m × 2 m × 2 m becomes 20 m × 2 m × 2 m. At present, the single fracturing section generally includes 3 ~ 8 perforation clusters. Therefore, 8 clusters of single section are taken for research. Based on the study of the angle between the fracture and the wellbore and the fracture width, the risk dip angle of hydraulic fractures is selected as 15°. Since the width of fractures formed by common large aperture perforations in the current oilfield is not more than 100 mm, so the model fracture width is set as 100 mm, so as to establish 8 cluster fracturing models. Based on the two cases of fracture dip angle of 90° (hydraulic fracture) and 15° (multi-cluster fracture uniform extrusion and non-uniform extrusion casing), the effects of rock elastic modulus, casing inner wall pressure, geostress field and cement sheath parameters on casing equivalent stress are simulated.

Elastic modulus of rock
The elastic modulus of rock directly affects the anti deformation ability of rock. The easier the rock is deformed, the stronger the compression strength of the surrounding rock near the wellbore against the casing, and the weaker the rock's ability to resist the far-field stress (Nie et al. 2019). Combined with the near wellbore induced stress field generated by fracturing fractures, the probability of casing damage is greatly improved. In order to study the influence of different rock elastic modulus on casing equivalent stress, the change range of rock elastic modulus is set as 5GPa ~ 55GPa. The nephogram of casing equivalent stress under different rock elastic modulus is shown in Fig. 12, and the curve of casing equivalent stress changing with rock elastic modulus is shown in Fig. 13.
It can be seen from Fig. 12 that casing equivalent stress is mainly concentrated on the inner wall of the perforated section, and stress concentration occurs on the outer wall of the far cluster perforated section. The stress dispersion  area between clusters is located near the perforation section. With the increase of the elastic modulus of rock, the range of stress concentration on the inner wall of the perforated section is reduced; the stress dispersion area extends to the entire area between clusters, and the stress concentration on the outer wall of the far shower section is relieved.
It can be seen from Fig. 13 that the casing equivalent stress decreases with the increase of the elastic modulus of the rock. This is because with the increase of the elastic modulus of the rock, the ability of the rock to resist the induced stress generated by cracks is enhanced, and it bears a part of the induced stress and far-field geostress that should have acted on the outer wall of the casing (Nie et al. 2019). At the same time, the lower the deformation degree of the formation under the same geostress condition, and the difficult deformation of the formation will effectively fix the geometry of the wellbore. The casing in the horizontal section of shale gas well is not easy to deform, so that it is not easy to produce large stress concentration on the casing. When the elastic modulus of the rock is the same, the extrusion of the fracture on the casing changes from uniform to non-uniform, which enhances the stress field induced by the near wellbore fracture. Therefore, when the elastic modulus of the rock is lower than 25GPa, the non-uniform extrusion of the fracture will make the casing equivalent stress reach the yield stress limit, and the casing begins to undergo plastic deformation. When the elastic modulus of rock reaches 5GPa, casing damage may occur when the casing is squeezed by cracks uniformly. Compared with non-uniform extrusion, uniform extrusion reduces the probability of plastic deformation of casing and plays a certain role in protecting casing. However, staged fracturing technology is easy to re fracturing shale in the same area, and multiple fracturing will sharply reduce the elastic modulus of rock. When the elastic modulus of rock continues to decrease, the influence on casing equivalent stress becomes greater and greater, especially when the elastic modulus of rock drops below 15 GPa, the trend of casing deformation increases sharply. This shows that shale gas well has experienced multiple fracturing, and the elastic modulus of formation rock decreases, which has a great impact on casing equivalent stress, that is, casing damage is easy to be caused by repeated fracturing (Zheng et al. 2021).

Casing inner wall pressure
Fracturing operation usually uses large construction pressure at the wellhead to break the rock with fracturing fluid to achieve the purpose of fracture formation, but it also generates huge internal wall pressure on the casing. In order to study the influence of casing inner wall pressure change on casing equivalent stress distribution under the condition of whether the fracture is evenly squeezed. According to the change of bottom hole pressure during fracturing of a well in Weiyuan Shale Gas Demonstration Area, the change range of casing inner wall pressure is set as 50 MPa ~ 100 MPa. The nephogram of casing equivalent stress under different casing wall pressures is shown in Fig. 14, and the curve of casing equivalent stress changing with casing wall pressure is shown in Fig. 15.
It can be seen from Fig. 14 that between clusters, the stress on the extrusion side of the cracks on the outer wall of the casing in the perforation section is large. The casing equivalent stress is mainly concentrated on the inner wall of the casing at the perforation location, and the inner wall of the casing between clusters forms a stress dispersion area. However, with the increase of the pressure on the inner wall of the casing, the stress concentration on the outer wall of the casing between clusters decreases. The stress concentration on the perforation part of the inner wall of the casing changes from uniform distribution to the stress concentration on the left and right perforations, and the stress concentration on the inner wall of the casing between clusters increases.
It can be seen from Fig. 15 that the casing equivalent stress increases with the increase of casing inner wall pressure, while whether the cracks squeeze the casing evenly does not affect the change of casing equivalent stress. When the pressure on the inner wall of the casing is the same, the casing equivalent stress increases with the increase of the non-uniform squeezing degree of the fracture. This is because the tilt of the fracture increases the non-uniformity of the induced stress field near the wellbore, which deepens the external squeezing force of the casing and its nonuniformity. The fluid entering the fracture is more likely to leak into the formation near the wellbore, which increases the formation load near the wellbore to a certain extent, and increases the additional load outside the casing. When the pressure on the inner wall of the casing increases to 80 MPa and the fracture compresses the casing unevenly, the casing equivalent stress exceeds the yield strength and the casing is damaged. When the fracture compresses the casing evenly, the casing equivalent stress concentration will be relieved and the casing is not easy to be damaged. This is because when the fracture compresses the casing evenly, the pressure of the liquid column on the inner wall of the casing in the perforation section will offset the load from the far field stress field and the near wellbore induced stress field to a certain extent ).

Change of geostress field
The geostress fields of different parts of horizontal wells are different. At the same time, the induced stress generated by fracturing fractures changes the original geostress field of 1 3 the formation near the wellbore. Therefore, it is necessary to study the stress change rule of casing under the effect of the induced stress field generated by multiple clusters of fractures.
Set the current vertical stress σ H varies from 30 to 50 MPa, and the minimum horizontal geostress σ min value is fixed as 30 MPa, and the maximum horizontal geostress σ max varies from 30 to 90 MPa, and the ground stress difference σ(σ = σ H − σ max ) is 0 MPa ~ 60 MPa. When σ H = 30Mpa, σ = 0 MPa or 20 MPa, the casing equivalent stress simulation results under different geostress differences are shown in Fig. 16. The three-dimensional change curve is as shown in Fig. 17, where the surface A is the plastic deformation surface of the casing, namely the casing damage surface.
It can be seen from Fig. 16 that when the difference of ground stress is small, the stress concentration range of the outer wall of the casing between clusters is large, while the stress concentration phenomenon occurs in the upper and lower parts of the inner wall of the casing between clusters, and the stress on the left and right sides is scattered. The stress concentration of the casing is mainly located in the inner wall of the casing at the perforation location. When the difference of ground stress increases to 20 MPa, the stress concentration range on the outer wall of the casing between clusters decreases; the stress concentration degree on the upper and lower inner wall of the casing decreases; and the stress concentration range at the perforation site decreases, which enhances the stress concentration degree of the casing and enhances the difference of stress distribution near the wellbore.
It can be seen from Fig. 17 that Under the influence of inner wall pressure of 70 MPa, in-situ stress field and multi-cluster fracture induced stress field, the casing equivalent stress increases with the increase of the geostress difference. That indicates that the increase of the geostress difference is conducive to enhancing the effect of multi-cluster crack induced stress field attached to the external load of the casing and the fluid pressure of the casing inner wall on the casing wall, thus enhancing the casing equivalent stress concentration. When casing is uniformly squeezed by multiple fractures, the change of ground stress difference has less influence on casing equivalent stress. While nonuniform extrusion has greater influence on casing equivalent stress, which shows that multiple induced stresses generated by uniform casing extrusion by multiple fractures can offset the damage of casing caused by far field stress to a certain extent; and the induced stress field generated by non-uniform extrusion instead strengthens the effect of near wellbore stress field, and the distribution of stress field is different. Finally, the casing equivalent stress is distributed unevenly, increasing the probability of casing damage.
When the geostress difference is constant and multiple fractures extrude the casing uniformly, the casing equivalent stress decreases with the increase of vertical stress. That indicates that the induced stress field generated by the uniform fracture extruding the casing can offset the vertical stress from the in-situ stress field, making the stress field around the wellbore weaken with the increase of vertical stress. Therefore, when the casing is uniformly squeezed, the influence of the change of vertical stress on the casing equivalent stress cannot be ignored.
However, when the casing is squeezed unevenly by multiple cracks, the casing equivalent stress basically remains stable with the increase of the vertical stress. When the vertical stress increases to 50 MPa, the casing equivalent stress changes significantly, which indicates that the vertical stress has little influence on the non-uniform crack induced stress field, and the crack induced stress field has the greatest influence on the casing equivalent stress at this time. When the difference of ground stress is 17 MPa and the vertical stress is 30 MPa, the casing equivalent stress reaches 758 MPa and the stress reaches the yield limit. The point where the difference of ground stress is 17 MPa, the vertical stress is 30 MPa, and the casing equivalent stress is 758 MPa is located on surface A of the casing damage surface in the three-dimensional coordinate system. If other points cross surface A, plastic deformation occurs in the casing fracturing section. Therefore, the well section corresponding to the appropriate geostress field should be selected according to the three-dimensional chart for multi-cluster perforation fracturing, so as to reduce the probability of casing damage (Lian et al. 2022).

Elastic modulus of cement sheath
As an important barrier for casing protection, the quality of cement sheath is crucial to the strength and safety of casing. The elastic modulus of cement sheath determines the softness and hardness of cement sheath. Based on the different effects of the softness and hardness of cement sheath under different working conditions, and considering the differences of construction parameters and geological conditions (Zhao et al. 2019). The influence of the change of elastic modulus of cement sheath on the equivalent stress of casing during the change of construction pressure and formation hardness is studied. Because it is ideal to keep the casing pressed uniformly by multi-cluster fractures, a multi-cluster fracture non-uniform squeeze wellbore model is selected for numerical simulation.

Construction pressure
The construction pressure only equals to wellhead pumping pressure. The pressure in the casing will decrease due to the reduction of wellhead pumping pressure during fracturing construction. For this reason, this paper has considered the more extreme working conditions in the calculation, assuming that there is no wellhead pumping pressure. The high or low pressure of the wellhead pumping pressure in the fracturing process has different effects on the casing, therefore, two conditions are set: 40 MPa (low pressure) and 90 MPa (high pressure). When there is no wellhead pumping pressure, the elastic modulus of cement sheath is 5GPa and 45GPa respectively, and the corresponding simulation results are shown in Fig. 18. The change curve of casing equivalent stress with elastic modulus of cement sheath under different construction pressures is shown in Fig. 19.
It can be seen from Fig. 18 that the increase of the elastic modulus of the cement sheath does not change the stress concentration position of the casing, that is, the casing equivalent stress is mainly concentrated in the upper and lower parts of the casing inner wall at the perforation position; however it will expand the stress concentration range of the inner wall. The stresses on the inner and outer walls of the casing between clusters are relatively dispersed. With the increase of the elastic modulus of the cement sheath, the stress dispersion range decreases.
It can be seen from Fig. 19 that when there is no wellhead pumping pressure, the casing equivalent stress will increase with the increase of the elastic modulus of the cement sheath. This is because with the increase of the elastic modulus of the cement sheath, although the cement sheath bears a certain far-field geostress field, its load transfer coefficient will increase due to the increase of the stiffness of the cement sheath, and the buffering effect of the cement sheath on the external load acting on the outer wall of the casing will decrease. The load on the casing caused by the induced stress field produced by the non-uniform extrusion of casing with multi-cluster fractures near the wellbore is increasing. When the construction pressure is 40 MPa and 90 MPa, the casing pressure decreases with the increase of the elastic modulus of the cement sheath; that indicates that there is pressure on the inner wall of the casing, which will make part of the external load offset with the pressure of the inner wall of the casing, thus reducing the stress on the casing. The elastic modulus of the cement sheath is higher than 25GPa, the casing equivalent stress decreases slowly, and the influence of the elastic modulus of the cement sheath on the change of the casing equivalent stress decreases. When the elastic modulus of the cement sheath is the same, the casing equivalent stress decreases first and then increases with the increase of the construction pressure. When the elastic modulus of the cement sheath is lower than 13 GPa and the construction pressure is 90 MPa, the casing equivalent stress exceeds the yield limit and the casing damage occurs. When the elastic modulus of the cement sheath is higher than 13 GPa and there is no wellhead pumping pressure on the inner wall of the casing, the casing begins to deform plastically. Therefore, after high construction pressure or no wellhead pumping pressure, the elastic modulus of the cement sheath 13 GPa is the risk value (Yin et al. 2018).

Coupling elastic modulus of rock with cement sheath
Assuming that multi-cluster fracturing is carried out in extremely soft formation, the elastic modulus of rock is set as 5GPa, and the elastic modulus of common tight reservoir  is mostly 20GPa ~ 40GPa. Therefore, the elastic modulus of rock is set as 35GPa. Considering the existence of extremely hard formation, the elastic modulus of rock is set as 65GPa. Under the above three rock conditions, the impact of elastic modulus of cement sheath on casing equivalent stress is studied. Under different construction pressures, the curve of casing equivalent stress changing with the elastic modulus of cement sheath is shown in Fig. 20.
It can be seen from Fig. 20 that when the elastic modulus of rock is 5GPa, the casing equivalent stress increases with the increase of the elastic modulus of cement sheath, which indicates that the formation is too soft to resist the far-field stress. The stress field coupled with multiple clusters of non-uniform fractures causes excessive load to act on the casing, so the softer cement sheath is required to play a role in buffering the load. When the elastic modulus of rock is 35GPa or 65GPa, the casing equivalent stress tends to decrease first and then increase. When the elastic modulus of rock is equal to the elastic modulus of cement sheath, it is the inflection point of the trend. This is because when the elastic modulus of cement sheath is lower than the elastic modulus of rock, the rock with high stiffness will bear most of the far-field stress. Increasing the elastic modulus of cement sheath can play a supporting role for the casing. When the elastic modulus of cement sheath is higher than the elastic modulus of rock. Then the stiffness of cement sheath is too large to play a buffer role, and the influence of multi-cluster fracture induced stress field on near wellbore stress is amplified; however the increase of casing equivalent stress is too small, which indicates that under the condition of non-uniform oiling casing with multi-cluster fractures. When the elastic modulus of rock has a certain stiffness, the impact of too large elastic modulus of cement sheath on casing equivalent stress is not obvious.
To sum up, the increase of the elastic modulus of the cement sheath after the fracturing pump shutdown will make the casing equivalent stress increase and finally break the yield limit. When the construction pressure is high pressure or low pressure, the elastic modulus of the casing equivalent stress decreases with the increase of the cement sheath. When the elastic modulus of cement sheath is higher than that of rock, the casing equivalent stress increases with the increase of the elastic modulus of cement sheath, and decreases on the contrary (Guo et al. 2018).

Case demonstration
Taking a single well section and a cluster in Weirong Block, southern Sichuan Basin, China as an example, the casing equivalent stress distribution is analyzed when there is an induced stress field in the fracturing section. Perforation parameter settings for a cluster of horizontal well: perforation density is 6 holes/m, perforation angles are 60°, 90° and 120° respectively, and the length of each cluster of fracturing sections is 1 m. According to the field data, the maximum wellhead construction pressure can reach 60 MPa. Therefore, the bottom hole pressure, that is, the internal pressure of the casing, is set to 80 MPa. A three-dimensional finite element model of multi fracture-formation-cement sheathcasing is developed, as shown in Fig. 21. The stress distribution of casing is shown in Fig. 22.
It can be seen from the calculation results in Fig. 22 that the in-situ stress field near the wellbore of the fracturing section changes due to the existence of the fracture induced stress field. Under the action of the near wellbore stress field and the internal pressure of the construction, the casing equivalent stress is mainly concentrated on the inner wall of the perforated section. Due to the existence of multi-fracture induced stress field, the equivalent stress of the outer wall of the casing attached to the near perforation position is offset to a certain extent, making the equivalent stress value of the casing near the perforation position smaller. Because there are 6 perforations that cause cracks to be too dense, the cumulative effect of the fracture induced stress field on the casing in the fracturing section is large, making the casing equivalent stress close to the yield strength. Although the casing has not undergone plastic deformation, there are certain risks.
According to the calculation results, high strength casing shall be used as far as possible at the fracturing section of the casing, and good cementing quality shall be used together to reduce the risk of casing damage.

Conclusions
The three-dimensional finite element model of fracture-formation-cement sheath-casing is developed, and the feasibility of the model is verified by theoretical formula, so as to otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.