Thermal-mechanical coupling frozen soil deformation rapid prediction method based on material point method

By using an improved Takashi 3D frost heave model and the thermo-mechanical coupling of the material point method, the problems of local large deformation and the coupling of temperature field and stress field in frozen soil deformation prediction are solved, realizing rapid and accurate frozen soil deformation prediction, which is applicable to infrastructure design and construction in cold regions.

CN120995939APending Publication Date: 2025-11-21HARBIN INST OF TECH
View PDF 0 Cites 0 Cited by

Patent Information

Application Number
CN202511166188.1
Authority / Receiving Office
CN · China
Patent Type
Applications(China)
Current Assignee / Owner
Filing Date
2025-08-20
Publication Date
2025-11-21

AI Technical Summary

Technical Problem

Existing technologies cannot quickly and accurately reflect the large local deformation of frozen soil under the coupling effect of temperature changes and external loads, and do not fully consider the coupling between temperature field and stress field, resulting in low accuracy and low computational efficiency in frozen soil deformation prediction.

Method used

An improved Takashi three-dimensional frost heave model was adopted, anisotropic parameters were introduced, and the frost heave rate was calculated by normalized direction factor. The thermo-mechanical coupling was combined with the material point method to establish a rapid prediction method for thermo-mechanical coupled frozen soil deformation. The combined influence of temperature field and stress field was considered, and the frost heave strain was calculated using Navier equation and constitutive relation.

Benefits of technology

It achieves adaptability to consider complex stress states in multidimensional frozen soil deformation prediction, improves prediction accuracy and computational efficiency, and can quickly and accurately predict frozen soil frost heave deformation, especially large frost heave deformation caused by water migration under sufficient water supply conditions.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN120995939A_ABST
    Figure CN120995939A_ABST
Patent Text Reader

Abstract

The invention relates to the technical field of frozen soil deformation prediction, and particularly discloses a thermal-mechanical coupling frozen soil deformation rapid prediction method based on a material point method, and the method comprises the steps: building an improved Takashi three-dimensional frost heaving model; thermal-mechanical coupling frozen soil deformation calculation: carrying out initialized background grid boundary division on the predicted frozen soil position; calculating the temperature gradient according to the node temperature, calculating the heat flux density of the material point, and updating the node temperature and the temperature of the material point; judging whether the material point reaches a frost heaving condition or not according to the updated material point temperature, and if the material point reaches the frost heaving condition, performing frost heaving strain calculation; calculating the strain rate and the rotation rate of the substance point, updating the volume of the substance point by using the strain rate of the substance point, and updating the stress of the substance point by using the constitutive relation; calculating node force and updating node momentum, applying displacement boundary conditions again, mapping node acceleration and velocity back to the material point, and updating the velocity and position of the material point; and repeating the calculation of the next time step by using the undeformed grid.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This invention relates to the field of frozen soil deformation prediction technology, specifically to a rapid prediction method for frozen soil deformation based on the material point method and thermo-mechanical coupling. Background Technology

[0002] With the increasing demand for infrastructure construction in cold regions (such as railways, highways, and oil pipelines), engineering defects caused by permafrost deformation are gradually increasing during construction. Frost heave of permafrost is the core issue causing these defects. Uneven frost heave of permafrost subgrades may lead to pavement cracking, and frost heave of bridge foundations may cause structural instability. Therefore, accurately predicting the deformation characteristics of permafrost under the coupled action of temperature changes and external loads is crucial for the design, construction, and operation and maintenance of projects in cold regions.

