ABAQUS Numerical Simulation Study on the Shear Instability of a Wellbore Induced by a Slip of the Natural Gas Hydrate Layer

Abstract

To study the shear deformation and failure characteristics of a wellbore and the interaction mechanism with its surrounding rocks induced by a layer slip during natural gas hydrates (NGHs) extraction, this paper conducted a numerical simulation study of wellbore shear induced by a layer slip using ABAQUS software and carried out a laboratory experiment of wellbore shear to verify the accuracy of the numerical model. The results show that the shear force–displacement curves obtained from the laboratory experiments and numerical simulations are consistent with five stages, including the compaction stage, linear stage, plastic stage, strain-softening stage and residual stage. The wellbore shows a “Z”-shaped deformation characteristic after its shear breakage. The shear force of the wellbore is maximum at the shear surface, and it is distributed in an approximate “M” shape along the shear surface. The axial force of the wellbore is small and uniformly distributed in the initial stage of the shear. The wellbore bending moment is minimum at the shear surface, with a value of 0, and it is distributed in a skew–symmetric wave shape along the shear surface. During the shearing, the evolution of the wellbore axial force and shear force can be classified into the distribution pattern along the radial direction on the shear surface and the pattern along the axial direction. The combination of the wellbore axial force and shear force causes the tensile–shear compound failure of the wellbore. During shearing, the wellbore and rock body gradually enter the plastic state with the increase in the shear displacement. When the entire cross-section of the wellbore is in the plastic state, a “necking” phenomenon of the wellbore begins to appear. During the shearing, the frictional dissipation energy and plastic dissipation energy increase constantly. In addition, the elastic strain energy increases to a peak and then decreases to a certain value, which remains unchanged along with the work conducted by the shear force.

1. Introduction

With the huge consumption of oil and gas resources, human beings will face energy depletion in the late twenty-first century, and the abundant oil and gas resources in the ocean will naturally become a strategic target for energy development and utilization by many countries. Natural gas hydrates (NGHs) are a new type of energy discovered in the ocean and tundra over the last two decades, which have the advantages of wide distribution, large resource reserves, high energy density and being clean and pollution-free [1]. Figure 1 shows the location of the NGHs test extraction area in Shenhu Sea, South China Sea. Therefore, NGHs will be an important target for future energy development. Since NGHs are generally endowed in shallow buried reservoirs in the deep sea at low temperature and high pressure and are sensitive to temperature and pressure changes, the existing extraction methods mainly include decomposition [2,3], thermal recovery [4], replacement [5] and chemical reagents [6]. Decomposition is a process that changes the stability of natural gas hydrate and causes it to decompose to produce natural gas and water by lowering the pressure of the hydrate reservoir and making it lower than the phase equilibrium pressure under other conditions in the region. Thermal recovery is a process that produces methane gas by increasing the temperature of the hydrate layer to reach the phase equilibrium temperature of the hydrate layer. The replacement method is a process of extraction by injecting CO₂ or other fluids and displacing and collecting the methane gas produced by the decomposition of natural gas hydrates. The chemical reagents method controls the rate of hydrate decomposition by injecting the corresponding synthetic inhibitor (methanol, etc.) into the hydrate layer, which indirectly improves the phase equilibrium conditions for the hydrate layer and, thus, causes hydrate decomposition to produce methane gas. In addition, gas hydrates may cause severe geological hazards due to the adverse effects of global warming on the geomechanical stability of gas hydrate deposits in marine and permafrost environments [7,8,9,10,11]. This is because the geomechanical, geophysical and hydrological properties of gas-hydrate-bearing sediments are essentially controlled by the presence of hydrates. Gas hydrates are sensitive to changes in temperature, and this sensitivity greatly affects the response of gas-hydrate-bearing sediments in different hydrate applications. The release of methane from gas-hydrate-bearing sediments is still considered to be a slow and chronic process. However, it may be exacerbated by anthropogenic global warming and rising local ambient temperatures. In particular, the thawing of permafrost and retreat of ice in the Arctic and Antarctic regions may lead to the increased emissions of trapped methane into the atmosphere, amplifying global warming and threatening global security. Before drilling for hydrate extraction, the subsurface rock layer is in stress equilibrium. Drilling will trigger changes in the temperature and pressure in the surrounding layer, which will cause the decomposition of NGHs and disrupt the original stress equilibrium of the layer, thus causing a redistribution of stress around the wellbore. When this stress equilibrium cannot be established, the layer will break down and shear failure will occur in the surrounding rock of the wellbore because the stresses applied exceed its own strength (Fereidounpour et al. [12]; Ahmadi et al. [13]; Wang et al. [14]; Ning et al. [15]). Therefore, the change in the mechanical properties of the subsea layer induced by NGH extraction is likely to lead to a large area of subsea landslide, resulting in shear instability and failure of the wellbore. Figure 2 shows the schematic diagram of a wellbore and the layer location in NGHs extraction.

Theoretical analysis is an effective approach to study the layer stability and wellbore shear failure in NGHs extraction, which focuses on the establishment of the corresponding theoretical model for shear strength calculation. In the study of layer and wellbore stability in NGHs extraction, Reem Freij et al. [16] established a pore elastoplastic mechanics model to study the thermodynamic change pattern of NGH layers during extraction based on the Mohr–Coulomb failure criterion, and they analyzed the pore pressure change of the NGH layers using the Darcy flow model and then predicted the potential failure modes of the layers and wellbores in the extraction. C. P. Tan et al. [17] developed an integrated risk management framework integrating a coupled model of wellbore stability in NGH layers, an optimization method for drilling fluid design, and risk assessment by studying the physical, mechanical and thermodynamic properties of NGH layers.

