A method for determining the multi-field coupling exploitation dynamics of a hydrate cap shallow gas reservoir considering sand particle migration

Abstract

The application discloses a kind of hydrate cap shallow gas reservoir multi-field coupling exploitation dynamic determination methods considering sand migration.The method of the application comprises the following steps: according to the geological data of actual sea area shallow gas reservoir, establish three-dimensional mechanism model;Based on three-dimensional mechanism model, by introducing sand erosion, deposition reaction and hydrate decomposition kinetics, establish multi-phase multi-component numerical model;Production well is set in three-dimensional mechanism model, and according to target horizon, determine perforation horizon, determine reasonable production scheme, carry out CH4 pressure reduction production;According to the result of hydrate cap shallow gas reservoir depressurization exploitation multi-field coupling numerical simulation considering sand migration, obtain the dynamic data of gas production, water production and sand production in CH4 exploitation period, and comprehensively analyze the spatial distribution characteristics of three-dimensional mechanism model whole and wellbore reservoir multi-physical field and its time-space evolution law with time, that is, the multi-field coupling exploitation dynamic of hydrate cap shallow gas reservoir considering sand migration is determined.

Description

Technical Field

[0001] This invention belongs to the technical field of shallow gas reservoir development and sand production dynamic prediction in marine areas, and relates to a multi-field coupled dynamic measurement method for shallow gas reservoirs with hydrate caprock that takes into account sand grain migration. Background Technology

[0002] Shallow offshore gas reservoirs are typically characterized by loosely cemented formations and a high risk of sand production, often overlain by a thermodynamically stable hydrate caprock. During depressurization extraction, intense hydrodynamic scouring can induce formation particle stripping and migration, leading to wellbore blockage and damage to surface equipment. Existing numerical simulation methods primarily focus on the gas-water two-phase flow in loose sandstone, lacking consideration of the dynamic coupled migration characteristics of solid sand particles with the fluid under the overlying hydrate caprock. Furthermore, they struggle to integrate experimental criticality indicators for specific sites, resulting in significant biases in predicting production performance and risks. Therefore, developing a predictive method that integrates experimental data, couples the particle migration process influenced by the hydrate caprock, and reveals multi-physics evolution characteristics is crucial for the safe and efficient development of shallow offshore gas reservoirs with hydrate caprocks. Summary of the Invention

[0003] The purpose of this invention is to provide a method for dynamic determination of multi-field coupled exploitation of shallow gas reservoirs with hydrate caprock that takes into account sand grain migration.

[0004] The present invention provides a method for dynamic determination of multi-field coupled exploitation of shallow gas reservoirs with hydrate caprock considering sand grain migration, comprising the following steps: (1) establishing a three-dimensional mechanism model based on geological data of actual shallow gas reservoirs in the sea area; (2) Based on the three-dimensional mechanism model, by introducing sand grain erosion, sedimentation reaction and hydrate decomposition kinetics, a multiphase and multicomponent numerical model is established to realize multi-field coupled numerical simulation of depressurization mining of shallow gas reservoirs with hydrate caprock considering sand grain migration. (3) Set up production wells in the three-dimensional mechanism model, determine the perforation layer according to the target layer, determine a reasonable production plan, and carry out depressurization production of CH4; (4) Based on the results of the multi-field coupled numerical simulation of depressurization development of shallow gas reservoirs with hydrate caprock considering sand migration as described in step (2), obtain the dynamic data of gas production, water production and sand production during the CH4 development cycle, and comprehensively analyze the overall three-dimensional mechanism model and the spatial distribution characteristics of multi-physics fields of the well perimeter reservoir and their spatiotemporal evolution over time, that is, determine the multi-field coupled development dynamics of shallow gas reservoirs with hydrate caprock considering sand migration.

[0005] In the above method, in step (1), the geological data includes at least one basic parameter among the following: location, thickness, permeability, porosity, hydrate saturation, gas saturation, and water saturation of the overlying rock layer, CH4 hydrate layer, shallow gas layer, and underlying rock layer.

[0006] In the above method, in step (2), the multiphase multicomponent numerical model includes a formation property dynamic correction model, a hydrate reaction kinetic model, a sand grain erosion and deposition reaction model, a mass conservation equation, an energy conservation equation, and a gas-water relative permeability model. The following assumptions are made to establish the multiphase, multicomponent numerical model: ① The model considers four phases and six components. The four phases include gas phase, water phase, hydrate phase and solid sand phase; the six components include CH4, water, flowing sand, NaCl, natural gas hydrate and detachable sand. ②The hydrates and detachable sands in the reservoir are not flowing. Only the gas and water two-phase flow is considered. The seepage conforms to Darcy's law. ③ The flowing sand and NaCl exist in the aqueous phase as dispersed components, which do not affect the thermal conductivity of the aqueous phase. They migrate with the aqueous phase and are accompanied by sensible heat transfer. ④ Descalable sand and mobile sand have the same thermodynamic properties as the rock skeleton; the energy conservation equation considers the latent heat of phase change of hydrate decomposition and ignores the mechanical reaction heat of sand particles. ⑤ Consider the influence of salinity on hydrate phase equilibrium and the dynamic regulation of reservoir permeability by pore structure evolution.

