Multi-physics field deformation analysis method and system based on COMSOL platform and particle finite element

Abstract

The invention relates to a multi-physical-field deformation analysis method and system based on a COMSOL platform and a particle finite element. The method comprises the steps that a physical deformation analysis model is built, and an initial finite element grid is generated; solving the control equation by using a solver of the COMSOL platform to obtain a multi-physical-field state variable; carrying out particle space position updating on finite element grid nodes according to the multi-physical field state variables, and carrying out Delou inner triangulation and boundary identification based on an updated particle point set to obtain a free boundary; reconstructing the grid based on the free boundary, and mapping the state variable of the previous grid to the new finite element grid; according to the new finite element mesh, updating boundary conditions of the physical deformation analysis model, and iteratively updating the physical deformation analysis model in a COMSOL platform; and outputting a deformation analysis result based on the final physical deformation analysis model. According to the invention, stable and high-precision numerical analysis of large-deformation and multi-field coupling processes in energy geotechnical engineering is realized.

Description

Technical Field

[0001] This invention relates to the field of physical field deformation analysis technology, and in particular to a multiphysics deformation analysis method and system based on the COMSOL platform and particle finite element method. Background Technology

[0002] Numerical analysis methods based on the Finite Element Method (FEM) have become indispensable computational tools in geotechnical engineering. These methods can simulate complex boundary conditions and multiphysics coupling processes, and have been widely applied in energy geotechnical engineering fields such as energy geology, geothermal development, natural gas hydrate extraction, carbon dioxide sequestration, and energy pile heat exchange.

[0003] However, when dealing with significant changes in geometry or large deformations in soil, the traditional finite element method (FEM) is prone to severe mesh distortion, leading to decreased computational accuracy or even solution failure. To address this issue, existing techniques have proposed methods such as the arbitrary Lagrange-Euler method, the coupled Euler-Lagrange method, and the small-strain remesh interpolation method to improve mesh quality and computational stability to some extent. However, these methods still suffer from drawbacks such as computational instability, complex reconstruction, and low computational efficiency when handling large geometric deformations, free surface evolution, and strongly nonlinear material behavior.

[0004] Furthermore, while existing COMSOL multiphysics platforms can efficiently solve coupled thermal, hydrological, and mechanical problems, most are based on a fixed-mesh finite element framework, which struggles to maintain numerical stability under severe deformation. Consequently, they cannot accurately describe phenomena such as large displacements, free surface motion, and material flow, often requiring complex user subroutines or re-meshing algorithms, resulting in low computational efficiency and limited applicability. For example, in problems such as geothermal reservoir rupture, interaction between subsea pipelines and soil, or thermally induced landslides, traditional FEM models are prone to losing convergence due to mesh distortion.

[0005] Therefore, there is an urgent need for a multiphysics large deformation analysis scheme that can implement the particle finite element method in the COMSOL environment to improve its applicability and computational efficiency in energy geotechnical engineering. Summary of the Invention

[0006] (a) Technical problems to be solved

[0007] In view of the above-mentioned shortcomings and deficiencies of the prior art, the present invention provides a multiphysics deformation analysis method and system based on the COMSOL platform and particle finite element method, which solves the technical problems of severe mesh distortion, computational instability and complex implementation of the existing finite element method when dealing with large geometric deformation, free surface evolution and strongly coupled multiphysics process in energy geotechnical engineering.

[0008] (II) Technical Solution

[0009] To achieve the above objectives, the main technical solutions adopted by the present invention include:

[0010] In a first aspect, embodiments of the present invention provide a multiphysics deformation analysis method based on the COMSOL platform and particle finite element method, including:

[0011] Based on the multiphysics constraints of the object to be analyzed, a physical deformation analysis model is established in COMSOL, and an initial finite element mesh is generated.

[0012] The multiphysics control equations are solved at each time step using the solver of the COMSOL platform to obtain the multiphysics state variables at the current time step.

[0013] The nodes of the finite element mesh are used as a set of particle points. The spatial position of each particle is updated by combining the state variables of the multiphysics field. Based on the updated set of particle points, Delaunay triangulation is performed. The boundary of the obtained triangulated convex domain is identified to obtain the free boundary of the current computational domain.

[0014] The finite element mesh of the current time step is reconstructed based on the free boundary, and the multiphysics state variables of the initial finite element mesh or the finite element mesh of the previous time step are mapped onto the reconstructed finite element mesh.

[0015] The boundary conditions of the physical deformation analysis model are updated based on the reconstructed finite element mesh, and the physical deformation analysis model is iteratively updated in the COMSOL platform until the preset termination condition is reached.

[0016] Based on the final updated physical deformation analysis model, the visualized deformation analysis results under the multi-physics coupling effect of the target are output.

[0017] Optionally, based on the multiphysics constraints of the object to be analyzed, a physical deformation analysis model is established in COMSOL, and an initial finite element mesh is generated, including:

[0018] Obtain the multiphysics constraints of the object to be analyzed. The multiphysics constraints include the geometric domain, boundary conditions, material parameters, and multiphysics modules.

[0019] Based on the geometric domain and material parameters, a geometric model of the object to be analyzed is constructed in the COMSOL platform. Based on boundary conditions and the multiphysics module, physical field constraints are applied to the geometric model to generate a physical deformation analysis model.

[0020] The physical deformation analysis model is meshed using finite element methods to generate an initial finite element mesh. The attribute information of the nodes in the computational domain is then exported through a real-time data interaction interface. The attribute information includes node coordinates, element connection information, and initial state variables.

[0021] Optionally, the multiphysics governing equations are solved using the solver on the COMSOL platform at each time step to obtain the multiphysics state variables for the current time step, including:

[0022] Based on the constitutive relations and boundary conditions of the physical deformation analysis model, the multiphysics field control equations to be solved are constructed. The multiphysics field control equations include the momentum balance equation describing mechanical equilibrium, the displacement and strain relationship equation describing deformation motion, and the heat transfer equation describing thermal effects.

[0023] The displacement field variables are spatially discretized using quadratic shape functions, and the auxiliary field variables containing stress and temperature gradients are spatially discretized using linear shape functions.

[0024] While performing spatial discretization, backward difference formulas are used to discretize all field variables in time, so as to transform the multiphysics field control equations into a set of nonlinear equations that can be solved numerically.

[0025] The built-in nonlinear solver of the COMSOL platform is invoked to solve the nonlinear equation system at each time step until a solution that meets the convergence condition is obtained. This solution is a multi-physics state variable that includes displacement, stress, and temperature fields at the current time step.

[0026] Alternatively, the momentum balance equation is:

[0027] ;

[0028] In the formula, Here, σ is the gradient operator, b is the stress tensor, ρ is the body force vector, and a is the density.

[0029] The equation relating displacement and strain is:

[0030] ;

[0031] In the formula, ε is the strain tensor and u is the displacement;

[0032] The heat transfer equation is:

[0033] ;

[0034] In the formula, C p q is the specific heat capacity under constant stress, T is the absolute temperature, q is the conduction heat flux, and Q is the source term.

[0035] 5. The method as described in claim 1, characterized in that, the nodes of the finite element mesh are used as a set of particle points, the spatial positions of each particle are updated in conjunction with multiphysics state variables, and Delaunay triangulation is performed based on the updated set of particle points. Boundary identification is performed on the obtained triangulated convex domain to obtain the free boundary of the current computational domain, including:

[0036] Obtain the multiphysics state variables obtained from the current time step from the COMSOL platform, and extract the displacement field data of the nodes;

