CN114510775A - A 3D Space Curved Meshing Method for Complex Models - Google Patents

Abstract

The invention belongs to the technical field of aircraft design, aircraft numerical simulation and numerical calculation technology and grid division, and particularly relates to a complex model three-dimensional space curved grid division method. The invention provides a new control point and space curved grid point selection strategy, which can well meet the grid division requirement of a complex model under the condition of not improving the grid magnitude. The selection of control points and deformation points is reduced while the grid quality is ensured, a large amount of position calculation in the grid division process is reduced, and the grid division strategy is optimized. And a basis function selection mode based on the cosine value of the included angle is provided, on the premise of ensuring the grid quality, the calculation required by the deformation quantity can be greatly reduced by using the low-order polynomial basis function with a simple form, the time consumption of grid division is reduced, and the calculation force requirement is reduced. The problems of poor mesh division universality and poor mesh quality in the mesh division process of the existing complex model are effectively solved.

Description

A 3D Space Curved Meshing Method for Complex Models

technical field

The invention belongs to the technical fields of aircraft design, aircraft numerical simulation and numerical calculation, and grid division, and in particular relates to a three-dimensional curved grid division method for complex models.

Background technique

With the continuous development of computer technology, numerical simulation technology and numerical calculation technology are more and more widely used in the related research of aerospace vehicles. Compared with experimental methods such as wind tunnel test, numerical simulation technology has the advantages of low experimental cost, strong repeatability, and strong applicability. Through the numerical simulation of the aerospace vehicle, the flow characteristics, thermal characteristics and electromagnetic characteristics of the aerospace vehicle can be obtained, and then the design of the relevant parameters of the aerospace vehicle can be optimized according to the results of the numerical simulation. In the process of numerical simulation, various numerical calculation methods (such as finite element method, finite volume method, discontinuous Galerkin method, etc.) are widely used in aircraft design, and efficient and accurate meshing process is a key to these numerical calculation methods. step.

Meshing is the basis of high-precision numerical simulation and numerical calculation, and the quality of mesh directly affects the accuracy of numerical simulation and numerical calculation. And in the process of numerical simulation, the meshing process often occupies most of the workload of the entire numerical simulation process, so meshing, especially the meshing of complex models, is a very important research part in numerical simulation. There are higher requirements for meshing in the numerical simulation of complex models and complex conditions. Complex models usually contain high-curvature parts or locally tiny geometric adjacent features. In order to ensure the discrete accuracy and mesh quality, the linear tetrahedral mesh needs to continuously reduce the mesh size to match its geometric characteristics in such geometric feature regions, which further leads to the continuous increase of the mesh magnitude, which leads to an increase in the amount of calculation. and an increase in computation time. For complex shapes, poor mesh quality may also lead to divergence in the calculation process, resulting in calculation failure. Therefore, it is necessary to use some new technologies to improve the mesh quality under complex conditions such as supersonic speed or multi-scale, and curved mesh is a better way to solve the existing mesh quality problems.

The existing curved meshing methods basically include direct method and indirect method. The direct method is rarely used in numerical calculation due to the problems of complex calculation and poor adaptability. The indirect method can be regarded as a mesh deformation based on linear mesh division, and is a commonly used method at present. In the indirect method, mesh deformation mainly includes extrapolation method, spring method and interpolation method. The extrapolation method is difficult to apply in complex models; the spring method is more complicated, the calculation process is time-consuming and the stability is poor. The radial basis function method is an interpolation method that constructs an interpolation basis by interpolating the distance between points. Because of its good versatility and good quality of deformed grids, it has been widely used in dynamic grids and grid deformation in recent years. However, the radial basis function method usually requires a lot of calculations, and the solution process is very time-consuming and unstable.

SUMMARY OF THE INVENTION

In view of the above-mentioned problems or deficiencies, in order to solve the problems of poor universality of grid division and poor grid quality in the process of grid division of existing complex models, the present invention provides a three-dimensional curved grid of complex models. division method. By proposing a new control point and deformation point selection strategy, the selection of control points and deformation points is reduced while ensuring the quality of the mesh, reducing a large number of position calculations in the process of meshing, and optimizing the meshing strategy. Reduced meshing time and reduced computing power requirements.

