A meshless analysis method considering background mesh partition refinement and truncation optimization

Abstract

A meshless analysis method considering background mesh partition refinement and truncation optimization belongs to the technical field of numerical simulation of meshless method. The steps are as follows: S1, reading in the boundary node information, internal node information, material information and load information of the meshless model in sequence; S2, positioning the background mesh near the boundary of the meshless model, and performing local quartering refinement on the background mesh; S3, generating a plurality of Gaussian points in the background mesh after local quartering refinement, and calculating the integral area of each Gaussian point for subsequent numerical integration; S4, calculating the shape function at each Gaussian point, and discretizing and solving the meshless balance equation based on the integral area of each Gaussian point in S3, to complete the meshless analysis considering background mesh partition refinement and truncation optimization. Through local quartering refinement and truncation optimization of the background mesh, the meshless method can improve the representation ability of the meshless method for complex boundaries, improve the problem of the decline of the simulation accuracy of the meshless method at the boundary, and expand the application range of the meshless method.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of meshless numerical simulation, and relates to a meshless analysis method considering background mesh partition refinement and truncation optimization. BACKGROUND

[0002] The traditional finite element method (FEM) has been the mainstream method in the field of engineering simulation for a long time due to its mature mathematical theory and calculation framework. Compared with the traditional finite element method, the meshless method (MFM) has obvious advantages in dealing with large deformation, dynamic fracture and multi-field coupling problems because it does not need to be bound by elements in the process of interpolation function calculation or numerical integration. The concept of meshless method is broad, and it is classified differently according to different criteria. The most recognized classification is to divide the meshless method into meshless method based on particle method, meshless method based on global weak form and meshless method based on local weak form according to the different ways of discrete equations.

[0003] The MFM based on particle method adopts node integration, and the representative methods include smoothed particle hydrodynamics, generalized finite difference method and finite point method. The meshless method based on global weak form performs Gaussian integration by means of background mesh covering the entire solution domain, and the representative methods include element Galerkin method, radial basis interpolation method and reproducing kernel particle method. The meshless method based on local weak form generates independent background mesh near each node for Gaussian integration, and the representative methods include meshless local Petrov-Galerkin method and local point interpolation method. Among them, the meshless method based on global weak form is widely used in engineering due to its clear physical meaning, high stability and easy implementation, including problems such as rainfall and self-weight induced slope sliding, strong earthquake induced ground liquefaction and crack development under combined load.

[0004] The numerical integration process of the global weak form meshless method depends on the background mesh, and the background mesh needs to cover the area before and after the model deformation, so the numerical integration process cannot truly reflect the geometric boundary and material partition boundary of the model, and under complex boundary conditions (the boundary does not coincide with the background mesh), the precision will decrease, which undoubtedly seriously limits the application of the meshless method in practical engineering and large deformation analysis. For example, Chinese invention patent CN201910212828.6 discloses a GPU acceleration method for thin plate meshless Galerkin structure dynamic response analysis, which does not consider improving the simulation precision at the boundary of the meshless method, is suitable for simple geometric model problems such as thin plates, and is difficult to apply to practical engineering; Chinese invention patent CN202210028137.2 discloses a statics calculation method based on a new type of meshless stable distribution point, which combines the distribution point method and the boundary stiffness array technology to improve the simulation precision of the meshless boundary, but it is only suitable for linear elastic constitutive model for mechanical problems. The boundary to be solved in practical engineering is generally complex, and in the process of large deformation analysis, the boundary is also dynamic, and a complex soil constitutive model also needs to be introduced.

[0005] Therefore, it is necessary to propose a background mesh partition refinement and truncation optimization method to improve the simulation precision of MFM at the boundary in practical engineering. SUMMARY

[0006] In view of the defects of the meshless method in the precision decrease at the complex boundary, the present application provides a meshless analysis method considering background mesh partition refinement and truncation optimization. The present application adopts the strategy of local four-part refinement of the background mesh at the geometric boundary and material boundary and integration domain truncation optimization, improves the capture ability of MFM for dynamic boundary, and further improves the simulation precision, which has important significance for expanding its application in practical engineering.

[0007] To achieve the above purpose, the technical scheme adopted by the present application is as follows:

[0008] A meshless analysis method considering background mesh partition refinement and truncation optimization, comprising the following steps:

[0009] S1, read in the boundary node information, internal node information, material information and load information of the meshless model in sequence, specifically:

[0010] S1.1, read in the boundary node information and internal node information of the meshless model, wherein the boundary nodes include geometric boundary nodes and material boundary nodes;

[0011] S1.2, based on the boundary node information and the internal node information in S1.1, obtaining the average distance of the meshless nodes, calculating the size of the rectangular background mesh, the size of the background mesh in the present application is twice the average distance of the meshless nodes, the background mesh here is defined as the initial background mesh;

[0012] S1.3, then reading in the material information of the meshless model, including material number, material partition, material parameters;

[0013] S1.4, then reading in the load information of the meshless model, including body force, surface force, node force, displacement condition;

[0014] S2, positioning the background mesh near the boundary of the meshless model, and performing local quadtree refinement thereon, in particular:

[0015] S2.1, sequentially connecting the boundary nodes in S1 to form a polyline representing the model boundary;

[0016] S2.2, performing a cross test on each initial background mesh in S1, sequentially determining whether the polyline in S2.1 representing the model boundary intersects with the initial background mesh, when both formula (1) and formula (2) are satisfied, it is considered that the two intersect:

[0017] (1)

[0018] (2)

[0019] In the formula, (x b1 ,y b1 ) is the coordinates of a meshless node in the model boundary in S2.1; (x b2 ,y b2 ) is the coordinates of the next meshless node in the clockwise direction in the model boundary in S2.1; (X i ,Y i ) and (X i+1 ,Y i+1 ) represent the position coordinates of two nodes on each side of the rectangular background mesh;

[0020] S2.3, for the initial background mesh intersecting with the model boundary in S2.2, performing n times of local quadtree refinement according to the refinement number n, the size of the background mesh after local quadtree refinement is shown in formula (3); for the initial background mesh not intersecting with the model boundary in S2.2, the size remains unchanged.

[0021] (3)

[0022] In the formula, n represents the refinement number, represents the size of the initial background mesh in S1, represents the size of the refined background mesh in S2.

[0023] S3, generate a plurality of Gaussian points in the background mesh after local quarter subdivision in S2, and calculate the integral area of each Gaussian point for subsequent numerical integration, specifically:

[0024] S3.1, generate a plurality of Gaussian points in the background mesh after local quarter subdivision in S2, and calculate the integral area of each Gaussian point for subsequent numerical integration, specifically:

[0025] S3.2, then use formula (1) and formula (2) to determine whether the background mesh after local quarter subdivision intersects with the model boundary;

[0026] S3.3, for the background mesh that does not intersect with the model boundary, the integral area of the Gaussian point inside is the area of the background mesh; for the background mesh that intersects with the model boundary, the model boundary can truncate the background mesh into a polygon, then the integral area of the Gaussian point inside needs to be calculated according to formula (4) and formula (5);

[0027] (4)

[0028] (5)

[0029] In the formula, (xx i ,yy i ) represents the node coordinates of the truncated background mesh, i represents the i-th node in the truncated background mesh, =1, 2, 3, …, k; k represents the number of polygon edges, represents the area of the truncated background mesh.

[0030] S4, calculate the shape function at each Gaussian point, and based on the integral area of each Gaussian point in S3, discretize and solve the meshless balance equation, specifically:

[0031] S4.1, calculate the meshless nodes in the support domain of each Gaussian point, see formula (6), and solve the shape function of each Gaussian point based on the radial basis function of MQ, see formula (7);

[0032] (6)

[0033] In the formula, (x i , y i ) is the meshless node coordinate, (X k , Y k ) is the meshless node coordinate, and r dis is the radius of the Gaussian point support domain;

[0034] (7)

[0035] In the formula, represents the displacement vector of the Gauss point; m represents the number of meshless nodes in the support domain of the Gauss point respectively; represents the displacement vector of the meshless node in the support domain; represents the shape function vector of the Gauss point;

[0036] S4.2, if the meshless nodes in the support domain of the Gauss point in S4.1 belong to different material partitions, the shape function of the Gauss point in formula (7) needs to be optimized, so as to better simulate the interaction between different materials, see formula (8) and formula (9);

[0037] (8)

[0038] (9)

[0039] In the formula, m represents the total number of meshless nodes in the support domain of the Gauss point; n1 represents the number of meshless nodes in the first material partition; n2 represents the number of meshless nodes in the second material partition; represents the optimized shape function vector of the Gauss point; and represents the constant parameter of the optimized shape function of the Gauss point;