[0037] Each node of the initial finite element mesh is regarded as a particle, and the coordinate position of each particle in space is updated based on the displacement field data to form an updated particle point set;

[0038] The updated particle point set is subjected to Delaunay triangulation to generate a triangulated convex domain covering all particles. The triangulated convex domain is composed of a set of non-overlapping triangular units.

[0039] The Alpha Shape Algorithm is used to identify the boundaries of the triangulated convex domain, thereby obtaining the free boundaries of the current computational domain.

[0040] Optionally, the Alpha Shape Algorithm is used to identify the boundaries of the triangulated convex domain to obtain the free boundaries of the current computational domain, including:

[0041] Iterate through the circumcircles of each triangular unit generated after Delaunay triangulation and determine the relationship between the radius of each circumcircle and the preset radius threshold.

[0042] Triangular elements with a circumcircle radius greater than or equal to a radius threshold are classified as boundary elements, and triangular elements with a circumcircle radius less than a radius threshold are classified as internal elements.

[0043] The triangle edges of all boundary cells are filtered out, and the internal edges that are shared by two boundary cells are removed to obtain the boundary edges that belong to the outermost contour of the computational domain.

[0044] Connect the boundary edges to generate a closed and non-convex boundary polygon, which is the free boundary of the current computational domain.

[0045] Optionally, reconstructing the finite element mesh for the current time step based on the free boundary, and mapping the multiphysics state variables of the initial finite element mesh or the finite element mesh from the previous time step onto the reconstructed finite element mesh includes:

[0046] Using the free boundary as a geometric constraint, a distance function-driven mesh generation strategy is employed to reconstruct the finite element mesh for the current time step within the free boundary region.

[0047] By comparing the geometric relationships of the finite element meshes before and after reconstruction, a mapping relationship between the nodes of the old and new meshes is established.

[0048] Based on the mapping relationship, a unique element mapping strategy is adopted to map the multiphysics state variables of the initial finite element mesh or the finite element mesh of the previous time step to the finite element mesh reconstructed in the current time step.

[0049] After completing the mapping of multiphysics state variables, the updated multiphysics state variables are used as the initial conditions for solving the multiphysics control equations in the next time step and loaded into the physical deformation analysis model of the COMSOL platform.

[0050] Optionally, the boundary conditions of the physical deformation analysis model are updated based on the reconstructed finite element mesh, and the physical deformation analysis model is iteratively updated in the COMSOL platform until the preset termination conditions are met, including:

[0051] Based on the node coordinates of the reconstructed finite element mesh, obtain the location of the free boundary of the current computational domain;

[0052] The boundary conditions defined in the original physical deformation analysis model are reapplied based on the current geometric position of the free boundary. The boundary conditions include mechanical boundary conditions and thermal boundary conditions.

[0053] In the COMSOL platform, the physical deformation analysis model with reapplied boundary conditions and completed state variable mapping is updated to the solution model for the current time step;

[0054] Determine whether the current iteration state has reached or exceeded the preset termination condition;

[0055] If the preset termination condition is not met, the current time step is taken as the new initial state, the time step size is incremented, and the process returns to the step of solving the multiphysics control equations in each time step using the solver of the COMSOL platform, and then iterative calculation is performed for the next time step.

[0056] If the preset termination condition has been met or exceeded, the iteration process will terminate and the final updated physical deformation analysis model will be output.

[0057] Optionally, based on the final updated physical deformation analysis model, the output of visualized deformation analysis results under the multiphysics coupling of the target includes:

[0058] Obtain the final updated physical deformation analysis model from the COMSOL platform;

[0059] Post-processing of the physical deformation analysis model generates a result file that can intuitively reflect the deformation process and final state of the object under the coupling of multiple physics fields.

[0060] The results file is rendered and displayed using preset visualization tools, and a visualization chart of the deformation analysis results is output.

[0061] Secondly, embodiments of the present invention provide a multiphysics deformation analysis system based on the COMSOL platform and particle finite element method, comprising:

[0062] The model initialization and mesh generation module is used to establish a physical deformation analysis model in COMSOL and generate an initial finite element mesh based on the multiphysics constraints of the object to be analyzed.

[0063] The multiphysics calculation module is used to solve the multiphysics control equations at each time step using the solver of the COMSOL platform, and obtain the multiphysics state variables at the current time step.

[0064] The geometric boundary identification module is used to take the nodes of the finite element mesh as a set of particles, update the spatial position of each particle in combination with the multiphysics state variables, and perform Delaunay triangulation based on the updated set of particles. It then identifies the boundary of the obtained triangulated convex domain to obtain the free boundary of the current computational domain.

[0065] The mesh reconstruction and variable mapping module is used to reconstruct the finite element mesh of the current time step based on the free boundary, and map the multiphysics state variables of the initial finite element mesh or the finite element mesh of the previous time step onto the reconstructed finite element mesh.

[0066] The boundary condition update and model iteration module is used to update the boundary conditions of the physical deformation analysis model based on the reconstructed finite element mesh, and iteratively update the physical deformation analysis model in the COMSOL platform until the preset termination condition is reached.

[0067] The results visualization output module is used to output visualized deformation analysis results of the target under multi-physics coupling based on the final updated physical deformation analysis model.

[0068] (III) Beneficial Effects

[0069] This invention proposes a multiphysics deformation analysis method based on the COMSOL platform and particle finite element method, which effectively integrates the geometric adaptability of the particle method with the physical modeling advantages of the finite element method, demonstrating significant benefits in the field of energy and geotechnical engineering. This method treats finite element mesh nodes as moving particles. After solving the multiphysics governing equations at each time step, it updates the particle positions based on the displacement field and dynamically identifies geometric boundaries based on Delaunay triangulation, achieving a natural description of the evolution and large deformation processes of free surfaces. Compared to the traditional finite element method, this method avoids numerical instability caused by mesh distortion, eliminates the need for complex remeshing techniques or user-defined subroutines, and significantly improves computational robustness.

[0070] Secondly, by integrating the particle finite element framework into the COMSOL platform, the method described in this invention can directly call the platform's built-in multiphysics solver and material model, balancing computational accuracy and engineering practicality. Especially in strongly coupled scenarios, this method can effectively simulate large material displacements, free surface motion, and multiphysics coupling behavior, providing a reliable tool for deformation prediction in energy geotechnical engineering.

[0071] Furthermore, this invention improves numerical efficiency while maintaining computational continuity through a mesh reconstruction and state variable mapping strategy, providing a scalable solution for long-term deformation analysis and multi-field coupled simulation in practical engineering. Attached Figure Description

[0072] Figure 1 A flowchart illustrating a multiphysics deformation analysis method based on the COMSOL platform and particle finite element method, provided as an embodiment of the present invention;

[0073] Figure 2 This is a schematic diagram of the COMSOL platform and particle finite element calculation framework provided in an embodiment of the present invention;

[0074] Figure 3 This is a schematic diagram of a beam analysis numerical model provided in an embodiment of the present invention;

[0075] Figure 4 This is a schematic diagram illustrating the analysis results of beam deformation on a single COMSOL platform under different stress distributions at different times, according to an embodiment of the present invention.

[0076] Figure 5 This is a schematic diagram illustrating the analysis results of beam deformation under different stress distributions at different times according to an embodiment of the present invention.

[0077] Figure 6 A schematic diagram of parameter configuration for bearing capacity verification analysis of a strip foundation provided in an embodiment of the present invention;

[0078] Figure 7 This is a schematic diagram comparing the numerical solution provided by the COMSOL platform and the particle finite element framework with other numerical solutions according to an embodiment of the present invention.

