A particle method for simulating the process of landslides in sedimentary layers

Abstract

The present application discloses a particle method for simulating the landslide process of a sediment layer, belonging to the fields of ocean engineering and port and waterway. Aiming at the problems of "numerical simulation methods being unable to simulate the large deformation movement of sediment layers" and "existing SPH methods being unable to simulate the deformation process of sediment layers containing hydrate decomposition" in the prior art, the method of the present invention establishes a gas-liquid-solid three-phase particle model for hydrate decomposition based on the SPH method, uses the Drucker-Prager elastoplastic constitutive model to simulate the sediment layer soil body, and uses the hydrate decomposition kinetic equation to simulate the hydrate decomposition process of the sediment layer. This method can not only simulate the transport process of gas-liquid-solid in the sediment layer, but also simulate the phase change heat-fluid coupling problem of the sediment layer, providing an effective numerical analysis tool for predicting the mechanical properties and movement laws of the sediment layer, and having important engineering significance.

Description

Technical Field

[0001] The present invention belongs to the fields of ocean engineering and port and waterway, and particularly relates to a particle method for simulating the process of sediment layer landslide. Background Art

[0002] In the field of ocean engineering, traditional numerical simulations of submarine sediment layer landslides are usually based on a computational domain discretized by grids, and the classical calculation method is generally the finite element method (FEM). However, as a grid-based method, although the finite element method gives relatively accurate calculation results before the deformation of soil and water, after a large deformation landslide occurs in the submarine sediment layer, it is often restricted by grid distortion and produces results with large errors. As an emerging meshless method, the Smoothed Particle Hydrodynamics (SPH) method has natural advantages in solving large deformation problems and is suitable for numerical simulations of ocean engineering problems.

[0003] Natural gas hydrates in nature widely exist in submarine sediment layers. The decomposition characteristics of hydrates are a key factor in studying the migration process of natural gas hydrate sediment layers. Many domestic and foreign scientific researchers have conducted experimental studies on the decomposition behavior of hydrate samples. It has been found through research that hydrates will decompose under certain temperature and pressure conditions, which will lead to significant changes in the mechanical properties of the ocean sediment layer, and may cause landslides in the originally stable soil layer. At the same time, the decomposition of hydrates will cause changes in the temperature and pressure within the ocean sediment layer, which will further affect the deformation of the submarine sediment layer soil. Therefore, it is very necessary to establish a numerical method for simulating the landslide of sediment layers containing hydrate decomposition.

[0004] In the field of SPH, the research on soils containing hydrates has just started, and no scholar has proposed a particle method for simulating the landslide of sediment layers containing hydrates yet. Summary of the Invention

[0005] Aiming at the problems of "numerical simulation methods cannot simulate the large deformation movement of sediment layers" and "existing SPH methods cannot simulate the deformation process of sediment layers containing hydrate decomposition" in the prior art, the present invention proposes a particle method for simulating the process of sediment layer landslide, and the method includes the following steps:

[0006] 1) At the initial moment, according to the scale of the submarine sediment layer to be simulated, use SPH particles to discretize the submarine sediment layer containing hydrates, establish a numerical model of the submarine sediment layer to be simulated and set boundary conditions; wherein, the discretized submarine sediment layer particles include boundary particles and internal particles;

[0007] 2) Perform step-by-step iterative calculations and updates on the particles in the seabed sediment layer at each time step. At each time step, output the current landslide morphology of the seabed sediment layer based on the current information of the updated particles.

[0008] Step 2) specifically includes:

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

[0010] 2.2) For internal particles, each discrete internal particle carries the physical information of water, hydrate, methane, and soil. Use the continuity equation in the Navier-Stokes equations to solve for the saturation increments and pressure information of water, hydrate, and methane at the particle positions.

[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 to obtain the temperature increment of the internal particles.

[0013] 2.5) Update the soil mechanical parameters of the sediment layer according to the saturation increments of each component obtained in step 2.2). The soil mechanical parameters include shear modulus G, bulk modulus K, cohesion coefficient c, and dilatancy angle ψ.

[0014] 2.6) For each time step of iterative solution, after obtaining the velocity increment, temperature increment, pressure information, and saturation increments of each component of the internal particles according to steps 2.2), 2.3), and 2.4), update the corresponding physical information of the internal particles, and calculate the soil mechanical parameters of the sediment layer using step 2.5) according to the hydrate saturation; then enter the iterative solution of the next time step. Among them, at each time step, output the updated particle information of all discrete seabed sediment layer particles, and output the current landslide morphology of the sediment layer according to the particle information.

