A hydrate dynamic evolution simulation method based on temperature transfer hysteresis effect and related equipment

Abstract

The application discloses a hydrate dynamic evolution simulation method based on temperature transmission lag effect and related equipment, breaks through the traditional steady-state assumption, quantitatively depicts the control mechanism of the temperature transmission lag on the stable domain evolution, can accurately reveal the non-equilibrium nature of BSR, dynamically reproduces the hydrate system evolution path in the geological history period, can significantly improve the resource quantity estimation precision, the geological disaster prediction reliability and the safety of the drilling engineering, and can be widely applied to the technical field of data processing.

Description

Technical Field

[0001] This invention relates to the field of data processing technology, and in particular to a method and related equipment for simulating the dynamic evolution of hydrates based on the temperature transfer hysteresis effect. Background Technology

[0002] Natural gas hydrates are an important new type of clean energy, and their dynamic evolution is closely related to the global carbon cycle, climate change, and seafloor geological hazards. Bottom-like reflectors (BSRs) on seismic profiles are widely regarded as key indicators for identifying hydrates. Traditionally, BSRs are considered to correspond to the bottom boundary of the hydrate stability domain, i.e., the phase boundary between hydrates and free gas under thermodynamic equilibrium conditions. Under this steady-state equilibrium assumption, numerous studies have directly used BSR depth to invert seafloor geothermal gradients, estimate seafloor heat flows, and assess hydrate resource quantities.

[0003] However, hydrate systems are not always in a steady state during actual geological processes. Temperature and pressure are the core factors controlling the stability domain. Pressure changes can be transmitted to the reservoir instantaneously, while temperature transmission exhibits a significant lag effect due to the limitations of sediment thermal conductivity and heat diffusion rate, with lag times ranging from thousands to tens of thousands of years. Since the Last Glacial Maximum—approximately 18,000 years ago (18 ka BP)—rapid sea-level rise has caused an instantaneous response in seafloor pressure, while the heat from the warming seawater has been slowly conducted downwards. This results in a significant asynchrony in the adjustment of the hydrate stability domain, and the current BSR may not correspond to the system's equilibrium boundary.

[0004] Currently, most related technologies are based on steady-state assumptions, directly using BSR (Body Temperature Scale) for geothermal gradient inversion and heat flow calculation under non-equilibrium conditions, which easily leads to significant biases. Resource calculations all use equilibrium models, ignoring the dynamic processes of hydrate decomposition and regeneration, resulting in inaccurate reserve assessments. Furthermore, there is a lack of evolutionary predictions from a non-equilibrium perspective for geological hazards such as formation instability and submarine landslides induced by dynamic adjustments in the stable domain. The reservoir locations inferred from BSR in drilling projects deviate from the actual bottom boundary of the stable domain, increasing engineering risks. Overall, existing methods struggle to reveal the true evolutionary patterns of hydrate systems under non-equilibrium conditions, hindering accurate resource evaluation and engineering safety assurance. Summary of the Invention

[0005] The main objective of this invention is to propose a method, apparatus, electronic device, storage medium, and program product for simulating the dynamic evolution of hydrates based on the temperature transfer hysteresis effect, aiming to solve at least one problem in the prior art.

[0006] To achieve the above objectives, one aspect of this invention proposes a method for simulating the dynamic evolution of hydrates based on the temperature transfer hysteresis effect, the method comprising: Acquire basic survey data for the target study area, and extract key thermodynamic and geological parameters through standardized processing; Seismic data of the target study area were acquired, seafloor reflectors were identified and their actual burial depth below the seabed was quantified, and a depth dataset of seafloor reflectors was obtained. Based on key thermodynamic parameters and pre-defined boundary conditions for different geological periods, a transient heat conduction equation is constructed by coupling the temperature transfer hysteresis effect, and the transient stratum temperature field for different geological periods is output. Based on key geological parameters, a hydrate phase equilibrium equation is established. The hydrate phase equilibrium equation is then coupled with the transient formation temperature field time-sequentially to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain at different geological histories. The theoretical depth trajectory includes the depth result, adjustment direction, and adjustment rate of the bottom boundary of the hydrate stability domain. Identify depth offsets and assess imbalance states by comparing depth datasets with depth results; Based on the adjustment direction, adjustment rate, and non-equilibrium state predictions, the spatial occurrence and transformation mode of hydrates and free gas, as well as the engineering geological risk level, are obtained.

[0007] In some embodiments, basic survey data of the target study area are acquired, and key thermodynamic parameters and key geological parameters are extracted through standardization processing, including the following steps: Acquire multimodal basic survey data for the target study area; including multichannel seismic data, well logging data, sediment samples, in-situ seafloor temperature data, water depth data, sediment pore fluid salinity data, and hydrocarbon gas composition data; The basic survey data is preprocessed to obtain preprocessed data; the preprocessing includes format standardization, depth realignment, time-depth conversion calibration, and outlier removal. Based on the preprocessed data, the key thermodynamic parameters of the sediments are quantified using measured methods or empirical formulas. Among them, key thermodynamic parameters include sediment thermal conductivity, sediment density and specific heat capacity. The empirical formula method includes using the seismic layer velocity corresponding to multichannel seismic data or the porosity corresponding to the sediment sample to invert the sediment thermal conductivity through the Woodside equation or geometric mean model, and using the seismic velocity corresponding to multichannel seismic data to derive the sediment density and specific heat capacity through the Gardner equation and solid-liquid mixing law. Key geological parameters were determined based on the hydrate decomposition gas components corresponding to hydrocarbon gas components in the preprocessed data and the formation pore water salinity information corresponding to sediment pore fluid salinity data.

[0008] In some embodiments, identifying the seabed-like reflective layer and quantifying its true burial depth below the seabed includes the following steps: Based on seismic data, determine the seismic profile and seismic layer velocity model for the target study area; Based on seismic profiles, quasi-seabed reflective layers are identified by amplitude intensity, negative polarity, and reflective layers parallel to the seabed. Read the two-way reflection time of the seabed-like reflector layer from the seismic data; Using the seismic layer velocity model of the target study area, the two-way reflection time is converted into the actual burial depth below the seabed.

[0009] In some embodiments, key thermodynamic parameters include sediment thermal conductivity, sediment density, and specific heat capacity. Based on these key thermodynamic parameters and the boundary conditions of a predetermined geological period, a transient heat conduction equation is constructed by coupling the temperature transfer hysteresis effect, outputting the transient formation temperature field for different geological periods. This includes the following steps: The thermal diffusivity is obtained by using the sediment thermal conductivity as the numerator and the product of sediment density and specific heat capacity as the denominator. Obtain the boundary conditions for geological history periods; among which, the boundary conditions include the initial seafloor temperature at the last glacial maximum, the geothermal gradient, and the total change in seafloor temperature since the last glacial maximum. Based on the thermal diffusivity, the temperature transfer hysteresis effect is quantitatively characterized by the complementary error function, and then the analytical solution model of the one-dimensional transient heat conduction equation is constructed by combining the boundary conditions. The expression for the analytical solution model is: , In the formula, Indicates burial depth and time nodes The corresponding in-situ formation temperature, Indicates the initial temperature of the seabed. Represents the geothermal gradient. This represents the total change in seabed temperature. Represents the complementary error function. Indicates the thermal diffusivity; By solving the analytical model at preset time intervals, the function of formation temperature variation with depth at multiple time points since the Last Glacial Maximum is reconstructed, and the transient formation temperature field at different geological periods is obtained.

[0010] In some embodiments, key geological parameters include hydrate decomposition gas components and formation pore water salinity information. Based on these key geological parameters, a hydrate phase equilibrium equation is established. This equation is then coupled temporally with the transient formation temperature field to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain at different geological histories. This process includes the following steps: A hydrate phase equilibrium equation was established based on information on hydrate decomposition gas components and formation pore water salinity. Transient formation temperature fields from various geological periods are collected at preset time intervals. The corresponding transient formation temperature fields are intersected with the hydrate phase equilibrium equations to obtain the depth of the bottom boundary of the hydrate stability domain for each geological period. Based on the depth results at each time point in the time series, the adjustment direction and adjustment rate of the bottom boundary of the hydrate stability domain were statistically obtained.

[0011] In some embodiments, identifying depth offsets and evaluating imbalance states by comparing depth datasets with depth results includes the following steps: The depth dataset is spatially overlaid and compared with current depth results, and the depth difference between the two is quantified as the depth offset. When the depth offset is greater than or equal to the first offset threshold, the non-equilibrium state is determined to be an intense non-equilibrium state in the early stage. When the depth offset is less than or equal to the second offset threshold, the non-equilibrium state is determined to be the thermodynamic equilibrium state in the late stage. When the depth offset is between the first offset threshold and the second offset threshold, the unbalanced state is determined to be the dynamic adjustment unbalanced state in the intermediate stage.

[0012] In some embodiments, the spatial occurrence and transformation mode and engineering geological risk level of hydrates and free gas are predicted based on adjustment direction, adjustment rate and non-equilibrium state, including the following steps: Based on the adjustment direction of the bottom boundary of the hydrate stability domain, the spatial transformation relationship between hydrate and free gas is defined; specifically, the following steps are included: if the adjustment direction is upward migration, the area between the original bottom boundary and the new bottom boundary is indicated as the hydrate decomposition zone and the free gas generation zone; if the adjustment direction is downward migration, the newly covered area is indicated as the hydrate generation zone and the free gas consumption zone. By coupling the adjustment rate with the non-equilibrium state, the engineering geological risk level is qualitatively classified. The non-equilibrium state includes intense non-equilibrium, dynamically adjusting non-equilibrium, and thermodynamic equilibrium. The classification of engineering geological risk levels includes: when in intense non-equilibrium and the adjustment rate is greater than the first rate threshold, it is judged as extremely high risk, and the drilling is predicted to encounter a large-scale overpressure free gas layer generated by rapid decomposition; when in thermodynamic equilibrium and the adjustment rate is less than the second rate threshold, it is judged as low risk, and the predicted seafloor reflector basically coincides with the bottom boundary of the actual stable domain; when in thermodynamic equilibrium and the adjustment rate is between the first rate threshold and the second rate threshold, it is judged as medium risk, and there is a dynamic transformation zone between the predicted seafloor reflector and the bottom boundary of the actual stable domain.

[0013] To achieve the above objectives, another aspect of the present invention proposes a hydrate dynamic evolution simulation device based on the temperature transfer hysteresis effect, the device comprising: The first module is used to acquire basic survey data of the target study area and extract key thermodynamic parameters and key geological parameters through standardization processing. The second module is used to acquire seismic data of the target study area, identify the seafloor reflector and quantify its true burial depth below the seafloor, and obtain the depth dataset of the seafloor reflector. The third module is used to construct transient heat conduction equations based on key thermodynamic parameters and pre-set boundary conditions for geological history periods, and output transient stratum temperature fields for different geological history periods by coupling temperature transfer hysteresis effects. The fourth module is used to establish hydrate phase equilibrium equations based on key geological parameters, and to couple the hydrate phase equilibrium equations with the transient formation temperature field time-sequentially to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain at different geological histories. The theoretical depth trajectory includes the depth result of the bottom boundary of the hydrate stability domain, the adjustment direction, and the adjustment rate. The fifth module is used to identify depth offsets and evaluate non-equilibrium states by comparing the depth dataset with the depth results. The sixth module is used to predict the spatial occurrence and transformation mode of hydrates and free gas and the engineering geological risk level based on the adjustment direction, adjustment rate and non-equilibrium state.

[0014] To achieve the above objectives, another aspect of the present invention provides an electronic device, which includes a memory and a processor. The memory stores a computer program, and the processor executes the computer program to implement the aforementioned method.

