Numerical simulation method for multi-field coupling heat exchange mechanism of deep carbonate rock geothermal reservoir
Through the geological three-dimensional model and microseismic monitoring network combined with seismic event positioning algorithm and finite element method, the heat flow field changes of deep carbonate geothermal reservoirs are simulated, and the problem of inaccurate prediction of thermal flow field changes in the existing technology is solved, and the precise control of the heat flow of deep geothermal reservoirs and the safe and efficient development of geothermal resources are achieved.
Patent Information
- Application Number
- CN202510361225.8
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-03-26
- Publication Date
- 2025-04-25
- Estimated Expiration
- 2045-03-26
AI Technical Summary
The prior art is difficult to effectively predict and control the thermal flow field changes in deep carbonate geothermal reservoirs, resulting in low reliability and practicality of the thermal flow field changes in geothermal mining.
By acquiring geological three-dimensional model and microseismic monitoring network data, seismic event positioning algorithm and finite element method are used to construct a mathematical model of fissure expansion, combined with multi-field coupled models to simulate heat conduction, fluid flow, solid deformation and chemical reaction processes, evaluate the changes in the heat flow field and formulate management strategies.
It improves the real-time monitoring ability of dynamic changes of cracks, accurately simulates the direction of cracks and activity intensity, effectively integrates multiple field coupling effects, improves the control ability of the heat flow of deep carbonate geothermal reservoirs, reduces geological risks, and optimizes the development and utilization of geothermal resources.
Smart Images

