A simulation method for quantitatively predicting evolution of gas leakage and crater formation in a hydrate zone and related equipment

Abstract

The application discloses a simulation method for quantitatively predicting gas seepage and crater formation evolution process in a hydrate area and related equipment, and the method comprises the following steps: obtaining three-dimensional seismic data of a target area, and obtaining thickness data through layer velocity modeling and time-depth conversion; based on the thickness data and pre-acquired physical parameters, a first quantitative relationship between gas column growth rate and seepage velocity is quantified by using Darcy's law; based on the first quantitative relationship, a second quantitative relationship of gas migration time is obtained through integral method processing; the first quantitative relationship and the critical condition of quicksand are solved together to obtain a third quantitative relationship of crater depth, and then a crater space form constraint is constructed; and based on the crater space form constraint and the pre-acquired average methane flux, the total methane flux of the target area is simulated. The application can realize physical constraint quantitative prediction of deep fluid seepage process and crater formation evolution, 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 simulation method and related equipment for quantitatively predicting the gas leakage and pit formation evolution process in hydrate zones. Background Technology

[0002] Submarine pockmarks are widespread, depression-like landforms on the surface of deep-sea sediments, typically circular or elliptical in shape. Their formation is primarily controlled by the upward seepage of overpressurized gas-bearing fluids (such as methane) from deep within the sedimentary layer. In natural gas hydrate deposits, pockmark development is often closely related to the dynamic changes in the free gas layer beneath the hydrate stability zone. When overpressure arises from tectonic activity or hydrate decomposition, the gas migrates upward along seepage channels. If the gas pressure exceeds the capillary closure threshold of the overlying hydrate layer, it penetrates the stability zone and enters the shallow sediments. During this process, the gas-water two-phase flow alters the pore pressure distribution, reducing the effective stress in the sediments. When the shear strength is insufficient to maintain the stability of the original sedimentary strata, the sediments undergo localized liquefaction or disturbance and are eroded or carried away by bottom currents, ultimately forming pockmarks on the seabed. Therefore, pockmarks are not only a geomorphic manifestation of fluid seepage but also a direct surface response to the dynamic evolution of deep hydrate systems.

[0003] In resource exploration, the Mapeng area is often associated with active gas seepage systems, serving as an important indicator for identifying underlying free gas enrichment zones or hydrate decomposition zones. In marine engineering, the Mapeng area has loose sediment structure and low bearing capacity, and may contain unexplored shallow gas reservoirs or fluid channels, posing a potential threat to subsea pipeline laying and platform pile foundation stability. In terms of the ecological environment, the continuously seeping methane provides an energy source for chemosynthetic microorganisms and cold seep communities (such as tube worms and mussels), supporting the long-term evolution of this unique ecosystem. These multiple values ​​make Mapeng research significant in the fields of resources, engineering, and ecology.

[0004] Currently, the identification of submarine pits mainly relies on geological and geophysical methods such as multibeam bathymetry, side-scan sonar, shallow seismic profiling, and multichannel seismic data. These methods can effectively obtain information on the planar morphology, spatial distribution, and geometric relationship with deep structures (such as chimneys and faults). However, existing research mainly focuses on qualitative descriptions or statistical correlations: on the one hand, although geophysical data can delineate the location of pits, limitations in vertical resolution and interpretability make it difficult to accurately infer their formation depth, duration, and leakage intensity; on the other hand, the analysis of formation mechanisms is mostly based on empirical hypotheses proposed from field observations (such as "overpressure-induced instability" and "hydrate decomposition-driven"), lacking physical models that quantitatively correlate deep fluid transport processes with shallow geomorphological responses.

[0005] Furthermore, most existing methods only provide static characterization of existing pits. Reflecting on the physical mechanisms, pit formation involves the geological coupling of multiple processes, including gas-water-hydrate multiphase flow, pore pressure evolution, and sediment mechanical response. However, there is currently a lack of a unified computational framework that can integrate seismic inversion parameters, fluid dynamics equations, and sediment instability criteria. This makes it difficult to quantitatively predict key parameters (such as leakage rate, critical liquefaction depth, fluid movement time, and pit size) under physical constraints. Summary of the Invention

[0006] The main objective of this invention is to provide a simulation method, apparatus, electronic device, storage medium, and program product for quantitatively predicting gas leakage and pit formation evolution in hydrate zones, aiming to solve at least one problem in the prior art.

[0007] To achieve the above objectives, one aspect of this invention proposes a simulation method for quantitatively predicting gas leakage and pit formation evolution in hydrate zones, the method comprising: Acquire 3D seismic data of the target area, and obtain thickness data through layer velocity modeling and time-depth conversion; Based on thickness data and pre-collected physical parameters, the first quantitative relationship between the gas column growth rate and the seepage velocity is quantified using Darcy's law. Based on the first quantification relationship, a second quantification relationship for gas migration time is obtained through integration; the second quantification relationship is used to provide geological constraints for the onset time and duration of gas leakage. By combining the first quantitative relationship with the critical condition of quicksand formation, the third quantitative relationship of pit depth is obtained, and then the spatial morphological constraints of pit are constructed. The total methane flux in the target region was obtained by simulating the spatial morphology constraints of the pit and the pre-collected average methane flux.

[0008] In some embodiments, the thickness data includes the thickness of the free gas layer and the thickness of the stable domain from the bottom boundary of the hydrate stability domain to the seabed. The thickness data is obtained through layer velocity modeling and time-depth conversion, including the following steps: Extract the root mean square velocity and two-way travel time of each key layer from the 3D seismic data; Based on the root mean square velocity and two-way travel time, the Dix formula is used to convert them into layer velocity, and a layer velocity model is established. Based on the layer velocity model, the time-domain seismic data volume is converted into the depth-domain seismic data volume to obtain the absolute depth of key geological interfaces. The thickness of the stable domain and the thickness of the free gas layer are obtained by quantifying the difference between the absolute depths of each key geological interface.

[0009] In some embodiments, the thickness data includes the thickness of the free gas layer and the thickness of the stable domain from the bottom boundary of the hydrate stable domain to the seabed. Physical parameters include water phase density, gas phase density, gravitational acceleration, gas viscosity, gas column height, water viscosity, water flow column height, and absolute permeability of the sedimentary body. The gas column height is a time-varying variable, and the water flow column height is determined by the difference between the stable domain thickness and the gas column height. Based on the thickness data and pre-collected physical parameters, the first quantitative relationship between the gas column growth rate and the seepage velocity is quantified using Darcy's law, including the following steps: The fourth quantitative relationship of the driving pressure difference is obtained by quantifying the product of the difference between the water phase density and the gas phase density with the gravitational acceleration and the thickness of the free gas layer. The fifth quantitative relationship for gas column resistance is obtained by quantifying the product of the ratio of gas viscosity and absolute permeability of the sediment with gas column height and seepage velocity. The sixth quantitative relationship for water column resistance is obtained by quantifying the ratio of water viscosity to absolute permeability of sediment, multiplied by the height of the water column and the seepage velocity. The fifth and sixth quantitative relationships are combined to obtain the seventh quantitative relationship for total flow resistance; Based on the relationship that the driving force equals the total flow resistance under steady seepage conditions, the fourth and seventh quantitative relationships are combined to obtain the first quantitative relationship of the seepage velocity. Among them, since the height of the air column changes with time, the growth rate of the air column is equal to the seepage velocity.