The specific scheme is as follows, a method for dividing a complex model three-dimensional space curved mesh, including the following steps:

Step 1. Input the target 3D physical model and mesh parameters.

Step 2. Convert the target 3D physical model into 3D space geometric information, mark and store the computational domain boundary information.

Step 3. According to the three-dimensional spatial geometric information in step 2, T initial linear tetrahedral elements are generated based on the frontal advancement method and the three-dimensional constraint Delaunay method in turn, and the grid topology relationship of the T initial linear tetrahedral elements is constructed. The hierarchy of tetra (body) → face (face) → edge (edge) → node (point) stores mesh information.

Step 4. According to the mesh information of the initial linear tetrahedral unit generated in step 3 and the computational domain boundary information in step 2, the T initial linear tetrahedral units generated in step 3 are retrieved one by one.

If the t-th element is a boundary element, a boundary label is attached to this element and counted as the ѱ-th boundary element, t∈{1,2,3,…, T }. Traversing the T initial linear tetrahedral elements in step 3, a total of Ψ boundary elements can be obtained, ѱ∈{1,2,3,…,Ψ}.

Step 5. According to the Ψ boundary elements marked in step 4, retrieve the Ψ boundary elements one by one, and extract the boundary surface, boundary edge and boundary point information of the ѱth boundary element.

The boundary surface in the tetrahedral element contains three boundary edges and three vertices. The edges on all boundary surfaces in the Ψ boundary elements are truncated to generate boundary surface truncation points. The generated boundary surface truncation point is projected to the computational domain boundary in step 2 to obtain the projection point, and there is a one-to-one correspondence between the projection point and the boundary surface truncation point.

Step 6. According to the information of the Ψ boundary elements obtained in step 5, the surface, edge and point information of the non-boundary surfaces of the Ψ boundary elements are retrieved one by one.

Step 7. According to the information of the Ψ boundary units obtained in step 5, the adjacent units of the Ψ boundary units are retrieved one by one, and the information of the faces, edges and points of the non-adjacent faces in the adjacent units is extracted.

The adjacent units: two units in the tetrahedron unit that share one face are defined as adjacent units, there are four adjacent units of the non-boundary unit, and the adjacent units of the boundary unit are less than four but at least one.

Step 8. Traverse all boundary surface truncation points and their corresponding projection points, and calculate the cosine of the angle between the boundary surface truncation points and their corresponding projection points.

Step 9. Substitute all the truncation points obtained by the truncation in steps 5, 6 and 7 into the radial basis function calculation formula, and the value obtained by the radial basis function calculation formula is the truncation point deformation variable.

The control point is the truncation point of the boundary surface in step 5, that is, a subset of the truncation points of all boundary surfaces; the weight is first calculated by the control point and the projection point; The sequence of the local numbers of the loop traverses the truncation point of the boundary surface to calculate the displacement of the control point, and repeats the point to skip.

Step 10. Substitute the deformation variables calculated in step 9 into the cutoff point coordinates to calculate the spatial curvature grid coordinates (x _c , y _c , z _c ).

Step 11. The spatially curved grid points calculated in step 10 are used as new grid nodes, inserted into grid data, and the grid topology relationship is updated. Then, the new spatial curve grid data is output in the corresponding grid data file format, and the output grid data file is imported into the visualization.

Further, the specific process of the step 3 is as follows: first, the frontal advance method is used to generate a layered grid of W layer, where W≥1, and W can be selected by the user according to the actual situation (W The larger the space curved grid, the higher the accuracy, but the greater the amount of calculation, the value of W is W∈{1,2,3,…,10}).

Then in the remaining computational domain, the initial linear tetrahedral element is generated by using the three-dimensional constraint Delaunay method with the current hierarchical mesh surface as the initial constraint.

Further, when calculating the position of the truncation point in the step 5, the scale parameter λ takes 1, which is the midpoint of the edge.

Further, if a higher-precision curved mesh is to be constructed in the step 5, new truncation points are continued to be generated on the basis of the first truncation of each boundary element. Then, the new first-order truncation points obtained respectively are projected to the boundary of the computational domain to obtain the corresponding projection points, which are stored in the truncation point data class and the projection point data class.

