Three-dimensional finite element mesh rapid generation method based on octree

Abstract

A three-dimensional finite element grid rapid generation method based on an octree comprises the steps that an STL format file exported by a BIM model is read, a computational domain covering the whole domain and with the size being 2 < n > power is constructed, and an initial octree unit grid is generated; on the basis of spatial smoothness and spatial continuity constraints, identifying feature lines and grid boundary lines, and recursively dividing grids according to a 2: 1 balanced subdivision principle to obtain balanced octree grids; the boundary units are classified into boundary smoothing units and feature constraint units, and boundary grids are optimized through intersection point calculation, tolerance control node replacement and unit surface cutting in combination with a Laplace method; and identifying independent nodes of which the size ratio of adjacent units is 2: 1 in the internal grid, and realizing node sharing of the whole grid by newly adding a central node splitting unit. According to the method, efficient and high-quality grid automatic generation of the imperfect geometric model can be realized.

Description

Technical Field

[0001] This invention belongs to the technical field of finite element analysis, and specifically relates to a method for rapid generation of three-dimensional finite element meshes based on octrees. Background Technology

[0002] Finite element analysis (FEM) is a core technology in engineering numerical simulation, and mesh generation is a crucial preliminary step in FEM, its quality directly affecting the accuracy and efficiency of the analysis results. Among existing mesh generation techniques, the octree method is widely used due to its advantages such as regular structure and fast generation speed; however, traditional octree meshes have the following shortcomings: High requirements are placed on the integrity of the geometric model. For imperfect geometric models with defects such as gaps and overlaps that are common in engineering practice such as BIM, tedious geometric repair and topology correction are required in advance, resulting in low automation. The generated mesh boundaries are jagged and have poor fit with the boundary surfaces of the geometric model, requiring additional optimization. When the size ratio of adjacent elements is 2:1, independent nodes are easily generated, leading to inconsistent mesh topology and affecting the accuracy of finite element analysis. Distorted meshes are prone to occur during boundary cutting, and the element angle, aspect ratio and other indicators are difficult to meet the requirements of engineering analysis.

[0003] Existing open-source mesh generation tools such as TetGen and Gmsh can generate high-quality meshes, but they also rely on perfect geometric models and it is difficult to balance adaptability to complex geometric features with generation efficiency. Summary of the Invention

[0004] This invention provides a method for rapid generation of three-dimensional finite element meshes based on octrees to solve the problem of insufficient integrity in the generation of three-dimensional finite element meshes.

[0005] To solve the above-mentioned technical problems, the technical solution adopted by the present invention is as follows: A method for rapid generation of 3D finite element meshes based on octrees includes the following steps: Step 1: Read the STL format file exported from the BIM model, construct a computational domain that covers the entire area and has a size of 2 to the power of n, and generate the initial octree cell mesh. Step 2: Based on spatial smoothness and spatial continuity constraints, identify feature lines and grid boundary lines, and recursively divide the grid according to the 2:1 balanced partitioning principle to obtain a balanced octree grid. Step 3: Classify the boundary elements into smooth boundary elements and feature-constrained elements. Optimize the boundary mesh by calculating intersection points, replacing nodes with tolerance control, and cutting element surfaces, combined with the Laplace method. Step 4: Identify independent nodes in the internal mesh whose adjacent cell size ratio is 2:1, and split the cells by adding a central node to achieve a shared node across the entire mesh.

[0006] Furthermore, in step one, when constructing the computational domain, a computational domain covering the entire domain is constructed based on the three-dimensional coordinates of the triangular mesh. The size of the computational domain is set to a power of 2, and the calculation formula is as follows: ; in, Round the result up. , , These represent the maximum spatial span of the geometric model in the x, y, and z directions, respectively.

[0007] Furthermore, when generating the initial octree cell mesh in step one, the initial size is set. The computational domain is uniformly divided into cubes of equal size along the x, y, and z directions to form an initial octree cell grid. Each cell records its spatial coordinate range and information about its neighboring cells.