[0007] In the above method, the dynamic correction model of formation properties is designed to address the evolution of formation microstructure caused by hydrate decomposition, detachable sand erosion, and mobile sand deposition. The dynamic correction model of formation properties consists of three parts: a hierarchical porosity evolution model, a saturation equation, and a dynamic permeability adjustment model, in order to achieve dynamic coupling between physicochemical reactions and seepage characteristics. (a) Hierarchical porosity evolution model Equation (1); In equation (1), The absolute porosity of the formation is the pore space in which gas, water, and hydrates coexist without considering rock deformation and damage. Total porosity refers to the total pore space in which gas, water, hydrates, and detachable sand coexist. The occupancy rate of detachable sand is the volume percentage of detachable sand. Fluid porosity that truly allows fluid flow Then the space occupied by hydrates needs to be further deducted: Equation (2); In equation (2), For fluid porosity, The absolute porosity of the formation. In absolute porosity space hydrate saturation in; (b) Saturation equation Equation (3); In equation (3), S g , S w , S H These represent the saturation levels of the gas phase, aqueous phase, and hydrate phase, respectively, all based on absolute pore space. Calculations show that changes in the content of detachable sand alter the absolute pore space. The magnitude of the saturation affects the global saturation, while the sum of the saturations always equals 1. (c) Dynamic model of penetration rate absolute permeability of formation With fluid porosity The changes in fluid porosity are dynamically adjusted, and the formation permeability under different fluid porosities is calculated using the Carmen–Kozeny formula: Equation (4); In equation (4), k The absolute permeability of the formation. The initial absolute permeability of the formation without hydrates, in mD; For fluid porosity, The initial fluid porosity without hydrates, without units; n This is the penetration rate decline index, and can specifically be set to a value of 5.

[0008] In this invention, rock deformation and failure are not considered; the absolute porosity of the formation is defined. The pore space is defined as the coexistence of gas, water, and hydrate phases; while to characterize the modification of the stratigraphic framework by erosion and deposition processes, total porosity is introduced. As a spatial reference (i.e., the total pore space where gas, water, hydrates, and detachable sand coexist); detachable sand is considered as the solid sand phase component occupying the total pore space, which is the part of the sandstone skeleton that can be transformed into flowing sand by fluid erosion. Before being disturbed, it does not have fluidity, and its volume proportion is defined as the occupancy rate. The transformation between detachable sand and mobile sand is the core of describing the impact of sand erosion and deposition on reservoir modification. The model realizes the dynamic transformation of sand components between different phases through erosion and deposition reactions. Erosion drives detachable sand to transform from solid sand phase components into mobile sand components existing in the water phase; conversely, deposition causes the mobile sand grains to be transformed back into detachable sand components and remain in the reservoir. Due to the erosion and deposition reactions of sand components, the occupancy rate of detachable sand is affected. The absolute porosity of the formation will change as the simulation progresses. The dynamic evolution logic is shown in equations (1) and (2) above.

[0009] In the above method, the hydrate reaction kinetic model is the reaction equation for the formation or decomposition of natural gas hydrates, as shown in equation (5): Equation (5); In the formula, n is the hydration number of the natural gas hydrate molecule, which can specifically be 5.75; The enthalpy of the reaction for the formation of natural gas hydrates.