[0040] S4.3, based on the integral area of each Gauss point in S3, the meshless equilibrium equation is discretized by formula (10);

[0041] (10)

[0042] In the formula, Omega represents the integral domain of the Gauss point, K is the stiffness matrix, B is the displacement-strain matrix, D is the material constitutive matrix, and F is the load vector.

[0043] S4.4, formula (10) is solved and the calculation result is output, and then the meshless analysis considering the background mesh partition refinement and truncation optimization is completed.

[0044] Compared with the prior art, the beneficial effects of the present application are:

[0045] (1) The present application locally refines the background mesh at the boundary under the premise that the overall background mesh size remains unchanged, increases the calculation burden at a small cost, and realizes the refinement of the key area of the model boundary;

[0046] (2) In the present application, the background mesh can be truncated into a polygon, and the polygon is used as the integral area of the Gauss point, which enhances the representation ability of the meshless method for complex boundaries, and thus improves the simulation accuracy;

[0047] (3) When the meshless nodes in the support domain of the Gauss point belong to different material partitions, the shape function of the Gauss point is optimized,

[0048] Thus, the interaction between the two materials is better captured;

[0049] To sum up, the present application aims at the problem of the decline of simulation accuracy of the meshless method in the boundary area, and proposes an improved algorithm of background mesh partition refinement and truncation optimization, which more truly captures the dynamic boundary evolution of the model in the deformation process. That is, the present application improves the representation ability of the meshless method for complex boundaries by local four-part refinement and truncation optimization of the background mesh, improves the problem of the decline of simulation accuracy of the meshless method at the boundary, and effectively expands the application range of the meshless method. BRIEF DESCRIPTION OF DRAWINGS

[0050] Figure 1 It is a main flow diagram of the method of the present application;

[0051] Figure 2 It is a background mesh refinement diagram near the geometric boundary of the model (one refinement);

[0052] Figure 3 It is a background mesh refinement diagram near the material boundary of the model (one refinement);

[0053] Figure 4 It is a high point generation diagram in the refined background mesh;

[0054] Figure 5 It is a background mesh truncation diagram of the model boundary;

[0055] Figure 6 It is an improved Gaussian point function diagram;

[0056] Figure 7 It is a geometric model of the cantilever beam end stress;

[0057] Figure 8 It is a meshless node model of the cantilever beam end stress;

[0058] Figure 9 It is a cantilever beam node deformation diagram at step=500;

[0059] Figure 10 It is a cantilever beam x-direction normal stress at step=500;

[0060] Figure 11 It is a cantilever beam y-direction normal stress at step=500;

[0061] Figure 12 It is a cantilever beam xy-plane shear stress at step=500;

[0062] Figure 13 It is a geometric model of the foundation bearing capacity problem;

[0063] Figure 14A meshless node model for foundation bearing capacity problem

[0064] Figure 15 A foundation node deformation graph for step=500

[0065] Figure 16 A foundation node deformation graph for step=1000

[0066] Figure 17 A foundation node deformation graph for step=1500

[0067] Figure 18 A foundation node deformation graph for step=2000 DETAILED DESCRIPTION

[0068] The application will be further described below in conjunction with the accompanying drawings and specific embodiments, but the scope of protection of the application is not limited thereto.

[0069] Reference Figure 1 A meshless analysis method considering background mesh partition refinement and truncation optimization, comprising the following steps:

[0070] S1, sequentially read in the boundary node information, internal node information, material information and load information of the meshless model;

[0071] S2, locate the background mesh near the boundary of the meshless model, and perform four-part refinement thereon, see Figure 2 and Figure 3 ;

[0072] S3, generate a plurality of Gaussian points in the background mesh locally refined in S2, see Figure 4 , and calculate the integral area of each Gaussian point for subsequent numerical integration, see Figure 5 ;

[0073] S4, calculate the shape function at each Gaussian point, see Figure 6 , and based on the integral area of each Gaussian point in S3, discretize and solve the meshless equilibrium equation.

[0074] To prove the effectiveness of the technical solution, the technical solution is used to respectively simulate the cantilever beam end under concentrated load (Example 1) and the eccentric compression problem of foundation on slope (Example 2).

[0075] Example 1: Cantilever beam end under concentrated load problem