[0015] The method of the present invention adds the Kim-Bishnoi hydrate decomposition kinetic equation to the migration process of the marine sediment layer, and adds the hydrate decomposition process to the soil movement; at the same time, considers the migration processes of water and gas during the migration of the seabed sediment layer, and calculates the pressure and saturation changes of water and gas in the pores of the sediment layer, which is beneficial to better simulating the influence of hydrate decomposition on the landslide process of the seabed sediment layer.

[0016] The method of the present invention uses a formula to correlate the mechanical properties of the soil in the marine sediment layer with the hydrate concentration. Such treatment makes the sediment layer soil closer to the actual seabed soil, which is beneficial to better simulating the migration process of the hydrate sediment layer during the hydrate decomposition process.

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

[0018] Figure 1 Schematic diagram of each component of the hydrate-bearing sediment layer;

[0019] Figure 2 Phase equilibrium curve of hydrate decomposition;

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

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

[0022] Figure 5 Cumulative plastic strain distribution of the sediment layer under different pressure difference conditions. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0023] The present invention will be further described and illustrated below in conjunction with the specific embodiments. The embodiments are only examples of the present disclosure and do not delimit the scope of limitation. The technical features of each embodiment of the present invention can be combined correspondingly without conflict.

[0024] As Figure 1 shown, it is a schematic diagram of each component of the hydrate-bearing sediment layer studied by the present invention. The hydrate-bearing sediment layer mainly includes components such as natural gas, water, hydrate, and soil. 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, the hydrate can be decomposed from the solid hydrate form into natural gas and water, resulting in a change in the phase state. The phase equilibrium is related to temperature and pressure. Figure 2 Schematically shows the phase equilibrium curve of the hydrate.

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

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

[0027] 1) At the initial moment, according to the scale of the submarine sediment layer to be simulated, use SPH particles to discretize the submarine sediment layer containing hydrate, establish a numerical model of the submarine sediment layer to be simulated and set boundary conditions. The discretized submarine sediment layer particles include boundary particles and internal particles. Among them, establishing a numerical model of the submarine sediment layer to be simulated and setting boundary conditions specifically means: performing mathematical modeling on 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) Perform step-by-step iterative calculation and update on the submarine sediment layer particles according to time steps. At each time step, output the current landslide morphology of the submarine sediment layer according to the current information of the updated particles.

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

[0030] 2.1) For boundary particles, use Shepard interpolation to solve the pressure, stress, velocity, and temperature of the particles. Among them, the Shepard interpolation method is a well-known method in the art. Obtaining the pressure, stress, velocity, and temperature of the boundary particles in this step 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. Use the continuity equation in the Navier-Stokes equations 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 of hydrate:

[0035]

[0036] Continuity equation of water:

[0037]

[0038] Continuity equation of methane gas:

[0039]

[0040] where t is time; is the divergence symbol; are the water production rate, hydrate dissociation rate, and gas production rate, respectively; S h , S w and S g are the saturations of hydrate, water, and methane gas, respectively; u h , u w and u g are the velocities of hydrate, water, and methane gas, 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 hydrate dissociation rate is calculated according to the Kim - Bishnoi hydrate dissociation kinetic equation:

[0042]

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

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

[0045]

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

[0047] The velocities of water and gas are represented by the Darcy seepage velocities q w and q g respectively, and the calculation formulas are:

[0048]

[0049] where μ w and μg are the dynamic viscosities of water and methane gas respectively; P w and P g are the pressures of water and methane gas respectively, is the gradient operator, and are the hydraulic and gas 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 sorting and simplification, the formulas for calculating the saturation increments of water, hydrate and methane gas are respectively:

[0051] Formula for calculating the saturation increment of water:

[0052]

[0053] Formula for the saturation increment of hydrate:

[0054]

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

[0056]

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

[0058] The increment of the particle gas pressure is calculated by the continuity equation of the gas:

[0059]

[0060] Where, is the gas pressure increment.

[0061] 2.3) Calculate the velocity increment of the particle by using the momentum equation through the soil stress and viscous drag force The calculation formula is as follows:

[0062]

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

[0064] The formula for calculating the displacement increment of particles is:

[0065]

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

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

[0068]

[0069] where, is the temperature increment of the particles in the deposition layer, is the divergence symbol, is the temperature gradient, C T is the total heat capacity of the marine deposition layer, K T is the thermal conductivity of the marine deposition layer, Q h is the energy generated by the decomposition of hydrates.

[0070] 2.5) According to the saturation increment of each component obtained in step 2.2), update the soil mechanics parameters of the deposition layer. The soil mechanics parameters include the shear modulus G, the bulk modulus K, the cohesion coefficient c, and the 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 are adjustable coefficients; the superscript "0" represents the material parameter value when the hydrate concentration is zero.

