Saturated soil and structure large deformation analysis method based on semi-implicit substance point method

Abstract

The invention relates to a saturated soil and structure large deformation analysis method based on a semi-implicit substance point method. The method comprises the following steps: decoupling pore pressure and motion variables by adopting an incremental method; introducing a penalty function of an enhanced frictional contact algorithm; manual compressibility parameters are introduced into a pressure Poisson equation; the volume locking effect is reduced by adopting a B-bar technology; an affine particle cellular method scheme is adopted; modifying a calculation formula for an axial symmetry problem, and considering volume integration and gradient evaluation in a cylindrical coordinate system; and performing verification. The method has the advantages that the incremental method, the artificial fluid compressibility and the B-bar technology are adopted, numerical instability caused by the coupling problem in a traditional method is avoided, the convergence of the simulation process is ensured, and non-physical stress oscillation is avoided.

Description

Technical Field

[0001] The present invention relates to the technical field of geotechnical engineering numerical simulation, and more specifically, it relates to a large deformation analysis method for saturated soil and structure based on the semi-implicit material point method. Background Technique

[0002] In geotechnical engineering, the large deformation problems of saturated soil masses (such as pile foundation installation, slope instability, etc.) are highly non-linear and complex. The traditional finite element method (FEM) is prone to calculation failure due to mesh distortion when dealing with large deformations. As a meshless method, the material point method (MPM) can effectively handle large deformation problems. However, when dealing with the water-hydrodynamic coupling problem, explicit MPM has problems such as limited time step size and pore pressure oscillation. Although the existing semi-implicit MPM methods alleviate these problems to a certain extent, there are still problems such as insufficient stability and low calculation efficiency when dealing with soil-structure interaction, axisymmetric problems, and low-permeability materials.

[0003] Although the prior art discloses the use of polynomial pressure projection stabilization techniques in FEM analysis to alleviate the checkerboard stress pattern caused by volume locking, some non-physical oscillations still exist. Summary of the Invention

[0004] The purpose of the present invention is to address the deficiencies of the prior art and propose a large deformation analysis method for saturated soil and structure based on the semi-implicit material point method.

[0005] In the first aspect, a large deformation analysis method for saturated soil and structure based on the semi-implicit material point method is provided, including:

[0006] Step 1: Decouple the pore pressure and motion variables using the incremental step method;

[0007] Step 2: Introduce the penalty function of the enhanced friction contact algorithm, adjust the velocity correction of the contact nodes, and simplify the calculation of the contact force using the lumped mass matrix;

[0008] Step 3: Introduce an artificial compressibility parameter into the pressure Poisson equation to alleviate the pore pressure oscillation in low-permeability materials, and optimize the selection of the time step by adjusting the artificial compressibility parameter;

[0009] Step 4: Use the B-bar technique to reduce the volume locking effect, and calculate the particle strain increment and the internal force of the solid phase nodes by modifying the B matrix;

[0010] Step 5: Adopt the affine particle cell method (APIC) scheme to add a local affine description for each particle;

[0011] Step 6: Modify the calculation formula for axisymmetric problems, consider the volume integral and gradient evaluation in the cylindrical coordinate system, and expand the application scope of MPM in axisymmetric problems;

[0012] Step 7: Verify the effectiveness of the semi-implicit MPM in dealing with soils with different permeabilities, the accuracy of the coupled axisymmetric formula, and the stability and accuracy of the semi-implicit MPM in dealing with soil-structure interaction problems by analyzing one-dimensional consolidation problems, Cryer sphere problems, and strip foundation problems.

[0013] Preferably, in Step 1, it is divided into a prediction step and a correction step; in the prediction step, the pore pressure gradient term is ignored to solve the intermediate velocity, and in the correction step, the intermediate velocity is explicitly corrected by implicitly solving the pressure Poisson equation.

[0014] Preferably, in Step 2, the penalty function is used to simulate the contact between an impermeable rigid structure and a deformable saturated medium.

[0015] Preferably, in Step 6, the modification of the calculation formula includes: modifying the calculation formulas for particle volume, strain, pressure Poisson equation, and nodal internal force.

[0016] Preferably, in Step 6, the application scope of MPM in axisymmetric problems is expanded by introducing additional shape functions and gradient formulas.

