Particle method for simulating sedimentary layer landslide process

Abstract

The invention discloses a particle method for simulating a sedimentary layer landslide process, and belongs to the fields of ocean engineering and port and navigation. Aiming at the problems that a numerical simulation method cannot simulate large deformation motion of a sedimentary layer and an existing SPH method cannot simulate the deformation process of the sedimentary layer containing hydrate decomposition in the prior art, the method provided by the invention establishes a gas-liquid-soil three-phase particle model containing hydrate decomposition based on the SPH method; a Drucker-Prager elastic-plastic constitutive model is adopted to simulate a sedimentary layer soil body, a hydrate decomposition kinetic equation is adopted to simulate a sedimentary layer hydrate decomposition process, the method not only can simulate a sedimentary layer gas-liquid-soil transportation process, but also can simulate a sedimentary layer phase change heat-flow coupling problem, and the method has the advantages of being simple in structure and convenient to operate. An effective numerical analysis tool is provided for predicting the mechanical property and the motion law of the sedimentary layer, and the method has important engineering significance.

Description

Technical Field

[0001] The invention belongs to the field of marine engineering and port navigation, and specifically is a particle method for simulating the landslide process of a sediment layer. Background Art

[0002] In the field of marine engineering, traditional numerical simulation of marine sediment landslides is usually based on the computational domain of mesh discretization, and the classic calculation method is generally the finite element method (FEM). However, as a mesh-based method, the finite element method is more accurate in calculating the results before the deformation of soil and water, but after the large deformation landslide of the submarine sediment layer occurs, it is often limited by mesh distortion, resulting in results with large errors. As an emerging meshless method, the smoothed particle hydrodynamics method (SPH) has a natural advantage in solving large deformation problems and is suitable for the numerical simulation of marine engineering problems.

[0003] Natural gas hydrates are widely present in submarine sediments in nature. The decomposition characteristics of hydrates are a key factor in studying the migration process of natural gas hydrate sediments. Many researchers at home and abroad have conducted experimental studies on the decomposition behavior of hydrate samples. Studies have found that hydrates will decompose under certain temperature and pressure conditions, which will cause significant changes in the mechanical properties of marine sediments and may cause landslides in originally stable soil layers. At the same time, hydrate decomposition will cause changes in temperature and pressure in marine sediments, which will further affect the deformation of the submarine sediments. Therefore, it is very necessary to establish a numerical method for simulating landslides in sediments containing hydrate decomposition.

[0004] In the field of SPH, the study of hydrate-containing soils has just started, and no scholar has yet proposed a particle method to simulate landslides in hydrate-containing sediments. Summary of the invention

[0005] In view of the problems existing in the prior art that "numerical simulation methods cannot simulate large deformation movement of sedimentary layers" and "existing SPH methods cannot simulate the deformation process of sedimentary layers containing hydrate decomposition", the present invention proposes a particle method for simulating the landslide process of sedimentary layers, which comprises the following steps:

[0006] 1) At the initial moment, according to the scale of the submarine sediment layer to be simulated, the SPH particles are used to discretize the submarine sediment layer containing hydrate, a numerical model of the submarine sediment layer to be simulated is established and boundary conditions are set; wherein the discrete submarine sediment layer particles include boundary particles and internal particles;

[0007] 2) Iteratively calculate and update the seabed sediment particles step by step according to the time step, and output the current landslide morphology of the seabed sediment according to the current information of the updated particles at each time step;

[0008] The step 2) specifically includes:

[0009] 2.1) For boundary particles, use Shepard interpolation to solve the pressure, stress, velocity and temperature of the particles;

[0010] 2.2) For internal particles, each discretized internal particle carries the physical information of water, hydrate, methane and soil. The continuity equation in the Navier-Stokes equations is used to solve the saturation increment and pressure information of water, hydrate and methane at the particle position;

[0011] 2.3) Use the momentum equation to calculate parameters such as stress and viscous drag force to obtain the velocity increment and displacement increment of the internal particles;

[0012] 2.4) Use the energy equation to calculate the temperature change in the marine sediment layer and obtain the temperature increment of the internal particles;