[0003] Currently, the prediction of permafrost deformation mainly includes methods such as empirical estimation, physical experiments, and numerical simulation. Empirical estimation is based on field observation data or indoor test data, and establishes empirical formulas between deformation characteristics and influencing factors through regression analysis. However, this method has low prediction accuracy and a narrow scope of application. Physical experiment methods mainly simulate the deformation process of permafrost using equipment such as indoor permafrost triaxial apparatus and frost heave apparatus, and can directly obtain deformation data. However, the test cycle is long, the cost is high, and it is difficult to reproduce the temperature field distribution and stress transfer law at the engineering scale. Numerical simulation is the main method for studying permafrost deformation. It establishes governing equations and simulates the deformation of permafrost in a computer. Commonly used numerical calculation methods include the finite element method and the finite difference method. However, in the case of permafrost deformation, the coupling effect of multiple physical fields needs to be considered, and its computational efficiency is relatively low and it cannot solve the problem of large local deformation.

[0004] The material point method, developed from the mass-point grid method, employs both Lagrangian and Eulerian descriptions, discretizing the object as a set of mass points moving within a spatial grid. Each mass point carries all material information, such as mass, velocity, strain, and stress. The motion equations are solved on the spatial grid, avoiding grid distortion problems and making it suitable for extremely large deformations and flow problems, overcoming the limitations of the two methods. However, existing research on frozen soil deformation based on the material point method still has shortcomings: most studies only focus on mechanical response, failing to fully consider the coupling between the temperature and stress fields; and the computational efficiency of solving traditional frozen soil deformation using the material point method is low, failing to meet the "rapid prediction" requirements in engineering practice. Summary of the Invention

[0005] The purpose of this invention is to provide a rapid prediction method for thermo-coupling frozen soil deformation based on the material point method, so as to solve the technical problems that the existing technology cannot reflect large local deformation and does not fully consider the coupling of temperature field and stress field.

[0006] The objective of this invention can be achieved through the following technical solutions: A rapid prediction method for frozen soil deformation based on thermo-mechanical coupling using the material point method includes the following steps: An improved Takashi three-dimensional frost heave model was established: anisotropic parameters were introduced to achieve multidimensional frost heave prediction. The frost heave rate in different directions was allocated through parameters, and the parameters in different directions were calculated through normalized direction factors, taking into account the influence of material, stress and heat flow on the anisotropic coefficient. Calculation of thermally coupled frozen soil deformation: The background grid boundary is initialized for the predicted permafrost location, and the node velocity and node temperature of each grid node are calculated. The temperature gradient is calculated based on the node temperature, and the heat flux density of the material point is calculated, where the material point represents the actual position coordinates of the node. Update the node temperature and material point temperature based on the heat flow at the nodes; Determine whether the material point has reached the frost heave condition based on the updated material point temperature. If the frost heave condition has been reached, then perform frost heave strain calculation. The velocity gradient of the material point is obtained using the background mesh, and then the strain rate and curl of the material point are calculated. The volume of the material point is updated using the strain rate of the material point and the stress of the material point is updated using the constitutive relation. Calculate nodal forces and update nodal momentum. After applying displacement boundary conditions again, map nodal acceleration and velocity back to material points and update the velocity and position of material points. The calculation for the next time step is repeated using the undeformed mesh.

[0007] As a further aspect of the present invention: in the improved Takashi three-dimensional frost heave model, the frost heave rate in different directions Through parameters Distribute, and simultaneously in different directions The sum is 1. By normalizing the direction factor In the two-dimensional model, the specific calculation formula is as follows: ; ; ; ; Among them, parameters This represents the dominant role of heat flow in the directional distribution of frost heave, obtained by the ratio of the heat flow direction under isotropic compression to the lateral frost heave; the material parameter k controls the anisotropic parameters. The degree of variation with stress is used to consider the impact of constraint stress on the distribution of frost heave in different types of soil; It is the effective stress in the direction of heat flow. It is an effective constraint stress perpendicular to the heat flow direction. This represents the average value of the principal effective stress.