[0017] In a second aspect, a large deformation analysis system for saturated soil and structure based on the semi-implicit material point method is provided for performing any of the methods in the first aspect, including:

[0018] A decoupling module for decoupling pore pressure and motion variables using the incremental step method;

[0019] An adjustment module for introducing a penalty function of an enhanced friction contact algorithm to adjust the velocity correction of contact nodes and simplifying the calculation of contact forces using a lumped mass matrix;

[0020] A mitigation module for introducing an artificial compressibility parameter into the pressure Poisson equation to mitigate pore pressure oscillations in low-permeability materials and optimizing the selection of time steps by adjusting the artificial compressibility parameter;

[0021] A first modification module for reducing the volume locking effect using the B-bar technique and calculating the particle strain increment and the nodal internal force of the solid phase by modifying the B matrix;

[0022] An addition module for adding a local affine description to each particle using the affine particle cell method scheme;

[0023] A second modification module is used to modify the calculation formula for axisymmetric problems, consider volume integration and gradient evaluation in the cylindrical coordinate system, and expand the application scope of MPM in axisymmetric problems;

[0024] An analysis and verification module is used to verify the effectiveness of the semi-implicit MPM in dealing with soils with different permeabilities, the accuracy of the coupled axisymmetric formula, and the stability and accuracy of the semi-implicit MPM in dealing with soil-structure interaction problems by analyzing one-dimensional consolidation problems, Cryer sphere problems, and strip foundation problems.

[0025] In a third aspect, a computer storage medium is provided, and a computer program is stored in the computer storage medium; when the computer program runs on a computer, the computer is enabled to execute any method described in the first aspect.

[0026] In a fourth aspect, an electronic device is provided, including:

[0027] A memory for storing a computer program;

[0028] A processor for executing the computer program to implement any method described in the first aspect.

[0029] The beneficial effects of the present invention are as follows:

[0030] 1. The present invention adopts an incremental step method, artificial fluid compressibility, and B-bar technology. The incremental step method effectively decouples pore pressure and motion variables, avoiding numerical instability caused by coupling problems in traditional methods; the introduction of the artificial fluid compressibility parameter alleviates pore pressure oscillations in low-permeability materials, ensuring the convergence of the simulation process; the B-bar technology further enhances the stability of the simulation by reducing the volume locking effect, avoiding non-physical stress oscillations.

[0031] 2. The present invention adopts the APIC scheme and a lumped mass matrix. The APIC scheme enhances the momentum conservation property by adding a local affine description for each particle, reducing energy dissipation and numerical noise, and at the same time improving the information transfer accuracy between particles and the background grid; the use of the lumped mass matrix simplifies the calculation of contact forces, reduces the computational complexity, and significantly improves the computational efficiency of large-scale numerical simulations.

[0032] 3. The present invention reduces the checkerboard stress pattern and improves the accuracy of the simulation results through an enhanced friction contact algorithm and B-bar technology. The enhanced friction contact algorithm reduces false gaps and numerical fluctuations in the soil-structure contact area by introducing a penalty function and adjusting the velocity correction of contact nodes, ensuring the accurate calculation of contact forces.

[0033] 4. The present invention is applicable to a variety of geotechnical engineering problems, such as pile foundation installation, soil-structure interaction, etc. Description of the Drawings

[0034] Figure 1 Schematic diagram of the single - point MPM description for two - phase porous media modeling provided by the present invention;

[0035] Figure 2 Implementation process of the coupled semi - implicit MPM algorithm provided by the present invention;

[0036] Figure 3 Geometry and boundary conditions of the one - dimensional consolidation test provided by the present invention;

[0037] Figure 4 Pore - pressure dissipation provided by the present invention: (a) Comparison between MPM simulation and Terzaghi solution; (b) Evolution at point A under different hydraulic conductivities;

[0038] Figure 5 Cryer sphere problem provided by the present invention: (a) Geometry and boundary conditions; (b) Pore - pressure distribution at ν = 0.3 at different times;

[0039] Figure 6 Evolution of pore - pressure at the center of the sphere under different Poisson's ratios provided by the present invention;

[0040] Figure 7 Schematic diagram of the geometry and boundary conditions of the strip - footing model provided by the present invention;

[0041] Figure 8 Evolution of (a) contact force and (b) pore - pressure at point A during the penetration process of the strip - footing provided by the present invention;

