A method for analyzing large deformation of saturated porous media by considering cross iteration and variable time step mapping

Abstract

A method for large deformation analysis of saturated porous media using MFM (Meshless Dynamic Model) considering cross-iteration and variable time-step mapping is presented, belonging to the field of numerical simulation of large deformation in geotechnical engineering. The steps are as follows: S1, sequentially read in the node information, material information, and load information of the meshless model to generate a background mesh; S2, construct the radial basis and additional basis for each Gaussian point and calculate the shape function; S3, discretize the governing equations of the meshless saturated porous media; S4, perform variable time-step mapping using the relative energy error of the Gaussian points as the criterion; S5, optimize the constitutive internal variable vector after mapping in S4 and output the calculation results. This invention avoids the problems of ill-conditioned stiffness matrix and difficulty in convergence of iterative solutions. By using the relative energy error as the mapping threshold and performing variable time-step mapping, it improves computational efficiency and has significant engineering and academic value for efficient large deformation analysis of saturated porous media.

Description

Technical Field

[0001] This invention belongs to the field of numerical simulation of large deformation in geotechnical engineering, and is a method for analyzing large deformation of MFM saturated porous media that considers cross-iteration and variable time step mapping. Background Technology

[0002] As the importance and complexity of practical engineering projects increase, the demands on numerical simulation are also constantly rising. Currently, numerical simulation techniques for many projects are developing towards stronger material nonlinearity, large geometric deformation, and multi-field coupling, such as structural instability in coastal nuclear power plants, foundation weakening-liquefaction, earth-rock dam failure, and landslides and debris flows. It is evident that in the field of geotechnical engineering, the simulation of large deformation problems in saturated porous media has become a key and challenging area of ​​numerical simulation research. Therefore, capturing the spatiotemporal evolution of field variables in saturated soil under complex loads and achieving efficient large-scale analysis has significant engineering implications and academic value.

[0003] Meshless methods have been an important tool for large deformation analysis since their inception. They calculate shape functions based on the relative distances between meshless nodes and perform numerical integration through a background mesh, thus eliminating dependence on individual elements during the solution process. Scholars have also conducted research on improving the shape functions of meshless methods (MFM), optimizing the background mesh, multi-numerical coupling, multi-field coupling, and improving computational efficiency. However, the following problems still exist: 1) The stiffness matrix of saturated porous media is ill-conditioned in meshless methods, making iterative solutions difficult to converge, and direct solutions further increase the computational burden; 2) In the process of mapping Gaussian point field variables, frequent mapping processes may lead to error accumulation and significantly affect computational efficiency. For example, Chinese invention patent CN202310964065.7 discloses a method for analyzing large deformations in soil liquefaction based on meshless RBF mapping technology. This method uses a meshless approach to simulate the large deformation process of saturated soil, but it still employs a traditional direct solution framework and does not consider variable time step mapping, thus limiting computational efficiency and stability. Similarly, Chinese invention patent CN202211594488.6 discloses a numerical simulation method for dynamic large deformations of dams in earthquake-induced liquefaction sites. This simulation method uses the idea of ​​mesh re-subdivision and mapping combined with polygonal scale boundary finite element calculation to simulate the large deformation process of saturated soil, but it also does not employ a cross-iterative solution framework or post-mapping parameter correction. Due to these problems, MFM has been difficult to widely apply in practical engineering problems involving large deformations of saturated porous media.

[0004] To address the aforementioned issues, this invention introduces a cross-iterative calculation method for saturated soil within a meshless arbitrary Lagrange-Euler framework, avoiding the problems of ill-conditioned stiffness matrix and difficulty in convergence of iterative solutions. Simultaneously, using relative energy error as the mapping threshold, variable time-step mapping is performed to reduce the computational burden while maintaining accuracy. Furthermore, after the field variable mapping at each Gaussian point is completed, the mapped constitutive variables are corrected based on the stress state of the current calculation step, improving the robustness of large deformation solutions. This invention can provide strong technical support for the large deformation analysis of saturated porous media in practical engineering. Summary of the Invention

[0005] To address the challenge of low efficiency in MFM simulation of saturated porous media, this invention aims to provide a method for large deformation analysis of saturated porous media using MFM, incorporating cross-iteration and variable time-step mapping. This invention avoids the ill-conditioned stiffness matrix and convergence difficulties in iterative solutions by introducing a cross-iteration calculation method for saturated soil within a meshless arbitrary Lagrange-Eulerian framework. Furthermore, by using relative energy error as the mapping threshold and implementing variable time-step mapping, the computational efficiency is significantly improved. This method has significant engineering and academic value for efficient large deformation analysis of saturated porous media.