[0008] As a further aspect of the present invention: the process of initializing the background grid boundary division for the predicted permafrost location specifically includes: Initialize the background mesh boundary, map the mass, momentum, heat capacity, and heat of the material points onto the background mesh nodes through shape functions, and apply displacement boundary conditions and temperature boundary conditions to the nodes; The calculation formula is as follows: ; ; ; ; Where, m g Indicates node mass; m p Indicates the mass of a point substance; This represents the value of the shape function of node g at material point p; This represents the momentum of a node; C represents the velocity of a point mass; g c represents the heat capacity of a node. p Q represents the specific heat capacity of a point in a substance; g T represents the heat at the node; p Indicates the temperature of a point in the substance; n p The total number of discrete material points; subscript p represents the variable of each material point; subscript g represents the variable of each node.

[0009] As a further aspect of the present invention, the specific steps for calculating the node velocity and node temperature of each grid node include: The calculation formula is as follows: ; ; in, Indicates the node velocity; This indicates the junction temperature.

[0010] As a further aspect of the present invention: calculating the temperature gradient based on the node temperature and calculating the heat flux density of the material point specifically includes: ; Where q represents the heat flux density vector; Indicates thermal conductivity; Represents the updated value; This represents the Nabla operator.

[0011] As a further aspect of the present invention: updating the node temperature and material point temperature based on the heat flow at the node specifically includes: The calculation formula for nodal heat flux, including internal and external heat flux, is as follows: ; Among them, F g Indicates the effect of nodal heat flow; V p This represents the volume of the substance point p; Indicates the intensity of the internal heat source; The nodal temperature is updated using the nodal heat flux effect, and the temperature boundary condition is reapplied. The calculation formula is as follows: ; in, This indicates the updated node temperature; Indicates the time step; The nodal heat flux is mapped back to the material point and the material point temperature is updated using the following formula: ; in, This indicates the updated temperature of the material point.

[0012] As a further aspect of the present invention, the specific steps for calculating frost heave strain are as follows: The freezing rate U is calculated using the following formula: ; ; ; Where s represents the freezing distance; Indicates the duration of the time step; U represents the temperature gradient; U represents the freezing rate. The frost heave rate was calculated based on the improved Takashi formula. ; The Navier equations are used as the governing equations for the stress field, and the equations are as follows: ; Where C is the material stiffness; u is the displacement vector; and F is the volume force vector. Represents the Nabla operator; The strain-displacement relationship can be expressed as follows: ; in, For strain tensor; The frost heave strain tensor was calculated as follows: ; in, It is the frost heave strain tensor. It is elastic strain; It is frost heave strain; The stiffness matrix D of the soil is: .

[0013] As a further aspect of the present invention: the velocity gradient of the material point is obtained using the background mesh, thereby calculating the strain rate and curl of the material point. The specific steps for updating the volume of the material point using the strain rate and updating the stress of the material point using the constitutive relation are as follows: The velocity gradient of the material points is obtained using the background mesh, and the strain rate and curl of the material points are calculated using the following formulas: ; ; in, The subscripts i and j represent the strain rate at a material point; they also represent the spatial coordinate components. , Represents the velocity components of a substance point and Partial derivatives with respect to spatial coordinates; Indicates the spin rate; , Represents the velocity components of the background mesh nodes; Update the volume of a material point using its strain rate: ; in, This represents the updated volume of the material point; Indicates the time step; Represents the diagonal component of the strain rate; Represents the diagonal component of the thermal strain rate; Update the point stress of the material using constitutive relations; ; ; Where T is the operator; This represents the updated point stress of the material.

[0014] As a further aspect of the present invention: calculating nodal forces and updating nodal momentum, applying displacement boundary conditions again, mapping nodal accelerations and velocities back to material points, and updating the velocity and position of the material points, specifically includes: Update nodal momentum using nodal forces Apply displacement boundary conditions again: ; in, This represents the updated node momentum; Map the nodal accelerations and velocities back to the matter points and update the velocities and positions of the matter points: ; ; in, Indicates the updated position of the matter point; This indicates the updated velocity of the matter point.

[0015] The beneficial effects of this invention are: The model can consider the impact of complex stress states on frost heave, has strong adaptability to complex scenarios, and can make good predictions of multidimensional frost heave in permafrost. The improved Takashi model can effectively simulate frost heave characteristics under various temperature and stress boundary conditions, accurately predict frost heave deformation caused by water migration when water supply is sufficient, and further improve the prediction accuracy for complex scenarios.

