Carbonate geothermal reservoir acid fracturing evaluation method based on coupling model

Abstract

The application discloses a carbonate rock geothermal reservoir acid fracturing evaluation method based on a heat-water-force-chemical coupling model, and specifically relates to the following steps: collecting physical and mechanical parameters of a carbonate rock reservoir matrix and engineering parameters of a fracture fissure; according to the engineering parameters, a geometric model for simulating the carbonate rock reservoir is established, variables in the fractures in each field obtained are distributed to the fracture domain, and a calculation domain is obtained; a heat transfer field, a hydraulic field, a stress field and a chemical field are added to the calculation domain respectively, a coupling relationship among them is established, and a heat-water-force-chemical coupling model covering the deformation of the fractures is formed; and the separated coupling model is solved to perform acid fracturing evaluation. The application designs a full-coupling THMC model for the acid fracturing process of the carbonate rock reservoir, the model integrates the interaction among heat transfer, acidizing fluid flow, deformation of the rock matrix and the fractures, migration of reaction dilution material and chemical corrosion phenomena, and is used for quantifying the influence of the acid fracturing of the carbonate rock geothermal reservoir.

Description

Technical Field

[0001] This invention belongs to the field of geothermal extraction technology, specifically relating to an evaluation method for acid fracturing of carbonate geothermal reservoirs based on a coupled model. Background Technology

[0002] Among hydrothermal reservoir types, carbonate geothermal reservoirs are important, and their large-scale development and utilization has become a key focus of future geothermal exploration and development. Currently, the main fracturing techniques employed include hydraulic fracturing, acid fracturing, proppant fracturing, and high-energy gas fracturing. Of these techniques, acid fracturing is of great significance for enhancing carbonate geothermal systems. Acid fracturing involves injecting acid or acidic pre-fluid into the casing at a rate exceeding the reservoir capacity, thereby rapidly establishing wellbore pressure exceeding the formation compressive stress and rock tensile strength. This process causes formation fracturing and the formation of fractures. Continuous injection of acid can extend these fractures, forming acid-etched fractures, which enhance the conductivity of the formation compared to the original formation. This not only increases the production capacity of geothermal wells but also enhances their backfilling capacity.

[0003] Currently, acid fracturing technology has received widespread attention. Ma Jianjun and his research team proposed a thermo-hydraulic-chemical (THC) coupled model, using a unified pipeline network method to simulate the acid fracturing process, which neglects the interaction between stress and THC coupling. Researchers such as Mahmoodi have shown through parameter sensitivity analysis that the acid fracturing effect has a significant impact on the formation and expansion of boreholes. Existing technologies also establish coupled THMC models to appropriately consider the influence of heat generation during rock acid fracturing on the reaction rate, thereby affecting the chemical process. Furthermore, the THMC model has been applied to study the effects of chemical corrosion and precipitation, injection concentration, and injection rate on the evolution of fracture aperture. These studies demonstrate that exploring the THMC coupling mechanism is one of the key factors in optimizing geothermal resource development.

[0004] Discrete fracture networks (DFNs) are crucial for understanding geothermal reservoirs because they can identify geothermal fluid flows. Several researchers have explored the applications of DFNs in fluid flow. Enhanced fracture network connectivity obtained through DFNs can significantly improve solute concentration distribution. A thermo-hydraulic-mechanical (THM) coupled model based on embedded DFNs has also been proposed to simulate complex flow behavior and thermal properties during CO2 injection, thus validating the importance of DFNs in analyzing heat flow, fracture deformation, and flow characteristics. In studies of geothermal reservoirs containing natural fractures, DFNs are the preferred method for simulating the behavior and distribution of natural fractures. DFN models accurately capture the behavior and distribution of natural fractures, suggesting that they significantly influence stress changes and damage evolution. In carbonate reservoirs, the application of DFNs allows researchers to study the presence, pore size, distribution, and sensitivity of dissolution processes to acid injection rates, revealing the significant impact of fracture structure on reactive flow. However, the role of DFN characteristics in the acid fracturing process of carbonate geothermal rocks in THMC coupled models has rarely been comprehensively studied. Therefore, it is of great significance to propose a coupled THMC model that considers DFN characteristics and to study the chemical corrosion and mechanical deformation of fractures during acid fracturing in carbonate geothermal reservoirs. Summary of the Invention

[0005] The purpose of this invention is to provide an evaluation method for acid fracturing of carbonate geothermal reservoirs based on a coupled model. This model integrates the interactions between heat transfer, acidizing fluid flow, deformation of rock matrix and fractures, migration of reactive diluents, and chemical corrosion phenomena, and is used to quantify the impact of acid fracturing on carbonate geothermal reservoirs.

[0006] The technical solution adopted in this invention is an evaluation method for acid fracturing of carbonate geothermal reservoirs based on a coupled model, which is implemented according to the following steps:

[0007] Step 1: Collect the physical and mechanical parameters of the carbonate reservoir matrix and the engineering parameters of the fractures and fissures;

[0008] Step 2: Based on the engineering parameters, establish a geometric model of the simulated carbonate reservoir, and assign the variables in the fractures in each field to the fracture domain to obtain the computational domain;