[0076] S1, sequentially read in the boundary node information, internal node information, material information and load information of the cantilever beam model, specifically:

[0077] S1.1, see Figure 7For the geometric model of cantilever beam end stress, the model length L = 12 m, width H = 2 m, the distance between nodes is 0.5 m, the boundary nodes are 53, the internal nodes are 72, the total number of nodes is 125, see Figure 8 ;

[0078] S1.2, based on the average node spacing in S1.1, calculate the size of the rectangular background grid, the background grid size in the invention is 2 times the average distance between nodes, the initial background grid size is 1 m;

[0079] S1.3, the cantilever beam is simulated by using linear elastic model, the detailed parameters are shown in Table 1;

[0080] Table 1, cantilever beam linear elastic model parameters

[0081]

[0082] S1.4, the gravity of the cantilever beam is ignored, and a concentrated load is applied at the free end of the cantilever beam, the concentrated load is F1 = 240 kN, which is applied in 100 equal interval steps;

[0083] S2, positioning the background grid near the boundary of the cantilever beam model, and carrying out four subdivision, in detail:

[0084] S2.1, sequentially connect the boundary nodes in S1 to form a plurality of lines representing the model boundary;

[0085] S2.2, carry out cross test on each initial background grid in S1, sequentially determine whether the plurality of lines of the model boundary in S2.1 intersects with the initial background grid, when formula (1) and formula (2) are satisfied at the same time, it is considered that the two intersect:

[0086] (1)

[0087] (2)

[0088] In the formula, (x b1 ,y b1 ) is the coordinates of a node in the model boundary of S2.1; (x b2 ,y b2 ) is the coordinates of the next node in the clockwise direction in the model boundary of S2.1; (X i ,Y i ), (X i+1 ,Y i+1 ) represent the position coordinates of two nodes in the edges of the rectangular background grid;

[0089] S2.3, according to the number of refinement n, n times of local quadtree refinement is carried out on the initial background grid intersecting with the model boundary in S2.2, the size of the background grid after local quadtree refinement is seen from formula (3); for the initial background grid not intersecting with the model boundary in S2.2, the size is kept unchanged, in the embodiment n = 1.

[0090] (3)

[0091] In the formula, n represents the number of refinement, represents the size of the initial background grid in S1, represents the size of the background grid after refinement in S2.

[0092] S3, a plurality of Gaussian points are generated in the background grid after local quadtree refinement in S2, and the integral area of each Gaussian point is calculated for subsequent numerical integration, specifically:

[0093] S3.1, a plurality of Gaussian points are generated in the background grid after local quadtree refinement in S2, in the application, 2x2 Gaussian points are generated in each background grid;

[0094] S3.2, then formula (1) and formula (2) are used to judge whether the background grid after local quadtree refinement intersects with the model boundary;

[0095] S3.3, for the background grid not intersecting with the model boundary, the integral area of the internal Gaussian point is the area of the background grid; for the background grid intersecting with the model boundary, the model boundary can truncate the background grid into a polygon, then the integral area of the internal Gaussian point needs to be calculated according to formula (4) and formula (5). In order to better compare the refinement and truncation optimization effect of the background grid partition, four working conditions are simulated in the embodiment, working condition 1, the meshless method without considering boundary partition refinement and truncation optimization; working condition 2, the meshless method considering only boundary partition refinement; working condition 3, the meshless method considering only boundary truncation optimization; working condition 4, the meshless method considering boundary partition refinement and truncation optimization.

[0096] (4)

[0097] (5)

[0098] In the formula, (xx i ,yy i ) represents the node coordinates of the truncated background grid, i represents the i th node in the truncated background grid, = 1, 2, 3, …, k; k represents the number of polygon edges, represents the area of the truncated background grid.

[0099] S4, calculate the shape function at each Gauss point, and discretize and solve the meshless equilibrium equation based on the integral area of each Gauss point in S3; specifically:

[0100] S4.1, calculate the meshless nodes within the support domain of each Gauss point, see equation (6), and solve the shape function of each Gauss point based on the radial basis function of MQ, which can be represented by equation (11) by solving the constant matrix to be solved and , see equations (12) and (13), equation (11) can be re-expressed as equation (7);

[0101] (6)

[0102] where (x i ,y i ) is the meshless node coordinate, (X k ,Y k ) is the meshless node coordinate, and r dis is the radius of the Gauss point support domain;

