Carbonate geothermal reservoir acid fracturing evaluation method based on heat-water-force-chemical coupling model

Abstract

The invention discloses a carbonatite geothermal reservoir acid fracturing evaluation method based on a heat-water-force-chemical coupling model. The carbonatite geothermal reservoir acid fracturing evaluation method specifically comprises the steps that physical and mechanical parameters of a carbonatite reservoir matrix and engineering parameters of fracture fractures are collected; according to the engineering parameters, establishing a geometric model for simulating the carbonate reservoir, and distributing variables in cracks in each field to a crack domain to obtain a computational domain; a heat transfer field, a hydraulic field, a stress field and a chemical field are added to the computational domain, the mutual coupling relation is established, and a heat-water-force-chemical coupling model covering fracture deformation is formed; the separation coupling model is solved, and acid fracturing evaluation is carried out. According to the invention, the fully-coupled THMC model of the carbonate reservoir acid fracturing process is designed, and the model integrates heat transfer, acidizing fluid flow, rock matrix and crack deformation, reaction dilution substance migration and interaction among chemical corrosion phenomena, and is used for quantifying the influence of carbonate geothermal reservoir acid fracturing.

Description

Technical Field

[0001] The invention belongs to the technical field of geothermal exploitation, and in particular relates to an acid fracturing evaluation method for carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model. Background Art

[0002] Among the hydrothermal reservoir types, carbonate geothermal reservoirs are an important reservoir type, and their large-scale development and utilization have become the focus of future geothermal exploration and development. At present, hydraulic fracturing, acid fracturing, sand fracturing, high-energy gas fracturing and other fracturing processes are mainly used. Among these technologies, the acid fracturing method is of great significance for the strengthening of carbonate geothermal systems. Acid fracturing uses a displacement exceeding the reservoir capacity to inject acid or acid pre-liquid from the casing to quickly establish a wellbore pressure that exceeds the formation compression stress and the rock tensile strength. This process causes the formation to break and cracks to form. Continuous injection of acid can extend these cracks, thereby forming acid-etched fractures, which have enhanced conductivity compared to the original formation. This not only increases the production capacity of geothermal wells, but also enhances their backfilling capacity.

[0003] At present, acid fracturing technology has been receiving extensive attention. Ma Jianjun and his research team proposed a thermal-hydro-chemical (THC) coupling model, using a unified pipe network method to simulate the acid fracturing process, in which the interaction between stress and THC coupling was ignored. Researchers such as Mahmoodi showed through parameter sensitivity analysis that the acid fracturing effect has a significant impact on the formation and expansion of cavities. In the prior art, a coupled THMC model is established to appropriately consider the effect of heat generation on the reaction rate during the rock acid fracturing reaction, thereby affecting the chemical process. The THMC model is also used to study the effects of chemical corrosion and precipitation, injection concentration, injection rate and other factors on the evolution of fracture aperture. From the above research, it can be seen that exploring the THMC coupling mechanism is one of the key factors in optimizing the development of geothermal resources.

[0004] Discrete fracture networks (DFNs) are crucial in understanding geothermal reservoirs as they enable identification of geothermal fluid flows. The application of DFNs in fluid flows has been explored by several researchers. The enhanced connectivity of fracture networks obtained through DFNs can significantly improve the solute concentration profile. A coupled thermo-hydraulic-mechanical (THM) model based on embedded DFN has also been proposed in the prior art to simulate the complex flow behavior and thermal performance during CO2 injection, thus verifying the importance of DFN in analyzing heat flow, fracture deformation, and flow characteristics. In studies of geothermal reservoirs containing natural fractures, DFN is the preferred method for simulating the behavior and distribution of natural fractures. The DFN model accurately captures the behavior and distribution of natural fractures and believes that they significantly affect stress changes and damage evolution. In carbonate reservoirs, the application of DFN has enabled researchers to study the presence, pore size, distribution, and sensitivity of dissolution processes to acid injection rates, thereby revealing the significant impact of fracture structure on reactive flow. However, the role of DFN features in the acid fracturing process of carbonate geothermal rocks in the THMC coupled model has rarely been comprehensively studied. Therefore, it is of great significance to consider the DFN characteristics based on the coupled THMC model and study the chemical corrosion and mechanical deformation of fractures during acid fracturing in carbonate geothermal reservoirs. Summary of the invention