Further, in the step 9, the basis function in the radial basis function is a tightly supported radial basis function.

Further, in the step 9, the radial basis function uses a quadratic polynomial.

Further, in the step 11, the grid data file format output adopts CGNS format.

Further, in the step 6, the same method as in the step 5 is used to obtain the cutoff point of the non-boundary interface.

Further, in step 7, the same method as in step 5 is used to truncate each edge to obtain a truncation point of adjacent cells.

Further, in the step 7, for the multi-layer space curved grid, the adjacent units of the adjacent units are further cyclically searched.

The three-dimensional space curved meshing method for complex models proposed by the present invention is realized by truncating the original linear tetrahedron to generate curved tetrahedral elements, which can well meet the meshing requirements of complex models without increasing the magnitude of the mesh. , the initial linear mesh is based on the frontal advancement method and the 3D constrained Delaunay method. In the process of generating the space curved grid points, only the truncation point of the boundary surface in step 5 is selected as the control point, and the deformation point also only selects the truncation point of the boundary element and the adjacent element of the boundary element, which greatly reduces the number of control points and deformation points. The calculation related to the deformation points, that is, the spatially curved grid points, is greatly reduced. Compared with the prior art, step 8 innovatively proposes a basis function selection method based on the cosine value of the included angle; on the premise of ensuring the quality of the grid, the use of a simple low-order polynomial basis function can greatly reduce the amount of deformation. required calculation. The reduction of deformation point calculation and the simplification of the basis function form greatly improve the stability and speed of meshing.

To sum up, the present invention proposes a new control point and deformation point selection strategy, which reduces the selection of control points and deformation points while ensuring the quality of the mesh, and reduces a large number of position calculations in the process of mesh division. The grid division strategy is optimized, which reduces the time-consuming of grid division and reduces the computing power requirement. It can well meet the meshing requirements of complex models without increasing the mesh size.

Description of drawings

1 is a schematic diagram of a three-dimensional space curved grid of an embodiment target model;

Figure 2 is a comparison diagram of a blunt cone (10°) linear/curved grid;

Figure 3 is a schematic diagram of the cross-section (two-dimensional) of adjacent units;

Figure 4 is a schematic diagram of the X-51A linear grid;

Figure 5 is a schematic diagram of the X-51A space curved grid;

FIG. 6 is a flow chart of the present invention.

Detailed ways

The technical solutions of the present invention will be described in detail below with reference to the accompanying drawings and embodiments.

Referring to the accompanying drawings, a method for dividing a complex model three-dimensional curved mesh includes the following steps:

step 1. Enter the target 3D physical model and mesh parameters.

According to the 3D physical model (dimension parameters, formulas, drawings, physical objects, etc.), the digital model based on the ACIS kernel is drawn, and the grid parameter settings such as grid size and boundary type are read through the program.

Step 2. Convert the target 3D physical model into 3D space geometric information, mark and store the computational domain boundary information.

Create the Model class, import the three-dimensional spatial geometric information in the digital model into the Model class, and establish the boundary mark according to the set boundary type (the program is the member BCtype of the Model class).

Step 3. According to the three-dimensional spatial geometric information in step 2, T initial linear tetrahedral elements are generated based on the frontal advancement method and the three-dimensional constraint Delaunay method in turn, and the tetra body is constructed according to the mesh topology relationship of the T initial linear tetrahedral elements. The hierarchy of →face →edge →node points stores mesh information.

Read the grid parameters and the 3D space geometric information in the Model class, first use the frontal advancement method to advance from the surface of the calculation area to the inside to generate multi-layer (two layers in the spherical embodiment) grid unit, and then use the 3D constraint Delo The inner (Delaunay) method generates a mesh of the remaining area with the current array as the constraint. Construct tetraArray class, faceArray class, edgeArray class, nodeArray class at the same time to store unit information.

Step 4. According to the mesh information of the initial linear tetrahedral unit generated in step 3 and the computational domain boundary information in step 2, the T initial linear tetrahedral units generated in step 3 are retrieved one by one.