[0016] The model has relatively few parameters with clear physical meaning, resulting in high computational efficiency and lowering the application threshold. Since the Takashi model itself already encompasses the impact of moisture migration on frost heave deformation, the influence of the moisture field on frost heave deformation is simplified, eliminating the need for additional calculations. Only the total frost heave needs to be predicted using the improved Takashi equation. After inputting physical parameters such as temperature and stress, the calculation process is convenient and efficient, quickly yielding frozen soil deformation prediction results and enabling rapid calculation of frost heave. Attached Figure Description

[0017] The invention will now be further described with reference to the accompanying drawings.

[0018] Figure 1 This is a graph showing the variation of anisotropic parameters observed in an embodiment of the present invention; Figure 2 This is a flowchart illustrating the rapid prediction method for thermo-coupling frozen soil deformation based on the material point method of this invention. Detailed Implementation

[0019] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0020] This invention, based on the material point method for predicting frozen soil deformation, considers the influence of temperature and stress fields on frozen soil deformation, introduces an improved three-dimensional Takashi model, establishes a thermo-mechanical coupling analysis framework, and proposes a rapid prediction method for frozen soil deformation based on the material point method. This model can effectively capture the frost heave characteristics under different temperature and stress boundary conditions. Furthermore, the model can effectively predict large frost heave deformation caused by moisture migration under conditions of sufficient water supply. The good fit between the new method and experimental data verifies the rationality of this invention in multidimensional frost heave prediction. This method also has the advantages of fewer parameters and clear physical meanings for the parameters.

[0021] To address the problem of predicting frozen soil deformation, a rapid prediction method based on the material point method is proposed. This method can quickly predict frozen soil deformation under the combined influence of temperature and stress fields, solving the problems of slow calculation and inability to reflect large local deformations in traditional frozen soil deformation prediction.

[0022] refer to Figure 2 The present invention aims to achieve rapid prediction of frozen soil deformation based on the material point method using thermo-coupling, applicable to the design, construction, operation, and maintenance of infrastructure in cold regions. The following is an implementation method of this approach: 1. Establishment of the improved Takashi 3D frost heave model: Since the Takashi equation is a one-dimensional frost heave prediction model that only considers the frost heave rate in the direction of heat flow, an anisotropic parameter is introduced. To achieve multidimensional frost heave prediction. Frost heave rate in different directions. Through parameters Distribute, and simultaneously in different directions The sum is 1. The anisotropy coefficient is calculated using a normalized direction factor β, which takes into account the effects of material properties, stress, and heat flow. For the two-dimensional model, The parameter μ represents the dominant role of heat flow in the directional distribution of frost heave, and can be obtained by the ratio of the heat flow direction to the lateral frost heave under isotropic compression. The material parameter k controls the anisotropy parameter. The degree of variation with stress is used to account for the impact of constraint stress on the distribution of frost heave in different types of soil. It is the effective stress in the direction of heat flow. It is an effective constraint stress perpendicular to the heat flow direction. This represents the average value of the principal effective stress.

[0023] ; ; ; ; To verify the applicability of the proposed anisotropy parameters to actual multidimensional frost heave conditions, existing frost heave test results were compared and analyzed with the proposed calculation method. and The ratio is used as the x-axis. Figure 1 The results show that when k=3.2 and μ=5.2, the improved model accurately captures the anisotropic parameter changes observed in the experiment, and the fitting results verify the applicability and accuracy of the proposed model. The fitting of the experimental data also provides a reference range for the parameter values ​​of the numerical simulation model for specific frozen soils. Therefore, the improved model can effectively reflect the evolution law of frost heave anisotropy with the change of constraint stress under experimental conditions, and has the advantages of fewer parameters and clear physical meaning of each parameter.