[0010] In this invention, when the formation and decomposition of hydrates both follow the Kim-Bishnoi model, taking the decomposition reaction as an example, the hydrate decomposition rate is shown in equation (6): Equation (6); In equation (6), The molar concentration of natural gas hydrate, in mol / m³; The intrinsic frequency of the natural gas hydrate decomposition reaction, mol / (m 2 (kPa·d) Let m be the solid specific surface area of ​​the hydrate. 2 / m 3 ; is the activation energy of the reaction, J / mol; R is the gas constant, with a value of 8.314 J / (mol·K); To balance the pressure, kPa; K is the equilibrium constant of the hydrate phase, which is a function of pressure p and temperature T, and its value can be determined from the K value table.

[0011] In the above method, the sand grain erosion and deposition reaction model: on the one hand, regarding the destruction of the original stratigraphic framework and the process of particle stripping, the model defines an erosion reaction from detachable sand to flowing sand; physically, this process is controlled by the shear stress of the fluid on the rock surface, and particles will only detach and be initiated when the drag force generated by the pore fluid velocity overcomes the cementation strength and internal friction between particles; to this end, a velocity-dependent reaction rate control term is introduced in the numerical simulation settings, and a critical gas initiation velocity is set as the threshold for the reaction to occur; when the local gas velocity in the grid is lower than this critical value, the erosion reaction rate is zero, and the formation remains stable; once the velocity exceeds this threshold, the reaction is initiated immediately, and the reaction rate increases linearly with the increase of velocity; The liquefaction rate of mud and sand is shown in equation (7): Equation (7); In equation (7), The sand liquefaction rate; For detachable sand that remains solid per unit volume, m 3 / m 3 ; It is a detachable sand erosion reaction factor; The gas velocity can be obtained from Darcy's flow equation, in m / s; The critical flow velocity for sand discharge is obtained by fitting experimental data, in m / s; Equivalent to this in the general chemical reaction rate expression in numerical simulation is shown in equation (8): Equation (8); In equation (8), n s The liquefaction rate of mud and sand is expressed in mol / (m²). 3 ·d); The erosion reaction rate constant of detachable sand; For detachable sand in the total pore space Medium molar density, mol / m 3 ; The velocity index is used to characterize the nonlinear response relationship between erosion rate and flow velocity. It is dimensionless and obtained by fitting. for The reference gas flow rate for the value; On the other hand, suspended particles accompanying fluid transport, when flowing through narrow pore throats or tortuous channels, are affected by physical capture and adsorption, and will separate from the sand-water mixture (aqueous phase) and revert to detachable sand, i.e., particle deposition reaction. This process is mainly controlled by the geometry of the pore throat and can be indirectly characterized by permeability-related kinetic reactions: the lower the permeability (corresponding to the narrower the pore throat), the faster the deposition reaction rate. This mechanism follows the deep bed filtration law, and the deposition rate depends on the real-time molar fraction of sand particles in the sand-water mixture and the permeability-related reaction rate constant, which is equivalent to the general chemical reaction rate expression for numerical simulation, as shown in equation (9): Equation (9); In equation (9), n d The deposition rate is expressed in mol / (m²). 3 ·d), The deposition reaction rate constant; This refers to the water phase saturation. This represents the mole fraction of sand particles in the aqueous phase. The average molar mass of the aqueous phase is expressed in kg / mol. This is a multiple of the reaction rate controlled by formation permeability.

[0012] In this invention, the sand grain erosion and deposition reaction model also considers the secondary migration of sedimentary particles, that is, the detachable sand generated after deposition can be resuspended and transformed into flowing sand in the erosion reaction; similar to skeletal erosion, the captured sand grains are not permanently fixed, and when the subsequent gas reaches the critical flow rate, the detachable sand generated by deposition will be washed away again. The simulator calculates the occupancy rate of detachable sand in the total pore space in real time. This allows for the correction of formation permeability, thereby quantitatively characterizing reservoir blockage and reduced productivity caused by sand migration.

[0013] In the above method, the mass conservation equation and energy conservation equation are as follows: Based on the principle of mass conservation, and taking into account the effects of fluid flow, chemical reactions, sand particle mechanical reactions, and source-sink terms, the mass conservation equations for each component in the model, including CH4, water, NaCl, flowing sand, hydrates, and detachable sand, are as follows: CH4: Equation (10); water: Equation (11); NaCl: Equation (12); Flowing sand: Equation (13); Natural gas hydrates: Equation (14); Descalable sand: Equation (15); In the formula, 、 , , Densities of the gas phase, aqueous phase, hydrate phase, and solid sand phase, in kg / m³. 3 ; , , where represents the flow velocity of the gas phase and the water phase, in m / s; , , The saturation level represents the gas phase, aqueous phase, and hydrate phase; it has no unit. For detachable sand in the total pore space Volume fraction within, without units; , For gas and aqueous phases, the source and sink terms are kg / (m³). 3 ·s); , , The values ​​represent the mass fractions of water, NaCl, and flowing sand components in the aqueous phase; they are unitless. , , , The formation rates of CH4, water, flowing sand, and natural gas hydrates. The erosion reaction rate of the detachable sand component is expressed in kg / (m²). 3 ·s); Based on the law of conservation of energy, and considering the effects of heat conduction, chemical reactions, and source-sink terms, the energy conduction and transformation process in a unit volume of formation system can be described by the following equation: Equation (16); In this equation, and These are effective specific heat capacity and effective thermal conductivity, respectively. Equation (17); Equation (18); In the formula, T is the fluid temperature, K; 、 , , Densities of gas phase, aqueous phase, hydrate phase, and sandstone, in kg / m³ 3 ; , , , Specific heat capacity of gas phase, aqueous phase, hydrate phase, and sandstone, J / (kg·K); , , , , where W / (m·K) represents the thermal conductivity of the gas phase, water phase, hydrate phase, and sandstone. It is the heat released when one mole of natural gas hydrate is formed, in J / mol; It is the amount of natural gas hydrate formed per unit volume per unit time, expressed in mol / (m³). 3 ·s).

[0014] In the above method, the air-water relative permeability model is specifically as follows: The relative permeability curves of the gas-water two-phase flow in the porous medium in the model were calculated using the Brooks-Corey model, i.e., equations (19) and (20): Equation (19); Equation (20); In equations (19) and (20), , These are the relative permeabilities of the water phase and the gas phase in a two-phase flow, respectively, without units; , These represent the saturation levels of water and gas in the mobile phase, respectively; no units are provided. , These represent the bound water saturation and residual gas saturation in the formation pores, respectively, with specific values ​​of 0.1 and 0.03. , These are the fitting parameters, with specific values ​​of 5 and 3 respectively.

[0015] In the above method, in step (3), the production well includes a vertical well or a horizontal well; The target layer includes a hydrate layer or a shallow gas layer.

[0016] In the above method, in step (4), the multiphysics field includes at least one of temperature field, pressure field, porosity, permeability, hydrate saturation and gas saturation.

[0017] The present invention has the following beneficial effects: This invention provides a reliable numerical simulation tool for depressurization development of loose sandstone gas reservoirs with hydrate caprocks. It can quantitatively predict the production dynamics of gas and water, as well as the risk of sand production. Through analysis, it can further optimize the perforation site and production regime, reduce the wear of wellbore and surface equipment, and ensure the safe and efficient exploitation of gas reservoirs. It has clear engineering guidance significance and field application value. Attached Figure Description

[0018] Figure 1 This is a three-dimensional radial mechanism model and well location diagram of Example 1; Figure 2 The dynamic curves of gas production, water production, water-gas ratio, and sand production during the eight years of depressurization mining in Example 1 are shown. Figure 3 The pressure P, temperature T, and fluid porosity during the 8-year depressurization mining period of Example 1 are as follows: Horizontal penetration rate k h hydrate saturation S H Gas saturation S g Water saturation S w The spatiotemporal evolution; Figure 4 The hydrate saturation within a 10 m radial direction around the well during 8 years of depressurization extraction in Example 1. S H Fluid porosity Horizontal penetration rate k hMass fraction of sand particles in the aqueous phase w s The spatiotemporal evolution of the situation. Detailed Implementation