[0103] (11)

[0104] (12)

[0105] (13)

[0106] (7)

[0107] where, represents the displacement vector of the Gauss point, m respectively represents the number of meshless nodes within the support domain of the Gauss point, represents the displacement vector of the meshless node within the support domain, represents the shape function vector of the Gauss point, and are constant matrices to be solved, and are radial basis matrices and additional linear basis matrices.

[0108] S4.2, this embodiment only involves one material, and there is no interaction between two materials, so there is no need to optimize the shape function;

[0109] S4.3, based on the integral area of each Gauss point in S3, the meshless equilibrium equation is discretized using equation (10);

[0110] (10)

[0111] In the formula, Ω represents a Gaussian point integration domain, K is a stiffness matrix, B is a displacement-strain matrix, D is a material constitutive matrix, and F is a load vector.

[0112] S4.4, solving formula (10) and outputting the calculation result, thereby completing the meshless analysis considering the background mesh partition refinement and truncation optimization. First, the normalized vertical displacement of the cantilever beam under four working conditions is given, and the results of the four working conditions are close to the large deformation theoretical solution. Among them, working condition 1 does not consider the partition refinement and truncation optimization of the background mesh near the dynamic geometric boundary during the deformation process, so the calculation accuracy is poor. Although the accuracy of working conditions 2 and 3 is improved, it still cannot reach the accuracy of the present application (working condition 4), see the table, which verifies that the present application can effectively improve the simulation accuracy of the meshless method at the boundary.

[0113] Table 2, comparison of normalized vertical displacement of cantilever beam under four working conditions

[0114]

[0115] Figures 9-12 The node deformation diagram, x-direction normal stress, y-direction normal stress, and xy-plane shear stress of the cantilever beam under concentrated load at step=500 of the present application (working condition 4) are given respectively.

[0116] Example two: eccentric compression problem of foundation on slope;

[0117] S1, the boundary node information, internal node information, material information and load information of the slope foundation and foundation model are read in sequence, specifically:

[0118] S1.1, referring to Figure 13 , the geometric model of the eccentric compression problem of the foundation on the slope, the model has a slope top length L1=100m, a slope bottom length L1+L2=200m, a foundation length L3=50m, a load distance from the left end of the foundation L4=30m, a slope height L5=100m, and a meshless node spacing of 0.5m. In this embodiment, the model boundary includes a geometric boundary, a material boundary between the foundation and the foundation, and a meshless node average spacing of 7.5m, referring to Figure 14 ;

[0119] S1.2, based on the average node spacing in S1.1, the size of the rectangular background mesh is calculated, and the background mesh size in the present application is twice the average spacing of the meshless nodes, and the initial background mesh size is 15m;

[0120] S1.3, the slope foundation soil body is simulated by an ideal elastoplastic model of Mohr-Coulomb criterion, and the strip foundation is simulated by a linear elastic model, and the related parameters are shown in Tables 3 and 4.

[0121] Table 3, ideal elastoplastic model parameters of slope foundation soil

[0122]

[0123] Table 4, strip foundation baseline elastic model parameters

[0124]

[0125] S1.4, at the center of the strip foundation, a 10m vertical downward displacement condition is applied, a total of 2000 calculation steps are completed.

[0126] S2, locate the background mesh near the model boundary, including the geometric boundary and the material boundary, and perform quadtree subdivision, in particular:

[0127] S2.1, sequentially connect the boundary nodes in S1 to form a polyline representing the model boundary;

[0128] For each initial background mesh in S1, a cross-test is performed to determine whether the polyline in S2.1 representing the model boundary intersects the initial background mesh. When both formula (1) and formula (2) are satisfied, it is considered that the two intersect:

[0129] S2.3, for the initial background mesh in S2.2 that intersects the model boundary, perform n times of local quadtree subdivision according to the subdivision number n. The size of the background mesh after local quadtree subdivision is shown in formula (3). For the initial background mesh in S2.2 that does not intersect the model boundary, the size remains unchanged. In this embodiment, n = 1.

[0130] S3, generate multiple Gaussian points in the background mesh after local quadtree subdivision in S2, and calculate the integral area of each Gaussian point for subsequent numerical integration; in particular:

[0131] S3.1, generate multiple Gaussian points in the background mesh after local quadtree subdivision in S2. In this invention, 2x2 Gaussian points are generated in each background mesh;