[0024] 2. A thermo-coupling calculation method for frozen soil deformation based on the material point method: (1) Initialize the background mesh boundary, map the mass, momentum, heat capacity and heat of the material points onto the background mesh nodes through shape functions, and apply displacement boundary conditions and temperature boundary conditions to the nodes.

[0025] The calculation formula is as follows: ; ; ; ; Where, m g Indicates node mass; m p Indicates the mass of a point substance; This represents the value of the shape function of node g at material point p; This represents the momentum of a node; C represents the velocity of a point mass; g c represents the heat capacity of a node. p Q represents the specific heat capacity of a point in a substance; g T represents the heat at the node; p Indicates the temperature of a point in the substance; n p The total number of discrete material points; subscript p represents the variable of each material point; subscript g represents the variable of each node.

[0026] (2) Calculate the node velocity using the node mass and momentum, and calculate the node temperature using the node heat capacity and heat.

[0027] The calculation formula is as follows: ; ; in, Indicates the node velocity; This indicates the junction temperature.

[0028] (3) Calculate the temperature gradient using the nodal temperature and update the heat flux density of the material point.

[0029] The calculation formula is as follows: ; Where q represents the heat flux density vector; Indicates thermal conductivity; Represents the updated value; This represents the Nabla operator.

[0030] (4) Calculate the heat flow at the nodes, including the heat flow inside the nodes and the heat flow outside the nodes.

[0031] The calculation formula is as follows: ; Among them, F g Indicates the effect of nodal heat flow; V p This represents the volume of the substance point p; Indicates the intensity of the internal heat source; (5) Update the nodal temperature using the nodal heat flow and apply the temperature boundary conditions again.

[0032] The calculation formula is as follows: ; in, This indicates the updated node temperature; Indicates the time step; (6) Map the nodal heat flow back to the material point and update the material point temperature.

[0033] The calculation formula is as follows: ; in, This indicates the updated temperature of the material point.

[0034] (7) Determine whether the material point is frozen based on the updated material point temperature. If the freezing heave condition is met, calculate the freezing heave strain (thermal strain). First, calculate the temperature gradient to obtain its magnitude, and then calculate the freezing rate U according to the following formula.

[0035] ; ; ; Where s represents the freezing distance; Indicates the duration of the time step; U represents the temperature gradient; U represents the freezing rate. In the case of freezing, the projection of stress along the temperature gradient direction and the magnitude of stress perpendicular to the temperature gradient direction can be calculated to obtain the stress ratio, which is controlled within the range of 0.25-4.0.

[0036] The frost heave distribution coefficients in different directions are calculated using the following formulas: ; ; ; The frost heave rate is calculated based on the improved Takashi formula, as follows: ; The Navier equations are used as the governing equations for the stress field, and the equations are as follows: ; Where C is the material stiffness; u is the displacement vector; and F is the volume force vector. Represents the Nabla operator; The strain-displacement relationship can be expressed as follows: ; in, For strain tensor; The strain caused by frost heave should be considered as thermal strain. The simplified linear elastic constitutive equation can be derived as follows, and the frost heave strain tensor can be calculated: ; ; in, is the frost heave strain tensor; D is the soil stiffness matrix; It is elastic strain; It is frost heave strain. Calculations were performed based on the improved Takashi model, and simultaneously equal , which is the thermal strain rate in the calculation process below.

[0037] (8) Use the background grid to obtain the velocity gradient of the material point, and then calculate the strain rate and rotation rate of the material point. Thermal strain needs to be considered in the strain rate.

[0038] ; ; in, The subscripts i and j represent the strain rate at a material point; they also represent the spatial coordinate components. , Represents the velocity components of a substance point and Partial derivatives with respect to spatial coordinates; Indicates the spin rate; , Represents the velocity components of the background mesh nodes; (9) Update the volume of a material point using the strain rate of the material point.

[0039]