Laboratory experimentation is another important method to study wellbore shear deformation and failure in NGH layers. Yoneda et al. [18,19] conducted relevant experiments and concluded that the thickness of the shear band in the CH4–HBS appeared to be thinner than that in the hydrate-free sediment. Kato’s [20] research showed that the thickness of the shear band in the CH4–HBS appeared to have no dependence on the confining pressure. Kajiyama et al. [21] conducted a series of experiments and concluded that the increasing hydrate saturation leads to a steeper and narrower shear band for the CH4–HBS, and a greater fine component would delay the localized deformation development. Wu et al. [22] studied the localized deformation development of HBS and concluded that the shear band development showed a strong relationship with the stress–strain response and cementation failure. Wu et al. [23] mentioned that the localized deformation of the HFS was similar to that of dissociated HBS, presenting an obvious barrel-shaped deformation with X-shaped shear bands. Winters et al. [24] studied the differences in the acoustic and mechanical properties of NGH reservoirs with different fillings and different porosities by measuring the acoustic properties and mechanical properties of sediments filled with hydrate and ice in the pore space using an acoustic sonde and a triaxial experimental instrument, respectively. Georgia Tech conducted triaxial experiments on NGH-bearing rock samples made of three matrices—sand, clay and quartz powder—to determine the bulk modulus of elasticity, stress–strain curves, Poisson’s ratio and shear strength of these three NGH-bearing rock samples and to investigate the differences in the cohesion and friction angle between the quartz powder rock samples at 50% hydrate saturation and without hydrate based on the Mohr–Coulomb criterion. Hyodo et al. [25] investigated the mechanical properties of methane NGH core specimens at different saturations by indoor triaxial compression experiments, followed by experimental verification of the effects of temperature and effective confining pressure on the mechanical properties, and finally derived quantitative relationships for the variation of the mechanical properties of the NGH deposits with saturation, temperature and effective confining pressure. Clayton et al. [26] and Masui et al. [27] conducted triaxial compression experiments on four NGH core specimens drilled in the South Sea Trough area of Japan and synthetic NGH core specimens in the laboratory to compare the differences between in situ and synthetic specimens, and the results showed that the strength of the synthetic specimens was slightly stronger than that of in situ specimens, which may be related to the disturbance of in situ specimens during the acquisition or handling process.

The properties of the NGH layers obtained from the above experiments are of great significance for the selection of parameters in the numerical simulation modeling of wellbore shear in NGH layers. Due to the special stable environment for the presence and complex mechanical properties of NGHs, it is usually difficult to simulate and systematically study the actual working conditions when conducting relevant mechanical experiments. Numerical simulations can extend and enrich the results of laboratory experiments and obtain stress, strain and failure evolution patterns inside rocks and wellbores that are difficult to measure in laboratory experiments.

Ning et al. [28] analyzed the wellbore stability of drilling fluid intrusion into NGH layers during drilling by establishing a finite element model and concluded that shear instability is the main form of wellbore instability. Cheng [29], Shen [30] and Li et al. [31] studied the layer stability in the process of NGHs’ depressurization extraction based on the flow–solid coupled seepage model (FEPG, finite element platform), analyzed the yield stress variation law of an NGHs reservoir and the risk of layer yield instability caused by the borehole effect, and developed simulation and evaluation software for the “flow–solid coupled stability analysis of NGHs reservoir depressurization extraction” on this basis. Freij-Rutqvist et al. [32] proposed a numerical model for evaluating the wellbore stability of drilled NGH-bearing deposits and concluded that the heating effect of the drilling fluids on NGH-bearing layers can decompose NGHs, which greatly exacerbates the yielding of layers around the wellbore and causes wellbore failure. Freij-Ayoub R. et al. [17] found that the use of an overbalanced drilling method can reduce the extent of the yield zone and, thus, maintain the stability of the wellbore by using numerical simulations. Cheng et al. [33] proposed a nonisothermal, transient and two-phase fluid–solid coupled mathematical model to measure wellbore stability in NGH-bearing layers, and the results showed that decreasing the drilling fluid pressure and increasing the temperature can accelerate the decomposition of NGHs, which in turn can significantly increase the wellbore instability. Sun et al. [34] indicated that the free gas production from offshore gas NGH-bearing sediments could be controlled by setting the salinity of the drilling mud appropriately, which could avoid further destabilization of the wellbore. Wei et al. [35] revealed that the use of high-density, high-temperature drilling fluids can increase layer temperatures through thermal convection, exacerbating NGH decomposition and making it difficult to maintain wellbore stability. These studies have focused on the effects of drilling fluid density, temperature and salinity on wellbore stability in NGH-bearing layers when drilling in marine layers.

However, these studies did not further analyze the effect of the layer instability caused by the decomposition of NGHs on the stress distribution of the wellbore from a mechanical point of view. The effect of the formation of shear surface roughness was ignored, and no systematic analysis of the interaction and distribution of the plastic zone during the destruction of the wellbore and its surrounding rock was carried out. Therefore, this paper applied ABAQUS numerical simulation software to establish a wellbore shear model for the NGHs layer and its overlying layer, and it carried out numerical wellbore shear experiments to compare with laboratory experimental results for verification. The shear slip surface roughness of different layers is highly stochastic due to the variability in the composition of NGHs layers and their overlying layers, so shear surfaces with different roughness will affect the shear characteristics of the wellbore. Therefore, numerical simulations were carried out using the above model to analyze the shear characteristics of the wellbore and the interaction between the wellbore and the layer during the failure in depth considering factors such as the shear surface roughness of the layer. In addition, the failure mode and destabilization mechanism of the wellbore were further studied from the perspectives of the shear stress and energy dissipation of the shear system. The research results are expected to provide reference for layer stability and wellbore instability control during extraction of NGHs.