[0006] To achieve the above objectives, the technical solution adopted by the present invention is as follows: A method for large deformation analysis of MFM-saturated porous media considering cross-iteration and variable time-step mapping includes the following steps: S1 sequentially reads in the node information, material information, and load information of the meshless model, and generates the background mesh. Specifically: S1.1, Read in the node information of the meshless model, including node number and node coordinates; S1.2, then read in the material and load information of the meshless model; S1.3 generates a background mesh that can cover the meshless model in S1.1, and generates Gaussian points within it.

[0007] S2, construct the radial basis and additional basis for each Gaussian point in S1, and calculate the shape function accordingly. Specifically: S2.1, Construct the radial basis of the Gaussian point; The radial basis includes composite quadratic (MQ), Gaussian (EXP), thin plate spline (TPS), or logarithmic (LOG) types, and their formulas are as follows: (1) (2) (3) (4) In the formula, Radial basis; r The distance between Gaussian points and unmesh nodes; d s The average spacing between gridless nodes in the model; a c = 3, q= 1.03, η= 0.98 are all radial basis constant parameters.

[0008] S2.2, Construct additional basis for Gaussian points; The additional base can be either a primary additional base or a secondary additional base, and the formulas are as follows: (5) (6) In the formula, For additional bases; x , y These are the coordinates of the gridless nodes.

[0009] S2.3, based on the radial basis and additional basis in S2.1 and S2.2 respectively, obtain the radial basis matrix formed by the radial basis and the additional basis matrix formed by the additional basis, calculate the shape function of the Gaussian point in S1.3, see formula (7). (7) In the formula, Represents a Gaussian point-shaped function vector; and Let be the constant matrix to be determined; This is the radial basis matrix formed by the radial basis in S2.1; This is the additional basis matrix formed by the additional basis in S2.2.

[0010] S3, based on the shape function in S2, is the discretization of the governing equations for meshless saturated porous media. Specifically: S3.1, based on the incremental equilibrium equation of the meshless saturated porous medium as shown in formula (8), it is first assumed that the estimated values ​​of the meshless nodal displacement vector and the pore pressure increment vector are respectively Δu ’ and Δp ’ Substitute into formula (8) to solve for the displacement vector of the meshless nodes. Δu Velocity increment vector Δv and acceleration increment vector Δa .

[0011] (8) In the formula, M MFM , CMFM , K MFM , K MFM,liq , L MFM These represent the meshless mass matrix, meshless damping matrix, meshless solid stiffness matrix, meshless pore water stiffness matrix, and Laplace operator matrix, respectively. This represents the solid load increment vector.

[0012] S3.2, the velocity increment vector Δ obtained in S3.1 v Substitute the equation into the incremental meshless continuity equation shown in formula (9) and solve for the nodal pore pressure increment vector Δ. p .

[0013] (9) In the formula, S MFM , H MFM These represent the meshless seepage matrix and the meshless compression matrix, respectively; The velocity vector representing the pore pressure increment; Represents the fluid load increment vector; L MFM Represents the Laplace operator matrix; S3.3, the solution obtained in S3.2 Δu and Δp Substitute the equation back into the equilibrium equation shown in formula (8) and iterate repeatedly until the relative deviations of the two are both less than the set error value. error If both formulas (10) and (11) are satisfied, then the solution is considered to meet the iterative accuracy. In this invention... error The recommended value is 0 to 0.02.

[0014] The formula for calculating the relative deviation is as follows: (10) (11) S4, based on the calculation results in S3, uses the relative energy error of the Gaussian point as the criterion to perform variable time-step mapping, specifically: S4.1, Calculate the relative energy error at the Gaussian point, see formula (12): (12) In the formula, abs It is an absolute value function; u and p These are the current displacement vector and pore pressure vector of the node; Indicates increment; D andB Here are the constitutive matrix and the strain-displacement matrix; k For the set energy error threshold, in this invention k The recommended value is 0 to 0.05.

[0015] S4.2 When the relative energy error in S4.1 is less than the set energy error threshold, Gaussian point field variable mapping is carried out to realize the dynamic update of field variables in space. The stress vector and pore pressure vector mapping is shown in formula (13), the strain vector mapping is shown in formula (14), and the constitutive intra-variable vector mapping is shown in formula (15). (13) (14) (8)(15) In the formula, , , , These represent the stress vector, pore pressure vector, strain vector, and constitutive internal variable vector before mapping, respectively. , , , These represent the mapped stress vector, pore pressure vector, strain vector, and constitutive internal variable vector, respectively. It is a mapping function based on radial basis.