[0019] Unless otherwise specified, the experimental methods used in the following examples are conventional methods.

[0020] Unless otherwise specified, all materials and reagents used in the following examples are commercially available.

[0021] It should be noted that the following detailed description is illustrative and intended to provide further explanation of the invention. Unless otherwise specified, 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.

[0022] Example 1 (1) Referring to the geological data of a hydrate-bearing loose sandstone gas reservoir at a certain site in the Qiongdongnan Basin, a three-dimensional radial mechanism model of appropriate size was established based on the relative positions and effective thicknesses of the overlying layer, hydrate layer, shallow gas layer, and underlying layer. The three-dimensional radial mechanism model and well location settings are as follows: Figure 1 As shown, the radial (R direction) range of the model is 0~1500 m, the vertical (Z direction) thickness is 252 m, and the top boundary of the model is the seabed mudline. To balance computational accuracy and efficiency, the model is divided into 200 grids radially (R direction), with the grid step size increasing logarithmically from the model center towards R. This provides high-resolution near-wellbore flow field simulation for subsequent production wells in the model center, accurately capturing the sand particle transport dynamics under high near-wellbore flow conditions. The vertical (K direction) is further subdivided into 293 simulation layers. The overlying layer is 180 m thick (i.e., the distance from the seafloor mudline to the top surface of the hydrate layer), and is divided into 36 simulation layers at 5 m thicknesses per layer. The hydrate layer and shallow gas layer are 10 m and 14 m thick, respectively, and are divided into 100 and 140 simulation layers at 0.1 m thicknesses per layer, respectively. The underlying layer is 48 m thick (i.e., 48 m below the bottom boundary of the shallow gas layer, equivalent to adding a heat-transferring reservoir boundary, while considering the presence of bottom water to simulate its pressure support and heat conduction effect on depressurization extraction), and is divided into 16 simulation layers at 3 m thicknesses per layer.

[0023] The geological parameters used in the model are shown in Table 1: Table 1. Geological parameters used in the model

[0024] (2) Based on the aforementioned three-dimensional mechanism model, a multiphase, multi-component numerical model is established by introducing sand grain erosion, deposition reaction, and hydrate decomposition kinetics to achieve coupled simulation of seepage field, thermal field, chemical field, and sand grain transport. The basic assumptions of the model are as follows: ① The model considers four phases and six components. The four phases are: gas phase, water phase, hydrate phase, and solid sand phase; the six components are: CH4, water, flowing sand, NaCl, natural gas hydrate, and detachable sand. ②The hydrates and detachable sands in the reservoir are not flowing. Only the gas-water two-phase flow is considered. The seepage conforms to Darcy's law. ③ The flowing sand and NaCl exist in the aqueous phase as dispersed components, which do not affect the thermal conductivity of the aqueous phase. They migrate with the aqueous phase and are accompanied by sensible heat transfer. ④ Descalable sand and mobile sand have the same thermodynamic properties as the rock skeleton; the latent heat of phase change due to hydrate decomposition is considered in the energy equation, while the heat of mechanical reaction of sand particles is ignored. ⑤ Consider the influence of salinity on hydrate phase equilibrium and the dynamic regulation of reservoir permeability by pore structure evolution.

[0025] Based on the above assumptions, dynamic correction models for formation properties, hydrate reaction kinetics, sand erosion and deposition reaction, mass conservation equations, energy conservation equations, and gas-water relative permeability models were constructed. Details are as follows: Dynamic correction model for formation properties: To address the evolution of formation microstructure caused by hydrate decomposition, detachable sand erosion, and sand grain deposition, a dynamic correction model for formation properties was constructed, consisting of three parts: a hierarchical porosity evolution model, a saturation equation, and a dynamic permeability adjustment model, in order to achieve dynamic coupling between physicochemical reactions and seepage characteristics.

[0026] (a) Hierarchical porosity evolution model Defining absolute porosity of formation without considering rock deformation and failure. It refers to the pore space where gas, water, and hydrates coexist. To characterize the alteration of the stratigraphic framework by erosion and deposition processes, total porosity is introduced. As a spatial reference (i.e., the total pore space where gas, water, hydrates, and detachable sand coexist), detachable sand is considered as the solid sand phase component occupying the total pore space. It is the part of the sandstone skeleton that can be transformed into flowing sand by fluid erosion. Before being disturbed, it does not possess fluidity, and its volume proportion is defined as the occupancy rate. δ sThe transformation between detachable sand and mobile sand is the core of describing the impact of sand erosion and deposition on reservoir modification. The model realizes the dynamic transformation of sand components between different phases through erosion and deposition reactions: erosion drives detachable sand to transform from solid sand phase components into sand grain components existing in the water phase; conversely, deposition causes mobile sand to transform back into detachable sand components and remain in the reservoir.

[0027] Due to the erosion and deposition reactions of sand components, the occupancy rate of detachable sand is affected. δ s The absolute porosity will change as the simulation progresses. The dynamic evolution logic is as follows: Equation (1); Fluid porosity that truly allows fluid flow Then the space occupied by hydrates needs to be further deducted: Equation (2); In the formula, S H In absolute porosity space The hydrate saturation in the solution.