[0009] Step 3: Add heat transfer field, hydraulic field, stress field and chemical field to the computational domain respectively, establish the coupling relationship between them, and thus form a thermo-hydraulic-mechanical-chemical coupling model covering crack deformation;

[0010] Step 4: Solve the decoupled model and perform acid fracturing evaluation.

[0011] The invention is further characterized in that,

[0012] In step 1, the physical and mechanical parameters of the carbonate reservoir matrix include: Young's modulus, Poisson's ratio, density, Biot ratio, coefficient of thermal expansion, initial porosity, initial permeability, storage coefficient, thermal conductivity, Poisson's ratio, and specific heat capacity; the engineering parameters of the fractures include: fracture width, Young's fracture modulus of the fracture, storage coefficient of the fracture, and initial pore size of the fracture.

[0013] In step 2, the specific steps for establishing a geometric model of the simulated carbonate reservoir using COMSOL software are as follows: Based on the engineering geological characteristics of the carbonate reservoir, a computational geometric model is established; the geometry and size of the reservoir are defined on the canvas; the thickness of the fractures in the carbonate reservoir is determined according to the engineering parameters; and the variables in the fractures in each field are assigned to the fracture domain using functions in COMSOL to obtain the computational domain.

[0014] In step 3, the thermo-hydraulic-mechanical-chemical coupling model covering crack deformation includes: porous media flow equation, crack flow equation, transport equation of diluted species in porous media, transport equation of diluted species in crack, heat transfer equation in porous media, heat transfer equation in fracture, solid mechanics equation, thin elastic layer equation, boundary ordinary differential equation and partial differential algebraic equation, and neighborhood ordinary differential equation and partial differential algebraic equation.

[0015] The flow equation of the acidizing fluid in the porous rock matrix, i.e. the flow equation of the porous medium, is shown in equation (1):

[0016] (1);

[0017] The flow equation of the acidizing fluid in the crack, i.e. the crack flow equation, is shown in equation (2):

[0018] (2);

[0019] in, The density of the acidified fluid. m It is the storage coefficient of the rock matrix;

[0020] f The storage coefficient for cracks. It's pressure. It is time. It is the fluid velocity vector of the acidifying fluid in the porous rock matrix. It is the fluid velocity vector of the acidizing fluid in the crack. It is the Biot-Willis coefficient of the rock matrix. It is the Biot-Willis coefficient of the crack; It is the volumetric strain of the rock matrix. It is the volumetric strain of the crack. It refers to the size of the crack opening. It is a gradient operator;

[0021] It is the mass flux from the rock matrix to the fracture, as shown in equation (3):

[0022] (3);

[0023] in, m It is the permeability of the rock matrix. acid It refers to the viscosity of the acidified fluid. n u It is the normal vector pointing to the rock-crack interface on the side facing the rock. n d It is the normal vector pointing towards the crack side of the rock-crack interface;

[0024] Based on the variable of HCl, its transport in the porous matrix is ​​controlled by formula (5), which is the transport equation of dilution species in the porous medium;

[0025] (5)

[0026] in, It refers to the porosity of the rock matrix. C HCl It is the concentration of HCl. D e It is the effective diffusion coefficient, which can be obtained by... The tortuosity of the crack The rate of acid-rock chemical reactions was calculated using diffusion coefficients and given sufficient amounts of solid reactants. As shown in equation (6);

[0027] (6);

[0028] Among them, reaction order for The reaction rate constant is shown in equation (7):

[0029] (7);

[0030] in, A It is the Arrhenius constant. e It is the base of the natural logarithm. E It is activation energy. R g It is the molar gas constant. ，T It is absolute temperature.

[0031] The transport equation for dilute species in the crack is shown in equation (8):

[0032] (8);

[0033] in, This indicates the rate of acidification reaction at the fluid-solid interface. Q c The mass flux from the porous matrix into the crack is represented by Equation (9):

[0034] (9);

[0035] It can be determined based on mass transfer, as shown in equation (10):

[0036] (10);

[0037] in, It is the mass transfer coefficient, as shown in equation (11):

[0038] (11);

[0039] in, The solute diffusion coefficient in the fracture is... Let Reynolds number be 1. For the Schmitt number, ; That is, the interface area of ​​the crack per unit volume. It is the HCl concentration at the fluid-solid interface; It can be expressed as equation (12):

[0040] (12).

[0041] The heat transfer equation in porous media is shown in equation (13):

[0042] (13);

[0043] in, It is the effective volumetric specific heat capacity of the rock matrix. It is the effective thermal conductivity of the rock matrix. It is the effective density of the rock matrix. It is the specific heat capacity under constant pressure. It's temperature. It is a heat source generated by a chemical reaction, as shown in equation (14);

[0044] (14)

[0045] in, For reaction enthalpy;

[0046] The heat transfer equation during fracture is shown in equation (15).

[0047] (15);

[0048] in, The effective volumetric heat capacity of the crack. The effective thermal conductivity of the crack, The heat source for the reaction in the crack is as shown in equation (16);

[0049] (16);

[0050] It is the heat source term from the matrix to the fracture, as shown in equation (17);

[0051] (17);

[0052] The equation for solid mechanics is shown in equation (18);

[0053] (18)

[0054] in, and These are Lamé constants, It is volumetric strain, It is the displacement vector and It is the bulk modulus;