[0015] To achieve the above objectives, another aspect of the present invention provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the aforementioned method.

[0016] To achieve the above objectives, another aspect of the present invention provides a computer program product, including a computer program that, when executed by a processor, implements the aforementioned method.

[0017] The embodiments of this invention include at least the following beneficial effects: This invention provides a method, apparatus, electronic device, storage medium, and program product for simulating the dynamic evolution of hydrates based on the temperature transfer hysteresis effect. This scheme acquires basic survey data of the target study area, extracts key thermodynamic and geological parameters through standardization processing; acquires seismic data of the target study area, identifies seafloor-like reflectors and quantifies their actual burial depth below the seafloor, obtaining a depth dataset of the seafloor-like reflectors; based on the key thermodynamic parameters and the boundary conditions of preset geological periods, a transient heat conduction equation is constructed by coupling the temperature transfer hysteresis effect, outputting transient stratigraphic temperature fields for different geological periods; a hydrate phase equilibrium equation is established based on the key geological parameters, and the hydrate phase equilibrium equation is coupled with the transient stratigraphic temperature field time-sequentially to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain for different geological periods; wherein, the theoretical depth trajectory includes the depth result, adjustment direction, and adjustment rate of the bottom boundary of the hydrate stability domain; the depth offset is identified and the non-equilibrium state is assessed by comparing the depth dataset with the depth result; based on the adjustment direction, adjustment rate, and non-equilibrium state, the spatial occurrence and transformation mode of hydrates and free gas and the engineering geological risk level are predicted. This invention breaks through the traditional steady-state assumption, quantitatively characterizes the control mechanism of temperature transfer hysteresis on the evolution of the steady-state domain, can accurately reveal the non-equilibrium nature of BSR, dynamically reproduce the evolution path of hydrate systems in geological history, and can significantly improve the accuracy of resource estimation, the reliability of geological hazard prediction, and the safety of drilling projects. Attached Figure Description

[0018] Figure 1 This is a schematic diagram of an implementation environment for a method for simulating the dynamic evolution of hydrates based on the temperature transfer hysteresis effect, provided in an embodiment of the present invention. Figure 2 This is a schematic flowchart of the hydrate dynamic evolution simulation method based on the temperature transfer hysteresis effect provided in the embodiments of the present invention; Figure 3 This is a schematic diagram of the overall process of the hydrate dynamic evolution simulation method based on the temperature transfer hysteresis effect provided in the embodiments of the present invention; Figure 4 This is a schematic diagram illustrating an example of the evolution of the seafloor temperature field during geological history constrained by the temperature transfer hysteresis effect, as provided in an embodiment of the present invention. Figure 5 This is a schematic diagram illustrating an example of the dynamic evolution of the hydrate stability domain constrained by the temperature transfer hysteresis effect provided in an embodiment of the present invention; Figure 6 This is a schematic diagram illustrating a comprehensive comparison example of earthquake identification BSR and dynamically evolving hydrate stability domain bottom boundary provided in an embodiment of the present invention; Figure 7This is a schematic diagram illustrating an example of the dynamic evolution process of a hydrate-free gas system under the constraint of temperature transfer hysteresis provided in an embodiment of the present invention. Figure 8 This is a schematic diagram of the hydrate dynamic evolution simulation device based on the temperature transfer hysteresis effect provided in the embodiment of the present invention; Figure 9 This is a schematic diagram of the structure of the electronic device provided in an embodiment of the present invention. Detailed Implementation

[0019] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the invention and are not intended to limit the invention. In the following description, when referring to the accompanying drawings, unless otherwise indicated, the same numbers in different drawings represent the same or similar elements. The embodiments described in the following exemplary embodiments do not represent all embodiments consistent with the embodiments of this invention; they are merely examples of apparatuses and methods consistent with some aspects of the embodiments of this invention as detailed in the appended claims.

[0020] It is understood that the terms "first," "second," etc., used in this invention may be used to describe various concepts, but unless specifically stated otherwise, these concepts are not limited by these terms. These terms are only used to distinguish one concept from another. For example, without departing from the scope of embodiments of this invention, first information may also be referred to as second information, and similarly, second information may also be referred to as first information. Depending on the context, the words "if" or "when" as used herein may be interpreted as "when," "in response to determination," or "in the event of a determination."

[0021] The terms “at least one,” “multiple,” “each,” “any,” etc., used in this invention, “at least one” includes one, two, or more than two; “multiple” includes two or more than two; “each” refers to each of the corresponding multiple; and “any” refers to any one of the multiple.

[0022] Unless otherwise defined, all technical and scientific terms used in this invention have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains. The terminology used in this invention is for descriptive purposes only and is not intended to limit the invention.

[0023] To facilitate understanding of the technical solution of this invention, the technical terms of the proprietary technical features that may be involved in the technical solution of this invention will first be explained: Natural gas hydrates are ice-like crystalline compounds formed by natural gas (mainly methane gas) and water molecules under low temperature and high pressure conditions on the seabed. They are a huge unconventional clean energy source and are widely found in seabed sediments and permafrost areas along continental margins.

[0024] Hydrate phase equilibrium: used to describe the thermodynamic equilibrium relationship when multiple phases such as hydrate phase, free gas phase, and water phase coexist, and can accurately determine the critical temperature and pressure conditions for hydrate formation and decomposition.

[0025] Hydrate stability region: This refers to the spatial region where natural gas hydrates can exist stably when conditions such as temperature, pressure, fluid composition, and formation medium meet thermodynamic stability requirements. Its range is determined by the hydrate phase equilibrium conditions and is an important basis for judging the occurrence, distribution, accumulation potential, and resource quantity of hydrates.

[0026] Base of Gas Hydrate Stability Zone (BGHSZ): This refers to the maximum burial depth at which hydrates in sediments can stably exist while meeting the temperature and pressure equilibrium conditions. Above this interface, hydrates can exist stably, while below the interface, hydrates decompose into free gas.

[0027] BSR (Bottom Simulating Reflector) is a high-amplitude seismic reflection interface in marine seismic profiles that extends approximately parallel to the seabed, has a reflection polarity opposite to the seabed, and can obliquely intersect sedimentary bedding. This interface corresponds to the bottom boundary of the natural gas hydrate stability domain and is formed by the difference in wave impedance between the overlying hydrate-bearing sediments and the underlying free gas-bearing sediments. It is a key geophysical marker for identifying marine natural gas hydrates.

[0028] Temperature transfer lag effect: refers to the time lag between heat transfer and temperature field response, which causes the actual temperature change in the formation to lag behind the thermal disturbance input, making the hydrate phase equilibrium conditions asynchronous with the real-time temperature and pressure state. This is an important factor that triggers the non-equilibrium phase transformation of hydrates and thus induces the dynamic decomposition of hydrates.

[0029] Geothermal gradient: refers to the rate of temperature change with depth in the strata. It is a key thermodynamic parameter that controls the thickness of the hydrate stability domain, the depth of the top and bottom interface and its spatial distribution. It has a significant controlling effect on the formation, occurrence state and dynamic evolution of hydrates.

[0030] Steady state: refers to a stable state in which key physicochemical parameters such as temperature, pressure, and component concentration of a hydrate system do not change significantly over a relatively long time scale, and the rates of hydrate formation and decomposition tend to be in dynamic equilibrium.

[0031] Transient: refers to an unstable state process in a hydrate system where the temperature, pressure, composition, and phase change rapidly over time after geological activity, thermal disturbance, or abrupt changes in external conditions, before reaching thermodynamic equilibrium. This process is often accompanied by the rapid formation or decomposition of hydrates, reflecting the system's instantaneous response to external disturbances.

[0032] Porosity is a physical parameter characterizing the proportion of pore volume to total volume within sediments, reflecting the formation's reservoir and fluid transport capacity. In natural gas hydrate systems, porosity directly controls pore fluid distribution, hydrate formation space, and heat and mass transfer efficiency, making it a key reservoir parameter affecting hydrate formation, occurrence state, and dynamic evolution.

[0033] Thermal conductivity is a thermodynamic parameter that measures the ability of sediments to conduct heat, characterizing the heat flux per unit temperature gradient. Thermal conductivity is influenced by lithology, porosity, saturation, and the state of hydrate occurrence, and is a key thermodynamic parameter determining the formation heat transfer rate, temperature field evolution, and hydrate phase stability.

[0034] Specific heat capacity: refers to the amount of heat absorbed by a unit mass of sediment to increase its temperature by 1°C, characterizing the reservoir's heat storage capacity. It is closely related to sediment composition, pore fluid, and hydrate content, determining the magnitude of formation temperature changes and thermal inertia, and is a fundamental parameter for simulating the temperature field evolution and thermodynamic processes of hydrate systems.

[0035] Thermal diffusivity is a thermophysical property parameter that describes the rate of temperature change propagation in sediments, reflecting the diffusion rate of temperature disturbances within the formation. It is determined by thermal conductivity, density, and specific heat capacity, directly controlling the strength of the temperature transfer hysteresis effect and the response speed of the temperature field, and has a significant impact on the unsteady formation and dynamic decomposition processes of hydrates.

[0036] In related technologies, existing methods are insufficient to reveal the true evolution of hydrate systems under non-equilibrium conditions, which restricts accurate resource evaluation and engineering safety assurance.

[0037] In view of this, this invention provides a method and related equipment for simulating the dynamic evolution of hydrates based on the temperature transfer hysteresis effect. This method acquires basic survey data of the target study area and extracts key thermodynamic and geological parameters through standardization processing; acquires seismic data of the target study area, identifies the seafloor-like reflector layer and quantifies its actual burial depth below the seafloor, obtaining a depth dataset of the seafloor-like reflector layer; based on the key thermodynamic parameters and the pre-set boundary conditions of geological history periods, a transient heat conduction equation is constructed by coupling the temperature transfer hysteresis effect, outputting the transient stratigraphic temperature field of different geological history periods; based on the key geological parameters, a hydrate phase equilibrium equation is established, and the hydrate phase equilibrium equation is coupled with the transient stratigraphic temperature field time-sequentially to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain in different geological history periods; wherein, the theoretical depth trajectory includes the depth result, adjustment direction, and adjustment rate of the bottom boundary of the hydrate stability domain; the depth offset is identified and the non-equilibrium state is assessed by comparing the depth dataset with the depth result; based on the adjustment direction, adjustment rate, and non-equilibrium state, the spatial occurrence and transformation mode of hydrates and free gas and the engineering geological risk level are predicted. This invention breaks through the traditional steady-state assumption, quantitatively characterizes the control mechanism of temperature transfer hysteresis on the evolution of the steady-state domain, can accurately reveal the non-equilibrium nature of BSR, dynamically reproduce the evolution path of hydrate systems in geological history, and can significantly improve the accuracy of resource estimation, the reliability of geological hazard prediction, and the safety of drilling projects.

[0038] It is understood that the hydrate dynamic evolution simulation method based on the temperature transfer hysteresis effect provided by this invention can be applied to any computer device with data processing and computing capabilities, and this computer device can be various terminals or servers. When the computer device in the embodiment is a server, the server is an independent physical server, or a server cluster or distributed system composed of multiple physical servers, or a cloud server that provides basic cloud computing services such as cloud services, cloud databases, cloud computing, cloud functions, cloud storage, network services, cloud communication, middleware services, domain name services, security services, CDN (Content Delivery Network), and big data and artificial intelligence platforms. Optionally, the terminal can be a smartphone, tablet, laptop, or desktop computer, but it is not limited to these.

[0039] like Figure 1 The diagram shown is a schematic representation of an implementation environment provided by an embodiment of the present invention. (Refer to...) Figure 1 The implementation environment includes at least one terminal 102 and a server 101. The terminal 102 and the server 101 can be connected via a network, either wirelessly or via a wired connection, to complete data transmission and exchange.