[0040] in, This represents the updated volume of the material point; Indicates the time step; Represents the diagonal component of the strain rate; Represents the diagonal component of the thermal strain rate; (10) Update the point stress of the material using constitutive relations.

[0041] ; ; Where T is the operator; This represents the updated point stress of the material.

[0042] (11) Calculate nodal forces This includes internal nodal forces and external nodal forces.

[0043] ; in, Indicates physical strength; Represents surface force; (12) Update nodal momentum using nodal forces Then, apply displacement boundary conditions again.

[0044] ; in, This represents the updated node momentum; (13) Map the nodal acceleration and velocity back to the material point and update the velocity and position of the material point.

[0045] ; ; in, Indicates the updated position of the matter point; This indicates the updated velocity of the matter point.

[0046] (14) Repeat the calculation for the next time step using the undeformed mesh.

[0047] The foregoing has provided a detailed description of one embodiment of the present invention, but this description is merely a preferred embodiment and should not be construed as limiting the scope of the invention. All equivalent variations and modifications made within the scope of the present invention should still fall within the scope of the present invention.

Claims

1. A rapid prediction method for frozen soil deformation based on the material point method using thermo-mechanical coupling, characterized in that, Includes the following steps: An improved Takashi three-dimensional frost heave model was established: anisotropic parameters were introduced to achieve multidimensional frost heave prediction. The frost heave rate in different directions was allocated through parameters, and the parameters in different directions were calculated through normalized direction factors, taking into account the influence of material, stress and heat flow on the anisotropic coefficient. Calculation of thermally coupled frozen soil deformation: The background grid boundary is initialized for the predicted permafrost location, and the node velocity and node temperature of each grid node are calculated. The temperature gradient is calculated based on the node temperature, and the heat flux density of the material point is calculated, where the material point represents the actual position coordinates of the node. Update the node temperature and material point temperature based on the heat flow at the nodes; Determine whether the material point has reached the frost heave condition based on the updated material point temperature. If the frost heave condition has been reached, then perform frost heave strain calculation. The velocity gradient of the material point is obtained using the background mesh, and then the strain rate and curl of the material point are calculated. The volume of the material point is updated using the strain rate of the material point and the stress of the material point is updated using the constitutive relation. Calculate nodal forces and update nodal momentum. After applying displacement boundary conditions again, map nodal acceleration and velocity back to material points and update the velocity and position of material points. The calculation for the next time step is repeated using the undeformed mesh.

2. The rapid prediction method for thermo-coupling frozen soil deformation based on the material point method according to claim 1, characterized in that, In the improved Takashi three-dimensional frost heave model, the frost heave rate in different directions Through parameters Distribute, and simultaneously in different directions The sum is 1. By normalizing the direction factor In the two-dimensional model, the specific calculation formula is as follows: ; ; ; ; Here, parameter μ represents the dominant role of heat flow in the directional distribution of frost heave, and is obtained by the ratio of the heat flow direction under isotropic compression to the lateral frost heave; material parameter k controls the anisotropic parameters. The degree of variation with stress is used to consider the impact of constraint stress on the distribution of frost heave in different types of soil; It is the effective stress in the direction of heat flow. It is an effective constraint stress perpendicular to the heat flow direction. This represents the average value of the principal effective stress.

3. The rapid prediction method for thermo-coupling frozen soil deformation based on the material point method according to claim 1, characterized in that, The process of initializing the background mesh boundary for the predicted permafrost location specifically includes: Initialize the background mesh boundary, map the mass, momentum, heat capacity, and heat of the material points onto the background mesh nodes through shape functions, and apply displacement boundary conditions and temperature boundary conditions to the nodes; The calculation formula is as follows: ; ; ; ; Where, m g Indicates node mass; m p Indicates the mass of a point substance; This represents the value of the shape function of node g at material point p; This represents the momentum of a node; C represents the velocity of a point mass; g c represents the heat capacity of a node. p Q represents the specific heat capacity of a point in a substance; g T represents the heat at the node; p Indicates the temperature of a point in the substance; n p The total number of discrete material points; subscript p represents the variable of each material point; subscript g represents the variable of each node.