2. Numerical Modeling and Laboratory Experimental Validation

2.1. Laboratory Experiments

To verify the accuracy of the numerical model of wellbore shear, this paper carried out a laboratory experiment of wellbore shear considering the rough shear surface, obtained the shear–displacement curve and compared the experimental results with the simulation results to perform a preliminary analysis of the deformation and displacement of the wellbore.

2.1.1. Specimen Preparation

In order to study the effect of the surface topography on the shear characteristics of the wellbore, shear surfaces with different roughness were selected, and the surface topography of the shear surfaces were obtained by 3D laser scanning, as shown in Figure 3. The root mean square calculation method of the slope was used to calculate the roughness of the four shear surfaces based on the 3D laser scan data.

$$Z_2=\sqrt{\frac{1}{\left(n-1\right){(\Delta x)}^2}{{\displaystyle \sum }}_{i=1}^{n-1}{(Z_{i+1}-Z_i)}^2}$$

(1)

$$\mathrm{SRC}=32.69+32.98\log Z_2 \left(R=0.993\right)$$

(2)

where $Z_2$ is the root mean square of the slope, $\Delta x$ is the interval between data points, $n$ is the number of data points on each curve, $Z_i$ is the z-coordinate of each curve, R is the correlation coefficient, and SRC is the shear surface roughness coefficient.

Then, the average value of the roughness of 50 curves was taken as the roughness value of the shear surface. Finally, the roughness value of the shear surfaces was calculated as 6.46. The 3D modeling software Rhino was used to model the 3D printed model of the upper and lower plates of the shear specimen, as shown in Figure 3. The shear upper and lower disc specimens used in the test were mainly obtained by 3D printing technology. Resin sand and quartz sand were used as the base material for printing specimens, and furan resin was used as the binder to bond their powder particles together to make the shear upper and lower plates using the VX2000 3D printer.

The material of the metal tube used in the experiment specimens was 304 steel, whose material properties are similar to those of the casing material applied in practical NGHs extraction engineering [36]. In the process of specimen preparation, due to the limitation of the specimen size, it is difficult to use cement material to fill the space between the metal tube and the upper and lower plates of the shear specimen, i.e., the “cement sheath” part in the engineering. Therefore, high-strength gypsum, which has similar mechanical properties to cement, was used as the filling material, and the mass ratio of gypsum:water = 1:0.3 was used in the preparation of the specimens. In order to test the basic mechanical parameters of high-strength gypsum and cement, two sets of standard specimens were made from high-strength gypsum and cement. The specimens were fully maintained and then subjected to uniaxial compression experiments on a Shimadzu AG-X250 testing machine at a loading rate of 0.5 MPa/s. The uniaxial compressive strength and modulus of elasticity of the high-strength gypsum were tested to be 3.11 MPa and 18.55 GPa, respectively, while the uniaxial compressive strength and modulus of elasticity of the cement specimens were 2.96 MPa and 17.80 GPa, respectively. Firstly, the metal tube was placed into the lower shear plate, the gypsum was poured into the pore between the lower part of the metal tube and the lower shear plate, then the upper part of the tube body was placed into the upper shear plate, the gypsum pouring process was repeated, and the shear specimen was finished.

2.1.2. Experimental System

The shear tests were conducted on an RMT-150C servo shear test device, as shown in Figure 4, which can test the shear mechanical behavior of the test specimen under different boundary conditions, and its key technical indicators include three main units. (1) hydraulic servo control system: The system is mainly composed of two jacks, which can apply an almost uniform normal stress in the shear plane. The system’s normal force and shear force are applied through a servo-controlled hydraulic oil pump and then the hydraulic cylinder is controlled. The maximum force of the vertical hydraulic cylinder was 1000 kN, and its piston stroke was 50 mm. The maximum force of the horizontal hydraulic cylinder was 500 kN, its piston stroke was 50 mm, and the deformation rate was 0.0001~1 mm/s. (2) Sensor unit: This system mainly includes two types of sensors—load sensors and displacement sensors. Among them, three digital load sensors are used to measure the shear and normal forces, including a tensile–compression-type sensor that can carry a 200 kN normal load, a shear-type sensor that can carry a 400 kN shear load, and a linear displacement sensor used to monitor the displacement. (3) Shear box unit: The shear box is divided into two parts: upper and lower shear box. In the shearing process, the upper shear box is fixed, while the lower shear box moves laterally to shear the test object. The shear test’s specific operational steps are as follows. Pre-preparation: number the specimen and mark the direction of the specimen shearing, and take photos to record the complete specimen before shearing. Specimen installation: first, loosen the screws of the upper and lower shear boxes, then put the specimen into the shear box and tighten the screws to ensure that the edges of the upper and lower shear boxes are approximately 5 mm each from the shear interface. Shear box installation: put the installed specimen and shear box into the shear track, and fix the upper shear box to prevent it from deflecting during the test; adjust the four vertical LVDT displacement gauges fixed on the tester frame to the appropriate position, and finally, install the lateral baffle of the shear box and horizontal LVDT displacement gauges. Shearing: Firstly, open the operation platform, start the oil pump and zero the data record. Next, apply a normal load at a rate of 0.5 MPa/min using the force control method, and keep the normal force constant after applying the normal load to a predetermined value. Then, the displacement control method is used to apply a shear load at a rate of 0.5 mm/min, and the target displacement is 20 mm. The computer records the shear stress, shear displacement, vertical displacement and other data throughout the test. End of the test: The shear stops when the shear displacement reaches the target value. First, unload the normal force to zero, and close the operating system. Remove the lateral baffle, pull out the shear box, and finally, take pictures of the damage of the specimen and collect the damage debris of the specimen. A set of parallel tests was set up after each group of tests, and the group with better results was selected for analysis.