[0013] 2.5) updating the soil mechanics parameters of the sedimentary layer according to the saturation increments of each component obtained in step 2.2), wherein the soil mechanics parameters include shear modulus G, bulk modulus K, cohesion coefficient c and dilatancy angle ψ;

[0014] 2.6) For each iterative solution time step, after obtaining the velocity increment, temperature increment, pressure information of the internal particles and the saturation increment of each component according to steps 2.2), 2.3) and 2.4), the corresponding physical information of the internal particles is updated, and the soil mechanics parameters of the sedimentary layer are calculated according to step 2.5) based on the hydrate saturation; then the iterative solution of the next time step is entered; wherein, each time step outputs the updated particle information of all discretized submarine sedimentary layer particles, and outputs the current landslide morphology of the sedimentary layer based on the particle information.

[0015] The method of the present invention adds the Kim-Bishnoi hydrate decomposition kinetic equation to the migration process of marine sediment layers, and adds the hydrate decomposition process to the soil movement; at the same time, the migration process of water and gas during the migration of submarine sediment layers is considered, and the pressure and saturation changes of water and gas in the pores of the sediment layers are calculated, which is conducive to better simulating the influence of hydrate decomposition on the landslide process of submarine sediment layers.

[0016] The method of the present invention uses a formula to relate the mechanical properties of marine sedimentary soil to hydrate concentration. This treatment makes the sedimentary soil closer to the actual seabed soil, which is conducive to better simulating the migration process of the hydrate sediment during the hydrate decomposition process.

[0017] Compared with the prior art, the beneficial effects of the present invention include: in response to the problems existing in the prior art that "numerical simulation methods cannot simulate large deformation movement of sedimentary layers" and "existing SPH methods cannot simulate the deformation process of sedimentary layers containing hydrate decomposition", the present application proposes a particle method and system for simulating the landslide process of sedimentary layers, establishes a gas-liquid-soil three-phase SPH particle model, adopts the Drucker-Prager elastoplastic constitutive model to simulate the seabed soil, and adds a hydrate decomposition kinetic model to the soil model. This method can not only simulate the gas-liquid-soil transport process of the sedimentary layer, but also simulate the phase change heat-flow coupling problem of the sedimentary layer, providing an effective numerical analysis tool for predicting the mechanical properties and movement laws of the sedimentary layer. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] Figure 1 Schematic diagram of the components of hydrate-containing sedimentary layers;

[0019] Figure 2 Hydrate decomposition phase equilibrium curve;

[0020] Figure 3 Schematic diagram of soil landslide model during hydrate decomposition;

[0021] Figure 4 At t = 8s, the effective stress cloud map, displacement cloud map, air pressure cloud map, and water pressure cloud map of the sedimentary layer;

[0022] Figure 5 Cumulative plastic strain distribution of sedimentary layer under different pressure difference conditions. DETAILED DESCRIPTION

[0023] The present invention is further described and illustrated below in conjunction with specific embodiments. The embodiments are merely exemplary of the present disclosure and do not define the scope of limitation. The technical features of each embodiment of the present invention may be combined accordingly without conflicting with each other.

[0024] like Figure 1 As shown in the figure, it is a schematic diagram of the components of the hydrate-containing sediment layer studied in the present invention. The hydrate-containing sediment layer mainly includes natural gas, water, hydrate, soil and other components. In the subsequent parameter symbols, the subscript g represents natural gas, the subscript w represents water, the subscript h represents hydrate, and the subscript s represents soil. Among them, hydrate can be decomposed into natural gas and water in the form of solid hydrate, which leads to the change of phase state. The phase equilibrium is related to temperature and pressure. Figure 2 The phase equilibrium curve of hydrate is shown.

[0025] The core idea of ​​the particle method for simulating the landslide process of the sedimentary layer proposed in the present invention is to establish a gas-liquid-soil three-phase particle model containing hydrate decomposition based on the SPH method, adopt the Drucker-Prager elastoplastic constitutive model to simulate the sedimentary layer soil, adopt the hydrate decomposition kinetic equation to simulate the hydrate decomposition process of the sedimentary layer, and then simulate the gas-liquid-soil transport process of the sedimentary layer to realize the simulation of the landslide process of the sedimentary layer.