[0005] The purpose of the present invention is to provide an acid fracturing evaluation method for carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model, which integrates the interaction between heat transfer, acidizing fluid flow, deformation of rock matrix and fractures, migration of reaction dilution substances and chemical corrosion phenomena, and is used to quantify the impact of acid fracturing on carbonate geothermal reservoirs.

[0006] The technical solution adopted by the present invention is an acid fracturing evaluation method for carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model, which is specifically implemented in the following steps:

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

[0008] Step 2: According to the engineering parameters, a geometric model for simulating the carbonate reservoir is established, and the variables in the fractures in each field are assigned to the fracture domain to obtain the calculation domain;

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

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

[0011] The present invention is also 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, thermal expansion coefficient, initial porosity, initial permeability, storage coefficient, thermal conductivity, Poisson's ratio and specific heat capacity; the engineering parameters of the fracture fissure include: crack width, Young's modulus of rupture of the crack, storage coefficient of the crack, and initial pore size of the crack.

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

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

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

[0016]

[0017] The flow equation of the acidizing fluid in the fracture, that is, the fracture flow equation, is shown in equation (2):

[0018]

[0019] Among them, ρ acid is the density of the acidizing fluid, S m is the storage coefficient of the rock matrix;

[0020] S f is the storage coefficient of the crack, p is the pressure, t is the time, v m is the fluid velocity vector of the acidizing fluid in the porous rock matrix, v f is the fluid velocity vector of the acidizing fluid in the fracture, α Bm is the Biot-Willis coefficient of the rock matrix, α mf is the Biot-Willis coefficient of the crack; ε vm is the volumetric strain of the rock matrix, ε vf is the volume strain of the crack, d f is the size of the crack opening, is the gradient operator;

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

[0022]

[0023] Among them, k m is the permeability of the rock matrix, μ acid is the viscosity of the acidizing fluid, n u is the normal vector pointing to the rock-crack interface toward the rock side, n d is the normal vector pointing to the rock-crack interface toward the crack side;

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

[0025]

[0026] Among them, φ m is the porosity of the rock matrix, C HCl is the concentration of HCl, D e is the effective diffusion coefficient, which can be expressed by D e =φ f D f / τ f The tortuosity of the crack f =φ -13 and diffusion coefficient; in the presence of sufficient solid reactants, the rate R of the acid-rock chemical reaction m HCl , as shown in formula (6);

[0027]

[0028] Among them, the reaction order K is n a The reaction rate constant is shown in formula (7):

[0029]

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

[0031] The transport equation of diluted species in fractures is shown in equation (8):

[0032]

[0033] in, represents the acidification reaction rate at the fluid-solid interface, Q c represents the mass flux from the porous matrix into the crack, as shown in formula (9):

[0034]

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

[0036]

[0037] Among them, κ c is the mass transfer coefficient, as shown in formula (11):

[0038]

[0039] in, is the solute diffusion coefficient in the crack, N Re is the Reynolds number, N Sc is the Schmidt number, a vf =2 / d f , that is, the interfacial area of ​​the crack per unit volume, is the HCl concentration at the fluid-solid interface; It can be expressed as formula (12):

[0040]

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

[0042]

[0043] Among them, (ρc p ) eff,m =φ m ρ acid c p,acid +(1-φ m )ρ s c p,s is the effective volume specific heat capacity of the rock matrix, λ eff,m =φλ acid +(1-φ)λ s is the effective thermal conductivity of the rock matrix, ρ is the effective density of the rock matrix, c p is the specific heat at constant pressure, T is the temperature, is the heat source generated by the chemical reaction, as shown in formula (14);

[0044]

[0045] Where H is the reaction enthalpy;

[0046] Heat transfer equation in fracture: as shown in equation (15);

[0047]

[0048] Among them, (ρc p ) eff,m ≈ρ acid c p,acid is the effective volume heat capacity of the crack, λ eff,m ≈λ acid is the effective thermal conductivity of the crack, is the reaction heat source in the crack, as shown in formula (16);

[0049]

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

[0051]

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

[0053]

[0054] Where L and G are the Lamé constant, is the volume strain, u is the displacement vector and K v is the bulk modulus;

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

[0056] u u -u d -u0=-F / k A (19)

[0057] Where F represents force / unit area u u ,u d are the displacement vectors of the upper and lower boundaries of the fracture layer and the initial compensation displacement u0 respectively;

[0058] The stiffness can be calculated by the normal stiffness k A and tangential stiffness k n To obtain, they represent the deformation in two corresponding directions, as shown in formula (20);

[0059]