[0079] Figure 8 A schematic diagram showing the corresponding distribution of cumulative plastic strain invariants and shear stress according to an embodiment of the present invention;

[0080] Figure 9 This is a schematic diagram illustrating the analysis and demonstration of cohesive soil collapse verification according to an embodiment of the present invention.

[0081] Figure 10This is a schematic diagram showing the distribution of plastic strain at four moments during the cohesive soil collapse verification analysis process provided in an embodiment of the present invention.

[0082] Figure 11 This is a schematic diagram illustrating the thermal effect of a two-dimensional plate in the verification analysis of slope failure provided in an embodiment of the present invention.

[0083] Figure 12 This is a schematic diagram illustrating thermally induced slope instability in a slope failure verification analysis according to an embodiment of the present invention.

[0084] Figure 13 Numerical results of thermally induced degradation in the verification analysis of slope failure provided in an embodiment of the present invention;

[0085] Figure 14 This is a numerical result of the post-instability evolution in the verification analysis of slope failure provided in an embodiment of the present invention;

[0086] Figure 15 This is a schematic diagram illustrating the verification and analysis of foundation settlement caused by heat transfer through buried pipelines in clay, provided in an embodiment of the present invention.

[0087] Figure 16 This is a schematic diagram illustrating the evolution of active layer properties caused by a hot pipe, according to an embodiment of the present invention. Detailed Implementation

[0088] To better explain and facilitate understanding of the present invention, the present invention will be described in detail below with reference to the accompanying drawings and specific embodiments.

[0089] refer to Figures 1 to 16As shown in the embodiment of the present invention, a multiphysics deformation analysis method based on the COMSOL platform and particle finite element method is proposed. This method can be applied to the dynamic analysis of large deformation under multiphysics coupling conditions of soil and rock. The method includes: establishing a physical deformation analysis model in COMSOL according to the multiphysics constraints of the object to be analyzed, and generating an initial finite element mesh; using the solver of the COMSOL platform to solve the multiphysics control equations in each time step to obtain the multiphysics state variables of the current time step; using the nodes of the finite element mesh as a particle point set, updating the spatial position of each particle in combination with the multiphysics state variables, and based on... The updated particle point set is subjected to Delaunay triangulation, and the boundary of the obtained triangulated convex domain is identified to obtain the free boundary of the current computational domain. Based on the free boundary, the finite element mesh of the current time step is reconstructed, and the multiphysics state variables of the initial finite element mesh or the finite element mesh of the previous time step are mapped to the reconstructed finite element mesh. The boundary conditions of the physical deformation analysis model are updated according to the reconstructed finite element mesh, and the physical deformation analysis model is iteratively updated in the COMSOL platform until the preset termination condition is reached. Based on the finally updated physical deformation analysis model, the visualized deformation analysis results under the multiphysics coupling of the target are output.

[0090] This embodiment treats finite element mesh nodes as moving particles. After solving the multiphysics governing equations at each time step, it updates particle positions based on the displacement field and dynamically identifies geometric boundaries based on Delaunay triangulation, achieving a natural description of free surface evolution and large deformation processes. Compared to the traditional finite element method, this embodiment avoids numerical instability caused by mesh distortion and does not rely on complex remeshing techniques or user-defined subroutines, significantly improving computational robustness. Secondly, by integrating the particle finite element framework on the COMSOL platform, the method described in this embodiment can directly call the platform's built-in multiphysics solver and material models, balancing computational accuracy and engineering practicality. Especially in strongly coupled scenarios, this method can effectively simulate large material displacements, free surface motion, and multiphysics coupling behavior, providing a reliable tool for deformation prediction in energy and geotechnical engineering. Furthermore, this embodiment improves numerical efficiency while maintaining computational continuity through mesh reconstruction and state variable mapping strategies, providing a scalable solution for long-term deformation analysis and multi-field coupled simulation in practical engineering.

[0091] To better understand the above technical solutions, exemplary embodiments of the present invention will be described in more detail below with reference to the accompanying drawings. Although exemplary embodiments of the present invention are shown in the drawings, it should be understood that the present invention can be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided so that the present invention can be understood more clearly and thoroughly, and that the scope of the present invention can be fully conveyed to those skilled in the art.

[0092] Specifically, refer to Figure 1 As shown, this embodiment proposes a multiphysics deformation analysis method based on the COMSOL platform and particle finite element method. This method uses the COMSOL multiphysics simulation platform as the main computing environment, and combines it with the MATLAB interface to implement the particle finite element algorithm, forming a computational framework capable of performing large deformation dynamic analysis of soil and rock masses under multiphysics coupling conditions. The method may include the following steps S100 to S600:

[0093] S100. Based on the multiphysics constraints of the object to be analyzed, establish a physical deformation analysis model in COMSOL and generate an initial finite element mesh.

[0094] In this embodiment, step S100 may include the following sub-steps S110 to S130:

[0095] S110. Obtain the multiphysics constraints of the object to be analyzed. The multiphysics constraints include the geometric domain, boundary conditions, material parameters, and multiphysics modules.

[0096] S120. Based on the geometric domain and material parameters, construct the geometric model of the object to be analyzed in the COMSOL platform, and apply physical field constraints to the geometric model based on boundary conditions and multiphysics modules to generate a physical deformation analysis model.

[0097] S130. Perform finite element mesh generation on the physical deformation analysis model to generate an initial finite element mesh, and export the attribute information of the nodes in the computational domain through the real-time data interaction interface. The attribute information includes node coordinates, element connection information and initial state variables.

[0098] This embodiment defines the geometric domain, boundary conditions, material parameters, and multiphysics modules (such as solid mechanics, heat conduction, and fluid flow) of the analysis model in COMSOL, and generates an initial finite element mesh. After the model is defined, the computational domain node coordinates, element connection information, and initial state variables are exported through the COMSOL-MATLAB LiveLink interface.

[0099] S200. Using the solver of the COMSOL platform, solve the multiphysics control equations at each time step to obtain the multiphysics state variables at the current time step.

[0100] In this embodiment, step S200 may include the following sub-steps S210 to S240:

[0101] S210. Based on the constitutive relations and boundary conditions of the physical deformation analysis model, construct the multiphysics field control equations to be solved. The multiphysics field control equations include the momentum balance equation describing mechanical equilibrium, the displacement and strain relationship equation describing deformation motion, and the heat transfer equation describing thermal effects.

[0102] Furthermore, assuming the material follows an elastoplastic constitutive relation, its elastoplastic constitutive model is modeled in the solid mechanics interface: First, the elastic response of the material is described by the generalized Hooke's law, i.e., stress and strain are linearly proportional. Then, to characterize the plastic behavior of the material after reaching the yield limit, the Drucker-Plag criterion and the Mohr-Coulomb criterion, widely used in geotechnical mechanics, are selected from the constitutive models available in the solid mechanics interface to define the yielding and hardening behavior of the material. These two criteria effectively reflect the sensitivity of geotechnical materials to hydrostatic pressure and their shear strength characteristics.

[0103] Furthermore, by combining the aforementioned momentum balance equation, displacement-strain relationship equation, and heat transfer equation with the selected elastoplastic constitutive model and corresponding boundary conditions, a closed, complete thermo-hydraulic-mechanical coupled multiphysics governing equation is finally constructed, which can be used to solve the displacement field, stress field, strain field, and temperature field. The momentum balance equation is as follows:

[0104] (1)

[0105] In equation (1), Here, σ is the gradient operator, b is the stress tensor, ρ is the body force vector, and a is the density.

[0106] The equation relating displacement and strain is as follows.

