Phase field simulation method for quasi-static hydrofracture of thermoelastic porous medium under consideration of high-energy fracturing condition

Abstract

The invention belongs to the technical field of microwave rock breaking, and provides a phase field simulation method for quasi-static hydrofracture of a thermoelastic porous medium under the condition of high-energy fracturing, which comprises the following steps of: firstly, acquiring parameters of a target reservoir and parameters of fluid in the target reservoir; secondly, microwave parameters are set on the basis of target reservoir parameters and fluid parameters in a target reservoir, and a phase field model of quasi-static hydraulic fracturing of the thermoelastic porous medium under the microwave-assisted fracturing condition is established; then, the target reservoir is irradiated with the microwave parameters, and the phase field model is used for simulating the crack evolution process during hydraulic fracturing; according to the method, the crack initiation difficulty of hydraulic fracturing can be reduced through microwave-assisted heating.

Description

Technical Field

[0001] This invention relates to the field of microwave rock breaking technology, and in particular to a phase field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions. Background Technology

[0002] The global energy landscape is entering a phase where coal, oil, natural gas, and new energy sources will share the market. It is projected that by 2040, natural gas will account for 25%–28% of global primary energy consumption. China has also entered a crucial stage of energy structure adjustment and the development of multiple energy sources. In 2024, China's total natural gas production was 2463.7 × 10⁻⁶. 8 m 3 Shale gas and tight gas production were 254 × 10⁻⁶ respectively. 8 m 3 620×10 8 m 3 Shale gas and tight gas production account for over 30% of China's total natural gas production. Natural gas consumption will continue to grow rapidly as the best partner for new energy sources, with unconventional natural gas such as shale gas, tight gas, and coalbed methane remaining the main drivers of production growth. It is projected that by 2035, their production will be comparable to that of conventional gas. Therefore, the efficient development and utilization of unconventional natural gas plays a crucial role in ensuring China's energy supply. Currently, my country is gradually entering the stage of deep resource extraction. The extraction depth of conventional oil and gas resources such as oil and natural gas has exceeded 7,500 meters, while the extraction depth of unconventional natural gas resources such as shale gas and tight gas has exceeded 4,000 meters. Oil and gas exploration wells, represented by the Chuanke 1 and Take 1 wells, have even broken the 10,000-meter mark.

[0003] Deep unconventional oil and gas reservoirs are characterized by low reserves, poor physical properties, low natural productivity, high production costs, and significant challenges in increasing production. Most are difficult-to-develop reservoirs, requiring large-scale hydraulic fracturing to achieve industrial production capacity during exploration, evaluation, and development. However, the high-pressure, low-permeability geological characteristics pose significant difficulties for fracturing operations. Excessive wellhead pressure often exceeds the pressure limits of wellhead and downhole tools, making normal fracturing impossible. Therefore, effectively improving near-wellbore petrological properties and reducing formation fracturing pressure are crucial for safe fracturing. Currently, acid fracturing technology is commonly used for deep unconventional oil and gas reservoirs to reduce pressure by altering the rock's mechanical parameters and mineral composition. Its simplicity and rapid effectiveness have led to its widespread use. However, acid fracturing still has some drawbacks. For example, its application is typically limited to carbonate formations, and it carries the risk of environmental pollution.

[0004] The physical and mechanical essence of microwave-assisted deep unconventional energy reservoir transformation is mainly reflected in the following aspects: On the one hand, microwave heating can improve the desorption of unconventional natural gas and remove the water-locking effect in tight reservoirs, thereby promoting the connectivity of primary pores and making it feasible to use heat injection to assist in production enhancement; on the other hand, the thermal stress generated by microwave heating will degrade the mechanical properties of the rock itself, thereby reducing the fracturing pressure and assisting in the completion of hydraulic fracturing projects.

[0005] Current research on microwave heating-enhanced shale gas recovery mainly focuses on the ability of reservoir temperature to assist oil and gas migration. Only a few scholars have studied the characteristics of reservoir rock cracking due to microwave radiation. The microwave-assisted hydraulic fracturing process involves a complex electromagnetic-thermal-hydraulic-mechanical (ETHM) coupling mechanism. Therefore, even when using numerical simulation methods, there is a problem of low efficiency in predicting complex fracture propagation modes. Summary of the Invention

[0006] The purpose of this invention is to provide a phase field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions, which can reduce the difficulty of hydraulic fracturing initiation through microwave-assisted heating.

[0007] The technical solution adopted by this invention to solve its technical problem is as follows: A phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions includes the following steps: Obtain the target reservoir parameters and the fluid parameters in the target reservoir; Microwave parameters are set based on the target reservoir parameters and fluid parameters in the target reservoir, and a quasi-static hydraulic fracturing phase field model considering microwave-assisted fracturing of thermoelastic porous media is established. The target reservoir was irradiated with microwave parameters, and the crack evolution process during hydraulic fracturing was simulated using a phase field model.

[0008] In some embodiments, considering the energy of fluid pressure in the porous medium and the heat generated by microwave electromagnetic loss, a quasi-static hydraulic fracturing phase-field model for thermoelastic porous media under microwave-assisted fracturing conditions is established. This phase-field model is then represented by an energy functional, and an initial expression for the energy functional of this phase-field model is established, which is expressed as follows: in, This represents the total energy of the system in the phase-field model. For elastic strain, For displacement, For phase field variables, Relative temperature For the computational domain, For fracture domain, For elastic energy density, For biot coefficient, Pore ​​pressure, For the density of the medium, For specific heat capacity, For heat flux, The critical energy release rate. The traction force acting on the Newman boundary, For the computational domain boundary, For volume-uniformly distributed force, For computational domain volume.

