Numerical simulation method for frozen wall formation of underground engineering based on mps

Abstract

The application discloses a numerical simulation method for frozen wall formation of underground engineering based on MPS, and comprises the following steps: a meshless technique is used to disperse media in a simulation area stratum of underground engineering into a particle set with temperature, mass and energy attributes; an improved MPS multi-physical field control equation system for simulating a frozen wall is constructed, a completely implicit viscosity term is introduced for a high-viscosity phase-change fluid, a double-error compensation term is added in a source term of a pressure Poisson equation, and a piecewise enthalpy-temperature phase-change model is coupled with a particle framework; the improved MPS multi-physical field control equation system is solved by time-step iteration to obtain dynamic evolution data of the frozen wall; and the formation time and closed state of the frozen wall are predicted according to the dynamic evolution data. The application solves the problem of the limitation of traditional analysis methods for frozen wall in the excavation process of underground engineering under the conditions of inhomogeneous, anisotropic media and complex boundary conditions.

Description

Technical Field

[0001] This invention relates to the field of frozen wall technology, and more specifically to a numerical simulation method for the formation of frozen walls in underground engineering based on MPS. Background Technology

[0002] In the field of underground engineering, especially during the excavation of deep underground spaces such as coal mines, tunnels, and subways, artificial ground freezing technology is often used to ensure construction safety and structural stability. This technology involves injecting a refrigerant into freezing pipes, lowering the temperature of the surrounding soil below freezing point, causing pore water to freeze and form a frozen wall with a certain strength and impermeability. This method can effectively control groundwater seepage, enhance the mechanical properties of the soil, and provide necessary temporary support for underground space excavation.

[0003] However, the formation of frozen walls during underground engineering excavation involves complex physical processes, including heat conduction and phase change (liquid water transforming into solid ice) effects. These processes not only interact but are also influenced by geological conditions (such as soil type and moisture content), environmental factors (such as temperature changes), and engineering operating parameters (such as freezing rate and cooling system configuration). Traditional analytical methods, such as empirical formulas and simplified models, have limitations in handling heterogeneous, anisotropic media and complex boundary conditions, making it difficult to accurately predict the development process of frozen walls and their impact on the surrounding environment. Summary of the Invention

[0004] To overcome the shortcomings of existing technologies, a numerical simulation method for the formation of frozen walls in underground engineering based on MPS is provided to address the limitations of traditional analysis methods for frozen walls during underground engineering excavation under heterogeneous, anisotropic media and complex boundary conditions.

[0005] To achieve the above objectives, a numerical simulation method for the formation of frozen walls in underground engineering based on MPS is provided, comprising the following steps: The meshless technique is used to discretize the medium in the strata of the simulated area of ​​underground engineering into a collection of particles with temperature, mass, and energy properties; An improved MPS multiphysics field governing equation system was constructed for simulating frozen walls. A fully implicit viscosity term was introduced for high-viscosity phase change fluids, and a double error compensation term was added to the source term of the pressure Poisson equation. At the same time, the piecewise enthalpy-temperature phase change model was coupled with the particle frame. The improved MPS multiphysics control equation system is solved iteratively step by step to obtain the dynamic evolution data of the frozen wall; Based on the dynamic evolution data, the formation time and closure state of the frozen wall are predicted.

[0006] Furthermore, the medium includes soil, water, and ice.

[0007] Furthermore, the momentum equation of the improved MPS multiphysics governing equation system employs a fully implicit viscosity term: ; in, For particles i Temporary velocity vector at the prediction step; For particles i In the k The known velocity vector at the end of each time step; v The kinematic viscosity of the fluid; For time step; d The number of spatial dimensions; n 0 The initial standard particle number density; l These are the correction coefficients for discretization of the Laplace operator; w(˙) For kernel functions; For particles i With neighboring particles j In the k The spatial distance at the end of each time step.

[0008] Furthermore, the pressure Poisson equation for the improved MPS multiphysics governing equation system is: ; in, r For fluid density; P For pressure field; ( 2 P) K+1 For the first k+1 The pressure Laplace operator for each time step; n ** The temporary particle number density is the particle's movement after an explicit prediction step. For time step; β This is the first error compensation coefficient used to compensate for the truncation error when discretizing the time rate of change of particle number density; c This is the second error compensation coefficient used to compensate for the cumulative error in particle number density deviating from the standard value.