[0026] The method of the present invention is specifically implemented according to time steps, and the method comprises the following steps:

[0027] 1) At the initial moment, according to the scale of the submarine sediment layer to be simulated, the SPH particles are used to discretize the submarine sediment layer containing hydrates, a numerical model of the submarine sediment layer to be simulated is established, and boundary conditions are set. The discrete submarine sediment layer particles include boundary particles and internal particles. Among them, the numerical model of the submarine sediment layer to be simulated and the boundary conditions are established as follows: mathematical modeling of the submarine sediment layer to be simulated at the initial moment, and setting the pressure, temperature, and velocity data of the submarine sediment layer to be simulated at the initial moment.

[0028] 2) The seabed sediment particles are iteratively calculated and updated step by step according to the time step. At each time step, the current landslide morphology of the seabed sediment is output according to the current information after the particles are updated.

[0029] In a specific embodiment of the present invention, the step 2) specifically includes:

[0030] 2.1) For boundary particles, Shepard interpolation is used to solve the pressure, stress, velocity and temperature of the particles. The Shepard interpolation method is a well-known method in the art. This step obtains the pressure, stress, velocity and temperature of the boundary particles, which is to update the boundary conditions for the next time step.

[0031] 2.2) For internal particles, each discretized internal particle carries the physical information of water, hydrate, methane and soil. The continuity equation in the Navier-Stokes equations is used to solve the saturation increment and pressure information of water, hydrate and methane at the particle position;

[0032] The step 2.2) specifically includes:

[0033] The continuity equation in the Navier-Stokes equations includes:

[0034] Continuity equation for hydrate:

[0035]

[0036] Continuity equation for water:

[0037]

[0038] Continuity equation for methane gas:

[0039]

[0040] Where, t is time; is the divergence symbol; are the water production rate, hydrate decomposition rate and gas production rate respectively; S h , S w and S g are the saturations of hydrate, water and methane respectively; u h ,u w and u g are the velocities of hydrate, water and methane respectively; φ is the porosity of the marine sediment layer; ρ h ,ρ w and ρ g are the densities of hydrate, water, and methane gas, respectively;

[0041] The decomposition rate of hydrates is calculated according to the Kim-Bishnoi hydrate decomposition kinetic equation:

[0042]

[0043] In the formula, k d is the kinetic reaction rate, in mol / (m 2 ·Pa·s); M h is the molar mass of the hydrate; A s is the specific surface area of ​​hydrate decomposition per unit volume of sediment layer, in m -2 ;P g is the methane gas pressure, in Pa; P eq To balance the pressure;

[0044] The water and gas production rates are calculated using the mass conservation principle according to the corresponding hydrate decomposition chemical equations:

[0045]

[0046] Where M g and M w are the molar masses of methane and water, respectively;

[0047] The water and gas velocities are respectively based on the Darcy seepage velocity q w and q g The calculation formulas are:

[0048]

[0049] where μ w and μg are the dynamic viscosities of water and methane respectively; P w and P g are the pressures of water and methane gas, respectively. is the gradient operator, and are the water pressure and air pressure gradients; k is the absolute permeability of the sediment; k w and k g are the relative permeabilities of water and methane gas, respectively;

[0050] After arranging and simplifying, the saturation increment formulas for calculating water, hydrate and methane gas are:

[0051] The water saturation increment calculation formula is:

[0052]

[0053] Hydrate saturation increment formula:

[0054]

[0055] Because of the relationship Therefore, the methane saturation increment can be calculated by the following formula:

[0056]

[0057] in, is the water saturation increment, is the hydrate saturation increment, is the saturation increment of methane gas;

[0058] The increase in particle pressure is calculated using the continuity equation for gases:

[0059]

[0060] in, is the air pressure increment.

[0061] 2.3) Use the momentum equation to calculate the particle velocity increment through the soil stress and viscous drag force The calculation formula is as follows:

[0062]