[0040] Server 101 can be a standalone physical server, a server cluster or distributed system consisting of multiple physical servers, or a cloud server that provides basic cloud computing services such as cloud services, cloud databases, cloud computing, cloud functions, cloud storage, network services, cloud communication, middleware services, domain name services, security services, CDN (Content Delivery Network), and big data and artificial intelligence platforms.

[0041] Additionally, server 101 can also be a node server in a blockchain network. Blockchain is a novel application model of computer technologies such as distributed data storage, peer-to-peer transmission, consensus mechanisms, and encryption algorithms.

[0042] Terminal 102 can be a smartphone, tablet computer, laptop computer, desktop computer, smart speaker, smartwatch, etc., but is not limited to these. Terminal 102 and server 101 can be directly or indirectly connected via wired or wireless communication, and this embodiment of the invention does not impose any limitations.

[0043] For example, based on Figure 1 The implementation environment shown in this embodiment of the invention provides a hydrate dynamic evolution simulation method based on the temperature transfer hysteresis effect. The following description uses the application of this hydrate dynamic evolution simulation method based on the temperature transfer hysteresis effect in server 101 as an example. It can be understood that this hydrate dynamic evolution simulation method based on the temperature transfer hysteresis effect can also be applied to terminal 102.

[0044] Reference Figure 2 , Figure 2 This is an optional flowchart of a hydrate dynamic evolution simulation method based on temperature transfer hysteresis provided in an embodiment of the present invention. The execution subject of this hydrate dynamic evolution simulation method based on temperature transfer hysteresis can be any of the aforementioned computer devices (including servers or terminals). Figure 2 The method may include, but is not limited to, steps S100 to S600.

[0045] Step S100: Obtain basic survey data of the target study area, and extract key thermodynamic parameters and key geological parameters through standardization processing; It should be noted that in some embodiments, step S100 may include the following steps: acquiring multimodal basic survey data of the target study area; wherein, the basic survey data includes multichannel seismic data, well logging data, sediment samples, in-situ seafloor temperature data, water depth data, sediment pore fluid salinity data, and hydrocarbon gas composition data; preprocessing the basic survey data to obtain preprocessed data; wherein, the preprocessing includes format standardization, depth realignment, time-depth conversion calibration, and outlier removal; based on the preprocessed data, quantifying the key thermodynamic parameters of the sediments using measured methods or empirical formula methods. Key thermodynamic parameters include sediment thermal conductivity, sediment density, and specific heat capacity. Empirical formula methods include using seismic layer velocities corresponding to multichannel seismic data or porosity corresponding to sediment samples to invert sediment thermal conductivity using the Woodside equation or geometric mean model, and using seismic velocities corresponding to multichannel seismic data to derive sediment density and specific heat capacity using the Gardner equation and solid-liquid mixing law. Key geological parameters are determined based on hydrate decomposition gas components corresponding to hydrocarbon gas component data in preprocessed data and formation pore water salinity information corresponding to sediment pore fluid salinity data.

[0046] For example, in some specific implementations, basic survey data is first collected, organized, and standardized, which can be achieved as follows: Objective and significance: To collect and organize basic survey data of the study area, complete data preprocessing and obtain thermal property parameters, and provide data support for subsequent BSR identification, construction of unsteady seabed temperature field and dynamic evolution simulation of hydrate stability domain (GHSZ).

[0047] 1. Collect multichannel seismic survey data, well logging data, sediment samples, in-situ seafloor temperature measurement data, water depth data, sediment pore fluid salinity data, hydrocarbon gas composition data, etc. in the target area.

[0048] 2. Conduct data collection, including format standardization, depth realignment, time-depth conversion calibration, outlier removal, and quality control.

[0049] 3. Based on the actual measurements of hydrate drilling cores, or by using logging-while-drilling information and combining empirical formulas, obtain key thermodynamic parameters such as sediment density (ρ), thermal conductivity (λ), and specific heat capacity (c) in the study area.

[0050] In the absence of actual drilling and logging information, relevant parameters can be estimated using empirical formulas, with minimal impact on the quantitative judgment results. In practical applications, based on extensive research in marine geology and ocean drilling (IODP / DSDP), the thermal conductivity, density, and specific heat capacity of deep-sea sediments are highly correlated with porosity, acoustic velocity, and lithology. Specifically… Sediment thermal conductivity (λ): can be determined by seismic layer velocity (Vp) or porosity (λ / λ). High-precision inversion can be performed using the Woodside equation or the geometric mean model.

[0051] Sediment density (ρ) and specific heat capacity (c): The density can be estimated based on seismic velocity using the Gardner equation, and the specific heat capacity can be calculated according to the solid-liquid mixing law.

[0052] Step S200: Obtain seismic data of the target study area, identify the seafloor reflector and quantify its true burial depth below the seafloor to obtain the depth dataset of the seafloor reflector. It should be noted that, in some embodiments, step S200 may include the following steps: determining the seismic profile and seismic layer velocity model of the target study area based on seismic data; identifying a seafloor-like reflective layer based on the seismic profile by amplitude intensity, negative polarity, and a reflective layer parallel to the seafloor; reading the two-way reflection time of the seafloor-like reflective layer from the seismic data; and converting the two-way reflection time into the actual burial depth below the seafloor using the seismic layer velocity model of the target study area.

[0053] For example, in some specific implementations, BSR integrated recognition and deep computing can be achieved as follows: Objective and significance: To complete the identification, tracking and depth calculation of BSR in the study area, and to provide experimental evidence for subsequent comparison with the dynamic GHSZ evolution results.

[0054] 1. Conduct identification of BSR reflection features with strong amplitude, negative polarity, and near-parallelity to the seabed on high-precision multichannel seismic profiles.

[0055] 2. Read the two-way reflection time of BSR on the seismic profile to determine the distribution characteristics of BSR in the time domain.

[0056] 3. Using the seismic layer velocity model of the study area and combining the two-way travel time of the seabed and BSR reflections, the time-depth conversion calculation of the actual burial depth below the seabed at the BSR interface was completed.

[0057] 4. Obtain the depth profile and spatial distribution characteristics of BSR within the study area.

[0058] Step S300: Based on key thermodynamic parameters and the boundary conditions of the preset geological history period, a transient heat conduction equation is constructed by coupling the temperature transfer hysteresis effect, and the transient stratum temperature field of different geological history period is output. It should be noted that the key thermodynamic parameters include sediment thermal conductivity, sediment density, and specific heat capacity. In some embodiments, step S300 may include the following steps: using sediment thermal conductivity as the numerator and the product of sediment density and specific heat capacity as the denominator, the thermal diffusivity is obtained by ratio calculation; boundary conditions for geological history are obtained; wherein, the boundary conditions include the initial seafloor temperature at the last glacial maximum, the geothermal gradient, and the total change in seafloor temperature since the last glacial maximum; based on the thermal diffusivity, the temperature transfer hysteresis effect is quantitatively characterized using a complementary error function, and then an analytical solution model of the one-dimensional transient heat conduction equation is constructed in combination with the boundary conditions; the analytical solution model is solved by a preset time interval to reconstruct the function of formation temperature change with depth at multiple time points since the last glacial maximum, thereby obtaining the transient formation temperature field for different geological history periods.

[0059] For example, in some specific implementations, the construction of transient stratigraphic temperature fields during geological history can be achieved as follows: Purpose and Significance: Starting from the time when the seafloor temperature began to rise since the Last Glacial Maximum (18 ka BP), this study aims to quantitatively characterize the changes in temperature propagation to the strata below the seafloor and the lag effect of temperature transfer, construct a transient heat transfer model, calculate the variation function of in-situ strata temperature below the seafloor with burial depth during geological history, and output transient temperature field information.

[0060] 1. Based on the key thermodynamic parameters of sediments in the study area, such as density (ρ), thermal conductivity (λ), and specific heat capacity (c), calculate the sediment thermal diffusivity k = λ / (ρ). c).

[0061] 2. The starting point is set as the time since the last glacial maximum when the seafloor temperature has increased (18 ka BP). Initial boundary conditions are determined, and a time evolution sequence of formation temperature is established, with time intervals set sequentially as 18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, and 0 ka BP (i.e., the present).

[0062] 3. Construct a one-dimensional transient heat conduction equation (combined with a compensation function). By setting boundary and initial conditions, use numerical simulation methods to reconstruct the distribution function of stratum temperature with burial depth at each time point in geological history (18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, 0 ka BP (i.e., present)). This is the transient temperature field of sedimentary strata below the seabed that changes with time.

[0063] Specifically, a one-dimensional transient heat conduction model can be constructed. By setting initial and boundary conditions and inputting the aforementioned parameters, the one-dimensional transient heat conduction function can be solved. In a specific embodiment, the initial temperature distribution function in the sedimentary layers below the seabed during the geological history evolution since the Last Glacial Maximum is as follows: ; The boundary conditions are: ; The boundary value solution (i.e., the analytical solution model) is: ; in, Indicates burial depth and time nodes The corresponding in-situ formation temperature, This refers to the seabed temperature during the Last Glacial Maximum. The subsea temperature gradient; This refers to the depth of the strata buried below the seabed. This represents the change in seabed temperature since the Last Glacial Maximum, i.e., the magnitude of the temperature change when the system reaches a steady state. is the thermal diffusivity.

[0064] The erfc() function (compensation function) used is a compensation term added to a linear function. Its essence is a precise characterization and mathematical representation of the temperature field transmission hysteresis effect. The selection of this function form is based on the fact that the erfc() function can accurately describe the temperature response characteristics of a semi-infinite object when the boundary temperature undergoes a step change. In other words, it can effectively and accurately characterize the spatiotemporal variation law of the temperature field under the temperature field transmission hysteresis effect in this invention, and can more realistically reflect the temperature field characteristics during the dynamic evolution of the hydrate stability domain.

[0065] Step S400: Based on key geological parameters, establish the hydrate phase equilibrium equation, and couple the hydrate phase equilibrium equation with the transient formation temperature field time-sequentially to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain in different geological periods. Among them, the theoretical depth trajectory includes the depth result, adjustment direction and adjustment rate of the bottom boundary of the hydrate stability domain; It should be noted that the key geological parameters include hydrate decomposition gas components and formation pore water salinity information. In some embodiments, step S400 may include the following steps: establishing a hydrate phase equilibrium equation based on the hydrate decomposition gas components and formation pore water salinity information; collecting transient formation temperature fields at preset time intervals for each geological period, performing intersection analysis on the corresponding transient formation temperature fields and the hydrate phase equilibrium equation, and determining the depth result of the hydrate stability domain bottom boundary for each geological period through the depth corresponding to the intersection point; and statistically obtaining the adjustment direction and adjustment rate of the hydrate stability domain bottom boundary based on the depth results at each time node in the time sequence.

[0066] For example, in some specific embodiments, the dynamic evolution simulation of hydrate stability domains during geological history can be implemented as follows: Objective and Significance: To conduct hydrate phase equilibrium simulations and couple them with an established unsteady temperature field to complete the dynamic evolution simulation of the theoretical depth of the bottom boundary of the hydrate stability domain at different geological periods (18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, 0 ka BP (i.e., the present)) since the seafloor temperature began to rise (starting from 18 ka BP). This simulation will reconstruct the dynamic evolution process of the hydrate stability domain and determine the trajectory of the bottom boundary of the hydrate stability domain over time.

[0067] 1. Establish phase equilibrium equations: Based on gas composition and salinity, establish phase equilibrium equations (“pressure-P vs. temperature-T” curves) that can accurately calculate the critical temperature and pressure conditions for hydrate formation / decomposition.