[0042] Figure 9 Pore - pressure at the end of penetration provided by the present invention: (a) MPM proposed by the present invention and (b) Comparison chart of FEM results of Monforte et al. (2023);

[0043] Figure 10 Pore - pressure contour lines of the MPM results at the end of penetration provided by the present invention: (a - b) Displacement and deviatoric strain, (c - d) Effective vertical stress and mean stress;

[0044] Figure 11 Load - displacement curve of the large - penetration analysis provided by the present invention;

[0045] Figure 12 Simulation results of the strip - footing at the end of large - penetration provided by the present invention. Detailed Implementation Manner

[0046] The present invention will be further described below in conjunction with embodiments. The description of the following embodiments is only for helping to understand the present invention. It should be noted that for those of ordinary skill in the art, without departing from the principle of the present invention, several modifications can be made to the present invention, and these improvements and modifications also fall within the protection scope of the claims of the present invention.

[0047] Embodiment 1:

[0048] As an embodiment, as Figure 1 and Figure 2 shown, this large deformation analysis method of saturated soil and structure based on the semi - implicit material point method uses a single - layer material point to represent the solid phase and the liquid phase.

[0049] This large deformation numerical simulation method of saturated soil mass based on the semi - implicit material point method includes the following steps:

[0050] Step 1, incremental fractional step method: The pore pressure and motion variables are decoupled by the incremental fractional step method, which is divided into a prediction step and a correction step.

[0051] In the prediction step, the intermediate velocity is solved by ignoring the pore pressure gradient term; in the correction step, the intermediate velocity is explicitly corrected by implicitly solving the pressure Poisson equation, specifically as follows:

[0052] (1) According to the mixture theory, saturated soil is homogenized into a continuous medium composed of a solid phase and a liquid phase. The partial density of each component respectively representing the solid phase and the liquid phase) is determined by its volume fraction and the intrinsic density ρ α as:

[0053]

[0054] where n is the porosity, and ρ m is the total density of the mixture.

[0055] In saturated porous media, Terzaghi effective stress is used, expressed as:

[0056] σ = σ' - pI (2)

[0057] where σ is the total stress tensor (tensile is positive), σ' is the effective stress tensor, p is the pore pressure (compressive is positive), and I is the unit tensor. The partial stresses of the solid phase and the liquid phase are defined as:

[0058] σ s = σ' - (1 - n)pI, σ l = -npI (3)

[0059] Assuming that both solid particles and liquid are incompressible, the mass balance equations of the solid phase and the liquid phase are written as:

[0060]

[0061]

[0062] where v α is the velocity vector of the α-phase, representing the time derivative with respect to the α-phase. In the single-point two-phase MPM, the time derivative of the liquid phase with respect to the motion of the solid phase is described as:

[0063]

[0064] Combining equations (4a), (4b), and (5), and neglecting the porosity gradient (i.e., ), the mass balance equation of the mixture can be expressed as:

[0065]

[0066] The motion of the saturated medium is controlled by the momentum conservation of the involved solid and liquid phases:

[0067]

[0068] where b is the body force vector, and f d is the drag force due to the solid-liquid interaction. The present invention considers laminar flow, and the drag force satisfies Darcy's law:

[0069]

[0070] where g is the acceleration due to gravity (g = 9.81 m / s 2 ), and k is the hydraulic conductivity.

[0071] The sum of equations (7a) and (7b) gives the momentum balance equation of the mixture:

[0072]

[0073] (2) By discretizing the velocity variable using the forward Euler scheme, the time-discretized form of the momentum equation can be written as:

[0074]

[0075] where the superscript k represents the discrete approximation at time t k , and Δt is the time increment.

[0076] (3) Applying the incremental step algorithm to decouple the pore pressure and motion variables in the solid-liquid coupled system facilitates the equal-order interpolation of the displacement field and the pressure field. According to the incremental step procedure, the momentum equation is solved by introducing an intermediate velocity It is divided into a prediction step and a correction step. The drag term in Equation (10b) is formulated based on the intermediate velocity, which makes the time step independent of the permeability. The prediction step is expressed as:

[0077]

[0078] The correction step is expressed as:

[0079]

[0080] Substituting Equation (12) into Equation (6), the pressure Poisson equation can be derived from the mass balance equation of the mixture as follows:

[0081]

[0082] (4) In the present invention, the generalized interpolation material point (GIMP) method and B-spline-MPM are adopted to spatially discretize the momentum balance equation (Equation 11) in the prediction step as:

[0083]

[0084] where, v α,i , and respectively represent the nodal velocity, internal force and external force of phase α at node i. The nodal matrices related to mass and drag are respectively denoted as M α and Q:

[0085]

[0086] where, S iq is the shape function evaluated at material point p; the subscripts i and p represent the node and the material point respectively; N p is the total number of material points related to node i, and V p is the volume of each material point. It should be noted that although the drag term is implicitly calculated based on the intermediate velocity, the above momentum equation can be converted into a nodal explicit expression by replacing M α and Q. The calculation procedure is simplified and the calculation efficiency is improved by the nodal explicit expression, while maintaining the prediction accuracy.

[0087] The pressure Poisson equation in Equation (13) is spatially discretized as:

[0088]

[0089] The velocity is updated through an explicit correction procedure as shown in Equation (12), and its spatial discretization is:

[0090]

[0091] Among them, φ l,p = n p and φ s,p = 1 - n p respectively represent the volume fractions of the solid and liquid phases. A lumped mass matrix is adopted in the correction step.

[0092] Step 2: Enhanced frictional contact algorithm: Introduce a penalty function. By adjusting the velocity correction of the contact nodes, reduce the spurious gaps and numerical fluctuations in the soil-structure contact area. Adopt a lumped mass matrix to simplify the calculation of contact forces and improve the calculation efficiency.

[0093] The contact constraint on the interface nodes is determined by the impenetrability detection condition as:

[0094] (v rigid,i - v α,i )·n r,i > 0 (20)

[0095] Among them, n r,i is the unit outer normal vector of the rigid body at node i, determined by its mass gradient, and v rigid,i represents the velocity vector of the rigid body.

[0096] When it is detected that the interface nodes are in contact, the nodal velocities and accelerations of the deformable body in the prediction and correction steps should be adjusted to avoid mutual penetration. The correction of the nodal intermediate velocity and the corresponding contact force applied to phase α are calculated as:

[0097]

[0098] Among them, ω is the unit tangential vector of the contact surface, μ is the interface friction coefficient, and μ c is the Coulomb friction coefficient. The first term on the right side of equation (24) represents the normal velocity component constrained by the impenetrability condition and the velocity continuity assumption, while the second term represents the tangential velocity correction controlled by the Coulomb friction criterion. Therefore, for the liquid phase, the tangential correction related to the friction term should be discarded (i.e., μ = 0). For simplicity, a lumped mass matrix is adopted to allow explicit and nodal correction of nodal velocities.

[0099] The nodal velocity correction and contact force in the correction step are calculated as:

[0100]

[0101] Combining the above, the nodal contact force applied to the saturated medium is:

[0102]

[0103] Step 3, Artificial Fluid Compressibility: Introduce an artificial compressibility parameter into the pressure Poisson equation to effectively alleviate pore pressure oscillations in low-permeability materials. By adjusting the artificial compressibility parameter, optimize the selection of the time step to improve numerical stability.

[0104] Step 3 is specifically as follows:

[0105]

[0106] where ε is a sufficiently small parameter used to control the degree of artificial compressibility, expressed as:

[0107]

[0108] where the degree of artificial compressibility can be adjusted by the parameter β0, and the optimal solution is Δt≈Δt cr . The critical time step Δt cr is expressed as:

[0109]

[0110] where L min is the minimum element size, and K and G are the bulk modulus and shear modulus of the solid, respectively.

[0111] Step 4, B-bar Technique: Adopt the B-bar technique to reduce the volume locking effect. By modifying the B matrix to calculate the particle strain increment and the nodal internal force of the solid phase, reduce non-physical stress oscillations.

[0112] Among them, reduced integration of the volumetric component and full integration of the deviatoric component are adopted. The specific form of the B-bar matrix in the plane strain configuration is:

[0113]

[0114] where and the gradients of the shape functions evaluated at the material point p, in the x and y directions respectively, and are the gradients of the shape functions evaluated at the element center.