[0063] In the formula, is the particle's velocity increment, σ s is the effective stress of the soil, is the effective stress divergence of soil, P w is the water pressure, ρ s is the density of soil, f dis the viscous drag force, f g is the acceleration due to gravity;

[0064] The calculation formula of particle displacement increment is:

[0065]

[0066] in, is the displacement increment of the particle, and u is the velocity of the particle.

[0067] 2.4) The energy equation used is as follows to calculate the temperature increment of the particles in the deposition layer:

[0068]

[0069] in, is the temperature increment of the particles in the sediment layer, is the divergence symbol, is the temperature gradient, C T is the total heat capacity of the ocean sediments, K T is the thermal conductivity of the marine sediment layer, Q h It is the energy produced by the decomposition of hydrates.

[0070] 2.5) updating the soil mechanics parameters of the sedimentary layer according to the saturation increments of each component obtained in step 2.2), wherein the soil mechanics parameters include shear modulus G, bulk modulus K, cohesion coefficient c and dilatancy angle ψ;

[0071] K=α K S h +K0

[0072] G=α G S h +G0

[0073]

[0074]

[0075] S h is the hydrate concentration; α K , α G , α c and β c is an adjustable coefficient; the superscript “0” indicates the material parameter value when the hydrate concentration is zero.

[0076] 2.6) For each iterative solution time step, after obtaining the velocity increment, temperature increment, pressure information of the internal particles and the saturation increment of each component according to steps 2.2), 2.3) and 2.4), the corresponding physical information of the internal particles is updated, and the soil mechanics parameters of the sedimentary layer are calculated according to step 2.5) based on the hydrate saturation; then the iterative solution of the next time step is entered; wherein, each time step outputs the updated particle information of all discretized submarine sedimentary layer particles, and outputs the current landslide morphology of the sedimentary layer based on the particle information.

[0077] The method of the present invention is applied in combination with a specific example. This example takes soil landslide during hydrate decomposition as an example to illustrate the technical solution of the present invention. The simulation includes the following steps:

[0078] 1) Initially, establish Figure 3 The hydrate soil model shown has an upper base length of 1m, a lower base length of 2m, and a height of 1m. The initial temperature T of the hydrate soil is set to 275.45K. The left boundary I is set as the decompression boundary, and the air pressure P g0 The pressure is constant at 3.15MPa, and the other boundaries Ⅱ, Ⅲ, and Ⅳ are all non-penetrating closed boundaries. The temperature of the four sides is constant, and the temperature is 275.45K. In terms of motion constraints, the left boundary Ⅰ and the lower boundary Ⅳ are no-slip boundary conditions, and the motion of the boundaries Ⅱ and Ⅲ is unconstrained. When performing numerical simulation calculations, only the effect of air pressure on hydrate decomposition is considered, and the effect of air pressure on soil movement is not considered;

[0079] 2) If Figure 1 As shown in the figure, the hydrate-containing sediment layer is composed of water, hydrate, methane and soil. After discretization, each sediment layer particle carries the physical information of water, hydrate, methane and soil. The generation rate of hydrate, water and methane gas is calculated by the hydrate decomposition kinetic equation according to the sediment layer temperature, air pressure and other information at the current moment. The calculation formula is as follows:

[0080]

[0081]

[0082]

[0083] Where M h 、M g and M w are the molar masses of hydrate, methane and water, respectively, M h =0.124 kg / mol, M g =0.016 kg / mol and M w=0.018 kg / mol. The continuity equation in the aforementioned Navier-Stokes equations is used to solve the saturation S of water, hydrate and methane at the particle position. w , S h and S g The current air pressure is obtained through the above state equation;

[0084] 3) If a particle carries soil information, the soil stress and strain parameters are calculated according to the above momentum equation based on the current air pressure, water pressure and hydrate saturation concentration;

[0085] 4) For all particles, the temperature change is calculated using the aforementioned energy equation based on the decomposition rate, total specific heat capacity, and thermal conductivity of the hydrate;

[0086] 5) The mechanical properties of the marine sediment layer will change according to the hydrate saturation. The current hydrate concentration is obtained through step 2), and the mechanical parameters of the current soil are calculated according to the above formula;