[0028] (b) Saturation equation The saturation of the gas phase, aqueous phase, and hydrate phase are all defined based on the absolute pore space. Calculate and satisfy the following equation: Equation (3); Under this definition, changes in the content of detachable sand alter the absolute pore space. The magnitude of the saturation affects the overall saturation, while the sum of the saturation always equals 1.

[0029] (c) Dynamic model of penetration rate absolute permeability of formation k With fluid porosity The changes in fluid porosity are dynamically adjusted, and the formation permeability under different fluid porosities is calculated using the Carmen–Kozeny formula: Equation (4); In the formula, k 0 represents the initial absolute permeability of the formation without hydrates, in mD; The initial fluid porosity without hydrates, without units; n This is the penetration rate decline index, with a value of 5.

[0030] Hydrate reaction kinetic model: The reaction equation for the formation or decomposition of natural gas hydrates can be described as equation (5): Equation (5); In the formula, n is the hydration number of the natural gas hydrate molecule, which is taken as 5.75 in this paper; Δ H The enthalpy of the reaction for the formation of natural gas hydrates.

[0031] Assuming that the formation and decomposition of hydrates both follow the Kim-Bishnoi model, taking the decomposition reaction as an example, the hydrate decomposition rate is shown in equation (6): Equation (6); In the formula, c H The molar concentration of natural gas hydrate, in mol / m³; The intrinsic frequency of the natural gas hydrate decomposition reaction, mol / (m 2 ·kPa·d); A HS The solid specific surface area of ​​the hydrate. m 2 / m 3 ;Δ E is the activation energy of the reaction, J / mol; R is the gas constant, with a value of 8.314 J / (mol·K); p e To balance the pressure, kPa; K K is the equilibrium constant of the hydrate phase, which is a function of pressure p and temperature T, and its value can be determined from the K value table.

[0032] Sand erosion and deposition response model: On the one hand, regarding the destruction of the original stratigraphic framework and the process of particle stripping, the model defines an erosion reaction from detachable sand to flowing sand. Physically, this process is controlled by the shear stress of the fluid on the rock surface. Particles will only detach and be initiated when the drag force generated by the pore fluid velocity overcomes the cementation strength and internal friction between particles. Therefore, a velocity-dependent reaction rate control term is introduced in the numerical simulation settings, and the critical gas initiation velocity is set as the threshold for the reaction to occur. When the local gas velocity in the grid is lower than this critical value, the erosion reaction rate is zero and the formation remains stable; once the velocity exceeds this threshold, the reaction is initiated immediately, and the reaction rate increases linearly with the increase of velocity. The mud-sand liquefaction rate can be obtained by equation (7): Equation (7); In the formula, m si For detachable sand that remains solid per unit volume, m 3 / m 3 ;α m It is a detachable sand erosion reaction factor; u gThe gas velocity can be obtained from Darcy's flow equation, in m / s; u ci The critical flow velocity for sand discharge is obtained by fitting experimental data, in m / s.

[0033] Equivalent to this, the general chemical reaction rate expression in numerical simulation is shown in equation (8). Equation (8); In the formula, n s The liquefaction rate of mud and sand is expressed in mol / (m²). 3 ·d); Let be the rate constant for the erosion reaction of detachable sand, (s·mol). -1 c s For detachable sand in the total pore space Medium molar density, mol / m 3 ; α The velocity index is used to characterize the nonlinear response relationship between erosion rate and flow velocity. It is dimensionless and obtained by fitting. u ref The reference gas velocity for the value of ω.

[0034] On the other hand, suspended particles accompanying fluid transport, when flowing through narrow pore throats or tortuous channels, are affected by physical capture and adsorption, and will separate from the sand-water mixture (aqueous phase) and revert to detachable sand, i.e., particle deposition reaction. This process is mainly controlled by the geometry of the pore throat and can be indirectly characterized by permeability-related kinetic reactions: the lower the permeability (corresponding to the narrower the pore throat), the faster the deposition reaction rate. This mechanism follows the deep bed filtration law, and the deposition rate depends on the real-time molar fraction of sand particles in the sand-water mixture and the permeability-related reaction rate constant, which is equivalent to the general chemical reaction rate expression for numerical simulation, as shown in Equation (9).

[0035] Equation (9); In the formula, n d The deposition rate is expressed in mol / (m²). 3 ·d); β S is the deposition reaction rate constant; w This refers to the water phase saturation. x fines This represents the mole fraction of sand particles in the aqueous phase. M w The average molar mass of the aqueous phase is expressed in kg / mol. f k This is a multiple of the reaction rate controlled by formation permeability.

[0036] Furthermore, the model also considers the secondary transport of sedimentary particles, namely the erosion reaction in which detachable sand generated after deposition can be resuspended and transformed into flowing sand. Similar to skeletal erosion, the captured sand particles are not permanently fixed. When subsequent gas reaches a critical flow rate, the detachable sand generated after deposition will be washed away again.

[0037] The simulator calculates the occupancy rate of detachable sand in the total pore space in real time. δ s This allows for the correction of formation permeability, thereby quantitatively characterizing reservoir blockage and reduced productivity caused by sand migration.