2.1.3. Experimental Design

To verify the numerical simulation’s results, four sets of shear control tests under four normal stress conditions of 0.5, 1, 1.5 and 2 MPa were set up, and a set of parallel tests was set up under each test group. The group with the better results was selected for analysis. The specific experimental design is shown in

Table 1. Experimental design.

No.	Normal Stress (Mpa)	Test Result
Test group 1	0.5 MPa	Yes
Control group 1	0.5 MPa	No
Test group 2	1.0 MPa	Yes
Control group 2	1.0 MPa	Yes
Test group 3	1.5 MPa	No
Control group 3	1.5 MPa	Yes
Test group 4	2.0 MPa	Yes
Control group 4	2.0 MPa	Yes

In the shear test, the stress loading was performed in a stress-controlled mode with a loading rate of 1 kN/s. The normal force was applied in the normal direction of the specimen to a constant value, and then the tangential force was applied in the shear direction of the shear box at a loading rate of 1 mm/min to make the upper and lower discs of the specimen reach the maximum static friction (0.5 kN). Finally, the shear force was applied to the lower part of the shear box at a loading rate of 2 mm/min until the displacement of the specimen reached 20 mm or the specimen was damaged.

2.2. Numerical Modeling

Silica sand is a typical component of sandstone NGH deposits, and NGHs are endowed in the form of pore filling, bearing and cementation, forming a rough structure at the intersection of the NGHs layer and the overlying layer [37]. The consistency of the model settings in the numerical simulation with those in the shear test is a precondition to make the numerical simulation’s results validate each other with the shear test results. Therefore, the structural characteristics of the shear surface were obtained by 3D laser scanning of the shear surface of the shear test specimen, as shown in Figure 5. The 3D modeling software Rhino was used to model the 3D printed model of the upper and lower plates of the shear specimen, and the shear upper and lower plates with rough shear surfaces and the wellbore were imported into ABAQUS separately and assembled together using Boolean operations. Natural gas hydrate is generally endowed in the deep subsurface of terrestrial permafrost zones or in submarine strata at a certain depth from the seafloor with water depths greater than 300 m. Therefore, it is difficult to obtain the true shear strength of hydrate layers. Due to the specificity of hydrate fugitive conditions, the real shear strength in the engineering cannot be tested directly, and it needs to be scaled down and then modeled. However, regardless of the scaling method adopted, a scaling effect is bound to occur, i.e., the simulated results after scaling may not accurately reflect the true shear strength of the subsea wellbore [38]. In order to obtain more accurate simulation results, as close as possible to the engineering scale in hydrate mining, the initial numerical model was set at L × B × H = 200 m × 100 m × 140 m, the length of the borehole was 140 m and the borehole radius was 1.5 m. However, the number of meshes increased geometrically with the model size when meshing the model, which required too much computational power. Therefore, after preliminary trial calculations, the layer model was reduced to L × B × H = 2 m × 1 m × 1.4 m, the borehole length was shortened to 1.4 m and the borehole radius was reduced to 0.015 m. The model had three parts in total: shear upper plate, shear lower plate and borehole. The total number of units in the model was 11,246, of which the number of units in the shear upper plate was 5416, the number of units in the shear lower plate was 5712, the number of units in the wellbore was 136 and the computational model is shown in Figure 6.

Figure 6 is the surface morphology roughness scan of the shear surface, and the calculation model of the simulation is shown in Figure 7. The constitutive relation of the model is the Mohr–Coulomb strength criterion, as shown in Formulas (3). The model consists of solid units, and the tensile test results are adopted as the stress–strain relationship. Boundary conditions: shear displacement is constrained in the lower side of the model (X direction); normal displacement is constrained in the lower and upper sides of the model (Z direction); the bottom boundary of the model is fixed; and the top boundary is free. Loading: a load of 2.2 × 10⁶ kN was applied to the top of the model to simulate the subsea pressure of 11 MPa above the NGHs layer [16]; a load was applied in the shear direction (X direction) of the model at a rate of 1 mm/min to simulate the shear force of the shearing process. The reaction force on the velocity load acting on the surface during shear was monitored, and the horizontal displacement of the nodes on the loading surface was used to represent the shear displacement of the wellbore. Four observation points were set up at the load application point of the model’s shear lower plate and the upper, middle and lower parts of the wellbore to monitor the evolution of the stress–strain at different locations in the model at different times. The physical parameters of the wellbore and layer materials in the model are shown in

Table 2. Main parameters of the wellbore and rock layer in the numerical simulation [15,28,38].

Parameter	Value	Parameter	Value
Layer density, ρ/(Kg·m−3)	2600	Layer adhesive, strength/MPa	2
Poisson’s ratio of the layer, μ	0.5	Wellbore density, ρ/(Kg·m−3)	7800
Modulus of elasticity of the layer, EAO/MPa	3000	Wellbore elastic modulus, EAO/GPa	191.7
Friction angle of the layer, φ/(°)	30

. And the Figure 8 shows the simulation flow of the wellbore stability stages by ABAQUS analysis.

$$\sin \phi \le \frac{{\sigma }_1-{\sigma }_3}{{\sigma }_1+{\sigma }_3+2c\cot \phi }$$

(3)

where φ is the angle of the internal friction, c is the cohesive force, ${\sigma }_1$ is the maximum principal stress and ${\sigma }_3$ is the minimum principal stress.

2.3. Model Parameter