[0115] Step 5, Affine Particle Cell Method (APIC) Scheme: Enhance the momentum conservation property, reduce energy dissipation and numerical noise by adding a local affine description for each particle.

[0116] In addition, adopt the APIC scheme to improve the information transfer accuracy between particles and the background grid and enhance the stability of the calculation, specifically as follows:

[0117] To balance energy conservation and numerical stability, the present invention adopts the APIC scheme. In APIC, the velocity transfer operator enhances the momentum conservation property by adding a local affine description for each particle. The nodal mass m i and the nodal momentum are obtained by mapping the particle mass and the local affine velocity to the grid nodes, which are respectively:

[0118]

[0119] The local affine velocity field of the material point consists of a translational component and a linear component, and is specifically expressed as:

[0120]

[0121] Among them, C α,p represents the linear transformation of α at the material point p, v α,p represents the translational velocity, and x is the spatial position vector.

[0122] The affine matrix C α,p is calculated as:

[0123] C α,p = B α,p (D p ) -1 (33)

[0124]

[0125]

[0126] Among them, B α,p represents the affine state containing angular momentum, and D p is an inertia-like matrix derived by maintaining affine motion.

[0127] Step 6, Axisymmetric formula: For axisymmetric problems, modify the calculation formulas of particle volume, strain, pressure Poisson equation, and nodal internal force, and consider volume integration and gradient evaluation in the cylindrical coordinate system. By introducing additional shape functions and gradient formulas, expand the application range of MPM in axisymmetric problems.

[0128] The axisymmetric formula is characterized by the r-z plane (where r and z represent the radial and axial directions respectively). In the cylindrical coordinate system, the circumferential direction θ is taken as radian 1, and the particle volume in the wedge region of radian 1 is calculated as:

[0129]

[0130] Among them, χ pIt is a particle characteristic function in GIMP. The Dirac delta function is usually selected in the standard linear MPM or B-spline - MPM. A p and represent the area and radial position of particle p respectively.

[0131] In axisymmetric MPM, an additional shape function T ip is needed to represent the circumferential variable, denoted as T ip = S ip / r p , and the shape function S ip and the gradient are calculated according to the plane formula.

[0132] The components of the strain increment Δε in axisymmetric coordinates are calculated by the following formula:

[0133]

[0134] where, (v s,i ) r and (v s,i ) z are the radial and axial components of the nodal velocity respectively, and are the gradients of the shape function in the r and z directions respectively. Compared with the plane strain condition, an additional circumferential strain increment (Δε s,p ) θθ is considered.

[0135] For the pressure Poisson equation, the right - hand side term characterizing the volume change in formula (16) should be adjusted by considering the circumferential deformation, and the expression is as follows:

[0136]

[0137] The internal forces of the solid and liquid phases, adding the circumferential stress (σ α,p ) θθ are expressed as:

[0138]

[0139] The corresponding term in the correction step should also be modified by evaluating the contribution of the circumferential pore pressure increment as:

[0140]

[0141] Using the B - bar technique, the expression under axisymmetric conditions is adjusted to:

[0142]

[0143] Among them, variables with subscript p represent values evaluated at particle positions, while variables with subscript c represent values evaluated at cell centers.

[0144] Step 7, Application and Result Output: By analyzing one-dimensional consolidation problems, Cryer sphere problems, and strip foundation problems, verify the effectiveness of semi-implicit MPM in dealing with soils with different permeabilities, the accuracy of the coupled axisymmetric formula, and the stability and accuracy of semi-implicit MPM in dealing with soil-structure interaction problems.

[0145] Example 2:

[0146] Based on Example 1, Example 2 of this application verifies the effectiveness of semi-implicit MPM in dealing with soils with different permeabilities by simulating the classical one-dimensional consolidation problem through 2D plane strain MPM analysis, as follows:

[0147] Simulate the classical one-dimensional consolidation problem through 2D plane strain MPM analysis. As Figure 3 shown, a constant vertical stress of 10 kPa is applied at the top of a saturated soil column that is 1.0 m high and 0.02 m wide. Only the top surface is allowed to drain, rolling boundaries are applied on both sides, and the bottom is fully constrained. The computational domain is discretized by linear quadrilateral elements of size 0.02 m × 0.02 m, and each element contains 2 × 2 material points. The solid phase is considered linearly elastic, and the model parameters are set as follows: Young's modulus E = 1×10 7 Pa, Poisson's ratio ν = 0.3, soil particle density ρ s = 2650 kg / m 3 , pore fluid density ρ l = 1000 kg / m 3 , initial porosity n0 = 0.3, hydraulic conductivity k = 1×10 -3 m / s. Assume the initial pore pressure and effective stress are zero, and gravity is ignored. Take the time step as Δt = 1.0×10 -4 s.