[0038] mass conservation equation: Based on the principle of mass conservation, and taking into account the effects of fluid flow, chemical reactions, sand particle mechanical reactions, and source-sink terms, the mass conservation equations for each component in the model, including CH4, water, NaCl, flowing sand, hydrates, and detachable sand, are as follows: CH4: Equation (10); water: Equation (11); NaCl: Equation (12); Flowing sand: Equation (13); Natural gas hydrates: Equation (14); Descalable sand: Equation (15); In the formula, ρ g ρ w , ρ H , ρ s Densities of the gas phase, aqueous phase, hydrate phase, and solid sand phase, in kg / m³. 3 ;v g v w , where represents the flow velocity of the gas phase and the water phase, in m / s; S g , S w , S H The saturation level represents the gas phase, aqueous phase, and hydrate phase; it has no unit. δ s For detachable sand in the total pore space Volume fraction within, without units;q g , q w For gas and aqueous phase source and sink terms, kg / (m³) 3 ·s); , , The values ​​represent the mass fractions of water, NaCl, and flowing sand components in the aqueous phase; they are unitless. , , , The formation rates of CH4, water, flowing sand, and natural gas hydrates. The erosion reaction rate of the detachable sand component is expressed in kg / (m²). 3 ·s).

[0039] Energy conservation equation: Based on the law of conservation of energy, and considering the effects of heat conduction, chemical reactions, and source-sink terms, the energy conduction and transformation process in a unit volume of formation system can be described by the following equation: Equation (16); In this equation, and These are effective specific heat capacity and effective thermal conductivity, respectively.

[0040] Equation (17); Equation (18); In the formula, T is the fluid temperature, K; 、 , , Densities of gas phase, aqueous phase, hydrate phase, and sandstone, in kg / m³ 3 ; C g 、C w , C H , C s λ represents the specific heat capacity of the gas phase, aqueous phase, hydrate phase, and sandstone, in J / (kg·K); g 、 λ w , λ H , λ s The thermal conductivity of the gas phase, aqueous phase, hydrate phase, and sandstone is given in W / (m·K); Δ H H It is the heat released when one mole of natural gas hydrate is formed, in J / mol; n HIt is the amount of natural gas hydrate formed per unit volume per unit time, expressed in mol / (m³). 3 ·s).

[0041] Air-water relative permeability model: The relative permeability curves of the gas-water two-phase flow in the porous medium in the model were calculated using the Brooks-Corey model, i.e., equations (19) and (20): Equation (19); Equation (20); In the formula, k rw , k rg These are the relative permeabilities of the water phase and the gas phase in a two-phase flow, respectively, without units; , These represent the saturation levels of water and gas in the mobile phase, respectively; no units are provided. , These represent the bound water saturation and residual gas saturation in the formation pores, respectively, with values ​​of 0.1 and 0.03. n w , n g The fitting parameters are set to 5 and 3 respectively.

[0042] (3) In this embodiment, a vertical well is selected as the production well, such as Figure 1 As shown, the model was set at the center R=0 m of the three-dimensional mechanism model, and full perforation was carried out on the shallow gas layer. The perforation section was 14 m long and the wellbore radius was 0.1 m. Depressurization mining was carried out for 8 years.

[0043] (4) Based on the numerical simulation results in (2), obtain the dynamic data of gas production, water production, and sand production during the mining cycle, such as Figure 2 As shown, where, Figure 2 The four sub-graphs (a), (b), (c), and (d) are the dynamic curves of gas production, water production, water-to-gas ratio, and sand production, respectively.

[0044] Depend on Figure 2 As shown in Figure (a), the gas production rate peaks rapidly in the early stages of extraction and then declines exponentially until it stabilizes at approximately 26 × 10⁻⁶. 4 m 3 / d, the cumulative gas production increases rapidly in the first year and then slows down. After 8 years of extraction, the cumulative gas production is approximately 181,924 × 10⁶. 4 m 3 The gas production is quite considerable. For example... Figure 2As shown in (b) and (c), the water production rate and the water-to-gas ratio show a continuous upward trend throughout the entire extraction cycle, exhibiting a significant step-like increasing characteristic; by the end of the 8th year, the water production rate slowly increases to approximately 4.4 m. 3 / d, the production water-to-gas ratio simultaneously increased to approximately 0.118 m 3 / 10 4 m 3 Due to the continuously increasing water production rate, the cumulative water production curve showed an accelerated upward trend, and the cumulative water production over 8 years eventually reached approximately 4400 m³. 3 It is noteworthy that the water-to-gas ratio exhibits a significant step-like upward trend throughout the entire extraction cycle. This is because the initial gas saturation is high, resulting in a correspondingly low relative permeability of the water phase; as shallow gas is rapidly extracted, the gas saturation decreases, and the water phase flow capacity gradually recovers, thus leading to the step-like increase in the water-to-gas ratio. Figure 2 As shown in Figure (d), the sand production rate increased rapidly in the early stage of mining, reaching a peak of about 0.54 kg / d in about 0.5 years, and then dropped rapidly and significantly. In the middle and late stages of mining, the sand production rate remained at a low level, and the final cumulative sand production reached about 410 kg by the end of the 8th year.

[0045] Pressure P, temperature T, and fluid porosity during 8 years of depressurization mining Horizontal penetration rate hydrate saturation Gas saturation Water saturation The vertical evolution of the physical field is as follows Figure 3 As shown, the production well is located on the far left of the field map, and full-section perforation was performed on the shallow gas layer with a depth range of 190-204 mbsf. It can be seen that during production, the pressure drop funnel rapidly spreads throughout the entire reservoir, resulting in a significant and uniform decrease in overall pressure. The overall hydrate saturation does not change significantly, thus the temperature field remains generally stable, with only a slight throttling and cooling effect caused by the pressure drop and expansion of the shallow gas layer. The morphology and thickness of the shallow gas layer do not change significantly, but the gas saturation in the upper part of the shallow gas layer continuously decreases with production.