[0010] In some embodiments, based on the first quantization relationship, a second quantization relationship for gas transport time is obtained through integration, including the following steps: Based on the first quantization relation, the method of separation of variables is used for integration to obtain the eighth quantization relation; The expression for the first quantization relation is: ; The expression for the eighth quantization relation is: ; In the formula, Indicates the seepage velocity; Indicates the rate of air column growth. Represents the differential symbol; and For thickness data, Indicates the thickness of the free gas layer. This indicates the thickness of the stability domain from the bottom boundary of the hydrate stability domain to the seabed. , , , , , , For physical parameters, Indicates the absolute permeability of the sediment. This indicates the density of the aqueous phase. Indicates gas phase density, Represents gravitational acceleration. Indicates gas viscosity, Indicates the height of the air column. Indicates the viscosity of water. The height of the water column; Based on the eighth quantization relation, the left and right boundaries are integrated using preset integration limits, and then the integration results of the left and right boundaries are combined to obtain the second quantization relation of gas transport time.

[0011] In some embodiments, the first quantitative relationship is combined with the critical condition for quicksand formation to obtain the third quantitative relationship for the depth of the pit, including the following steps: The critical conditions for quicksand formation were constructed based on seepage velocity, water viscosity, relative water permeability, absolute permeability of sediment, sediment particle density, water phase density, and gravitational acceleration. Substituting the first quantitative relationship into the critical condition of quicksand formation, we can then transform and obtain the third quantitative relationship for the depth of the pit. The expression for the critical condition of quicksand formation is as follows: ; The expression for the first quantification relation is: ; The expression for the third quantification relation is: ; In the formula, Indicates the seepage velocity; and For thickness data, Indicates the thickness of the free gas layer. This indicates the thickness of the stability domain from the bottom boundary of the hydrate stability domain to the seabed. , , , , , , , , For physical parameters, Indicates the absolute permeability of the sediment. Indicates the relative permeability of water. Indicates the density of sediment particles. This indicates the density of the aqueous phase. Indicates gas phase density, Represents gravitational acceleration. Indicates gas viscosity, Indicates the height of the air column. Indicates the viscosity of water. The height of the water column. ; Indicates the depth of the pit.

[0012] In some embodiments, constructing the spatial morphological constraints of the pits includes the following steps: The planar distribution range of the pits was identified using underwater acoustic detection data; The projected area of ​​the pit is obtained by extracting the pit depth based on the planar display range and the third quantization relationship, and is used as the spatial morphological constraint of the pit.

[0013] In some embodiments, the spatial morphological constraint of the pits includes the projected area of ​​the pits. The total methane flux of the target region is simulated based on the spatial morphological constraint of the pits and the pre-collected average methane flux, including the following steps: Based on geological constraints, the average methane flux is predicted by conducting spatial representativeness analysis through preset observation points. The total methane flux in the target region is obtained by simulating the product of the average methane flux and the projected area of ​​the pit.

[0014] To achieve the above objectives, another aspect of this invention proposes a simulation device for quantitatively predicting gas leakage and pit formation evolution in hydrate zones. The device includes: The first module is used to acquire three-dimensional seismic data of the target area and obtain thickness data through layer velocity modeling and time-depth conversion. The second module is used to quantify the first quantitative relationship between the gas column growth rate and the seepage velocity based on thickness data and pre-collected physical parameters using Darcy's law. The third module is used to obtain a second quantification relationship for gas migration time based on the first quantification relationship through an integral method; wherein, the second quantification relationship is used to provide geological constraints for the start time and duration of gas leakage. The fourth module is used to combine the first quantitative relationship with the critical condition of quicksand formation to obtain the third quantitative relationship of pit depth, and then construct the spatial morphological constraint of pit. The fifth module is used to simulate the total methane flux in the target area based on the spatial morphology constraints of the pits and the pre-collected average methane flux.

[0015] 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.

[0016] 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.

[0017] 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.

[0018] The embodiments of this invention include at least the following beneficial effects: This invention provides a simulation method, apparatus, electronic device, storage medium, and program product for quantitatively predicting gas leakage and pit formation evolution in hydrate areas. This scheme acquires three-dimensional seismic data of the target area, obtains thickness data through layer velocity modeling and time-depth conversion; based on the thickness data and pre-acquired physical parameters, it uses Darcy's law to quantify the first quantitative relationship between the gas column growth rate and seepage velocity; based on the first quantitative relationship, it obtains a second quantitative relationship for gas migration time through integration; wherein, the second quantitative relationship provides geological constraints for the onset time and duration of gas leakage; combining the first quantitative relationship with the critical condition of quicksand formation, it obtains a third quantitative relationship for pit depth, thereby constructing a spatial morphological constraint for pits; based on the spatial morphological constraint of pits and the pre-acquired average methane flux, it simulates the total methane flux of the target area. The embodiments of this invention, by integrating three-dimensional seismic data, Darcy's law, and the critical condition of quicksand formation, establish a series of quantitative relationships between the gas column growth rate, seepage velocity, gas migration time, and pit depth, realizing quantitative prediction of the physical constraints of deep fluid leakage processes and pit formation evolution. This invention provides the first quantitative correlation between deep fluid transport and shallow geomorphological response, overcoming the ambiguity of geophysical interpretation. Furthermore, by coupling multiphase flow, pore pressure evolution, and sediment mechanical response within a unified framework, this invention can accurately predict leakage rates, critical liquefaction depths, transport times, and pit sizes. Specifically, this invention can provide a reliable basis for free gas identification in resource exploration, marine engineering hazard assessment, and ecological environment methane flux estimation, significantly enhancing the scientific rigor and application value of pit research. Attached Figure Description

[0019] Figure 1 This is a schematic diagram of an implementation environment for a simulation method provided in this invention for quantitatively predicting gas leakage and pit formation evolution in hydrate zones. Figure 2 This is a schematic flowchart of a simulation method for quantitatively predicting gas leakage and pit formation evolution in hydrate zones, provided by an embodiment of the present invention. Figure 3 This is a schematic diagram illustrating an example of a geological model for the formation of submarine pits by gas seepage in a hydrated region, as provided in an embodiment of the present invention. Figure 4 This is a schematic diagram of the overall process of the simulation method for quantitatively predicting gas leakage and pit formation evolution in hydrate zones, provided in an embodiment of the present invention. Figure 5 This is a schematic diagram illustrating the variation of gas leakage time with gas column height under different simulation conditions provided in the embodiments of the present invention; Figure 6 This is a schematic diagram illustrating an example of a numerical model of the gas leakage process from a hydrated region to a shallow layer, provided in an embodiment of the present invention. Figure 7 This is a schematic diagram illustrating the distribution of submarine pockmark depth under different simulation conditions provided by embodiments of the present invention, showing the relationship between water flow impedance and effective formation stress. Figure 8 This is a schematic diagram illustrating an example of fluid leakage monitoring and gas flux estimation in a submarine pit area provided in an embodiment of the present invention; Figure 9 This is a schematic diagram of the structure of a simulation device for quantitatively predicting the gas leakage and pit formation evolution process in hydrate zones, provided in an embodiment of the present invention. Figure 10 This is a schematic diagram of the structure of an electronic device provided in an embodiment of the present invention. Detailed Implementation

[0020] 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.

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

[0022] 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.

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

[0024] To facilitate understanding of the technical solution of this invention, the technical skills that may be involved in the technical solution of this invention will first be explained: Hydrates are cage-like solid compounds formed by gases such as methane and water under low temperature and high pressure conditions on the seabed. Their formation occupies pore space, reduces permeability, and alters seepage paths. Their caprock sealing effect can even trigger local overpressure, making them a key geological factor in regulating fluid leakage behavior and pit formation.