[0068] 2. Coupling of temperature and pressure fields: The formation temperature-depth curve (Tz) at each time point is combined with the formation pressure-depth curve (Pz) obtained from seawater depth to construct the formation temperature and pressure field at that time point.

[0069] 3. Solving for the bottom boundary of the stability domain: The "phase equilibrium curve (PT)" and "formation temperature and pressure field curve (PT)" at each time point are intersected. The depth corresponding to the intersection of the two curves is the theoretical depth of the bottom boundary of the hydrate stability domain (BGHSZ) for that geological period.

[0070] 4. Dynamic Evolution Reconstruction: Repeat the above coupling and solution process according to the time series (18 ka BP, 16 ka BP... 0 ka BP) to obtain a series of BGHSZ (Base of Gas Hydrate Stability Zone) depth values, thereby reconstructing its dynamic migration trajectory over time.

[0071] Step S500: Identify the depth offset and evaluate the unbalanced state by comparing the depth dataset with the depth results; It should be noted that in some embodiments, step S500 may include the following steps: spatially overlaying and comparing the depth dataset with the current depth results, quantifying the depth difference between the two as a depth offset; when the depth offset is greater than or equal to a first offset threshold, determining the non-equilibrium state as an early stage of intense non-equilibrium; when the depth offset is less than or equal to a second offset threshold, determining the non-equilibrium state as a late stage of thermodynamic equilibrium; when the depth offset is between the first offset threshold and the second offset threshold, determining the non-equilibrium state as a mid-stage of dynamically adjusting non-equilibrium.

[0072] For example, in some specific implementations, the comparison between the seismic identification BSR and the theoretical stability domain bottom boundary can be achieved as follows: Objective and significance: To conduct a comparative analysis of the measured depth of the BSR identified by seismic identification (seismic or well logging identification) and the current bottom boundary of the hydrate stability domain obtained based on thermodynamic simulation, calculate the offset between the two, and determine the gap between the hydrate stability domain and the final equilibrium state.

[0073] Spatial matching, overlay, and comparative analysis were conducted between the measured depth of the seismically identified bottom seismic revelation (BSR) and the theoretical bottom boundary depth of the hydrate stability domain obtained from dynamic simulation. The offset between the BSR and the theoretical bottom boundary depth of the hydrate stability domain in the seismic profile was calculated. Based on the quantitative characterization results of the temperature transfer hysteresis effect, combined with the dynamic evolution trajectory of the stability domain bottom boundary, the adjusted duration of the hydrate stability domain affected by temperature changes in the study area and the difference from the final equilibrium state of the system were determined.

[0074] The core focus of the offset analysis is the impact of the temperature field transmission lag effect on the dynamic evolution of the hydrate stability domain. Introducing too many factors such as stability domain offset, decomposition rate, and pressure changes would shift the focus of the analysis towards pressure correlation, deviating to some extent from the original research objective. Therefore, the subsequent step S600 does not involve precise quantitative calculations but serves only as an inference and guiding indicator for assessing geological risk levels.

[0075] Therefore, this invention plans to conduct a qualitative analysis by combining the stability domain adjustment time and offset. Here, adjustment time refers to the duration of adjustment the hydrate stability domain has undergone to reach its current state from approximately 18 ka BP to the present, which can be divided into early, middle, and late stages. Specifically, the time division can be roughly set as follows: 18 ka BP to 12 ka BP as the early stage, 12 ka BP to 6 ka BP as the middle stage, and 6 ka BP to the present as the late stage. Adjustments can also be made according to actual circumstances.

[0076] The overall qualitative understanding is as follows: In the late stage, the stable domain basically tends to be in dynamic equilibrium, the actual bottom boundary and the BSR identified on the seismic profile tend to be consistent, the reliability and certainty of drilling are high, and the geological risk level is low; in the early and middle stages, due to the lag effect of the temperature field, there is a significant difference between the actual stable domain distribution and the seismic profile identification results, and the uncertainty of geological risk during the drilling process will increase significantly.

[0077] Therefore, in general, the subsequent steps of S600 are positioned as trend-based and guiding risk prediction, only making qualitative risk classification judgments, without going into too much quantitative analysis, and without focusing on in-depth discussions on offset, gas production rate and pressure changes.

[0078] Step S600: Based on the adjustment direction, adjustment rate and non-equilibrium state prediction, the spatial occurrence and transformation mode of hydrate and free gas and the engineering geological risk level are obtained. It should be noted that in some embodiments, step S600 may include the following steps: defining the spatial transformation relationship between hydrates and free gas based on the adjustment direction of the bottom boundary of the hydrate stability domain; specifically including the following steps: if the adjustment direction is upward migration, the area between the original bottom boundary and the new bottom boundary is indicated as the hydrate decomposition zone and the free gas generation zone; if the adjustment direction is downward migration, the newly covered area is indicated as the hydrate generation zone and the free gas consumption zone; coupling the adjustment rate with the non-equilibrium state to qualitatively classify the engineering geological risk level; wherein, the non-equilibrium state includes the intense non-equilibrium state and the dynamically adjusted non-equilibrium state. Based on thermodynamic equilibrium, the classification of engineering geological risk levels includes: when in a state of intense non-equilibrium and the adjustment rate is greater than the first rate threshold, it is judged as extremely high risk, and drilling is predicted to encounter a large-scale overpressure free gas layer generated by rapid decomposition; when in a state of thermodynamic equilibrium and the adjustment rate is less than the second rate threshold, it is judged as low risk, and the predicted seafloor-like reflector layer basically coincides with the bottom boundary of the actual stable domain; when in a state of thermodynamic equilibrium and the adjustment rate is between the first rate threshold and the second rate threshold, it is judged as medium risk, and there is a dynamic transformation zone between the predicted seafloor-like reflector layer and the bottom boundary of the actual stable domain.

[0079] For example, in some specific embodiments, the prediction of the evolution trend of the hydrate-free gas system can be achieved as follows: Objective and Significance: Based on the dynamic adjustment direction and evolution law of hydrate stability domain, combined with the effect of temperature transfer hysteresis, this study aims to reveal the potential distribution regions of hydrates and free gas, as well as their transformation laws. The core reasoning chain is as follows: 1. Based on the adjustment direction, define the transformation space of "hydrate-free gas": The direction of migration is used to define the transformation region: the vertical migration direction of the bottom boundary of the hydrate stability domain directly indicates "who is replacing whom." If the bottom boundary migrates upward (stability domain contraction): this indicates that the region between the original bottom boundary and the new bottom boundary has changed from a thermodynamically stable region to an unstable region. The solid hydrate in this region is decomposing or has already decomposed, releasing free gas and water. This region is the hydrate decomposition region and also the free gas generation region (gas-water mixing and transport region). If the bottom boundary migrates downward (stability domain expansion): this indicates that the newly covered region has moved from the unstable region to the stable region, and the free gas and water in the pores will recombine into solid hydrates here. This region is the hydrate generation region and also the free gas consumption region.

[0080] 2. Based on the adjustment rate and the non-equilibrium state, qualitatively classify the risk level: The non-equilibrium state of the system directly determines the severity and duration of the transformation process, thus affecting the drilling project risks. Based on the assessment conclusions in step five, the risks can be categorized into three qualitative levels: (1) Extremely high risk (corresponding to intense non-equilibrium state / early stage): Transition mode characteristics: The system is rapidly adjusting from a glacial period to an interglacial period, with the migration rate (adjustment rate) of the bottom boundary of the stable domain reaching its peak. Large-scale hydrates are rapidly decomposing, generating large amounts of free gas and water, leading to a sharp increase in formation pore pressure.

[0081] Engineering geological risk derivation: Drilling in this area is highly likely to encounter overpressured free gas layers or unstable sections containing a mixture of gas, water, and hydrates. This could trigger serious geological disasters such as well kicks, wellbore instability, or even submarine landslides. The BSR, as a "fossil" boundary, cannot reflect the true free gas top boundary, causing serious misleading predictions of drilling fluid density windows.

[0082] (2) Medium risk (corresponding to the late / mid-stage of dynamic adjustment): Characteristics of the conversion mode: The main adjustments of the system have been basically completed, and the migration rate at the bottom boundary of the stability region has slowed significantly. The decomposition and gas production process is nearing its end, and most of the overpressure has been partially released or readjusted through fluid transport and other means.

[0083] Engineering geological risk derivation: The direct risk of hydrate decomposition is reduced, but a dynamic transformation zone with a discernible width still exists between the identified BSR and the true bottom boundary. Low-pressure anomalies or small amounts of free gas may remain within this zone. If the true depth and properties of this zone are not accurately estimated during drilling, drilling fluid loss or trace gas intrusion may still occur, representing a controllable but significant risk.

[0084] (3) Low risk (approaching or reaching thermodynamic equilibrium / late stage): Transformation mode characteristics: The migration rate at the bottom boundary of the stability region is almost stagnant. The hydrate decomposition and free gas transport processes have essentially ended, and the system reaches a new dynamic equilibrium at or very close to this state at the BSR.

[0085] Engineering geological risk derivation: The seismically identified BSR (Bottom Surface Range) largely coincides with the current actual stable domain bottom boundary; therefore, the BSR represents a usable and reliable geological interface. Following traditional hydrate drilling approaches for drilling fluid density and wellbore structure design is safe and carries low risk. The main risk regression is to conventional hydrate cuttings blockage, etc.

[0086] 3. Combining the overall model and risks to arrive at the final prediction result: Combining the above two steps, specific prediction results can be generated, providing guidance for engineering prevention and control. Taking the upward migration of the bottom boundary and the non-equilibrium state as the "middle to late stage" with a displacement of 6.8 meters as an example, the prediction results are as follows: Spatial occurrence and transformation model: It is predicted that the 6.8-meter interval above the current BSR (270 m) to the theoretical bottom boundary (263.2 m) is a dynamic transformation zone from hydrate to free gas. This zone may contain hydrate decomposition residues, trace amounts of newly generated free gas, and pore water.

[0087] Engineering geological risk level: classified as medium risk. The core engineering recommendation is that, in the selection of drilling target areas, the seismic BSR (270 m) should not be simply used as the upper limit of free gas to design the drilling fluid density window. Its lower limit must be increased to cover the slight overpressure or free gas intrusion risk that may exist above the theoretical lower limit (263.2 m).

[0088] To explain in detail the principle of the technical solution of the present invention, the overall process of the present invention will be described below with reference to some specific embodiments. It is easy to understand that the following is an explanation of the technical principle of the present invention and should not be regarded as a limitation of the present invention.

[0089] First, it should be noted that existing technologies mostly revolve around the identification of hydrate stability regions or multiphysics evolution simulation. Their core approach remains limited to the steady-state equilibrium assumption of hydrate systems, focusing primarily on determining stability region boundaries under specific scenarios or simulating instantaneous responses during extraction. They lack quantitative analysis of the lag effect of temperature transfer over geological history and fail to incorporate this lag effect as a key variable into the dynamic evolution simulation framework of hydrate systems for comprehensive consideration. Overall, the shortcomings of existing technologies are mainly reflected in the following aspects: 1. Existing theoretical models have assumption biases and do not fully consider the impact of temperature transfer hysteresis.

[0090] Existing technologies, when constructing computational models for hydrate stability domains, assume that the seafloor geothermal field and sedimentary temperature field are in a long-term stable equilibrium state, failing to consider the temperature transfer lag effect that is prevalent in geological history—that is, the objective fact that sedimentary temperature changes lag behind thermal disturbances in the geological environment, and hydrate phase transitions lag behind dynamic adjustments in the temperature field. This idealized steady-state assumption, differing from actual geological processes, causes the model to fail to reflect the true phase evolution characteristics of the hydrate system under tectonic thermal disturbances and transient temperature changes. This results in a cognitive bias between the BSR and the boundary of the hydrate stability domain, failing to reveal the non-equilibrium nature of the BSR.