Figure CN119882047A_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the field of geothermal development technology, and in particular to a numerical simulation method for multi-field coupled heat exchange mechanism of deep carbonate geothermal reservoir. Background Art
[0002] Deep carbonate geothermal reservoirs have rich heat storage capacity and have become an important target for the development of geothermal resources. However, traditional geothermal development technologies are mainly focused on shallow geothermal systems. Deep carbonate geothermal reservoirs are difficult to develop due to their complex geological conditions and multi-physical field coupling characteristics. Existing technologies often use simplified models to describe the thermal flow behavior of geothermal reservoirs, but because they fail to fully consider the fracture expansion, chemical reactions and multi-field coupling effects in the formation, the model is usually insufficiently accurate and cannot effectively predict and control the changes in the thermal flow field during geothermal mining. This limits the efficient utilization and further development of deep carbonate geothermal resources.
[0003] In the prior art, the publication number is CN115310329A, and the name is a multi-field coupled numerical simulation method for geothermal reservoirs; it includes the following steps: step 1, basic geological research; step 2, heat reservoir characteristics research; step 3, three-dimensional geological model construction; step 4, initial key parameter determination; step 5, mathematical model establishment; step 6, thermal reservoir development and utilization plan optimization. This method is based on the geological structure of the project area and the geological model established by the drilling data, and builds a hydrothermal development heat reservoir multi-field coupling mathematical model to simulate the dynamic response of various parameters under large-scale extraction and irrigation conditions, thereby optimizing the geothermal development extraction and irrigation operation parameters in the project area, and providing guiding suggestions for geothermal development in other regions.
[0004] For deep carbonate geothermal reservoirs, the existing technologies have the following shortcomings: lack of comprehensive consideration of multi-field coupling effects, especially in terms of fracture expansion and heat flow coupling changes, most studies fail to provide accurate simulation methods. In addition, for the dynamic monitoring of fracture expansion and the assessment of geological risks, traditional methods often rely on discrete monitoring data and simple simulation methods, resulting in low reliability and practicality of prediction results; The above information disclosed in the above Background section is only for enhancement of understanding of the background of the present disclosure and therefore it may contain information that does not form the prior art that is already known to one of ordinary skill in the art. Summary of the invention
[0005] The purpose of the present invention is to provide a numerical simulation method for multi-field coupled heat exchange mechanism of deep carbonate geothermal reservoirs to solve the problems raised in the above background technology.
[0006] To achieve the above object, the present invention provides the following technical solutions: A numerical simulation method for multi-field coupled heat transfer mechanism of deep carbonate geothermal reservoirs, comprising the following specific steps: Step S1: obtaining a geological three-dimensional model of a target geothermal area, and marking a number of target fractures in the geological three-dimensional model; Based on several target fractures on the geological three-dimensional model, a microseismic monitoring network is deployed in the target geothermal area to monitor several microseismic events occurring in the target geothermal area during geothermal development in real time, and collect relevant preliminary microseismic event data, including but not limited to time estimation, spatial coordinate estimation point, magnitude estimation value and focal mechanism estimation information of each microseismic event; Step S2: using a seismic event location algorithm to process the collected preliminary microseismic event data, thereby accurately determining the time, spatial coordinates, focal mechanism and magnitude of each microseismic event; Based on the spatial coordinates of each microseismic event, determining the associated impact area of each microseismic event, and determining the target fracture within the associated impact area; Step S3: extracting the seismic phase characteristics of each microseismic event, the seismic phase characteristics including the propagation time, amplitude and frequency of the P wave and the S wave; Based on the seismic phase characteristics of each microseismic event, a mathematical model of crack extension is constructed using the finite element method. This mathematical model is used to analyze the extension direction and activity intensity of several target cracks in the target geothermal area, thereby obtaining the predicted crack size of each target crack. Step S4: obtaining a multi-field coupling model of the target geothermal area, the multi-field coupling model is used to simulate the heat conduction, fluid flow, solid deformation and chemical reaction processes in the deep carbonate geothermal reservoir, and finally obtain the thermal flow field data of each associated influence area in the target geothermal area; Step S5: obtaining the thermal flow field data of each associated influence area in the target geothermal area, and importing the predicted fractures of each target fracture in the associated influence area into the multi-field coupling model for separate simulation, obtaining the thermal flow field change data of each associated influence area, and analyzing the thermal flow field change data of each associated influence area, and calculating the influence coefficient, which is used to evaluate the degree of change of the thermal flow field under different fracture extension situations, and finally obtaining the thermal flow field change evaluation result of each associated influence area; Step S6: Based on the calculated magnitude in each associated impact area and the predicted fracture of each target fracture, and through the thermal flow field change assessment results, corresponding management and adjustment strategies are formulated for the geological risk and thermal flow field change state of each target fracture to effectively control the heat flow of the geothermal reservoir.
[0007] Compared with the prior art, the beneficial effects of the present invention are as follows: by using geological three-dimensional models and microseismic monitoring networks, the real-time monitoring capability of dynamic changes of fractures is improved; by using seismic event location algorithms and finite element methods, the direction of fracture extension and activity intensity can be accurately simulated, and changes in the thermal flow field can be evaluated; this method effectively integrates the multi-field coupling effects of heat conduction, fluid flow, solid deformation and chemical reactions, and greatly improves the ability to control heat flow in deep carbonate geothermal reservoirs; in addition, by calculating influence coefficients and formulating management strategies, the potential geological risks caused by geothermal development can be significantly reduced, and the level of geothermal resource development and utilization can be optimized, which has important scientific value and practical significance for the safe, efficient and economic development of deep carbonate geothermal resources. BRIEF DESCRIPTION OF THE DRAWINGS
[0008] Figure 1 It is a schematic diagram of the overall method flow of the present invention. DETAILED DESCRIPTION
[0009] In order to make the objectives, technical solutions and advantages of the present invention more clearly understood, the present invention is further described in detail below in conjunction with specific embodiments.
[0010] It should be noted that, unless otherwise defined, the technical terms or scientific terms used in the present invention should be understood by people with ordinary skills in the field to which the present invention belongs. The words "first", "second" and similar words used in the present invention do not indicate any order, quantity or importance, but are only used to distinguish different components. "Include" or "comprise" and similar words mean that the elements or objects appearing before the word include the elements or objects listed after the word and their equivalents, without excluding other elements or objects. "Connect" or "connected" and similar words are not limited to physical or mechanical connections, but may include electrical connections, whether direct or indirect. "Up", "down", "left", "right" and the like are only used to indicate relative positional relationships. When the absolute position of the described object changes, the relative positional relationship may also change accordingly.
[0011] Embodiment 1: See also Figure 1 , the present invention provides a technical solution: A numerical simulation method for multi-field coupled heat transfer mechanism of deep carbonate geothermal reservoirs is applied in the development process of deep carbonate geothermal reservoirs, and the deep carbonate geothermal reservoirs are taken as target geothermal areas. The specific steps include: Step S1: obtaining a geological three-dimensional model of a target geothermal area, and marking a number of target fractures in the geological three-dimensional model; Based on several target fractures on the geological three-dimensional model, a microseismic monitoring network is deployed in the target geothermal area to monitor several microseismic events occurring in the target geothermal area during geothermal development in real time, and collect relevant preliminary microseismic event data, including but not limited to time estimation, spatial coordinate estimation point, magnitude estimation value and focal mechanism estimation information of each microseismic event; Step S2: using a seismic event location algorithm to process the collected preliminary microseismic event data, thereby accurately determining the time, spatial coordinates, focal mechanism and magnitude of each microseismic event; Based on the spatial coordinates of each microseismic event, determining the associated impact area of each microseismic event, and determining the target fracture within the associated impact area; Step S3: extracting the seismic phase characteristics of each microseismic event through the spatial coordinates, focal mechanism and magnitude determined for each microseismic event, the seismic phase characteristics including the propagation time, amplitude and frequency of the P wave and the S wave; Based on the seismic phase characteristics of each microseismic event, a mathematical model of crack extension is constructed using the finite element method. This mathematical model is used to analyze the extension direction and activity intensity of several target cracks in the target geothermal area, thereby obtaining the predicted crack size of each target crack. Step S4: obtaining a multi-field coupling model of the target geothermal area, the multi-field coupling model is used to simulate the heat conduction, fluid flow, solid deformation and chemical reaction processes in the deep carbonate geothermal reservoir, and finally obtain the thermal flow field data of each associated influence area in the target geothermal area; Step S5: obtaining the thermal flow field data of each associated influence area in the target geothermal area, and importing the predicted fractures of each target fracture in the associated influence area into the multi-field coupling model for separate simulation, obtaining the thermal flow field change data of each associated influence area, and analyzing the thermal flow field change data of each associated influence area, and calculating the influence coefficient, which is used to evaluate the degree of change of the thermal flow field under different fracture extension situations, and finally obtaining the thermal flow field change evaluation result of each associated influence area; Step S6: Based on the calculated magnitude in each associated impact area and the predicted fracture of each target fracture, and through the thermal flow field change assessment results, corresponding management and adjustment strategies are formulated for the geological risk and thermal flow field change state of each target fracture to effectively control the heat flow of the geothermal reservoir.
[0012] Further explanation: marking several target fractures in the geological three-dimensional model specifically includes: 1.1) Develop minimum requirements based on current heat transfer mechanisms and identify target geothermal areas: 1.11) Obtain the geological exploration data, historical seismic activity data and current geothermal development status of the current area to define the scope of the target geothermal area; the principles for defining the scope of the target geothermal area are as follows: Geophysical exploration: using geological radar, gravity anomaly detection and resistivity imaging to analyze the distribution and boundaries of underground thermal reservoirs; GIS technical analysis: input the collected geological, earthquake, ground temperature and other information into the Geographic Information System (GIS), and determine the target range through the overlay analysis model; 1.12) Ensure that the target geothermal area has good heat storage conditions, including heat flow intensity, geothermal temperature, and that the thermal reservoir meets the current minimum requirements for the development of heat exchange mechanisms.
[0013] 1.2) Construct a 3D geological model of the target geothermal area, including the following: 1.21) Obtain geological data from geological survey reports, drilling data, seismic data and geothermal development data of the target geothermal area; geological data types include but are not limited to stratigraphic information, lithological characteristics, the orientation and distribution of each fracture, and rock mechanical properties; ensure that the spatial resolution of the data meets the following modeling requirements; 1.22) For geological 3D model: Select Petrel or Leapfrog geological modeling software, import the collected geological data into the selected geological modeling software, perform data format conversion and verification, and ensure data integrity; 1.23) Based on the geological data, use the stratigraphic construction tool of the geological modeling software to construct a three-dimensional stratigraphic model to represent the spatial distribution of various types of rocks and rock formations; 1.24) Construct the three-dimensional location of each fracture in the target geothermal area; 1.25) Use historical seismic data to verify whether the fracture distribution in the geological 3D model is consistent with the actual situation; If there is a deviation in the data, the parameters are optimized according to the historical earthquake activities, and the fracture parameters and stratum interfaces in the model are adjusted; the specific implementation method of parameter optimization is a general technical means in the construction of geological three-dimensional models, which will not be described in detail; Determine the three-dimensional position of each fracture in the target geothermal area, as well as the fracture characteristics of each fracture, including fracture size, fracture density, and fracture distribution pattern; Analyze the fracture characteristics of each fracture in the geological three-dimensional model, including fracture size, fracture density, and fracture distribution law; The impact of each fracture characteristic on geothermal development is evaluated. The total score = size score + density score + distribution score. The larger the evaluation value, the greater the corresponding impact. The evaluation strategy is as follows: Each fracture feature is assigned a scoring value, and the impact of the fracture on development is determined based on the actual situation. The scoring rules in this embodiment are as follows: Crack Size Rating: 1 point: The crack size is very small, with a width of less than 1 mm, a depth of less than 2 m, and a length of less than 10 m.
[0014] 2 points: The crack size is small, with a width of less than 1mm, a depth of 2-4m, and a length of 10-20m.
[0015] 3 points: The crack size is medium, with a width of 1-2 mm, a depth of 4-6 m, and a length of 20-30 m.
[0016] 4 points: The crack size is relatively large, with a width of 2-3mm, a depth of 6-8m, and a length of 30-40m.
[0017] 5 points: The crack size is very large, with a width of ≥3 mm, a depth of ≥8 m, and a length of ≥40 m.
[0018] Fracture Density Score: 1 point: The density is very low, and the number of cracks is less than 5 per cubic meter.
[0019] 2 points: The density is low, and the number of cracks is 5-10 per cubic meter.
[0020] 3 points: medium density, with 10-15 cracks per cubic meter.
[0021] 4 points: The density is relatively high, with the number of cracks being 15-20 per cubic meter.
[0022] 5 points: The density is very high, with the number of cracks ≥ 20 / cubic meter.
[0023] Fracture distribution law score: 1 point: The distribution of cracks is very regular, with more than 80% of the cracks having the same direction, and more than 85% of the crack length having uniform spacing.
[0024] 2 points: The distribution of cracks is relatively regular, the directions of the cracks (70%, 80%) are consistent, and the spacing of the crack lengths (70%, 85%) is uniform.
[0025] 3 points: The cracks are distributed in a medium (50%, 70%) area with consistent directions, and the crack lengths (40%, 70%) are evenly spaced.
[0026] 4 points: The distribution of cracks is relatively irregular, the directions of the cracks in the (30%, 50%] region are consistent, and the spacing of the cracks in the (20%, 40%] length is uniform.
[0027] 5 points: The distribution of cracks is very irregular. The directions of cracks in the (0%, 30%] region are consistent, and the spacing of crack lengths (0%, 20%] is uniform.
[0028] According to these evaluation numerical results, these cracks are prioritized in descending order, and according to the size of the total score value, these cracks are divided into high priority cracks, medium priority cracks and low priority cracks, and the cracks that meet the high priority cracks are used as target cracks; High priority fractures: Comprehensive score of 13 points or above. These fractures have significant impacts in size, density and distribution pattern, and have a positive effect on geothermal development.
[0029] Medium priority cracks: The overall score is 9 to 12. These cracks have good local characteristics and are suitable as candidate targets.
[0030] Low priority fissures: Comprehensive scores are 8 or less. These fissures have little impact in all aspects and are generally not preferred. Target fissure selection is based on actual needs, and high priority fissures are selected as target fissures; In the geological three-dimensional model, each target fracture is assigned a number to form {1, 2, ..., i, ..., n}, where i represents the index of the target fracture i and n represents the total number of target fractures; and in the corresponding index of each target fracture, the corresponding fracture feature is marked; A microseismic monitoring network is deployed in the target geothermal area, specifically: The microseismic monitoring network is composed of seismographs and seismic sensors; this embodiment uses three-component seismographs and broadband seismic sensors for deep deployment; Deployment plan: Seismometers and seismic sensors are placed in boreholes at different depths in the target geothermal area, and combined with surface and underground equipment to form a multi-level microseismic monitoring network; and each monitoring location in the microseismic monitoring network is marked as a monitoring point to ensure monitoring of all levels of the fracture system; Selection basis: The densely deployed and high-precision seismic sensor array can improve the spatial resolution of microseismic activity detection and can be applied to environments with complex regional geology and multi-level distribution of fracture systems, thereby more accurately reflecting the changes in fractures; Selection and testing of seismic sensors: Three-component broadband seismic sensor: Select broadband sensors that meet the requirements of geothermal monitoring, such as GeospaceGS-30 or Geosense3C. The reason for choosing such sensors is that they can cover a wider frequency range and ensure sensitivity to small seismic activities. Selection of seismograph: Choose a digital seismograph with high precision to ensure that it can accurately capture the occurrence of microseismic events and can withstand high temperature and high pressure environments; Seismic sensor function test: Pre-install seismic sensors in a laboratory environment and test their frequency response, noise suppression and stability; Test the depth adaptability of seismic sensors to confirm that they can work stably in the high temperature and high pressure environment of geothermal wells; For microseismic monitoring network design: Depth distribution design: Based on the geological three-dimensional model of the target geothermal area, determine the high-frequency occurrence areas of microseismic activities; based on the distribution of fractures, design the layout of seismic sensors at different depths to ensure that microseismic activities at all levels of the fracture system can be monitored; For deep fracture systems, seismic sensors are placed at deeper well locations; for shallow fracture systems, seismic sensors are placed at shallower well locations; It is recommended to deploy at least one set of seismic sensors above, below, and on both sides of each target fracture to ensure efficient coverage; Surface and underground equipment coordination: At least three high-sensitivity seismometers are deployed at the installation location on the surface as a supplement to ensure that data is fully collected in different spatial dimensions; Combine underground seismic sensors with surface monitoring equipment to form a dual monitoring network on the surface and underground; The equipment is installed as follows: For seismic sensor deployment: Use standardized installation tools at the wellhead to install three-component seismic sensors at a specific depth in the well, ensuring that each seismic sensor is stably fixed and connected to the data transmission system; After the seismic sensor is installed, a functional test is performed to confirm the working status of each seismic sensor; Data transmission system construction: Configure a data transmission system to ensure that the underground seismic sensor and the surface data receiving device are connected via a stable communication link; the communication link in this embodiment includes optical fiber and wireless transmission; Install data receiving and processing equipment on the surface to ensure that data from underground seismic sensors can be received in real time; After the microseismic monitoring network is established, the preliminary microseismic event data sent back by the seismic sensors are collected in real time. It should be noted that the time estimation, spatial coordinate estimation point, magnitude estimation value and focal mechanism estimation information collected by the preliminary microseismic event data are all rough measurement results of the seismic sensors; Initial data processing: Use data analysis software to preliminarily process real-time data, remove noise and perform filtering; Use filtering algorithms or Fourier transforms to reduce noise in vibration signals to ensure accurate identification of microseismic events; Real-time monitoring system establishment: A real-time monitoring system is established at the data receiving center, using visualization tools such as charts or 3D geological displays to display the distribution of current microseismic activity.
[0031] Further explanation: The collected preliminary microseismic event data are processed using a seismic event location algorithm to accurately determine the time, spatial coordinates, focal mechanism and magnitude of each microseismic event, including: The preliminary microseismic event data collected in the microseismic monitoring network also includes: timestamp, arrival time difference and seismic wave information, including: The timestamp is the monitoring time point when each microseismic event occurs; The arrival time difference represents the propagation time difference from the earthquake source to each monitoring point; Seismic wave information is composed of the arrival time and amplitude of P and S waves; The above data are preprocessed to remove noise, correct errors, and ensure the accuracy of timestamps and epicenter locations. Signal decoding is then performed to determine the preliminary characteristics of the microseismic event, including the fault type and source intensity corresponding to the microseismic event. Seismic event location algorithm, including: using 3D seismic inversion technology, based on the microseismic monitoring network and the geological 3D model, deriving and correcting the preliminary microseismic event data obtained for each microseismic event, so as to accurately determine the time, spatial coordinates, focal mechanism and magnitude of each microseismic event; The specific steps of 3D seismic inversion technology are as follows: 1. Construction of seismic wave propagation velocity model: Construct a 3D seismic wave propagation velocity model based on the physical properties of different rock layers in the 3D geological model; the physical properties include rock density, elastic modulus, P-wave and S-wave velocities; The three-dimensional seismic wave propagation velocity model is used to simulate the propagation path of seismic waves in underground media to determine the arrival time difference of seismic waves at each monitoring point; 2. Inversion of seismic wave propagation path: An inversion algorithm is used to derive the spatial coordinates of each microseismic event, the spatial coordinates including the epicenter position and focal depth; the inversion algorithm includes the least squares method, genetic algorithm or simulated annealing algorithm; Based on the arrival time difference data and the three-dimensional seismic wave propagation velocity model, the propagation path of the seismic wave is inferred through reverse calculation, and finally the accurate location of each microseismic event is obtained; 3. Focal mechanism analysis: By analyzing the source mechanism of microseismic events, it is possible to determine which type of fault slip caused the microseismic event. Focal mechanisms include normal faults, reverse faults, and translation faults. They also include non-tectonic focal mechanisms, which include volume expansion type, vertical fracture movement type, and liquid induced type. Isotropic Source: Earthquakes caused by underground volume expansion or compression (such as underground reservoir explosion, steam expansion, etc.). Generally, they are not accompanied by obvious fault slip.
[0032] For vertical crack movement (Tensile Crack Opening / Closing): It is a microearthquake caused by the expansion or closure of vertical cracks, such as volcanic earthquakes caused by magma intrusion leading to crack expansion or microearthquakes caused by hydraulic fracturing.
[0033] For fluid-induced seismicity: earthquakes caused by the movement of fissures or faults due to geothermal extraction, oil and gas extraction, waste liquid injection or reservoir induction (such as reservoir earthquakes).
[0034] The fractures produced by normal faults are assumed to be tensile fractures: such fractures are conducive to the flow of thermal fluids and are ideal geothermal reservoir structures.
[0035] The fissures generated by reverse faults are assumed to be compression structures: the permeability of the fissures is poor, which is not conducive to the diffusion of thermal fluids, but there is a certain probability of forming high-pressure areas.
[0036] The direction, sliding mode and stress state of the earthquake source are calculated by using the waveform inversion method or the full-wave inversion technology, so as to determine the earthquake source mechanism. The technical means for determining the earthquake source mechanism is common knowledge and will not be elaborated on. 4. Magnitude calculation: According to the arrival time and amplitude information of the P wave and the S wave corresponding to the seismic wave information, the magnitude of each microseismic event is calculated using a seismological formula; the seismological formula uses the Richter magnitude formula or the moment magnitude formula; the specific calculation formula of the Richter magnitude formula or the moment magnitude formula is the conventional formula used in existing magnitude calculations and will not be described in detail; The Richter scale estimates the size of an earthquake by analyzing the amplitude of seismic waves recorded on a seismograph. The larger the amplitude of the seismic waves sensed by the seismograph, the more energy was released by the earthquake. The Richter scale also takes into account the distance from the epicenter to the monitoring point, because seismic waves gradually decay as they propagate, so the greater the distance, the smaller the amplitude. When calculating the Richter scale, the amplitude is corrected for distance, so that earthquakes in different locations and at different distances can be compared.
[0037] Moment magnitude is an assessment of the size of an earthquake based on the physical properties of the event itself. It takes into account the area of the fault rupture, the distance the fault has slipped, and the mechanical properties of the underground rock. The moment magnitude uses these factors to calculate the total energy released by the earthquake, and therefore more accurately reflects the true intensity of earthquakes of different sizes.
[0038] The calculation of earthquake magnitude requires combining the distance of the monitoring point and the amplitude of the received P and S waves to quantify the energy release of the microearthquake; Derivation and correction of preliminary microseismic event data: After the initial positioning results are obtained through inversion, the time estimate, spatial coordinate estimation point, magnitude estimation value and source mechanism estimation information are corrected in combination with the existing microseismic monitoring network data and geological three-dimensional model; the correction process can also be carried out by introducing correction factors, optimization algorithms or data fusion technology, with the aim of eliminating errors caused by factors such as the seismic wave propagation path and seismograph accuracy, thereby obtaining more accurate microseismic event positioning data. The above correction process adopts a mature application in existing seismological research and microseismic positioning technology, which will not be elaborated on; After obtaining the spatial coordinates of each microseismic event determined by the seismic event location algorithm, these spatial coordinates are mapped to the established geological three-dimensional model; and the microseismic event is associated with the spatial position of the target fracture. The association method includes: The spatial coordinates determined by each microseismic event are obtained, an influence radius RD is set with the spatial coordinates, and a circular area is constructed with the influence radius RD to obtain an associated influence area, and the target crack in the associated influence area is associated with the microseismic event. Specifically: Each microseismic event is marked to form an event sequence set {1, 2, ..., j, ..., M}, where j represents the index of the jth microseismic event and M represents the total number of microseismic events; and the target fracture in the associated impact area corresponding to the jth microseismic event is represented as ; where i1∈{1,2,…,n}; Define the target crack based on the Richter magnitude formula or moment magnitude formula The magnitudes of the corresponding j-th microseismic events are or ; Earthquakes with magnitudes between 2 and 6 are considered small to medium sized: using the Richter scale.
[0039] Earthquakes with a magnitude greater than 6 are considered large earthquakes: the moment magnitude is used because it more accurately reflects the energy and size of the earthquake.
[0040] Using the Richter magnitude formula: If the target rift The associated jth microseismic event is a small to medium-sized earthquake, and the amplitude has been measured near the epicenter. The magnitude of the event is calculated using the Richter magnitude formula: ; in, is the maximum amplitude of the jth microseismic event, is the distance from the epicenter of the microseismic event j to the observation point, is an empirical constant, which is taken as 1.0 in this embodiment; Using the moment magnitude formula: For large earthquakes or when studying deep rift faults, the seismic moment is considered, so the moment magnitude formula is used to calculate the magnitude of the jth microseismic event: ; in, is the seismic moment of the jth microseismic event, is the fault area, is the sliding displacement, is the stiffness of the rock, is the reference seismic moment constant.
[0041] The influence radius RD is determined as follows: The influence radius RD is determined based on the Inversion-Based Radius Adjustment Algorithm, and the specific contents are as follows: The algorithm uses microseismic monitoring data and the spatial distribution of known target fractures, combined with inversion technology, to dynamically optimize the influence radius RD of microseismic events, so as to more accurately associate microseismic events with target fractures. The core of the inversion process is to minimize the spatial error between the associated influence area of the microseismic event and the fracture, so that the influence radius is more matched with the actual underground environment. The main goal of the inversion algorithm is to obtain the optimal influence range that matches the actual geological conditions by iteratively optimizing the influence radius according to the spatial coordinates of different microseismic events, focal mechanisms, and spatial distribution of target fractures. This process includes the setting of the preliminary influence radius, the application of the inversion algorithm, the optimization of the objective function, and the verification and correction of the final results. The algorithm steps are: 1. Initial influence radius setting: Based on the magnitude, focal depth and geological conditions of the microseismic event, a preliminary impact radius is set; the impact radius is preliminarily estimated based on the magnitude, or roughly set based on the geological attributes in the geological 3D model; 2. Establish the objective function: Define an objective function, which is used to measure the matching degree between the associated impact area of the microseismic event and the target fracture. The basic form of the objective function is to optimize the value of the impact radius by calculating the spatial distance error between the microseismic event and the fracture, so that the associated impact area covers the target fracture. 3. Select the inversion algorithm: The inversion algorithm includes: Least squares method: Applicable to simple optimization problems, the optimal solution of the influence radius is directly obtained by minimizing the objective function; Genetic algorithm: Applicable to nonlinear and multi-peak optimization problems, and it can find the global optimal solution by simulating the process of biological evolution; Simulated annealing algorithm: Applicable to complex global optimization problems. It avoids local optimal solutions by simulating the physical annealing process and finds the global optimal influence radius. 4. Perform inversion adjustment: Taking the initial impact radius as the starting point, the objective function is optimized through the inversion algorithm to obtain the optimal impact radius; the inversion process will continuously adjust the value of the impact radius and evaluate the matching degree between the microseismic event impact area and the fracture; 5. Verification and correction: After the inversion adjustment is completed, the optimized influence radius is applied to the geological 3D model to verify the spatial correlation between the impact area of the microseismic event and the target fracture; if the matching degree is not ideal, the inversion adjustment can be performed again to further optimize the influence radius; 6. Feedback mechanism and iterative optimization After each new microseismic event data and fracture information is added, the inversion adjustment process can be re-executed to gradually optimize the influence radius; as more microseismic data and fracture information are accumulated, the adjustment of the influence radius will become more and more accurate.
[0042] Further explanation: A mathematical model of crack propagation is constructed to obtain the predicted cracks of each target crack, including: Fracture properties: Determine the initial fracture size from the fracture characteristics of each target fracture, and determine the physical properties of each target fracture; the initial fracture size includes length, width, and inclination, and the physical properties include elastic modulus and fracture toughness; External conditions: The temperature field, stress field and fluid pressure of the target geothermal area are taken as external conditions and used as driving factors of the mathematical model; and based on the stress field and fluid pressure, the stress state and fluid pressure state of each target fracture are determined; Stress state is defined as: Stress state refers to the distribution of internal and external forces to which the target fracture is subjected under specific conditions. It is expressed in the form of a tensor, including three principal stresses: maximum principal stress, intermediate principal stress, and minimum principal stress. In three-dimensional space, these stresses help describe the degree and direction of the force on the target fracture in different directions.
[0043] Function: Stress state affects the opening, closing and sliding of cracks.
[0044] The fluid pressure state is defined as: The fluid pressure state refers to the pressure state exerted by the fluid in the target fracture, which affects the stability and water conductivity of the fracture.
[0045] Function: Fluid pressure affects the water filling capacity, fluidity and response behavior of the fracture under stress field.
[0046] Fusion phase characteristics: The extracted seismic phase characteristics are combined with the stress state and fluid pressure state of the target fracture as input parameters of the fracture propagation mathematical model to describe the initial and boundary conditions of fracture propagation. Data weight allocation: Through parameter sensitivity analysis and statistical modeling, determine the impact of seismic phase characteristics, stress state and fluid pressure state on crack propagation behavior, assign appropriate weights, and ensure the balance and accuracy of the numerical model; Model solving and optimization: Numerical simulation: The finite element software ABAQUS or COMSOL Multiphysics is used to perform numerical simulation of crack extension. Based on the mathematical model of crack extension, the extension direction and activity intensity of several target cracks in the target geothermal area are analyzed to obtain the predicted cracks of each target crack. The predicted cracks are crack characteristics including crack size, crack density, and crack distribution law. The predicted cracks of the target cracks in the associated impact area corresponding to the jth microseismic event are expressed as The mathematical models of crack extension in this embodiment include but are not limited to: extended finite element method (XFEM, eXtended FiniteElement Method) and nonlinear fracture mechanics (EPFM, Elastic-Plastic Fracture Mechanics); EPFM combines elastic and plastic behaviors and is suitable for handling the plastic zone during fracture propagation. This is particularly important for fracture propagation in deep carbonate rocks under high temperature and high pressure environments, because under these conditions, rocks usually undergo significant plastic deformation.
[0047] XFEM is an efficient numerical simulation method that can simulate the expansion of fractures without re-meshing, and is particularly suitable for complex fracture expansion paths. For fracture expansion in deep carbonate reservoirs, XFEM can handle stress-fluid-heat-deformation coupling problems with great flexibility.
[0048] In terms of parameter optimization: Genetic algorithm or particle swarm optimization algorithm is used to optimize the parameters of the mathematical model of crack extension to minimize the error between the simulation results and the actual observation data; Model Validation: Historical data comparison: Compare the model prediction results with historical microseismic event data to verify the accuracy of the model; Field verification: Through drilling or other geological survey methods, the fracture extension predicted by the model is verified on site to further calibrate the model; The steps for analyzing the crack extension direction and activity intensity are as follows: 1. Analysis of expansion direction: Principal stress direction: Combine geological 3D model and stress field analysis to determine the main expansion direction of the target fracture; Seismic phase feature correlation: using seismic phase propagation time and amplitude changes to identify the potential expansion direction of the crack and adjust the model parameters to reflect the actual situation; 2. Activity intensity assessment: Amplitude and frequency correlation: By analyzing the amplitude and frequency of P and S waves, the intensity of rift activity and energy release can be evaluated; Dynamic monitoring: real-time monitoring of microseismic events, dynamic update of activity intensity assessment, and timely reflection of dynamic changes in fractures; 3. Identification of active and stable areas: Active areas: Based on activity intensity assessment, high activity areas of target fractures are identified as key areas for geothermal development; Stable Zones: Identify low-activity areas of the fracture to serve as areas for auxiliary development and safety monitoring; The division method of active area and stable area is as follows: The intensity of activity is expressed by indicators consisting of displacement rate, stress-deformation response, fluid pressure change, and seismic phase activity / microseismic signals; Displacement rate: The opening of the crack surface or the size of the shear displacement per unit time, which is used to measure the dynamics of crack expansion.
[0049] Stress-deformation response: the stress change rate or energy release rate in the crack corresponding area.
[0050] Stress concentration and rapid release often indicate that the crack is highly active.
[0051] Fluid pressure change: The magnitude of the fluid pressure gradient within the fracture, specifically the diffusion rate of high-pressure fluid, can reflect the conductivity of the fracture and the strength of fluid coupling.
[0052] Seismic phase activity / microseismic signals: Record tiny signal activities through seismic monitoring and analyze the dynamic change intensity and frequency of cracks.
[0053] Data statistics and analysis: By collecting long-term monitoring data, applying statistical analysis and numerical simulation, we determine the specific thresholds for each indicator and use them to divide active areas and stable areas.
[0054] Model simulation and verification: Numerical simulation tools are used to conduct sensitivity analysis on the numerical model to verify the rationality of the threshold and accurately divide different sections.
[0055] Historical data comparison: Combined with historical activity data, comprehensively evaluate the stability and changing trends of fracture activities in different regions to determine regional boundaries.
[0056] It is further explained that the multi-field coupling model of the target geothermal area is determined by using Huang Yonghui's "Multi-field Coupling Simulation Technology and Its Application in Geothermal Exploitation Safety Assessment" or "A Multi-field Coupling Numerical Simulation Method for Geothermal Reservoirs" with publication number CN115310329A; The implementation steps of the multi-field coupling model are as follows: 1. Data preparation and preliminary analysis: Collect basic data on geology, geophysics, and fluid properties of the target geothermal area.
[0057] Conduct preliminary analysis of carbonate rock layers and reservoir characteristics in the region to determine the physical and chemical parameters required for simulation. Reservoir characteristics include but are not limited to porosity, fracture density, and fracture size; Determine model parameters related to thermal, fluid, mechanical and chemical reactions of geothermal reservoirs, such as thermal conductivity, porosity, permeability, reaction rates, etc.
[0058] 2. Define the coupling equations of the model: Heat conduction equation: The heat conduction equation is established based on the thermal conductivity of the rock, the temperature of the fluid, etc. The heat conduction equation is determined based on Fourier's law; Fluid flow equation: The flow model of Darcy's law is used to describe the flow process of geothermal fluid.
[0059] Solid deformation equation: Through the stress-strain relationship, a mechanical model of solid deformation is established, taking into account the plastic, elastic or fracture behavior of rocks during geothermal mining.
[0060] Chemical reaction equation: Based on the reaction process of geothermal fluid and carbonate rock, a chemical reaction kinetic model is established, taking into account the reaction rate and product generation, etc.
[0061] 3. Determine boundary conditions and initial conditions: According to the geological structure and mining history of the region, the boundary conditions of the simulation area are set, including heat flow, temperature, pressure, and fluid injection.
[0062] Determine the initial conditions, including the initial temperature field, pressure field, chemical composition, etc., to simulate the state of the reservoir at the beginning.
[0063] 4. Construct numerical grid and discretization: The simulation area is meshed and the appropriate meshing strategy is selected based on the reservoir geometry, fracture distribution, physical properties and other information.
[0064] Discretize heat conduction, fluid flow, solid deformation and chemical reaction equations using finite element method, finite difference method or finite volume method.
[0065] 5. Multi-field coupling calculation and solution: The numerical solution is performed by coupling the four fields of heat, flow, force and chemical reaction. The time evolution is performed using explicit or implicit time stepping algorithms.
[0066] An iterative solver is used to perform multi-field interactions to account for the interaction between heat flow, fluid flow, and solid deformation.
[0067] 6. Result analysis and verification: Extract simulation results and analyze the distribution of various physical fields such as thermal field, flow field, stress field and chemical reaction field.
[0068] Verify the accuracy and reliability of the model. Compare it with field experimental data, data in existing literature or theoretical models to confirm the rationality of the simulation results.
[0069] The sensitivity of the model is analyzed to study the impact of different parameter changes on the simulation results.
[0070] 7. Application and evaluation of results: Based on the results of the multi-field coupling model, the thermal flow field data of each associated influence area in the target geothermal area are evaluated; The heat flow field data include heat flux density, thermal conductivity, heat reserve, and heat recovery rate. The heat flux density, thermal conductivity, heat reserve, and heat recovery rate of the associated impact area corresponding to the jth microseismic event are defined as , , , .
[0071] Heat flux is the heat flow per unit area, measured in watts per square meter (W / m²), which indicates the intensity of heat transfer. It is calculated from thermal conductivity and temperature gradient: ; in, is the heat flux density, is the thermal conductivity, is the temperature gradient.
[0072] Thermal conductivity is the change in thermal conductivity with changes in temperature, pressure and rock type; the calculation formula is as follows: ; is the thermal conductivity (unit: W / m·K); is the temperature (unit: K); is the pressure (unit: Pa); is the porosity; Thermal Reservoir Capacity is expressed as Heat Capacity: it describes the thermal energy that the reservoir can store; the calculation formula is: ; is the heat reserve of the reservoir (in J, joule); is the density of the rock (in kg / m³); is the specific heat capacity of the rock (unit: J / kg·K); is the volume of the reservoir (in m³); is the temperature change; this formula calculates the temperature change of the rock on the target fracture The amount of stored heat under the given conditions takes into account the mass, specific heat capacity and volume of the rock.
[0073] Heat Recovery Rate: By simulating the mining process, the heat recovery rate under different mining conditions is evaluated. This parameter determines the economy and sustainability of geothermal mining; the calculation formula is: ; is the heat recovery rate (unit: W, Watt); It is the mass flow rate of the fluid (unit: kg / s); it represents the mass of the fluid flowing through the geothermal reservoir per unit time and is obtained by measuring the flow rate.
[0074] is the specific heat capacity of the fluid (unit: J / kg·K), which represents the amount of heat required to change the temperature of a unit mass of fluid; is the temperature change of the fluid (in K).
[0075] Further explanation: the predicted fractures of each target fracture are imported into the multi-field coupling model for simulation, and the heat flow field change data of each associated impact area are obtained. The change values of heat flux density, thermal conductivity, heat reserve and heat recovery rate of the associated impact area corresponding to the jth microseismic event are expressed as , , , ; The heat flux density change value is defined as the change in heat flow transport intensity caused by predicted fracture introduction and geothermal production activities.
[0076] Variation factors: Changes in fractures are predicted to affect the path and intensity of heat flow; The propagation of cracks is predicted to increase or decrease the local heat flux; The calculation formula is ,in and They represent the new heat flux density after the introduction of the predicted crack and the heat flux density in the initial state respectively.
[0077] The thermal conductivity change value is defined to reflect the change in heat transfer capacity within the target geothermal area due to the impact of the predicted fractures.
[0078] Variation factors: Cracks are predicted to increase overall thermal conductivity by providing more paths for heat conduction.
[0079] Conversely, cracks are predicted to reduce local thermal conductivity due to increased thermal resistance; The calculation formula is: ,in and It is the thermal conductivity after crack introduction and in the initial state.
[0080] Definition Thermal reserve change The change in storable thermal energy due to extraction and changes in reservoir conditions.
[0081] Variation factors: Predicted fractures will change the volume and structure of the reservoir, affecting the total heat storage.
[0082] The predicted increase in fractures improves the thermal conductivity of the reservoir, potentially increasing heat storage.
[0083] The calculation formula is ,in and It is to predict the heat storage after fracture introduction and in the initial state.
[0084] Definition: Change in thermal recovery rate is the adjustment of thermal recovery efficiency due to fracture introduction and mining conditions; Thermal recovery rate refers to the ratio of thermal energy recovered from geothermal reservoirs to the thermal energy initially injected or stored during geothermal energy development; Variation factors: Fractures are predicted to change the flow path of the fluid, affecting the effective recovery of thermal energy.
[0085] Changes in heat exchange efficiency can result in an increase or decrease in recovery rate.
[0086] The calculation formula is ,in and It is to predict the heat recovery rate after fracture introduction and in the initial state.
[0087] To further illustrate, the comprehensive contribution of the associated impact area corresponding to the target fracture i1 of the jth microseismic event to the change of the thermal flow field is marked as , and As the influence coefficient, The specific calculation formula is as follows: ; in, is the target fracture number of the jth microseismic event Instead of predicting the crack number The absolute difference is used to indicate the crack number offset; The larger the value, the greater the deviation between the predicted crack and the actual situation; The meaning of the crack number offset is: It is the difference between the target fracture numbers, indicating the difference between the fracture number that actually changes and the predicted fracture number when a microseismic event occurs; this offset reflects the error or deviation between the actual fracture and the predicted fracture; The larger the offset, the greater the difference between the actual target fracture number and the predicted fracture number, which means that the range of fracture extension is more complex, or the deviation between the predicted model and the actual situation is greater.
[0088] The smaller the offset, the smaller the difference between the actual crack and the predicted crack, the simpler the crack expansion, and the more accurate the prediction result.
[0089] is the weight factor, corresponding to the influence weights of heat flux, thermal conductivity, heat reserve and heat recovery rate, with a value range of (0,1) and satisfying ; The contribution of different influencing factors in the target geothermal area experiment is statistically fitted and the weight ratio is adjusted. Specifically, the weight determination method combining the entropy weight method and the fuzzy analytic hierarchy process (FAHP) is used to determine the weight ratio. The weight factor is used to determine the proportion; In order to Conduct detailed qualification analysis and provide quantitative parameters and descriptions based on different scenarios in actual applications, analyzing two typical scenarios: Scenario 1 is a scenario with a large impact; Scenario 2 is a scenario with a small impact; Large impact scenarios: When any parameter among heat flux, thermal conductivity, heat reserve or heat recovery rate has a great influence on the geothermal reservoir, the corresponding weight coefficient will occupy a large proportion.
[0090] Scenarios with small impact: When the impact of any parameter among heat flux, thermal conductivity, heat reserve or heat recovery rate is small, the weight coefficient will be correspondingly smaller.
[0091] Description of scenes with great impact: If, during the development of geothermal reservoirs, the heat flux density has a significant effect on the heat flow field, while the heat recovery rate has a small effect; the weight coefficient is set to ; It means that the influence of heat flux is significant and the influence ratio is large; It means that thermal conductivity has the second greatest influence on the thermal flow field; It indicates that the influence of heat storage is small, but still an important factor; It indicates that the influence of heat recovery rate is minimal; In this scenario, the change of heat flux density has a stronger dominant effect on the change of heat flow field, so it is given a higher weight coefficient , the weights of other factors are relatively low.
[0092] Description of the scenario with small impact: If the heat flux changes slightly, and the change in heat recovery rate plays an important role in reservoir stability; in this case, the weight coefficient is ; It means that the heat flux density has little influence; It indicates that the influence of thermal conductivity is also small; It indicates that heat reserve is more important to reservoir temperature control and is given a higher weight; It indicates that the heat recovery rate has a great influence on the effect of geothermal development; In this scenario, the influence of heat reserves and heat recovery rate on the thermal flow field is more prominent, so the corresponding weight coefficient is higher.
[0093] The weight coefficients in the scenarios with large impact are determined based on the critical role of the actual heat flux changes on the reservoir temperature and energy recovery. In this implementation, through simulation experiments or historical data, it is determined that the change in heat flux has the most significant effect on the thermal management of the reservoir, so a higher weight is given.
[0094] The weight coefficients in scenarios with small impacts are set based on the needs of technical optimization. For example, in the early stages of reservoir development, the heat flux is relatively stable and the impact of heat recovery on the reservoir is more important, so in this scenario, the weight value of heat flux is relatively low.
[0095] The optimization strategies for large-scale scenarios are as follows: When the heat flux density has a greater influence on the thermal flow field, increasing The value of can more effectively evaluate the changes in the thermal flow field and ensure its dominant position in geothermal development.
[0096] Accordingly, reduce , which means that when the heat recovery rate is relatively unimportant, there is no need to rely too much on the effect of the heat recovery rate on the change of the thermal flow field, thereby reducing the computational complexity.
[0097] The optimization strategies for small-scale scenarios are as follows: When the heat flux changes slightly, the and The weights can fully reflect the impact of heat reserves and heat recovery rate and improve the accuracy of the evaluation results, which is very important in optimizing geothermal recovery and reserve management.
[0098] Reduce and The weight of heat flux and thermal conductivity is reduced, which makes the evaluation results focus more on the changes in heat reserves and recovery rate.
[0099] is a nonlinear adjustment factor, corresponding to the nonlinear weight of each change, and its value range is [0.5,2]. The value determination rules are as follows: When the nonlinear adjustment factor is greater than 1: the relative change will be amplified, indicating that the magnitude of the change in the physical quantity is of great importance to the overall influence coefficient.
[0100] When the nonlinear adjustment factor is in the interval (0,1): the relative change will be compressed, indicating that the weight of the change in the physical quantity in the overall influence coefficient is small; Example of technical significance: By adjusting the value of the nonlinear adjustment factor, reasonable weights are given to the relative changes of different physical quantities, reflecting the sensitivity of different physical quantities to the changes in the thermal flow field during the crack expansion process.
[0101] For example, if the heat recovery rate is a key indicator for evaluating the impact of fractures on the thermal flow field, then >1; The effect of thermal conductivity is not as significant as that of heat storage, so it is set <1; The determination method is as follows: According to historical experimental data or actual monitoring data, the change range of heat flux density, thermal conductivity, heat reserve, heat recovery rate and the corresponding heat flow field change results are calculated respectively, and the optimal least square method of nonlinear fitting is used to find value.
[0102] Conduct multiple fracture extension simulation experiments, record the change trend of heat flux density, thermal conductivity, heat reserve, heat recovery rate and corresponding heat flow field influence coefficient in the associated impact area corresponding to the jth microseismic event; use the fitting nonlinear formula Compare the actual impact coefficient results and adjust Select the value that makes the formula fit best.
[0103] is a normalization constant used to standardize the range of the influence coefficient so that The output value range is limited to (0,1); is the distance decay rate of the target crack number difference, which is a dimensionless constant and has a value range of [1,10]; The determination method is as follows: Experimental data analysis: In laboratory or field observations, the heat flow variation between different cracks is measured. Through statistical analysis, the relationship between the crack spacing and the heat flow variation is evaluated to determine the effective distance of the crack influence range.
[0104] By fitting the experimental data, a mathematical model that can describe the variation of heat flow with distance is selected, and appropriate attenuation parameters are extracted from it.
[0105] Field experience and expert opinions: Combined with practical experience in geothermal development and mining engineering, communicated with experts in the field, and determined the reasonable value.
[0106] It means that the influence caused by the offset between the target crack and the predicted crack gradually decays. The larger the distance, The smaller; The calculation logic analysis between the parameters in the formula is as follows: 1.1) The relationship between the calculation logic in the formula and the evaluation results of thermal flow field changes: The following parameter physical quantities are defined including heat flux, thermal conductivity, heat reserve, and heat recovery rate; The relative ratio of each physical quantity change By corresponding weight coefficient With nonlinear factors By combining the above, we can get the comprehensive contribution of each physical quantity to the change of thermal flow field. This part reflects the influence of crack expansion on the thermal flow field, and by adjusting the weights and nonlinear factors, we can flexibly simulate the change of thermal flow under different crack expansion situations. When the change amplitude of a physical quantity is large, such as When it increases, the corresponding relative change value Increases, and thus the influence coefficient will also increase; by weighting the relative changes of various physical quantities and adjusting them through nonlinear factors, the influence of each target crack on the thermal flow field can be accurately evaluated; Crack Number Offset It will decay according to the change of the crack position; the larger the crack number offset, the smaller the value of the exponential decay term, resulting in a smaller influence coefficient. Conversely, the closer the crack position is, the larger the influence coefficient is. When the crack position changes, and As the distance between the cracks increases, the effect of the cracks on the thermal flow field is significantly weakened, and this change is effectively reflected by the attenuation function. This design ensures that the change of the thermal flow field will not be affected by the cracks with excessive position changes, making the model more stable. when When a certain factor of increases, the influence of the relevant physical quantity is strengthened, and the degree of change of the thermal flow field will also increase; this is because the contribution of the physical quantity to the change of the thermal flow field becomes larger, so the influence coefficient increases, which helps to accurately predict the change of thermal flow caused by crack expansion; Increase the nonlinear factor: Increase nonlinear factor The value of can amplify the influence of the corresponding physical quantity change, making the formula more sensitive; by reasonably adjusting the value of the nonlinear factor, the change of the thermal flow field under different crack expansion situations can be more finely adjusted; The increase in the crack number offset leads to a decrease in the attenuation factor and ultimately a decrease in the influence coefficient; this design ensures that when the crack offset is too large, the model's prediction of the thermal flow field changes tends to be stable, thus avoiding overly drastic prediction results; The design goal of the formula is to accurately evaluate the changes in the thermal flow field under different crack extension situations, so its applicable scenarios are mainly based on the following conditions: Fracture propagation in geothermal reservoirs: The formula provides a quantitative assessment of the changes in the heat flow field during fracture propagation and is applicable in geothermal production and reservoir management, especially when predicting the impact of fracture propagation on heat flow; Microseismic events: The formula is used to analyze the impact of crack expansion on the thermal flow field based on the monitoring data of microseismic events. It is applicable to areas with active seismic activity and helps to determine the thermal flow impact of crack expansion.
[0107] like The increase in the value of any factor in the equation indicates that the corresponding heat flux, thermal conductivity, heat reserve, and heat recovery rate are more important, and have an impact on the coefficient The contribution of Will The value range is limited to the interval (0,1), and the evaluation results of the thermal flow field change are as follows: The closer it is to 0, the closer it is to the target crack. The smaller the impact on the change of thermal flow field; The closer it is to 1, the more the target rift is. The greater the impact on the thermal flow field changes.
[0108] When all parameters In the interval (0,0.2), It approaches 0, indicating that the thermal flow field changes little and the reservoir is in a stable state; when In the interval (0.9, 1), it indicates that the thermal flow field changes greatly and the reservoir is in an unstable state.
[0109] Further explanation: According to the geological risk and thermal flow field change status of each target fracture, corresponding management and adjustment strategies are formulated, including: For each target rift The influence coefficient value is combined with the magnitude database to analyze geological risks: like and , then the target crack is determined The corresponding area is of high geological risk; like and , it is judged as low geological risk; in is the risk classification threshold; The risk classification thresholds are determined through historical data analysis, experience accumulation, simulation experiments or model calculations; these thresholds play a key role in risk management, and are used to define what degree of changes or events should be regarded as "risk intervals" and provide a basis for further management measures and adjustment strategies.
[0110] The specific rules for risk classification thresholds are as follows: is the change threshold of the thermal flow field; it indicates the influence of the target crack on the thermal flow field, and and The calculation rules are as follows: The change state of the thermal flow field is expressed by the heat recovery rate, and the following calculation formula is extracted from the thermal flow field change data: The calculation formula is ,in and It is to predict the heat recovery rate after fracture introduction and in the initial state; Classification of thermal flow field effects: According to different The value range divides the thermal flow field change state into different levels: Low impact range: and ; Medium impact range: and ; High impact range: and ; The value is set by simulation results or historical monitoring data of different crack effects. The common value range is: If the impact of the target crack on the thermal flow field is less than 10%, it is considered a low-risk crack area; If the impact of the target crack on the thermal flow field is between 10% and 30%, it is considered a medium-risk crack area.
[0111] If the impact of the target crack on the thermal flow field exceeds 30%, it is considered a high-risk crack area.
[0112] In this embodiment, Set to 0.7; It is the threshold of magnitude change, and the magnitude is related to the expansion of the crack and geological stability; The magnitude classification is as follows: Low magnitude range: ; Moderate magnitude range: ; High magnitude range: ; The value is defined based on the magnitude of microseismic activity and the response characteristics of crack expansion. Through model analysis of the relationship between magnitude change and geological risk, the following thresholds are determined: If the magnitude is less than 2.0, it means that the microseismic activity has little impact on the reservoir and is considered a low-risk area.
[0113] If the magnitude is between 2.0-6.0, it means that the expansion of the cracks may have a certain impact on the thermal flow field and is considered a medium-risk area.
[0114] If the magnitude exceeds 6.0, it means that the cracks will expand on a large scale, accompanied by significant changes in the thermal flow field, and emergency adjustment measures will need to be taken, and it will be considered a high-risk area.
[0115] In this embodiment, Set to 0.3; It is the safety magnitude threshold, which indicates the maximum magnitude value that can be accepted for microseismic activity in different fracture areas. Exceeding this value will cause greater geological risks. Its determination rules are as follows: It means that during a certain crack expansion process, when the micro-earthquake magnitude exceeds this value, it will cause significant damage to the geothermal reservoir and affect the balance of the heat flow field.
[0116] The value of is set based on the seismic resistance of the geological region, historical magnitude data, and model predictions; For example, for relatively stable carbonate formations, based on the systematic experience of the expert group, The value is set to a magnitude of 4.0. If the magnitude is greater than this value, it may cause reservoir structure damage or severe imbalance of thermal flow field.
[0117] For areas with complex geological conditions and easy expansion of cracks, The value can be appropriately lowered based on actual geological exploration results, such as setting it to magnitude 3.0.
[0118] Comprehensive risk classification rules: Management and adjustment strategies for low-risk areas: maintain the existing injection and production plan, continue monitoring, and conduct regular assessments; Management and adjustment strategies for medium-risk areas: adjust injection and production flow rates and increase the frequency of dynamic monitoring.
[0119] Management and adjustment strategies for high-risk areas: Immediately adjust the injection and production strategy, reduce the injection and production pressure, and increase fracture sealing measures.
[0120] According to the distribution of influence coefficients and the changing state of thermal flow field, identify the local high-risk crack areas and draw the thermal flow field risk assessment map; Based on the assessment results of high-risk fracture areas, a thermal flow field optimization adjustment strategy is formulated, specifically: Adjust the injection and production ratio of geothermal fluids in the target area; Optimize injection pressure and temperature gradient in the target fracture direction; Apply physical compensation measures to areas where thermal flow field changes significantly; Continuously monitor the geological risks and thermal flow field changes of each target fracture and adjust strategies dynamically.
[0121] The following is a detailed description. The geological risk judgment rules are as follows: according to and Jointly determine the geological risk level using the following rules: ; Interaction rules and parameter impact analysis: Impact of changes on geological risks: Rules: Increase or decrease The value of directly affects the geological risk level, but does not change the magnitude .
[0122] Specific description: ; Impact of changes on geological risks: Rule: Increase or decrease magnitude , directly affects the risk level, but does not affect .
[0123] Specific description: ; Develop management and thermal flow field adjustment strategies: According to the geological risk assessment results and the changing state of the thermal flow field, quantitative management methods and adjustment strategies are adopted to optimize the distribution of the thermal flow field, improve the thermal energy utilization efficiency of the reservoir and reduce the potential geological disaster risks in high-risk fracture areas.
[0124] Regional priority sorting: The priority of the reservoir area is determined by the geological risk level results, which are divided into high priority, medium priority and low priority: High priority area: the impact of the crack on the thermal flow field is >70%, and the magnitude is ≥6.0; Medium priority area: The impact of the crack on the thermal flow field is 30%-70%, the magnitude ; Low priority areas: fractures have an impact on the thermal flow field of <30% and a magnitude of <2.0.
[0125] Thermal flow field control strategy design, specific suggestions are as follows: For the management and adjustment of high priority areas: Strategy 1: Non-uniform injection and production regulation: Increase the regional thermal fluid injection ratio by 50%, but reduce the recovery ratio to 60% of the current level, while reducing the injected fluid temperature gradient by 25%.
[0126] Technical effect: Enhance the crack closure trend and reduce crack expansion; reduce thermal flow field fluctuations by 30%-40%.
[0127] Strategy 2: Optimizing injection and production pressure by fracture orientation: According to the target crack number According to the position and direction of the pump, adjust the initial injection pressure to reduce it by 20%-30%, and control the pressure fluctuation range to ±5%.
[0128] Technical effect: Significantly reduce the damage of high-priority fractures to local reservoir stability, and reduce heat flow deviation by 15%-20%.
[0129] Management and adjustment of medium priority areas; Strategy 1: Uniform injection and production flow regulation: Adjust the injection and production liquid flow ratio in the target area to 1:1, and control the temperature gradient variation range to <15%.
[0130] Technical effect: Maintain fracture expansion within a controllable range and reduce the impact on reservoir temperature field by 20%-30%.
[0131] Strategy 2: Dynamic crack monitoring and adjustment: For earthquake magnitude In the target fracture area, increase the dynamic monitoring frequency by 30%-40% to reduce excessive injection and production.
[0132] Technical effect: Significantly reduce the fracture shear deformation rate and improve the overall reservoir stability by 15%-20%.
[0133] Management and adjustment of low priority areas: The strategy is to maintain the current injection and production strategy: The injection and production pressure, flow rate, and temperature gradient in the current fracture area remain stable and do not require dynamic adjustment.
[0134] Technical effect: Ensure maximum reservoir efficiency while reducing management costs by 10%-15%.
[0135] The quantitative effects of the interaction rules and adjustment strategies are as follows: Defining interaction rules is to adjust a parameter change in the strategy, and the technical effect on the thermal flow field or risk level must be clear; the following are examples: ; The beneficial effects of this embodiment are as follows: Through the technical solution of this application, the development of deep carbonate geothermal reservoirs will be more precise and stable, realizing the continuous and effective development of geothermal resources, while ensuring the sustainable use of the environment and resources; this solution not only provides a scientific basis and technical support for the development of geothermal resources, but also provides an important reference for the development of deep geothermal resources in other similar areas, and promotes the advancement and application of renewable energy utilization technology.
[0136] The above formulas are all dimensionless and numerical calculations. The formula is a formula for the most recent real situation obtained by collecting a large amount of data and performing software simulation. The preset parameters in the formula are set by technicians in this field according to actual conditions.
[0137] The above embodiments may be implemented in whole or in part by software, hardware, firmware or any other combination thereof. When implemented by software, the above embodiments may be implemented in whole or in part in the form of a computer program product. Those skilled in the art may appreciate that the units and algorithm steps of each example described in conjunction with the embodiments disclosed herein may be implemented by electronic hardware, or a combination of computer software and electronic hardware. Whether these functions are performed by hardware or software methods depends on the specific application and design constraints of the technical solution.
[0138] The units described as separate components may or may not be physically separated, and the components shown as units may or may not be physical units, and may be located in one place or distributed on multiple network units. Some or all of the units may be selected according to actual needs to achieve the purpose of the solution of this embodiment.
[0139] The above description is only a specific implementation manner of the present application, but the protection scope of the present application is not limited thereto. Any technician familiar with the technical field can easily think of changes or substitutions within the technical scope disclosed in the present application, which should be included in the protection scope of the present application.
Claims
1. A numerical simulation method for multi-field coupled heat transfer mechanism of deep carbonate geothermal reservoirs, characterized in that: The specific steps include: Step S1: obtaining a geological three-dimensional model of a target geothermal area, and marking a number of target fractures in the geological three-dimensional model; Based on several target fractures on the geological three-dimensional model, a microseismic monitoring network is deployed in the target geothermal area to monitor several microseismic events occurring in the target geothermal area during geothermal development in real time, and collect relevant preliminary microseismic event data, including but not limited to time estimation, spatial coordinate estimation point, magnitude estimation value and focal mechanism estimation information of each microseismic event; Step S2: Processing the collected preliminary microseismic event data using a seismic event location algorithm to accurately determine the time, spatial coordinates, focal mechanism and magnitude of each microseismic event; Based on the spatial coordinates of each microseismic event, determining the associated impact area of each microseismic event, and determining the target fracture within the associated impact area; Step S3: extracting the seismic phase characteristics of each microseismic event, the seismic phase characteristics including the propagation time, amplitude and frequency of the P wave and the S wave; Based on the seismic phase characteristics of each microseismic event, a mathematical model of crack extension is constructed using the finite element method. This mathematical model is used to analyze the extension direction and activity intensity of several target cracks in the target geothermal area, thereby obtaining the predicted crack size of each target crack. Step S4: obtaining a multi-field coupling model of the target geothermal area, the multi-field coupling model is used to simulate the heat conduction, fluid flow, solid deformation and chemical reaction processes in the deep carbonate geothermal reservoir, and finally obtain the thermal flow field data of each associated influence area in the target geothermal area; Step S5: obtaining the thermal flow field data of each associated influence area in the target geothermal area, and importing the predicted fractures of each target fracture in the associated influence area into the multi-field coupling model for separate simulation, obtaining the thermal flow field change data of each associated influence area, and analyzing the thermal flow field change data of each associated influence area, and calculating the influence coefficient, which is used to evaluate the degree of change of the thermal flow field under different fracture extension situations, and finally obtaining the thermal flow field change evaluation result of each associated influence area; Step S6: Based on the calculated magnitude in each associated impact area and the predicted fracture of each target fracture, and through the thermal flow field change assessment results, corresponding management and adjustment strategies are formulated for the geological risk and thermal flow field change state of each target fracture to effectively control the heat flow of the geothermal reservoir.
2. A method for numerical simulation of multi-field coupled heat exchange mechanism of deep carbonate geothermal reservoir according to claim 1, characterized in that: Mark several target fractures in the geological 3D model, including: Develop minimum requirements based on current heat transfer mechanisms and identify target geothermal areas: Construct a geological 3D model of the target geothermal area; Determine the three-dimensional position of each fracture in the target geothermal area, as well as the fracture characteristics of each fracture, including fracture size, fracture density, and fracture distribution pattern; The impact of each fracture characteristic on geothermal development is evaluated, and the total score = size score + density score + distribution score; the larger the evaluation value, the greater the corresponding impact; According to these evaluation numerical results, these cracks are prioritized in descending order, and according to the size of the total score value, these cracks are divided into high priority cracks, medium priority cracks and low priority cracks, and the cracks that meet the high priority cracks are used as target cracks; In the geological three-dimensional model, each target fracture is assigned a number to form {1, 2, ..., i, ..., n}, where i represents the index of the target fracture and n represents the total number of target fractures; A microseismic monitoring network is deployed in the target geothermal area, specifically: The microseismic monitoring network consists of seismometers and seismic sensors; Seismographs and seismic sensors are placed in boreholes at different depths in the target geothermal area, and combined with surface and underground equipment to form a multi-level microseismic monitoring network; and each monitoring location in the microseismic monitoring network is marked as a monitoring point.
3. A method for numerical simulation of multi-field coupled heat exchange mechanism of deep carbonate geothermal reservoir according to claim 2, characterized in that: The collected preliminary microseismic event data are processed using the earthquake event location algorithm to accurately determine the time, spatial coordinates, focal mechanism and magnitude of each microseismic event, including: The preliminary microseismic event data collected in the microseismic monitoring network also includes: timestamp, arrival time difference and seismic wave information, including: The timestamp is the monitoring time point when each microseismic event occurs; The arrival time difference represents the propagation time difference from the earthquake source to each monitoring point; Seismic wave information is composed of the arrival time and amplitude of P and S waves; Seismic event location algorithm, including: using 3D seismic inversion technology, based on the microseismic monitoring network and the geological 3D model, deriving and correcting the preliminary microseismic event data obtained for each microseismic event, so as to accurately determine the time, spatial coordinates, focal mechanism and magnitude of each microseismic event; Focal mechanisms include normal faults, reverse faults, and translation faults. They also include non-tectonic focal mechanisms, which include volume expansion type, vertical fracture movement type, and liquid induced type. According to the arrival time and amplitude information of the P wave and S wave corresponding to the seismic wave information, the magnitude of each microseismic event is calculated using the seismological formula; the seismological formula uses the Richter magnitude formula or the moment magnitude formula; After obtaining the spatial coordinates of each microseismic event determined by the seismic event location algorithm, these spatial coordinates are mapped to the established geological three-dimensional model; and the microseismic event is associated with the spatial position of the target fracture; the association method includes: The spatial coordinates determined by each microseismic event are obtained, an influence radius RD is set with the spatial coordinates, and a circular area is constructed with the influence radius RD to obtain an associated influence area, and the target crack in the associated influence area is associated with the microseismic event. Specifically: Each microseismic event is marked to form an event sequence set {1, 2, ..., j, ..., M}, where j represents the index of the jth microseismic event and M represents the total number of microseismic events; and the target fracture in the associated impact area corresponding to the jth microseismic event is represented as ; where i1∈{1,2,…,n}; Based on the earthquake severity, choose the Richter magnitude formula or the moment magnitude formula to calculate the magnitude. When calculating the magnitude based on the Richter magnitude formula, define the target fracture The magnitude of the corresponding j-th microseismic event is , when calculating the magnitude based on the moment magnitude formula, define the target fracture The magnitude of the corresponding j-th microseismic event is .
4. A method for numerical simulation of multi-field coupled heat exchange mechanism of deep carbonate geothermal reservoir according to claim 3, characterized in that: A mathematical model of crack propagation is constructed to obtain the predicted cracks of each target crack, including: determining an initial fracture size from the fracture characteristics of each target fracture, and determining a physical property of each target fracture; The temperature field, stress field and fluid pressure of the target geothermal area are taken as external conditions, and the external conditions are taken as driving factors of the mathematical model; and based on the stress field and fluid pressure, the stress state and fluid pressure state of each target fracture are determined.
5. A method for numerical simulation of multi-field coupled heat exchange mechanism of deep carbonate geothermal reservoir according to claim 4, characterized in that: The extracted seismic phase characteristics are combined with the stress state and fluid pressure state of the target fracture as input parameters of the fracture propagation mathematical model to describe the initial and boundary conditions of fracture propagation. Finite element software ABAQUS or COMSOL Multiphysics is used to perform numerical simulation of crack extension. Based on the mathematical model of crack extension, the extension direction and activity intensity of several target cracks in the target geothermal area are analyzed to obtain the predicted cracks of each target crack. The predicted cracks are crack characteristics including crack size, crack density, and crack distribution law. The predicted cracks of the target cracks in the associated impact area corresponding to the jth microseismic event are expressed as .
6. A method for numerical simulation of multi-field coupled heat exchange mechanism of deep carbonate geothermal reservoir according to claim 5, characterized in that: Based on the results of the multi-field coupling model, the thermal flow field data of each associated influence area in the target geothermal area are evaluated; The heat flow field data include heat flux density, thermal conductivity, heat reserve, and heat recovery rate; the heat flux density, thermal conductivity, heat reserve, and heat recovery rate of the associated impact area corresponding to the jth microseismic event are defined as , , , .
7. A method for numerical simulation of multi-field coupled heat exchange mechanism of deep carbonate geothermal reservoir according to claim 6, characterized in that: The predicted fractures of each target fracture are imported into the multi-field coupling model for simulation to obtain the heat flow field change data of each associated impact area. The change values of heat flux density, thermal conductivity, heat reserve and heat recovery rate of the associated impact area corresponding to the jth microseismic event are expressed as , , , ; The comprehensive contribution of the associated impact area corresponding to the target crack i1 of the jth microseismic event to the change of the thermal flow field is marked as , and As the influence coefficient, The specific calculation formula is as follows: ; in, is the target fracture number of the jth microseismic event Instead of predicting the crack number The absolute difference is used to indicate the crack number offset; The larger the value, the greater the deviation between the predicted crack and the actual situation; is the weight factor, corresponding to the influence weights of heat flux, thermal conductivity, heat reserve and heat recovery rate, with a value range of (0,1) and satisfying ; is the nonlinear adjustment factor, corresponding to the nonlinear weight of each change, and its value range is [0.5,2]; is a normalization constant used to standardize the range of the influence coefficient so that The output value range is limited to (0,1); is the distance decay rate of the target crack number difference, which is a dimensionless constant and has a value range of [1,10]; like The increase in the value of any factor in the equation indicates that the corresponding heat flux, thermal conductivity, heat reserve, and heat recovery rate are more important, and have an impact on the coefficient The contribution of Will The value range is limited to the interval (0,1), and the evaluation results of the thermal flow field change are as follows: The closer it is to 0, the closer it is to the target crack. The smaller the impact on the change of thermal flow field; The closer it is to 1, the more the target rift is. The greater the impact on the thermal flow field changes.
8. A method for numerical simulation of multi-field coupled heat exchange mechanism of deep carbonate geothermal reservoir according to claim 7, characterized in that: According to the geological risk and thermal flow field change status of each target fracture, corresponding management and adjustment strategies are formulated, including: For each target rift The influence coefficient value is combined with the magnitude database to analyze geological risks: like and , then the target crack is determined The corresponding area is of high geological risk; like and , it is judged as low geological risk; in is the risk classification threshold; According to the distribution of influence coefficients and the changing state of thermal flow field, identify the local high-risk crack areas and draw the thermal flow field risk assessment map; Based on the assessment results of high-risk fracture areas, a thermal flow field optimization adjustment strategy is formulated, specifically: Adjust the injection and production ratio of geothermal fluids in the target area; Optimize injection pressure and temperature gradient in the target fracture direction; Apply physical compensation measures to areas where thermal flow field changes significantly; Continuously monitor the geological risks and thermal flow field changes of each target fracture and adjust strategies dynamically.
Citation Information
Patent Citations
Multifunctional multi-field coupled seepage experiment device and testing method
CN107462508A
Geothermal reservoir multi-field coupling numerical simulation method
CN115310329A
Safe and efficient development simulation system for myriameter deep geothermal resources
CN116163721A
Complex anisotropic mechanical property conglomerate oil reservoir fracture propagation simulation method
CN117892565A
Geothermal resource exploration method based on multi-model combined decision
CN118378933A
Cited By
Layer fracture dynamic monitoring method and system based on three-dimensional resistivity imaging
CN120652554A
Geomechanical constraint-based geothermal reservoir fracture network dynamic modeling method
CN120781708A
A geomechanically-constrained dynamic modeling method for geothermal reservoir fracture networks
CN120781708B