[0148] Figure 4 (a) shows the spatial distribution of pore pressure along the column depth at different dimensionless times T (i.e., T = c v t / H 2 , where c v is the consolidation coefficient and H is the column height). It can be seen that the numerical results of the coupled MPM are in good agreement with Terzaghi's analytical solution. By fixing the time step Δt = 1.0×10 -4 s, a series of consolidation tests with different hydraulic conductivities are simulated. Figure 4(b) Plots the simulated pore pressure of the selected material point (i.e., point A at h = 0.5 m) and the corresponding analytical solution. The numerical results are close to the theoretical solution. Although there are some slight deviations in the case of lower permeability (k = 1.0×10 -6 m / s), they can also be solved by reducing the time step. When considering the consolidation of low-permeability soils such as clay, a smaller time step can be adopted.

[0149] It should be noted that the same or similar parts in this embodiment and Embodiment 1 can be referred to each other and will not be elaborated in this application.

[0150] Embodiment 3:

[0151] Based on Embodiment 1, Embodiment 3 of this application verifies the accuracy of the coupled axisymmetric formula by simulating the Cryer sphere problem using axisymmetric MPM. Specifically as follows:

[0152] The semi-implicit MPM and the simulation results are compared with the theoretical solution to verify the coupled axisymmetric formula. Figure 5 Shows how to establish an axisymmetric two-dimensional (2D) model in the cylindrical coordinate system (r, z, θ), where θ = 1 radian and the radius R = 0.5 m. Free drainage is allowed on the outer surface, and an impermeable boundary is assigned on the vertical axis of the sphere. An isotropic pressure P0 = 10 kPa is applied on the sphere surface. The background grid is constructed from 0.05 m × 0.05 m elements, and each element contains 4 material points. To improve the accuracy of the external pressure application, the material points near the outer surface are uniformly repositioned, allowing the calculation of particle traction based on a constant curvature. Let the initial pore pressure and stress be zero. The material properties of the isotropic elastic sphere are: E = 1×10 7 Pa, ρ s = 2143 kg / m 3 , ρ l = 1000 kg / m 3 , n0 = 0.3, k = 1.0×10 -3 m / s. Considering various Poisson's ratios, the time increment is 2.0×10 -4 s.

[0153] Figure 6 Shows the variation of the pore pressure at the center of the sphere with the dimensionless time T (i.e., T = c v t / R 2 ) under different Poisson's ratios. The MPM results are in good agreement with Cryer's analytical results. Figure 5 (b) Shows the pore pressure distribution inside the sphere for ν = 0.3 at T = 0.05, 0.15, and 0.20 (i.e., t = 0.0091 s, 0.0182 s, 0.0273 s).

[0154] It should be noted that the parts that are the same as or similar to those in Embodiment 1 in this embodiment can be referred to each other and will not be elaborated in this application.

[0155] Embodiment 4:

[0156] Based on Embodiment 1, in Embodiment 4 of this application, through small deformation and large deformation analyses, the stability and accuracy of semi - implicit MPM in dealing with soil - structure interaction problems are verified.

[0157] Using the coupled MPM proposed by the present invention for plane - strain strip - footing analysis, Figure 7 Figure shows the schematic diagram of the geometry and boundary conditions of the strip - footing model. The foundation width B = 1m, the width and height of the computational domain are both set to W = H = 4B for small - deformation analysis and extended to W = 12B and H = 10B for large - deformation analysis. The side edges are sliding supports and the bottom edge is a fixed support. Drainage is allowed on the top surface of the foundation, and the soil on both sides of the foundation is subjected to a vertical uniform load q = 10 kPa. The rigid foundation is vertically pushed into the saturated MCC soil at a specified constant speed. Gravity is ignored, the initial pore pressure is set to zero, and the initial effective vertical and horizontal stresses are set to 10 kPa. The slope of the unloading - reloading line κ = 0.01, the slope of the normal compression line λ = 0.1, the slope of the critical - state line M = 1, the Poisson's ratio ν = 0.3, the initial pre - consolidation pressure p co = 15 kPa, the soil - particle density ρ s = 2650 kg / m 3 , the fluid density ρ l = 1000 kg / m 3 , the initial porosity n0 = 0.45, the permeability k = 1.0×10 -10 m / s. The B - bar GIMP scheme and artificial compressibility (β0 = 3) are adopted to achieve better modeling stability.