[0091] 2. The understanding of the dynamic evolution of hydrate systems is one-sided and limited to the scope of steady-state (equilibrium) research.

[0092] Existing research focuses only on the distribution and evolution of hydrates under steady-state equilibrium conditions, neglecting the formation mechanism and evolution of non-equilibrium (transient) hydrate systems. Since the Last Glacial Maximum, under the influence of temperature transfer lag effects, deep hydrate systems have maintained non-equilibrium dynamic adjustment characteristics for a long period. Traditional techniques cannot capture this key feature, resulting in limitations in our understanding of the transient processes of hydrate formation, decomposition, and reformation, and a complete theoretical system for the dynamic evolution of hydrates has not yet been established.

[0093] 3. Discrepancies exist between resource assessment and exploration application, which restricts the scientific nature of hydrate drilling and extraction practices.

[0094] Based on cognitive biases in the steady-state equilibrium theory of hydrates, existing technologies often mistakenly treat the bottom boundary of the hydrate system's non-equilibrium state (BSR) as the bottom boundary of the stability domain during the calculation of hydrate resources and the spatial identification of hydrates and free gas. This easily leads to overestimation or underestimation of the hydrate occurrence range and resource reserves. Furthermore, it fails to accurately characterize the dynamic occurrence relationship between hydrates and free gas, resulting in biased judgments about the spatial distribution patterns of hydrates and free gas. These limitations and deficiencies in geological judgment directly restrict the selection of hydrate exploration targets, the determination of drilling and production targets, and the scientific formulation of pilot production plans. This can easily lead to problems such as low exploration efficiency, increased drilling and production risks, and poor resource development benefits, thus hindering the industrialization process of natural gas hydrate exploration and development.

[0095] In view of the shortcomings of existing technologies, this invention addresses the problems of biased steady-state assumptions, insufficient understanding of non-equilibrium evolution mechanisms, and limited accuracy of resource evaluation in existing hydrate stability domain simulations. It aims to establish a dynamic evolution simulation method and technical process for hydrate stability domains based on the temperature transfer lag effect. The main objective of this invention is to overcome the limitations of traditional steady-state thermodynamic models by introducing transient heat conduction theory and a temperature transfer lag compensation function. It constructs a transient evolution model of the temperature field, coupling multiple parameters such as sediment thermal conductivity, thermal diffusivity, and specific heat capacity. By quantifying the lag effect of temperature transfer during geological history and the time delay characteristics of hydrate phase transitions, it reveals the essence of transient evolution of hydrate systems from the perspective of temperature transfer thermodynamics, clarifies the differences between the BSR (Boundary Sequence Regulator) and the true stability domain boundary under non-equilibrium conditions, and solves the problem of misjudgment of resource occurrence range caused by neglecting the temperature transfer lag effect in traditional methods. This provides a more scientific and reliable theoretical basis and technical support for accurate evaluation of hydrate resources, selection of drilling and production targets, and engineering risk prevention and control.

[0096] This invention addresses the problem that existing natural gas hydrate (BSR) stability domain prediction technologies generally rely on steady-state models and neglect the temperature transfer lag effect. This leads to a limited understanding of the BSR formation mechanism, confined to static equilibrium, and unable to characterize the dynamic evolution of the stability domain. Consequently, it causes systematic biases in resource assessment and engineering prediction. Figure 3 As shown, the technical solution of this invention is based on high-resolution seismic data, logging-while-drilling (BSR) data, geothermal gradient measurements, sediment thermophysical parameters, and temperature field evolution information from geological history. It establishes a complete technical sequence of "basic data collection - BSR identification - transient temperature field reconstruction - stable domain calculation - dynamic evolution simulation of hydrate systems" to quantitatively measure the lag process of temperature field transmission during geological history, thereby achieving unsteady-state and dynamic simulation and prediction of the spatiotemporal evolution of hydrate stable domains. Specifically, the embodiments of this invention can be implemented through the following steps: Step 1: Collection, organization, and standardization of basic survey data; Objective and significance: To collect and organize basic survey data of the study area, complete data preprocessing and obtain thermal property parameters, and provide data support for subsequent BSR identification, construction of unsteady seabed temperature field and dynamic evolution simulation of hydrate stability domain (GHSZ).

[0097] 1. Collect multichannel seismic survey data, well logging data, sediment samples, in-situ seafloor temperature measurement data, water depth data, sediment pore fluid salinity data, hydrocarbon gas composition data, etc. in the target area.

[0098] 2. Conduct data collection, including format standardization, depth realignment, time-depth conversion calibration, outlier removal, and quality control.

[0099] 3. Based on the actual measurements of hydrate drilling cores, or by using logging-while-drilling information and combining empirical formulas, obtain key thermodynamic parameters such as sediment density (ρ), thermal conductivity (λ), and specific heat capacity (c) in the study area.

[0100] In the absence of actual drilling and logging information, relevant parameters can be estimated using empirical formulas, with minimal impact on the quantitative judgment results. In practical applications, based on extensive research in marine geology and ocean drilling (IODP / DSDP), the thermal conductivity, density, and specific heat capacity of deep-sea sediments are highly correlated with porosity, acoustic velocity, and lithology. Specifically… Sediment thermal conductivity (λ): can be determined by seismic layer velocity (Vp) or porosity (λ / λ). High-precision inversion can be performed using the Woodside equation or the geometric mean model.

[0101] Sediment density (ρ) and specific heat capacity (c): The density can be estimated based on seismic velocity using the Gardner equation, and the specific heat capacity can be calculated according to the solid-liquid mixing law.

[0102] Step 2: BSR Comprehensive Recognition and Depth Calculation; Objective and significance: To complete the identification, tracking and depth calculation of BSR in the study area, and to provide experimental evidence for subsequent comparison with the dynamic GHSZ evolution results.

[0103] 1. Conduct identification of BSR reflection features with strong amplitude, negative polarity, and near-parallelity to the seabed on high-precision multichannel seismic profiles.

[0104] 2. Read the two-way reflection time of BSR on the seismic profile to determine the distribution characteristics of BSR in the time domain.

[0105] 3. Using the seismic layer velocity model of the study area and combining the two-way travel time of the seabed and BSR reflections, the time-depth conversion calculation of the actual burial depth below the seabed at the BSR interface was completed.

[0106] 4. Obtain the depth profile and spatial distribution characteristics of BSR within the study area.

[0107] Step 3: Constructing the transient stratigraphic temperature field during geological history; Purpose and Significance: Starting from the time when the seafloor temperature began to rise since the Last Glacial Maximum (18 ka BP), this study aims to quantitatively characterize the changes in temperature propagation to the strata below the seafloor and the lag effect of temperature transfer, construct a transient heat transfer model, calculate the variation function of in-situ strata temperature below the seafloor with burial depth during geological history, and output transient temperature field information.

[0108] 1. Based on the key thermodynamic parameters of sediments obtained in step one, such as density (ρ), thermal conductivity (λ), and specific heat capacity (c), calculate the sediment thermal diffusivity k = λ / (ρ). c).

[0109] 2. The starting point is set as the time since the last glacial maximum when the seafloor temperature has increased (18 ka BP). Initial boundary conditions are determined, and a time evolution sequence of formation temperature is established, with time intervals set sequentially as 18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, and 0 ka BP (i.e., the present).

[0110] 3. Construct a one-dimensional transient heat conduction equation (combined with a compensation function). By setting boundary and initial conditions, use numerical simulation methods to reconstruct the distribution function of stratum temperature with burial depth at various time points in geological history (18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, 0 ka BP (i.e., present)). Figure 4 As shown, this represents the transient temperature field of sedimentary strata below the seabed that varies over time.

[0111] The erfc() function (compensation function) used is a compensation term set on the basis of a linear function. Its essence is a fine characterization and mathematical description of the hysteresis effect of temperature field transmission.

[0112] The basis for choosing this function form is that the erfc() function can accurately describe the temperature response characteristics of a semi-infinite object when the boundary temperature undergoes a step change. In other words, it can accurately depict the spatiotemporal variation law of the temperature field under the temperature field transmission hysteresis effect in this invention, and can more realistically reflect the temperature field characteristics in the dynamic evolution process of the hydrate stability domain.