[0076] 2.6) For each time step solved iteratively, after obtaining the velocity increment, temperature increment, pressure information, and the increment of the saturation of each component of the internal particles according to steps 2.2), 2.3), and 2.4), update the corresponding physical information of the internal particles, and calculate the geotechnical parameters of the sediment layer using step 2.5) according to the hydrate saturation; then enter the iterative solution of the next time step; among them, for each time step, output the updated particle information of all discretized seabed sediment layer particles, and output the current landslide form of the sediment layer according to the particle information.

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

[0078] 1) At the initial moment, establish a hydrate-bearing soil model as shown in the appendix Figure 3 , with the upper base length of 1 m, the lower base length of 2 m, and the height of 1 m. Set the initial temperature T of the hydrate soil to 275.45 K. Set the left boundary Ⅰ as the pressure-reducing boundary, with the air pressure P g0 constant at 3.15 MPa, and the remaining boundaries Ⅱ, Ⅲ, and Ⅳ are all non-penetrating closed boundaries. The temperatures of all four sides are constant, with the temperature of 275.45 K. In terms of motion constraints, the left boundary Ⅰ and the lower boundary Ⅳ are no-slip boundary conditions, and the motions of the Ⅱ and Ⅲ boundaries are unconstrained. When performing numerical simulation calculations, only consider the influence of air pressure on hydrate decomposition, and do not consider the influence of air pressure on soil body movement;

[0079] 2) As shown in Figure 1 , the hydrate-bearing sediment layer is composed of water, hydrate, methane, and soil. Each sediment layer particle after discretization carries the physical information of water, hydrate, methane, and soil. Calculate the generation rates of hydrate, water, and methane gas according to the hydrate decomposition kinetic equation based on the information such as the temperature and air pressure of the sediment layer at the current moment. The calculation formula is as follows:

[0080]

[0081]

[0082]

[0083] In the formula, M h , M g and M w are the molar masses of hydrate, methane, and water respectively, and are M h = 0.124 kg / mol, M g = 0.016 kg / mol, and M w= 0.018 kg / mol. Solve for the change in saturation S of water, hydrate, and methane at the position of this particle using the continuity equation in the aforementioned Navier-Stokes equations. w , S h and S g . Obtain the current air pressure through the aforementioned equation of state;

[0084] 3) If a certain particle carries soil information, calculate parameters such as the stress and strain of the soil according to the current air pressure, water pressure, and hydrate saturation concentration using the aforementioned momentum equation;

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

[0086] 5) The mechanical properties of the marine sediment layer will change according to the change in hydrate saturation. Obtain the current hydrate concentration through step 2), and calculate the mechanical parameters of the current soil mass according to the aforementioned formula;

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

[0088] 7) For all fluid particles, after obtaining the velocity increment, temperature increment, pressure increment, and saturation increment of each component of the particles according to steps 2), 3), and 4), update the corresponding physical information of the particles, and calculate the soil mechanical parameters of the sediment layer according to the hydrate saturation through step 5), and enter the solution of the next time step.

[0089] Figure 4 (a) in gives the effective stress distribution diagram of the sediment layer at t = 8 s. σ is the effective stress of the sediment 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 left side of the sediment layer is subjected to the greatest force and is most likely to deform. Figure 4 (b) in is the displacement nephogram of the sediment layer. Δx is the displacement of the sediment layer. There is a maximum displacement to the left near the upper left corner of the sediment layer, and a maximum displacement to the right on the right side of the sediment layer, which is consistent with the strain nephogram of the sediment layer. Figure 4 (c) and (d) in are the air pressure and water pressure nephograms of the sediment layer. P g is the air pressure, and P w is the water pressure. It can be found that the water pressure and air pressure at the boundary are the smallest, and the low pressure is transmitted from the left side to the right side. As can be seen from the figure, this method is a fully coupled method, and the stress, displacement, air pressure, water pressure, etc. in the seabed sediment layer at any moment can be extracted.

[0090] Figure 5 shows the cumulative plastic strain distribution of the sediment layer under different pressure difference conditions. εdev is the cumulative plastic strain of the sediment layer, representing the magnitude of the deformation of the sediment layer. Among them, (a), (b), and (c) are the strain diagrams after calculating to the stable state with the initial differential pressures of 50 kPa, 150 kPa, and 200 kPa respectively. From Figure 5 it can be seen that when the differential pressure is small, the force generated by the pressure gradient is less than the cohesive force of the sediment layer and is not sufficient to cause the sediment layer to deform. As the differential pressure gradually increases, the deformation near the upper left corner of the sediment layer becomes larger and larger, and large deformation occurs when the differential pressure is 200 kPa.