[0009] Furthermore, the segmented enthalpy-temperature phase transition model coupled with the particle framework includes each particle independently maintaining its enthalpy and temperature, and determining the phase transition stage in real time through an inverse relationship, which is: Liquid phase T > 0℃ ; Freezing stage T=0℃, ; Solid state stage T < 0℃ ; in, T The current temperature of the particle; T initial The initial temperature of the particle; h This is the current enthalpy value of the particle; h initial For particles at the initial temperature T initial The reference enthalpy value below; m The mass of the particle; C water is the specific heat capacity of liquid water; x The freezing ratio; h water The enthalpy of pure liquid water at 0°C; L The latent heat of phase transition of water; h ice The enthalpy of pure solid ice at 0°C; C ice is the specific heat capacity of solid ice.

[0010] The beneficial effects of this invention are that the MPS-based numerical simulation method for the formation of frozen walls in underground engineering overcomes the limitations of mesh distortion, achieving accurate simulation of large deformation freezing processes. This invention employs a meshless MPS method to discretize continuous strata into particle sets. Unlike the traditional mesh-dependent finite element method, the MPS method does not require a fixed mesh structure and does not suffer from mesh distortion. During the freezing process, when the soil undergoes frost heave deformation, the particles can move freely without re-meshing, flexibly adapting to changes in complex geometries and large deformation situations. Simultaneously, the MPS method has unique advantages in capturing free surface motion, phase changes, and solid deformation, accurately tracking the movement trajectory of the freezing front, overcoming the fundamental defect of traditional mesh methods that lead to simulation failure due to mesh distortion when dealing with large deformation problems. Attached Figure Description

[0011] Other features, objects, and advantages of this application will become more apparent from the following detailed description of non-limiting embodiments with reference to the accompanying drawings: Figure 1 This is a flowchart illustrating the numerical simulation method for the formation of frozen walls in underground engineering based on MPS, according to an embodiment of the present invention.

[0012] Figure 2 This diagram illustrates the three stages involved in the freezing process of water. Detailed Implementation

[0013] The present application will now be described in further detail with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the invention and not intended to limit it. Furthermore, it should be noted that, for ease of description, only the parts relevant to the invention are shown in the accompanying drawings.

[0014] It should be noted that, unless otherwise specified, the embodiments and features described in this application can be combined with each other. This application will now be described in detail with reference to the accompanying drawings and embodiments.

[0015] This invention provides a numerical simulation method for the formation of frozen walls in underground engineering based on MPS, comprising the following steps: S1. Meshless technology is used to discretize the medium in the strata of the underground engineering simulation area into a set of particles with temperature, mass and energy properties.

[0016] In this embodiment, the medium includes soil, water, and ice.

[0017] Meshless technology simplifies the computation process by using particle number density instead of fluid density in incompressible fluids, making the algorithm more flexible and efficient.

[0018] The continuous strata are discretized into a set of particles with definite spatial coordinates, and the position vector of each particle is determined. r i And the initial distance between particles.

[0019] S2. Construct an improved MPS multiphysics field governing equation system for simulating frozen walls. Introduce a fully implicit viscosity term for high-viscosity phase change fluids and add a double error compensation term to the source term of the pressure Poisson equation. At the same time, couple the piecewise enthalpy-temperature phase change model with the particle frame.

[0020] Specifically, the momentum equation of the improved MPS multiphysics governing equation system adopts a fully implicit viscosity term: ; in, For particlesi Temporary velocity vector at the prediction step; For particles i In the k The known velocity vector at the end of each time step; v The kinematic viscosity of the fluid; For time step; d The number of spatial dimensions; n 0 The initial standard particle number density; l These are the correction coefficients for discretization of the Laplace operator; w(˙) For kernel functions; For particles i With neighboring particles j In the k The spatial distance at the end of each time step.

[0021] The pressure Poisson equation for the improved MPS multiphysics governing equation system is: ; in, r For fluid density; P For pressure field; ( 2 P) K+1 For the first k+1 The pressure Laplace operator for each time step; n ** The temporary particle number density is the particle's movement after an explicit prediction step. For time step; β This is the first error compensation coefficient used to compensate for the truncation error when discretizing the time rate of change of particle number density; c This is the second error compensation coefficient used to compensate for the cumulative error in particle number density deviating from the standard value.