If the t-th element is a boundary element, attach a boundary label to the element and count it as the ѱ-th boundary element, t∈{1,2,3,…, T }; traverse the T initial linear tetrahedral elements in step 3 , a total of Ψ boundary elements can be obtained, ѱ∈{1,2,3,…,Ψ}.

In this embodiment: insert the boundary marker BCtype into the tetraArray class, and establish a corresponding relationship. Then insert the BCtype in the tetraArray class (the member BCtype of the tetraArray class in the program) into the faceArray class, edgeArray class, nodeArray class, and mark the boundary surface, boundary edge and boundary point in the unit.

Step 5. According to the Ψ boundary elements marked in step 4, retrieve the Ψ boundary elements one by one, and extract the boundary surface, boundary edge and boundary point information of the ѱth boundary element;

A boundary surface in a tetrahedral element contains three boundary edges and vertices P ₁ (x ₁ , y ₁ , z ₁ ), P ₂ (x ₂ , y ₂ , z ₂ ) and P ₃ (x ₃ , y ₃ , z _{3 )} ), x, y and z respectively represent the coordinate values of the point in three dimensions in the three-dimensional Cartesian coordinate system; the first edge corresponds to vertices P ₁ and P ₂ , the second edge corresponds to P ₁ and P ₃ , and the third edge corresponds to P 1 and P 3 . The edges correspond to P ₂ and P ₃ . Truncate edges on all boundary faces in Ψ boundary elements (tetrahedral element faces contain three edges) to generate boundary face truncation points.

The position of the truncation point of the first boundary edge (similar to the second and third boundary edges) is calculated by vertices P ₁ (x ₁ , y ₁ , z ₁ ) and P ₂ (x ₂ , y ₂ , z ₂ ) , its calculation formula is:

and project the generated boundary surface truncation point P _H (x _n,1 ,y _n,1 ,z _n,1 ) ( H represents the boundary surface truncation point label) to the computational domain boundary in step 2 to obtain the projected point P _sub ( x _n,1 ,y _n,1 ,z _n,1 ) ( sub represents the label of the projection point, and there is a one-to-one correspondence between the projection point and the truncation point of the boundary surface), λ is the scale parameter (in most cases, the calculation step is omitted λ=1 can be directly taken as the midpoint of the edge. In special complex models such as high curvature, it can be adjusted based on the physical model curvature and radial basis function, and the adjustment range is λ∈[0.5,2].

Every time a truncation point is generated, the truncation point P _H (x _n,1 ,y _n,1 ,z _n,1 ) is stored in the truncation point data class P _H (x _N,1 ,y _N,1 ,z _{N, 1} ), and count (n in P _H represents the nth truncation point, and the maximum count value is N ). At the same time, all the projection points P _sub (x _n,1 ,y _n,1 ,z _n,1 ) are obtained by projecting them to the boundary surface and also stored in the projected point data class P _sub (x _N,1 ,y _N,1 ,z _{N ) ,1} ), and count (n in P _sub represents the nth projection point, and the maximum count value is also N ).

Further, if a higher-precision curved mesh is to be constructed, new truncation points are continued to be generated on the basis of the first truncation of each boundary element. Take the points P _H (x _n,1 ,y _n,1 ,z _n,1 ) and P ₁ (x ₁ ,y ₁ ,z ₁ ), and P _H (x _n,1 ,y _n,1 ,z _{n , 1} ) and P ₂ (x ₂ , y ₂ , z ₂ ) are the endpoints, respectively, according to

Generate a truncation point, at this time (x _i , y _i , z _i ) is a truncation point P _H (x _n,1 ,y _n,1 ,z _n,1 ), (x _j ,y _j ,z _j ) is are the values of vertices P ₁ (x ₁ , y ₁ , z ₁ ) and P ₂ (x ₂ , y ₂ , z ₂ ), respectively.