[0158] (1) Small - deformation analysis

[0159] The small - deformation analysis of foundation penetration involves vertically displacing the rigid strip - footing to a maximum penetration depth of 0.05 m at a speed of 0.05 m / s, with a time step Δt = 2×10 -4 s and an interface friction coefficient μ = 1.

[0160] Figure 8 shows the evolution of the contact force and pore pressure at the bottom of the foundation (point A in Figure 7 ) during penetration, which is in good agreement with the results obtained by Monforte et al. (2023) using nodal - smoothed finite elements with variable element types and sizes. The pore - pressure contour lines of the MPM results at the end of penetration are as shown in Figure 9As shown, when presented together with the FEM results of Monforte et al. (2023), good agreement is shown between the two. The displacement field, deviatoric strain, and effective stress field at the end of the simulation are as Figure 10 shown. The MPM results are also very similar to the FEM results given by Monforte et al. (2023). It is worth noting that although a polynomial pressure projection stabilization technique was adopted in the FEM analysis to alleviate the checkerboard stress pattern caused by volumetric locking, some non-physical oscillations still exist, while these oscillations are avoided in the semi-implicit MPM analysis proposed in the present invention.

[0161] (2) Large deformation analysis

[0162] The model setup is the same as that for the small deformation analysis, except for the extended computational domain. Additionally, to improve the computational efficiency, the penetration speed was doubled to 0.10 m / s, and the final penetration depth was extended to 2 m. A parametric analysis was conducted to explore the influence of the interface friction coefficient μ on the contact force. Figure 10 Shows the load-displacement curves obtained from the coupled B-bar GIMP simulation with different interface friction coefficients. While the small deformation analysis stabilizes to a finite limit, the large deformation resistance steadily increases with the increase in penetration depth, almost at a constant rate, and is insensitive to μ.

[0163] To alleviate the contact oscillations and spurious gaps that are common during large deformations due to premature contact, the analysis adopted an enhanced contact algorithm, Figure 11 and the contact force curves in

[0164] Figure 12 verify the effectiveness of the enhanced contact algorithm.

[0165] It should be noted that the same or similar parts in this embodiment and Embodiment 1 can be referred to each other and will not be elaborated in this application.

[0166] Embodiment 5:

[0167] Based on Embodiment 1, Embodiment 5 of the present application provides a large deformation analysis system for saturated soil and structure based on the semi-implicit material point method, including:

[0168] A decoupling module for decoupling the pore pressure and kinematic variables using the incremental step method;

[0169] Adjustment module, which is used to introduce the penalty function of the enhanced frictional contact algorithm, adjust the velocity correction of the contact nodes, and simplify the calculation of the contact force by using the lumped mass matrix;

[0170] Mitigation module, which is used to introduce an artificial compressibility parameter into the pressure Poisson equation, mitigate the pore pressure oscillation in low-permeability materials, and optimize the selection of the time step by adjusting the artificial compressibility parameter;

[0171] First modification module, which is used to reduce the volume locking effect by using the B-bar technique, and calculate the particle strain increment and the internal force of the nodes in the solid phase by modifying the B matrix;

[0172] Addition module, which is used to adopt the affine particle cell method scheme to add a local affine description for each particle;

[0173] Second modification module, which is used to modify the calculation formula for axisymmetric problems, consider the volume integral and gradient evaluation in the cylindrical coordinate system, and expand the application range of MPM in axisymmetric problems;

[0174] Analysis and verification module, which is used to verify the effectiveness of the semi-implicit MPM in dealing with soils with different permeabilities, the accuracy of the coupled axisymmetric formula, and the stability and accuracy of the semi-implicit MPM in dealing with soil-structure interaction problems by analyzing the one-dimensional consolidation problem, the Cryer sphere problem, and the strip foundation problem.