[0008] Furthermore, in step two, the method for determining spatial smoothness is as follows: calculate the angle between the normal vectors of each triangular mesh within the mesh subdomain. If the included angles do not exceed the preset threshold If the condition is met, the unit is determined to satisfy spatial smoothness; otherwise, it is not.

[0009] Furthermore, in step two, the method for determining spatial continuity is as follows: calculate the boundaries of each subdomain, and obtain the shared boundary lines between subdomains and their intersection with the unit surface; If there are continuous shared boundary lines between subdomains, and the cell surface intersects only one shared boundary line, then the cell is considered to satisfy spatial continuity; otherwise, it is not.

[0010] Furthermore, in step two, the recursive mesh generation process includes: S21. Set initial dimensions Based on the computational domain, an initial mesh model is generated, and step S22 is executed for each element in the mesh model. S22: Determine whether the current cell size and its adjacent cell sizes satisfy the 2:1 partitioning principle of a balanced octree. If yes, proceed to step S23; otherwise, proceed to step S25. S23: Determine whether the current unit size represents the geometric features of the geometric model. If yes, proceed to step S24; otherwise, proceed to step S25. S24: Determine if the current unit size is greater than the maximum size. If yes, then proceed to step S25; otherwise, backtrack to the previous level of recursion. S25: Determine if the current unit size is greater than the minimum size. If so, the current mesh is divided into 8 equal-sized sub-cells, and step S22 is executed on each of the 8 sub-cells; otherwise, the process backtracks to the previous level of recursion.

[0011] Furthermore, in step three, the classification rules for boundary units are as follows: boundary units without significant feature lines and whose subdomains satisfy spatial smoothness are boundary smooth units; boundary units containing feature lines, smooth within subdomains, and continuous between subdomains are feature-constrained units.

[0012] Furthermore, in step three, "intersection point calculation" specifically includes: for each boundary element, solving for the intersection points of the element edge line with the geometric boundary and the feature line with the element surface, and obtaining the intersection points by solving the simultaneous formulas and parameters. Three-dimensional coordinates; "Node replacement" specifically includes: for each node coordinate axis direction Independently determine that the size ratio of adjacent units in this direction is 2:1, and then define an independent tolerance benchmark for each direction. .

[0013] Furthermore, in step three, when performing the unit surface cutting process, the boundary smooth unit surface is divided into two parts by a single cutting line to form a three-node or four-node unit surface; the feature-constrained unit surface is divided into multiple regions based on the number of cutting lines with the feature point as the center; the polygonal unit surface is converted into a three-node surface by adding a center node, and the collinear distorted quadrilateral is divided into a three-node surface.

[0014] Furthermore, after the unit surface cutting process in step three, the boundary mesh generation is performed by the unit surface cutting process, which specifically includes: cutting lines forming cutting surfaces, classifying and combining three-node surfaces and four-node surfaces, directly constructing tetrahedral, pyramid, wedge and hexahedral elements by boundary smoothing elements, and constructing tetrahedral or pyramidal elements by feature constraint elements with feature points as vertices.

[0015] The present invention can achieve the following beneficial effects: 1. This method overcomes the bottleneck of traditional octree mesh generation technology's reliance on perfect geometric models. Through precise computational domain construction, balanced partitioning constraints, boundary classification optimization, and internal common node handling, it eliminates the need for tedious geometric repair and topological correction of imperfect BIM models with gaps, overlaps, or other flaws, enabling direct automated mesh generation. This method retains the advantages of octree structure regularity and fast generation speed while overcoming the shortcomings of traditional methods through multi-stage constraints and optimizations, thus improving its applicability and efficiency in practical engineering scenarios.