[0046] Furthermore, the hydrate saturation within a 10 m radial radius around the well was analyzed during the 8 years of depressurization mining. Fluid porosity Horizontal penetration rate Mass fraction of sand particles in the aqueous phase Evolution over time ( Figure 4The vertical extent of the field map includes the hydrate layer and the shallow gas layer. It can be seen that throughout the entire 8-year mining period, due to the relatively small production pressure differential, the main hydrate layer maintained a stable state with high saturation. Hydrate decomposition was mainly limited to the area near the interface between the hydrate layer and the shallow gas layer, exhibiting significant radial spatial differences: in the region above the interface near the wellbore perforation point, the hydrate underwent complete decomposition due to the most severe pressure drop effect; however, as the radial distance increased, the pressure drop gradually weakened, and the degree of decomposition decreased accordingly, manifested as a slight decrease in hydrate saturation at the far end of this interface. With the local dissociation of hydrates, the pore space occupied by hydrates was released, leading to a significant local increase in fluid porosity and horizontal permeability in the corresponding area. Observation of the water saturation profile reveals a local increase in water saturation near the interface between the hydrate layer and the underlying gas layer, and a decrease in water saturation from the far well to the near well as mining progresses. This is because the decomposition of hydrates produces a large amount of water, which migrates towards the perforation point during mining. Furthermore, as depressurization extraction continued (0-8 years), the gas saturation at the bottom of the underlying gas layer gradually decreased while the water saturation increased significantly, indicating that during the gas extraction process, not only did water from the decomposition of hydrates above seep downwards, but deep bottom water also intruded upwards.

[0047] Based on the above analysis of the production dynamics and spatiotemporal evolution of key physical fields during depressurization mining, a certain reference basis can be provided for the design of development schemes for actual shallow gas reservoirs with hydrate caprocks.

[0048] The above description is merely a preferred embodiment of the present invention and is not intended to limit the invention. Various modifications and variations can be made to the present invention by those skilled in the art. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A method for dynamic determination of multi-field coupled development in shallow gas reservoirs with hydrate caprocks considering sand grain migration, characterized in that, The steps include: (1) Establishing a three-dimensional mechanism model based on the geological data of the actual shallow gas reservoir in the sea area; (2) Based on the three-dimensional mechanism model, by introducing sand grain erosion, sedimentation reaction and hydrate decomposition kinetics, a multiphase and multicomponent numerical model is established to realize multi-field coupled numerical simulation of depressurization mining of shallow gas reservoirs with hydrate caprock considering sand grain migration. (3) Set up production wells in the three-dimensional mechanism model, determine the perforation layer according to the target layer (hydrate layer or shallow gas layer), determine a reasonable production plan, and carry out depressurization production of CH4; (4) Based on the results of the multi-field coupled numerical simulation of depressurization development of shallow gas reservoirs with hydrate caprock considering sand migration as described in step (2), obtain the dynamic data of gas production, water production and sand production during the CH4 development cycle, and comprehensively analyze the overall three-dimensional mechanism model and the spatial distribution characteristics of multi-physics fields of the well perimeter reservoir and their spatiotemporal evolution over time, that is, determine the multi-field coupled development dynamics of shallow gas reservoirs with hydrate caprock considering sand migration.

2. The method according to claim 1, characterized in that, In step (1), the geological data includes at least one of the following basic parameters: location, thickness, permeability, porosity, hydrate saturation, gas saturation, and water saturation of the overlying rock layer, CH4 hydrate layer, shallow gas layer, and underlying rock layer.

3. The method according to claim 1 or 2, characterized in that, In step (2), the multiphase and multicomponent numerical model includes a formation property dynamic correction model, a hydrate reaction kinetic model, a sand grain erosion and deposition reaction model, a mass conservation equation, an energy conservation equation, and a gas-water relative permeability model. The following assumptions are made to establish the multiphase, multicomponent numerical model: ① The model considers four phases and six components. The four phases include gas phase, water phase, hydrate phase and solid sand phase; the six components include CH4, water, flowing sand, NaCl, natural gas hydrate and detachable sand. ②The hydrates and detachable sands in the reservoir are not flowing. Only the gas and water two-phase flow is considered. The seepage conforms to Darcy's law. ③ The flowing sand and NaCl exist in the aqueous phase as dispersed components, which do not affect the thermal conductivity of the aqueous phase. They migrate with the aqueous phase and are accompanied by sensible heat transfer. ④ Descalable sand and mobile sand have the same thermodynamic properties as the rock skeleton; the energy conservation equation considers the latent heat of phase change of hydrate decomposition and ignores the mechanical reaction heat of sand particles. ⑤ Consider the influence of salinity on hydrate phase equilibrium and the dynamic regulation of reservoir permeability by pore structure evolution.

4. The method according to claim 3, characterized in that, The formation property dynamic correction model is designed to address the evolution of formation microstructure caused by hydrate decomposition, detachable sand erosion, and mobile sand deposition. The dynamic correction model consists of three parts: a hierarchical porosity evolution model, a saturation equation, and a permeability dynamic adjustment model, in order to achieve dynamic coupling between physicochemical reactions and seepage characteristics. (a) Hierarchical porosity evolution model Equation (1); In equation (1), The absolute porosity of the formation is the pore space in which gas, water, and hydrates coexist without considering rock deformation and damage. Total porosity refers to the total pore space in which gas, water, hydrates, and detachable sand coexist. The occupancy rate of detachable sand is the volume percentage of detachable sand. Equation (2); In equation (2), For fluid porosity, This represents the absolute porosity of the formation. In absolute pore space hydrate saturation in; (b) Saturation equation Equation (3); In equation (3), S g , S w , S H These represent the saturation levels of the gas phase, aqueous phase, and hydrate phase, respectively, all based on absolute pore space. Calculations show that changes in the content of detachable sand alter the absolute pore space. The magnitude of the saturation affects the global saturation, while the sum of the saturations always equals 1. (c) Dynamic model of penetration rate Equation (4); In equation (4), k The absolute permeability of the formation. The initial absolute permeability of the formation without hydrates, in mD; For fluid porosity, The initial fluid porosity without hydrates; n This is the penetration rate decline index.

5. The method according to claim 3, characterized in that, The hydrate reaction kinetic model is the reaction equation for the formation or decomposition of natural gas hydrates, as shown in equation (5): Equation (5); In the formula, n is the hydration number of a natural gas hydrate molecule; The enthalpy of the reaction for the formation of natural gas hydrates.

6. The method according to claim 3, characterized in that, The sand grain erosion and deposition reaction model: The liquefaction rate of mud and sand is shown in equation (7). Equation (7); In equation (7), The sand liquefaction rate; For detachable sand that remains solid per unit volume, m 3 / m 3 ; It is a detachable sand erosion reaction factor; The gas velocity can be obtained from Darcy's flow equation, in m / s; The critical flow velocity for sand discharge is obtained by fitting experimental data, in m / s; Equivalent to this in the general chemical reaction rate expression in numerical simulation is shown in equation (8): Equation (8); In equation (8), n s The liquefaction rate of mud and sand is expressed in mol / (m). 3 ·d); The erosion reaction rate constant of detachable sand; For detachable sand in the total pore space Medium molar density, mol / m 3 ; The velocity index is used to characterize the nonlinear response relationship between erosion rate and flow velocity. It is dimensionless and obtained by fitting. for The reference gas flow rate for the value; Equation (9); In equation (9), n d The deposition rate is expressed in mol / (m²). 3 ·d); The deposition reaction rate constant; This refers to the water phase saturation. This represents the mole fraction of sand particles in the aqueous phase. The average molar mass of the aqueous phase is expressed in kg / mol. This is a multiple of the reaction rate controlled by formation permeability.

7. The method according to claim 3, characterized in that, The mass conservation equation and energy conservation equation are as follows: Based on the principle of mass conservation, and taking into account the effects of fluid flow, chemical reactions, sand particle mechanical reactions, and source-sink terms, the mass conservation equations for each component in the model, including CH4, water, NaCl, flowing sand, hydrates, and detachable sand, are as follows: CH4: Equation (10); water: Equation (11); NaCl: Equation (12); Flowing sand: Equation (13); Natural gas hydrates: Equation (14); Descalable sand: Equation (15); In the formula, 、 , , Densities of the gas phase, aqueous phase, hydrate phase, and solid sand phase, in kg / m³. 3 ; , , where represents the flow velocity of the gas phase and the water phase, in m / s; , , The saturation level represents the gas phase, aqueous phase, and hydrate phase; it has no unit. For detachable sand in the total pore space Volume fraction within, without units; , For gas and aqueous phase source and sink terms, kg / (m³) 3 ·s); , , The values ​​represent the mass fractions of water, NaCl, and flowing sand components in the aqueous phase; they are unitless. , , , The formation rates of CH4, water, flowing sand, and natural gas hydrates. The erosion reaction rate of the detachable sand component is expressed in kg / (m²). 3 ·s); Based on the law of conservation of energy, and considering the effects of heat conduction, chemical reactions, and source-sink terms, the energy conduction and transformation process in a unit volume of formation system can be described by the following equation: Equation (16); In this equation, and These are effective specific heat capacity and effective thermal conductivity, respectively. Equation (17); Equation (18); In the formula, T is the fluid temperature, K; 、 , , Densities of gas phase, aqueous phase, hydrate phase, and sandstone, in kg / m³ 3 ; , , , , where J / (kg·K) represents the specific heat capacity of gas phase, aqueous phase, hydrate phase, and sandstone. , , , , where W / (m·K) represents the thermal conductivity of the gas phase, water phase, hydrate phase, and sandstone. It is the heat released when one mole of natural gas hydrate is formed, in J / mol; It is the amount of natural gas hydrate formed per unit volume per unit time, expressed in mol / (m³). 3 ·s).

8. The method according to claim 3, characterized in that, The specific air-water relative permeability model is as follows: The relative permeability curves of the gas-water two-phase flow in the porous medium in the model were calculated using the Brooks-Corey model, i.e., equations (19) and (20): Equation (19); Equation (20); In equations (19) and (20), , These are the relative permeabilities of the water phase and the gas phase in a two-phase flow, respectively, without units; , These represent the saturation levels of water and gas in the mobile phase, respectively; no units are provided. , These represent the bound water saturation and residual gas saturation in the formation pores, respectively. , These are the fitting parameters.

9. The method according to claim 1 or 2, characterized in that, In step (3), the production well includes a vertical well or a horizontal well; The target layer includes a hydrate layer or a shallow gas layer.

10. The method according to claim 1 or 2, characterized in that, In step (4), the multi-physics field includes at least one of temperature field, pressure field, porosity, permeability, hydrate saturation and gas saturation.

Images