The lithology of overlying layers of NGH reservoirs is usually determined by the local geological properties. The NGH sediments are mainly muddy and sandy lithologies, and shear often occurs in sandy NGH reservoirs and leads to wellbore failure [28]. In actual oil and gas pumping projects, corrosion-resistant 304 steel is generally used as wellbore casing for production [36]. The main stress field calculation parameters for the layer and wellbore in the numerical simulation are as follows.

The density of the layer in the parameter table is the density of the shear upper and lower plates in the numerical simulation, and the relevant value of the hard sandstone in the reference was adopted, which was 2600 ρ/(kg·m⁻³). The Poisson’s ratio is the ratio of transverse positive strain to the axial positive strain of the shear upper and lower plates strata in the numerical model, which reflects the transverse deformation elastic constants of the layer in the model. The modulus of elasticity of the layer in the table is the ratio of the stress and the stress of the shear plates of the model in the elastic phase, which are proportional to each other. The angle of the internal friction of the friction angle of the layer in the table reflects the magnitude of the internal friction among the particles inside the shear plates of the model, which is an indicator of its shear strength. The layer adhesive strength in the table is the bonding force that the upper and lower plates can withstand per unit bond surface in the model. The wellbore density in the table is the density parameter of the wellbore part of the numerical model, using the relevant parameters of 304 steel. The modulus of elasticity of the wellbore in the table reflects the proportional relationship between the stress and strain in the elastic phase of the wellbore part of the numerical model.

3. Results Analysis

3.1. Model Validation

To verify the accuracy of the above numerical model of wellbore shear, the shear–displacement curves and wellbore failure characteristics of the laboratory shear experiments and numerical simulation were compared.

A comparison of the shear–displacement curves from the laboratory experiments and numerical simulations is shown in Figure 9. It can be seen that the trends of the shear–displacement curves obtained from the experiments and numerical simulations are basically the same and can be generally divided into five stages: Stage I (compaction stage), Stage II (linear stage), Stage III (plastic stage), Stage IV (strain-softening stage), and Stage V (residual stage). Stage I (compaction phase): this stage starts from the initial point of shear to T_a, and the curve of this stage is slightly concave, reflecting the closure phenomenon of the shear surface during the shear. Stage II (linear stage): this phase, from T_a to T_b, is an approximately linear segment that reflects the elastic behavior of the shear surface under compression–shear loading. Stage III (plastic stage): the curve is concave down to the peak shear point T_p, reflecting the process of local microfracture generation and development until damage occurs in the contact area of the shear upper and lower plates and wellbore during shear. Stage IV (strain-softening stage): the shear load gradually decreases with an increasing shear displacement, exhibiting an obvious softening phenomenon and reflecting the brief post-peak process of the wellbore and the shear surface. Stage V (residual stage): the curve in this stage is approximately horizontal and reflects the final or residual strength of the shear surface.

In the first two stages, the mechanical behavior is approximately elastic and exhibits a slight hysteresis, and neither loading nor unloading of shear stresses produces irreversible changes in the structure of the shear surface. The T_b point is the yield point, which corresponds to the yield stress of the wellbore and shear surface. After the yield point, irreversible plastic deformation will occur in the wellbore and shear surface.

The shear deformation and breakage characteristics of the wellbore are shown in Figure 10. As can be seen from the figure, the wellbore is subjected to tensile and shear effects, and a plastic hinge with nearly oblique symmetric distribution is formed on each side of the shear surface of the upper and lower plates. The deformation characteristics of the wellbore after shear breakage were all “Z” shaped. In addition, the deformation of the wellbore at the shear surface and the shear breakage displacement are consistent with the wellbore shear–displacement curve in Figure 9. The shear stress of the wellbore disappears due to the shear breakage of the wellbore, causing the curve to sink in stages and remain around a value until the end of the shearing process. Different degrees and characteristics of the deformation and breakage of the wellbore at the shear surface occur because the wellbore extrudes the surrounding rock in the shearing process to cause damage. The location where a large degree of irregular extrusion deformation of the wellbore occurs may be at the maximum bending moment (i.e., at the plastic hinge) or at the intersection of the wellbore and the shear surface. As can be seen from Figure 10, the location of the major deformation of the wellbore occurs at the plastic hinge, and the degree of deformation decreases gradually from the ends of the plastic hinge to the shear surface.

On the one hand, the normal stress acting on the shear surface will cause the deformed and broken wellbore to be partially embedded in the surrounding rock, and the wellbore displacement will be limited. On the other hand, a large tensile stress will appear in a side of the plastic hinge of the wellbore. The combined effect of the two aspects causes the flexural failure of the wellbore. The position of the wellbore breakage occurs at the junction between the wellbore and the shear surface, which is due to the large degree of deformation of the wellbore gradually concentrated near the shear surface under the action of normal stress, and the wellbore gradually tends to undergo shear failure.

From the deformation and failure characteristics of the wellbore shown in Figure 10, the deformation and shear breakage of the wellbore are mainly influenced by the length of the wellbore between the plastic hinges and the degree of deformation of the wellbore in the normal direction of the shear surface (related to the dilatancy deformation of the shear plane). On the one hand, since rocks have the property that their compressive strength increases as the confining pressure grows, the plastic hinge distance of the wellbore decreases as the compressive strength of the rock increases with the normal stress. On the other hand, the dilatancy deformation of the shear surface of the rock layer decreases with the application of normal stress. In this paper, the mechanism of deformation and the shear breakage displacement law of the wellbore are explained from the above viewpoint based on the shear displacement evolution of the wellbore, shown in Figure 11 as follows.