4. The rapid prediction method for frozen soil deformation based on the material point method using thermo-mechanical coupling as described in claim 3, characterized in that, The specific steps for calculating the nodal velocity and nodal temperature of each grid node include: The calculation formula is as follows: ; ; in, Indicates the node velocity; This indicates the junction temperature.

5. The rapid prediction method for thermo-coupling frozen soil deformation based on the material point method according to claim 4, characterized in that, Calculating the temperature gradient based on the nodal temperatures and the heat flux density at the material points specifically includes: ; Where q represents the heat flux density vector; Indicates thermal conductivity; Represents the updated value; This represents the Nabla operator.

6. The rapid prediction method for frozen soil deformation based on the material point method using thermo-coupling as described in claim 5, characterized in that, The updating of nodal temperatures and material point temperatures based on the heat flow at the nodes specifically includes: The calculation formula for nodal heat flux, including internal and external heat flux, is as follows: ; Among them, F g Indicates the effect of nodal heat flow; V p This represents the volume of the substance point p; Indicates the intensity of the internal heat source; The nodal temperature is updated using the nodal heat flux effect, and the temperature boundary condition is reapplied. The calculation formula is as follows: ; in, This indicates the updated node temperature; Indicates the time step; The nodal heat flux is mapped back to the material point and the material point temperature is updated using the following formula: ; in, This indicates the updated temperature of the material point.

7. The rapid prediction method for thermo-coupling frozen soil deformation based on the material point method according to claim 6, characterized in that, The specific steps for calculating frost heave strain are as follows: The freezing rate U is calculated using the following formula: ; ; ; Where s represents the freezing distance; Indicates the duration of the time step; U represents the temperature gradient; U represents the freezing rate. The frost heave rate was calculated based on the improved Takashi formula. ; The Navier equations are used as the governing equations for the stress field, and the equations are as follows: ; Where C is the material stiffness; u is the displacement vector; and F is the volume force vector. Represents the Nabla operator; The strain-displacement relationship can be expressed as follows: ; in, For strain tensor; The frost heave strain tensor was calculated as follows: ; in, It is the frost heave strain tensor. It is elastic strain; It is frost heave strain; The stiffness matrix D of the soil is: 。 8. The rapid prediction method for thermo-coupling frozen soil deformation based on the material point method according to claim 7, characterized in that, The velocity gradient of a material point is obtained using the background mesh, and then the strain rate and curl of the material point are calculated. The specific steps for updating the volume of the material point using the strain rate and updating the stress of the material point using the constitutive relation are as follows: The velocity gradient of the material points is obtained using the background mesh, and the strain rate and curl of the material points are calculated using the following formulas: ; ; in, The subscripts i and j represent the strain rate at a material point; they also represent the spatial coordinate components. , Represents the velocity components of a substance point and Partial derivatives with respect to spatial coordinates; Indicates the spin rate; , Represents the velocity components of the background mesh nodes; Update the volume of a material point using its strain rate: ; in, This represents the updated volume of the material point; Indicates the time step; Represents the diagonal component of the strain rate; Represents the diagonal component of the thermal strain rate; Update the point stress of the material using constitutive relations; ; ; Where T is the operator; This represents the updated point stress of the material.

9. The rapid prediction method for thermo-coupling frozen soil deformation based on the material point method according to claim 3, characterized in that, Calculate nodal forces and update nodal momentum. After reapplying displacement boundary conditions, map nodal accelerations and velocities back to material points and update the velocities and positions of the material points. Specifically, this includes: Calculate nodal forces : ; in, Indicates physical strength; Represents surface force; Update nodal momentum using nodal forces Apply displacement boundary conditions again: ; in, This represents the updated node momentum; Map the nodal accelerations and velocities back to the matter points and update the velocities and positions of the matter points: ; ; in, Indicates the updated position of the matter point; This indicates the updated velocity of the matter point.