[0091] The above-described embodiments only express one feasible implementation manner of the present invention, and its description is relatively specific and detailed, but it should not be construed as a limitation on the scope of the present invention patent. It should be noted that for those of ordinary skill in the art, without departing from the concept of the present invention, several deformations and improvements can still be made, and these all belong to the protection scope of the present invention. Therefore, the protection scope of the present invention shall be subject to the appended claims.

Claims

1. A particle method for simulating the landslide process of a sediment layer, characterized in that, It includes the following steps: 1) At the initial moment, according to the scale of the seabed sediment layer to be simulated, use SPH particles to discretize the seabed sediment layer containing hydrates, establish a numerical model of the seabed sediment layer to be simulated and set boundary conditions; among them, the discretized seabed sediment layer particles include boundary particles and internal particles; 2) Perform step-by-step iterative calculations and updates on the seabed sediment layer particles at each time step. At each time step, output the current landslide form of the seabed sediment layer according to the current information of the updated particles; The specific steps of step 2) include: 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, hydrates, methane, and soil. The continuity equation in the Navier-Stokes equations is used to solve the saturation increment and pressure information of water, hydrates, and methane at the particle positions; 2.3) Use the momentum equation to calculate the stress and viscous drag force parameters 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 to obtain the temperature increment of the internal particles; 2.5) Update the soil mechanics parameters of the sediment layer according to the saturation increments of each component obtained in step 2.2). The soil mechanics parameters include shear modulus G, bulk modulus K, cohesion coefficient c, and dilatancy angle ψ; 2.6) For each time step of iterative solution, after obtaining the velocity increment, temperature increment, pressure information, and saturation increments of each component of the internal particles according to steps 2.2), 2.3), and 2.4), update the corresponding physical information of the internal particles, and calculate the soil mechanics parameters of the sediment layer using step 2.5) according to the hydrate saturation; then enter the iterative solution of the next time step; among them, at each time step, output the updated particle information of all discretized seabed sediment layer particles, and output the current landslide form of the sediment layer according to the particle information.

2. The particle method for simulating the landslide process of a sedimentary layer according to claim 1, characterized in that, In step 1), the establishment of the numerical model of the seabed sediment layer to be simulated and the setting of boundary conditions include: performing mathematical modeling on 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 the landslide process of a sedimentary layer according to claim 1, wherein, The specific steps of step 2.2) include: The continuity equation in the Navier-Stokes equations includes: The continuity equation of hydrates: The continuity equation of water: The continuity equation of methane gas: where t is time; is the divergence symbol; are the water production rate, the hydrate dissociation rate, and the gas production rate, respectively; S h , S w and S g are the saturations of hydrate, water, and methane gas, respectively; u h , u w and u g are the velocities of hydrate, water, and methane gas, 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: where k d is the kinetic reaction rate, with the unit of mol / (m 2 ·Pa·s); M h is the molar mass of the hydrate; A s is the specific surface area of hydrate decomposition in the unit volume of the sediment layer, with the unit of m -2 ; P g is the methane gas pressure, with the unit of Pa; P eq is the equilibrium pressure. The water production and gas production rates are calculated using mass conservation according to the corresponding hydrate decomposition chemical equations: Where, M g and M w are the molar masses of methane and water, respectively; The velocities of water and gas are represented by the Darcy seepage velocities q w and q g respectively, and the calculation formulas are as follows: where μ w and μ g are the dynamic viscosities of water and methane gas, respectively; P w and P g are the pressures of water and methane gas, respectively. is the gradient operator, and are the gradients of water pressure and gas pressure; k is the absolute permeability of the sediment; k w and k g are the relative permeabilities of water and methane gas, respectively. After sorting and simplifying, the formulas for calculating the saturation increments of water, hydrates, and methane gas are respectively: The calculation formula for the saturation increment of water: The formula for the saturation increment of hydrates: Due to the existence of a relationship Therefore, the increment of methane gas saturation is calculated by the following formula: wherein, is the increment of water saturation, is the increment of hydrate saturation, is the increment of methane gas saturation; The increment of particle air pressure is calculated by the continuity equation of gas: wherein, is the air pressure increment.

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

5. A particle method for simulating the landslide process of a sedimentary layer, characterized in that, The energy equation formula adopted in step 2.4) is as follows and is used to calculate the temperature increment of the particles in the deposition layer Among them, is the temperature increment of the sediment layer particles, is the divergence symbol, is the temperature gradient, C T is the total heat capacity of the ocean sediment layer, K T is the thermal conductivity of the ocean sediment layer, Q h is the energy generated by hydrate decomposition.

6. A particle method for simulating the landslide process of a sedimentary layer 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 = α G S h + G0 S h is the hydrate saturation; α K , α G , α c and β c are adjustable coefficients; the superscript "0" represents the value of the material parameter when the hydrate saturation is zero.

Images