[0016] 2. By classifying and matching boundary elements with different geometric features, a tolerance-controlled node replacement strategy ensures the rationality of intersection positions. Targeted element face cutting methods prevent the generation of distorted elements. Furthermore, the Laplace method is used to adjust node coordinates, enabling the boundary mesh to accurately fit the original geometric surface and achieve continuous boundary curvature. The final generated boundary mesh's angles, aspect ratios, and other parameters all meet engineering requirements, significantly improving the accuracy of boundary condition application and the reliability of calculation results in finite element analysis.

[0017] 3. The spatial continuity constraint and internal mesh shared-node handling scheme solve the topological inconsistency problem caused by the 2:1 size ratio of adjacent elements in traditional octree meshes. The spatial continuity constraint avoids the generation of independent nodes from the source, while the internal mesh splitting strategy eliminates existing independent nodes, achieving topological consistency across the entire mesh with shared nodes. Simultaneously, the spatial smoothness constraint optimizes the element size distribution, resulting in uniform mesh element quality, with over 98% of elements having an angle range of 15° to 135° and a minimum Jacobian value greater than 0.2. This effectively ensures the accuracy of finite element analysis. Engineering examples verify that the displacement and stress field calculation results are highly consistent with those of manually generated high-quality hexahedral meshes. Attached Figure Description

[0018] The present invention will be further described below with reference to the accompanying drawings and embodiments: Figure 1 This is a flowchart of a rapid three-dimensional finite element mesh generation method based on an octree according to the present invention; Figure 2 This refers to the triangle faces and their ASCII representations in the STL format file of this invention; Figure 3 This is a schematic diagram of the octree mesh division of the 1 / 4 cylinder of the present invention, (a) is the initial mesh, (b) is the balanced octree mesh, and (c) is the octree mesh; Figure 4 This is a schematic diagram of the spatial smoothness of the discrimination unit of the present invention, (a) a unit that satisfies spatial smoothness, and (b) a unit that does not satisfy spatial smoothness; Figure 5 This is a schematic diagram of the spatial continuity of the discrimination unit in this invention, (a) a unit that satisfies spatial continuity, and (b) a unit that does not satisfy spatial continuity; Figure 6 The diagrams show the boundary element cutting types of the present invention: (a) a boundary smoothing element diagram, and (b) a feature constraint element diagram. Figure 7 This is a schematic diagram of node replacement intersection on the unit surface of the present invention, (a) boundary smooth unit surface, (b) feature constrained unit surface; Figure 8 This is the boundary smoothing unit surface processing method of the present invention; Figure 9 This is the feature constraint unit surface processing method of the present invention; Figure 10 The diagram illustrates the generation of boundary elements in this invention: (a) boundary element generation method at smooth boundary points, (b) boundary element generation method with feature constraints (feature lines), and (c) boundary element generation method with feature constraints (corner points). Figure 11 This is a schematic diagram of the node offset topology of the present invention, showing (a) boundary elements and (b) to (d) node replacement and offset processes; Figure 12 This is a schematic diagram of the discrete element surface of the present invention ((a) internal element, (b) pyramid element, (c) discrete element surface). Figure 13 This is the common node cutting method for the unit surface of the present invention; Figure 14 This is the non-closed geometric model with gaps in the verification results of this invention; Figure 15 This is a schematic diagram of the boundary mesh generation at the gap in this invention, showing (a) the replacement of intersection points of unit nodes, (b) the cutting and optimization of the boundary mesh, and (c) the boundary mesh generation result. Figure 16 The diagram shows the boundary mesh division in the verification results of this invention: (a) boundary background mesh, (b) boundary mesh after cutting, and (c) boundary mesh after optimization. Figure 17 These are the mesh generation results for different mesh sizes in the verification results of this invention; Figure 18 The geometric models of the underground cavern engineering of this invention are: (a) the core area model of the underground cavern, and (b) the excavation area model of the underground cavern engineering. Figure 19 The following is a perspective view of the underground cavern engineering mesh model of the present invention: (a) hexahedral mesh model, (b) algorithm-generated mesh model; Figure 20 The images show a comparison of the excavation portions of the two mesh models of the present invention: (a) the excavation portion of the hexahedral mesh model, and (b) the excavation portion of the algorithm-generated mesh model. Figure 21 The maximum principal stress comparison cloud map of section 1 of the present invention is shown in (a) hexahedral mesh model and (b) algorithm-generated mesh model. Figure 22 The minimum principal stress comparison cloud map of section 1 of the present invention is shown in (a) hexahedral mesh model and (b) algorithm-generated mesh model. Figure 23 The maximum principal stress comparison cloud map of section 2 of the present invention is shown in (a) hexahedral mesh model and (b) algorithm-generated mesh model. Figure 24The minimum principal stress comparison cloud map of section 2 of the present invention is shown in (a) hexahedral mesh model and (b) algorithm-generated mesh model. Figure 25 The horizontal displacement comparison cloud map of section 2 of the present invention is shown in (a) hexahedral mesh model and (b) algorithm-generated mesh model. Figure 26 The vertical displacement comparison cloud map of section 2 of this invention is shown in (a) hexahedral mesh model and (b) algorithm-generated mesh model. Detailed Implementation