[0025] Submarine pits are unique geomorphic structures formed on the seabed surface by the continuous upward seepage of gas-bearing fluids from deep within the seabed to shallower layers. This seepage causes localized liquefaction or instability of sediments due to increased pore pressure and decreased effective stress. The size (depth and area) typically approximates the intensity and duration of the seepage.

[0026] Pressure gradient: refers to the rate of change of pressure per unit depth in the vertical direction of the seabed. It is mainly driven by the difference between the overpressure generated by the accumulation of free gas in the deep sea and the static rock pressure in the shallow sea. It is the direct driving force for gas-bearing fluids to overcome formation resistance and migrate upward.

[0027] Absolute permeability: refers to the inherent seepage capacity of shallow seabed sediments unaffected by fluid type. It depends only on pore structure (such as pore size distribution, connectivity, tortuosity, etc.) and is a fundamental physical property parameter for assessing whether gas-bearing fluids can migrate upwards from deep layers.

[0028] Relative permeability: When multiple phases such as free gas, water, and hydrates coexist in shallow seabed pores, the ratio of the effective permeability of a particular phase (e.g., gas) to the absolute permeability. Relative permeability typically decreases significantly with increasing hydrate saturation, directly controlling the actual leakage rate of gas-containing fluids.

[0029] Darcy's Law: A fundamental principle describing the linear seepage behavior of fluids in porous media (such as seabed sediments) under low-velocity laminar flow conditions. In porous media systems, the fluid seepage velocity is positively correlated with the media permeability and the pressure gradient along the seepage path, and negatively correlated with the fluid viscosity. It is a key basis for analyzing fluid transport behavior in seabed sediments.

[0030] Capillary force: The additional pressure generated by the surface tension at the gas-water interface in the tiny pores of seafloor sediments, directed from the water phase to the gas phase, and its magnitude is inversely proportional to the pore radius. In fine-grained sediments, capillary force can effectively block the rise of free gas, forming a microscopic resistance to leakage.

[0031] Fluid resistance refers to the combined resistance exerted by shallow seabed sediments and their internal multiphase fluids (gas, water, and hydrates) on the upward seepage of gas-containing fluids. When hydrates crystallize and fill pores or reduce pore throat size, permeability decreases, and fluid resistance increases significantly, thereby inhibiting the intensity of leakage. It is an important indicator for judging whether leakage can continue.

[0032] Static pressure: The vertical stress generated by the gravity of the overlying seabed sediments, which increases linearly with depth. It is the main background pressure that gas-bearing fluids must overcome during their migration from deep to shallow layers, and together with pore fluid pressure, it determines the effective stress state and stability of the sediments.

[0033] Among related technologies, there is a lack of a unified computational framework that can integrate seismic inversion parameters, fluid dynamics equations and sediment instability criteria, which makes it difficult to quantitatively predict key parameters (such as leakage rate, critical liquefaction depth, fluid movement time, and pit size) under physical constraints.

[0034] In view of this, this invention provides a simulation method and related equipment for quantitatively predicting gas leakage and pit formation evolution in hydrate areas. This method acquires three-dimensional seismic data of the target area, obtains thickness data through layer velocity modeling and time-depth conversion; based on the thickness data and pre-acquired physical parameters, it uses Darcy's law to quantify the first quantitative relationship between the gas column growth rate and seepage velocity; based on the first quantitative relationship, it uses an integral method to obtain a second quantitative relationship for gas migration time; wherein, the second quantitative relationship provides geological constraints for the onset time and duration of gas leakage; the first quantitative relationship is combined with the critical condition for quicksand formation to obtain a third quantitative relationship for pit depth, thereby constructing a spatial morphological constraint for the pits; based on the spatial morphological constraint of the pits and the pre-acquired average methane flux, the total methane flux of the target area is simulated. This invention, by integrating three-dimensional seismic data, Darcy's law, and the critical condition for quicksand formation, establishes a series of quantitative relationships between the gas column growth rate, seepage velocity, gas migration time, and pit depth, realizing quantitative prediction of the physical constraints of deep fluid leakage processes and pit formation evolution. This invention provides the first quantitative correlation between deep fluid transport and shallow geomorphological response, overcoming the ambiguity of geophysical interpretation. Furthermore, by coupling multiphase flow, pore pressure evolution, and sediment mechanical response within a unified framework, this invention can accurately predict leakage rates, critical liquefaction depths, transport times, and pit sizes. Specifically, this invention can provide a reliable basis for free gas identification in resource exploration, marine engineering hazard assessment, and ecological environment methane flux estimation, significantly enhancing the scientific rigor and application value of pit research.

[0035] It is understood that the simulation method for quantitatively predicting gas leakage and pit formation evolution in hydrate zones provided by this invention can be applied to any computer device with data processing and computing capabilities, and this computer device can be various types of terminals or servers. When the computer device in the embodiments 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.

[0036] 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 1The 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.

[0037] 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.

[0038] 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.

[0039] 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.

[0040] For example, based on Figure 1 The implementation environment shown in this embodiment of the invention provides a simulation method for quantitatively predicting the gas leakage and pit formation evolution process in hydrate areas. The following description uses the application of this simulation method for quantitatively predicting the gas leakage and pit formation evolution process in hydrate areas to server 101 as an example. It can be understood that this simulation method for quantitatively predicting the gas leakage and pit formation evolution process in hydrate areas can also be applied to terminal 102.

[0041] Reference Figure 2 , Figure 2 This is an optional flowchart of a simulation method for quantitatively predicting the gas leakage and pit formation evolution process in hydrate zones, provided in an embodiment of the present invention. The subject executing this simulation method for quantitatively predicting the gas leakage and pit formation evolution process in hydrate zones 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 S500.

[0042] Step S100: Obtain three-dimensional seismic data of the target area, and obtain thickness data through layer velocity modeling and time-depth conversion; It should be noted that the thickness data includes the thickness of the free gas layer and the thickness of the stable domain from the bottom boundary of the hydrate stable domain to the seafloor. In some embodiments, the thickness data is obtained through layer velocity modeling and time-depth conversion, which may include the following steps: extracting the root mean square velocity and two-way travel time of each key layer from the 3D seismic data; converting the root mean square velocity and two-way travel time into layer velocity using the Dix formula based on the root mean square velocity and two-way travel time to establish a layer velocity model; converting the time-domain seismic data volume to the depth-domain seismic data volume based on the layer velocity model to obtain the absolute depth of the key geological interface; and quantizing the difference between the absolute depths of each key geological interface to obtain the stable domain thickness and the free gas layer thickness.

[0043] For example, in some specific embodiments, the present invention aims to construct a high-precision three-dimensional geological model to provide an accurate geometric framework and basic parameters for subsequent quantitative calculations. An accurate depth domain model is crucial for eliminating the influence of seismic wave velocity differences and achieving "depth positioning," providing an important data foundation for subsequent depth identification of the bottom boundary of the hydrate stability domain (BGHSZ), seepage simulation, and prediction of pit size. Specifically, this can be achieved through the following steps: 1. Basic parameter extraction and layer velocity model construction: Three-dimensional seismic data of the target area were acquired, and the root mean square velocity (RMS) and two-way travel time (TWT) of each key layer (such as the seafloor, free gas cap, and gas-water interface) were extracted. Based on the extracted basic parameters, the root mean square velocity (RMS) in the time-domain seismic data was converted into layer velocity using the DIX formula to establish a preliminary layer velocity model.

[0044] (1) in: Root mean square velocity (m / s) is obtained from seismic velocity spectrum analysis and is usually extracted through stacking velocity analysis or pre-stack depth migration.