[0113] Step 4: Dynamic evolution simulation of hydrate stability domains during geological history; Objective and Significance: To conduct hydrate phase equilibrium simulations and couple them with the unsteady temperature field established in step three, to complete the dynamic evolution simulation of the theoretical depth of the bottom boundary of the hydrate stability domain at different geological periods (18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, 0 ka BP (i.e., the present) since the seafloor temperature began to rise (starting from 18 ka BP), to reconstruct the dynamic evolution process of the hydrate stability domain, and to determine the trajectory of the bottom boundary of the hydrate stability domain over time.

[0114] 1. Based on the information on hydrate decomposition gas components and formation pore water salinity obtained in Step 1, a hydrate phase equilibrium equation is established. This equation is then combined with the transient temperature field function of the subsea strata at different geological periods (18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, 0 ka BP (i.e., present)) established in the above steps since the seafloor temperature began to rise (starting from 18 ka BP). For example... Figure 5 As shown, the intersection of the curves represents the theoretical depth of the bottom boundary of the hydrate stability domain in the corresponding geological history period.

[0115] 2. Calculations were performed at set time intervals (18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, 0 ka BP (i.e., the present)) to reconstruct the dynamic evolution of the bottom boundary of the hydrate stability domain since the Last Glacial Maximum, i.e. after the seafloor temperature began to rise. The results revealed the trajectory of the bottom boundary of the hydrate stability domain over time and the influence of the temperature transfer lag effect on the migration of the hydrate stability domain.

[0116] Step 5: Compare the earthquake-identified BSR with the theoretical stability domain bottom boundary; Objective and significance: To conduct a comparative analysis of the measured depth of the BSR identified by seismic identification (seismic or well logging identification) and the current bottom boundary of the hydrate stability domain obtained based on thermodynamic simulation, calculate the offset between the two, and determine the gap between the hydrate stability domain and the final equilibrium state.

[0117] like Figure 6As shown, spatial matching, superposition, and comparative analysis were conducted between the measured depth of the seismically identified BSR and the theoretical bottom boundary depth of the hydrate stability domain obtained from dynamic simulation. The offset between the BSR and the theoretical bottom boundary depth of the hydrate stability domain in the seismic profile was calculated. Based on the quantitative characterization results of the temperature transfer hysteresis effect, combined with the dynamic evolution trajectory of the stability domain bottom boundary, the adjusted duration of the hydrate stability domain affected by temperature changes in the study area and the difference from the final equilibrium state of the system were determined.

[0118] The core focus of the study on offset is the impact of temperature field transmission lag on the dynamic evolution of hydrate stability domains. Introducing too many factors such as stability domain offset, decomposition rate, and pressure changes would shift the focus of the study towards pressure-related analysis, deviating to some extent from the original research objectives. Therefore, subsequent step six will not involve detailed quantitative calculations, but will only serve as an extension and guiding measure for assessing geological risk levels.

[0119] Therefore, this invention plans to conduct a qualitative analysis by combining the stability domain adjustment time and offset. Here, adjustment time refers to the duration of adjustment the hydrate stability domain has undergone to reach its current state from approximately 18 ka BP to the present, which can be divided into early, middle, and late stages. Specifically, the time division can be roughly set as follows: 18 ka BP to 12 ka BP as the early stage, 12 ka BP to 6 ka BP as the middle stage, and 6 ka BP to the present as the late stage. Adjustments can also be made according to actual circumstances.

[0120] The overall qualitative understanding is as follows: In the late stage, the stable domain basically tends to be in dynamic equilibrium, the actual bottom boundary and the BSR identified on the seismic profile tend to be consistent, the reliability and certainty of drilling are high, and the geological risk level is low; in the early and middle stages, due to the lag effect of the temperature field, there is a significant difference between the actual stable domain distribution and the seismic profile identification results, and the uncertainty of geological risk during the drilling process will increase significantly.

[0121] Therefore, in general, the subsequent step six is ​​positioned as a trend-based and guiding risk prediction, which only makes qualitative risk classification judgments, without going into too much quantitative analysis, and without focusing on in-depth discussions on topics such as offset, gas production rate and pressure changes.

[0122] Step Six: Prediction of the Evolution Trend of the Hydrate-Free Gas System; Objective and significance: Based on the dynamic adjustment direction and evolution law of hydrate stability domain, combined with the effect of temperature transfer hysteresis, this study aims to reveal the potential distribution areas of hydrates and free gas and their transformation laws.

[0123] 1. For example Figure 7As shown, based on the direction and rate of dynamic adjustment in the stable domain and the BSR indication characteristics in the seismic profile, the evolution stages and future adjustment trends of the hydrate system are clarified, the potential stable distribution area and unstable decomposition area of ​​hydrates are determined, the spatial occurrence and transformation relationship of hydrates and free gas are clarified, and the distribution evolution mode of hydrates and free gas under the control of temperature transfer lag effect is summarized.

[0124] 2. Based on the above evolution model, the risks of formation pore pressure fluctuations and free gas release caused by hydrate decomposition during drilling can be further predicted, providing specific reference for engineering control measures such as drilling mud ratio optimization.

[0125] To more clearly illustrate the specific implementation process of the technical solution of this invention, this invention uses a key sea area as a case study area, and combines the actual geological and geophysical data of the area to demonstrate the complete operation steps of the method in detail. This embodiment aims to verify the feasibility and operability of the method in practical applications by conducting dynamic evolution simulation of the hydrate stability domain based on the temperature transfer hysteresis effect in this sea area, and to predict the distribution and evolution trend of hydrates and associated free gas in the study area.

[0126] Step 1: Collection, organization, and standardization of basic data for the study area; Purpose and significance: To complete the collection, collation, preprocessing, and measurement and calculation of thermal physical parameters of basic data for a key sea area, providing data support for subsequent simulation work and clarifying important basic geological parameters of the study area.

[0127] 1. Collect basic geological survey data, including multichannel seismic data, well logging data, in-situ temperature measurements, seawater depth, formation salinity, alkane gas composition, and bottom water temperature, for key sea areas.

[0128] 2. Conduct unified data processing, in-depth correction, and outlier removal. Determine the key geological parameters for hydrate phase equilibrium simulation and transient temperature field establishment in the study area: alkane gas composition: hydrate decomposition gas (CH4: 99.11%, C2H6: 0.88%, C3H8: 0.01%), pore water salinity of 35‰, geothermal gradient of 50℃ / km, and seawater depth of 1240m.

[0129] 3. Using measured data or well logging data, and in conjunction with empirical formulas, obtain the sediment density ρ, thermal conductivity λ, and specific heat capacity c of the study area. Based on this, calculate the thermal diffusivity. =λ / (ρ c).

[0130] Step 2: BSR Comprehensive Recognition and Depth Calculation; Objective and significance: To identify, track, and calculate the depth of the BSR interface in key marine areas, obtain BSR depth information of the study area, and provide a basis for subsequent comprehensive comparative analysis with hydrate stability domains.

[0131] 1. Conduct detailed interpretation of high-precision multichannel seismic profiles in key sea areas and complete the identification and tracking of BSR interfaces.

[0132] 2. Read the two-way reflection time of the BSR and use the regional seismic layer velocity model to complete the time-depth conversion calculation of the actual burial depth of the BSR below the seabed.

[0133] 3. Determine the burial depth and spatial distribution characteristics of seismic identification BSRs in the target study area.

[0134] 4. In this specific embodiment, through BSR identification in the seismic profile and time-depth conversion calculation, the seismic identification BSR depth in the study area is obtained as 270 m below the seabed.

[0135] Step 3: Constructing the transient stratigraphic temperature field during geological history; Objective and Significance: Starting from the time when the seafloor temperature began to rise since the Last Glacial Maximum (18 ka BP), a transient heat transfer model is constructed to establish the evolution sequence of formation temperature over time under the influence of the temperature transfer lag effect, and the in-situ formation temperature at different time points is calculated. ) with burial depth ( The function of ) is used to output transient temperature field information for different geological periods.

[0136] 1. In this specific embodiment, the time point of seabed temperature rise since the Last Glacial Maximum (18kaBP) is set as the starting point.

[0137] 2. Construct a one-dimensional transient heat conduction model. By setting initial and boundary conditions and inputting the above parameters, solve for the one-dimensional transient heat conduction function. In a specific embodiment, the initial temperature distribution function in the sedimentary layers below the seabed during the geological history evolution since the Last Glacial Maximum is as follows: (1) The boundary conditions are: (2) The boundary value solution is: (3) in, This refers to the seabed temperature during the Last Glacial Maximum. The subsea temperature gradient; This refers to the depth of the strata buried below the seabed. This represents the change in seabed temperature since the Last Glacial Maximum, i.e., the magnitude of the temperature change when the system reaches a steady state. This represents the thermal diffusivity. In a specific embodiment, the initial temperature is set by combining global sea-level changes and measured scientific parameters obtained from hydrate drilling and ocean drilling expeditions in key sea areas. =1.2℃; geothermal gradient =50℃ / km; Seabed temperature change =2℃; Thermal diffusivity =5.9*10 -7 m 2 / s.