[0175] It should be noted that the system provided in this embodiment is the system corresponding to the method provided in Embodiment 1. Therefore, for the parts that are the same or similar in this embodiment and Embodiment 1, they can be referred to each other and will not be elaborated in this application.

Claims

1. A large deformation analysis method for saturated soil and structure based on semi-implicit material point method, characterized in that: include: Step 1: Decouple pore pressure and motion variables using incremental step method; Step 2: Introduce the penalty function of the enhanced friction contact algorithm, adjust the velocity correction of the contact node, and use a concentrated mass matrix to simplify the calculation of the contact force; Step 3: Introduce an artificial compressibility parameter into the pressure Poisson equation to alleviate the pore pressure oscillation in the low permeability material, and optimize the selection of the time step by adjusting the artificial compressibility parameter; Step 4: Use B-bar technology to reduce the volume locking effect and calculate the particle strain increment and the nodal internal force of the solid phase by modifying the B matrix; Step 5: Use the affine particle cell method to add a local affine description to each particle; Step 6: Modify the calculation formula for axisymmetric problems, consider the volume integral and gradient evaluation in the cylindrical coordinate system, and expand the application scope of MPM in axisymmetric problems; Step 7. By analyzing the one-dimensional consolidation problem, Cryer sphere problem, and strip foundation problem, verify the effectiveness of the semi-implicit MPM in dealing with soils of different permeability, the accuracy of the coupled axisymmetric formula, and the stability and accuracy of the semi-implicit MPM in dealing with soil-structure interaction problems.

2. The method for large deformation analysis of saturated soil and structure based on semi-implicit material point method according to claim 1 is characterized in that: In step 1, there are prediction step and correction step; in the prediction step, the pore pressure gradient term is ignored to solve the intermediate velocity, and in the correction step, the intermediate velocity is explicitly corrected by implicitly solving the pressure Poisson equation.

3. The method for large deformation analysis of saturated soil and structure based on semi-implicit material point method according to claim 2 is characterized in that: In step 2, the penalty function is used to simulate the contact between the impermeable rigid structure and the deformable saturated medium.

4. The method for large deformation analysis of saturated soil and structure based on semi-implicit material point method according to claim 3 is characterized in that: In step 6, the modification of the calculation formula includes: modifying the calculation formulas of particle volume, strain, pressure Poisson's equation and node internal force.

5. The method for large deformation analysis of saturated soil and structure based on semi-implicit material point method according to claim 3 is characterized in that: In step 6, the application scope of MPM is extended to axisymmetric problems by introducing additional shape functions and gradient formulations.

6. A saturated soil and structure large deformation analysis system based on semi-implicit material point method, characterized in that: The method for executing any one of claims 1 to 5 comprises: Decoupling module for decoupling pore pressure and kinematic variables using an incremental step method; A tuning module for introducing penalty functions for the enhanced friction contact algorithm, adjusting velocity corrections for contact nodes, and simplifying the calculation of contact forces using a lumped mass matrix; Mitigation module, which is used to introduce artificial compressibility parameters into the pressure Poisson equation to mitigate pore pressure oscillations in low permeability materials and optimize the choice of time step by adjusting the artificial compressibility parameters; The first modification module is used to reduce the volume locking effect by using the B-bar technique, and calculate the particle strain increment and the nodal internal force of the solid phase by modifying the B matrix; Added a module for adding a local affine description to each particle using the affine particle cellular method; The second modification module is used to modify the calculation formula for axisymmetric problems, consider the volume integral and gradient evaluation in the cylindrical coordinate system, and expand the application scope of MPM in axisymmetric problems; The analysis and verification module is used to verify the effectiveness of the semi-implicit MPM in dealing with soils of different permeability, the accuracy of the coupled axisymmetric formula, and the stability and accuracy of the semi-implicit MPM in dealing with soil-structure interaction problems by analyzing one-dimensional consolidation problems, Cryer sphere problems, and strip foundation problems.

7. A computer storage medium, characterized in that: The computer storage medium stores a computer program; when the computer program is executed on a computer, the computer executes any one of the methods described in claims 1 to 5.

8. An electronic device, characterized in that: include: Memory, used to store computer programs; A processor, configured to execute the computer program to implement the method according to any one of claims 1 to 5.

Images