[0132] S3.2, subsequently, formula (1) and formula (2) are used to determine whether the background mesh after local quadtree subdivision intersects the model boundary;

[0133] S3.3, for the background mesh that does not intersect the model boundary, the integral area of the internal Gaussian points is the area of the background mesh. For the background mesh that intersects the model boundary, the model boundary can truncate the background mesh into a polygon, and the integral area of the internal Gaussian points needs to be calculated according to formula (4) and formula (5).

[0134] S4, calculate the shape function at each Gaussian point, and based on the integral area of each Gaussian point in S3, discretize and solve the meshless balance equation; in particular:

[0135] S4.1, calculate the meshless nodes in each Gaussian point support domain, see formula (6), and solve the shape function of each Gaussian point based on the MQ-based radial basis function, see formula (7), formula (11) to formula (13);

[0136] S4.2, in the embodiment, there are meshless nodes in the support domain of some Gaussian points, and the meshless nodes are in different material partitions, so the shape function of the Gaussian point needs to be optimized, see formula (7), so as to better simulate the interaction between different materials, see formula (8) and formula (9); in order to better compare the effect of shape function optimization, two working conditions are given in the embodiment, including working condition 1 without considering shape function optimization, and working condition 2 considering shape function optimization, that is, the present application.

[0137] (8)

[0138] (9)

[0139] In the formula, m represents the total number of meshless nodes in the support domain of the Gaussian point, n1 represents the number of meshless nodes in the first material partition, n2 represents the number of meshless nodes in the second material partition, represents the optimized Gaussian point shape function vector; and represents the constant parameter of the optimized Gaussian point shape function;

[0140] S4.3, based on the integral area of each Gaussian point in S3, the meshless equilibrium equation can be discretized as, see formula (10):

[0141] S4.4, based on the above equation, the calculation results are solved and output, and the meshless analysis considering the background mesh partition refinement and truncation optimization is completed. First, the basic reaction force when the displacement condition reaches 10m is given, and the results in the foregoing literature (KHOEIAR, TABARRAIEAR, GHAREHBAGHISA. H-adaptive mesh refinement for shear band localization in elasto-plasticity Cosserat continuum [J]. Communications in Nonlinear Science and Numerical Simulation, 2005, 10(3): 253-286) are selected as reference values, see table 5. It can be seen that the calculation results considering the optimization of the shape function are more close to the literature, which verifies that the present application can better reflect the interaction between different materials.

[0142] Table 5, comparison of basic reaction forces under different working conditions

[0143]

[0144] Figures 15-18 The node deformation diagrams of the foundation under the bearing capacity at different typical time steps in the working condition 2, i.e., step=500, step=1000, step=1500 and step=2000, are given respectively. With the increase of displacement, the strip foundation gradually deforms in rotation, which further causes the right end of the foundation to gradually penetrate into the slope foundation, and the extrusion trend inside the soil slope pointing to the slope direction.

[0145] The above-described embodiments are only expressions of the embodiments of the present application, but should not be construed as limiting the scope of the present application. It should be noted that, for those skilled in the art, without departing from the concept of the present application, a number of modifications and improvements can be made, which are all within the scope of the present application.

Claims

1. A meshless analysis method considering background mesh partitioning refinement and truncation optimization, characterized in that, The meshless analysis method includes the following steps: S1, sequentially read in the boundary node information, internal node information, material information and load information of the meshless model; S2, locate the background mesh near the boundary of the meshless model and refine it locally into four subdivisions; S3 generates multiple Gaussian points within the background mesh after local quartic refinement in S2, and calculates the integral area of ​​each Gaussian point for subsequent numerical integration. S4 calculates the shape function at each Gaussian point, and based on the integral area of ​​each Gaussian point in S3, discretizes and solves the meshless equilibrium equation, completing the meshless analysis considering background mesh partitioning refinement and truncation optimization.

2. The meshless analysis method considering background mesh partitioning refinement and truncation optimization according to claim 1, characterized in that, Specifically, S1 refers to: S1.1, Read in the boundary node information and internal node information of the meshless model, where the boundary nodes include geometric boundary nodes and material boundary nodes; S1.2, based on the boundary node information and internal node information in S1.1, obtain the average spacing of the gridless nodes and calculate the size of the rectangular background grid. Here, the background grid is defined as the initial background grid. S1.3, then read in the material information of the meshless model, including material number, material partition, and material parameters; S1.4, then read in the load information of the meshless model, including body forces, surface forces, nodal forces, and displacement conditions.