[0138] 3. In this specific embodiment, the formation temperature ( ) is calculated at each 2-ka interval corresponding to the geological history time points (18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, 0 ka BP (i.e., present)) since the Last Glacial Maximum (18 ka BP). ) with burial depth ( The function is calculated using a continuously changing transient temperature field function.

[0139] 4. In this specific embodiment, as... Figure 4 As shown, the curves depicting the temperature variations of subsea strata at different geological time points (including 18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, and 0 ka BP (present)) at 2-ka intervals since 18 ka BP are constructed for the study area. These curves characterize the temperature distribution of subsea strata at different geological historical stages since the Last Glacial Maximum. For shallow seafloor strata (e.g., within 0-50 m below the seafloor), the strata temperature can reach equilibrium within a relatively short geological period; however, the lag in temperature transfer becomes particularly pronounced with increasing burial depth. For deep strata (e.g., hydrate reservoirs 200-300 m below the seafloor), Figure 4 In the area marked by b in the middle, the temperature difference of the strata at different geological time points is still quite obvious. The strata are significantly affected by the temperature transmission lag effect. There is a significant time lag in the transmission of seafloor temperature signals to the deep strata, which verifies the characteristics of the unsteady temperature field.

[0140] 5. In this specific embodiment, the formation temperature varies significantly during the 18-8 ka BP period. As geological time progresses, the formation temperature variation within the same time interval gradually decreases, the formation temperature difference caused by temperature transfer gradually narrows, the temperature transfer lag effect weakens, and the formation temperature field gradually stabilizes.

[0141] Step 4: Dynamic evolution simulation of hydrate stability domains during geological history; Objective and Significance: By coupling the hydrate phase equilibrium curve with the unsteady formation temperature field established in step three, dynamic evolution simulations of the theoretical depth of the bottom boundary of the hydrate stability domain are completed for different geological periods since 18 ka BP (18 ka BP, 16 ka BP, 14 ka BP, 12 ka BP, 10 ka BP, 8 ka BP, 6 ka BP, 4 ka BP, 2 ka BP, 0 ka BP (i.e., present)). This will reconstruct the dynamic evolution process of the hydrate stability domain and determine its migration and evolution trajectory over time.

[0142] 1. In this specific implementation process, based on the decomposition gas composition (CH4: 99.11%, C2H6: 0.88%, C3H8: 0.01%) and formation pore water salinity (35‰) information obtained in Step 1 from the pressurized hydrate samples in the study area, a hydrate phase equilibrium equation was established. Based on this, it was combined with the transient functions of the temperature field from different geological periods established in Step 3, such as... Figure 5 As shown, the depth at the intersection point is the bottom boundary depth of the hydrate stability domain in the corresponding geological history period.

[0143] 2. In this specific embodiment, as follows: Figure 5 As shown, from 18 ka BP to the present (0 ka BP), the bottom boundary of the hydrate stability domain has risen from 283.5 m below the seabed in geological history to 263.2 m below the seabed today. Specifically, during the geological history period from 18 to 10 ka BP, the hydrate stability domain exhibited a relatively rapid upward adjustment; from 8 ka BP to the present, the rate of upward adjustment has slowed significantly. Within the resolvable accuracy range of seismic data, it can be approximated that the stability domain has approached equilibrium during this period. This change reveals the migration and evolution of the hydrate stability domain during its dynamic adjustment process since the last glacial period, influenced by the lag effect of temperature transmission.

[0144] Step 5: Compare the earthquake-identified BSR with the theoretical stability domain bottom boundary; Purpose and significance: To conduct a comparative analysis of the measured depth of BSR (based on seismic or well logging identification) and the current bottom boundary of the hydrate stability domain obtained through thermodynamic simulation, calculate the offset between the two, and determine the gap between the hydrate stability domain and the final equilibrium state.

[0145] 1. Conduct spatial matching, superposition, and comparative analysis of the measured depth of the seismically identified BSR and the theoretical bottom boundary depth of the hydrate stability domain obtained from dynamic simulation, and calculate the offset between the BSR and the current theoretical bottom boundary depth of the hydrate stability domain in the seismic profile. Based on the quantitative characterization results of the temperature transfer hysteresis effect, combined with the dynamic evolution trajectory of the stability domain bottom boundary (e.g., Figure 6 As shown in the figure, the adjusted duration of the hydrate stability domain in the study area affected by temperature changes and the difference from the final equilibrium state of the system were determined.

[0146] 2. In this specific embodiment, as Figure 6 As shown, through comprehensive analysis of seismic identification and well logging information, the current seismically identified BSR depth is approximately 270 m below the seabed, which differs from the current theoretical bottom boundary depth of the hydrate stability domain of 263.2 m by 6.8 m. The theoretical bottom boundary depth of the hydrate stability domain is shallower than the current BSR depth, indicating that the hydrate system has not yet reached a complete equilibrium state and is still in a dynamic evolution process of continuous upward migration, that is, the bottom boundary of the stability domain is adjusting to the shallower layers, in order to gradually approach the final equilibrium state of the system.

[0147] Step Six: Prediction of the Evolution Trend of the Hydrate-Free Gas System; Objective and significance: Based on the dynamic adjustment direction and evolution law of hydrate stability domain, combined with the effect of temperature transfer hysteresis, this study aims to reveal the potential distribution areas of hydrates and free gas and their transformation laws.

[0148] 1. Based on the direction and rate of dynamic adjustment in the stable region and the BSR indication characteristics in seismic profiles, clarify the evolutionary stages and future adjustment trends of the hydrate system, determine the potential stable distribution area and unstable decomposition region of hydrates, and clarify the spatial transformation relationship between hydrates and free gas (e.g., Figure 7 (As shown). Based on this, we can further predict the risks of formation pore pressure fluctuations and free gas release caused by hydrate decomposition during drilling, providing specific references for engineering control measures such as drilling mud ratio optimization.

[0149] 2. In this specific embodiment, combined with Figure 7 The study area hydrates were described Dynamic transformation of the free gas system: Hydrate decomposition occurs in the region above the BSR seismic profile (e.g. Figure 7(The portion marked 'a' in the middle) The fluids and free gases produced by the decomposition of hydrates migrate towards the shallow seafloor until the system tends to stabilize; when the hydrate system reaches equilibrium, the BSR on the seismic profile should coincide with the current bottom boundary of the hydrate stability domain (e.g., ...). Figure 7 (The part marked with b in the middle).

[0150] 3. In this specific embodiment, it is important to note during the selection of target areas and resource evaluation for hydrate drilling: the upper part of the seismic BSR is not entirely composed of solid hydrates; its response may only reflect the non-equilibrium dynamic evolution process of hydrate transformation into free gas. Therefore, in engineering geological exploration, hydrate drilling deployment, and resource evaluation, especially in avoiding drilling geological risks and conducting hydrate exploration... When conducting a comprehensive evaluation of free gas, the potential impact of non-equilibrium evolution in the stability region should be fully considered.

[0151] In summary, this invention, based on exploration practices regarding the non-equilibrium evolution of hydrate systems caused by the temperature transfer lag effect, takes transient temperature field reconstruction as its starting point and achieves dynamic evolution simulation of the steady-state domain by coupling multiple physical parameters. Its innovation lies in the quantitative characterization of the temperature transfer lag effect and the construction of transient simulation techniques. Specifically, it achieves its purpose through techniques such as transient heat conduction equations modified by compensation functions, reconstruction of temperature field evolution sequences at multiple time nodes, and quantification of the contribution rate of the temperature transfer lag effect. Specifically, the technical solution of this invention proposes the following technical innovations: Innovation Point 1: A nonsteady-state temperature field model coupled with temperature transfer hysteresis effect was constructed.

[0152] This invention overcomes the limitations of existing technologies that generally employ steady-state heat conduction models by quantitatively introducing the temperature transfer lag effect caused by seafloor temperature changes since the Last Glacial Maximum into the simulation of hydrate system evolution. By establishing a dynamic evolution model of the temperature field based on transient heat conduction equations (combined with compensation functions), and integrating key thermophysical parameters such as sediment porosity, thermal conductivity, specific heat capacity, and thermal diffusivity, the non-equilibrium and transient evolution process of the seafloor temperature field is reconstructed. This model serves as the core theoretical foundation of this method, providing accurate and reliable temperature field input conditions for subsequent dynamic simulation of hydrate systems.

[0153] Innovation Point 2: Achieved quantitative simulation of non-equilibrium evolution of hydrate stability domain based on dynamic temperature-pressure coupling.

[0154] This invention overcomes the technical limitations of traditional static equilibrium simulation. It takes the unsteady temperature field coupled with the temperature transfer hysteresis effect as input, and simultaneously substitutes it into the hydrate phase equilibrium equation. It iteratively calculates the dynamic migration trajectory of the bottom boundary of the hydrate stability domain in a time sequence, realizing a fully quantitative simulation of the evolution of the hydrate stability domain from static equilibrium to dynamic non-equilibrium. It reveals the intrinsic mechanism of the deviation between the BSR and the steady-state theoretical stability domain from the essential level of the temperature transfer hysteresis effect, and realizes the accurate prediction of the spatiotemporal evolution law of the stability domain.

[0155] Innovation Point 3: It can effectively serve resource quantity assessment and engineering geological early warning.

[0156] This invention innovates a technical approach for hydrate resource evaluation and drilling geological risk prevention. Based on the revealed dynamic evolution law of hydrate stability domain under the temperature transfer lag effect and the changing characteristics of the hydrate-free gas system, it can effectively serve geological resource evaluation and engineering geological early warning. This technology can accurately delineate hydrate distribution, assess dynamic changes in resource quantity, predict formation pressure fluctuations, free gas release, and related geological risks during drilling, providing guidance for drilling implementation and engineering control (such as mud ratio optimization) to avoid drilling geological risks.

[0157] Compared with the prior art, the embodiments of the present invention have at least the following beneficial effects: To address the issue that existing natural gas hydrate stability domain prediction technologies generally rely on steady-state models and neglect the temperature transfer lag effect, resulting in a limited understanding of the BSR formation mechanism and difficulty in characterizing the dynamic evolution process of the stability domain, which in turn leads to a series of geological and engineering problems such as systematic biases in resource assessment and engineering prediction, this invention proposes a dynamic evolution simulation method for hydrate stability domains based on the temperature transfer lag effect.

[0158] Compared with traditional steady-state models, the significant advantage of this invention lies in its quantitative incorporation of the key dynamic process of temperature transfer lag effect into the simulation of hydrate system evolution. By constructing a transient heat conduction model, the lag process of seafloor temperature changes being transmitted to deep strata is accurately characterized, thereby achieving dynamic, non-equilibrium simulation of the evolution of the hydrate stability domain over geological time. This method not only theoretically explains the discrepancy between observed BSR depths and steady-state theoretical predictions, but also more realistically reflects the response of the hydrate system to changes in paleoenvironment (especially paleotemperature) during geological history. At the application level, by providing information on the dynamic evolution history and current non-equilibrium state of the stability domain, this method can improve the accuracy of hydrate resource assessment, enhance the understanding of the dynamic occurrence of free gas, and provide more reliable scientific support for drilling target selection, exploitation scheme design, and early warning of seafloor geological hazards, overcoming the limitations of traditional steady-state methods in dynamic process and time-sensitive prediction.

[0159] Specifically, the dynamic evolution of the natural gas hydrate stability domain is the result of the combined effects of multiple factors, including thermal disturbance, changes in temperature and pressure conditions, and sediment thermal properties. Among these factors, the temperature transfer hysteresis effect, as a key dynamic mechanism, directly leads to the deviation between the seismically identified bottom seismic hydrate (BSR) and the theoretical bottom boundary of the hydrate stability domain, which has a significant impact on resource assessment and exploration deployment.

[0160] This invention overcomes the shortcomings of traditional techniques, which are limited to steady-state equilibrium assumptions, lack consideration of temperature transfer lag effects, and have insufficient multi-parameter coupling. Based on the intrinsic relationship between global environmental changes and hydrate system responses since the Last Glacial Maximum (18 ka BP), it takes the quantitative characterization of temperature transfer lag effects as a starting point and establishes a complete technical sequence of "basic data collection - BSR identification - transient temperature field reconstruction - steady-state domain calculation - dynamic evolution simulation of hydrate systems" to achieve the technical transformation from static equilibrium simulation to dynamic non-equilibrium simulation.

[0161] This invention effectively overcomes the technical bottlenecks of existing technologies, such as insufficient understanding of the non-equilibrium evolution mechanism of hydrate systems, low simulation accuracy, and disconnection from actual exploration needs. It explores and forms a scientific, accurate, and operable dynamic evolution simulation method, providing a new technical approach for the accurate exploration and effective development of highly enriched hydrates, and has important practical value for promoting the commercial exploitation of natural gas hydrates.

[0162] The technical challenge of this invention lies in the quantitative characterization of the temperature transfer hysteresis effect and the solution of the transient heat conduction equation. By introducing a compensation function, clarifying the parameter calculation model, and establishing a standardized process, the above problems are effectively solved, ensuring the reliability and practicality of the simulation results.

[0163] like Figure 8 As shown, this embodiment of the invention also provides a hydrate dynamic evolution simulation device 900 based on the temperature transfer hysteresis effect, which can implement the above-described method. This device may include: The first module 901 is used to acquire basic survey data of the target study area and extract key thermodynamic parameters and key geological parameters through standardization processing. The second module 902 is used to acquire seismic data of the target study area, identify the seafloor reflector and quantify its true burial depth below the seafloor, and obtain the depth dataset of the seafloor reflector. The third module 903 is used to construct transient heat conduction equations based on key thermodynamic parameters and pre-set boundary conditions for geological history periods, and output transient stratum temperature fields for different geological history periods. Module 4, 904, is used to establish hydrate phase equilibrium equations based on key geological parameters. The hydrate phase equilibrium equations are coupled with transient formation temperature fields over time to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain at different geological histories. The theoretical depth trajectory includes the depth result, adjustment direction, and adjustment rate of the bottom boundary of the hydrate stability domain. The fifth module 905 is used to identify depth offsets and evaluate unbalanced states by comparing the depth dataset with the depth results. Module 6, 906, is used to predict the spatial occurrence and transformation mode of hydrates and free gas, as well as the engineering geological risk level, based on the adjustment direction, adjustment rate, and non-equilibrium state.

[0164] It is understood that the content of the above method embodiments is applicable to the present device embodiments. The specific functions implemented by the present device embodiments are the same as those of the above method embodiments, and the beneficial effects achieved are also the same as those achieved by the above method embodiments.

[0165] This invention also provides an electronic device, which includes a memory and a processor. The memory stores a computer program, and the processor executes the computer program to implement the method described above. This electronic device can be any smart terminal, including tablet computers, in-vehicle computers, etc.

[0166] It is understood that the content of the above method embodiments is applicable to this device embodiment. The specific functions implemented by this device embodiment are the same as those of the above method embodiments, and the beneficial effects achieved are also the same as those achieved by the above method embodiments.

[0167] like Figure 9 As shown, Figure 9 The hardware structure of an electronic device 1000 according to another embodiment is illustrated. The electronic device 1000 includes: The processor 1001 can be implemented using a general-purpose CPU (Central Processing Unit), microprocessor, application-specific integrated circuit (aSIC), or one or more integrated circuits, and is used to execute relevant programs to implement the technical solutions provided in the embodiments of the present invention. The memory 1002 can be implemented as a read-only memory (ROM), a static storage device, a dynamic storage device, or a random access memory (RaM). The memory 1002 can store the operating system and other application programs. When the technical solutions provided in the embodiments of this specification are implemented through software or firmware, the relevant program code is stored in the memory 1002 and is called and executed by the processor 1001. Input / output interface 1003 is used to implement information input and output; The communication interface 1004 is used to enable communication and interaction between this device and other devices. Communication can be achieved through wired means (such as USB, network cable, etc.) or wireless means (such as mobile network, WIFI, Bluetooth, etc.). Bus 1005 transmits information between various components of the device (e.g., processor 1001, memory 1002, input / output interface 1003, and communication interface 1004); The processor 1001, memory 1002, input / output interface 1003 and communication interface 1004 are connected to each other within the device via bus 1005.

[0168] The electronic device embodiments described above are merely illustrative. The units described as separate components may or may not be physically separate; that is, they may be located in one place or distributed across multiple network units. Some or all of the modules can be selected to achieve the purpose of this embodiment according to actual needs.

[0169] This invention also provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the above-described method.

[0170] It is understood that the content of the above method embodiments is applicable to this storage medium embodiment. The specific functions implemented in this storage medium embodiment are the same as those in the above method embodiments, and the beneficial effects achieved are also the same as those achieved in the above method embodiments.

[0171] This invention also provides a computer program product, including a computer program that, when executed by a processor, implements the above-described method.

[0172] It is understood that the content of the above method embodiments is applicable to the embodiments of this program product. The specific functions implemented by the embodiments of this program product are the same as those of the above method embodiments, and the beneficial effects achieved are also the same as those achieved by the above method embodiments.

[0173] Memory, as a non-transitory computer-readable storage medium, can be used to store non-transitory software programs and non-transitory computer-executable programs. Furthermore, memory may include high-speed random access memory, and may also include non-transitory memory, such as at least one disk storage device, flash memory device, or other non-transitory solid-state storage device. In some embodiments, memory may optionally include memory remotely located relative to the processor, and these remote memories can be connected to the processor via a network. Examples of such networks include, but are not limited to, the Internet, intranets, local area networks, mobile communication networks, and combinations thereof.

[0174] The present invention provides a method, apparatus, electronic device, storage medium, and program product for simulating the dynamic evolution of hydrates based on the temperature transfer hysteresis effect. This method acquires basic survey data of the target study area and extracts key thermodynamic and geological parameters through standardization processing. It then acquires seismic data of the target study area, identifies a seafloor-like reflector layer, and quantifies its actual burial depth below the seafloor, obtaining a depth dataset of the seafloor-like reflector layer. Based on the key thermodynamic parameters and pre-defined boundary conditions for geological periods, a transient heat conduction equation is constructed by coupling the temperature transfer hysteresis effect, outputting transient stratigraphic temperature fields for different geological periods. A hydrate phase equilibrium equation is established based on the key geological parameters, and the hydrate phase equilibrium equation is coupled temporally with the transient stratigraphic temperature field to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain for different geological periods. The theoretical depth trajectory includes the depth result, adjustment direction, and adjustment rate of the bottom boundary of the hydrate stability domain. The depth offset is identified and non-equilibrium states are assessed by comparing the depth dataset with the depth result. Based on the adjustment direction, adjustment rate, and non-equilibrium states, the spatial occurrence and transformation mode of hydrates and free gas, as well as the engineering geological risk level, are predicted. This invention breaks through the traditional steady-state assumption, quantitatively characterizes the control mechanism of temperature transfer hysteresis on the evolution of the steady-state domain, can accurately reveal the non-equilibrium nature of BSR, dynamically reproduce the evolution path of hydrate systems in geological history, and can significantly improve the accuracy of resource estimation, the reliability of geological hazard prediction, and the safety of drilling projects.

[0175] The embodiments described in this invention are for the purpose of more clearly illustrating the technical solutions of the embodiments of this invention, and do not constitute a limitation on the technical solutions provided by the embodiments of this invention. As those skilled in the art will know, with the evolution of technology and the emergence of new application scenarios, the technical solutions provided by the embodiments of this invention are also applicable to similar technical problems.

[0176] Those skilled in the art will understand that the technical solutions shown in the figures do not constitute a limitation on the embodiments of the present invention, and may include more or fewer steps than shown, or combine certain steps, or different steps.

[0177] The device embodiments described above are merely illustrative. The units described as separate components may or may not be physically separate; that is, they may be located in one place or distributed across multiple network units. Some or all of the modules can be selected to achieve the purpose of this embodiment according to actual needs.

[0178] Those skilled in the art will understand that all or some of the steps in the methods disclosed above, as well as the functional modules / units in the systems and devices, can be implemented as software, firmware, hardware, or suitable combinations thereof.

[0179] The preferred embodiments of the present invention have been described above with reference to the accompanying drawings, but this does not limit the scope of the claims of the present invention. Any modifications, equivalent substitutions, and improvements made by those skilled in the art without departing from the scope and spirit of the present invention should be within the scope of the claims of the present invention.

Claims

1. A method for simulating the dynamic evolution of hydrates based on the temperature transfer hysteresis effect, characterized in that, The method includes the following steps: Acquire basic survey data for the target study area, and extract key thermodynamic and geological parameters through standardized processing; Seismic data of the target study area were acquired, the seafloor reflector was identified and its actual burial depth below the seafloor was quantified, and a depth dataset of the seafloor reflector was obtained. Based on the key thermodynamic parameters and the boundary conditions of the preset geological history period, a transient heat conduction equation is constructed by coupling the temperature transfer hysteresis effect, and the transient stratum temperature field of different geological history period is output. Based on the key geological parameters, a hydrate phase equilibrium equation is established. The hydrate phase equilibrium equation is then coupled with the transient formation temperature field time-sequentially to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain at different geological histories. The theoretical depth trajectory includes the depth result, adjustment direction, and adjustment rate of the bottom boundary of the hydrate stability domain. The depth offset is identified and the unbalanced state is evaluated by comparing the depth dataset with the depth results. Based on the adjustment direction, the adjustment rate, and the non-equilibrium state, the spatial occurrence and transformation mode of hydrates and free gas, as well as the engineering geological risk level, are predicted.

2. The method according to claim 1, characterized in that, The acquisition of basic survey data for the target study area, followed by the extraction of key thermodynamic and geological parameters through standardization processing, includes the following steps: Acquire multimodal basic survey data for the target study area; wherein, the basic survey data includes multichannel seismic data, well logging data, sediment samples, in-situ seafloor temperature data, water depth data, sediment pore fluid salinity data, and hydrocarbon gas composition data; The basic survey data is preprocessed to obtain preprocessed data; wherein, the preprocessing includes format unification, depth realignment, time-depth conversion calibration, and outlier removal. Based on the preprocessed data, the key thermodynamic parameters of the sediments are quantified using either experimental or empirical formula methods. The key thermodynamic parameters include sediment thermal conductivity, sediment density, and specific heat capacity. The empirical formula method includes using the seismic layer velocity corresponding to the multichannel seismic data or the porosity corresponding to the sediment sample to invert the sediment thermal conductivity through the Woodside equation or geometric mean model, and using the seismic velocity corresponding to the multichannel seismic data to derive the sediment density and specific heat capacity through the Gardner equation and solid-liquid mixing law. The key geological parameters are determined based on the hydrate decomposition gas components corresponding to the hydrocarbon gas component data in the preprocessed data and the formation pore water salinity information corresponding to the sediment pore fluid salinity data.

3. The method according to claim 1, characterized in that, The process of identifying the apparent seabed reflective layer and quantifying its actual burial depth below the seabed includes the following steps: Based on the earthquake data, determine the seismic profile and seismic layer velocity model for the target study area; Based on the seismic profile, the seabed-like reflective layer is identified by amplitude intensity, negative polarity, and reflective layer parallel to the seabed. Read the two-way reflection time of the seafloor-like reflector from the seismic data; Using the seismic layer velocity model of the target study area, the two-way reflection time is converted into the actual burial depth below the seabed.

4. The method according to claim 1, characterized in that, The key thermodynamic parameters include sediment thermal conductivity, sediment density, and specific heat capacity. Based on these key thermodynamic parameters and pre-defined boundary conditions for different geological periods, a transient heat conduction equation is constructed by coupling the temperature transfer hysteresis effect, outputting the transient formation temperature field for different geological periods. This includes the following steps: The thermal diffusivity is obtained by taking the thermal conductivity of the sediment as the numerator and the product of the sediment density and the specific heat capacity as the denominator, and then calculating the ratio. Obtain boundary conditions for geological historical periods; wherein, the boundary conditions include the initial seafloor temperature at the last glacial maximum, the geothermal gradient, and the total change in seafloor temperature since the last glacial maximum. Based on the thermal diffusivity, the temperature transfer hysteresis effect is quantitatively characterized using complementary error functions, and then an analytical solution model of the one-dimensional transient heat conduction equation is constructed by combining the boundary conditions. The expression for the analytical solution model is as follows: ， In the formula, Indicates burial depth and time nodes The corresponding in-situ formation temperature, Indicates the initial temperature of the seabed. Represents the geothermal gradient. This represents the total change in seabed temperature. Represents the complementary error function. Indicates the thermal diffusivity; By solving the analytical solution model at preset time intervals, the function of formation temperature variation with depth at multiple time points since the Last Glacial Maximum is reconstructed, and the transient formation temperature field at different geological historical periods is obtained.

5. The method according to claim 1, characterized in that, The key geological parameters include hydrate decomposition gas components and formation pore water salinity information. The process of establishing a hydrate phase equilibrium equation based on these key geological parameters, and then coupling the hydrate phase equilibrium equation with the transient formation temperature field time-series to obtain the theoretical depth trajectory of the hydrate stability domain floor at different geological histories, includes the following steps: A hydrate phase equilibrium equation is established based on the hydrate decomposition gas components and the formation pore water salinity information. Transient formation temperature fields from various geological periods are collected at preset time intervals. The corresponding transient formation temperature fields are intersected with the hydrate phase equilibrium equations to obtain the depth of the bottom boundary of the hydrate stability domain for each geological period. Based on the depth results at each time point in the time sequence, the adjustment direction and adjustment rate of the bottom boundary of the hydrate stability domain are statistically obtained.

6. The method according to claim 1, characterized in that, The step of identifying depth offset and evaluating imbalance by comparing the depth dataset with the depth results includes the following steps: The depth dataset is spatially superimposed and compared with the current depth results, and the depth difference between the two is quantified as the depth offset. When the depth offset is greater than or equal to the first offset threshold, the non-equilibrium state is determined to be an intense non-equilibrium state in the early stage. When the depth offset is less than or equal to the second offset threshold, the non-equilibrium state is determined to be a thermodynamic equilibrium state in the late stage. When the depth offset is between the first offset threshold and the second offset threshold, the unbalanced state is determined to be a dynamic adjustment unbalanced state in the intermediate stage.

7. The method according to claim 1, characterized in that, The method for predicting the spatial occurrence and transformation mode and engineering geological risk level of hydrates and free gas based on the adjustment direction, the adjustment rate, and the non-equilibrium state includes the following steps: Based on the adjustment direction of the bottom boundary of the hydrate stability domain, the spatial transformation relationship between hydrate and free gas is defined; specifically, the following steps are included: if the adjustment direction is upward migration, the area between the original bottom boundary and the new bottom boundary is indicated as the hydrate decomposition zone and the free gas generation zone; if the adjustment direction is downward migration, the newly covered area is indicated as the hydrate generation zone and the free gas consumption zone. The adjustment rate is coupled with the non-equilibrium state to qualitatively classify the engineering geological risk level; The non-equilibrium state includes a violent non-equilibrium state, a dynamically adjusted non-equilibrium state, and a thermodynamic equilibrium state. The classification of engineering geological risk levels includes: when in the violent non-equilibrium state and the adjustment rate is greater than a first rate threshold, it is judged as extremely high risk, and the drilling is predicted to encounter a large-scale overpressure free gas layer generated by rapid decomposition; when in the thermodynamic equilibrium state and the adjustment rate is less than a second rate threshold, it is judged as low risk, and the predicted seafloor reflector layer basically coincides with the bottom boundary of the actual stable domain; when in the thermodynamic equilibrium state and the adjustment rate is between the first rate threshold and the second rate threshold, it is judged as medium risk, and there is a dynamic transformation zone between the predicted seafloor reflector layer and the bottom boundary of the actual stable domain.

8. A hydrate dynamic evolution simulation device based on the temperature transfer hysteresis effect, characterized in that, The device includes: The first module is used to acquire basic survey data of the target study area and extract key thermodynamic parameters and key geological parameters through standardization processing. The second module is used to acquire seismic data of the target study area, identify the seafloor reflector and quantify its true burial depth below the seafloor, and obtain the depth dataset of the seafloor reflector. The third module is used to construct transient heat conduction equations based on the key thermodynamic parameters and the boundary conditions of the preset geological history period, by coupling the temperature transfer hysteresis effect, and output transient stratum temperature fields of different geological history periods. The fourth module is used to establish a hydrate phase equilibrium equation based on the key geological parameters, and to couple the hydrate phase equilibrium equation with the transient formation temperature field time-sequentially to obtain the theoretical depth trajectory of the bottom boundary of the hydrate stability domain at different geological histories; wherein, the theoretical depth trajectory includes the depth result of the bottom boundary of the hydrate stability domain, the adjustment direction and the adjustment rate. The fifth module is used to identify the depth offset and evaluate the unbalanced state by comparing the depth dataset with the depth results; The sixth module is used to predict the spatial occurrence and transformation mode of hydrates and free gas and the engineering geological risk level based on the adjustment direction, the adjustment rate and the non-equilibrium state.

9. An electronic device, characterized in that, The electronic device includes a memory and a processor, the memory storing a computer program, and the processor executing the computer program to implement the method according to any one of claims 1 to 7.

10. A computer program product, characterized in that, The computer program product includes a computer program that, when executed by a processor, implements the method according to any one of claims 1 to 7.

Images