Abstract
To solve the difficult problems in the field of unconventional oil and gas extraction and hard rock excavation in urban underground spaces, this paper proposes a rock-breaking technique with high-temperature and high-pressure water under thermally driven conditions and establishes a coupled thermal-fluid-solid model with COMSOL Multiphysics (COMSOL Co., Ltd. Shanghai, China). Different simulation groups are established by controlling variables to explore the effects of the surrounding rock load, heat source power, and Biot coefficient on the damage evolution during thermally driven rock breaking. To make the results relevant to practical engineering, the damage evolution results under the maximum normal stress criterion, maximum normal strain criterion, and Coulomb-Navier damage criterion are considered, and a comparative analysis is performed. The results of this study show that an increase in unilateral load and heat source power accelerates the damage evolution rate, while an increase in bilateral load and Biot coefficient has the opposite effect. The damage evolution rate controlled by the maximum normal stress criterion is the fastest under general conditions. Finally, the advantages in rock breaking provided by the established method are verified by a comparison of results from the proposed model and a conventional hydraulic fracturing model.
1. Introduction
With the development of science and technology, the exploration of the Earth’s interior is progressing deeper. The problem of breaking hard rock has been paid more and more attention. Effective rock-breaking technology developing important at the time, especially in the processes of low-permeability oil and gas resource exploitation and fracturing of deep rock mass, and in the area of complex and sensitive subsurface environments such as urban area.
Currently, the widely used rock-breaking methods include tunnelling with blasting, TBMs, and various new rock-breaking technologies (such as water jets [1], lasers, microwaves [2], particle impacts [3], and dry ice powder pneumatic rock breaking [4, 5]). Of these, blasting is the earliest proposed rock-breaking technology with the highest construction efficiency, but its application is limited in many scenarios due to dangerousness. From the rapid development of technology by Robbins in 1956, TBM technology, which is now relatively mature, is still not sufficiently adaptable to complex geological conditions and protective against geological hazards in the actual engineering applications undertaken in China. Specially, in ultrahard rock formations, the cost of time and economy is too high for changing worn-out cutter [6]. In terms of new rock-breaking technologies, water jets, lasers, and other technologies began to be more commonly used in the second half of the last century, but they are still in the laboratory stage [1] [7, 8]. Furthermore, the shale gas revolution in the United States driven by the emergence of hydraulic fracturing [9] has prompted many scholars to conduct related research on crustal stress [10, 11], natural cracks [12, 13], fracturing fluid viscosity [14], bedding orientation [15], and other factors affecting the fracturing effect. In recent years, supercritical carbon dioxide rock breaking has attracted increasing attention because it can prevent the occurrence of the water blocking effect, promote oil and gas production, lower the pressure threshold [16] and easily form crack networks [17], where Du at el. [18]. proved the advantages of supercritical carbon dioxide jets over high-pressure water jets in deep rock breaking through experiments and obtained the relationship between supercritical carbon dioxide and nozzle diameter, pressure threshold, jet pressure, rock compressive strength, and jet temperature. Wang et al. [19] achieved the industrialization of supercritical carbon dioxide fracturing by combining laboratory and field tests. Li et al. [16] obtained the influence of the temperature field and rock parameters on the resulting damage by establishing a thermal-fluid-solid coupling model, while Zou et al. [20] and Song et al. [21] explored the influence of supercritical carbon dioxide on rock parameters. Although most supercritical CO2 fracking cases are not at ultra-high temperatures, but it is worth mentioning that the increase in temperature has a great impact on the increase in permeability and porosity [22]. Compared with hydraulic fracturing, supercritical carbon dioxide fracturing has a great advantage in that it has a lower viscosity and higher expansion rate than conventional fracturing fluids. However, it can also be achieved during hydraulic fracturing by changing the temperature and pressure conditions in situ. The influence of temperature on rock fracturing is significant. Therefore, the effect of high temperature and pressure on hydraulic fracturing is worth looking forward to. Due to the difficulty of performing experiments in this field and relevant lack of data, this paper will investigate the rock-breaking mechanism by a numerical simulation method.
At present, there have been a large number of numerical simulation studies on the factors affecting hydraulic fracturing, including rock permeability, loading rate (including water injection rate and pressurization rate), fluid viscosity, and rock anisotropy [17] [23–25]. In addition, Yi et al. [26] considered the influence of damage evolution on pore elastic parameters when studying hydraulic fracturing based on the phase field method. Nguyen et al. [27] proposed a simulation method for hydraulic fracturing with good convergence and low calculation cost. Furthermore, the impacts of the thermal effect on hydraulic fracturing have been studied through numerical simulation, including the impacts of the natural ambient temperature of the rock mass [28], the effect of artificial heat treatment [29], the improvement of relevant numerical methods [30], and even the surrounding thermal field of supercritical water thermal jets [31]. Therefore, in terms of numerical simulation, the foundation is laid through the previous research, and the rock-breaking problem of high-temperature and high-pressure fluid under thermal driving conditions is studied on this basis.
At the moment, there is no effective rock-breaking method to supplement and replace the TBM and blasting while the novel thermally driven rock-breaking method involving high-temperature and high-pressure water presented in this work is possible to achieve it. This method takes advantage of the low viscosity, high expansion rate, and high diffusivity of high-temperature fluids and gets satisfactory fracturing results. Until then, some scholars have studied the influence of temperature on rock breaking, but these studies are limited to cracking caused by crude oil viscosity reductions [32] or temperature differences [29]. Therefore, on the basis of relevant previous work, this paper studies the mechanism of the rapid rock-breaking process by high-temperature and high-pressure fluid under thermal driving through numerical simulation. The study will further promote the development of rock-breaking technology, thus driving innovation and breakthroughs in deep development and unconventional oil and gas resource development.
2. Simulation Details
This paper develops a study of rock breaking by high-temperature and high-pressure thermally driven water through a multiphysical field coupling model of thermal-fluid-solid damage. The 2D model is established according to a liquid-filled cylinder and concrete test blocks with holes and energy-accumulating agent mixed with fluid which provide thermal energy after being excited. The exterior of the model is rock, the interior is water, and the center is a heat source. The volume of the energy-accumulating agent is negligible; thus, the heat source is set as a fluid. In addition, the model simplifies the pipe wall as a boundary to study the direct action of water on rock.
2.1. Model Geometry
The size of the 2D model is 200 mm200 mm, with an axis of symmetry in the middle. The horizontal axis represents the horizontal direction, and the vertical axis represents the vertical direction (as shown in Figure 1).
2.2. Governing Equations
This simulation includes water phase change, heat and mass transfer, nonlinear flow in fractures, rock damage, and other processes, involving multiple physical fields such as flow, seepage, temperature, and solid mechanics. The main governing equations and coupling equations are as follows: (1)Correlation equation of elastic-brittle mechanics
In this model, rock damage studied is assumed to be elastic-brittle damage, and Hooke’s law is used as the instanton equation: where is the prestress and is the elastic matrix.
In the damage model established with COMSOL, the material damage is determined by the exponential strain-softening relationship, in which the independent variable describing the damage is equivalent strain. When the strain threshold is reached, the damage begins to evolve. In addition, the effects of fracture energy, induced fracture characteristic size, rock tensile strength, and rock elastic modulus on damage are considered by strain-softening parameters, which are expressed as follows: where is the fracture energy per unit area of the material, is the tensile strength of the material, and is the characteristic size of the induced fracture.
Damage is defined as follows: where is the equivalent strain, is the strain threshold, and eps is equal to .
In addition, to determine the stress and strain of the solid, the equilibrium differential equation needs to be added. After considering the inertia term, it is as follows: where is the solid density, is the solid displacement, is the unit volume force, and is the time. (2)Fluid equation
The central fluid domain of the model is controlled by the N-S equation: where is the fluid density, is the fluid velocity, is the fluid pressure, is the viscous stress tensor, is the volume force, and is the gravity term. The flow of fluid in rock pores and fractures also has strong nonlinearity. Therefore, Darcy’s law is not adopted, but the porosity and permeability terms are considered in the N-S equation, and the combined expression is as follows: where is the porosity, is the permeability, is isothermal compressibility, and is the source term.
Considering the changes in rock porosity and permeability in the process of damage, these two parameters are defined by damage [33]. Due to the rapid nonlinear increase of both during fracture propagation, the linear exponential function expression of the natural constant is adopted, and the permeability is as follows: where is the coefficient, which is set to 5 in this paper.
According to the famous K-C equation [34–36]: where is the specific surface area of the solid phase. When is approximately equal to 1, , and the porosity is as follows:
After combining damage and flow, coupling between the internal fluid and the rock is still not achieved. This is achieved by using the fluid pressure at the inner surface of the rock as a load.
Finally, the temperature field is introduced to more closely reflect the real thermally driven rock-breaking process. The flow under nonisothermal conditions is represented by the energy equation, which is expressed as follows: where is the temperature, is the thermal conductivity, and is the power of the heat source.
The above equations are used to connect the temperature field, stress field, damage, flow field, and seepage field to realize the transient calculation of the whole process of hydrothermal shock-driven rock breaking during the simulation.
2.3. Initial and Boundary Conditions
The model is initially under a static state at room temperature. Therefore, the initial temperature is set to 293.15 K, the pressure of the flow field and seepage field is 1 atm, and the displacement and velocity of the solid structure are 0.
For the plane boundary conditions of the solid mechanics field, the vertical load and horizontal load are set as and , fixed constraints are set at the bottom, and the fluid pressure is set at the inner boundary. For the flow and seepage field, the outlet pressure boundary is set at the boundary with a pressure of 1 atm. In the temperature field, since the temperature variation of the solid mechanics field is small, the outer boundary is set as an adiabatic boundary.
2.4. Material Parameters and Meshing
The specific values of the water parameters are from NIST. The rock parameters adopted are as shown in Table 1.
The model adopts a quadrilateral mesh (as shown in Figure 2), which is more suitable for highly nonlinear problems and has stronger convergence [37]. In addition, to obtain better initial conditions and consider the prestrain caused by the boundary load in two directions, the model is used for both steady-state and transient research. The steady-state research calculates the change in only the solid mechanical field under loading.
3. Results
This chapter describes the role of various factors in damage evolution by changing the parameters of the established model and using the control variable method.
3.1. Simulation Set-up
The model simulates the effects of the confining pressure, heat source power, and Biot coefficient on damage evolution under different parameters by changing the boundary load, source term, coupling expression, and material parameters. The set simulation group results are shown in Table 2.
In all the above-mentioned simulation groups, the failure criterion is the maximum tensile stress criterion. The damage evolution of group O is shown in Figure 3.
(a) =2 s
(b) =5 s
(c) =10s
(d) =15 s
The damage extends radially along the initial damage direction, which is consistent with the diagonal direction when the confining pressures on both sides are equal (see Figure 4) The damage development path for group A3 with different loads on either side is roughly consistent with the corresponding test results of Liu et al. [38], which indicates the reliability of the model.
According to the damage curve of group O (Figure 5), the evolution process is divided into the pressure rise stage, damage initiation stage, damage stagnation stage, and damage rise stage. Among them, the damage stagnation stage is caused by the decrease in pressure after damage.
3.2. Influence of Surrounding Rock Pressure on Damage Evolution
When the confining pressure condition is changed on only one side, the results shown in Figure 6 are obtained. The curve shows that when the load on one side is changed, the damage growth rate will increase with the increase in the larger (or equal) load on one side. In this process, the length of the damage stagnation stage decreases with increasing load, indicating that the time from initial crack formation to crack re-expansion will shorten with increasing confining pressure on one side. This result is independent of whether the pressure is applied horizontally or vertically but is related to the magnitude of pressure.
It can be seen from Figure 7 that under the condition of an equal on both sides but simultaneously changing loads, due to the relatively uniform distribution of stress, a totally different result is observed in the fluid domain; that is, the greater the load is, the higher the fracturing difficulty and the slower the damage evolution speed. This result is also consistent with the experimental results of Ge et al. [39] obtained under a low confining pressure. In addition, similarly, when the pressure on both sides is less than 1 MPa, the damage stagnation stage extends with increasing confining pressure, that is, a result opposite to that observed under the condition of a dominant loading side.
3.3. Influence of Heat Source Power on Damage Evolution
When the heat source power is greater, due to the rapid rise in fluid temperature, the volume expansion and pressure rise are also faster, and the damage degree at the same time is greater. Therefore, this paper studies the influence of heat source power itself by comparing the damage caused by different powers at the same temperature. As shown in Figure 8, when the heat source power changes, the damage rate will increase with increasing power. Even at the same temperature, the damage driven by a higher heat source power is greater. The damage curve of B1 is quite different from the damage curves of the other groups. This is further investigated, and it is found that when the heat source power is reduced to a certain level, damage does not occur. In addition, group B3 and group B4 show a more pronounced temperature drop, as well as a larger extent of damage, during damage development than the other groups do, suggesting that the damage would develop more rapidly under larger temperature fluctuations.
Figure 9 shows that there is an obvious relationship between the changes in temperature and flow rate due to the treatment and that the temperature fluctuation occurs slightly earlier than the flow rate fluctuation, indicating that the rise in temperature causes the acceleration of the flow rate. Figure 10 shows that the first temperature drop corresponds to the occurrence of initial damage, reflecting that the damage propagation causes the decrease in temperature. However, it is uncertain whether the expansion rate is positively correlated with the temperature drop.
3.4. Effect of Biot Coefficient on Damage Evolution
The Biot coefficient, also known as the effective stress coefficient, is related to the pore size of the rock. It is between zero and one. It reflects the influence of the pore fluid pressure on the solid skeleton. The expression of the effective stress considering the Biot coefficient is as follows:
In this equation, is the effective stress, is the total stress, is the Biot coefficient, and is the pore pressure.
The model is used to investigate the damage evolution results with different Biot coefficient values for the same rock under the same confining pressure conditions. It can be clearly seen from Figure 11 that both the time of initial damage occurrence and the damage growth rate of the damage rising stage increase with the Biot coefficient. When the Biot coefficient decreases to 0.7, the damage evolution rate decreases significantly. This can be explained by the Mohr-Coulomb circle (Figure 12). When the Biot coefficient is larger, the threshold of effective stress is smaller, and damage is more likely to be caused.
4. Analysis and Contrast
The third section shows the simulation results of rock breaking by high-temperature and high-pressure water under thermal driving, but the causes and internal mechanism are not discussed. This section further analyzes the effects of the surrounding pressure conditions, heat source power, and Biot coefficient and compares the damage evolution results under different damage criteria. In addition, thermally driven fracturing is compared with conventional hydraulic fracturing.
4.1. Analysis of Different Factors Affecting the Results
The length of the damage arrest stage essentially reflects the damage evolution rate at that stage, so both the overall damage growth rate and the length of the damage arrest segment are compared here. In the study of the effect of the pressure conditions, the damage growth rate increases with the increase in the large (or equal) unilateral load. Here, it is considered that the damage growth rate increases due to different stress fields caused by the surrounding rock pressure. At the same time, the distance from the hole to the boundary (75 mm) is greater than or equal to 1.5 times the hole diameter (50 mm). Therefore, the calculation method of the small orifice problem in elasticity is used to analyze the inner boundary stress. Since there is no radial stress or shear force on the inner hole boundary, the circumferential pressure is as follows:
Let ; the circumferential stress at the inner boundary is as follows:
When , corresponds to the minimum pressure . The greater the load on one side (or both sides), the smaller the minimum stress is; thus, when the value of is greater, the rock is more prone to damage in the process of thermal shock. When , corresponds to the minimum pressure . The greater the load on one side (or both sides), the smaller the value is, the lower the minimum stress is, and the more likely damage is to occur.
When the loads on both sides are equal, the circumferential pressure at the inner boundary is as follows:
Therefore, the increase in load will inhibit the evolution of damage.
In the research on the influence of heat source power, we should first pay attention to the second phenomenon, that is, when the heat source power is lower than a certain value, damage does not occur. Research in the field of hydraulic fracturing shows that when the water injection rate is lower than the damage threshold, continuous water injection will not make the rock crack [24]. The same is true for thermally driven rock breaking. The key is that the volume expansion and seepage of the fluid balance so that the pressure cannot rise further. This also explains the first phenomenon; that is, at the same temperature, the damage evolution rate driven by higher heat source power is faster. On the one hand, this is due to the rapid expansion brought on by the high power. On the other hand, it is also due to the smaller amount of seepage in a shorter time.
The results suggest that the temperature is not only the driving force of damage evolution but also subject to the counteraction of damage evolution. When a fracture further expands, on the one hand, the fluid expands to do work to release heat; on the other hand, it contacts cold rock and transfers heat, which reduces the fluid temperature. This relationship between temperature and velocity is also supported by molecular thermodynamics.
The results on the influence of the Biot coefficient can be explained by the Mohr-Coulomb circle. When the Biot coefficient is large and the Mohr circle moves to the left under the same fluid pressure, the pressure value required to achieve damage is smaller, that is, the rock is more prone to damage.
4.2. Damage Evolution of Rock under Different Failure Criteria
In this section, three failure criteria for brittle rocks, tensile stress failure, tensile strain failure, and Coulomb-Never criterion (compression-shear failure), are compared and analyzed for damage evolution.
As stated in Section 2, the prerequisite for the initiation of damage is that the equivalent strain is greater than the strain threshold. Under the criterion of tensile stress failure, the equivalent strain is expressed as follows:
The strain threshold is as follows: where is the tensile strength.
Under the tensile strain failure criterion, the stresses in the lateral directions are considered, and the equivalent strain is expressed as follows:
Equation (17) considers the influence of the principal stresses in the other two directions. Since the 3D principal stress state during damage is tensile, as shown in Figure 13, the damage is minimal.
Compressive shear damage considers both shear and compressive (positive) stresses, and the shear rupture surface of the rock is considered to be the most unfavourable surface for the combination of shear and positive stresses. Under the Coulomb-Neville fracture criterion, the shear failure condition can be expressed as follows:
Since the 3D principal stress state in the model is expressed as tensile and the Mohr circle is located on the left side of the -axis when approaching the damage envelope, the equation is transformed. After transformation, the following expression can be obtained: where is the first principal stress, is the third principal stress, is the internal friction angle, and c is the cohesion.
Therefore, the equivalent strain εq and strain thresholds ε0 are set as follows:
Substituting into equation (17) results in the following expression:
Comparing equations (15) and (21) and combining the simulation results (Figure 14), when the control parameters are equal, the damage rate under the Coulomb-Neville criterion is faster. However, its control parameter is usually much greater than the tensile strength. Therefore, the failure under this condition follows the tensile stress failure criterion.
4.3. Comparison of Damage Evolution of Conventional Hydraulic Fracturing
To ensure the comparability of conventional hydraulic fracturing and thermal drive fracturing, this paper selects the comparison group by considering the rate of increase in pressure and makes a comparative analysis mainly based on the stress intensity factor and pressure variation law. The difference between the two models is that the hydraulic fracturing group is pressurized by injecting fluid into the hole, while the thermal drive fracturing is pressurized by adjusting the heat source power. After adjusting various flow injection rate parameters, it is found that the rate of increase in the pressure of the hydraulic fracturing group with an injection flow of 200 ml/min is similar to that of thermally driven fracturing with a thermal power of 5e108; the former is greater than the 171.43 ml/min converted from the contact area in other models [40], and the latter easily occurs in actual projects. The pressure changes at the center point of the model Figure 15. In this paper, the stress intensity factors of the two are analyzed by the fracture mechanics method [41], and the stress intensity at the crack end is considered the superposition of the stress intensity factors caused by each stress: where is the stress intensity factor of the mode crack and () and () are the stress intensity factors caused by the horizontal external load and vertical external load, respectively. is the pressure acting on the pore, and indicates the fluid pressure from to , , , and .
The following can be deduced from the two models where and correspond to hydraulic fracturing and and correspond to thermally driven fracturing. Therefore, only the stress intensity factor caused by fluid pressure should be considered. The fluid pressure in intact rock is as follows: where is the hole wall pressure, and the approximate analytical expression of dimensionless constant is as follows: where.
The expression of the stress intensity factor caused by pressure in the crack is as follows:
In the formula, has different expressions according to different pressure gradients along the crack, and when the pressure drop is linear, the expression is as follows:
The cracking state is assumed to be , and the crack length and the corresponding stress intensity factor and their distributions are thus analyzed. The crack at a certain moment is divided into four sections, and the corresponding stress intensity factors are determined according to the parameters of the five corresponding nodes. For hydraulic fracturing, the crack length is roughly divided into five sections, and the parameters of the corresponding moment are taken to derive the results. For thermally driven fracturing, since the crack length is longer, the parameters of the five moments are taken at positions similar to the former positions on the basis of the five equal sections. The final stress intensity factor curves are shown in Figure 16, and it can be found that the two patterns are basically the same, while the intensity factor of thermally driven fracturing is slightly greater for the same crack length.
This explains why the thermally driven fractures are longer but does not explain why they are so much longer. Therefore, comparing the curves of the rates of increase in pressure, although the rates are similar, the pressure rise of thermally driven fracturing is faster and then slower compared with that of hydraulic fracturing, and there are two pressure rise segments. To make this phenomenon more obvious, this paper uses multiple simulations, increasing the rock tensile strength step by step, and the same two patterns mentioned above can be found in Figures 17 and 18. When the tensile strength of the rock continues to increase above 40.1 MPa, Figure 19 shows that hydraulic fracturing is unable to cause damage to the rock, while thermally driven fracturing is unable to cause damage until the tensile strength is 80 MPa, and the upper limit of fracturing is high and shows obvious multistage pressure increase characteristics, which is conducive to the expansion of cracks [41].
Clearly, it is rare to encounter rock with such a high tensile strength, and high ground stress conditions are more commonly encountered, so this paper also simulates and analyzes the fracturing situation under different surrounding pressure conditions (Figure 20). As shown in Figure 21, the water injection flow rate of this hydraulic fracturing simulation reaches 500 ml/min at the laboratory scale, which is difficult to achieve in actual engineering; a similar effect can be achieved with thermally driven fracturing.
(a) Thermally driven damage
(b) Conventional hydraulic damage
Advice on Tables
5. Conclusion
This paper investigates the effect of several factors on the effectiveness of high-temperature and high-pressure water breaking under a thermal drive by varying the surrounding rock load, heat source power, and Biot coefficient under the same damage criterion. Additionally, to increase the practicability of the research, the damage evolution results under different failure criteria are studied and finally compared with common rock breaking by water injection; the results show the advantages of this rock-breaking method.
In the study of establishing a thermal-fluid solid coupling model for the process of high-temperature and high-pressure water impacting rock, the effects of three factors on the damage rate are studied according to the control variable method. The results show that under the condition of plane stress, an increase in the unilateral load will accelerate the damage evolution rate, but under the condition of a bilateral loading, an increase in the bilateral load will change the distribution of stress and reduce the damage evolution rate. Regarding heat source power, when the heat source power increases, even at the same temperature, the damage evolution rate driven by the heat source power is faster; an equilibrium point exists between power and strength. When the heat source power is lower than this threshold, there will be no damage, and when the heat source power approaches it, the temperature at initial damage will rise rapidly. Furthermore, there is a positive correlation between the influence of the pore water pressure on the rock skeleton, that is, the Biot coefficient, and damage growth rate. When the Biot coefficient of rock is smaller, the influence of the change in Biot coefficient is more significant.
In addition, it is found that under the condition of this model, the damage evolution rate controlled by the maximum normal stress failure criterion is the fastest among the different failure criteria considered. The maximum normal stress failure criterion depends on the relationship between rock tensile strength and cohesion; the smaller of the two determines the failure criterion to which the damage follows.
In the comparison between conventional hydraulic fracturing and thermally driven fracturing, compared to the former, the latter can achieve higher peak pressures and exhibit umps in the peak pressure, both of which facilitate fracture expansion. In addition, thermally driven fracturing does not require continuous water injection, reducing water consumption.
Data Availability
Fluid parameters are from the https://webbook.nist.gov/chemistry/fluid/.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This study was supported by The Surface Project of the Natural Science Foundation of Jiangsu Province (No. BK20201313), The State Key Laboratory Open Fund Project (No. HKLBEF202004), and Science and Technology Program of the Ministry of Housing and Urban Rural Development (NO. 2021-K-087).