3. The meshless analysis method considering background mesh partitioning refinement and truncation optimization according to claim 2, characterized in that, In S1.2, the background grid size is twice the average spacing of the gridless nodes.

4. The meshless analysis method considering background mesh partitioning refinement and truncation optimization according to claim 2, characterized in that, Specifically, S2 is: S2.1 connects the boundary nodes in S1 in sequence to form a polyline representing the boundary of the model; S2.2, perform a cross-test on each initial background mesh in S1, and sequentially determine whether the polyline of the model boundary in S2.1 intersects with the initial background mesh. When both formula (1) and formula (2) are satisfied, they are considered to intersect: （1）； （2）； In the formula, (x b1 ,y b1 (x) represents the coordinates of a meshless node in the boundary of model S2.1; b2 ,y b2 (X) represents the coordinates of the next unmesh node in a clockwise direction on the boundary of model S2.1; i ,Y i ), (X) i+1 ,Y i+1 () represents the position coordinates of two nodes on each side of the rectangular background grid; S2.3 For the initial background mesh that intersects with the model boundary in S2.2, perform n local quartering refinements according to the refinement number n. The size of the background mesh after local quartering refinement is shown in formula (3). For the initial background mesh that does not intersect with the model boundary in S2.2, keep the size unchanged. （3）； In the formula, n represents the refinement number. This represents the initial background mesh size in S1. This indicates the size of the refined background grid in S2.

5. The meshless analysis method considering background mesh partitioning refinement and truncation optimization according to claim 4, characterized in that, Specifically, S3 is: S3.1 generates multiple Gaussian points within the background mesh after local quartic refinement in S2; S3.2, then use formula (1) and formula (2) to determine whether the background mesh after local quadrature refinement intersects with the model boundary; S3.3 For a background mesh that does not intersect with the model boundary, the area of ​​Gaussian points inside it is the area of ​​the background mesh; for a background mesh that intersects with the model boundary, the model boundary cuts the background mesh into polygons, and the area of ​​Gaussian points inside it needs to be calculated according to formula (4) and formula (5). （4）； （5）； In the formula, (xx) i yy i () represents the coordinates of a node in the truncated background mesh, where i represents the i-th node in the truncated background mesh, i = 1, 2, 3, ..., k; k represents the number of sides of the polygon. This represents the area of ​​the background grid after truncation.

6. The meshless analysis method considering background mesh partitioning refinement and truncation optimization according to claim 5, characterized in that, In S3.1, 2×2 Gaussian points are generated within each background grid.

7. A meshless analysis method considering background mesh partitioning refinement and truncation optimization according to claim 5, characterized in that, Specifically, S4 is: S4.1 Calculate the meshless nodes in the support domain of each Gaussian point, see formula (6), and solve the shape function of each Gaussian point based on the radial basis function of MQ, see formula (7). （6）； In the formula, (x i ,y i ) represents the coordinates of the unmesh node, (X) k ,Y k ) represents the coordinates of the unmesh node, r dis Let the radius of the supporting domain be the Gaussian point. （7）； In the formula, represents the displacement vector of the Gaussian point; m represents the number of meshless nodes in the support domain of the Gaussian point. This represents a displacement vector that supports no mesh nodes within the domain; Represents a Gaussian point-shaped function vector; S4.2 If the unmesh nodes in the Gaussian point support domain in S4.1 belong to different material partitions, then the shape function of the Gaussian point in formula (7) needs to be optimized, see formula (8) and formula (9). （8）； （9）； In the formula, m represents the total number of meshless nodes in the Gaussian point support domain; n1 represents the number of meshless nodes in the first material partition; and n2 represents the number of meshless nodes in the second material partition. This represents the optimized Gaussian point-shape function vector; and This represents the constant parameters of the Gaussian point shape function after optimization; S4.3, based on the integral area of ​​each Gaussian point in S3, the meshless equilibrium equation is discretized using formula (10); （10）； In the formula, Ω represents the Gaussian point integral domain, K is the stiffness matrix, B is the displacement-strain matrix, D is the material constitutive matrix, and F is the load vector. S4.4 Solve formula (10) and output the calculation results to complete the meshless analysis considering background mesh partitioning refinement and truncation optimization.

Images