Then, the obtained truncation points P _H (x _n,21 ,y _n,21 ,z _n,21 ) and P _H (x _n,22 ,y _n,22 ,z _n,22 ) are projected to the computational domain boundary Get P _sub (x _n,21 ,y _n,21 ,z _n,21 ) and P _sub (x _n,22 ,y _n,22 ,z _n,22 ) and store them in the corresponding truncation point data class and projection point data class. In this way, more high-order truncation points and high-order projection points are obtained. Generally, each edge is truncated once to obtain a first-order truncation point or truncates twice to obtain three second-order truncation points, which can meet the grid accuracy requirements. Although more truncation points can be cyclically intercepted to improve the grid accuracy, the more truncation times, the greater the amount of computation, the slower the computation speed, and the greater the computing power requirement. Therefore, it is not recommended to use a calculation scheme higher than the second-order truncation point in meshing.

In this embodiment, the boundary marker BCtype is retrieved in the tetraArray class, and when the unit is a boundary unit, the edgeArray class and nodeArray class corresponding to the unit are retrieved. Further search for the boundary marker BCtype of the edge. When the boundary edge is retrieved, the truncation point P _H (x _n,1 ,y _n,1 ,z _n,1 ) and the projection point P are obtained through the point coordinates and formula corresponding to the boundary edge _sub (x _n,1 ,y _n,1 ,z _n,1 ). After the truncation point P _H (x _n,1 ,y _n,1 ,z _n,1 ) and the projection point P _sub (x _n,1 ,y _n,1 ,z _n,1 ) are obtained, they are stored in the truncation point data class P _H (x _{N ,1} ,y _{N ,1} ,z _{N ,1} ) and P _sub (x _{N ,1} ,y _{N ,1} ,z _{N ,1} ), and count (the maximum count value of n is N ).

Step 6. According to the information of the Ψ boundary elements obtained in step 5, retrieve the surface, edge and point information of the non-boundary planes of the Ψ boundary elements one by one; adopt the same method in step 5 to obtain the non-boundary plane truncation point P _inr (x _{m ,1} ,y _m,1 ,z _m,1 ), stored in the truncation point data class P _inr (x _M,1 ,y _M,1 ,z _M,1 ), and counted, inr represents the non-boundary interface truncation point label , the maximum count value of m in Pinr is _denoted as M.

Step 7. According to the information of the Ψ boundary units obtained in step 5, the adjacent units of the Ψ boundary units are retrieved one by one, and the information of the faces, edges and points of the non-adjacent surfaces in the adjacent units is taken out. Further loops to retrieve adjacent cells of adjacent cells.

The adjacent units: two units in the tetrahedron unit that share one face are defined as adjacent units, there are four adjacent units of the non-boundary unit, and the adjacent units of the boundary unit are less than four but at least one.

Use the same method in step 5 to truncate each edge to obtain the adjacent unit truncation point P _adj (x _k,1 ,y _k,1 ,z _k,1 ) and store it in the truncation point data class P _adj (x _{K, 1} , y _K,1 ,z _K,1 ), and count, the maximum count value of k in Padj is K, and _adj represents the label of the truncation point of the adjacent unit.

In this embodiment, the class member tetra[A,B] of the faceArray class is retrieved, and A and B are the global numbers of the two units to which the face belongs. Then according to the global number A or B, the adjacent units are retrieved in the tetraArray class, and the information of the edges and points of the non-adjacent faces is taken out.

Step 8. Traverse all boundary surface truncation point data class P _H (x _N,1 ,y _N,1 ,z _N,1 ) and projected point data class P _sub (x _N,1 ,y _N,1 ,z _N,1 ) point, calculate the boundary surface truncation point P _H (x _n,1 ,y _n,1 ,z _n,1 ) and its corresponding projection point P _sub (x _n,1 ,y _n,1 ,z _n,1 ) between The cosine value of the included angle, denoted as cos( θ ) _H-sub , is stored in the corresponding position of each point in the nodeArray class.

where x _n,H ,y _n,H ,z _n,H are the three-dimensional coordinate values of the nth boundary surface truncation point, representing x _{n, sub} ,y _{n, sub} ,z _{n, sub} for the nth projection point 3D coordinate value.

Step 9. Substitute all truncation points obtained by truncation in steps 5, 6 and 7 into the radial basis function calculation formula, and the value of f ( x, y, z ) obtained by the radial basis function calculation formula is the truncation point shape variable; into the temporary variable mesh-distancecurved (total N+M+K ). The radial basis function is calculated as:

为权重，I _N 为控制点总数且I _N =N；

为径向基函数，

表示全部截断点到对应控制点的欧式距离，控制点

为步骤5中边界面截断点P _H （x_N,1,y_N,1,z_N,1），即全部边界面截断点的子集；where α is the coefficient,

is the weight, I _N is the total number of control points and I _N = N ;

is the radial basis function,

Indicates the Euclidean distance from all truncation points to the corresponding control point, the control point

is the boundary interface truncation point PH (x _{N ,1} , y _N,1 , z _N,1 ) in step 5, that is, a subset of all boundary interface truncation points;

In the embodiment, based on the maximum angle cosine value of the truncation point (that is, the space curved grid point), the radial basis function calculation formula used in the present invention is:

Radial basis functions only use quadratic polynomials in most cases, and the use of low-order polynomial basis functions with simple forms can greatly reduce the calculation required for deformation variables; only in complex cases such as high curvature in complex models, high-order polynomials are used. Meet mesh quality requirements.

，计算权重时

中X代表边界面截断点（即为步骤5中边界面截断点P _H （x_N,1,y_N,1,z_N,1），所以总数为I _N =N；计算权重时f (x,y,z)为已知量，代表控制点位移量，f (x,y,z)取同一边界面上截断点到投影点的欧式距离的最小值，以避免有平凡解的情况出现；First calculate the weights through the control points (boundary interface truncation points) and projection points

, when calculating the weights