[0060] Where n is the unit normal vector of the crack; the normal and tangential stiffness of the thin elastic layer can be calculated by its material parameters, as shown in Equation (21);

[0061]

[0062] Among them, E f is the elastic modulus, ν f is Poisson's ratio, G f is the shear modulus.

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

[0064]

[0065] The ordinary differential equations and partial differential algebraic equations in the field, i.e., the change in the crack opening due to mechanical deformation can be calculated by using the parameters of the thin elastic layer as well as the normal strain of the layer and the fluid pressure in the crack, as shown in Eq. (27);

[0066]

[0067] Under the thermo-hydraulic-chemical (THMC) coupling condition, the crack opening is finally obtained by integrating equations (26) and (27), which is equation (28);

[0068]

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

[0070]

[0071] Based on the above, a thermal-hydraulic-mechanical-chemical coupling model of covering crack deformation was established.

[0072] In step 4, specifically: using a separation coupling solver to solve the separation coupling, performing numerical calculation and sensitivity analysis on the coupling model, and evaluating the acid fracturing effect.

[0073] Beneficial effects of the present invention:

[0074] In the development of carbonate geothermal reservoirs, the implementation of acid fracturing technology plays a vital role in improving reservoir permeability. The present invention simulates the complex phenomena involved in the acid fracturing process, including solute migration, acid-rock reaction, heat conduction, and formation deformation, by constructing a comprehensive thermal-hydro-mechanical-chemical coupling model. The model is specifically targeted at field-scale discrete fracture networks, comprehensively considering the chemical corrosion and mechanical deformation of the fractures, as well as the chemical corrosion effects of the bedrock and fracture surfaces, and the model is able to track the changes in fracture aperture and matrix porosity over time. BRIEF DESCRIPTION OF THE DRAWINGS

[0075] Figure 1 It is a flow chart of the carbonate geothermal reservoir acid fracturing evaluation method based on the thermal-hydraulic-mechanical-chemical coupling model of the present invention;

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

[0077] Figure 3 is a diagram of the geometry of natural fracture networks with different numbers and sizes of fractures;

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

[0079] Figure 5 It is a schematic diagram of the pressure distribution of DFN at different times of the present invention;

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

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

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

[0083] Fig. 9 This is the spatiotemporal distribution diagram of the chemical aperture of the DFN observation line of the present invention (II). DETAILED DESCRIPTION

[0084] The present invention is described in detail below with reference to the accompanying drawings and specific embodiments.

[0085] Example 1

[0086] The present invention discloses an acid fracturing evaluation method for carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model, which specifically comprises the following steps: collecting physical and mechanical parameters of a carbonate reservoir matrix and engineering parameters of fracture fissures; using COMSOL software to establish a geometric model simulating a carbonate reservoir according to the engineering parameters, and assigning variables in the fractures in each field obtained to a fracture domain through functions in COMSOL to obtain a calculation domain; adding a heat transfer field, a hydraulic field, a stress field and a chemical field to the calculation domain to establish a coupling relationship between them, thereby forming a thermal-hydraulic-mechanical-chemical coupling model covering fracture deformation; allocating initial conditions and boundaries, creating and optimizing a mesh, solving a separation coupling model, and performing acid fracturing evaluation;

[0087] Example 2

[0088] The present invention provides an acid fracturing evaluation method for carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model. Figure 1 As shown, the specific implementation steps are as follows:

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

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

[0091] The engineering parameters of fracture fissures include: crack width, Young's rupture modulus, storage coefficient, and initial aperture of the crack.

[0092] Step 2: Based on the engineering parameters, COMSOL software is used to establish a geometric model for simulating carbonate reservoirs. Through the functions in COMSOL, the variables in the fractures in each field are assigned to the fracture domain to obtain the calculation domain;

[0093] The specific steps of using COMSOL software to establish the geometric model of the carbonate reservoir to be simulated are as follows: Based on the engineering geological characteristics of the carbonate reservoir, a computational geometry model is established. The geometric shape and size of the reservoir are defined on the canvas, and the thickness of the fractures in the carbonate reservoir is determined based on the engineering parameters. Through the functions in COMSOL, the variables in the fractures in each field are assigned to the fracture domain to obtain the computational domain;

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

[0095] The coupled thermo-hydraulic-mechanical-chemical model covering fracture deformation includes: porous media flow equations, fracture flow equations, transport equations for diluted species in porous media, transport equations for diluted species in fractures, heat transfer equations in porous media, heat transfer equations in fractures, solid mechanics equations, thin elastic layer (crack) equations, boundary ordinary differential equations and partial differential algebraic equations (fracture chemical aperture), field 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, that is, the porous medium flow equation, is shown in formula (1):