[0016] S5 optimizes the constitutive intrastructure variable vector after mapping in S4, then solves and outputs the calculation results, completing the large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step. Specifically: The S5.1 and S4 mapping processes do not consider the compatibility between the constitutive variables at Gaussian points, which may lead to a decrease in solution accuracy or even interruption of calculation. Therefore, based on the stress history, the constitutive variable vector is optimized, as shown in formula (16): (16) In the formula k , k These represent the constitutive internal variable vectors before and after optimization, respectively; Represents stress history; It represents the history of adaptation.

[0017] S5.2, loop through S2 to S5.1, solve and output the results, and thus complete the large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step.

[0018] Compared with the prior art, the beneficial effects of the present invention are as follows: (1) This invention proposes a cross-iteration calculation method for meshless saturated porous media, which avoids the problem of difficulty in convergence of iterative solution caused by ill-conditioning of stiffness matrix and improves calculation efficiency; (2) This invention proposes to carry out variable time step mapping with relative energy error as the mapping threshold, which can reduce unnecessary computational burden while ensuring the accuracy of meshless large deformation simulation; (3) Based on the constitutive internal variables after Gaussian point stress state correction mapping, this invention improves the applicability and robustness of the meshless large deformation method to complex constitutive models. In summary, the meshless method for large deformation analysis of saturated porous media considering cross-iteration and variable time-step mapping can significantly improve computational efficiency while maintaining accuracy. Furthermore, by correcting the constitutive intrinsic variables based on the Gaussian point stress state, the robustness of the solution is enhanced. This invention provides strong technical support for large deformation analysis of saturated porous media in practical engineering. Attached Figure Description

[0019] Figure 1 This is a schematic diagram of the main process of the method of the present invention; Figure 2 This is a geometric model of a typical caisson wharf; Figure 3 This is a meshless node model of a typical caisson wharf; Figure 4 The input earthquake time history curve ( x direction); Figure 5 The input earthquake time history curve ( y direction); Figure 6 The time history curve of the vertical displacement of observation point A of the caisson; Figure 7 The time history curve of the horizontal displacement of observation point A of the caisson; Figure 8 for t= 5.00s Nodal deformation diagram of the caisson wharf; Figure 9 for t= Node deformation diagram of the caisson wharf at 10:00s; Figure 10 for t= 15.00s Nodal deformation diagram of the caisson wharf; Figure 11 for t= Node deformation diagram of the caisson wharf at 25.00s. Detailed Implementation

[0020] The present invention will be further described below with reference to the accompanying drawings and specific embodiments, but the scope of protection of the present invention is not limited thereto.

[0021] See Figure 1 A method for large deformation analysis of MFM-saturated porous media considering cross-iteration and variable time-step mapping includes the following steps: S1, sequentially read in the node information, material information, and load information of the meshless model, and generate the background mesh; S2, construct the radial basis and additional basis for each Gaussian point in S1, and calculate the shape function accordingly; S3, based on the shape function in S2, is the discretized governing equation for meshless saturated porous media. S4, based on the calculation results in S3, uses the relative energy error of the Gaussian point as the criterion to perform variable time step mapping; S5 optimizes the constitutive variables mapped from S4, then solves and outputs the calculation results, completing the large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step.

[0022] To demonstrate the effectiveness of this technical solution, an analysis of the liquefaction deformation problem of a typical caisson wharf was conducted using this technical solution.

[0023] S1, sequentially read in the node information, material information, and load information of the caisson wharf model without mesh, and generate the background mesh, reference. Figure 2 Seabed clay layer height L 1 = 22m, height of clay layer L 2 = 10m x Directional model length L 3 = 185m, length of the replacement sand layer L 4 = 35m, caisson width L 5 = 12m, caisson height L 6 = 18m, specifically: S1.1, Reference Figure 3 Read in the node information of the meshless model, including node number, node coordinates, and the average spacing of the meshless model is 2m; S1.2, then the material and load information of the meshless model is read in. The caisson wharf can be divided into the caisson area, the replacement sand area, the backfill gravelly soil area, the seabed clay area, the clay area, and the backfill gravel area. Among them, the caisson is simulated using linear elastic constitutive modeling, and the caisson density is... ρ Take 2500 kg / m 3 Elastic modulus E = 25500 MPa, Poisson's ratio v=0.167, and observation point A was set at the upper left corner of the caisson. The generalized elastoplastic model was used for the sand replacement area, backfilled gravelly soil area, seabed clay area, clay area, and the back gravel area; parameters are shown in Table 1. The sand replacement area and backfilled gravelly soil area may experience liquefaction, therefore they are considered saturated porous media with permeability coefficients of 0.0001 m / s and 0.0005 m / s, respectively. Other areas are considered single-phase media. Gravity was applied to the model as a body force, and... x Direction input Figure 4 Seismic waves, in y Direction input Figure 5 seismic waves.