[0107] (2)

[0108] In equation (2), ε is the strain tensor and u is the displacement.

[0109] The heat transfer equation is:

[0110] (3)

[0111] In equation (3), C p q is the specific heat capacity under constant stress, T is the absolute temperature, q is the conduction heat flux, and Q is the source term.

[0112] S220. Spatial discretization of displacement field variables is performed using quadratic shape functions, and spatial discretization of auxiliary field variables containing stress and temperature gradients is performed using linear shape functions.

[0113] Furthermore, to enhance the comprehensiveness of the numerical calculation, this embodiment selects a hybrid element of quadratic and linear shape functions to co-discrete the displacement, stress, and temperature field variables. The expression for this hybrid space discretization strategy can be uniformly stated as:

[0114] (4)

[0115] In equation (4), N q Let N be a matrix containing quadratic shape functions. 1 Let be a matrix containing linear shape functions, and let “^” represent the discrete values ​​at the nodes.

[0116] S230. While performing spatial discretization, the backward difference formula is used to discretize all field variables in time, so as to transform the multiphysics field control equations into a set of nonlinear equations that can be solved numerically.

[0117] Furthermore, in the COMSOL platform configuration, for time-varying problems, this embodiment selects an implicit time-dependent solver. Specifically, backward difference formula (BDF) is used for time discretization, and the variable order and step size are automatically adjusted. BDF constructs an interpolation polynomial based on the numerical solutions of several previous time steps, and uses this polynomial to approximate the derivative of the current time step. Then, this derivative approximation is used to solve the differential equation. The solution for the field variable x can be approximated by the following formula:

[0118] (5)

[0119] In equation (5), Δt is the time step, n is the index of the current time, s is the BDF order (from order 1 to order 5), and α k The coefficients are determined by the order of the BDF.

[0120] S240. Call the built-in nonlinear solver of the COMSOL platform to solve the nonlinear equation system in each time step until a solution that meets the convergence condition is obtained. The solution is a multi-physics state variable containing displacement field, stress field and temperature field in the current time step.

[0121] In this embodiment, within each time step, the built-in nonlinear solver of COMSOL is used. The nonlinear solver solves the multiphysics control equations using the deployed Newton-Raphson iterative method and BDF time integration algorithm, obtaining the results of multiphysics state variables such as nodal displacement, stress, strain, and temperature. This embodiment constructs a fully coupled thermo-hydraulic-mechanical control equation and employs a hybrid element discretization strategy, effectively controlling the computational scale while ensuring the accuracy of displacement field calculations. Furthermore, the implicit BDF time integration method and the Newton-Raphson iterative solver significantly enhance the stability and efficiency of solving nonlinear coupled problems. This method systematically balances numerical accuracy and computational cost, providing a highly reliable simulation solution for the response analysis of geotechnical materials under multiphysics.

[0122] S300. The nodes of the finite element mesh are used as a set of particles. The spatial position of each particle is updated by combining the multiphysics state variables. Based on the updated set of particles, Delaunay triangulation is performed. The boundary of the obtained triangulated convex domain is identified to obtain the free boundary of the current computational domain.

[0123] In this embodiment, reference Figure 2 As shown in neutron diagram (a), firstly, the nodes of the finite element mesh are treated as a set of particle points, thus treating the computational domain as a point cloud. Then, a triangulated convex domain covering all particles is constructed using Delaunay triangulation; this triangulated convex domain consists of a set of non-overlapping triangular elements. Finally, the Alpha shape algorithm is used to identify the boundaries of the triangulated convex domain, obtaining the free boundaries of the current computational domain. Figure 2 The red lines in 'a' represent free boundaries. This embodiment uses a combination of Delaunay triangulation and the Alpha shape algorithm, which can significantly reduce the number of iterations in the boundary identification process while ensuring the accuracy of the geometric boundaries. Specifically, step S300 may include the following sub-steps S310 to S340:

[0124] S310. Obtain the multiphysics state variables obtained from the current time step from the COMSOL platform, and extract the displacement field data of the nodes.

[0125] S320. Treat each node of the initial finite element mesh as a particle, and update the coordinate position of each particle in space based on the displacement field data to form an updated particle point set.

[0126] S330. Perform Delaunay triangulation on the updated particle point set to generate a triangulated convex domain covering all particles. The triangulated convex domain consists of a set of non-overlapping triangular units.

[0127] S340. Use the Alpha Shape Algorithm to identify the boundaries of the triangulated convex domain and obtain the free boundaries of the current computational domain.

[0128] Further, step S340 may include the following sub-steps S341 to S344:

[0129] S341. Traverse the circumcircle of each triangular unit generated after Delaunay triangulation, and determine the relationship between the radius of each circumcircle and the preset radius threshold.

[0130] S342. Determine triangular elements with a circumscribed circle radius greater than or equal to the radius threshold as boundary elements, and determine triangular elements with a circumscribed circle radius less than the radius threshold as internal elements.

[0131] S343. Filter out the triangle edges of all boundary elements, remove the internal edges that are shared by two boundary elements at the same time, and obtain the boundary edges that belong to the outermost contour of the computational domain.

[0132] S344. Connect the boundary edges to generate a closed and non-convex boundary polygon, which is the free boundary of the current computational domain.

[0133] S400: Reconstruct the finite element mesh of the current time step based on the free boundary, and map the multiphysics state variables of the initial finite element mesh or the finite element mesh of the previous time step onto the reconstructed finite element mesh.

[0134] In this embodiment, step S400 may include the following sub-steps S410 to S440:

[0135] S410. Using the free boundary as a geometric constraint, a distance function-driven mesh generation strategy is adopted to reconstruct the finite element mesh of the current time step within the region of the free boundary.

[0136] Furthermore, during the process of identifying free boundaries using the Alpha Shape Algorithm, a new mesh is generated as a byproduct. However, this mesh often fails to directly meet the mesh quality requirements of finite element analysis, potentially exhibiting issues such as poor element shape, uneven size distribution, or insufficient boundary fitting accuracy. Therefore, this embodiment uses the identified free boundaries as geometric constraints to restart the mesh generator, generating a high-quality mesh suitable for the current time step of the finite element analysis. The distance function-driven mesh generation strategy in the mesh generator effectively controls the size distribution of mesh elements and the boundary fitting accuracy, thereby generating a high-precision, high-quality mesh near the free boundaries while ensuring a smooth transition of the mesh in the internal regions, improving element shape quality and overall topological continuity.

[0137] S420. By comparing the geometric relationships of the finite element mesh before and after reconstruction, a mapping relationship is established between the nodes of the old and new meshes.

[0138] S430. Based on the mapping relationship, a unique element mapping strategy is adopted to map the multiphysics state variables of the initial finite element mesh or the finite element mesh of the previous time step to the finite element mesh reconstructed in the current time step.