[0097]

[0098] The flow equation of the acidizing fluid in the fracture, that is, the fracture flow equation, is shown in equation (2):

[0099]

[0100] Among them, ρ acid is the density of the acidizing fluid, S m is the storage coefficient of the rock matrix;

[0101] S f is the storage coefficient of the crack, p is the pressure, t is the time, v m is the fluid velocity vector of the acidizing fluid in the porous rock matrix, v f is the fluid velocity vector of the acidizing fluid in the fracture, α Bm is the Biot-Willis coefficient of the rock matrix, α mf is the Biot-Willis coefficient of the crack; ε vm is the volumetric strain of the rock matrix, ε vf is the volume strain of the crack, d f is the size of the crack opening, is the gradient operator;

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

[0103]

[0104] Among them, k m is the permeability of the rock matrix, μ acid is the viscosity of the acidizing fluid, n u is the normal vector pointing to the rock-crack interface toward the rock side, n d is the normal vector pointing to the rock-crack interface towards the crack side.

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

[0106] CaCO3+2HCl→CaCl2+H2O+CO2 (4);

[0107] According to formula (4), it is proved that 1 mol HCl can dissolve 0.5 mol CaCO3. In order to describe the chemical reaction, based on the variable of HCl, its transport in the porous matrix is ​​controlled by formula (5), which is the transport equation of dilute species in porous media;

[0108]

[0109] Among them, φ m is the porosity of the rock matrix, CHCl is the concentration of HCl, D e is the effective diffusion coefficient, which can be expressed by D e =φ f D f / τ f The tortuosity of the crack f =φ -13 The rate of the acid-rock chemical reaction in the presence of sufficient solid reactants is As shown in formula (6);

[0110]

[0111] Among them, the reaction order K is n a The reaction rate constant is shown in formula (7):

[0112]

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

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

[0115]

[0116] in, represents the acidification reaction rate at the fluid-solid interface, Q c represents the mass flux from the porous matrix into the crack, as shown in formula (9):

[0117]

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

[0119]

[0120] Among them, κ c is the mass transfer coefficient, as shown in formula (11):

[0121]

[0122] Among them, D i f is the solute diffusion coefficient in the crack (i = HCl or CaCl2), N Re is the Reynolds number, N Sc is the Schmidt number, a vf =2 / d f , that is, the interfacial area of ​​the crack per unit volume, is the HCl concentration at the fluid-solid interface. For a first-order kinetic reaction (k reaction rate constant), It can be expressed as formula (12):

[0123]

[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 is also adopted, that is, it is assumed that the temperature of water is the same as that of nearby rocks. 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]

[0128] Among them, (ρc p ) eff,m =φ m ρ acid c p,acid +(1-φ m )ρ s c p,s is the effective volume specific heat capacity of the rock matrix, λ eff,m =φλ acid +(1-φ)λ s is the effective thermal conductivity of the rock matrix, ρ is the effective density of the rock matrix, c p is the specific heat at constant pressure, T is the temperature, is the heat source generated by the chemical reaction, as shown in formula (14);

[0129]

[0130] Where H is the reaction enthalpy, which is a negative value because the acidification reaction is an exothermic reaction.

[0131] Heat transfer equation in fracture: as shown in equation (15);

[0132]

[0133] Among them, (ρc p ) eff,m ≈ρ acid cp,acid is the effective volume heat capacity of the crack, λ eff,m ≈λ acid is the effective thermal conductivity of the crack, is the reaction heat source in the crack, as shown in formula (16);

[0134]

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

[0136]

[0137] For rock materials, linear elasticity (Hooke's law) and small deformation theory are applied to derive the equilibrium equation of rock matrix deformation under pure gravity load conditions, and the solid mechanics equation expression is shown in Equation (18);

[0138]

[0139] Where L and G are the Lamé constant, is the volume strain, u is the displacement vector and K v is the bulk modulus;

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

[0141] u u -u d -u0=-F / k A (19)

[0142] Where F represents force / unit area u u ,u d are the displacement vectors of the upper and lower boundaries of the crack layer and the initial compensation displacement u0 respectively.

[0143] The stiffness can be expressed by the normal stiffness k A and tangential stiffness k n To obtain, they represent the deformation in two corresponding directions, as shown in formula (20);

[0144]

[0145] Where n is 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]