[0024] Table 1. Material parameters of the generalized elastoplastic model

[0025] S1.3, Generate a background mesh that can cover the meshless model in S1.1. In this embodiment, the background mesh size is 4m, and Gaussian points are generated inside it. The number of Gaussian points in each background mesh is 2×2.

[0026] S2.1, Construct the radial basis of the Gaussian point; The radial basis includes composite quadratic (MQ), Gaussian (EXP), thin plate spline (TPS), or logarithmic (LOG) types, and their formulas are as follows: (1) (2) (3) (4) In the formula, Radial basis; r The distance between Gaussian points and unmesh nodes; d s The average spacing between gridless nodes in the model; a c = 3, q= 1.03, η= 0.98 are all radial basis constant parameters. In this implementation scheme, a composite quadratic (MQ) radial basis is selected.

[0027] S2.2, Construct additional basis for Gaussian points; The additional base can be either a primary additional base or a secondary additional base, and the formulas are as follows: (5) (6) In the formula, For additional bases; x, y These are the coordinates of the gridless nodes.

[0028] S2.3, based on the radial basis and additional basis in S2.1 and S2.2 respectively, obtain the radial basis matrix formed by the radial basis and the additional basis matrix formed by the additional basis, calculate the shape function of the Gaussian point in S1.3, see formula (7). (7) In the formula, Represents a Gaussian point-shaped function vector; and Let be the constant matrix to be determined; This is the radial basis matrix formed by the radial basis in S2.1; This is the additional basis matrix formed by the additional basis in S2.2.

[0029] S3, based on the shape function in S2, is the discretization of the governing equations for meshless saturated porous media. Specifically: S3.1, based on the incremental equilibrium equation of the meshless saturated porous medium as shown in formula (8), it is first assumed that the estimated values ​​of the meshless nodal displacement vector and the pore pressure increment vector are respectively Δu ’ and Δp ’ Substitute into formula (8) to solve for the displacement vector of the meshless nodes. Δu Velocity increment vector Δv and acceleration increment vector Δa .

[0030] (8) In the formula, M MFM , C MFM , K MFM , K MFM,liq , L MFM These represent the meshless mass matrix, meshless damping matrix, meshless solid stiffness matrix, meshless pore water stiffness matrix, and Laplace operator matrix, respectively. This represents the solid load increment vector.

[0031] S3.2, the velocity increment vector Δ obtained in S3.1 v Substitute the equation into the incremental meshless continuity equation shown in formula (9) and solve for the nodal pore pressure increment vector Δ. p .

[0032] (9) In the formula, SMFM , H MFM These represent the meshless seepage matrix and the meshless compression matrix, respectively; The velocity vector representing the pore pressure increment; Represents the fluid load increment vector; L MFM Represents the Laplace operator matrix; S3.3, the solution obtained in S3.2 Δu and Δp Substitute the equation back into the equilibrium equation shown in formula (8) and iterate repeatedly until the relative deviations of the two are both less than the set error value. error If both formulas (10) and (11) are satisfied, the solution is considered to meet the iteration accuracy. To study the influence of iteration error on the calculation results, the error values ​​in this embodiment are 0.005, 0.01, 0.02, 0.05, and 0.1, respectively. Table 2 compares the relative deviation of the horizontal displacement of point A under each working condition. The results show that when the error is less than 0.02, the relative deviation can be controlled within 2%.

[0033] The formula for calculating the relative deviation is as follows: (10) (11) Table 2. Impact of iteration error on calculation results