In the initial compaction stage of shear, the shear surface starts to close under the normal stress, and the shear displacement is basically in the range of 0.9 mm. As the shear reaches the elastic stage, tensile yielding occurs at the junction of the wellbore and the shear surface, and the wellbore deforms to a large extent, with the shear displacement increasing from 0.9 mm to 2.4 mm. At this time, the normal stress in the shear surface remains basically constant with the increase in the shear displacement. With the same shear displacement, the plastic hinge distance and shear surface expansion deformation are relatively larger under small normal stress conditions compared to larger normal stress conditions. The deflection of the plastic hinge part of the wellbore along the shear direction of the shear surface is conducive to mobilizing the tensile-bearing capacity of the wellbore to resist the external shear load on the shear surface to cause a larger shear displacement of the wellbore.

When the shear reaches the plastic stage, the shear displacement of the wellbore increases from 2.4 mm to 3.4 mm, and the normal stress at the shear surface increases with the increase in the boundary normal stress. The distance between the plastic hinge and the shear surface decreases, resulting in less dilatancy deformation at the shear surface with the same shear displacement, which means that the wellbore deformation along the shear surface in the normal direction decreases and the wellbore requires a larger shear displacement to reach its ultimate breakage load. Therefore, the shear breakage displacement of the wellbore at this stage tends to increase gradually with the increase in the normal stress, and the dilatancy deformation of the shear surface dominates the changing trend of the wellbore deformation at this stage.

The strain-softening stage of the shear begins at the peak point, T_p. The normal stress at the shear surface is greater at this point, and the wellbore displacement increases from 3.4 mm to 11.9 mm during this stage, which is reflected in the rapid increase in wellbore breakage and shear displacement. The distance of the plastic hinge of the wellbore dominates the trend of the shear breakage displacement of the wellbore. As the normal stress increases, the bearing capacity of the layer is strengthened, which in turn causes an increase in the forces applied to the wellbore, ultimately leading to a decrease in the distance between the plastic hinge and the shear surface. The wellbore yields in bending at the plastic hinge, which is manifested as bending failure. Eventually, the shear failure of the wellbore occurs at the shear surface. The last stage of shear is the residual stage, where the shear displacement increases from 11.9 mm to 20 mm.

3.2. Force Distribution and Evolution of the Wellbore

The evolution of the axial force, shear force and bending moment in the shearing process of wellbore in the NGHs layer can well reflect the force characteristics and deformation patterns of wellbore. In this paper, the wellbore shear at the shear surface, with a roughness value of 9.52, was used as an example, and the shear force, axial force and bending moment of the wellbore at different shear displacements were plotted to obtain the force distribution and evolution curve during the shear process, as shown in Figure 12, where the positive and negative values represent the direction of each physical quantity.

From Figure 12a, it can be seen that the wellbore shear force was basically symmetrically distributed along the shear surface in an approximate “M” shape. The wellbore shear force reached its maximum at the shear surface, gradually decreased in the direction away from the shear surface and reached a minimum value of 0 at the plastic hinge yield point. In addition, the shear force distribution showed a trend of increasing and then decreasing. With the increase in the shear displacement, the overall shear force of the wellbore showed a tendency to increase, and the shape of the shear force distribution also developed from a “flat M” type to “vertical M” type. When the shear displacement was less than 3.4 mm, the shear force near the shear surface grew rapidly, but the shear force outside the yield point grew slowly. After the shear displacement was greater than 3.4 mm, the wellbore shear force at the shear surface was maintained at a maximum value of 25.6 kN and did not change significantly due to the yielding of the wellbore. However, the wellbore shear force in the direction away from the yield point of the shear surface maintained a high growth rate, and there was an obvious peak shear force. At this point, the distribution of the shear force was more uneven, and a change in direction occurred at the two end positions of the wellbore. In fact, the evolution of the shear force distribution of the wellbore is closely related to the evolution of the plastic hinge during the shear process. After the generation of the plastic hinge, the overall deflection trend of the wellbore caused a change in the direction of the shear force at the two end positions of the wellbore.

As shown in Figure 12b, the wellbore axial force was more uniformly distributed in the initial stage of shear, and its curve shape was close to horizontal. During this stage, the value of the wellbore axial force was small, and the maximum value was only 2.4 kN at a shear displacement of 2.5 mm. However, as the shear displacement increased, the axial force at the shear surface of the wellbore started to be gradually higher than the axial force at both ends of the wellbore, and the distribution of the axial force along the length of the wellbore was in the form of “high in the middle and low on both sides”. In addition, the larger the shear displacement, the more significant the phenomenon. This is due to the fact that the shear force between the wellbore and layer interface will counteract part of the axial force when the wellbore axial force is transmitted from the shear surface to the ends. Moreover, the wellbore axial force increased with increasing shear displacement until it suddenly dropped to 0 after the wellbore breakage. The axial force of the wellbore at the shear surface was 4.8 kN, 5.9 kN, 7.1 kN, 8.2 kN and 9.2 kN for shear displacements of 3 mm, 3.4 mm, 6 mm, 8 mm and 12 mm, respectively.

From Figure 12c, it can be seen that the bending moment of the wellbore was basically distributed along the shear surface, with a skew–symmetric wave shape. The bending moment of the wellbore had a minimum value of 0 near the shear surface, increased gradually along the direction away from the shear surface and reached a maximum at the yield point. After, the distribution of the bending moment of the wellbore showed a decreasing trend when the shear displacement was less than 3.4 mm, while the direction changed when the shear displacement was greater than 3.4 mm. In fact, the variation pattern of the bending moment distribution was consistent with that of the shear force, which is due to the fact that the bending moment of the wellbore is mainly determined by the shear force to which it is subjected. In particular, the bending moment of the wellbore grew rapidly near the yield point when the shear displacement was less than 3.4 mm, but the bending moment grew slowly near the two end positions of the wellbore. In addition, after the shear displacement was greater than 3.4 mm, the maximum value of the wellbore bending moment near the yield point was kept at 1.92 kN·m. However, the wellbore bending moment in the direction away from the yield point at the shear surface will keep a higher growth rate until the peak in the reverse direction appears, and the maximum reverse peak reached 1.06 kN·m when the shear displacement was 12 mm.