[0147] Among them, E f is the elastic modulus, ν f is Poisson's ratio, G f is the shear modulus.

[0148] The acidification reaction leads to an increase in the porosity of the rock matrix and an increase in the fracture opening, thereby increasing the permeability of the carbonate geothermal reservoir. For the rock matrix, the relationship between the porosity change and the reaction can be described by the continuity equation of the rock skeleton, as shown in Equation (22);

[0149]

[0150] in, is the matrix porosity caused by chemical reactions, ω = ζ / M is the number of moles of calcium carbonate (CaCO3) per unit mass of rock, ζ and M represent the content and molar mass of calcium carbonate minerals in the reservoir, respectively;

[0151] In addition, the porosity change due to the thermo-poroelastic effect can be expressed as follows, as shown in Equation (23);

[0152]

[0153] Among them, ε v represents volumetric strain, K s represents the tangential stiffness.

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

[0155]

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

[0157]

[0158] Among them, k m0 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 aperture change rate caused by chemical corrosion, are shown in Equation (26);

[0160]

[0161] Furthermore, the ordinary differential equations and partial differential algebraic equations of the domain, i.e., the change in the crack opening due to mechanical deformation can be calculated by using the parameters of the thin elastic layer as well as the normal strain of the layer and the fluid pressure in the crack, as shown in Eq. (27);

[0162]

[0163] Under the thermo-hydraulic-chemical (THMC) coupling condition, the crack opening is finally obtained by integrating equations (26) and (27), which is equation (28);

[0164]

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

[0166]

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

[0168] Step 4: Assign initial conditions and boundaries, create and optimize the mesh, solve the separation coupling model, and perform acid fracturing evaluation;