[0019] To facilitate understanding of this application, a more complete description will be provided below with reference to the accompanying drawings, which illustrate embodiments of the present application. However, the present application can be implemented in many different forms and is not limited to the embodiments described herein. Rather, these embodiments are provided so that the disclosure of this application will be thorough and complete.

[0020] A method for rapid generation of 3D finite element meshes based on octrees, such as Figure 1 As shown, it includes the following steps: Step 1: Initial Octree Mesh Generation Input the STL format file exported from the BIM model, read the vertex coordinates and normal vector data of the triangular mesh from the STL file, store it as a triangular mesh dataset, and extract the surface contour and key feature regions of the geometric model. The STL format triangular mesh data is as follows: Figure 2 As shown.

[0021] The maximum and minimum coordinates of the computational geometric model in the x, y, and z directions are determined to obtain the maximum spatial span in each direction. A computational domain covering the entire area is constructed based on the three-dimensional coordinates of a triangular mesh. The size of the computational domain is set to a power of 2, and the calculation formula is as follows: ,in Round the result up. , , These represent the maximum spatial spans of the geometric model in the x, y, and z directions, respectively, to ensure that the computational domain completely covers the geometric model.

[0022] Set initial size The computational domain is uniformly divided into cubes of equal size along the x, y, and z directions to form an initial octree cell grid. Each cell records its spatial coordinate range and information about its neighboring cells.

[0023] Step 2: Recursively partition the balanced octree mesh The unit scale distribution is constrained by spatial smoothness and spatial continuity. Units that intersect with the triangular mesh are set as boundary units. Mesh edges where the normal vector changes abruptly are identified as feature lines, and isolated mesh edges are identified as mesh boundary lines. Triangular meshes are classified into different subdomains according to feature lines and mesh boundary lines.

[0024] Set the threshold for the angle between the normal vectors. Maximum size of boundary element and minimum size Based on the 2:1 balanced partitioning principle, the mesh elements are recursively cut, such as... Figure 3 As shown.

[0025] Specifically, the initial mesh cells are traversed, and it is determined whether a cell intersects with a triangular mesh. Intersecting cells are marked as boundary cells. The normal vectors of the triangular mesh are calculated, and cells with an angle exceeding the boundary value are marked as boundary cells. The mesh edges are marked as feature lines, and mesh edges without adjacent triangular faces are determined as mesh boundary lines.

[0026] Based on the feature lines and grid boundary lines, the triangular meshes intersecting the boundary cells are classified into different grid subdomains. The angle between the normal vectors of the triangular meshes in each subdomain is calculated to determine the spatial smoothness. The continuity of the shared boundary lines between subdomains and their intersection with the cell faces are calculated to determine the spatial continuity.