[0034] S4, based on the calculation results in S3, uses the relative energy error of the Gaussian point as the criterion to perform variable time-step mapping, specifically: S4.1, Calculate the relative energy error at the Gaussian point, see formula (12): (12) In the formula, abs It is an absolute value function; u and p These are the current displacement vector and pore pressure vector of the node; Indicates increment; D and B Here are the constitutive matrix and the strain-displacement matrix; kThe energy error threshold is set. In this embodiment, to better study the impact of the relative energy error threshold on the calculation results, six working conditions were studied (error threshold = 0%, 1%, 2%, 5%, 10%, 20%). The working condition with an error threshold of 0% was used as the benchmark, and the relative deviation of the horizontal displacement at point A and the normalized time were compared. The results are shown in Table 3. The results show that when the relative energy error threshold is less than 5%, the simulation relative deviation can be controlled within 3%, while improving the calculation efficiency by 15% while ensuring accuracy.

[0035] Table 3. Impact of relative energy error threshold on calculation results

[0036] S4.2 When the relative energy error in S4.1 is less than the set energy error threshold, Gaussian point field variable mapping is carried out to realize the dynamic update of field variables in space. The stress vector and pore pressure vector mapping is shown in formula (13), the strain vector mapping is shown in formula (14), and the constitutive intra-variable vector mapping is shown in formula (15). (13) (14) (8)(15) In the formula, , , , These represent the stress vector, pore pressure vector, strain vector, and constitutive internal variable vector before mapping, respectively. , , , These represent the mapped stress vector, pore pressure vector, strain vector, and constitutive internal variable vector, respectively. It is a mapping function based on radial basis.

[0037] S5 optimizes the constitutive intrastructure variable vector after mapping in S4, then solves and outputs the calculation results, completing the large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step. Specifically: The S5.1 and S4 mapping processes do not consider the compatibility between the constitutive variables at Gaussian points, which may lead to a decrease in solution accuracy or even interruption of calculation. Therefore, based on the stress history, the constitutive variable vector is optimized, as shown in formula (16): (16) In the formula k , k These represent the constitutive internal variable vectors before and after optimization, respectively; Represents stress history; It represents the history of adaptation.

[0038] S5.2 iterates through S2 to S5.1, solving and outputting the results, thus completing the large deformation analysis of MFM-saturated porous media considering cross-iteration and variable time steps. To compare the simulation accuracy and efficiency of the traditional solution framework and the cross-iteration solution framework for meshless saturated porous media, Figure 6 and Figure 7 The time history curves of horizontal and vertical displacements at observation point A during the earthquake are presented. The results show that the calculation results based on cross-iterative solution are consistent with the traditional solution framework, with both methods showing vertical and horizontal displacements of 1.83m and 4.17m respectively. These figures are also largely consistent with the maximum horizontal displacement of 5m and residual settlement of 1-2m at the top of the caisson observed in the earthquake damage survey, verifying the accuracy of the invention. Furthermore, the iterative method used in this invention to solve the governing equations improves computational efficiency by 23% compared to the traditional direct solution (traditional solution: 14.7h; cross-iterative solution: 11.3h).

[0039] Figures 8-11 The typical moments in the earthquake process calculated using this invention are given. t =5 s , t =10 s,t =15 s,t =25 s The nodal deformation diagram of the caisson wharf shows that the present invention can effectively reproduce the development process of wharf destruction and evolution.

[0040] The above embodiments are merely illustrative of the implementation methods of the present invention, but should not be construed as limiting the scope of the present invention. It should be noted that those skilled in the art can make various modifications and improvements without departing from the concept of the present invention, and these modifications and improvements all fall within the protection scope of the present invention.

Claims

1. A method for analyzing large deformations of MFM-saturated porous media considering cross-iteration and variable time-step mapping, characterized in that, The method for analyzing large deformations of MFM-saturated porous media includes the following steps: S1, sequentially read in the node information, material information and load information of the meshless model, and generate the background mesh; S2, construct the radial basis and additional basis for each Gaussian point in S1, and calculate the shape function accordingly; S3, based on the shape function in S2, is the discretized governing equation for meshless saturated porous media. S4, based on the calculation results in S3, uses the relative energy error of the Gaussian point as the criterion to perform variable time step mapping; S5 optimizes the constitutive variable vector after mapping in S4, and then solves and outputs the calculation results to complete the large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step.

2. The method for large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step mapping according to claim 1, characterized in that, Specifically, S1 is: S1.1, Read in the node information of the meshless model, including node number and node coordinates; S1.2, then read in the material and load information of the meshless model; S1.3 generates a background mesh that can cover the meshless model in S1.1, and generates Gaussian points within it.