[0022] The segmented enthalpy-temperature phase transition model coupled with the particle framework includes each particle independently maintaining its enthalpy and temperature, and determining the phase transition stage in real time through an inverse relationship, which is: Liquid phase T > 0℃ ; Freezing stage T=0℃, ; Solid state stage T < 0℃ ; in, T The current temperature of the particle; T initial The initial temperature of the particle; h This is the current enthalpy value of the particle; h initial For particles at the initial temperature T initial The reference enthalpy value below; m The mass of the particle; C water is the specific heat capacity of liquid water; x The freezing ratio; h water The enthalpy of pure liquid water at 0°C; L The latent heat of phase transition of water; h ice The enthalpy of pure solid ice at 0°C; C ice is the specific heat capacity of solid ice.

[0023] In this embodiment, the governing equations of the basic MPS method are the mass, momentum, and energy conservation equations: ; in, r , t , u , P , v , F , h , K , T , Q These are the fluid densities (kg / m³) 3 Time, velocity vector (m / s), pressure (Pa), kinematic viscosity (m³ / s) 2 / s), external force per unit mass (such as gravity, m / s) 2 The parameters are: enthalpy (J / kg), thermal conductivity, thermal diffusivity (m² / s), temperature (K or °C), and heat source (low-temperature freezing tube heat source (W / kg), simulating the cold injection or latent heat release of the freezing tube). In the MPS method, the convection term is directly incorporated into the calculation via moving particles. All differential operators in the governing equations are replaced by the following particle interaction model.

[0024] Common kernel functions used for discretizing differential operators are described as follows: ; in, r e It is the cutoff radius used to limit interactions. r It is the distance between two particles. (In position) r i At that point, the particle i The particle number density is defined as the sum of the weighting functions: ; n i For particles i The particle number density at that location is used to replace the fluid density.

[0025] n 0 The initial standard particle number density (calculated from the initial uniform arrangement).

[0026] Fluid density is proportional to particle number density. For incompressible fluids, the standard values ​​evaluated under the initial particle configuration are used as a substitute in the calculations. n i .

[0027] arbitrary scalar The gradient operator is defined as a particle. i Rather than at the cutoff radius r e Local weighted average of gradient vectors among neighboring particles within the range: ; in, For any scalar field (such as pressure, temperature). d It refers to the number of spatial dimensions. The Laplace operator considers it as a diffusion process towards neighboring particles: ; Among them, the coefficient is introduced l This makes the increase in the variance of the Laplace operator equal to the variance of the analytical solution: ; Substituting the above particle interaction model into the governing equations generates the particle dynamics equations to simulate fluid flow, ensuring that the variance of the discrete Laplace operator matches the analytical solution. No mesh is required in this discretization process.

[0028] A semi-implicit algorithm is employed in the MPS method. The mass and momentum equations are solved in two sub-steps at each time step. In the first step, the temporal velocity and position are obtained by explicitly solving the momentum equation (excluding the pressure gradient term). However, after completing the first step, the incompressibility of the fluid is not satisfied because the particle number density deviates from the standard value calculated from the initial particle configuration. In the second step, the following Poisson equation for pressure is derived from the pressure gradient in the mass and momentum equations to adjust the particle number density to the standard value: ; in, This represents the time increment. The pressure field is updated by implicitly solving the Poisson equation, and then the particle's velocity and position are corrected by solving the pressure gradient term. In the algorithm, the particle's temperature is explicitly calculated using the discretized energy equation.

[0029] In the MPS method, particles that satisfy the following conditions are considered free surface particles and are given Dirichlet boundaries: n i < αn 0 ; in, α It is a constant value < 1.0.

[0030] For high-viscosity fluids, such as materials in the phase transition stage, the original MPS method cannot be effectively used for simulation because the time increment is determined not only by the Courant stability criterion but also by the diffusion stability criterion. To satisfy the diffusion stability criterion, the time increment must be very small due to the high fluid viscosity. However, excessively small time increments lead to very high computational costs and cause numerical instability when solving the Poisson equation for pressure. To address this issue, it is necessary to implicitly calculate the viscosity term. In this embodiment, the following fully implicit viscosity term calculation technique is employed: ; in, The velocity at the current time step.

[0031] The temporary velocity for the predicted step (after explicit movement) is the unknown quantity to be solved.

[0032] , This is the temporary position after an explicit move.