[0027] The grid is divided according to the following recursive logic, and the recursive division steps are as follows: S21. Set initial dimensions Based on the computational domain, an initial mesh model is generated, and step S22 is executed on each element in the mesh model.

[0028] S22: Determine whether the current cell size and its adjacent cell sizes satisfy the 2:1 partitioning principle of a balanced octree. If yes, proceed to step S23; otherwise, proceed to step S25.

[0029] S23: Determine whether the current unit size represents the geometric features of the geometric model. If yes, proceed to step S24; otherwise, proceed to step S25.

[0030] S24: Determine if the current unit size is greater than the maximum size. If so, proceed to step S25; otherwise, backtrack to the previous level of recursion.

[0031] S25: Determine if the current unit size is greater than the minimum size. If so, the current mesh is divided into 8 equal-sized sub-cells, and step S22 is executed on each of the 8 sub-cells; otherwise, the process backtracks to the previous level of recursion.

[0032] After the recursion terminates, a balanced octree mesh covering the geometric model and highlighting its features is obtained.

[0033] The specific method for determining spatial smoothness is as follows: calculate the angle between the normal vectors of each triangular face within the grid subdomain. If the included angles do not exceed the preset threshold If the element satisfies spatial smoothness, then it is determined that the element does not satisfy it; otherwise, it does not. Figure 4 As shown.

[0034] The specific method for determining spatial continuity is as follows: calculate the boundaries of each subdomain, obtain the shared boundary lines between subdomains and their intersection with the element surface. If there are continuous shared boundary lines between subdomains, and the element surface intersects only one shared boundary line, then the element is deemed to satisfy spatial continuity; otherwise, it is not satisfied. Figure 5 As shown.

[0035] Step 3: Boundary Mesh Cutting Optimization Based on geometric characteristics, boundary elements are divided into smooth boundary elements and feature-constrained elements; the intersection points of the boundary elements and the geometric boundary are solved, and relative tolerances are used. The controlled node replacement strategy replaces intersection points, generates cutting lines and cutting surfaces in node order, and processes the element surface into three-node or four-node element surfaces; the node coordinates are adjusted using the Laplace method to fit the original geometric surface. Specifically, the following steps are included: S31, Boundary Unit Classification The specific rules for classifying boundary units include traversing the boundary units; those without significant feature lines and whose subdomains satisfy spatial smoothness are classified as boundary smooth units; those containing feature lines, smooth within subdomains, and continuous between subdomains are classified as feature-constrained units, such as... Figure 6 As shown.

[0036] S32, Intersection Calculation and Node Replacement The specific rules for intersection point calculation include finding the intersection points of each boundary element with the geometric boundary and the feature lines with the element faces. The intersection points are obtained by solving simultaneous formulas and resolving parameters. The three-dimensional coordinates are given by the following simultaneous formulas: ; ; In the formula, , , The parameters to be solved are... The three vertices of the triangle face. For the endpoints of the line segment, The intersection points are to be solved; like Figure 7 As shown, the specific rules for node replacement include for each node. coordinate axis direction Independently determine that the size ratio of adjacent units in this direction is 2:1, and then define an independent tolerance benchmark for each direction. The formula for determining the value is as follows: ; in, ∈{ }, The node is Adjacent elements in a direction have the smallest element size.

[0037] For the intersection point There are adjacent nodes The criteria for determining intersection points by replacing nodes are as follows: ; If the intersection point satisfies the condition, then the adjacent node replaces the intersection point; if not, the projection point of the adjacent element edge replaces the intersection point. The determination formula is as follows: ; In the formula, For all adjacent nodes.

[0038] S33, Unit surface cutting process: For a smooth boundary element surface, it is divided into two parts by a single cutting line to form a three-node or four-node element surface, such as... Figure 8 As shown. For the feature-constrained unit surface, it is divided into multiple regions centered on the feature points and according to the number of cutting lines, such as... Figure 9 As shown, a new center node is added to a polygonal unit surface to transform it into a three-node surface, and a collinear distorted quadrilateral is divided into a three-node surface.