[0045] :depth The seismic wave velocity (m / s) at a given location is the layer velocity.

[0046] Two-way travel time (s): The TWT value of the target layer is read directly from the seismic profile.

[0047] 2. Time-depth conversion calculation: By employing layered velocity models (such as seawater velocity and sedimentary velocity) and two-way travel time (TWT) data, the conversion of time-domain seismic data volumes to depth-domain seismic data volumes is completed, and the absolute depth of key geological interfaces is obtained.

[0048] (2) in: Depth (m): The calculated result, used for all subsequent depth-related analyses.

[0049] Root mean square velocity (m / s) is derived from the earthquake velocity model.

[0050] Two-way travel time (s) is derived from seismic profile interpretation results.

[0051] 3. Error Analysis and Verification: The key layer depths after time-depth conversion are compared with the logging depths of adjacent wells (such as LWD data) to calculate the depth error and ensure that the conversion accuracy meets engineering requirements (usually controlled within 5%). If necessary, the velocity model is iteratively optimized using drilling control points.

[0052] This step will obtain the precise depth of the depth-domain seismic data volume and key interfaces (such as the seabed, BGHSZ, and gas-water interface), which is the basis for subsequent calculations of hydrate stability domains, determination of free gas layer thickness, and identification of gas chimney top boundaries.

[0053] Step S200: Based on the thickness data and pre-collected physical parameters, the first quantitative relationship between the gas column growth rate and the seepage velocity is quantified using Darcy's law. It should be noted that the thickness data includes the thickness of the free gas layer and the thickness of the stable domain from the bottom boundary of the hydrate stable domain to the seabed. Physical parameters include water phase density, gas phase density, gravitational acceleration, gas viscosity, gas column height, water viscosity, water flow column height, and absolute permeability of the sediment. The gas column height is a time-varying variable, and the water flow column height is determined by the difference between the stable domain thickness and the gas column height. In some embodiments, step S200 may include the following steps: quantifying the driving pressure difference based on the product of the difference between water phase density and gas phase density and gravitational acceleration and free gas layer thickness to obtain a fourth quantification relationship; based on gas viscosity, sediment... The fifth quantitative relationship for air column resistance is obtained by quantifying the ratio of absolute permeability of the volume to the product of air column height and seepage velocity. The sixth quantitative relationship for water column resistance is obtained by quantifying the ratio of water viscosity to absolute permeability of the sediment to the product of water column height and seepage velocity. The seventh quantitative relationship for total flow resistance is obtained by summing the fifth and sixth quantitative relationships. According to the relationship that the driving force equals the total flow resistance under steady seepage conditions, the fourth and seventh quantitative relationships are combined to obtain the first quantitative relationship for seepage velocity. Among them, since the air column height changes with time, the air column growth rate is equal to the seepage velocity.

[0054] Exemplary examples, in some specific embodiments, the present invention aims to establish a dynamic mathematical model of the gas leakage process based on the geological model provided in the foregoing steps, quantitatively describing the instantaneous rate of free gas migration from the bottom boundary of the hydrate stability domain to the seabed. This rate is a key dynamic parameter connecting the deep gas source and the shallow geomorphological response (pits), providing crucial dynamic parameters for calculating the leakage time of gas-bearing fluids and assessing the drag force of deep fluid activity on sediments. Specifically, this can be achieved through the following steps: Based on Darcy's law, this model derives the gas-water two-phase flow velocity expression during the upward movement of deep gas by establishing a balance between driving force and flow resistance.

[0055] (1) Darcy's Law: Darcy's law quantitatively describes the seepage behavior of fluids in porous media, and its basic form is: (3) in: The seepage velocity is (m / s). Permeability (m²); Fluid viscosity (Pa·s); Pressure difference (Pa); The seepage path length is in meters (m).

[0056] (2) Fluid driving force: The core driving force for upward gas migration originates from the hydrostatic pressure difference (i.e., buoyancy) generated by the density difference between the gas column and the water column. In practice, it is necessary to constrain the gas saturation of the effective gas layer thickness using well logging data and set a fluid migration saturation threshold of more than 20%—only when the gas saturation reaches this threshold value can the conditions required for gas flow be met, thereby forming the fluid pressure driving migration.

[0057] (4) The driving pressure difference (Pa) represents the buoyancy that pushes the underlying gas upward and displaces the pore fluid.

[0058] Aqueous phase density (kg / m³) was determined by CTD profile measurement combined with pore water salinity analysis.

[0059] Gas phase density (kg / m³) is calculated using the gas state equation. The input parameters are gas composition and formation temperature and pressure conditions, or it can be obtained through drilling and logging data and experimental tests.

[0060] : Gravitational acceleration (m / s²), taken as 9.8 m / s².

[0061] The thickness of the free gas layer (m) was obtained by seismic impedance inversion technology and verified by seismic-well logging joint calibration.

[0062] (3) Flow resistance: As the gas rises, it must overcome the viscous resistance from both the air column and the water column, as shown below: Resistance caused by the air column: (5) in: : Gas partial pressure (Pa) generated by the gas column in the system; Gas viscosity (Pa·s) is calculated based on the gas chromatography analysis results of the decomposition gas from the core of the pressure-maintaining drilling, combined with the formation temperature and pressure conditions.

[0063] : Air column height (m), a variable that changes over time.

[0064] Absolute permeability of sedimentary bodies (m²) is determined by core experiments (such as pulse decay method) on pressure-controlled core samples or estimated by well logging curves.

[0065] Resistance caused by the water column: (6) in: : The partial pressure (Pa) of the fluid generated by the water column in the system; The viscosity of water (Pa·s) was obtained by laboratory testing of pore water samples using a rotational viscometer. The height of the water column (m) is equal to ; It is determined by the difference between the bottom boundary of the hydrate stability domain (BGHSZ) on the depth-domain seismic profile and the seabed depth.

[0066] The total flow resistance is the sum of the resistances of the air column and the water column: (7) in The total flow resistance is expressed in Pa.

[0067] (4) Analysis of the balance between driving force and resistance: Under steady-state seepage conditions, the driving force equals the total flow resistance, that is: (8) Solution of flow rate : (9) because ,but: (10) Due to the height of the air column Over time Change, its growth rate Equal to seepage velocity Therefore: (11) This step, based on Darcy's law, establishes a balance between the driving force and multi-stage flow resistance under two-phase fluid activity, providing a theoretical basis for the quantitative description of gas transport processes. This step relies on the information provided in step S100. and As input parameters, its output parameter is the seepage velocity. (i.e., air column growth rate) This is the direct input for subsequent calculations of fluid leakage time and a key basis for judging the conditions of fluidization.

[0068] Step S300: Based on the first quantization relationship, the second quantization relationship of gas transport time is obtained by processing it through an integration method; The second quantification relationship is used to provide geological constraints for the onset time and duration of gas leakage. It should be noted that in some embodiments, step S300 may include the following steps: based on the first quantization relationship, integration is performed using the method of separation of variables to obtain the eighth quantization relationship; based on the eighth quantization relationship, integration is performed on the left and right boundaries using a preset integration limit, and then the integration results of the left and right boundaries are combined and transformed to obtain the second quantization relationship of gas transport time.