[0033] superscript kThe asterisks and * represent the time step in the calculation and a time substep within a time step, respectively. Furthermore, the above equation can be considered a positive definite matrix, with the particle velocity of the last time step as the source term. This matrix is ​​solved using the incomplete Cholesky conjugate gradient method.

[0034] The following is an explicit calculation of enthalpy using the MPS Laplace model: .

[0035] While the MPS method has certain advantages in simulating free surfaces and multiphase flows, spurious pressure oscillations are an inherent drawback that can significantly affect fluid properties. In this embodiment, the pressure Poisson equation is improved using the following formula, where two error compensation parts are added to the source terms of the pressure Poisson equation to compensate for the truncation error generated during discretization: ; in n ** This represents the particle number density after shifting the particle over time. n k , n k-1 : The current and previous particle number density.

[0036] n 0 Initial standard particle number density.

[0037] Δt : Time step (s).

[0038] β and c These are coefficients used to adjust the compensation portion; the relationship between them should be maintained at 0 ≤ c ≤ β ≤ 1。

[0039] Another cause of pressure oscillations is the singularity in the kernel function used in the original MPS at the target particle. To address this issue, a similar kernel function without singularities was employed: .

[0040] The density difference between solids and liquids is small compared to the density difference between gases and liquids. Therefore, the phase transition is simply simulated using only two particle types with different properties, and volume changes are not considered during solidification / melting. Temperature is a function of enthalpy per unit volume, as shown below: ; in, T mIt is the melting point. C p It is specific heat. h 0 and h 1 These are the enthalpy at the start and end of melting, respectively. (Subscript) s and l These represent solid and liquid states, respectively.