[0039] S34. Unit surface cutting process for boundary mesh generation: Cutting lines form cutting surfaces, classifying and combining three-node and four-node surfaces. Smooth boundary elements directly construct tetrahedral, pyramidal, wedge-shaped, and hexahedral elements. Feature-constrained elements construct tetrahedral or pyramidal elements using feature points as vertices. Figure 10 As shown.

[0040] S35, Boundary Optimization: Calculate the node offset coordinates using the Laplace method, such as... Figure 11 As shown, the boundary element nodes are migrated to the fitted surface position to eliminate jagged distortion and make the boundary curvature continuous.

[0041] Step 4: Internal Mesh Cutting Traverse the internal cells of the balanced octree mesh, identifying cells with a size ratio of 2:1 and their corresponding independent nodes. Add a new node at the center of the cell containing the independent node, dividing the cell into 6 pyramidal cells. Discretize the cell faces, as shown below. Figure 12As shown. For the case where there are independent nodes at the midpoint of the element surface edge, according to... Figure 13 The pyramid element is cut in the manner shown to form a combination of tetrahedral and pyramid elements, eliminating independent nodes and achieving full mesh sharing of nodes.

[0042] The verification results of the octree-based rapid generation method for three-dimensional finite element meshes proposed in this application are as follows: For non-closed geometric models with gaps, where the model features are non-closed geometric models with gaps on the surface, no geometric repair is required before direct input, such as... Figure 14 As shown; The generated result is as follows Figure 15 , Figure 16 As shown, the mesh spans geometric gaps and is unaffected by defects. By adjusting the mesh size or angle tolerance, detailed features of the geometric surface can be precisely captured, such as... Figure 17 As shown. Compared with traditional tools, the octree-based rapid generation method for 3D finite element meshes proposed in this application eliminates the need for tedious geometric corrections, significantly improving the level of automation.

[0043] For a certain underground cavern excavation model, the excavation area model is as follows: Figure 18 As shown in (a). To facilitate comparison and verification of the computational feasibility and accuracy of the target algorithm's mesh generation, a manually generated full hexahedral mesh model is used as a reference model, as shown in [image]. Figure 19 As shown in (a), the grid of its excavation area is as follows: Figure 20 As shown in (a), the initial model size is preset to 10, the maximum boundary element size is 2, the minimum boundary element size is 0.5, and the preset included angle threshold is 20°. The algorithm automatically generates the mesh model, as shown in (a). Figure 19 As shown in (b), the grid of its excavation area is as follows: Figure 20 As shown in (b).

[0044] The algorithm generates 1,317,008 mesh cells with 416,233 cell nodes; over 98% of the cells have angles between 15° and 135°, aspect ratios ≤10, minimum Jacobian value >0.2, average Jacobian value >0.9, and no distorted cells.

[0045] Therefore, the mesh quality meets the requirements of practical engineering finite element analysis, such as... Figures 21-26 As shown, the displacement and stress field distributions of the two models are highly consistent. The mesh quality generated by the method presented in this paper meets the requirements of finite element analysis, and the calculation results have good reliability. The computational feasibility and accuracy meet the standards.

[0046] The above description is only a preferred embodiment of this application and is not intended to limit this application. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of this application should be included within the protection scope of this application.

Claims

1. A method for fast generation of three-dimensional finite element mesh based on octree, characterized in that, The method comprises the following steps: Step one, read the STL format file exported by the BIM model, build a calculation domain covering the whole domain and the size of 2n, and generate an initial octree unit grid; Step two, based on the spatial smoothness and spatial continuity constraints, identify feature lines and grid boundary lines, recursively divide the grid according to the 2:1 balance subdivision principle, and obtain a balanced octree grid; Step three, classify the boundary elements into boundary smooth elements and feature constraint elements, and optimize the boundary grid through intersection point calculation, tolerance control node replacement and element face cutting combined with Laplace method; Step four, identify independent nodes in the internal grid whose adjacent element size ratio is 2:1, split the elements by adding a center node, and realize the full grid node sharing.