[0009] In some embodiments, after establishing the initial expression of the energy functional of the phase-field model, the total strain of the thermoelastic porous medium is determined. The total strain is the sum of the elastic strain and the thermal strain, wherein the thermal strain is expressed as... Indicates total strain. express; The total strain The expression is: thermal strain The expression is: ,in, The coefficient of thermal expansion is... It is a second-order unit tensor; After determining the total strain of the thermoelastic porous medium, the final expression of the energy functional of the phase field model is established.

[0010] After determining the total strain of the thermoelastic porous medium, the final expression of the energy functional of the phase field model is established.

[0011] In some embodiments, establishing the final expression of the energy functional of the phase-field model includes the following steps: Assuming crack propagation is driven solely by tension, and the elastic strain is decomposed into tensile and compressive strain tensors such that cracks only occur under tension, the decomposition of elastic strain into tensile and compressive strain tensors is expressed as: in, These are the tensile and compressive elastic strain tensors, respectively. It is the principal elastic strain. It is the direction of the principal elastic strain, and the operation symbol. Defined as ; The elastic energy density is decomposed into tensile and compressive components, expressed as: in, and Let represent the elastic energy density under tension and compression, respectively, and their expressions are: in, These are the first principal elastic strain, the second principal elastic strain, and the third principal elastic strain, respectively. and Let be the Lamé constant, and its calculation formulas are as follows: In the formula, E and These are Young's modulus and Poisson's ratio, respectively. Assuming that the compressive portion of the elastic energy density does not affect crack propagation, the elastic energy can be rewritten as follows: in, It is a degenerate function, in d Monotonically decreasing over [1, 0], this function describes the degradation of elastic energy due to the evolution of material properties, and satisfies the following conditions: Select The specific form, namely ,in k= 10 -9 It is a stability coefficient used to prevent when d Numerical singularity occurs when =0; The fracture energy during crack propagation can be approximately rewritten in the following form: Substituting the rewritten expressions for elastic energy and fracture energy into the initial expression of the energy functional of the phase-field model, we obtain the final expression of the energy functional of the phase-field model, which is expressed as:

[0012] In some embodiments, after obtaining the final expression of the energy functional of the phase-field model, it is subjected to variational processing to make the energy functional Minimize, and the total variational form is 0, i.e. , These correspond to three types of coupled control equations: phase field crack evolution equation, displacement equilibrium equation, and heat conduction equation.

[0013] In some embodiments, the derivation process of the phase field crack evolution equation includes the following steps: Applying energy functionals to phase field variables d variational form Setting it to 0 means: in, for The derivative function, This is a test function; Applying the Gaussian divergence theorem to the third term, it transforms into the following form: in, For the Laplace operator of the crack phase field, Since it is an arbitrary test function, the original governing equations of the phase field are obtained: Introducing a historical field As an irreversible condition in the calculation process, it is used to ensure that the phase field value increases monotonically. The history field is defined as follows: use The original governing equations replacing the phase field and will Substituting the specific form into the equation, we obtain the strong form of the phase field, which is the phase field crack evolution equation:

[0014] In some embodiments, the derivation of the displacement equilibrium equation includes the following steps: Applying the energy functional to the displacement variable variational form Setting it to 0 means: The displacement equilibrium equation is derived as follows: Among them, stress The specific form is as follows: in, It is a fourth-order elastic tensor matrix, which deteriorates with increasing phase field value, and is obtained by differentiating the stress tensor with respect to the strain tensor: in, It is a heaviside function, when hour, ,when hour, , It is a fourth-order tensor. , and It's the Kronecker symbol. and It is a fourth-order tensor, and its components are shown below: in, Representing vectors The Each component.

[0015] In some embodiments, the derivation of the heat conduction equation includes the following steps: Based on the initial expression of the energy functional of this phase-field model, the governing equation of the temperature field, namely the heat conduction equation, is derived as follows: Considering thermal convection between the boundary temperature and the ambient temperature, add boundary conditions, assuming the thermal conductivity is... The convection coefficient is Boundary temperature The external temperature is The following boundary conditions are obtained: The constitutive relation of matter under an electromagnetic field is described using Maxwell's equations, as shown below: The time-domain harmonic electromagnetic field problem is solved using the frequency domain method, and Maxwell's equations are simplified to Helmholtz vector equations: in, denoted as the relative magnetic permeability of the mineral. For free space wavenumber, is the relative permittivity of the mineral. For electrical conductivity, Angular frequency, It is the vacuum permittivity; Write the relative permittivity in complex form: in, The dielectric constant of the material, The dielectric loss coefficient of the material; Microwave energy is lost when propagating through a lossy dielectric. The lost microwave energy is absorbed by the dielectric material as a heat source and converted into heat energy, which is then dissipated as power. The microwave energy loss is calculated using the following formula:

[0016] In some embodiments, the process of irradiating the target reservoir with microwave parameters and simulating the crack evolution during hydraulic fracturing using a phase field model includes: calculating the fluid flow in the thermoelastic porous medium; The calculation of fluid flow in the thermoelastic porous medium includes the following steps: The thermoelastic porous medium is divided into three different regions based on the phase field value: when When it is a storage area, This is a transitional zone, when The time was a fault zone. and Thresholds between different regions; Linear interpolation of the reservoir and fracture domains is used to calculate the hydraulic parameters in the transition domain. The interpolation function is defined as follows: The effective permeability of the thermoelastic porous medium as a function of phase field value was calculated using an interpolation function. and effective biot coefficient The calculation formula is: in, and These represent the permeability and biot coefficient of the fracture, respectively. , These represent the reservoir's permeability and biot coefficient, respectively. Assume that the hydraulic fracturing fluid obeys Darcy's law and satisfies the following fluid mass conservation equation: in, For fluid density, The water storage coefficient, For flow rate, For the source term of the flow, the formula for calculating the flow velocity is as follows: in, The viscosity of the fluid; Water storage coefficient The specific form is as follows: in, Porosity is the porosity of the porous medium. The compressibility of the fluid. Bulk modulus; Assuming the initial fluid pressure is The fluid flux at the boundary is The boundary conditions satisfying the Dirichlet and Neumann conditions are given as follows:

[0017] In some embodiments, after calculating the fluid flow in the thermoelastic porous medium, the method further includes: The phase-field model is discretized into the microwave electromagnetic field generation process and the thermo-hydraulic-mechanical coupling process. The process of microwave generating electromagnetic fields is treated as a frequency domain problem, and the thermal-water-mechanical coupling process is treated as a transient problem. The electromagnetic field generation process of microwaves is solved, and electromagnetic loss is transferred as a heat source to the transient problem for solution, realizing the calculation of electromagnetic-thermal coupling and obtaining the phase field value.

[0018] The beneficial effects of this invention are as follows: Taking the much-discussed microwave-assisted unconventional natural gas extraction as its research background and the problem of excessively high wellhead pressure faced during conventional fracturing of deep unconventional oil and gas tight reservoirs as its starting point, this invention innovatively focuses on the assisting efficiency of microwave heating for hydraulic fracturing from the perspective of rock fragmentation. Utilizing the advantages of the phase-field method in simulating crack initiation and propagation, a quasi-static hydraulic fracturing phase-field model considering thermoelastic porous media under microwave-assisted fracturing conditions is established. The target reservoir is irradiated with microwave parameters, and the phase-field model is used to efficiently simulate the crack evolution process during hydraulic fracturing. Therefore, the difficulty of initiating hydraulic fracturing can be reduced by using microwave heating assistance. Attached Figure Description

[0019] Figure 1 This is a flowchart of the phase field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions in Embodiment 1 of the present invention. Figure 2 This is a schematic diagram of the sharp crack regularization treatment in Embodiment 1 of the present invention; Figure 3 This is a schematic diagram of the microwave-heated rock test equipment and numerical simulation geometric model in Embodiment 3 of the present invention; Figure 4 This is a schematic diagram comparing the heating curve of microwave-heated rock with the surface temperature field in Embodiment 3 of the present invention; Figure 5 This is a schematic diagram of the geometric model and boundary conditions for the numerical simulation experiment of ceramic quenching in Embodiment 3 of the present invention; Figure 6 This is a schematic diagram comparing the crack propagation morphology when comparing the simulation results and experimental results of ceramic quenching phase field method in Embodiment 3 of the present invention. Figure 7 This is a schematic diagram comparing the crack length when comparing the simulation results and experimental results of ceramic quenching phase field method in Embodiment 3 of the present invention; Figure 8 This is a schematic diagram of the phase field distribution in the hydraulic fracturing KGD problem in Embodiment 3 of the present invention; Figure 9 This is a schematic diagram comparing the numerical solution of fracture opening under hydraulic fracturing with the analytical solution of the phase field method and the classical analytical solution of fracture opening, fracture length, and fracture pressure in the KGD problem in Embodiment 3 of the present invention. Figure 10 This is a schematic diagram comparing the numerical solution of fracture length under hydraulic fracturing with the analytical solution when comparing the phase field method and the classical analytical solution of fracture opening, fracture length and fracture pressure in the KGD problem in Embodiment 3 of the present invention. Figure 11 This diagram illustrates the comparison between the numerical solution and the analytical solution of fracture pressure under hydraulic fracturing when comparing the phase field method for fracture opening, fracture length, and fracture pressure in the KGD problem in Embodiment 3 of the present invention. Detailed Implementation

[0020] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. The components of the embodiments of the present invention described and shown in the accompanying drawings can generally be arranged and designed in various different configurations.

[0021] Example 1

[0022] This embodiment provides a phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions. The flowchart is shown below. Figure 1 The method may include the following steps: S1. Obtain the target reservoir parameters and the fluid parameters in the target reservoir; S2. Set microwave parameters based on target reservoir parameters and fluid parameters in the target reservoir, and establish a phase field model for quasi-static hydraulic fracturing of thermoelastic porous media under microwave-assisted fracturing conditions. S3. Irradiate the target reservoir with microwave parameters and use a phase field model to simulate the crack evolution process during hydraulic fracturing.

[0023] In practical applications, generally speaking, the phase-field method can be used to detect sharp cracks. Cracks dispersed into finite width bands ,and And introduce the crack phase field Describe its state. Indicates materials that are intact and undamaged. And it indicates a complete break.

[0024] According to image segmentation theory, the crack area can be expressed as the following integral formula (1): In equation (1), the first integral represents the area of ​​the sharp crack. The second integral represents the regularized crack area. ,in The crack area density is the crack phase field. and its gradient The function is shown in the specific form of equation (2): In the formula, This represents the intrinsic length scale parameter, therefore it can be proven that: when the crack scale... When the crack approaches zero, equation (1) gives the crack area. Satisfying convergence means The above regularization process disperses sharp cracks into crack bands, which greatly reduces the difficulty of minimizing the energy functional when using variational methods. A schematic diagram of the sharp crack regularization process can be found here. Figure 2 .

[0025] It should be noted that the fracture phase-field method originated from the Griffith energy criterion and is a variational method based on Griffith theory. In Griffith fracture theory, it is assumed that the energy required to generate a fracture surface per unit area is equal to the critical fracture energy density. For brittle fracture, the energy functional of the phase-field model has the form shown in equation (3): In the formula, Represents the total energy in the system. For the elastic energy of the system, For fracture energy, For external potential energy, the energy contribution of fluid pressure in porous media also needs to be considered for hydraulic fracturing. Therefore, a new energy functional considering the influence of pore pressure can be given (4): In the formula, For elastic energy density, For biot coefficient, Pore ​​pressure, For displacement, The critical energy release rate. It is the traction force acting on the Newman boundary. It is a volume-uniformly distributed force.