[0055] The equation for the thin elastic layer is shown in equation (19);

[0056] (19)

[0057] in, Force / unit area , The displacement vectors at the upper and lower boundaries of the crack layer are respectively... Initial compensation displacement;

[0058] Stiffness can be measured by normal stiffness and tangential stiffness They are obtained by means of two corresponding deformations, as shown in equation (20);

[0059] (20)

[0060] in, The unit normal vector of the crack; the normal and tangential stiffness of the thin elastic layer can be calculated using its material parameters, as shown in equation (21);

[0061] (twenty one);

[0062] in, It is the elastic modulus. It is Poisson's ratio. It is the shear modulus.

[0063] For cracks, the boundary ordinary differential equation and partial differential algebraic equation, i.e. the crack pore size change rate caused by chemical corrosion, are shown in equation (26).

[0064] (26);

[0065] The ordinary differential equations and partial differential algebraic equations of the field, namely the crack opening change caused by mechanical deformation, can be calculated by using the parameters of the thin elastic layer and the normal strain of the layer and the fluid pressure in the crack, as shown in equation (27).

[0066] (27);

[0067] Under the coupled thermo-hydraulic-mechanical (THMC) conditions, the crack opening is finally obtained by integrating equations (26) and (27), which is equation (28).

[0068] (28);

[0069] The permeability of the crack can be calculated using a parallel plate flow model, as shown in equation (29);

[0070] (29);

[0071] A thermo-hydraulic-mechanical-chemical coupling model covering crack deformation was established based on the above.

[0072] Step 4 specifically involves: using a decoupled solver to solve the decoupled model, performing numerical calculations and sensitivity analysis on the coupled model, and evaluating the acid fracturing effect.

[0073] The beneficial effects of this invention are:

[0074] In the development of carbonate geothermal reservoirs, acid fracturing technology plays a crucial role in improving reservoir permeability. This invention constructs a comprehensive thermo-hydraulic-mechanical-chemical coupled model to simulate the complex phenomena involved in acid fracturing, including solute migration, acid-rock reactions, heat conduction, and formation deformation. This model is specifically designed for discrete fracture networks at the field scale, comprehensively considering the chemical corrosion and mechanical deformation of fractures, as well as the chemical corrosion effects on the bedrock and fracture surfaces. Furthermore, the model can track changes in fracture aperture and matrix porosity over time. Attached Figure Description

[0075] Figure 1 This is a flowchart of the acid fracturing evaluation method for carbonate geothermal reservoirs based on a coupled model, as described in this invention.

[0076] Figure 2 This is a schematic diagram of the boundary conditions of the THMC field in this invention;

[0077] Figure 3 It is a geometric diagram of a natural crack network with different numbers and sizes of cracks;

[0078] Figure 4 This is a schematic diagram of the THMC coupling model of the present invention;

[0079] Figure 5 This is a schematic diagram of the pressure distribution of the DFN at different times in this invention;

[0080] Figure 6 This is a spatial and temporal distribution diagram of the HCl concentration of the reactant at different times on the observation line of the DFN of this invention;

[0081] Figure 7 This is a spatial and temporal distribution diagram of the CaCl2 concentration of the product at different times on the observation line of the DFN of this invention;

[0082] Figure 8 This is a spatiotemporal distribution diagram of the chemical aperture of the DFN observation line of this invention (I);

[0083] Figure 9 This is the spatiotemporal distribution map (II) of the chemical aperture of the DFN observation line of this invention. Detailed Implementation

[0084] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments.

[0085] Example 1

[0086] This invention discloses an acid fracturing evaluation method for carbonate geothermal reservoirs based on a coupled model. Specifically, the method involves: collecting the physical and mechanical parameters of the carbonate reservoir matrix and the engineering parameters of the fractures; establishing a geometric model of the carbonate reservoir using COMSOL software based on the engineering parameters; allocating the variables in the fractures within each field to the fracture domain using functions in COMSOL to obtain the computational domain; adding heat transfer, hydraulic, stress, and chemical fields to the computational domain to establish their coupling relationships, thus forming a thermo-hydraulic-mechanical-chemical coupled model covering fracture deformation; allocating initial conditions and boundaries, creating and optimizing the mesh, solving the decoupled coupled model, and performing acid fracturing evaluation.

[0087] Example 2

[0088] This invention relates to a method for evaluating acid fracturing in carbonate geothermal reservoirs based on a coupled model, such as... Figure 1 As shown, please follow these steps:

[0089] Step 1: Collect the physical and mechanical parameters of the carbonate reservoir matrix and the engineering parameters of the fractures and fissures;

[0090] The physical and mechanical parameters of carbonate reservoir matrix include: Young's modulus, Poisson's ratio, density, Biot ratio, coefficient of thermal expansion, initial porosity, initial permeability, storage coefficient, thermal conductivity, Poisson's ratio, and specific heat capacity.

[0091] The engineering parameters of fractures include: fracture width, Young's modulus of fracture, fracture storage coefficient, and initial pore size of fracture.

[0092] Step 2: Based on the engineering parameters, use COMSOL software to establish a geometric model of the simulated carbonate reservoir. Through the functions in COMSOL, the variables in the fractures of each field are assigned to the fracture domain to obtain the computational domain.