[0087] 6) If a particle is a boundary particle, use Shepard interpolation to solve the pressure, stress, velocity and temperature of the particle;

[0088] 7) For all fluid particles, after obtaining the particle velocity increment, temperature increment, pressure increment and saturation increment of each component according to steps 2), 3) and 4), the corresponding physical information of the particles is updated, and the soil mechanics parameters of the sedimentary layer are calculated according to the hydrate saturation through step 5), and the solution of the next time step is entered.

[0089] Figure 4 (a) in the figure gives the effective stress distribution diagram of the sedimentary layer at t=8s, where σ is the effective stress of the sedimentary layer. It can be seen from the figure that under the action of gravity and water pressure, the horizontal effective stress is the largest on the left side, and the sedimentary layer is subjected to the largest force on the left side and is most susceptible to deformation. Figure 4 (b) is the displacement cloud map of the sedimentary layer, Δx is the displacement of the sedimentary layer, the maximum displacement to the left is near the upper left corner of the sedimentary layer, and the maximum displacement to the right is on the right side of the sedimentary layer, which is consistent with the strain cloud map of the sedimentary layer. Figure 4 (c) and (d) are the cloud diagrams of air pressure and water pressure in the sedimentary layer. g is the air pressure, P w As shown in the figure, this method is a fully coupled method that can extract stress, displacement, air pressure, and water pressure in the seafloor sediment layer at any time.

[0090] Figure 5 The cumulative plastic strain distribution of the sedimentary layer under different pressure difference conditions is shown.dev is the accumulated plastic strain of the sedimentary layer, indicating the deformation of the sedimentary layer. (a), (b) and (c) are the strain diagrams after the initial pressure difference is 50kPa, 150kPa and 200kPa respectively, calculated to a stable state. Figure 5 It can be seen that when the pressure difference is small, the force generated by the pressure gradient is smaller than the cohesive force of the sediment layer and is insufficient to cause the sediment layer to deform. As the pressure difference gradually increases, the deformation near the upper left corner of the sediment layer becomes larger and larger, and a large deformation occurs when the pressure difference is 200 kPa.

[0091] The above-described embodiment only expresses a feasible implementation of the present invention, and its description is relatively specific and detailed, but it cannot be understood as limiting the scope of the present invention. It should be pointed out that for ordinary technicians in this field, several modifications and improvements can be made without departing from the concept of the present invention, which all belong to the protection scope of the present invention. Therefore, the protection scope of the present invention shall be based on the attached claims.

Claims

1. A particle method for simulating a sedimentary landslide process, characterized in that: The steps include: 1) At the initial moment, according to the scale of the submarine sediment layer to be simulated, the SPH particles are used to discretize the submarine sediment layer containing hydrate, a numerical model of the submarine sediment layer to be simulated is established and boundary conditions are set; wherein the discrete submarine sediment layer particles include boundary particles and internal particles; 2) Iteratively calculate and update the seabed sediment particles step by step according to the time step, and output the current landslide morphology of the seabed sediment according to the current information of the updated particles at each time step; The step 2) specifically includes: 2.1) For boundary particles, use Shepard interpolation to solve the pressure, stress, velocity and temperature of the particles; 2.2) For internal particles, each discretized internal particle carries the physical information of water, hydrate, methane and soil. The continuity equation in the Navier-Stokes equations is used to solve the saturation increment and pressure information of water, hydrate and methane at the particle position; 2.3) Use the momentum equation to calculate parameters such as stress and viscous drag force to obtain the velocity increment and displacement increment of the internal particles; 2.4) Use the energy equation to calculate the temperature change in the marine sediment layer and obtain the temperature increment of the internal particles; 2.5) updating the soil mechanics parameters of the sedimentary layer according to the saturation increments of each component obtained in step 2.2), wherein the soil mechanics parameters include shear modulus G, bulk modulus K, cohesion coefficient c and dilatancy angle ψ; 2.6) For each iterative solution time step, after obtaining the velocity increment, temperature increment, pressure information of the internal particles and the saturation increment of each component according to steps 2.2), 2.3) and 2.4), the corresponding physical information of the internal particles is updated, and the soil mechanics parameters of the sedimentary layer are calculated according to step 2.5) based on the hydrate saturation; then the iterative solution of the next time step is entered; wherein, each time step outputs the updated particle information of all discretized submarine sedimentary layer particles, and outputs the current landslide morphology of the sedimentary layer based on the particle information.