[0026] Since the problem in this embodiment is the assisting effect of microwave heating on hydraulic fracturing, it is also necessary to incorporate the energy input of the system by microwave radiation into the expression of the energy functional. In practical applications, microwaves are electromagnetic waves with a frequency range of 300 MHz to 300 GHz. The principle of microwave heating is that under the action of a high-frequency alternating electromagnetic field, the internal molecules of the material are polarized and their direction of motion changes continuously with the change of the positive and negative poles of the electric field, thereby generating heat. The ability of the material to absorb microwaves depends on its dielectric properties. The principle of microwave-assisted rock fracturing is that microwave heating causes unevenly distributed thermal stress inside the material, which leads to cracks in the rock. This is essentially a thermoelastic mechanical problem. Therefore, the energy equation provided by the heat source can be given first, as shown in equation (5): Therefore, based on the above formulas (1)-(5), this embodiment establishes a quasi-static hydraulic fracturing phase field model for thermoelastic porous media under microwave-assisted fracturing conditions, and uses an energy functional to represent the phase field model, and establishes the initial expression of the energy functional of the phase field model, that is: a new energy functional is established for the phase field model of hydraulic fracturing under microwave assistance, as shown in formula (6): in, This represents the total energy of the system in the phase-field model. For elastic strain, For displacement, For phase field variables, Relative temperature For the computational domain, For fracture domain, For elastic energy density, For biot coefficient, Pore ​​pressure, For the density of the medium, For specific heat capacity, For heat flux, The critical energy release rate. The traction force acting on the Newman boundary, For the computational domain boundary, For volume-uniformly distributed force, For the computational domain. In formula (6), the first and third terms together reflect the influence of the heat source on the energy in the system. In the problem studied in this embodiment, the heat source is the heat generated by microwave electromagnetic loss. The specific control equation will be explained later.

[0027] Since this embodiment discusses thermoelastic porous media, after establishing the initial expression of the energy functional of the phase-field model, the total strain of the thermoelastic porous media is determined. The total strain is the sum of the elastic strain and the thermal strain, where the thermal strain is expressed as... Indicates total strain. The expressions are shown in formulas (7) and (8) respectively: in, The coefficient of thermal expansion is... As a second-order unit tensor, after determining the total strain of the thermoelastic porous medium, the final expression of the energy functional of this phase field model is established.

[0028] Specifically, in this embodiment, establishing the final expression of the energy functional of the phase-field model may include the following steps: In this embodiment, it can be assumed that crack propagation is driven only by tension. Therefore, the elastic energy needs to be decomposed into two parts: tension and compression, to ensure that cracks are generated only under tension. Following the strain spectrum decomposition method, the elastic strain tensor can be decomposed into the following form, as shown in equation (9): in, These are the tensile and compressive elastic strain tensors, respectively. It is the principal elastic strain. It is the direction of the principal elastic strain, and the operation symbol. Defined as ; The elastic energy density is decomposed into tensile and compressive components, as shown in equation (10): In the formula and Let represent the elastic energy density under tension and compression, respectively, as shown in equation (11): in, These are the first principal elastic strain, the second principal elastic strain, and the third principal elastic strain, respectively. and Lamé constant can be calculated using equation (12): In the formula, E and These are Young's modulus and Poisson's ratio, respectively. Assuming that the compressive portion of the elastic energy density does not affect crack propagation, the elastic energy is rewritten as shown in equation (13): In the formula It is a degenerate function, in d Monotonically decreasing over [1, 0], this function describes the degradation of elastic energy due to the evolution of material properties, and satisfies the following condition: for The specific form, namely ,in k= 10 -9 It is a stability coefficient used to prevent when d Numerical singularity occurs when =0.

[0029] To minimize the energy functional using the variational method, the fracture energy also needs to be regularized. Based on the crack regularization mentioned above in this embodiment, assuming that the intrinsic length scale is much larger than the pore size of the porous medium, the fracture energy can be approximately rewritten in the following form, as shown in Equation (15): Substituting the elastic energy (13) and fracture energy (15) into equation (6), the energy functional can be rewritten as equation (16), which gives the final expression of the energy functional of the phase-field model: In order to make the energy functional Minimization requires variational processing, and the total variational form is 0, i.e. , These correspond to three types of coupled control equations: phase field crack evolution equation, displacement equilibrium equation, and heat conduction equation.

[0030] Specifically, in this embodiment, the derivation process of the phase field crack evolution equation includes the following steps: Applying energy functionals to phase field variables d variational form Set to 0, as shown in equation (17): in, for The derivative function, To test the function, applying the Gaussian divergence theorem to the third term, equation (17) can be easily transformed into the form shown in equation (18): in, For the Laplace operator of the crack phase field, furthermore, because It is any test function, and it must have: The above set of equations constitutes the original governing equations of the phase field. Since crack healing is not considered, an irreversible condition is required during the calculation to ensure that the phase field value is monotonically increasing. Therefore, a history field needs to be introduced. The specific definition is shown in equation (20): use Replace (19-1) and will Substituting the specific form into the equation, we obtain the strong form of the phase field, which is the phase field crack evolution equation, as shown in equation (21): Specifically, the derivation of the displacement equilibrium equation includes the following steps: Applying the energy functional to the displacement variable variational form Set to 0, as shown in equation (22): The control equation for displacement can be derived from equation (22), as shown in equation (23): Stress in the formula The specific form is shown in equation (24): in, It is a fourth-order elastic tensor matrix, which deteriorates with increasing phase field value and can be obtained by differentiating the stress tensor with respect to the strain tensor. in, It is a heaviside function, when hour, ,when hour, , It is a fourth-order tensor. , and It's the Kronecker symbol. and It is a fourth-order tensor, and its specific component forms are shown in equations (26) and (27): in, Representing vectors The Each component.