where X represents the boundary surface cut-off point (that is, the boundary surface cut-off point P _H (x _N,1 ,y _N,1 ,z _N,1 ) in step 5, so the total number is I _N = N ; when calculating the weight, f ( x , y, z ) is a known quantity, representing the displacement of the control point, f ( x, y, z ) takes the minimum value of the Euclidean distance from the truncation point to the projection point on the same boundary surface to avoid the occurrence of trivial solutions;

f ( x,y,z )=min{║ X _b - X _sub ║}

In the case of the first-order truncation point, X _b ,b∈{1,2,3} are the three truncation points corresponding to the three boundary edges on the same boundary surface, and X _sub is the coordinate value of the projection point corresponding to the truncation point; then According to the order of the global numbering of the outer loop of the boundary surface, then the local number of the inner loop of the boundary surface, traverse the truncation point of the boundary surface to calculate the displacement of the control point, and repeat the point skipping;

后代入截断点计算f (x,y,z)，此时

为已知量，f (x,y,z)为待求量；计算网格点位移时

中X=X _j ，X _j 代表步骤5中边界面截断点、步骤6中非边界面截断点和步骤7中相邻单元截断点求得的全部截断点，截断点总数I _SUM =N+M+K；get weight

The descendants enter the truncation point to calculate f ( x, y, z ), at this time

is a known quantity, f ( x, y, z ) is the quantity to be determined; when calculating the displacement of grid points

where X = X _j , X _j represents the boundary surface truncation point in step 5, the non-boundary surface truncation point in step 6 and the truncation point of adjacent units in step 7. The total number of truncation points is I _SUM = N+ M + K;

Step 10. Substitute the deformation variable f ( x,y,z ) calculated in step 9 into the cut-off point coordinates (x _J,1 ,y _J,1 ,z _J,1 ) to calculate the space curved grid coordinates (x _c ,y _c , z _c ), its calculation formula is:

(x _c , y _c , z _c ) = (x _J,1 ,y _J,1 ,z _J,1 ) + f ( x,y,z ),J∈{ N } , {M} , {K}

Traverse the temporary variable mesh-distancecurved to take out the deformation variable corresponding to each truncation point, and add it to the initial position of the truncation point to obtain the position of the truncation point after deformation (the position of the added truncation point is the space curve grid point), which is Coordinates of the space curved grid points. Then save the coordinates of the obtained space curve grid points into the truncation point class to replace the original truncation point coordinates.

Step 11. The high-order grid points (x _c , y _c , z _c ) calculated in step 10 are used as new grid nodes, that is, the coordinates of the spatially curved grid points calculated in step 10 are read from the truncation point class, as the newly added grid points. The grid nodes node5 to node10 are inserted into the grid data, and the grid topology relationship is updated; then the new spatial curve grid data is output in the corresponding grid data file format (such as CGNS format), and the output grid data File import visualization.

In this embodiment, three target objects are used as examples to verify the effectiveness of the present invention; the spatial curved grid results of final grid division in each embodiment are shown in Figures 1 and 2 on the right, and Figure 5 .

Figure 1 is a sphere as the target object, the left is a global surface view of the space curved grid, the middle is a cross-sectional view of the space curved grid, and the right is a partial enlarged view.

Figure 2 is the result of a blunt cone (10°) space curved grid. The left picture is the initial linear tetrahedral mesh obtained in step 3, and the right picture is the schematic diagram of the completed space curved mesh.

Figures 4 and 5 are results plots for the X-51A aircraft. FIG. 4 is the initial linear tetrahedral mesh obtained in step 3, and FIG. 5 is a schematic diagram of the completed space curved mesh.

It can be seen from the comparison between the left picture and the right picture in Fig. 2 and the comparison between Fig. 4 and Fig. 5 that the space curved grid division method proposed by the present invention realizes the grid quality and the basic topology structure unchanged. A high-precision space curve grid is obtained. The existing curved grid technology basically only realizes the object surface curved grid. Although the object surface curved grid improves the calculation accuracy at the boundary, the numerical calculation accuracy is often limited by the accuracy of the internal grid. The invention adopts a brand-new selection strategy of control points and deformation points (that is, space curved grid points), and extends the curved grid to the inner space grid to realize high-precision spatial curved grid division of complex models. To sum up, it can be seen that the present invention is practical and effective.

To sum up, it can be seen that, by proposing a new control point and deformation point selection strategy, the present invention can well meet the grid division requirements of complex models without increasing the grid magnitude. Compared with the prior art: 1. The present invention realizes the division of space curved grids, and the prior art is basically an object surface curved grid; 2. A new control point and deformation point selection strategy and calculation method are adopted, and the grid is not improved by adopting a new selection strategy and calculation method. Under the circumstance of the order of magnitude, it can well meet the meshing requirements of complex models; 3. The radial basis function is selected according to the cosine value of the included angle. The invention reduces the selection of control points and deformation points while ensuring the quality of the grid, reduces a large number of position calculations in the process of grid division, optimizes the grid division strategy, reduces the time consumption of grid division and reduces the computing power need. The reduction of deformation point calculation and the simplification of the basis function form greatly improve the stability and speed of meshing.

Claims (10)

1. A complex model three-dimensional space curved grid division method is characterized by comprising the following steps:

step 1, inputting a target three-dimensional physical model and grid parameters;

step 2, converting the target three-dimensional physical model into three-dimensional space geometric information, and marking and storing the boundary information of the calculation domain;

and 3, generating the three-dimensional space geometrical information according to the three-dimensional space geometrical information in the step 2 on the basis of a wavefront advancing method and a three-dimensional constrained Delaunay methodTAn initial linear tetrahedral unit, according toTConstructing a tetra body → face surface → edge → node point hierarchical structure to store grid information by the grid topological relation of the initial linear tetrahedron units;

step 4, according to the grid information of the initial linear tetrahedron units generated in the step 3 and the calculation domain boundary information in the step 2, the grid information generated in the step 3 is searched one by oneTAn initial linear tetrahedral unit;

if the t-th cell is a border cell, then the cell is given a border label and is counted as the ѱ -th border cell, t e { 1,2,3, …,T｝；in step 3 of traversalTObtaining a total of psi boundary units by using an initial linear tetrahedron unit, wherein ѱ is equal to { 1,2,3, …, psi };

step 5, retrieving psi boundary units one by one according to the psi boundary units marked in the step 4, and extracting boundary surfaces, boundary edges and boundary point information of ѱ th boundary units;

the boundary surface in the tetrahedron unit comprises three boundary edges and three vertexes, and edges on all the boundary surfaces in the psi boundary units are truncated to generate boundary surface truncation points; projecting the generated boundary surface truncation point to the calculation domain boundary in the step 2 to obtain a projection point, wherein the projection point and the boundary surface truncation point have a one-to-one correspondence relationship;

step 6, retrieving the surface, edge and point information of the non-boundary surface of the psi boundary units one by one according to the psi boundary unit information obtained in the step 5;

step 7, retrieving adjacent units of the psi boundary units one by one according to the psi boundary unit information obtained in the step 5, and extracting the information of surfaces, edges and points of non-adjacent surfaces in the adjacent units;

the adjacent unit: two units sharing a face in tetrahedral units are defined as adjacent units, four adjacent units are not boundary units, and less than four adjacent units are boundary units but at least one;

step 8, traversing all the boundary surface truncation points and the corresponding projection points thereof, and calculating the cosine value of the included angle between the boundary surface truncation points and the corresponding projection points thereof;

step 9, substituting all the truncation points obtained by truncation in the steps 5, 6 and 7 into a radial basis function calculation formula, wherein the value obtained through the radial basis function calculation formula is the deformation amount of the truncation point;

the control points are the boundary surface truncation points in the step 5, namely the subsets of all the boundary surface truncation points; firstly, calculating the weight through the control point and the projection point; then, traversing the boundary surface interception points according to the sequence of the global serial numbers of the outer layer cycle of the boundary surface and the sequence of the local serial numbers of the inner layer cycle of the boundary surface to calculate the displacement of the control points, and skipping the repeated points;

step 10, substituting the deformation amount calculated in the step 9 into the coordinate of the truncation point to calculate the coordinate of the space curved grid;

step 11, the high-order grid points obtained by calculation in the step 10 are used as newly added grid nodes, inserted into grid data and updated in the grid topological relation; and then outputting the new space curved grid data by using a corresponding grid data file format, and importing the output grid data file into visualization.

2. The method for dividing the three-dimensional space curved grid of the complex model according to claim 1, wherein in the step 3: firstly, a wavefront advancing method is adopted to advance the calculation domain surface serving as an initial wavefront to the inside to generate a hierarchical grid of a W layer, wherein W is larger than or equal to 1, the spatial curve grid precision is higher when W is larger, but the calculation amount is larger, and W belongs to { 1,2,3, …,10 };

and then, generating initial linear tetrahedral units by adopting a three-dimensional constrained Delaunay method in the residual calculation domain by taking the current hierarchical grid surface as an initial constraint condition.

3. The method for dividing the three-dimensional space curved grid of the complex model according to claim 1, characterized in that: and (5) when the position of the truncation point is calculated in the step (5), taking 1 as the proportional parameter lambda, namely the middle point of the edge.

4. The method for dividing the three-dimensional space curved grid of the complex model according to claim 1, characterized in that:

if a curved grid with higher precision is to be constructed in the step 5, continuously generating a new truncation point on the basis of the first truncation of each boundary unit;

and then projecting the obtained new first-order truncation points to the boundary of the calculation domain to obtain corresponding projection points, and storing the corresponding projection points into a truncation point data class and a projection point data class, and circularly obtaining more high-order truncation points and high-order projection points.

5. The method for dividing the three-dimensional space curved grid of the complex model according to claim 1, characterized in that: and in the step 9, the basis function in the radial basis function is a tight-branch radial basis function.

6. The method for dividing the three-dimensional space curved grid of the complex model according to claim 1, characterized in that: the radial basis function in step 9 uses a quadratic polynomial.

7. The method for dividing the three-dimensional space curved grid of the complex model according to claim 1, characterized in that: in the step 11, the format output of the grid data file adopts a CGNS format.

8. The method for dividing the three-dimensional space curved grid of the complex model according to claim 1, characterized in that: and step 6, calculating the non-boundary surface truncation point by the same method in step 5.

9. The method for dividing the three-dimensional space curved grid of the complex model according to claim 1, characterized in that: and 7, truncating each edge by adopting the same method in the step 5 to obtain a truncation point of the adjacent unit.

10. The method for dividing the three-dimensional space curved grid of the complex model according to claim 1, characterized in that: in step 7, for the multi-layer space curved grid, the neighboring cells of the neighboring cells are further retrieved in a loop.

Images