[0093] The specific steps for establishing a geometric model of a simulated carbonate reservoir using COMSOL software are as follows: First, establish a computational geometric model based on the engineering geological characteristics of the carbonate reservoir. Second, define the reservoir's geometry and dimensions on the canvas, and determine the thickness of fractures in the carbonate reservoir based on engineering parameters. Third, use COMSOL functions to assign variables from the fractures in each field to the fracture domain, thus obtaining the computational domain.

[0094] Step 3: Add heat transfer field, hydraulic field, stress field and chemical field to the computational domain respectively, establish the coupling relationship between them, and thus form a thermo-hydraulic-mechanical-chemical coupling model covering crack deformation;

[0095] The thermo-hydraulic-mechanical-chemical coupling model covering crack deformation includes: flow equations in porous media, flow equations in cracks, transport equations for diluted species in porous media, transport equations for diluted species in cracks, heat transfer equations in porous media, heat transfer equations in fractures, solid mechanics equations, equations for thin elastic layers (cracks), boundary ordinary differential equations and partial differential algebraic equations (crack chemical aperture), and neighborhood ordinary differential equations and partial differential algebraic equations (chemical and mechanical porosity).

[0096] The flow equation of the acidizing fluid in the porous rock matrix, i.e. the flow equation of the porous medium, is shown in equation (1):

[0097] (1);

[0098] The flow equation of the acidizing fluid in the crack, i.e. the crack flow equation, is shown in equation (2):

[0099] (2);

[0100] in, The density of the acidified fluid. m It is the storage coefficient of the rock matrix;

[0101] f The storage coefficient for cracks. It's pressure. It is time. It is the fluid velocity vector of the acidifying fluid in the porous rock matrix. It is the fluid velocity vector of the acidizing fluid in the crack. It is the Biot-Willis coefficient of the rock matrix. It is the Biot-Willis coefficient of the crack; It is the volumetric strain of the rock matrix. It is the volumetric strain of the crack. It refers to the size of the crack opening. It is a gradient operator;

[0102] It is the mass flux from the rock matrix to the fracture, as shown in equation (3):

[0103] (3);

[0104] in, m It is the permeability of the rock matrix. acid It refers to the viscosity of the acidified fluid. n u It is the normal vector pointing to the rock-crack interface on the side facing the rock. n d It is the normal vector pointing towards the crack side of the rock-crack interface.

[0105] In carbonate geothermal reservoirs (such as limestone), the main component of the rock is calcium carbonate (CaCO3). It can react with hydrochloric acid (HCl), which is commonly used in fracturing processes, as shown in equation (4):

[0106] (4);

[0107] Equation (4) proves that 1 mol of HCl can dissolve 0.5 mol of CaCO3. To describe the chemical reaction, based on the variable of HCl, its transport in the porous matrix is ​​controlled by Equation (5), which is the transport equation of dilution species in the porous medium;

[0108] (5)

[0109] in, It refers to the porosity of the rock matrix. C HCl It is the concentration of HCl. D e It is the effective diffusion coefficient, which can be obtained by... tortuosity of the crack The rate of acid-rock chemical reactions was calculated using diffusion coefficients. This was determined under conditions of sufficient solid reactants. As shown in equation (6);

[0110] (6);

[0111] Among them, reaction order for The reaction rate constant is shown in equation (7):

[0112] (7);

[0113] in, A It is the Arrhenius constant. e It is the base of the natural logarithm. E It is activation energy. R g It is the molar gas constant. ，T It is absolute temperature.

[0114] In addition to the chemical reaction process, the transport of reactants and products in the cracks should also be considered. The transport equation for dilute species in the cracks is shown in equation (8):

[0115] (8);

[0116] in, This indicates the rate of acidification reaction at the fluid-solid interface. Q c The mass flux from the porous matrix into the crack is represented by Equation (9):

[0117] (9);

[0118] It can be determined based on mass transfer, as shown in equation (10):

[0119] (10);

[0120] in, It is the mass transfer coefficient, as shown in equation (11):

[0121] (11);