2. The method of claim 1, wherein: In step one, a calculation domain covering the whole domain is constructed based on the three-dimensional coordinates of the triangular face grid, and the size of the calculation domain is set to 2n, and the calculation formula is: ； wherein is the result rounded up, , , are the maximum spatial span of the geometric model in x, y, z direction, respectively.

3. The method of claim 1, wherein: When generating the initial octree unit grid in step one, set the initial size The calculation domain is uniformly divided into equal cubes along the x, y, and z directions to form an initial octree unit grid, and each unit records its spatial coordinate range and adjacent unit information.

4. The method of claim 1, wherein: In step two, the way to judge the spatial smoothness is to calculate the included angle of the normal vector of each triangular mesh in the grid sub-domain If the included angle is not more than the preset threshold , the judging unit is determined to satisfy the spatial smoothness; otherwise, it does not satisfy the spatial smoothness.

5. The method of claim 1, wherein: In step two, the judgment method of spatial continuity is to calculate the boundary of each subdomain, obtain the shared boundary line between subdomains and the intersection state with the element face; If there is a continuous shared boundary line between subdomains, and the element face only intersects with one shared boundary line, it is determined that the element satisfies the spatial continuity, otherwise it does not satisfy.

6. The method of claim 1, wherein: In step two, the recursive division step includes: S21, setting initial size Based on the calculation domain, an initial mesh model is generated, and each cell in the mesh model is executed step S22. S22: Determine whether the size of the current element and its adjacent element satisfies the 2:1 subdivision principle of balanced octree, if yes, execute step S23, otherwise execute step S25; S23: Determine whether the size of the current element represents the geometric feature of the geometric model, if yes, execute step S24, otherwise execute step S25; S24: determine whether the current cell size is larger than the maximum size If yes, proceed to step S25, otherwise backtrack to the previous recursion level. S25: Determine whether the current cell size is greater than the minimum size If yes, split the current grid into 8 sub-cells of equal size and perform step S22 for each of the 8 sub-cells, otherwise backtrack to the previous recursion level.

7. The method of claim 1, wherein: In step three, the classification rule of boundary elements is that the boundary element without significant feature line and the subdomain satisfying spatial smoothness is a boundary smooth element; the boundary element containing feature line, smooth in subdomain and continuous between subdomains is a feature constraint element.

8. The method of claim 1, wherein: In step three, the "intersection calculation" specifically includes: for each boundary element, solving the intersection of the element edge line and the geometric boundary, the intersection of the characteristic line and the element face, and obtaining the intersection point by solving the parameters through the simultaneous equations Three-dimensional coordinates; "Node replacement" specifically includes: for each node Coordinate axis direction , independently determine the adjacent unit size ratio of the direction to be 2:1, and then define an independent tolerance reference for each direction .

9. The method of claim 8, wherein: In step three, when performing element face cutting processing, the boundary smooth element face is divided into two parts according to a single cutting line, forming a three-node or four-node element face; the feature constraint element face is divided into multiple regions according to the number of cutting lines with the feature point as the center; the multi-edge element face adds a center node to become a three-node face, and the collinear abnormal quadrilateral is divided into a three-node face.

10. The method of claim 1 or 9, wherein: After the element face cutting processing in step three, the boundary grid generation is performed, which specifically includes: cutting line forms a cutting face, three-node and four-node faces are classified and combined, boundary smooth elements directly construct tetrahedron, pyramid, wedge and hexahedron elements, and feature constraint elements construct tetrahedron or pyramid elements with feature points as vertices.

Images