[0139] Furthermore, this embodiment employs the unique element method for variable mapping. This method exhibits high robustness and accuracy in analyzing various geotechnical engineering problems, and is particularly suitable for complex conditions involving large deformations, material nonlinearity, and mesh reconstruction. The unique element method is based on the finite element shape function and achieves the transfer of state variables by establishing logical connections between old and new mesh nodes. Specifically, for each new Gaussian point (e.g., ... Figure 2 As shown by the "red cross" marking in sub-figure b, the corresponding natural coordinates are first determined based on its spatial position relative to the nearest old element node. Then, interpolation or extrapolation calculations are performed using the state variable values ​​stored at the original Gaussian point ("black cross" marking) in the old element to complete the high-precision mapping of physical field variables. The element-unique mapping strategy effectively maintains the continuity and balance of field variables at the element level, ensuring the continuity and conservation of stress, strain, temperature, and other field variables.

[0140] S440. After completing the mapping of multiphysics state variables, the updated multiphysics state variables are used as the initial conditions for solving the multiphysics control equations in the next time step, and loaded into the physical deformation analysis model of the COMSOL platform.

[0141] S500: Update the boundary conditions of the physical deformation analysis model based on the reconstructed finite element mesh, and iteratively update the physical deformation analysis model in the COMSOL platform until the preset termination condition is reached.

[0142] In this embodiment, step S500 may include the following sub-steps S510 to S540:

[0143] S510. Based on the node coordinates of the reconstructed finite element mesh, obtain the free boundary position of the current computational domain.

[0144] Furthermore, after the finite element mesh is reconstructed, the topology and geometry of its computational domain have changed. Therefore, it is first necessary to identify the mesh boundaries that constitute the current computational domain. This embodiment achieves this by comparing the internal element surfaces and the external contour surfaces, ensuring that subsequent boundary conditions can be accurately applied to the deformed actual physical boundaries.

[0145] S520. Reapply the boundary conditions defined in the original physical deformation analysis model based on the current geometric position of the free boundary. The boundary conditions include mechanical boundary conditions and thermal boundary conditions.

[0146] S530. In the COMSOL platform, the physical deformation analysis model that has been re-applied with boundary conditions and completed state variable mapping is updated to the solution model for the current time step.

[0147] S540. Determine whether the current iteration state has reached or exceeded the preset termination condition.

[0148] S550a If the preset termination condition is not met, the current time step is taken as the new initial state, the time step size is incremented, and the process returns to the step of solving the multiphysics control equations in each time step using the solver of the COMSOL platform, and then iterative calculations are performed for the next time step.

[0149] S550b: If the preset termination condition has been met or exceeded, the iteration process is terminated, and the final updated physical deformation analysis model is output.

[0150] S600, based on the final updated physical deformation analysis model, outputs visualized deformation analysis results under the multi-physics coupling effect of the target.

[0151] In this embodiment, step S600 may include the following sub-steps S610 to S640:

[0152] S610. Obtain the latest physical deformation analysis model from the COMSOL platform.

[0153] S620. Post-process the physical deformation analysis model to generate a result file that can intuitively reflect the deformation process and final state of the object under analysis under the coupling of multiple physics fields.

[0154] Furthermore, the post-processing operations include: extracting data on key physical quantities, such as the displacement field, stress field, strain field, and distribution data of related field variables; generating a series of result sequences of time steps or parameter scan steps based on the distribution data to dynamically display the deformation process; and finally, packaging the result data into a result file for output, wherein the result file may include a geometric deformation animation cloud map, the distribution curves of key physical quantities along the path, and the historical curves of the physical quantities at specific points changing over time.

[0155] S630. Render and display the result file using preset visualization tools, and output a visualization chart of the deformation analysis results.

[0156] Furthermore, this embodiment also proposes a multiphysics deformation analysis system based on the COMSOL platform and particle finite element method, which includes:

[0157] The model initialization and mesh generation module is used to establish a physical deformation analysis model in COMSOL and generate an initial finite element mesh based on the multiphysics constraints of the object to be analyzed.

[0158] The multiphysics calculation module is used to solve the multiphysics control equations at each time step using the solver of the COMSOL platform, and obtain the multiphysics state variables at the current time step.

[0159] The geometric boundary identification module is used to treat the nodes of the finite element mesh as a set of particles, update the spatial position of each particle in combination with multiphysics state variables, and perform Delaunay triangulation based on the updated set of particles. The module then identifies the boundaries of the obtained triangulated convex domain to obtain the free boundary of the current computational domain.

[0160] The Mesh Reconstruction and Variable Mapping module is used to reconstruct the finite element mesh of the current time step based on the free boundary, and to map the multiphysics state variables of the initial finite element mesh or the finite element mesh of the previous time step onto the reconstructed finite element mesh.

[0161] The boundary condition update and model iteration module is used to update the boundary conditions of the physical deformation analysis model based on the reconstructed finite element mesh, and iteratively update the physical deformation analysis model in the COMSOL platform until the preset termination condition is reached.

[0162] The results visualization output module is used to output visualized deformation analysis results of the target under multi-physics coupling based on the final updated physical deformation analysis model.

[0163] Based on this, this embodiment uses the following five specific examples to demonstrate the beneficial effects of the proposed multiphysics deformation analysis method based on the COMSOL platform and particle finite element method. Each example verifies the advantages of this method in terms of computational accuracy, convergence efficiency, and multi-field coupling analysis capability, starting from different physical field coupling conditions, different material models, and different boundary conditions:

[0164] Example 1: Comparative Verification Analysis of Beam Bending

[0165] This embodiment uses a single COMSOL platform and the method proposed in this embodiment to perform two independent beam bending analyses, such as... Figure 3 As shown, the bending behavior of an elastic cantilever beam, 20m long and 1m wide, under a bending moment M at its right end is analyzed. The beam's material properties are elastic modulus E = 5 GPa and Poisson's ratio v = 0. In the simulation, the left end of the beam is fixed, while the top and bottom surfaces are assumed to be free. The model is discretized into 2036 elements and 1124 mesh nodes. Theoretical analysis shows that the beam can bend into a circle when the bending moment reaches M = 2πEI / l, where E is the elastic modulus, I is the area of ​​inertia, and l is the beam length. In both analyses, the applied bending moment gradually increases from zero in increments of 0.01 Ma.

[0166] refer to Figure 4As shown, it displays the results obtained from a single COMSOL platform at four time points, with bending moments of 0.25 Ma, 0.5 Ma, 0.75 Ma, and 1.0 Ma. Because large deformations are involved as the bending moment increases, a single COMSOL platform cannot capture the accurate evolution of the beam, specifically because similar final profiles are obtained when subjected to different magnitudes of bending moments.

[0167] In contrast, the method proposed in this embodiment is able to capture, such as Figure 5 The evolution of a beam with extreme deformation is shown. When the applied bending moment increases to 1.04 Ma, a near-perfect circle is obtained, with the deviation from the analytical solution being only about 4%, further confirming the accuracy of the multiphysics deformation analysis method based on the COMSOL platform and particle finite element method in the elastic analysis model.

[0168] Example 2: Bearing Capacity Verification Analysis of Strip Foundations

[0169] This embodiment performs a bearing capacity verification analysis on a strip foundation on a uniform, weightless soil surface. For example... Figure 6 As shown, the following material parameters are used: Young's modulus E = 100 MPa, Poisson's ratio v = 0.49, and undrained shear strength S. u =100kPa, and internal friction angle φ=0°. The width of the rigid strip foundation is set to B=1m, and a rolling boundary condition is applied on the left side to account for symmetry. The entire computational domain is selected as a square with a side length of 10B. The boundary conditions are set as follows: Figure 6 Fixed boundary conditions are applied to the right and bottom of the region shown in Figure a, and the load is applied by defining incremental vertical displacements of the nodes within the loading region. The incremental vertical displacement is set to a constant ∆u. y = -0.005m, the simulation ends when the total vertical displacement reaches 1m. The mesh generation scheme is as follows: Figure 6 As shown in b, the entire computational domain is divided into three sub-regions (I, II, III), with corresponding grid side lengths of 0.025m, 0.15m, and 0.625m, respectively.

[0170] The curves of normalized vertical reaction force versus normalized vertical displacement are as follows: Figure 7 As shown. The results of COMSOL-PFEM (the computational framework of the method described in this embodiment) fall within the range of Prandtl's analytical solution, specifically 5.14 for small deformation and 8.28 for deep penetration, consistent with the theoretical expectation of Meyerhof's ultimate load. Furthermore, it can be seen that the normalized vertical reaction increases with increasing penetration depth, consistent with the foundation bearing capacity theory. (Reference) Figure 7As shown, when the results of COMSOL-PFEM are compared with other independent simulation results, it can be seen that the curve of COMSOL-PFEM is closer to the curve calculated by using rigid-plastic analysis and smooth particle finite element analysis, and the arbitrary Lagrange and Euler method (ALE) failed to provide sufficiently accurate results in the bearing capacity verification of strip foundations.

[0171] in u y The corresponding distribution of cumulative plastic strain invariants and shear stress at point -1m is as follows: Figure 8 As shown. Plastic strain distribution ( Figure 8 a): A wedge-shaped region with predominantly elastic behavior forms near the foundation. A distinct transition zone exists below and adjacent to the loading area, within which significant plastic strain concentration is observed. Stress distribution characteristics ( Figure 8 b): The simulation results show a smooth and continuous stress field distribution with no obvious stress oscillation phenomenon, which verifies the correctness and stability of the elastoplastic analysis framework implemented in this paper.

[0172] Example 3: Verification and Analysis of Cohesive Soil Collapse

[0173] This embodiment analyzes the collapse problem of a rectangular region composed of cohesive soil, specifically simulating particle flow and landslide issues. The setup is as follows: Figure 9 As shown in figure a, its model parameters are: Young's modulus E = 1.8 MPa, Poisson's ratio ν = 0.2, density ρ = 1850 kg / m³ 2 The friction angle φ = 25° and the cohesion c = 5 kPa. The bottom and left boundaries of the soil are constrained by fixed boundary conditions, while the top and right boundaries are assumed to be free boundary conditions. Figure 9 b shows the final sediment profile output by the method described in this embodiment, and compares it with the reference simulation results (SPH results, a simulation based on smooth particle hydrodynamics). The two result curves show a high degree of agreement. Furthermore, referring to... Figure 10 As shown, the simulation results of this embodiment show good consistency with the equivalent plastic strain fields of the reference simulation results at four time points. Specifically, in the initial deformation stage (t=0.6s), the equivalent plastic strain is mainly concentrated in the slope toe region, and the distribution pattern is basically consistent with the reference simulation. As the deformation develops (t=0.88s), the plastic zone gradually expands upward, and the results of this embodiment accurately capture the formation process of the sliding surface. In the sliding intensification stage (t=1.28s), the plastic strain further extends to the rear edge of the slope, which highly coincides with the strain concentration area in the reference simulation results. Finally, in the stable stage (t=2.0s), the plastic strain distribution tends to stabilize, and the strain field shape, magnitude, and spatial distribution simulated in this embodiment are consistent with the reference results.

[0174] The above comparison results show that the method proposed in this embodiment can effectively simulate the particle flow behavior, final accumulation morphology and plastic strain development process of cohesive soil during the collapse process, verifying its accuracy and applicability in landslide simulation.

[0175] Example 4: Verification and Analysis of Soil Slope Failure under Thermal Load

[0176] This embodiment demonstrates the applicability of the COMSOL-PFEM framework, a multiphysics deformation analysis method based on the COMSOL platform and particle finite element method. It uses the COMSOL-PFEM framework to perform multiphysics deformation analysis on energy geotechnical engineering, specifically the dynamic response of soil around energy piles and methane hydrate-containing sediments under temperature variations. This embodiment first verifies the correctness of the thermal module, and then further simulates the progressive slope failure under thermal loading conditions.

[0177] Thermal effect of a two-dimensional plate:

[0178] First, construct a 1×1 meter rectangular board, such as... Figure 11 As shown in figure a, the initial temperature is set to T0 = 0℃. Assume the rectangular plate is an isotropic material and follows a linear thermoelastic relationship. The mechanical parameters are selected as follows: material density ρ = 1 kg / m³. 3 Young's modulus E = 1.0 Pa, Poisson's ratio ν = 0.25, specific heat capacity ν = 0.25, thermal conductivity k = 0.1 W / (m∙K), linear thermal expansion coefficient β = 1 × 10⁻⁶ -6 / ℃. Roller boundary conditions are set along the left and lower boundaries, while the upper and right boundaries are assumed to be free boundaries. At the start of the simulation, the temperature of the left and lower boundaries is set to T1 = 1℃. The locations of three monitoring points A, B, and C are marked. The simulation uses a fixed time step ∆t = 0.0001s and ends at t = 0.5s.

[0179] refer to Figure 11 As shown in Figure b, the temperature changes at three monitoring points A, B, and C within the COMSOL-PFEM framework were compared with the temperature changes calculated using a single COMSOL platform, and these simulation results showed a high degree of similarity. Furthermore, in Figure 11 c and Figure 11 The temperature distribution at t=0.5s was compared in section d, and the two results matched very well. Therefore, the correctness of the thermal module was verified.

[0180] Thermally induced slope instability:

[0181] refer to Figure 12 The slope model shown in figure a, where a 6.0-meter-high slope is discretized into 4223 elements, has boundary conditions as follows. Figure 12As shown in b: the top and slope surfaces are designated as free boundaries, the left boundary as a roller boundary, and the bottom as a fixed boundary; a constant temperature of 50°C is set along the left boundary of the slope, and the remaining parts are set as adiabatic. To account for soil weakening behavior caused by temperature changes and plastic strain accumulation, this embodiment uses the following relationship to update the relationship between temperature T and plastic deviation strain. Influenced cohesion c and internal friction angle φ:

[0182] (6)

[0183] (7)

[0184] In equations (6) and (7), c0 and φ0 represent the cohesive force and friction angle at the reference temperature T0, respectively. res and φ res These are the residual values ​​of c and φ, respectively. η and η T This is a material constant used to control the effects of strain softening and thermal softening.

[0185] The model parameters used are: Young's modulus E = 20 MPa, Poisson's ratio ν = 0.3, and density ρ = 1850 kg / m³. 3 The initial soil cohesion c0 = 30 kPa, and the residual soil cohesion c res =5kPa, initial friction angle φ0=25°, residual friction angle φ res =15°, expansion angle ψ=0°, specific heat capacity of soil C s =2080J / (kgK), coefficient of thermal expansion C T =10 -6 The initial temperature T0 = 0℃, η = 100 and η T =1. To speed up the simulation, the specific heat capacity of the soil has been amplified by a factor of 1.728 × 10⁵ days / second, making 1 second in the simulation equivalent to 2 days in reality.

[0186] To better illustrate the evolutionary behavior of slopes, thermo-induced degradation was demonstrated. Figure 12 ) and post-instability evolution ( Figure 13 Numerical results are presented to highlight the formation of the overall slip surface and the transient slip evolution process after failure.

[0187] Thermotropic degradation:

[0188] refer to Figure 13 Figure a shows the temperature distribution within the slope at five representative time points (t=1, 2, 3, 4, and 5 days), clearly depicting the process of heat propagation from the left boundary inwards. The thickness of the perturbed region (i.e., T>0°C) gradually increases, and a significant geometric change is observed at t=5 days. Reference Figure 13As shown in b, it illustrates the formation process of the overall sliding surface, in which plastic yielding begins first at the left boundary (see t<3 days) and then extends further to the upper and bottom surfaces as temperature spreads (t=4 and 5 days).

[0189] Post-instability evolution:

[0190] After the formation of the global sliding surface, the stable slope evolves into an unstable state, which is related to the acceleration of the sliding body and the significant evolution of the free surface. Considering the difference in time scale between the pre-failure and post-failure stages, this embodiment defines t as follows. post To indicate, for example Figure 14 The time stage shown is after t=5 days. To illustrate the post-failure evolution, the velocity, equivalent plastic strain, and cohesion distributions at four time points were plotted. From the velocity distribution ( Figure 14 a) It can be seen that the formed sliding body consists of two main blocks, with a maximum velocity greater than 0.5 m / s, and at t post The process stopped abruptly at 3.7 s. The equivalent plastic strain is as follows: Figure 14 As shown in b, it illustrates the block differences separated by two sliding surfaces (i.e., the two red curves). The relevant distribution of cohesion is as follows: Figure 14 As shown in c, mechanical degradation is mainly concentrated in the regions with high temperature and high equivalent plastic strain.

[0191] Example 5: Verification and Analysis of Foundation Settlement Caused by Heat Transfer from Buried Pipelines in Clay

[0192] This embodiment aims to further verify the application of the COMSOL-PFEM framework in energy geotechnical engineering by conducting numerical simulations of foundation settlement caused by soil strength weakening under pipeline heat transfer conditions in clay.

[0193] refer to Figure 15 As shown in figure a, it illustrates the geometry of a model with a discrete mesh, which contains 3305 triangular elements. Figure 15 As shown in b, the model contains two objects: the foundation and the soil (the pipe is simplified to a fixed boundary condition with a radius R = 0.5 m). Based on whether the soil strength is affected by temperature, the soil is divided into an active layer (2 m) and a passive layer (4 m). The boundary conditions are set as follows: the top of the soil is set as a free boundary, the left and right boundaries are set as rolling boundaries, and the bottom is set as a fixed boundary; a pressure of 40 kPa is applied to the top of the foundation, and the sides of the foundation are set as free boundaries; a constant temperature of 50°C is set along the pipe boundary, and the other parts are set as adiabatic. In this model, the foundation is considered a rigid body. Both the active and passive layers are clay, and their mechanical behavior is described by the Drucker-Prag yield criterion, while the mechanical parameters of the active layer are updated by equations (6) and (7).

[0194] As computation time increases, reference Figure 16As shown in Figure a, it illustrates the evolution of the soil temperature field in the active layer under the influence of pipeline heat transfer. The temperature gradually increases from the middle section above the pipeline. At t=40 days, the temperature of the soil surrounding the foundation is also affected and rises, resulting in a significant tilting of the foundation. (Reference) Figure 16 As shown in b, it illustrates the evolution of soil cohesion within the active layer. The weakening trend of soil strength is closely consistent with the rising trend of temperature, which directly reveals the cause of foundation settlement.

[0195] In summary, this invention proposes a multiphysics deformation analysis method and system based on the COMSOL platform and particle finite element method. By coupling the particle finite element method (PFEM) with the COMSOL multiphysics platform, stable and high-precision numerical analysis of large deformation and multi-field coupling processes in energy geotechnical engineering is achieved. Compared with the traditional finite element method, this invention has significant effects and advantages in the following aspects:

[0196] 1. Improved computational stability: In cases of large geometric deformation, the traditional COMSOL finite element model suffers from mesh distortion, leading to solution failure. However, the method of this invention reconstructs the mesh through Delaunay triangulation and the Alpha shape algorithm, ensuring computational stability under extreme deformation and enabling stable convergence throughout the entire computation process.

[0197] 2. Improved calculation accuracy: In the case of bending of elastic cantilever beam, the maximum bending radius error obtained by this invention is less than 4%, which is highly consistent with the analytical solution; in the case of foundation bearing capacity, the result is between Prandtl's theoretical solution and Meyerhoff's ultimate load limit solution, with an error of less than 5%.

[0198] 3. Enhanced multi-field coupling performance: By integrating with COMSOL's native thermo-hydraulic-mechanical module, this invention can continuously simulate the interaction between temperature field, seepage field and stress field, accurately capturing the entire process of thermally induced landslides and pipeline foundation settlement.

[0199] 4. Improved efficiency and ease of operation: Compared with traditional self-written particle finite element programs, this invention is based on the COMSOL platform and can directly call its graphical modeling and post-processing functions.

[0200] 5. High scalability: This invention can expand the large deformation calculation function without modifying the COMSOL core solver, and is suitable for various engineering scenarios such as energy pile heat exchange, geothermal extraction, and carbon sequestration. It has good versatility and engineering application prospects.

[0201] Since the systems / devices described in the above embodiments of the present invention are systems / devices used to implement the methods of the above embodiments of the present invention, those skilled in the art can understand the specific structure and modifications of the systems / devices based on the methods described in the above embodiments of the present invention, and therefore will not be repeated here. All systems / devices used in the methods of the above embodiments of the present invention fall within the scope of protection of the present invention.

[0202] Those skilled in the art will understand that embodiments of the present invention can be provided as methods, systems, or computer program products. Therefore, the present invention can take the form of a completely hardware embodiment, a completely software embodiment, or an embodiment combining software and hardware aspects. Furthermore, the present invention can take the form of a computer program product embodied on one or more computer-usable storage media (including, but not limited to, disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0203] This invention is described with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the invention. It will be understood that each block of the flowchart illustrations and / or block diagrams, as well as combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions.

[0204] It should be noted that in the description of this invention, the word "a" or "an" preceding a component does not exclude the existence of multiple such components. This invention can be implemented by means of hardware comprising several different components and by means of a suitably programmed computer. The use of terms such as first, second, third, etc., is merely for convenience and does not indicate any order. These terms can be understood as part of the component names.

[0205] Furthermore, it should be noted that in the description of this specification, the terms "one embodiment," "some embodiments," "embodiment," "example," "specific example," or "some examples," etc., refer to specific features, structures, materials, or characteristics described in connection with that embodiment or example, which are included in at least one embodiment or example of the present invention. In this specification, the illustrative expressions of the above terms do not necessarily refer to the same embodiment or example. Moreover, the specific features, structures, materials, or characteristics described may be combined in any suitable manner in one or more embodiments or examples. Furthermore, without contradiction, those skilled in the art can combine and integrate the different embodiments or examples described in this specification, as well as the features of different embodiments or examples.

[0206] Although preferred embodiments of the invention have been described, those skilled in the art, upon learning of the basic inventive concept, can make other changes and modifications to these embodiments.

[0207] Obviously, those skilled in the art can make various modifications and variations to this invention without departing from the spirit and scope of the invention.

Claims

1. A COMSOL platform and particle finite element based multi-physics deformation analysis method, characterized in that, The method comprises the following steps: establishing a physical deformation analysis model in COMSOL according to the multi-physical field constraint conditions of an object to be analyzed, and generating an initial finite element grid; solving the multi-physical field control equations at each time step by using a solver of the COMSOL platform to obtain the multi-physical field state variables at the current time step; updating the spatial positions of the particles by using the nodes of the finite element grid as a particle point set and combining the multi-physical field state variables, and performing Delaunay triangulation on the updated particle point set to obtain a triangulated convex domain, and identifying the boundary of the triangulated convex domain to obtain the free boundary of the current calculation domain; reconstructing the finite element grid at the current time step based on the free boundary, and mapping the multi-physical field state variables of the initial finite element grid or the finite element grid at the previous time step to the reconstructed finite element grid; updating the boundary conditions of the physical deformation analysis model according to the reconstructed finite element grid, and iteratively updating the physical deformation analysis model in the COMSOL platform until a preset termination condition is reached; outputting a visual deformation analysis result under the coupling action of the target multi-physical field based on the finally updated physical deformation analysis model.

2. The method of claim 1, wherein, The method comprises the following steps: obtaining the multi-physical field constraint conditions of the object to be analyzed, wherein the multi-physical field constraint conditions comprise a geometric domain, boundary conditions, material parameters, and a multi-physical field module; constructing a geometric model of the object to be analyzed in the COMSOL platform according to the geometric domain and the material parameters, and applying physical field constraints to the geometric model based on the boundary conditions and the multi-physical field module to generate a physical deformation analysis model; dividing the physical deformation analysis model into finite element grids to generate an initial finite element grid, and exporting attribute information of the nodes in the calculation domain through a real-time data interaction interface, wherein the attribute information comprises node coordinates, element connection information, and initial state variables.

3. The method of claim 1, wherein, The method comprises the following steps: constructing the multi-physical field control equations to be solved according to the constitutive relation of the physical deformation analysis model and the boundary conditions, wherein the multi-physical field control equations comprise a momentum balance equation describing mechanical equilibrium, a displacement and strain relationship equation describing deformation motion, and a heat transfer equation describing thermal effects; spatially discretizing the displacement field variables by using a quadratic shape function, and spatially discretizing auxiliary field variables including stress and temperature gradient by using a linear shape function; spatially discretizing and time-discretizing all field variables by using a backward difference formula to convert the multi-physical field control equations into a nonlinear equation set that can be used for numerical solution; calling a nonlinear solver built in the COMSOL platform to solve the nonlinear equation set at each time step until a solution satisfying a convergence condition is obtained, and the solution is the multi-physical field state variables including displacement field, stress field, and temperature field at the current time step.

4. The method of claim 2, wherein the momentum balance equation is: the displacement and strain relationship equation is: ； wherein is the gradient operator, σ is the stress tensor, b is the body force vector, p is the density, and a is the acceleration; wherein ε is a strain tensor and u is a displacement; ； the heat transfer equation is: ​ ； In the formula, C p is the specific heat capacity under constant stress, T is the absolute temperature, q is the heat flux of conduction, and Q is the source term.

5. The method of claim 1, wherein, The nodes of the finite element mesh are taken as a particle point set, the spatial positions of each particle are updated in combination with multi-physical field state variables, and based on the updated particle point set, Delaunay triangulation is performed, the boundary of the obtained triangulated convex domain is identified, and the free boundary of the current calculation domain is obtained, including: The multi-physical field state variables solved at the current time step are obtained from the COMSOL platform, and the displacement field data of the nodes are extracted; Each node of the initial finite element mesh is regarded as a particle, and the coordinate position of each particle in space is updated based on the displacement field data to form an updated particle point set; Delaunay triangulation is performed on the updated particle point set to generate a triangulated convex domain covering all particles, and the triangulated convex domain is composed of a group of non-overlapping triangular elements; The alpha shape algorithm is used to identify the boundary of the triangulated convex domain to obtain the free boundary of the current calculation domain.

6. The method of claim 5, wherein, The alpha shape algorithm is used to identify the boundary of the triangulated convex domain to obtain the free boundary of the current calculation domain, including: The circumscribed circle of each triangular element generated after Delaunay triangulation is traversed, and the size relationship between the radius of each circumscribed circle and the preset radius threshold is judged; Triangular elements with a circumscribed circle radius greater than or equal to the radius threshold are determined as boundary elements, and triangular elements with a circumscribed circle radius less than the radius threshold are determined as internal elements; The triangular edges of all boundary elements are filtered to remove internal edges shared by two boundary elements to obtain boundary edges belonging to the outermost contour of the calculation domain; The boundary edges are connected to generate a closed and non-convex boundary polygon, which is the free boundary of the current calculation domain.

7. The method of claim 1, wherein, The finite element mesh of the current time step is reconstructed based on the free boundary, and the multi-physical field state variables of the initial finite element mesh or the finite element mesh of the previous time step are mapped to the reconstructed finite element mesh, including: The free boundary is taken as a geometric constraint, and a distance function driven mesh generation strategy is used to reconstruct the finite element mesh of the current time step in the area of the free boundary; By comparing the geometric relationship of the reconstructed finite element mesh with the previous finite element mesh, a mapping relationship between the nodes of the new and old meshes is established; Based on the mapping relationship, an element unique mapping strategy is used to map the multi-physical field state variables of the initial finite element mesh or the finite element mesh of the previous time step to the finite element mesh reconstructed at the current time step; After completing the mapping of the multi-physical field state variables, the updated multi-physical field state variables are taken as the initial conditions for solving the multi-physical field control equations at the next time step, and are loaded into the physical deformation analysis model in the COMSOL platform.

8. The method of claim 1, wherein, The boundary conditions of the physical deformation analysis model are updated according to the reconstructed finite element mesh, and the physical deformation analysis model is iteratively updated in the COMSOL platform until a preset termination condition is reached, including: The free boundary position of the current calculation domain is obtained according to the node coordinates of the reconstructed finite element mesh; The boundary conditions defined in the original physical deformation analysis model are re-applied according to the geometric position of the current free boundary, and the boundary conditions include mechanical boundary conditions and thermal boundary conditions; In the COMSOL platform, the physical deformation analysis model with the re-imposed boundary conditions and the completed state variable mapping is updated as the solution model of the current time step; It is judged whether the current iteration state reaches or exceeds the preset termination condition; If the preset termination condition is not reached, the current time step is taken as a new initial state, the time step length is incremented, and the step of solving the multi-physics control equation in each time step by using the solver of the COMSOL platform is returned to perform the iteration calculation of the next time step; If the preset termination condition is reached or exceeded, the iteration process is terminated, and the final updated physical deformation analysis model is output.

9. The method of claim 1, wherein, Based on the final updated physical deformation analysis model, the visualized deformation analysis result under the target multi-physics coupling effect is output, including: The final updated physical deformation analysis model is obtained from the COMSOL platform; The physical deformation analysis model is post-processed to generate a result file capable of directly reflecting the deformation process and final state of the object to be analyzed under the multi-physics coupling effect; The result file is rendered and displayed by using the preset visualization tool, and the visualized chart of the deformation analysis result is output.

10. A COMSOL platform and particle finite element based multi-physics deformation analysis system, characterized in that, It includes: A model initialization and mesh generation module for establishing a physical deformation analysis model in COMSOL according to the multi-physics constraint conditions of the object to be analyzed, and generating an initial finite element mesh; A multi-physics calculation module for solving the multi-physics control equation in each time step by using the solver of the COMSOL platform to obtain the multi-physics state variable of the current time step; A geometric boundary identification module for taking the nodes of the finite element mesh as a particle point set, updating the spatial position of each particle in combination with the multi-physics state variable, performing Delaunay triangulation based on the updated particle point set, identifying the boundary of the obtained triangulation convex domain, and obtaining the free boundary of the current calculation domain; A mesh reconstruction and variable mapping module for reconstructing the finite element mesh of the current time step based on the free boundary, and mapping the multi-physics state variable of the initial finite element mesh or the finite element mesh of the previous time step to the reconstructed finite element mesh; A boundary condition updating and model iteration module for updating the boundary conditions of the physical deformation analysis model according to the reconstructed finite element mesh, and iteratively updating the physical deformation analysis model in the COMSOL platform until the preset termination condition is reached; A result visualization output module for outputting the visualized deformation analysis result under the target multi-physics coupling effect based on the final updated physical deformation analysis model.

Images