[0169] Specifically, the numerical calculation of THMC (thermal-hydraulic-chemical) modeling was performed using the finite element method using COMSOL software, and the influence of discrete fracture networks (DFNs) in acid fracturing was further considered. These fractures are randomly generated, their lengths follow a power law distribution, and their directions follow a Fisher distribution. The "Darcy's Law" physics module is used to calculate the acid flow in the matrix (Formula (1)) and the fracture (Formula (2)). The "Transport of Diluted Species in Porous Media" module is used to calculate the solute concentration in the matrix (Formula (5)) and the fracture (Formula (8)). In order to determine the temperature distribution in the matrix (Formula (13)) and the fracture (Formula (15)), the "Heat Transfer in Porous Media" module is used. In addition, the “Solid Mechanics” module was applied to obtain the deformation of the rock matrix (Eq. (18). Subsequently, in order to obtain the change in matrix porosity due to chemical corrosion (Eq. (22)), the domain ordinary differential equations (ODEs) and differential algebraic equations (DAEs) modules were enabled. For the change in fracture aperture (Eq. (26)), the boundary ODEs and DAEs modules were activated. For modeling, the computational domain was discretized using a free triangular mesh. The calculations started with the default initial step and the variables were updated using the implicit backward difference formula. In order to improve the computational convergence of the proposed model, a segregated coupled solver was used to solve the segregated coupled model, and the coupled THMC model was numerically calculated and sensitivity analyzed, and the acid fracturing effect was evaluated;

[0170] Through the method of the present invention, it is revealed that the structure of the fracture network, the acid injection rate, the acid concentration and other factors have a significant impact on the fracture porosity change, reaction rate and pressure evolution. In particular, the connectivity of the fracture plays a key role in fluid flow and pressure distribution.

[0171] The model of the present invention simulates the chemical corrosion and mechanical deformation of the fracture in detail based on the assumption of thin elastic layer and fracture unit, including the expansion and closure of the fracture under chemical action, and the resulting change in matrix porosity. In addition, the model also considers the dynamic evolution of fracture aperture and matrix porosity, which is crucial for understanding the long-term changes in reservoir permeability. By simulating the effect of acid fracturing under different reservoir conditions and operating parameters, the influence of fracture location and fracture network type on the difference in fracture aperture is revealed. The study found that the injection rate is a key factor affecting the average chemical aperture, while the acid concentration of HCl has a positive effect on improving the average chemical aperture. In addition, there is a negative correlation between the initial aperture of the fracture and the peak average reaction rate, while it is positively correlated with the average chemical aperture. The initial temperature of the reservoir also has a significant effect on the change in chemical reaction rate and fracture aperture. Chemical factors play an important role in fracture extension, while mechanical factors impose significant constraints on the fracture extension process. In contrast, the effect of temperature on fracture extension is relatively small. These findings provide a valuable scientific basis for optimizing acid fracturing in carbonate reservoirs and help guide practical geothermal development work. Through the study of multiple actual cases, the practicability and effectiveness of the proposed model are further verified.

[0172] Example 3

[0173] The acid fracturing evaluation method for carbonate geothermal reservoirs based on the thermal-hydraulic-mechanical-chemical coupling model of this embodiment is specifically as follows:

[0174] Taking the coupled THMC model for numerical calculation and sensitivity analysis as an example, a plane strain square rock with a dimension of 200 m in one direction was used in the study. Figure 2The boundary conditions of 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 calculation steps are taken before injection to achieve an equilibrium initial state of ground stress. For the chemical field, the 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's ratio, thermal expansion coefficient, initial porosity, initial permeability, storage coefficient, thermal conductivity, Poisson's ratio, and specific heat capacity; the engineering parameters of the fracture fissure include: the width of the fracture, the Young's modulus of rupture of the fracture, the storage coefficient of the fracture, and the initial aperture of the fracture.

[0175] Table 1 Parameters used in the simulation

[0176]

[0177]

[0178] According to the engineering parameters, COMSOL software was used to establish a geometric model for simulating carbonate reservoirs. The formation of the fracture network was simulated by COMSOL to obtain the variables in the computational domain and the fractures. The heat transfer field was added to the computational domain, and the coupling relationship between the hydraulic field, stress field and chemical field was established, thus forming a thermal-hydraulic-mechanical-chemical coupling model covering the fracture deformation.

[0179] A fracture network is added to the computational domain. The length of the fracture follows an exponential distribution law, and the inclination of the fracture is considered to be random and uniformly distributed along the fracture. Figure 3 As shown, a discrete fracture network is generated. Figure 3 The discrete fracture network (DFN) shown in the figure is composed of 40 fractures with lengths ranging 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 in the fracture network. The formation of the fracture network is simulated by the interpolation function in COMSOL to obtain the variables in the computational domain and the fractures. Then, the heat transfer field is added to the computational domain, and the coupling relationship between the hydraulic field, the stress field, and the chemical field is established, thereby forming a thermal-hydraulic-mechanical-chemical coupling model covering the fracture deformation (such as Figure 4 Assign initial conditions and boundaries, create and optimize the mesh, solve for the separation coupling, and perform acid fracturing evaluation.

[0180] Example 4

[0181] Figure 5 The distribution of pressure, acid concentration, temperature and chemical aperture width (fracture width created by chemical reactions) at four different moments in DFN, namely 30 minutes, 60 minutes, 90 minutes and 120 minutes, is shown. The results change over time as the acid injection progresses. At 30 minutes, the pressure in most areas of the site exceeded 30MPa, and the highest local pressure reached about 45MPa. As time goes by, the pressure gradually decreases, and the overall value is close to the initial pressure of 10MPa. Regarding the concentration of acid, at 30 minutes of injection, the concentration is higher in the fractures close to the inflow boundary, and can reach 8750mol / m 3 . Subsequently, the acid concentration in the internal cracks connected to the cracks near the inflow area gradually increases with the diffusion of solutes. On the other hand, for the unconnected cracks, the acid concentration changes little. It can be clearly seen from the temperature map 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 the cracks improves the liquid flow in the network, which increases the reaction rate. In this case, the crack corrosion becomes faster and enhances the convection, so the crack temperature drops rapidly. In addition, the chemical crack width gradually increases over time, and some of them can reach about 5 mm. The identification of these cracks with large chemical crack widths shows their basic consistency with the low-temperature cracks observed in the temperature map, which is caused by chemical corrosion. The change in crack width is affected not only by the crack width, but also by factors such as pore pressure and temperature. Among these factors, the effect of temperature is smaller because the temperature exposure time is short and the heat transfer is limited. In addition, the thermal expansion of the matrix is ​​smaller than the deformation of the crack.

[0182] Example 5

[0183] In order to more accurately describe the concentration changes of the reactant HCl and the 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 certain rate for 20 minutes. Figure 6In the DFN, as the reactant HCl advances along the crack and undergoes consumption, its concentration decreases. In the initial stage of DFN, the concentration of the reactant decreases significantly, marked by the rapid appearance of an inflection point, after which the change tends to stabilize. As the acidification reaction becomes more and more complete, the HCl zero concentration point advances forward. 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 start of the acidification reaction, and the decrease is due to the gradual outflow of the product from the crack. At the observed crack, 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, slowly rises to a peak at approximately 2.5 meters, and then slowly decreases. At the same time, the peak value increases with time. As Figure 7 The results show that the overall trend of CaCl2 concentration in DFN is to increase first and then decrease. This increase is attributed to the start of the acidification reaction, and the subsequent decrease is the result of the gradual outflow of products from the fractures. Specifically, along the observed fractures, the CaCl2 concentration gradually increased in the first 20 minutes and then dropped to zero at 40 minutes. Subsequently, in the period from 40 to 120 minutes, the product concentration increased rapidly, reaching a peak of approximately 2.5 meters, and then gradually decreased. At the same time, the peak value showed an upward trend over time.

[0184] Example 6

[0185] The chemical aperture width of the crack is related to the extent of the chemical reaction. Figure 8 and Fig. 9 The spatial and temporal distribution of chemical aperture widths of observed cracks in the crack network is shown. For cracks such as Figure 8 As shown in Figure 2, in DFN, the chemical aperture width generally shows a downward trend and drops to zero at about 30 meters. It is worth noting that the chemical aperture width at the far left end of the fracture increases with time. Fig. 9 As shown in Figure 1, in the DFN, the chemical aperture width increases very slowly within 40 meters. Subsequently, the chemical aperture width grows rapidly near the fracture intersection. As time goes by, the magnitude of the change becomes more obvious, and the peak reaches a higher value.

[0186] The present invention designs a fully coupled THMC model of carbonate reservoir acid fracturing process, which integrates the interaction among heat transfer, acidizing fluid flow, deformation of rock matrix and fractures, migration of reaction dilution substances 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 of carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model, characterized in that: Follow the steps below to implement it: Step 1: Collect the physical and mechanical parameters of carbonate reservoir matrix and engineering parameters of fractures and fissures; Step 2: According to the engineering parameters, a geometric model simulating the carbonate reservoir is established, and the variables in the fractures in each field are assigned to the fracture domain to obtain the calculation domain; Step 3: Add heat transfer field, hydraulic field, stress field and chemical field to the calculation domain respectively, establish the coupling relationship between them, and then form a thermal-hydraulic-mechanical-chemical coupling model covering the crack deformation; Step 4: Solve the separation coupling model and perform acid fracturing evaluation.

2. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model according to 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, thermal expansion coefficient, initial porosity, initial permeability, storage coefficient, thermal conductivity, Poisson's ratio and specific heat capacity; the engineering parameters of the fracture fissure include: crack width, Young's modulus of rupture of the crack, storage coefficient of the crack, and initial pore size of the crack.

3. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model according to claim 2, characterized in that: In step 2, the specific steps of using COMSOL software to establish a geometric model of the carbonate reservoir to be simulated are: based on the engineering geological characteristics of the carbonate reservoir, a computational geometry model is established, the geometric shape and size of the reservoir are defined on the canvas, and the thickness of the cracks generated in the carbonate reservoir is determined according to the engineering parameters; through the function in COMSOL, the variables in the cracks in each field are assigned to the crack domain to obtain the computational domain.

4. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model according to claim 3, characterized in that: In step 3, the thermal-hydraulic-mechanical-chemical coupling model covering crack deformation includes: porous medium flow equation, crack flow equation, transport equation of diluted species in porous medium, transport equation of diluted species in crack, heat transfer equation in porous medium, heat transfer equation in fracture, solid mechanics equation, thin elastic layer equation, boundary ordinary differential equation and partial differential algebraic equation, field ordinary differential equation and partial differential algebraic equation.

5. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model according to claim 4, characterized in that: The flow equation of the acidizing fluid in the porous rock matrix, that is, the porous medium flow equation, is shown in formula (1): The flow equation of the acidizing fluid in the fracture, that is, the fracture flow equation, is shown in equation (2): Among them, ρ acid is the density of the acidizing fluid, S m is the storage coefficient of the rock matrix; S f is the storage coefficient of the crack, p is the pressure, t is the time, v m is the fluid velocity vector of the acidizing fluid in the porous rock matrix, v f is the fluid velocity vector of the acidizing fluid in the fracture, α Bm is the Biot-Willis coefficient of the rock matrix, α mf is the Biot-Willis coefficient of the crack; ε vm is the volumetric strain of the rock matrix, ε vf is the volume strain of the crack, d f is the size of the crack opening, is the gradient operator; Q p is the mass flux from the rock matrix to the fracture, as shown in formula (3): Among them, k m is the permeability of the rock matrix, μ acid is the viscosity of the acidizing fluid, n u is the normal vector pointing to the rock-crack interface toward the rock side, n d is the normal vector pointing to the rock-crack interface toward the crack side; Based on the variable of HCl, its transport in the porous matrix is ​​controlled by equation (5), which is the transport equation for dilute species in porous media; Among them, φ m is the porosity of the rock matrix, C HCl is the concentration of HCl, D e is the effective diffusion coefficient, which can be expressed by D e =φ f D f / τ f The tortuosity of the crack f =φ -13 and diffusion coefficient; in the presence of sufficient solid reactants, the rate of acid-rock chemical reaction As shown in formula (6); Among them, the reaction order K is n a The reaction rate constant is shown in formula (7): Where A is the Arrhenius constant, e is the base of the natural logarithm, E is the activation energy, and R g is the molar gas constant, T is the absolute temperature; The transport equation of diluted species in fractures is shown in equation (8): in, represents the acidification reaction rate at the fluid-solid interface, Q c represents the mass flux from the porous matrix into the crack, as shown in formula (9): It can be determined based on mass transfer, as shown in formula (10): Among them, κ c is the mass transfer coefficient, as shown in formula (11): in, is the solute diffusion coefficient in the crack, N Re is the Reynolds number, N Sc is the Schmidt number, a vf =2 / d f , that is, the interfacial area of ​​the crack per unit volume, is the HCl concentration at the fluid-solid interface; It can be expressed as formula (12):

6. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model according to claim 5, characterized in that: The heat transfer equation in porous media is shown in equation (13): Among them, (ρc p ) eff,m =φ m ρ acid c p,acid +(1-φ m )ρ s c p,s is the effective volume specific heat capacity of the rock matrix, λ eff,m =φλ acid +(1-φ)λ s is the effective thermal conductivity of the rock matrix, ρ is the effective density of the rock matrix, c p is the specific heat at constant pressure, T is the temperature, is the heat source generated by the chemical reaction, as shown in formula (14); Where H is the reaction enthalpy; Heat transfer equation in fracture: as shown in equation (15); Among them, (ρc p ) eff,m ≈ρ acid c p,acid is the effective volume heat capacity of the crack, λ eff,m ≈λ acid is the effective thermal conductivity of the crack, is the reaction heat source in the crack, as shown in formula (16); Q T is the heat source term from the matrix to the fracture, as shown in equation (17); The solid mechanics equation expression is shown in formula (18); Where L and G are the Lamé constant, is the volume strain, u is the displacement vector and K v is the bulk modulus; The thin elastic layer equation is shown in equation (19); u u -u d -u0=-F / k A (19) Where F represents force / unit area u u ,u d are the displacement vectors of the upper and lower boundaries of the fracture layer and u0 initial compensation displacement respectively; The stiffness can be expressed by the normal stiffness k A and tangential stiffness k n To obtain, they represent the deformation in two corresponding directions, as shown in formula (20); Where n is the unit normal vector of the crack; the normal and tangential stiffness of the thin elastic layer can be calculated by its material parameters, as shown in Equation (21); Among them, E f is the elastic modulus, ν f is Poisson's ratio, G f is the shear modulus.

7. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model according to claim 6, characterized in that: For cracks, the boundary ordinary differential equation and partial differential algebraic equation, i.e., the crack aperture change rate caused by chemical corrosion, are shown in Equation (26); The ordinary differential equations and partial differential algebraic equations in the field, i.e., the change in the crack opening due to mechanical deformation can be calculated by using the parameters of the thin elastic layer as well as the normal strain of the layer and the fluid pressure in the crack, as shown in Eq. (27); Under the coupled conditions of thermal-hydraulic-mechanical-chemical, the crack opening is finally obtained by integrating equations (26) and (27), which is equation (28); The permeability of the fracture can be calculated using the parallel plate flow model, as shown in Equation (29); Based on the above, a thermal-hydraulic-mechanical-chemical coupling model of covering crack deformation was established.

8. The method for evaluating acid fracturing of carbonate geothermal reservoirs based on a thermal-hydraulic-mechanical-chemical coupling model according to claim 7, characterized in that: In the step 4, specifically: using a separation coupling solver to solve the separation coupling, performing numerical calculation and sensitivity analysis on the coupling model, and evaluating the acid fracturing effect.

Images