The variation of the axial and shear forces of the wellbore during shear can be explained by two laws. Law 1 shows the variation of the axial force and shear force in a single section along the radial direction of the wellbore as the shear displacement increases, as shown in Figure 13. Before the shear starts, the axial force of the wellbore is uniformly distributed along the axial direction with a value of 5 kN. When the shear displacement was 0~4 mm, the axial force of the wellbore was basically the same as when the displacement was zero, without significant change. When the shear displacement reached approximately 3.4 mm, the axial force of the wellbore near the shear surface started to increase gradually with the increase in the shear displacement, and the axial force at the cross-section increased from 5 kN to 14 kN. The wellbore axial force at the shear surface was slightly greater than that at other sections, and there was the characteristic of the uneven distribution of the wellbore axial force along the axial direction at this time. However, this phenomenon was not yet obvious, and the overall distribution was still uniform. The axial force of the wellbore tended to decrease from this point, but the axial force near the shear surface was still greater than at other locations, which indicates that the wellbore had actually yielded and damaged and entered a failure state at this point. This change of wellbore axial force can be used as an important indicator of wellbore failure criterion. Subsequently, the wellbore axial force slowly decreased after reaching a peak until the end of the shearing process. As shown in Figure 13, the shear force near the shear surface of the wellbore was 0 at the start of shear. The shear force near the shear surface increased rapidly from 0 to a peak of 25.6 kN during the compaction stage, elastic and plastic phases of shear, decreased to 15 kN during the strain-softening stage, and then remained near this value from the residual stage to the end of the shear process. The main reason is that the wellbore is subjected to active normal thrust from one side of the rock (the lower left of the four orientations at the intersection of the wellbore and the shear surface) and passive normal resistance from the other side of the rock (the upper right of the four orientations at the intersection of the wellbore and the sliding surface) during the shear process. After a downward trend of shear force, it showed fluctuation in the strain-softening stage. The main reason is that the elongation variable of the wellbore is increasing rapidly when the shear force rises again, and the cross-section of the wellbore has moved away from the shear slip surface into the separation gap between the wellbore and the rock body as the large axial deformation of the wellbore is generated, which makes the normal resistance of the rock body on the right side of the wellbore decrease.

Law 2 is the variation law of the axial force and shear force along the axial direction of the wellbore, as shown in Figure 14. Although the axial forces of the wellbore are uniformly distributed in the axial direction before the start of shear, the wellbore is deformed and subjected to the compression of the shear surface, i.e., the upper and lower plates of the shear, which generates shear forces in the cross-section of the wellbore and causes the wellbore to elongate under the shear strain. At the same time, the component of the normal force in the vertical direction causes the axial force in the wellbore to expand along the axial direction, even in the unbonded state, resulting in an uneven distribution of the axial force along the axial direction. In fact, the surrounding rock of the wellbore has been damaged by extrusion, so there will be friction with the wellbore, which will generate shear stresses on the surface of the wellbore along the axial direction and cause an uneven distribution of the wellbore axial force along the axial direction. In summary, in the initial stage of shear, the original stress (i.e., subsea ground stress) and the action of the axial force exert shear resistance. In addition, in the later stage of shear, the axial force increment exerts the main shear resistance until the wellbore enters a plastic state and finally breaks down.

3.3. Wellbore–Layer Interactions during Shear

In the shearing process, the wellbore and the surrounding rock will enter a plastic state, and the expansion of the plastic zone can reflect the internal stress redistribution process of the wellbore, which is the “gestation” stage before wellbore failure. At this stage, microfractures start to develop inside the wellbore until the wellbore microfracture penetration is destroyed. Therefore, the evolution pattern of the plastic zone of the wellbore can well reflect the catastrophic process of the wellbore.

As can be seen from Figure 15, with the increase in the shear displacement, the wellbore gradually enters the plastic state, and the initial plastic zone is between the two reverse bending points above and below the shear surface, i.e., the plastic hinge. With the further increase in the shear displacement, the plastic zone of the wellbore penetrates from the surface to the interior and then expands from the shear surface to both sides. However, its change range is limited, and it starts to expand slowly when it reaches the exterior of the plastic hinge, which is mainly due to the fact that the rock deformation is smaller at the farther side of the shear surface, the shear force is smaller and the wellbore is mainly in a tensile state. When the entire cross-section of the wellbore is in a plastic state, the phenomenon of “necking” of the wellbore begins to occur, accompanied by the failure of the wellbore, which can be used as a critical indicator of wellbore failure.

As can be seen from the Figure 15, the plastic zone first appeared at the junction of the shear surface and the wellbore, where the stresses were concentrated and the plastic breakage zone appeared first. In addition, the large deformation increased the contact area between the wellbore and the rock and dispersed the interaction force between them. On the other side of the wellbore, the contact surface between the wellbore and the surrounding rock separated so that the force between the two was zero, and no plastic zone damage was generated. Since the stiffness of the wellbore was less than that of the rock, the compression deformation of the wellbore was larger, and the plastic zone expanded faster than that of the rock. The deformation of the rock body was small, and the plastic zone expanded slowly and to a small extent. The rock material was sensitive to the deformation of the wellbore, and the degree of the “Z-shaped” deformation of the wellbore determined the proportional relationship between the axial and shear forces of the internal force components of the wellbore, i.e., the failure mode of the wellbore.