[0031] Similarly, the derivation of the heat conduction equation includes the following steps: Similarly, the governing equations for the temperature field can be derived from equation (6), as shown in equation (28): Considering the thermal convection between the boundary temperature and the ambient temperature, boundary conditions need to be added. Assume the thermal conductivity is... The convection coefficient is Boundary temperature The external temperature is The boundary conditions can be obtained as shown in equation (29): The microwave excitation problem can be described by Maxwell's equations, which describe the constitutive relationship of matter under electromagnetic fields, as shown in equation (30): The frequency domain method is typically used to solve time-domain harmonic electromagnetic field problems. This method can simplify Maxwell's equations into Helmholtz vector equations. in, The relative magnetic permeability of the mineral, N / A 2 , For free space wavenumber, is the relative permittivity of the mineral. The conductivity is expressed in S / m. Angular frequency, rad / s The vacuum permittivity is F / m; The relative permittivity can be written in complex form: in, The dielectric constant of the material, The dielectric loss coefficient of the material; Microwave energy is lost when propagating through a lossy dielectric. The lost microwave energy is absorbed by the dielectric material as a heat source and converted into heat energy, which can be dissipated through power loss. (W / m 3 Calculation of microwave energy loss: In practical applications, it is also necessary to calculate fluid flow during crack evolution. Specifically, the process of irradiating the target reservoir with microwave parameters and simulating crack evolution during hydraulic fracturing using a phase-field model includes calculating fluid flow in the thermoelastic porous medium. This calculation may include the following steps: To facilitate fluid flow calculations, the thermoelastic porous medium can be divided into three distinct regions based on the phase field value: when... When it is a storage area, This is a transitional zone, when The time was a fault zone. and Thresholds between different regions; In this embodiment, linear interpolation of the reservoir and fracture domain is used to calculate the hydraulic parameters in the transition domain, and the interpolation function is defined as shown in equation (34): The effective permeability of thermoelastic porous media as a function of phase field can be calculated using the interpolation function described above. and effective biot coefficient As shown in equation (35): In the formula, and These represent the permeability and biot coefficient of the fracture, respectively. , These represent the reservoir's permeability and biot coefficient, respectively. Assume that the hydraulic fracturing fluid obeys Darcy's law and satisfies the following fluid mass conservation equation: in, For fluid density, The water storage coefficient, For flow rate, For the source term of the flow, the formula for calculating the flow velocity is as follows: In the formula, The viscosity of the fluid; The specific form is shown in equation (38): in, Porosity is the porosity of the porous medium. The compressibility of the fluid. Bulk modulus; Meanwhile, assuming the initial fluid pressure is The fluid flux at the boundary is The boundary conditions satisfying the Dirichlet and Neumann conditions are given as follows: Finally, the phase field value can be calculated based on the phase field model of quasi-static hydraulic fracturing of thermoelastic porous media under microwave-assisted fracturing conditions. Therefore, in this embodiment, after calculating the fluid flow in the thermoelastic porous media, the following steps are also included: The phase-field model is discretized into the microwave electromagnetic field generation process and the thermo-hydraulic-mechanical coupling process. The process of microwave generating electromagnetic fields is treated as a frequency domain problem, and the thermal-water-mechanical coupling process is treated as a transient problem. The electromagnetic field generation process of microwaves is solved, and electromagnetic loss is transferred as a heat source to the transient problem for solution, realizing the calculation of electromagnetic-thermal coupling and obtaining the phase field value.

[0032] Example 2

[0033] Based on Example 1, this example provides a detailed explanation of the calculation process for the phase field value.

[0034] The phase-field model in this embodiment involves two processes: frequency domain and transient. The process of microwaves generating electromagnetic fields is a frequency domain problem, while the subsequent thermo-water-mechanical coupling problem is a transient problem.

[0035] First, we introduce the discretization and solution of frequency domain problems. In this embodiment, we use the weighted residual method to process the Helmholtz equation, and a complex-valued trial function is introduced here. And perform volume integration over the entire computational domain: It should be noted that, because this is a complex domain problem, therefore It is a complex conjugate. This is the conjugate transpose. By performing partial integration over the first term, and assuming the existence of excitations at the boundary, we obtain the complete weak form: In the formula, To contribute to microwave excitation, the weak form described above, using shape function interpolation, can be rewritten as follows: In the formula, For vector shape functions, , For the electric field, represents the degree of freedom.

[0036] Equation (42) is simplified to matrix form, and the residuals are constructed: The damped Newton method is used to solve the above-mentioned weak form of nonlinear electromagnetic field, and the iterative scheme is as follows: In the formula, and These are index symbols used in nonlinear iterative processes to indicate the current iteration step and the next step in the damped Newton's method. For Jacobian matrices, The damping coefficient is used for adaptive contraction / expansion step size to improve computational convergence.

[0037] After obtaining the solution for the electromagnetic field, the electromagnetic loss is transferred as a heat source to the transient problem for solution, thus realizing the electromagnetic-thermal coupling calculation.

[0038] For transient problems, the nonlinear partial differential governing equations were discretized using the finite element method, including the stress balance equation, fluid continuity equation, phase field equation, and heat transfer equation.

[0039] Multiply the left side of equation (23-1) by a weight function For the entire region Integrating yields the weighted residual form of the equilibrium differential equation, as shown in equation (45): Using the divergence theorem, we can solve the above equation... Perform partial integration and substitute the boundary conditions. Finally, substituting equation (24) into the equation, we obtain the complete expression in the weak form with boundary terms: Similar to the above operation, the phase field governing equation Eq.(21), the temperature field governing equation Eq.(28), and the fluid continuity equation Eq.(36) are multiplied by weight functions respectively. , and By integrating over the computational domain and then applying the divergence theorem in conjunction with boundary conditions, we obtain the equivalent weak integral forms of each governing equation: Within the finite element framework, having The displacement, phase field, temperature, pressure, and weighting function in each node element are interpolated by shape functions, as shown in equations (50) and (51): In the formula, For the spatial location within the unit, These are vectors representing nodal displacement, phase field, temperature, and pressure, respectively. and These are the interpolation shape functions for displacement, phase field, temperature, and pressure, respectively. , .

[0040] The gradients of displacement, phase field, temperature, pressure, and weighting function are shown in equations (52) and (53): In the formula, and These are the gradients of the interpolation shape functions for displacement, phase field, temperature, and pressure, respectively. Where: , .

[0041] Combining equations (46)-(49) with equations (50)-(53), we construct the residuals: Since electromagnetic-thermal-hydraulic-mechanical-phase-field coupling is a complex problem that exhibits high nonlinearity in both space and time, an interleaved solution scheme is adopted to solve the problem in order to reduce computational costs.

[0042] Specifically, the transient problem is divided into two parts, and the classic Newton's method is used to solve it. Assuming that the phase field value of each element node is fixed in the s-th iteration step, the iterative format of the electromagnetic-thermal-hydraulic-mechanical coupling equations is as follows: In the formula For Jacobian matrices, n and m Superscript symbol, number In each iteration step, the displacement, temperature, and pressure are respectively... Then, the displacement, temperature, and pressure in this iteration step are fixed, and the solution is obtained for the 1st iteration. The phase field values ​​in each iteration step are given by the following iteration format: In the formula It is a Jacobian matrix.

[0043] Example 3

[0044] Based on Examples 1 and 2, since there are currently no experimental reports on microwave-assisted hydraulic fracturing, it is impossible to directly compare numerical simulation results with experimental results to verify the effectiveness of the simulation method. Therefore, to verify the rationality of the numerical simulation method, this example adopts a step-by-step verification approach, verifying it from three aspects: electromagnetic-thermal coupling, thermo-mechanical coupling, and hydraulic-mechanical coupling. The effectiveness of electromagnetic-thermal coupling and thermo-mechanical coupling is verified by comparing indoor experimental results with numerical simulation results, while the effectiveness of hydraulic fracturing simulation is verified by comparing analytical solutions with numerical solutions.

[0045] I. Numerical simulation verification of microwave-heated rocks: like Figure 3As shown in the figure, the indoor experimental equipment and numerical simulation geometric model for microwave-heated rock are illustrated. The heating cavity is a 400mm × 400mm × 400mm all-copper cubic cavity. Two BJ-26 rectangular waveguides are positioned at the horizontal and vertical ends of the heating cavity, respectively. The dielectric material inside the heating cavity and waveguides is air. The microwave source excites the heating cavity through the rectangular waveguides at a frequency of 2.45 GHz, operating in TE10 mode. The sample-bearing base is a boron nitride cylinder with a diameter of 80mm and a height of 150mm. The sample, a standard 50mm × 100mm cylindrical specimen, is positioned above the boron nitride base. A pre-fabricated crack is created at the center of the sample using Boolean operations; the crack has a diameter of 5mm, a length of 20mm, and an opening of 2mm.

[0046] Experimental and numerical simulation results of pre-fractured sandstone after 270s of microwave irradiation with a power of 3000W are as follows: Figure 4 As shown, Figure 4 The image shows the temperature rise curve of the highest temperature on the rock surface and the temperature field distribution on the rock surface. Figure 4 As can be seen, the heating curves obtained from the experiment and numerical simulation show good matching. Due to the uneven distribution of electric field intensity in space, obvious temperature zoning phenomenon is formed on the sample surface after microwave radiation. The simulation results show that the positions of the high temperature zone and low temperature zone on the sample surface correspond well with the experimental results, and the simulation results reproduce the experimental results well.

[0047] II. Numerical simulation verification of ceramic quenching test: Quenching tests on ceramic materials have been widely used by scholars to verify the reliability of phase-field thermo-mechanical coupling models. This embodiment also verifies the thermo-mechanical coupling reliability of the model through a numerical simulation case of a ceramic quenching test. Indoor test results were selected for comparison with the simulation results of this model. A ceramic plate with dimensions of 1.0 mm wide, 9.8 mm thick, and 50 mm long was heated to a specified temperature and then cooled in a water bath. After cooling, the sample was impregnated with blue dye to observe the cracks formed during thermal shock. This embodiment selects the case of ΔT = 380 K for simulation, such as... Figure 5 As shown, the geometric dimensions of the ceramic plate in the numerical simulation are the same as those in the experiment, and other physical and mechanical parameters are shown in Table 1. In this embodiment, after introducing phase field variables to describe the crack, crack initiation is related not only to the model geometry and thermal load, but also to the parameters. Regarding this, to ensure a sufficient number of microcracks initiate, take Grid size set to The mesh type is a free quadrilateral mesh.

[0048] Table 1 Physical and mechanical parameters of ceramic plates Figure 6 and Figure 7 The comparison between the numerical simulation results and the final experimental crack morphology is shown. Figure 6 It can be seen that the simulated crack distribution morphology is highly consistent with the experimental results, demonstrating excellent reproduction capabilities for both long and short cracks. Figure 7 As shown, by combining the experimental and numerical simulation results, it can be observed that the crack lengths are almost identical. Furthermore, cracks exceeding one-quarter of the plate width are defined as long cracks; the experimental results show 24 long cracks, while the numerical simulation results show 25. The numerical simulation results are close to the experimental results both qualitatively and semi-quantitatively, thus demonstrating the good reliability of the phase-field model in describing thermo-mechanical coupling problems in this embodiment.

[0049] III. KGD problem in hydraulic crack propagation: Hydraulic fracture propagation in infinitely porous elastic media can be used to verify the reliability of the phase-field method hydraulic fracturing model. The hydraulic fracture propagation problem in infinitely porous elastic media can be described by the KGD (Khristianovic–Geertsma–deKlerk) model. This model is one of the classic hydraulic fracturing theoretical models, used to describe the propagation behavior of two-dimensional planar fractures. This embodiment also verifies the reliability of the phase-field model by comparing a numerical solution and an analytical solution for hydraulic fracture propagation. The geometric model is as follows: Figure 8 As shown, due to the significant difference between the geometric models of the cracks and the medium, the porous elastic medium can be considered as an infinitely porous elastic medium. The parameters of the numerical model are shown in Table 2.

[0050] Table 2 Physical and mechanical parameters of KGD Since there are analytical solutions for crack orifice displacement (CMOD), crack length (CL), and crack orifice pressure (CMP), as shown in equations (60)-(62): In the formula, A , B , C For fitting parameters, A =1.87, B =0.68, C= 1.38, A =2.14, B =0.65, C= 1.97, q Normal velocity, Let Lamé constant be . For fluid viscosity, Poisson's ratio, S It is the in-situ stress perpendicular to the crack propagation direction. t It is the fluid injection time.

[0051] Furthermore, a formula for calculating the width of hydraulic fracturing fractures based on strain is proposed, as shown in equation (63). This formula can effectively calculate the width of fractures in the phase field model: In the formula, The element size at the fracture region, For the body to adapt to strain.

[0052] A comparison of the numerical and classical analytical solutions for fracture opening, fracture length, and fracture pressure using the phase-field method is shown below. Figure 9 , Figure 10 and Figure 11 As shown. From Figure 9-11 As can be seen, COMD and CL increase with time, and the rate of increase gradually decreases. The numerical solution curve of COMD falls within the envelope of the two analytical solutions, and the numerical solution of CL is very close to both analytical solutions. CMP first increases and then decreases with time, and the curve of the numerical solution also falls within the envelope of the two analytical solutions. Therefore, the phase-field model used in this embodiment also has good reliability in hydraulic fracturing simulation.

[0053] Therefore, this embodiment verifies the feasibility of the model from three aspects: electromagnetic-thermal, thermo-mechanical, and hydro-mechanical coupling. It studies the auxiliary effect of microwave heating on hydraulic fracturing under different characteristic parameters, reveals the sensitivity of unconventional energy reservoirs with different mechanical properties to microwave-assisted hydraulic fracturing, and explores the effectiveness of microwave heating in reducing rock formation fracture pressure under different confining pressures. Therefore, this embodiment aims to provide evidence for the effectiveness of this technology from the perspective of reducing the difficulty of hydraulic fracturing initiation by microwave heating, which is different from the previous research on microwave heating to help oil and gas migration. This can help microwave-assisted unconventional natural gas extraction.

[0054] The above are merely preferred embodiments of the present invention and are not intended to limit the present invention. Various modifications and variations can be made to the present invention by those skilled in the art. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions, characterized in that, Includes the following steps: Obtain the target reservoir parameters and the fluid parameters in the target reservoir; Microwave parameters are set based on the target reservoir parameters and fluid parameters in the target reservoir, and a quasi-static hydraulic fracturing phase field model considering microwave-assisted fracturing of thermoelastic porous media is established. The target reservoir was irradiated with microwave parameters, and the crack evolution process during hydraulic fracturing was simulated using a phase field model.

2. The phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions according to claim 1, characterized in that, Considering the energy of fluid pressure in porous media and the heat generated by microwave electromagnetic losses, a quasi-static hydraulic fracturing phase-field model for thermoelastic porous media under microwave-assisted fracturing conditions is established. This phase-field model is represented by an energy functional, and an initial expression for the energy functional of this phase-field model is established, which is expressed as follows: in, This represents the total energy of the system in the phase-field model. For elastic strain, For displacement, For phase field variables, Relative temperature For the computational domain, For fracture domain, For elastic energy density, For biot coefficient, Pore ​​pressure, For the density of the medium, For specific heat capacity, For heat flux, The critical energy release rate. The traction force acting on the Newman boundary, For the computational domain boundary, For volume-uniformly distributed force, For computational domains.

3. The phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions according to claim 2, characterized in that, After establishing the initial expression for the energy functional of the phase-field model, the total strain of the thermoelastic porous medium is determined. The total strain is the sum of the elastic strain and the thermal strain, where the thermal strain is expressed as... Indicates total strain. express; The total strain The expression is: thermal strain The expression is: ,in, The coefficient of thermal expansion is... It is a second-order unit tensor; After determining the total strain of the thermoelastic porous medium, the final expression of the energy functional of the phase field model is established.

4. The phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions according to claim 3, characterized in that, The final expression for the energy functional of the phase-field model is established by the following steps: Assuming crack propagation is driven solely by tension, and the elastic strain is decomposed into tensile and compressive strain tensors such that cracks only occur under tension, the decomposition of elastic strain into tensile and compressive strain tensors is expressed as: in, These are the tensile and compressive elastic strain tensors, respectively. It is the principal elastic strain. It is the direction of the principal elastic strain, and the operation symbol. Defined as ; The elastic energy density is decomposed into tensile and compressive components, expressed as: in, and Let represent the elastic energy density under tension and compression, respectively, and their expressions are: in, These are the first principal elastic strain, the second principal elastic strain, and the third principal elastic strain, respectively. and Let be the Lamé constant, and its calculation formulas are as follows: In the formula, E and These are Young's modulus and Poisson's ratio, respectively. Assuming that the compressive portion of the elastic energy density does not affect crack propagation, the elastic energy can be rewritten as follows: in, It is a degenerate function, in d Monotonically decreasing over [1, 0], this function describes the degradation of elastic energy due to the evolution of material properties, and satisfies the following conditions: Select The specific form, namely ,in k= 10 -9 It is a stability coefficient used to prevent when d Numerical singularity occurs when =0; The fracture energy during crack propagation can be approximately rewritten in the following form: Substituting the rewritten expressions for elastic energy and fracture energy into the initial expression of the energy functional of the phase-field model, we obtain the final expression of the energy functional of the phase-field model, which is expressed as:

5. The phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions according to claim 3, characterized in that, After obtaining the final expression of the energy functional of the phase-field model, a variational process is applied to it to make the energy functional... Minimize, and the total variational form is 0, i.e. , These correspond to three types of coupled control equations: phase field crack evolution equation, displacement equilibrium equation, and heat conduction equation.

6. The phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions according to claim 5, characterized in that, The derivation process of the phase field crack evolution equation includes the following steps: Applying energy functionals to phase field variables d variational form Setting it to 0 means: in, for The derivative function, This is a test function; Applying the Gaussian divergence theorem to the third term, it transforms into the following form: in, For the Laplace operator of the crack phase field, Since it is an arbitrary test function, the original governing equations of the phase field are obtained: Introducing a historical field As an irreversible condition in the calculation process, it is used to ensure that the phase field value increases monotonically. The history field is defined as follows: use The original governing equations that replace the phase field and will Substituting the specific form into the equation, we obtain the strong form of the phase field, which is the phase field crack evolution equation:

7. The phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions according to claim 5, characterized in that, The derivation of the displacement equilibrium equation includes the following steps: Applying the energy functional to the displacement variable variational form Setting it to 0 means: The displacement equilibrium equation is derived as follows: Among them, stress The specific form is as follows: in, It is a fourth-order elastic tensor matrix, which deteriorates with increasing phase field value, and is obtained by differentiating the stress tensor with respect to the strain tensor: in, It is a heaviside function, when hour, ,when hour, , It is a fourth-order tensor. , and It's the Kronecker symbol. and It is a fourth-order tensor, and its components are shown below: in, Representing vectors The Each component.

8. The phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions according to claim 5, characterized in that, The derivation of the heat conduction equation includes the following steps: Based on the initial expression of the energy functional of this phase-field model, the governing equation of the temperature field, namely the heat conduction equation, is derived as follows: Considering thermal convection between the boundary temperature and the ambient temperature, add boundary conditions, assuming the thermal conductivity is... The convection coefficient is Boundary temperature The external temperature is The following boundary conditions are obtained: The constitutive relation of matter under an electromagnetic field is described using Maxwell's equations, as shown below: The time-domain harmonic electromagnetic field problem is solved using the frequency domain method, and Maxwell's equations are simplified to Helmholtz vector equations: in, denoted as the relative magnetic permeability of the mineral. For free space wavenumber, is the relative permittivity of the mineral. For electrical conductivity, Angular frequency, It is the vacuum permittivity; Write the relative permittivity in complex form: in, The dielectric constant of the material, The dielectric loss coefficient of the material; Microwave energy is lost when propagating through a lossy dielectric. The lost microwave energy is absorbed by the dielectric material as a heat source and converted into heat energy, which is then dissipated as power. The microwave energy loss is calculated using the following formula:

9. The phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions according to any one of claims 1-8, characterized in that, The process of irradiating the target reservoir with microwave parameters and simulating the crack evolution during hydraulic fracturing using a phase field model includes: calculating the fluid flow in the thermoelastic porous medium. The calculation of fluid flow in the thermoelastic porous medium includes the following steps: The thermoelastic porous medium is divided into three different regions based on the phase field value: when When it is a storage area, This is a transitional zone, when The time was a fault zone. and The threshold between different regions; Linear interpolation of the reservoir and fracture domains is used to calculate the hydraulic parameters in the transition domain. The interpolation function is defined as follows: The effective permeability of the thermoelastic porous medium as a function of phase field value was calculated using an interpolation function. and effective biot coefficient The calculation formula is: in, and These represent the permeability and biot coefficient of the fracture, respectively. , These represent the reservoir's permeability and biot coefficient, respectively. Assume that the hydraulic fracturing fluid obeys Darcy's law and satisfies the following fluid mass conservation equation: in, For fluid density, The water storage coefficient, For flow rate, For the source term of the flow, the formula for calculating the flow velocity is as follows: in, The viscosity of the fluid; Water storage coefficient The specific form is as follows: in, Porosity is the porosity of the porous medium. The compressibility of the fluid. Bulk modulus; Assuming the initial fluid pressure is The fluid flux at the boundary is The boundary conditions satisfying the Dirichlet and Neumann conditions are given as follows:

10. The phase-field simulation method for quasi-static hydraulic fracturing of thermoelastic porous media under high-energy fracturing conditions according to claim 9, characterized in that, After calculating the fluid flow in thermoelastic porous media, the following is also included: The phase-field model is discretized into the microwave electromagnetic field generation process and the thermo-hydraulic-mechanical coupling process. The process of microwave generating electromagnetic fields is treated as a frequency domain problem, and the thermal-water-mechanical coupling process is treated as a transient problem. The electromagnetic field generation process of microwaves is solved, and electromagnetic loss is transferred as a heat source to the transient problem for solution, realizing the calculation of electromagnetic-thermal coupling and obtaining the phase field value.

Images