[0069] For example, in some specific embodiments, based on the seepage rate model established in the aforementioned steps, this embodiment of the invention aims to calculate, through an integral method, the time *t* that deep gas takes to seep from the bottom boundary of the hydrate stability domain until it finally reaches the seabed. This timescale is crucial for understanding the periodicity of seepage events and assessing the dynamic evolution rate of the hydrate system. Further, in specific embodiments, the dating results of cold seep carbonates obtained through drilling can be combined, using their formation time as an age constraint, and superimposed with the migration time *t*, to deduce the initiation time of fluid seepage. It should be noted that because there is a time difference between carbonate precipitation and the initial fluid seepage, the derived result is an approximation. The specific implementation process is as follows: Based on the differential equation for the air column growth rate obtained in the preceding steps, this step aims to solve for the air column's growth from an initial height (assumed to be 0) to a target height by integrating in the time and depth domains. Time required In the specific implementation process, the method of separation of variables is used for integration. The following equation can be obtained by transforming formula (11): (12) Set a points limit, allowing points to accumulate from 0 to [amount missing] over time. (The desired time); Integrate the height of the air column from 0 to... (Target Altitude): (13) in and It is the integral variable.

[0070] The result of the left boundary integral is: (14) The result of the right boundary integral is: (15) Make the left and right sides equal: (16) Solve (17) When the air column grows to (Right now Substituting this into the above formula, we can obtain the total time required for the gas to reach the seabed: (18) In the formula: The viscosity of water (Pa·s) is obtained by testing pore water samples in the laboratory using a rotational viscometer, and the test temperature must be consistent with the formation temperature.

[0071] Gas viscosity (Pa·s) is calculated based on the gas chromatography analysis results of the decomposition gas from the core of the pressure-maintaining drilling, combined with the formation temperature and pressure conditions.

[0072] It is determined by the difference between the bottom boundary of the hydrate stability zone (BGHSZ) on the depth-domain seismic profile and the seabed depth.

[0073] The thickness of the free gas layer (m) was obtained by seismic impedance inversion technology and verified by seismic-well logging joint calibration.

[0074] Absolute permeability of sedimentary bodies (m²) is determined by core experiments (such as pulse decay method) on pressure-controlled core samples or estimated by well logging curves.

[0075] , The densities of the aqueous phase and the gas phase (kg / m³) were obtained by CTD profile measurements combined with salinity analysis and calculations based on the gas state equation, respectively.

[0076] This step involves measuring the instantaneous seepage velocity. Integrating the spatial path, we can convert it into the total time required for deep gas to migrate from the bottom boundary of the hydrate layer to the seabed. This allows for a quantitative characterization of the seepage process over time. Furthermore, by combining the dating results of submarine cold seep carbonate rocks (formed by gaseous fluids seeping to the seabed, inducing anaerobic methane oxidation (AOM) and subsequent carbonate precipitation), geological constraints can be provided on the onset and duration of fluid activity.

[0077] Step S400: Combine the first quantification relationship with the critical condition of quicksand formation to obtain the third quantification relationship of pit depth, and then construct the spatial morphological constraint of pit. It should be noted that in some embodiments, the third quantitative relationship of pit depth is obtained by combining the first quantitative relationship with the critical condition of quicksand formation. This may include the following steps: constructing the critical condition of quicksand formation based on seepage velocity, water viscosity, relative water permeability, absolute permeability of sediment, sediment particle density, water phase density, and gravitational acceleration; substituting the first quantitative relationship into the critical condition of quicksand formation, and then transforming it to obtain the third quantitative relationship of pit depth.

[0078] For example, in some specific embodiments, based on the geological model, seepage velocity, and time evolution information formed in the aforementioned steps, the embodiments of the present invention aim to establish a quantitative relationship model between the depth of seafloor pit formation and deep gas reservoir parameters. By linking the deep seepage process with the shallow geomorphological response, a direct theoretical basis and practical means are provided for predicting the potential of deep oil and gas resources through the characteristics of seafloor pits. Specifically, this can be achieved as follows: 1. Critical conditions for quicksand formation in pockmarked areas: When gas-bearing fluids seep from deep strata to shallow seafloor layers, the water flow exerts an upward drag force (i.e., flow resistance) on the shallow sedimentary skeleton. When this drag force exceeds the effective stress of the sediments, the sedimentary skeleton becomes unstable, liquefies, and becomes quicksand, eventually forming pitted structures on the seafloor.

[0079] The critical condition for mechanical equilibrium in this process is: water drag force = effective stress of sediment; According to fluid mechanics and Darcy's law, the drag force per unit area exerted by water flow on the sedimentary skeleton is numerically equal to the pressure gradient during seepage. Therefore, the critical condition for sediment shifting can be expressed as: (19) in: The fluid seepage velocity (m / s) was calculated from the model in step two.

[0080] The viscosity of water (Pa·s) is obtained from the same steps described above.

[0081] : Relative permeability of water (dimensionless), usually taken as 1 when the water phase is dominant.

[0082] Absolute permeability of the sediment (m²), sourced from the same steps described above.

[0083] Sediment particle density (kg / m³) is determined by obtaining seabed sediment samples through hydrate drilling cores and performing gamma-ray density testing or drying and weighing in the laboratory.

[0084] : Aqueous phase density (kg / m³), sourced from the same steps described above.

[0085] : Gravitational acceleration (m / s²), taken as 9.8 m / s².

[0086] 2. Seepage velocity and gas layer overpressure: Under critical fluidization conditions, the seepage velocity The critical value Vcrit is reached. According to Darcy's law and the principle of pressure balance, the seepage velocity can be expressed as: (20) The height of the water column (m) is equal to ; : Identified and determined by the bottom boundary of the hydrate stability domain (BGHSZ) on the depth domain seismic profile; The thickness of the free gas layer (m) is obtained from the seismic inversion in step one. Gas phase density (kg / m³), calculated using the gas law; 3. Establishment of a pit depth prediction model: At the critical moment of sand shifting, the seepage velocity The critical condition is just met, at which point the water column height is... This refers to the thickness of the sediment that is about to be fluidized and removed, i.e., the pit depth. ; By refining the above formulas, a predictive model for the depth of pockmarks can be established. =1 (simplification process), resulting in: (twenty one) After simplification, the final numerical model for calculating the depth of the pit is obtained.

[0087] (twenty two) in: The predicted depth of the pit (m) is the target value for this step.

[0088] The difference between the bottom boundary of the hydrate stability zone (BGHSZ) and the seabed depth.

[0089] : Thickness of free gas layer (m).

[0090] : Aqueous phase density (kg / m³), sourced from the same steps described above.

[0091] : Gas phase density (kg / m³), calculated by the gas law.

[0092] Sediment particle density (kg / m³) was obtained from field sampling and laboratory testing.

[0093] The viscosity of water (Pa·s) is obtained from the same steps described above.

[0094] Gas viscosity (Pa·s) was calculated based on the analysis of gas components from the core of the pressure-maintaining drilling and the formation temperature and pressure conditions.

[0095] It should be noted that, in some embodiments, the construction of the spatial morphological constraint of the pit may include the following steps: using seabed acoustic detection data to identify the planar distribution range of the pit; and extracting the pit projection area as the spatial morphological constraint of the pit based on the pit depth obtained by combining the planar distribution range with the third quantization relationship.

[0096] For example, in some specific implementations, based on the establishment of a vertical depth prediction model for the pits, the planar distribution range of the pits is identified using seabed acoustic detection data, thereby achieving a preliminary characterization of their lateral scale and providing geometric constraints for constructing the spatial morphology of the pits.

[0097] Submarine pockmarks typically appear as circular or elliptical depressions in their topography, and their boundaries can be identified using a high-resolution digital terrain model (DTM) acquired by a multibeam echo sounder. The projected area of ​​the pockmarks can then be extracted by combining side-scan sonar or shallow seismic profile data. Or the equivalent diameter, used as a representative parameter of its lateral dimension. This parameter is related to the vertical depth obtained through the above steps. Together, they form the basic spatial framework of the pit, reflecting the extent of deep fluid leakage on the seabed surface and providing a spatial integration domain for subsequent flux calculations.

[0098] Step S500: Based on the spatial morphology constraints of the pits and the pre-collected average methane flux, the total methane flux of the target area is simulated. It should be noted that the spatial morphological constraints of the pits include the projected area of ​​the pits. In some embodiments, step S500 may include the following steps: based on geological constraints, perform spatial representativeness analysis through preset observation points to predict the average methane flux; simulate the total methane flux of the target area based on the product of the average methane flux and the projected area of ​​the pits.

[0099] For example, in some specific embodiments, based on the preliminary clarification of the spatial morphology of the Makeng area, the embodiments of the present invention combine in-situ methane observation data to make a rough estimate of the methane release flux in the Makeng area, providing basic parameters for assessing its role in the regional carbon cycle.

[0100] In practice, methane concentration or flux monitoring devices are deployed in the center, edges, and surrounding areas of Makeng to obtain local release intensity. By conducting spatial representativeness analysis on a limited number of monitoring points, the average methane flux in the Makeng area can be estimated. (Unit: mol / (m²·s)). The projected area of ​​the pit. For the spatial integration range, the total methane flux can be approximately expressed as:

[0101] It should be noted that this invention is based on the geological phenomenon that there is a physical relationship between deep gas-bearing fluid seepage and the morphology of seabed pits. Taking the establishment of a quantitative relationship between seepage process, sediment instability and geomorphological response as the starting point, this invention achieves the coordinated prediction of the vertical depth, lateral range and seepage duration of pits by coupling gas-water two-phase seepage calculation, judgment of critical conditions for sand flow and time integration inversion.

[0102] The key technology of this invention lies in utilizing the relationship between seepage rate and sediment liquefaction threshold, combined with deep gas reservoir parameters inverted from seismic data and seafloor topography data, to construct a predictive process from deep fluid release to shallow pockmark formation, providing a basis for pockmark formation analysis and seepage history reconstruction. It should be noted that if other forms of seepage equations (e.g., considering non-Darcy flows), different sediment instability criteria (e.g., based on shear strength or pore pressure ratio), alternative time inversion methods (e.g., based on age data fitting or numerical inversion), or seafloor topographic data from different sources (e.g., side-scan sonar, AUV micro-topography, etc.) are used to identify pockmark boundaries, as long as the basic idea remains the same—constraining pockmark size and seepage time through the coupling relationship between fluid seepage dynamics and shallow sedimentary response—then its technical essence is consistent with this invention and should be considered an equivalent implementation of this invention.

[0103] 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.

[0104] First, it should be noted that current research on the formation and evolution of submarine pockmarks is mainly qualitative, focusing on static characteristics, and lacks a clear explanation of the relevant mechanisms. This makes it difficult to accurately characterize the intrinsic relationship between fluid kinematics (such as leakage time and migration paths) and dynamics (such as pressure transmission and fluid-sediment interactions) during pockmark formation. This limits a deeper understanding of the formation mechanisms of pockmarks and restricts the ability to invert their historical evolution and predict future development trends. Consequently, it fails to meet practical needs such as high-precision resource target area selection, dynamic risk assessment of submarine engineering, and semi-quantitative estimation of cold seep carbon flux. Therefore, there is an urgent need to develop a quantitative prediction method that integrates deep geological parameters (such as hydrate layer and free gas layer thickness), seepage dynamics models, and critical conditions for shallow sedimentary material desertification. This method would effectively link the formation and evolution of pockmarks with deep fluid activity processes, enabling a shift from "phenomenon identification" to "process reconstruction and parameter prediction," providing more reliable technical support for research on marine natural gas hydrate leakage systems and related engineering applications.

[0105] Specifically, existing technologies either focus on simulating local phenomena under controlled conditions or are limited to passive monitoring of shallow topography, generally lacking the ability to dynamically correlate deep gas-bearing fluid activity with shallow pit responses through physical mechanisms, especially lacking quantitative characterization of the temporal information of the seepage process. To address these shortcomings, this invention aims to provide a quantitative analysis method that integrates geological modeling, seepage dynamics, and sediment instability criteria. This method constructs a physical model coupling gas-water two-phase seepage and critical conditions for fluidized bed formation based on hydrate layer and free gas distribution inverted from actual seismic data, quantitatively estimating the vertical development depth of pits; it defines the lateral extent of pits using multibeam echo acoustic topographic data, and dynamically reconstructs the formation history of pits by inverting the duration of seepage activity through time integration; further, using the projected area of ​​the pit as the integration domain, combined with methane flux observations at a limited number of points, it conducts a semi-quantitative assessment of the regional carbon release potential. Figure 3 The diagram shown is a schematic representation of a geological model where gas seepage from a hydrated area forms submarine pits.

[0106] In some specific implementations, such as Figure 4 As shown, this invention systematically integrates multi-source data such as seismic, drilling, well logging, experimental testing, temperature, and pressure data to construct a quantitative prediction chain from "geological model → seepage rate → leakage time → pit size". The steps are logically rigorous, with preceding steps providing necessary input parameters for subsequent steps. The data sources and engineering acquisition methods for key parameters in each formula are clearly defined, ensuring the operability of the technical process and the reliability of the results.

[0107] The embodiments and technical methods of this invention provide important technical means for understanding the evolution of seafloor geomorphology under the background of deep gas-bearing fluids, and also provide important basis for inverting the potential of deep oil and gas resources through the characteristics of seafloor pockmarks, which has good theoretical and applied value. The main technical steps, specific content and implementation process of this invention are described below.

[0108] Step 1: Identification and Model Building of Key Geological Interfaces This step aims to construct a high-precision three-dimensional geological model, providing an accurate geometric framework and basic parameters for subsequent quantitative calculations. An accurate depth domain model is crucial for eliminating the influence of seismic wave velocity differences and achieving "depth positioning," providing an important data foundation for subsequent depth identification of the bottom boundary of the hydrate stability domain (BGHSZ), seepage simulation, and prediction of pit size.

[0109] Step 2: Calculation of the flow rate of gas-containing fluids: This step aims to establish a dynamic mathematical model of the gas leakage process based on the geological model provided in the previous steps, quantitatively describing the instantaneous rate of free gas migration from the bottom boundary of the hydrate stability zone to the seabed. This rate is a key dynamic parameter connecting the deep gas source and the shallow geomorphological response (pits), providing crucial dynamic parameters for calculating the leakage time of gas-bearing fluids and assessing the drag force of deep fluid activity on sediments.

[0110] Step 3: Estimation of the time it takes for fluid to seep to the seabed: Based on the seepage rate model established in step two, this step aims to calculate, using an integral method, the time *t* it takes for deep gas to seep from the bottom boundary of the hydrate stability domain until it finally reaches the seabed. This timescale is crucial for understanding the periodicity of seepage events and assessing the dynamic evolution rate of the hydrate system. Further, in a specific embodiment, the dating results of cold seep carbonates obtained through drilling can be combined, using their formation time as an age constraint, and superimposed with the migration time *t*, to deduce the initiation time of fluid seepage. It should be noted that because there is a time difference between carbonate precipitation and the initial fluid seepage, the derived result is an approximation. Figure 5 The diagram shown illustrates the variation of gas leakage time with gas column height under different simulation conditions (D: gas layer thickness; h). g (Air column height).

[0111] Step 4: Prediction of the size of the seabed pits: Building upon the geological model, seepage velocity, and temporal evolution information established in the first three steps, this step aims to establish a quantitative relationship model between the depth of submarine pit formation and deep gas reservoir parameters. By linking the deep seepage process with the shallow geomorphological response, this provides a direct theoretical basis and practical means for predicting deep oil and gas resource potential through submarine pit characteristics. Figure 6 The diagram shown is a schematic representation of a numerical model example of the gas leakage process from the hydrate zone to the shallow layer. Figure 7 The diagram shown illustrates the depth distribution of submarine pockmarks under different simulation conditions, reflecting the relationship between water flow impedance and effective formation stress. (D: gas layer thickness; h: stability region thickness; h...) pm : Depth of the pit; h w (Water column / water flow column height).

[0112] Step 5: Methane flux estimation under pitting constraint: Based on the preliminary clarification of the spatial morphology of Makeng, this step combines in-situ methane observation data to make a rough estimate of the methane release flux in the Makeng area, providing basic parameters for assessing its role in the regional carbon cycle.

[0113] In practice, methane concentration or flux monitoring devices are deployed in the center, edges, and surrounding areas of Makeng to obtain local release intensity. By conducting spatial representativeness analysis on a limited number of monitoring points, the average methane flux in the Makeng area can be estimated. (Unit: mol / (m²·s)). The projected area of ​​the pit. For the spatial integration range, the total methane flux can be approximately expressed as:

[0114] In a specific embodiment, such as Figure 8 As shown, taking a gas chimney and its associated submarine pits discovered during exploration in a basin as an example, the target pit is approximately circular, with its leakage intensity gradually decreasing from the center to both sides. By deploying methane monitoring probes at different spatial locations and within different annular zones of the pit, the methane flux corresponding to each annulus is calculated using an integral method. This flux value can be converted into carbon flux, and combined with the leakage duration or period obtained in step three, the cumulative methane release can be further estimated. Combined with the dating results of cold seep carbonate rocks, long-term methane release history can be cross-constrained, providing a preliminary basis for quantifying the contribution of the submarine cold seep system to regional carbon output.

[0115] It should be noted that the specific logical principles and formula calculation processes involved in the above steps are described in the aforementioned specific implementation method description, and will not be repeated here.

[0116] In summary, this invention addresses the uncertainty in predicting the shallow geomorphological response during fluid seepage in natural gas hydrate zones by proposing a method for predicting the size of submarine pockmarks. Based on geological parameters such as the thickness of the hydrate layer and free gas layer obtained through seismic inversion, and combined with a gas-water two-phase flow model and the critical conditions for shallow sedimentary fluid siltation, the method quantitatively predicts the duration of fluid seepage activity and the vertical depth of the pockmarks. Furthermore, it integrates submarine acoustic topographic data to identify the lateral distribution range of the pockmarks and preliminarily constructs their spatial morphology. Finally, it combines in-situ methane observations to provide a generalized estimate of the methane release flux in the pockmark area. This method achieves a cross-scale correlation between deep seepage drivers and shallow pockmark geomorphology, providing technical support for the identification of submarine cold seep systems and the assessment of carbon release potential.

[0117] Compared with the prior art, the embodiments of the present invention have at least the following beneficial effects: 1. A coupled seepage-siltation model was constructed to support the physical and quantitative prediction of the vertical depth of pits: By constructing a gas-water two-phase flow model and introducing critical conditions for shallow sedimentary fluid desertification, this method correlates deep geological parameters such as hydrate layer thickness, free gas layer thickness, sediment density, and fluid viscosity with shallow mechanical responses, enabling quantitative prediction of the vertical depth of submarine pockmarks. Based on Darcy flow and the effective stress principle and its correlation, this method establishes a functional relationship between pockmark depth and deep gas reservoir structure, fluid properties, and sedimentary mechanical parameters by solving the imbalance threshold between pore pressure and sediment desertification intensity driven by fluid overpressure. This overcomes the limitations of traditional methods that rely on geomorphological inversion or empirical analogy to effectively constrain pockmark depth.

[0118] 2. Effectively integrates the temporal evolution of seepage with geological record constraints, overcoming the limitations of static identification and dynamic reconstruction in pit research: By combining fluid leakage dynamics with time integration, the complete geological history of gas-bearing fluid release from deep layers to pockmark formation is inverted, and future pockmark development trends are predicted. This system uses leakage rate as a time variable and incorporates seismic stratigraphic framework and cold seep carbonate dating data to constrain the model, ensuring the inversion results are geologically interpretable. This effectively solves the problem that traditional methods can only identify existing pockmarks but cannot obtain temporal information about fluid activity.

[0119] 3. Based on the constraint of pit morphology, methane flux estimation is achieved, enabling the assessment of carbon release from seafloor cold seeps to move from qualitative to semi-quantitative methods: By combining seafloor acoustic topography data to identify the lateral distribution of pits and integrating it with in-situ methane observations, the regional average flux is estimated using the projected area of ​​the pits as the spatial integration domain. This method treats pits as basic units of methane release, delineates boundaries using multibeam or side-scan sonar, obtains the average flux based on observations at a limited number of points, and combines this with the identified pit size to achieve a rough quantification of methane fluid flux, providing effective geological constraints for regional carbon cycle research.

[0120] like Figure 9 As shown, this embodiment of the invention also provides a simulation device 900 for quantitatively predicting the gas leakage and pit formation evolution process in hydrate zones, which can implement the above-mentioned method. This device may include: The first module 910 is used to acquire three-dimensional seismic data of the target area and obtain thickness data through layer velocity modeling and time-depth conversion. The second module 920 is used to quantify the first quantitative relationship between the gas column growth rate and the seepage velocity based on thickness data and pre-acquired physical parameters using Darcy's law. The third module 930 is used to obtain a second quantification relationship for gas migration time by using an integral method based on the first quantification relationship; wherein, the second quantification relationship is used to provide geological constraints for the start time and duration of gas leakage; The fourth module 940 is used to combine the first quantitative relationship with the critical condition of quicksand formation to obtain the third quantitative relationship of pit depth, and then construct the spatial morphological constraint of pit. The fifth module, 950, is used to simulate the total methane flux in the target area based on the spatial morphology constraints of the pit and the pre-collected average methane flux.

[0121] 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.

[0122] 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.

[0123] 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.

[0124] like Figure 10 As shown, Figure 10 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.

[0125] 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.

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

[0127] 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.

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

[0129] 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.

[0130] 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.

[0131] The present invention provides a simulation method, apparatus, electronic device, storage medium, and program product for quantitatively predicting gas leakage and pit formation evolution in hydrate zones. It acquires three-dimensional seismic data of the target area, obtains thickness data through layer velocity modeling and time-depth conversion; based on the thickness data and pre-acquired physical parameters, it quantifies the first quantitative relationship between gas column growth rate and seepage velocity using Darcy's law; based on the first quantitative relationship, it obtains a second quantitative relationship for gas migration time through integration; the second quantitative relationship provides geological constraints for the onset time and duration of gas leakage; it combines the first quantitative relationship with the critical condition of quicksand formation to obtain a third quantitative relationship for pit depth, thereby constructing a spatial morphological constraint for pits; based on the spatial morphological constraint of pits and the pre-acquired average methane flux, it simulates the total methane flux of the target area. The present invention integrates three-dimensional seismic data, Darcy's law, and the critical condition of quicksand formation to establish a series of quantitative relationships between gas column growth rate, seepage velocity, gas migration time, and pit depth, achieving quantitative prediction of the physical constraints of deep fluid leakage processes and pit formation evolution. This invention provides the first quantitative correlation between deep fluid transport and shallow geomorphological response, overcoming the ambiguity of geophysical interpretation. Furthermore, by coupling multiphase flow, pore pressure evolution, and sediment mechanical response within a unified framework, this invention can accurately predict leakage rates, critical liquefaction depths, transport times, and pit sizes. Specifically, this invention can provide a reliable basis for free gas identification in resource exploration, marine engineering hazard assessment, and ecological environment methane flux estimation, significantly enhancing the scientific rigor and application value of pit research.

[0132] 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.

[0133] 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.

[0134] 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.

[0135] 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.

[0136] 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 simulation method for quantitatively predicting gas leakage and pit formation evolution in hydrate zones, characterized in that, The method includes the following steps: Acquire 3D seismic data of the target area, and obtain thickness data through layer velocity modeling and time-depth conversion; Based on the thickness data and pre-collected physical parameters, the first quantitative relationship between the gas column growth rate and the seepage velocity is quantified using Darcy's law. Based on the first quantification relationship, a second quantification relationship for gas migration time is obtained through integration; wherein, the second quantification relationship is used to provide geological constraints for the start time and duration of gas leakage. By combining the first quantitative relationship with the critical condition of quicksand formation, a third quantitative relationship for the depth of the pit is obtained, thereby constructing the spatial morphological constraint of the pit. The total methane flux in the target region is simulated based on the spatial morphological constraints of the pit and the pre-collected average methane flux.

2. The method according to claim 1, characterized in that, The thickness data includes the thickness of the free gas layer and the thickness of the stable domain from the bottom boundary of the hydrate stability domain to the seabed. The thickness data is obtained through layer velocity modeling and time-depth conversion, including the following steps: The root mean square velocity and two-way travel time of each key layer are extracted from the three-dimensional seismic data. Based on the root mean square velocity and the two-way travel time, the Dix formula is used to convert them into layer velocity, and a layer velocity model is established. Based on the layer velocity model, the time-domain seismic data volume is converted into the depth-domain seismic data volume to obtain the absolute depth of the key geological interface. The thickness of the stable region and the thickness of the free gas layer are obtained by quantifying the difference between the absolute depths of each of the key geological interfaces.

3. The method according to claim 1, characterized in that, The thickness data includes the thickness of the free gas layer and the thickness of the stable domain from the bottom boundary of the hydrate stable domain to the seabed. The physical parameters include water phase density, gas phase density, gravitational acceleration, gas viscosity, gas column height, water viscosity, water flow column height, and absolute permeability of the sediment. The gas column height is a time-varying variable, and the water flow column height is determined by the difference between the stable domain thickness and the gas column height. The first quantitative relationship between the gas column growth rate and the seepage velocity, based on the thickness data and the pre-collected physical parameters, using Darcy's law, includes the following steps: The fourth quantitative relationship of the driving pressure difference is obtained by quantifying the product of the difference between the water phase density and the gas phase density with the gravitational acceleration and the thickness of the free gas layer. The fifth quantitative relationship of gas column resistance is obtained by quantifying the product of the ratio of gas viscosity and absolute permeability of the sediment with gas column height and seepage velocity. The sixth quantitative relationship of water column resistance is obtained by quantifying the product of the ratio of the viscosity of the water and the absolute permeability of the sediment with the height of the water column and the seepage velocity. The fifth and sixth quantitative relationships are combined to obtain the seventh quantitative relationship for total flow resistance; Based on the relationship that the driving force equals the total flow resistance under steady seepage conditions, the fourth quantitative relationship is combined with the seventh quantitative relationship to obtain the first quantitative relationship of the seepage velocity. Wherein, since the height of the air column changes with time, the growth rate of the air column is equal to the seepage velocity.

4. The method according to claim 1, characterized in that, The process of obtaining a second quantification relationship for gas transport time based on the first quantification relationship through an integration method includes the following steps: Based on the first quantization relation, the method of separation of variables is used for integration to obtain the eighth quantization relation; The expression for the first quantization relationship is: ； The expression for the eighth quantization relation is: ； In the formula, Indicates the seepage velocity; Indicates the rate of air column growth. Represents the differential symbol; and For thickness data, Indicates the thickness of the free gas layer. This indicates the thickness of the stability domain from the bottom boundary of the hydrate stability domain to the seabed. , , , , , , For physical parameters, Indicates the absolute permeability of the sediment. This indicates the density of the aqueous phase. Indicates gas phase density, Represents gravitational acceleration. Indicates gas viscosity, Indicates the height of the air column. Indicates the viscosity of water. The height of the water column; Based on the eighth quantization relation, the left and right boundaries are integrated using a preset integration limit, and then the integration results of the left and right boundaries are combined and transformed to obtain the second quantization relation of the gas transport time.

5. The method according to claim 1, characterized in that, The process of combining the first quantitative relationship with the critical condition for quicksand formation to obtain the third quantitative relationship for the depth of the pit includes the following steps: The critical conditions for quicksand formation are constructed based on the seepage velocity, water viscosity, relative water permeability, absolute permeability of sediment, sediment particle density, water phase density, and gravitational acceleration. Substituting the first quantitative relationship into the critical condition of quicksand formation, and then transforming it to obtain the third quantitative relationship of the pit depth; The expression for the critical condition of sandification is as follows: ； The expression for the first quantization relation is: ； The expression for the third quantization relation is: ； In the formula, Indicates the seepage velocity; and For thickness data, Indicates the thickness of the free gas layer. This indicates the thickness of the stability domain from the bottom boundary of the hydrate stability domain to the seabed. , , , , , , , , For physical parameters, Indicates the absolute permeability of the sediment. Indicates the relative permeability of water. Indicates the density of sediment particles. This indicates the density of the aqueous phase. Indicates gas phase density, Represents gravitational acceleration. Indicates gas viscosity, Indicates the height of the air column. Indicates the viscosity of water. The height of the water column. ; Indicates the depth of the pit.

6. The method according to claim 1, characterized in that, The construction of the pit space shape constraint includes the following steps: The planar distribution range of the pits was identified using underwater acoustic detection data; The depth of the pit, obtained by combining the planar distribution range with the third quantization relationship, is used to extract the projected area of ​​the pit as the spatial morphological constraint of the pit.

7. The method according to any one of claims 1 to 6, characterized in that, The spatial morphological constraint of the pits includes the projected area of ​​the pits. The process of simulating the total methane flux of the target region based on the spatial morphological constraint of the pits and the pre-collected average methane flux includes the following steps: Based on the geological constraints, the average methane flux is predicted by performing spatial representativeness analysis through preset observation points. The total methane flux in the target region is simulated based on the product of the average methane flux and the projected area of ​​the pit.

8. A simulation device for quantitatively predicting gas leakage and pit formation evolution in hydrate zones, characterized in that, The device includes: The first module is used to acquire three-dimensional seismic data of the target area and obtain thickness data through layer velocity modeling and time-depth conversion. The second module is used to quantify the first quantitative relationship between the gas column growth rate and the seepage velocity based on the thickness data and the pre-collected physical parameters using Darcy's law. The third module is used to obtain a second quantification relationship for gas migration time based on the first quantification relationship through an integral method; wherein the second quantification relationship is used to provide geological constraints for the start time and duration of gas leakage. The fourth module is used to combine the first quantitative relationship with the critical condition of quicksand formation to obtain the third quantitative relationship of pit depth, and then construct the spatial morphological constraint of pit. The fifth module is used to simulate the total methane flux of the target area based on the spatial morphological constraints of the pit and the pre-collected average methane flux.

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, comprising a computer program, characterized in that, When the computer program is executed by a processor, it implements the method of any one of claims 1 to 7.

Images