2. A particle method for simulating a sedimentary landslide process according to claim 1, characterized in that: In step 1), the numerical model of the seabed sediment layer to be simulated is established and boundary conditions are set, including: mathematically modeling the seabed sediment layer to be simulated at the initial moment, and setting the pressure, temperature, and velocity data of the seabed sediment layer to be simulated at the initial moment.

3. A particle method for simulating a sedimentary landslide process according to claim 1, characterized in that: The step 2.2) specifically includes: The continuity equation in the Navier-Stokes equations includes: Continuity equation for hydrate: Continuity equation for water: Continuity equation for methane gas: Where, t is time; is the divergence symbol; are the water production rate, hydrate decomposition rate and gas production rate respectively; S h , S w and S g are the saturations of hydrate, water and methane respectively; u h ,u w and u g are the velocities of hydrate, water and methane respectively; φ is the porosity of the marine sediment layer; ρ h ,ρ w and ρ g are the densities of hydrate, water, and methane gas, respectively; The decomposition rate of hydrates is calculated according to the Kim-Bishnoi hydrate decomposition kinetic equation: In the formula, k d is the kinetic reaction rate, in mol / (m 2 ·Pa·s); M h is the molar mass of the hydrate; A s is the specific surface area of ​​hydrate decomposition per unit volume of sediment layer, in m -2 ;P g is the methane gas pressure, in Pa; P eq To balance the pressure; The water and gas production rates are calculated using the mass conservation principle according to the corresponding hydrate decomposition chemical equations: Where M g and M w are the molar masses of methane and water, respectively; The water and gas velocities are respectively based on the Darcy seepage velocity q w and q g The calculation formulas are: where μ w and μ g are the dynamic viscosities of water and methane respectively; P w and P g are the pressures of water and methane gas, respectively. is the gradient operator, and are the water pressure and air pressure gradients; k is the absolute permeability of the sediment; k w and k g are the relative permeabilities of water and methane gas, respectively; After arranging and simplifying, the saturation increment formulas for calculating water, hydrate and methane gas are: The water saturation increment calculation formula is: Hydrate saturation increment formula: Because of the relationship Therefore, the methane saturation increment can be calculated by the following formula: in, is the water saturation increment, is the hydrate saturation increment, is the saturation increment of methane gas; The increase in particle pressure is calculated using the continuity equation for gases: in, is the air pressure increment.

4. A particle method for simulating a sedimentary landslide process according to claim 1, characterized in that: In step 2.3), the momentum equation is used to calculate the particle velocity increment through the soil stress and viscous drag force The calculation formula is as follows: In the formula, is the particle's velocity increment, σ s is the effective stress of the soil, is the effective stress divergence of soil, P w is the water pressure, ρ s is the density of soil, f d is the viscous drag force, f g is the acceleration due to gravity; The calculation formula of particle displacement increment is: in, is the displacement increment of the particle, and u is the velocity of the particle.

5. The particle method for simulating a sedimentary landslide process according to claim 1, characterized in that: The energy equation used in step 2.4) is as follows, which is used to calculate the temperature increment of the particles in the deposition layer: in, is the temperature increment of the particles in the sediment layer, is the divergence symbol, is the temperature gradient, C T is the total heat capacity of the ocean sediments, K T is the thermal conductivity of the marine sediment layer, Q h It is the energy produced by the decomposition of hydrates.

6. A particle method and system for simulating a sedimentary landslide process according to claim 1, characterized in that: In step 2.5), the following formulas are used to calculate the shear modulus G, bulk modulus K, cohesion coefficient c and dilatancy angle ψ: K=α K S h +K0 G=a G S h +G0 S h is the hydrate concentration; α K , α G , α c and β c is an adjustable coefficient; the superscript "0" indicates the material parameter value when the hydrate concentration is zero.

Images