3.4. Energy Analysis during Shear

As is known from the laws of thermodynamics, energy transformation is the essential feature of the physical process of matter, and matter destruction is a state instability driven by energy. During the shearing process, the movement of the shear upper and lower plates and the deformation and failure of the wellbore are always accompanied by different energy input, transmission and transformation, so the study of the energy evolution during the shearing process can well reveal the shear deformation and failure characteristics of the wellbore and its surrounding rock. In order to study the energy variation pattern of the system during shearing, three types of energy closely related to the shear characteristics, namely, elastic strain energy, frictional dissipation energy and plastic dissipation energy, were tracked and recorded in the numerical simulation, as shown in Figure 16. According to the characteristics of the shear force and energy curves, the evolution process is divided into three stages: Stage I (elastic stage), Stage II (pre-strain softening stage) and Stage III (late strain-softening residual stage). Stage I (compaction, elastic and plastic stage): This stage starts from the initial point of shear to the first peak value of the shear force, T_p. At this stage, the shear upper and lower plates and the wellbore begin to deform under the action of the shear force, and the work performed by the shear force on the system is mainly converted into the elastic strain energy and stored, while the frictional and plastic dissipation energy remain basically unchanged. Stage II (pre-strain softening stage): This stage is from the first peak of the shear force to the first breakage of the wellbore, and the shear force at the break of the wellbore is T_f. Under the action of the shear force, the elastic strain energy in the system continues to rise to its peak. It is known from the laws of thermodynamics that energy dissipation is the essential property of rock and wellbore deformation and failure, which reflects the process of continuous development, weakening and leading to the eventual loss of the strength of the micro-defects inside the rock and wellbore. Therefore, energy dissipation is directly related to the damage and strength loss of the wellbore, and the dissipation energy reflects the degree of original strength decay. As can be seen from Figure 12, the frictional dissipation energy and plastic dissipation energy increase rapidly in this stage under the action of the shear force, reflecting the continuous development of the internal defects and overall strength loss in the shear upper and lower plates and the wellbore during this process. Finally, the wellbore undergoes shear breakage, and the elastic strain energy is released at point T_f. Stage III (late strain-softening residual stage): The shear force–displacement curve in this stage shows a trend of bench sinking, and the shear force drops to a lower value and then basically remains around a value until the end of the shear process. This is due to the fact that when the shear deformation reaches a certain value, the wellbore deformation exceeds its limit, leading to its breakage, and the shear resistance in the wellbore disappears, which in turn leads to a stepwise slow sinking of the shear force. After, the shear strength of the shear surface is mainly maintained by the shear force generated by the shear surface itself, so it will remain around a value until the end of the shear process. During this process, the strain energy continues to decrease and reaches its lowest point at the T_e, and the growth rate of the frictional and plastic dissipation energy also begins to slow down until the end of the shearing process.

4. Conclusions

In this paper, ABAQUS finite element numerical simulation software was applied to carry out a numerical simulation of wellbore shear induced by the slip of deep sea NGHs layers, and the numerical simulation model was validated by laboratory wellbore shear experiments. Most of the research of numerical models related to wellbore stability studies in existing references focus on the effects of factors causing hydrate decomposition, such as temperature and liquid density on wellbore stability, without further analyzing the shear deformation characteristics of the wellbore from a mechanical perspective. The model in this paper considered the effect of the layer shear surface roughness on the wellbore shear characteristics and focused on the stress distribution around the wellbore, the interaction between the wellbore and the surrounding rock mass and the development of the elastoplastic zone during shear while monitoring the energy evolution of the shear system. The research in this paper can provide some reference for the study of wellbore deformation damage during gas hydrate extraction and its prevention and control. The following main conclusions were obtained.

(1)

The shear–displacement curves obtained from the laboratory experiments and numerical simulations were basically consistent, with the five stages of the compaction stage, linear stage, plastic stage, strain-softening stage and residual stage, in general, and the wellbore showed “Z” deformation characteristics after shear breakage.

(2)

The wellbore shear force was the maximum at the shear surface, with an approximate “M”-shaped distribution along the shear surface and an overall increasing trend with the increase in the shear displacement, and the distribution shape developed from a “flat M” shape to a “vertical M” shape. The wellbore axial force was small and uniformly distributed in the initial stage of shear, and its distribution gradually showed the form of “high in the middle and low on both sides” with the increase in the shear displacement. The minimum value of the wellbore bending moment was 0 at the shear surface, which was distributed in an obliquely symmetric wave shape along the shear surface. The evolution of the wellbore axial and shear forces during shear can be divided into the distribution patterns along the radial direction at the shear surface and along the axial direction of the wellbore. The combination of the wellbore axial force and shear force caused the wellbore’s tensile–shear compound failure.

(3)

In the shearing process, both the wellbore and the rock body gradually entered the plastic state and formed a plastic hinge until failure with the increase in the shear displacement. When the entire cross-section of the wellbore is in a plastic state, the wellbore begins “necking”, which is accompanied by the failure of the wellbore. The plastic zone of the rock and the wellbore first appeared at the junction of the shear surface and the wellbore, where the rock stress was concentrated, forming a plastic crushing zone and generating large deformation.

(4)

With the change of shear stress, the elastic strain energy first increased, and the frictional dissipation energy and plastic dissipation energy started to increase rapidly after the shear force reached the peak, T_p. The internal defects in the shear upper and lower plates and the wellbore continued to develop until the wellbore broke at T_f, where the strain energy reached a peak and then released rapidly. The strain energy continued to decrease and reached its lowest point at the T_e point. At this time, the growth rate of the frictional and plastic dissipation energy began to slow down until the end of the shearing process.