[0122] in, The solute diffusion coefficient in the fracture ( i = HCl or CaCl2). The Reynolds number is... For the Schmitt number, ; That is, the interface area of ​​the crack per unit volume. This refers to the HCl concentration at the fluid-solid interface. Because... For first-order kinetic reactions ( (reaction rate constant) It can be expressed as equation (12):

[0123] (12);

[0124] As shown in equation (12) in the current process, this indicates that the acid is consumed immediately upon transfer to the liquid-solid interface. Since the reaction rate is limited by the mass transfer rate here, this reaction can be classified as a mass transfer-controlled reaction.

[0125] The water-rock thermal equilibrium assumption was also adopted, which assumes that the water and the nearby rocks have the same temperature. Considering heat convection, heat conduction, and the heat source generated by the reaction, the law of conservation of energy can be obtained.

[0126] The heat transfer equation in porous media is shown in equation (13):

[0127] (13);

[0128] in, It is the effective volumetric specific heat capacity of the rock matrix. It is the effective thermal conductivity of the rock matrix. It is the effective density of the rock matrix. It is the specific heat capacity under constant pressure. It's temperature. It is a heat source generated by a chemical reaction, as shown in equation (14);

[0129] (14)

[0130] in, This is the enthalpy of the reaction; since the acidification reaction is exothermic, it is a negative value.

[0131] The heat transfer equation during fracture is shown in equation (15).

[0132] (15);

[0133] in, The effective volumetric heat capacity of the crack. The effective thermal conductivity of the crack, The heat source for the reaction in the crack is as shown in equation (16);

[0134] (16);

[0135] It is the heat source term from the matrix to the fracture, as shown in equation (17);

[0136] (17);

[0137] For rock materials, applying linear elasticity (Hooke's Law) and small deformation theory, the equilibrium equation for the deformation of the rock matrix under pure gravity load is derived, and the solid mechanics equation expression is shown in equation (18).

[0138] (18)

[0139] in, and These are Lamé constants, It is volumetric strain, It is the displacement vector and It is the bulk modulus;

[0140] For cracks, since their size is significantly smaller compared to the size of the surrounding matrix, in numerical simulations, cracks can be simplified to a low-dimensional thin elastic layer with a specific nominal thickness. Figure 1 The crack layer follows Hooke's law, therefore the equation for the thin elastic layer (crack) is as shown in equation (19);

[0141] (19)

[0142] in, Force / unit area , The displacement vectors at the upper and lower boundaries of the crack layer are respectively... Initial compensation displacement.

[0143] Stiffness can be measured by normal stiffness and tangential stiffness They are obtained by means of two corresponding deformations, as shown in equation (20);

[0144] (20)

[0145] in, Let be the unit normal vector of the crack. The normal and tangential stiffness of the thin elastic layer can be calculated using its material parameters, as shown in equation (21);

[0146] (twenty one);

[0147] in, It is the elastic modulus. It is Poisson's ratio. It is the shear modulus.

[0148] Acidification reactions lead to increased porosity and fracture aperture in the rock matrix, thereby increasing the permeability of carbonate geothermal reservoirs. For the rock matrix, the relationship between porosity changes and reactions can be described by the continuity equation of the rock skeleton, as shown in equation (22).

[0149] (twenty two);

[0150] in, It is matrix porosity caused by chemical reactions. It is the number of moles of calcium carbonate (CaCO3) per unit mass of rock. and These represent the content and molar mass of calcium carbonate minerals in the reservoir, respectively.

[0151] Furthermore, the porosity change caused by the thermo-porous elastic effect can be expressed as follows, as shown in equation (23);

[0152] (twenty three);

[0153] in, Indicates volumetric strain. This indicates tangential stiffness.

[0154] Combining equation (22) and matrix (23), the porosity obtained by coupling the THMC model is shown in equation (24);

[0155] (twenty four);

[0156] The permeability of the rock matrix at a certain time point can be determined by the following formula, as shown in Equation (25);

[0157] (25)

[0158] in, It is the initial permeability of the rock matrix.

[0159] For cracks, the boundary ordinary differential equation and partial differential algebraic equation, i.e. the crack pore size change rate caused by chemical corrosion, are shown in equation (26).

[0160] (26);

[0161] In addition, the ordinary differential equations and partial differential algebraic equations of the field, namely the crack opening change caused by mechanical deformation, can be calculated by using the parameters of the thin elastic layer and the normal strain of the layer and the fluid pressure in the crack, as shown in Equation (27).

[0162] (27);

[0163] Under the coupled thermo-hydraulic-mechanical (THMC) conditions, the crack opening is finally obtained by integrating equations (26) and (27), which is equation (28).

[0164] (28);

[0165] The permeability of the crack can be calculated using a parallel plate flow model, as shown in equation (29);

[0166] (29);

[0167] Based on the above, a coupled THMC model containing multiple mechanisms was established.

[0168] Step 4: Assign initial conditions and boundaries, create and optimize the mesh, solve the decoupled model, and evaluate acid fracturing;

[0169] Specifically, numerical calculations for THMC (thermal-hydraulic-mechanical-chemical) modeling were performed using the finite element method in COMSOL software, further considering the influence of discrete fracture networks (DFNs) in acid fracturing. These fractures are randomly generated, with their lengths following a power-law distribution and their directions following a Fisher distribution. The "Darcy's Law" physics module was used to calculate the acid flow rate in the matrix (Equation (1)) and fractures (Equation (2)). The "Transport of Dilute Species in Porous Media" module was used to calculate the solute concentration in the matrix (Equation (5)) and fractures (Equation (8)). The "Heat Transfer in Porous Media" module was used to determine the temperature distribution in the matrix (Equation (13)) and fractures (Equation (15)). In addition, the "Solid Mechanics" module was used to obtain the deformation of the rock matrix (Equation (18)). Subsequently, in order to obtain the change in matrix porosity due to chemical corrosion (Equation (22)), the Domain Ordinary Differential Equations (ODEs) and Differential Algebraic Equations (DAEs) modules were activated. For the change in crack aperture (Equation (26)), the Boundary ODEs and DAEs modules were activated. For modeling, the computational domain was discretized using a free triangular mesh. The calculation started from the default initial step, and the variables were updated using implicit backward difference formulas. In order to improve the computational convergence of the proposed model, a separate coupled solver was used to solve the separate coupled model, perform numerical calculations and sensitivity analysis on the coupled THMC model, and evaluate the effect of acid fracturing.

[0170] The method of this invention reveals the significant influence of factors such as fracture network structure, acid injection rate, and acid concentration on fracture porosity changes, reaction rates, and pressure evolution. In particular, fracture connectivity plays a crucial role in fluid flow and pressure distribution.

[0171] This invention's model, based on the assumptions of a thin elastic layer and fracture units, provides a detailed simulation of the chemical corrosion and mechanical deformation of fractures, including fracture propagation and closure under chemical action, and the resulting changes in matrix porosity. Furthermore, the model considers the dynamic evolution of fracture aperture and matrix porosity, crucial for understanding long-term changes in reservoir permeability. Simulations of the effects of acid fracturing under different reservoir conditions and operating parameters reveal the influence of fracture location and fracture network type on fracture aperture differences. The study found that the injection rate is a key factor affecting the average chemical aperture, while the HCl concentration has a positive effect on increasing the average chemical aperture. Moreover, there is a negative correlation between the initial fracture aperture and the peak average reaction rate, and a positive correlation with the average chemical aperture. The initial reservoir temperature also significantly affects the chemical reaction rate and fracture aperture. Chemical factors play a crucial role in fracture propagation, while mechanical factors significantly constrain the process. In contrast, temperature has a relatively small impact on fracture propagation. These findings provide valuable scientific evidence for optimizing acid fracturing in carbonate reservoirs and can help guide practical geothermal development. The practicality and effectiveness of the proposed model were further verified through the study of multiple real-world cases.

[0172] Example 3

[0173] The acid fracturing evaluation method for carbonate geothermal reservoirs based on a coupled model in this embodiment is as follows:

[0174] Taking the coupled THMC model for numerical calculation and sensitivity analysis as an example, the study uses a plane strain square rock with a dimension of 200m in one direction. Figure 2 Boundary conditions for the thermal module (T), hydraulic module (H), mechanical module (M), and chemical module (C) are shown. The temperature of the left boundary is kept constant at 80°C, while the other three boundaries are isolated. The left boundary is given a constant inflow rate, and the right boundary is under a constant pressure of 10 MPa. The left and right boundaries are rollers, and the lower boundary is fixed. In addition, the upper boundary is subjected to a uniformly distributed external load. Note that sufficient computational steps were taken to reach the initial state of geostress equilibrium before injection. For the chemical field, acid is injected from the left boundary at a constant concentration, while the other three boundaries are impermeable. The physical and mechanical parameters used in the simulation are shown in Table 1. The physical and mechanical parameters of the carbonate rock matrix include: Young's modulus, Poisson's ratio, density, Biot ratio, coefficient of thermal expansion, initial porosity, initial permeability, storage coefficient, thermal conductivity, Poisson's ratio, and specific heat capacity; the engineering parameters of the fracture include: fracture width, Young's fracture modulus of the fracture, fracture storage coefficient, and initial pore size of the fracture.

[0175] Table 1. Parameters used in the simulation

[0176]

[0177] Based on engineering parameters, a geometric model simulating a carbonate reservoir was established using COMSOL software. The formation of the fracture network was simulated using COMSOL to obtain the variables in the computational domain and fractures. Heat transfer, hydraulic, stress, and chemical fields were added to the computational domain, and their coupling relationships were established, thus forming a thermo-hydraulic-mechanical-chemical coupled model covering fracture deformation.

[0178] A crack network was added to the computational domain. The crack lengths followed an exponential distribution, and the crack dip angles were considered random and uniformly distributed along the cracks. Figure 3 As shown, a discrete crack network was generated. Figure 3 The image shows a Discrete Fracturing Network (DFN) consisting of 40 fractures ranging in length from 50 to 100 meters. In this case, the computational domain is meshed using 82,734 free triangle elements. Six observation points and one observation seam are established within the fracture network. The formation of the fracture network is simulated using interpolation functions in COMSOL, yielding variables within the computational domain and fractures. Subsequently, heat transfer, hydraulic, stress, and chemical fields are added to the computational domain, establishing coupling relationships between them, thus forming a thermo-hydraulic-mechanical-chemical coupled model covering fracture deformation (e.g., ...). Figure 4 (As shown). Initial conditions and boundaries are assigned, mesh creation and optimization are performed, decoupling is solved, and acid fracturing evaluation is conducted.

[0179] Example 4

[0180] Figure 5 The distribution of pressure, acid concentration, temperature, and chemical pore width (the width of cracks resulting from the chemical reaction) at four different time points—30, 60, 90, and 120 minutes—is shown. The results change over time as acid injection progresses. At 30 minutes, the pressure in most areas of the site exceeds 30 MPa, with the highest local pressure reaching approximately 45 MPa. The pressure gradually decreases over time, with the overall value approaching the initial pressure of 10 MPa. Regarding acid concentration, at 30 minutes of injection, the concentration is higher in the cracks near the inflow boundary, reaching up to 8750 mol / m³. 3Subsequently, as the solute diffuses, the acid concentration in the internal cracks connected to the cracks near the inflow area gradually increases. On the other hand, for unconnected cracks, the acid concentration changes very little. It is clear from the temperature diagram that the temperature of most cracks is close to the initial injection temperature of 80°C, while for a few cracks, the temperature is lower due to the injection of low-temperature acid. The interconnection between cracks improves the liquid flow in the network, which increases the reaction rate. In this case, crack corrosion becomes faster, convection is enhanced, and therefore the crack temperature drops rapidly. Furthermore, the width of the chemical cracks gradually increases over time, with some reaching approximately 5 mm. Identifying these cracks with large chemical crack widths shows their basic consistency with the low-temperature cracks observed in the temperature diagram, which is due to chemical corrosion. The variation in crack width is influenced not only by crack width but also by factors such as pore pressure and temperature. Among these factors, temperature has a smaller effect because the temperature exposure time is short and heat transfer is limited. Furthermore, the thermal expansion of the matrix is ​​less than the deformation of the crack.

[0181] Example 5

[0182] To more accurately describe the concentration changes of reactant HCl and product CaCl2, different observation fractures (i.e., fractures) were selected from the fracture network for detailed examination. 26.8% HCl was continuously injected at a controlled rate for 20 minutes. Figure 6 In the DFN, the concentration of HCl decreases as the reactant advances along the crack and is consumed. In the initial stages of DFN, the reactant concentration decreases significantly, marked by a rapid inflection point, after which the change tends to plateau. As the acidification reaction becomes more complete, the zero concentration point of HCl advances further. Specifically, at 20 minutes, the zero concentration point of the reactant is located at approximately 40 meters, and by 120 minutes, the HCl concentration drops to zero at approximately 30 meters. The CaCl2 concentration in DFN generally increases first and then decreases; the increase is due to the onset of the acidification reaction, while the decrease is due to the gradual outflow of products from the crack. At the observed crack location, the CaCl2 concentration increases slowly at 20 minutes but drops to zero at approximately 40 meters. From 40 minutes to 120 minutes, the product concentration increases rapidly, rises slowly to a peak at approximately 2.5 meters, and then decreases slowly. Simultaneously, the peak value increases over time. Figure 7 The CaCl2 concentration in DFN generally showed an initial increase followed by a decrease. This increase was attributed to the onset of the acidification reaction, while the subsequent decrease was a result of the gradual outflow of products from the crack. Specifically, along the observed crack, the CaCl2 concentration gradually increased over the first 20 minutes, then dropped to zero at 40 minutes. Subsequently, during the period from 40 to 120 minutes, the product concentration rapidly increased, reaching a peak at approximately 2.5 m, before gradually decreasing. Simultaneously, the peak value showed an upward trend over time.

[0183] Example 6

[0184] The chemical pore width of a crack is related to the extent of the chemical reaction. Figure 8 and Figure 9 The spatial and temporal distributions of the observed chemical pore widths of fractures in the fracture network are presented. For fractures, such as... Figure 8 As shown, in DFN, the chemical pore width generally exhibits a decreasing trend, dropping to zero at approximately 30 meters. Notably, the chemical pore width increases over time at the leftmost end of the fracture. Figure 9 As shown, in DFN, the chemical pore width increases very slowly within 40 meters. Subsequently, the chemical pore width increases rapidly near the fracture confluence. Over time, the magnitude of the change becomes more pronounced, with the peak value reaching a higher level.

[0185] This invention designs a fully coupled THMC model for acid fracturing processes in carbonate reservoirs. This model integrates the interactions between heat transfer, acidizing fluid flow, deformation of rock matrix and fractures, migration of reactive diluents, and chemical corrosion phenomena, and is used to quantify the impact of acid fracturing on carbonate geothermal reservoirs.

Claims

1. A method for evaluating acid fracturing in carbonate geothermal reservoirs based on a thermo-hydraulic-mechanical-chemical coupling model, characterized in that, The specific steps are as follows: Step 1: Collect the physical and mechanical parameters of the carbonate reservoir matrix and the engineering parameters of the fractures and fissures; Step 2: Based on the engineering parameters, establish a geometric model of the simulated carbonate reservoir, and assign the variables in the fractures in each field to the fracture domain to obtain the computational domain; Step 3: Add heat transfer field, hydraulic field, stress field and chemical field to the computational domain respectively, establish the coupling relationship between them, and thus form a thermo-hydraulic-mechanical-chemical coupling model covering crack deformation; The thermo-hydraulic-mechanical-chemical coupling model covering crack deformation includes: porous media flow equation, crack flow equation, transport equation of diluted species in porous media, transport equation of diluted species in cracks, heat transfer equation in porous media, heat transfer equation in fracture, solid mechanics equation, thin elastic layer equation, boundary ordinary differential equation and partial differential algebraic equation, and neighborhood ordinary differential equation and partial differential algebraic equation. The flow equation of the acidizing fluid in the porous rock matrix, i.e. the flow equation of the porous medium, is shown in equation (1): （1）； The flow equation of the acidizing fluid in the crack, i.e. the crack flow equation, is shown in equation (2): （2）； in, The density of the acidified fluid. m It is the storage coefficient of the rock matrix; f The storage coefficient for cracks. It's pressure. It is time. It is the fluid velocity vector of the acidifying fluid in the porous rock matrix. It is the fluid velocity vector of the acidizing fluid in the crack. It is the Biot-Willis coefficient of the rock matrix. It is the Biot-Willis coefficient of the crack; It is the volumetric strain of the rock matrix. It is the volumetric strain of the crack. It refers to the size of the crack opening. It is a gradient operator; It is the mass flux from the rock matrix to the fracture, as shown in equation (3): （3）； in, m It is the permeability of the rock matrix. acid It refers to the viscosity of the acidified fluid. n u It is the normal vector pointing to the rock-crack interface on the side facing the rock. n d It is the normal vector pointing towards the crack side of the rock-crack interface; Based on the variable of HCl, its transport in the porous matrix is ​​controlled by formula (5), which is the transport equation of dilution species in the porous medium; （5） in, It refers to the porosity of the rock matrix. C HCl It is the concentration of HCl. D e It is the effective diffusion coefficient, which can be obtained by... The tortuosity of the crack The rate of acid-rock chemical reactions was calculated using diffusion coefficients and other parameters, given a sufficient amount of solid reactants. As shown in equation (6); （6）； Among them, reaction order for The reaction rate constant is shown in equation (7): （7）； in, A It is the Arrhenius constant. e It is the base of the natural logarithm. E It is activation energy. R g It is the molar gas constant. ，T It is absolute temperature; The transport equation for dilute species in the crack is shown in equation (8): （8）； in, This indicates the rate of acidification reaction at the fluid-solid interface. Q c The mass flux from the porous matrix into the crack is represented by Equation (9): （9）； It can be determined based on mass transfer, as shown in equation (10): （10）； in, It is the mass transfer coefficient, as shown in equation (11): （11）； in, The solute diffusion coefficient in the fracture is... The Reynolds number is... For the Schmitt number, ; That is, the interface area of ​​the crack per unit volume. It is the HCl concentration at the fluid-solid interface; It can be expressed as equation (12): （12） The heat transfer equation in porous media is shown in equation (13): （13）； in, It is the effective volumetric specific heat capacity of the rock matrix. It is the effective thermal conductivity of the rock matrix. It is the effective density of the rock matrix. It is the specific heat capacity under constant pressure. It's temperature. It is a heat source generated by a chemical reaction, as shown in equation (14); （14） in, For reaction enthalpy; The heat transfer equation during fracture is shown in equation (15). （15）； in, The effective volumetric heat capacity of the crack. The effective thermal conductivity of the crack, The heat source for the reaction in the crack is as shown in equation (16); （16）； It is the heat source term from the matrix to the fracture, as shown in equation (17); （17）； The equation for solid mechanics is shown in equation (18); （18） in, and These are Lamé constants, It is volumetric strain, It is the displacement vector and It is the bulk modulus; The equation for the thin elastic layer is shown in equation (19); （19） in, Force / unit area , The displacement vectors at the upper and lower boundaries of the crack layer are respectively... Initial compensation displacement; Stiffness can be measured by normal stiffness and tangential stiffness They are obtained by means of two corresponding deformation directions, as shown in equation (20); （20） in, The unit normal vector of the crack; the normal and tangential stiffness of the thin elastic layer can be calculated using its material parameters, as shown in equation (21); （21）； in, It is the elastic modulus. It is Poisson's ratio. It is the shear modulus; For cracks, the boundary ordinary differential equation and partial differential algebraic equation, i.e. the crack pore size change rate caused by chemical corrosion, are shown in equation (26). （26）； The ordinary differential equations and partial differential algebraic equations of the field, namely the crack opening change caused by mechanical deformation, can be calculated by using the parameters of the thin elastic layer and the normal strain of the layer and the fluid pressure in the crack, as shown in equation (27). （27）； Under the coupled conditions of thermo-hydraulic-mechanical-chemical, the crack opening is finally obtained by integrating equations (26) and (27), which is equation (28). （28）； The permeability of the crack can be calculated using a parallel plate flow model, as shown in equation (29); （29）； A thermo-hydraulic-mechanical-chemical coupling model covering crack deformation was established based on the above. Step 4: Solve the decoupled model and perform acid fracturing evaluation.

2. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermo-hydraulic-mechanical-chemical coupling model as described in claim 1, characterized in that, In step 1, the physical and mechanical parameters of the carbonate reservoir matrix include: Young's modulus, Poisson's ratio, density, Biot ratio, coefficient of thermal expansion, initial porosity, initial permeability, storage coefficient, thermal conductivity, Poisson's ratio, and specific heat capacity; the engineering parameters of the fractures include: fracture width, Young's fracture modulus of the fracture, storage coefficient of the fracture, and initial pore size of the fracture.

3. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermo-hydraulic-mechanical-chemical coupling model as described in claim 2, characterized in that, In step 2, the specific steps for establishing a geometric model of the simulated carbonate reservoir using COMSOL software are as follows: Based on the engineering geological characteristics of the carbonate reservoir, a computational geometric model is established; the geometry and size of the reservoir are defined on the canvas; the thickness of the fractures in the carbonate reservoir is determined according to engineering parameters; and the variables in the fractures in each field are assigned to the fracture domain through functions in COMSOL to obtain the computational domain.

4. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermo-hydraulic-mechanical-chemical coupling model as described in claim 3, characterized in that, In step 4, specifically: a decoupled solver is used to solve the decoupled model, numerical calculations and sensitivity analysis are performed on the coupled model, and the acid fracturing effect is evaluated.

Images