3. The method for large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step mapping according to claim 2, characterized in that, Specifically, S2 is: S2.1, Construct the radial basis of the Gaussian point; The radial basis is selected from composite quadratic, Gaussian, thin-plate spline, or logarithmic types, and their formulas are as follows: （1） （2） （3） （4） In the formula, Radial basis; r The distance between Gaussian points and unmesh nodes; d s The average spacing between gridless nodes in the model; a c = 3, q= 1.03, η= 0.98 represents radial basis constant parameters; S2.2, Construct additional basis for Gaussian points; The additional base can be either a primary additional base or a secondary additional base, and the corresponding formulas are as follows: （5） （6） In the formula, For additional bases; x , y The coordinates of the gridless nodes; S2.3, based on the radial basis and additional basis in S2.1 and S2.2 respectively, obtain the radial basis matrix formed by the radial basis and the additional basis matrix formed by the additional basis, calculate the shape function of the Gaussian point in S1.3, see formula (7). （7） In the formula, Represents a Gaussian point-shaped function vector; and Let be the constant matrix to be determined; This is the radial basis matrix formed by the radial basis in S2.1; This is the additional basis matrix formed by the additional basis in S2.

2.

4. The method for large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step mapping according to claim 3, characterized in that, Specifically, S3 is: S3.1, based on the incremental equilibrium equation of the meshless saturated porous medium as shown in formula (8), it is first assumed that the estimated values ​​of the meshless nodal displacement vector and the pore pressure increment vector are respectively Δu ’ and Δp ’ Substitute into formula (8) to solve for the displacement vector of the meshless nodes. Δu Velocity increment vector Δv and acceleration increment vector Δa ; （8） In the formula, M MFM , C MFM , K MFM , K MFM,liq , L MFM These represent the meshless mass matrix, meshless damping matrix, meshless solid stiffness matrix, meshless pore water stiffness matrix, and Laplace operator matrix, respectively. Represents the solid load increment vector; S3.2, the velocity increment vector Δ obtained in S3.1 v Substitute into the incremental meshless continuity equation shown in formula (9) and solve for the nodal pore pressure increment vector Δ. p ; （9） In the formula, S MFM , H MFM These represent the meshless seepage matrix and the meshless compression matrix, respectively; The velocity vector representing the pore pressure increment; Represents the fluid load increment vector; L MFM Represents the Laplace operator matrix; S3.3, the solution obtained in S3.2 Δu and Δp Substitute the equation back into the equilibrium equation shown in formula (8) and iterate repeatedly until the relative deviations of the two are both less than the set error value. error If both formula (10) and formula (11) are satisfied, then the solution is considered to meet the iteration accuracy. The formula for calculating the relative deviation is as follows: （10） （11）。 5. The method for large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step mapping according to claim 4, characterized in that, In S3.3 error The value ranges from 0 to 0.

02.

6. The method for large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step mapping according to claim 4, characterized in that, Specifically, S4 is: S4.1, Calculate the relative energy error at the Gaussian point, see formula (12): （12） In the formula, abs It is an absolute value function; u and p These are the current displacement vector and pore pressure vector of the node; Indicates increment; D and B Here are the constitutive matrix and the strain-displacement matrix; k The set energy error threshold; S4.2 When the relative energy error in S4.1 is less than the set energy error threshold, Gaussian point field variable mapping is carried out to realize the dynamic update of field variables in space.

7. The method for large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step mapping according to claim 6, characterized in that, In S4.1 k The value ranges from 0 to 0.

05.

8. The method for large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step mapping according to claim 6, characterized in that, In S4.2, the mapping of stress vector and pore pressure vector is shown in Equation (13), the mapping of strain vector is shown in Equation (14), and the mapping of constitutive variable vector is shown in Equation (15): （13） （14） （8）（15） In the formula, , , , These represent the stress vector, pore pressure vector, strain vector, and constitutive internal variable vector before mapping, respectively. , , , These represent the mapped stress vector, pore pressure vector, strain vector, and constitutive internal variable vector, respectively. It is a mapping function based on radial basis.

9. The method for large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step mapping according to claim 6, characterized in that, Specifically, S5 is: S5.1, Based on the stress history, optimize the constitutive intermolecular variable vector, as shown in formula (16): （16） In the formula k , k These represent the constitutive internal variable vectors before and after optimization, respectively; Represents stress history; Represents the history of response; S5.2, loop through S2 to S5.1, solve and output the results, completing the large deformation analysis of MFM saturated porous media considering cross-iteration and variable time step.

Images