[0041] In phase transition modeling, the solid fraction is introduced. c To identify the state of particles, such as Figure 2 As shown, the calculation method is as follows: ; solid particles ( γ=1 Particles that remain stationary while undergoing a phase transition (0 < 0) c < 1 Then, by solving the Navier-Stokes equations, it moves in the same way as other fluid particles.

[0042] The freezing process of water involves three stages: liquid water cooling to 0°C, phase change freezing at 0°C, and continued cooling of the ice. The relationship between enthalpy h and temperature T needs to be described in stages, and the specific formula is as follows: (1) Liquid water cools down from high temperature to 0°C (liquid stage).

[0043] When the temperature of water is above 0°C, its enthalpy change is determined by its specific heat capacity. Assume the initial temperature is... T initial ＞0℃ The mass is m, and the specific heat capacity is C water enthalpy h With temperature T The relationship is: h(T)=h initial -m·C water (T initial -T) ; Formula meaning: Temperature from T initial Down to T At that time, the decrease in enthalpy is: m·C water · T ; Parameter value: Specific heat capacity of water C water ≈4.1868KJ / (kg˙K); h initial This is the reference enthalpy value at the initial temperature.

[0044] (2) Phase change freezing of water at 0°C (liquid-solid coexistence stage).

[0045] Once water reaches 0°C, the temperature remains constant, but the enthalpy decreases with the degree of freezing. Let the freezing ratio be denoted as . ( x =0 indicates that it is completely water. x =1 indicates that it is completely ice, then the enthalpy is h(x)=h water -m·L·x .

[0046] Formula meaning: h water It is the enthalpy of water at 0°C (usually taken as...). m ˙ C water˙ T initial ).

[0047] L It is the latent heat of phase change of water (the heat released per unit mass). .

[0048] When x=1, the enthalpy reaches the value after complete freezing. h ice =h water -m·L .

[0049] (3) Ice cools down from 0°C to low temperature (solid stage).

[0050] If frozen ice continues to be cooled (below 0°C), the enthalpy change is determined by the specific heat capacity of the ice. Let the temperature be: T At temperatures below 0°C, the specific heat capacity of ice is... C ice The total enthalpy is h(T)=h ice -m˙C ice ˙(0-T) .

[0051] Parameter value: Specific heat capacity of ice: C ice ≈2.1KJ / (kg˙K) .

[0052] h ice It is the enthalpy of ice at 0°C (i.e., h water -mL ).

[0053] The above stages are combined into a single piecewise function: ; h initial Initial temperature T initial Reference enthalpy (J) below.

[0054] h water The enthalpy of liquid water at 0°C h water = h initial -mC water ( T initial - 0).

[0055] h ice The enthalpy of ice at 0°C h ice = h water -mL .

[0056] x For the frozen ratio ( x =0 total water, x =1 (all ice), during the phase transition stage T It remains at 0℃.

[0057] In this embodiment, temperature T With enthalpy h The inverse relationship, if enthalpy is known. h Temperature can be determined by reverse calculation. T Stage at which: 1. Liquid stage T >0°C, ; 2. Freezing Phase T =0°C, ; 3 solid-state stage T <0°C, .

[0058] The above components together constitute the complete final physical and mathematical model. The mass, momentum, and energy equations define the macroscopic physical laws; the MPS discrete operator transforms the partial differential equations into algebraic operations between particles; the improved PPE and implicit viscosity treatment ensure the numerical stability and feasibility of high-viscosity phase change flow with large time steps; and the piecewise phase change constitutive model accurately characterizes the release of latent heat of freezing and temperature evolution.

[0059] S3. The dynamic evolution data of the frozen wall is obtained by iteratively solving the improved MPS multiphysics field control equation system step by step.

[0060] The transient evolution of the frozen wall is simulated using a semi-implicit time-stepping algorithm, specifically including: Explicit prediction step: Without considering pressure gradient, explicitly calculate temporary values ​​of particle velocity and position (based on viscous force, gravity, etc.).

[0061] Implicit correction step: Calculate the deviation between the temporary particle number density and the initial standard density.

[0062] Solve the improved pressure Poisson equation to obtain the pressure field.

[0063] The incompressibility is restored by using pressure gradients to correct particle velocity and position.

[0064] Energy update: Explicitly calculate the energy equation, update the enthalpy and temperature of each particle, and update the solid fraction and physical properties (such as viscosity and density) according to the phase transition model.

[0065] The movement of the freezing front, the evolution of the temperature field, and the advancement of the phase transition interface are simulated by calculating the gradient and Laplace operator through inter-particle interactions (kernel function weighting).

[0066] S4. Based on dynamic evolution data, predict the formation time and closure state of the frozen wall.

[0067] The thickness, closure time, and integrity of the frozen wall are determined based on the particle temperature field and solid fraction distribution.

[0068] Determine whether the frozen wall meets the design requirements for strength and impermeability.

[0069] It provides quantitative basis for optimizing the layout of freezing pipes, predicting construction period and safety control, and is especially suitable for heterogeneous soil layers and complex boundary conditions.

[0070] See Figure 1 As shown, the specific process of the MPS-based numerical simulation method for the formation of frozen walls in underground engineering according to the present invention is as follows: a. Input initial parameters ( , , , ).

[0071] The continuous strata are discretized into a set of particles with temperature, mass, and energy properties, and the initial position of each particle is defined. ,speed ,pressure ,temperature and enthalpy .

[0072] Calculate the initial standard particle number density based on the initial uniform distribution: This value serves as the benchmark for subsequent incompressibility criteria.

[0073] Set calculation parameters: time step Δt Total simulation duration, cutoff radius r e Kernel function form, physical properties (specific heat, latent heat, viscosity, etc.).

[0074] b. Calculate the energy equation and phase transition: h k+1 =h k + (k 2 T k +Q) .

[0075] At the beginning of each time step, the energy conservation equation is explicitly solved using the temperature and enthalpy fields from the previous time step, updating the enthalpy and temperature of each particle. Simultaneously, the solid fraction is updated according to the piecewise phase transition constitutive relation. c .

[0076] Based on the updated enthalpy value, the new temperature is determined using the inverse relationship. and solid fraction γ i (γ=1 represents the solid state, γ=0 represents the liquid state, and 0<γ<1 represents the phase transition state). Solid particles will no longer participate in the solution of the momentum equation.

[0077] c. Implicit calculation of viscosity term solution: .

[0078] For high-viscosity fluids (especially phase transition states), a fully implicit scheme is used to solve for the viscous force contribution, obtaining a temporary velocity that does not contain pressure gradients or external forces.

[0079] d. Explicitly calculate surface tension and gravity terms: .

[0080] Based on the temporary velocity obtained in step c, the effects of gravity g and surface tension F are explicitly added to obtain the final predicted velocity.

[0081] e. Estimate particle position: .

[0082] Using the predicted velocity obtained in step d, the particle position is explicitly advanced to obtain temporary coordinates. .

[0083] f. Implicit calculation of the Poisson equation for pressure: .

[0084] Calculation of temporary particle number density based on temporary coordinates n ** And combined with the current step density n k and the density of the previous step n k−1 Substituting the pressure Poisson equation, which contains double error compensation terms, the new pressure field is implicitly solved. P k+1 .

[0085] g. Update instantaneous velocity and position: .

[0086] Using the pressure gradient obtained in step f, the predicted velocity in step d and the predicted position in step e are implicitly corrected to obtain a new velocity that satisfies incompressibility. and new location .

[0087] in, P Discrete calculations are performed using the MPS gradient operator.

[0088] h. Determine the termination condition. If the condition is met, the process ends; otherwise, return to step b.

[0089] Check if the preset total number of simulation time steps has been reached or if the frozen wall thickness meets the closure criteria. If so, terminate the loop and proceed to the result output; otherwise, update the particle state ( u k+1 ,r k+1 ,T k+1 ,c Use this as the initial value for the next time step, and return to step b to continue iterating.

[0090] After the cycle is completed, based on the final saved particle temperature field, solid fraction distribution and historical evolution data, key engineering indicators such as frozen wall thickness, closure time and freezing front movement trajectory are quantitatively extracted to provide a basis for design optimization and construction safety.

[0091] This invention presents a numerical simulation method for the formation of frozen walls in underground engineering based on MPS, overcoming the limitations of mesh distortion and achieving accurate simulation of large deformation freezing processes. This invention employs a meshless MPS method to discretize continuous strata into a collection of particles. Unlike the traditional mesh-dependent finite element method, the MPS method does not require a fixed mesh structure and does not suffer from mesh distortion. During the freezing process, when the soil undergoes frost heave deformation, the particles can move freely without re-meshing, flexibly adapting to changes in complex geometries and large deformation situations. Simultaneously, the MPS method has unique advantages in capturing free surface motion, phase changes, and solid deformation, accurately tracking the movement trajectory of the freezing front, and overcoming the fundamental defect of traditional mesh methods that lead to simulation failure due to mesh distortion when dealing with large deformation problems.

[0092] This invention provides a numerical simulation method for the formation of frozen walls in underground engineering based on MPS, overcoming time step limitations and significantly improving the computational efficiency for high-viscosity phase change fluids. This invention employs a fully implicit viscous term handling technique in the momentum equation. The standard MPS method uses an explicit scheme for the viscous term, and its time step is strictly limited by the diffusion stability criterion. Δt≤(Δx) 2 / (2 n) During freezing, the water-ice phase transition causes a sharp increase in fluid viscosity. If an explicit scheme is used, the time step is forced to the microsecond level, making engineering-scale freezing simulations (lasting days to months) computationally infeasible. This invention constructs the viscosity term as a positive definite matrix and solves it using the incomplete Cholesky conjugate gradient method, freeing the time step from viscosity limitations. This significantly improves the computational efficiency of simulating high-viscosity phase transition fluids, making efficient numerical simulation of the entire freezing wall formation process possible.

[0093] The numerical simulation method for frozen wall formation in underground engineering based on MPS in this invention effectively suppresses spurious pressure oscillations, significantly improving simulation stability and reliability. This invention innovatively adds two error compensation terms (coefficients β and γ, satisfying 0 ≤ γ ≤ β ≤ 1) to the source term of the pressure Poisson equation. Spurious pressure oscillations are an inherent defect of the MPS method, significantly affecting fluid characteristics and contaminating the prediction accuracy of the temperature field and phase transition interface. This invention, by introducing a dual error compensation mechanism into the source term, effectively compensates for the truncation error generated during discretization, significantly suppresses non-physical high-frequency pressure oscillations, obtains a smoother and more stable pressure field, and ensures the reliability of the predicted frozen wall thickness and closure time.

[0094] This invention presents a numerical simulation method for the formation of frozen walls in underground engineering based on MPS, achieving real-time linkage between the thermodynamic and mechanical states at the particle level and accurately characterizing the entire phase transition process. This invention deeply couples a segmented enthalpy-temperature phase transition constitutive model with the MPS particle framework. Each particle independently maintains its enthalpy and temperature, and its phase transition stage is determined in real time through an inverse relationship. A solid fraction γ is introduced to achieve linkage between the thermodynamic and mechanical states of the particles: solid particles with γ=1 remain stationary, while phase transition particles with 0<γ<1 move according to the Navier-Stokes equations. This parallel computational architecture, where "each particle independently judges the phase transition progress," allows the spatial evolution of the frozen front to be accurately tracked at the particle level, rather than relying on interpolation estimation as in grid methods, thus achieving precise characterization of the release of latent heat of freezing and temperature evolution.

[0095] The MPS-based numerical simulation method for the formation of frozen walls in underground engineering, as presented in this invention, has a wide range of applications, providing a reliable scientific basis for the design of frozen engineering projects. This invention can be applied to various complex geological conditions and engineering scenarios. Regardless of the soil medium with different thermophysical properties or different freezing pipe layout schemes, this method can effectively simulate and analyze them. It can accurately predict key engineering indicators such as frozen wall thickness, closure time, temperature field distribution, and the movement trajectory of the freezing front, providing quantitative basis for optimizing freezing pipe layout schemes, predicting construction periods, and ensuring safety control. This effectively solves the technical problem of limitations in traditional analysis methods when dealing with heterogeneous, anisotropic media and complex boundary conditions.

[0096] The above description is merely a preferred embodiment of this application and an explanation of the technical principles employed. Those skilled in the art should understand that the scope of the invention involved in this application is not limited to technical solutions formed by specific combinations of the above-described technical features, but should also cover other technical solutions formed by arbitrary combinations of the above-described technical features or their equivalents without departing from the inventive concept. For example, technical solutions formed by substituting the above features with (but not limited to) technical features with similar functions disclosed in this application.

Claims

1. A numerical simulation method for the formation of frozen walls in underground engineering based on MPS, characterized in that, Includes the following steps: The meshless technique is used to discretize the medium in the strata of the simulated area of ​​underground engineering into a collection of particles with temperature, mass, and energy properties; An improved MPS multiphysics field governing equation system was constructed for simulating frozen walls. A fully implicit viscosity term was introduced for high-viscosity phase change fluids, and a double error compensation term was added to the source term of the pressure Poisson equation. At the same time, the piecewise enthalpy-temperature phase change model was coupled with the particle frame. The improved MPS multiphysics control equation system is solved iteratively step by step to obtain the dynamic evolution data of the frozen wall; Based on the dynamic evolution data, the formation time and closure state of the frozen wall are predicted.

2. The numerical simulation method for the formation of frozen walls in underground engineering based on MPS according to claim 1, characterized in that, The medium includes soil, water, and ice.

3. The numerical simulation method for the formation of frozen walls in underground engineering based on MPS according to claim 1, characterized in that, The momentum equation of the improved MPS multiphysics governing equation system adopts a fully implicit viscosity term: ； in, For particles i Temporary velocity vector at the prediction step; For particles i In the k The known velocity vector at the end of each time step; v The kinematic viscosity of the fluid; ∆t For time step; d The number of spatial dimensions; n 0 The initial standard particle number density; λ These are the correction coefficients for discretization of the Laplace operator; w(˙) For kernel functions; For particles i With neighboring particles j In the k The spatial distance at the end of each time step.

4. The numerical simulation method for the formation of frozen walls in underground engineering based on MPS according to claim 1, characterized in that, The pressure Poisson equation for the improved MPS multiphysics control system is: ； in, ρ For fluid density; P For pressure field; (∇ 2 P) K+1 For the first k+1 The pressure Laplace operator for each time step; n ** The temporary particle number density is the particle's movement after an explicit prediction step. ∆t For time step; β This is the first error compensation coefficient used to compensate for the truncation error when discretizing the time rate of change of particle number density; γ This is the second error compensation coefficient used to compensate for the cumulative error in particle number density deviating from the standard value.

5. The numerical simulation method for the formation of frozen walls in underground engineering based on MPS according to claim 1, characterized in that, The segmented enthalpy-temperature phase transition model coupled with the particle framework includes each particle independently maintaining its enthalpy and temperature, and determining the phase transition stage in real time through an inverse relationship, which is: Liquid phase T > 0℃ ; Freezing stage T=0℃, ; Solid state stage T < 0℃ ; in, T The current temperature of the particle; T initial The initial temperature of the particle; h This is the current enthalpy value of the particle; h initial For particles at the initial temperature T initial The reference enthalpy value below; m The mass of the particle; C water is the specific heat capacity of liquid water; x The freezing ratio; h water The enthalpy of pure liquid water at 0°C; L The latent heat of phase transition of water; h ice The enthalpy of pure solid ice at 0°C; C ice is the specific heat capacity of solid ice.

Images