CN114510775A - Complex model three-dimensional space curved grid division method - Google Patents
Complex model three-dimensional space curved grid division method Download PDFInfo
- Publication number
- CN114510775A CN114510775A CN202111637505.5A CN202111637505A CN114510775A CN 114510775 A CN114510775 A CN 114510775A CN 202111637505 A CN202111637505 A CN 202111637505A CN 114510775 A CN114510775 A CN 114510775A
- Authority
- CN
- China
- Prior art keywords
- grid
- boundary
- point
- points
- truncation
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/15—Vehicle, aircraft or watercraft design
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
- G06T17/20—Finite element generation, e.g. wire-frame surface description, tesselation
- G06T17/205—Re-meshing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Geometry (AREA)
- General Physics & Mathematics (AREA)
- General Engineering & Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Aviation & Aerospace Engineering (AREA)
- Automation & Control Theory (AREA)
- Computer Graphics (AREA)
- Software Systems (AREA)
- Image Generation (AREA)
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
Technical Field
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.
Background
With the continuous development of computer technology, the application of numerical simulation technology and numerical calculation technology in the relevant research of aerospace aircrafts is more and more extensive. Compared with experimental methods such as wind tunnel tests and the like, the numerical simulation technology has the advantages of low experimental cost, strong repeatability, strong applicability and the like. The flow characteristics, the thermal characteristics and the electromagnetic characteristics of the aerospace craft can be obtained through the numerical simulation of the aerospace craft, and then the design of relevant parameters of the aerospace craft is optimized according to the numerical simulation result. Various numerical calculation methods (such as a finite element method, a finite volume method, a discontinuous Galerkin method and the like) are widely applied to aircraft design in a numerical simulation process, and an efficient and accurate grid division process is a key step of the numerical calculation methods.
The grid division is the basis for carrying out high-precision numerical simulation and numerical calculation, and the quality of the grid directly influences the precision of the numerical simulation and the numerical calculation. In the numerical simulation process, the meshing process usually occupies most of workload of the whole numerical simulation process, so that meshing, especially meshing of a complex model, is an important research part in the numerical simulation. The grid division has higher requirements in numerical simulation under complex models and complex conditions. Complex models often contain portions of high curvature or locally small geometrically neighboring features. In order to ensure the discrete precision and the grid quality, the linear tetrahedral grid needs to continuously reduce the grid size to match the geometric features in such geometric feature regions, which further leads to continuous increase of the grid magnitude, thereby leading to increase of the calculation amount and increase of the calculation time. For complex shapes, poor mesh quality may also lead to divergence in the computation process, which in turn leads to computation failures. Therefore, under complex conditions such as supersonic velocity or multi-scale conditions, it is necessary to improve the grid quality through some new techniques, and curved grids are a better way to solve the existing grid quality problem.
The existing curved grid division method basically comprises a direct method and an indirect method. The direct method is less used in numerical calculation due to problems of complicated calculation, poor adaptability and the like. The indirect method can be regarded as mesh deformation based on linear mesh division, and is a commonly used method at present. In the indirect method, the mesh deformation mainly includes an extrapolation method, a spring method, an interpolation method, and the like. The extrapolation method is difficult to be applied to a complex model; the theory of the spring method is complex, the calculation process consumes long time, and the stability is poor. The radial basis function method is an interpolation method for constructing an interpolation basis through the distance between interpolation points, and is widely applied to the fields of dynamic grids, grid deformation and the like in recent years due to good universality and good quality of deformed grids. However, the radial basis function method usually requires a large amount of calculation, and the solving process is very time-consuming and unstable.
Disclosure of Invention
Aiming at the problems or the defects, the invention provides a complex model three-dimensional space curved grid division method for solving the problems of poor grid division universality and poor grid quality in the existing complex model grid division process. By providing a new control point and deformation point selection strategy, 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, the grid division strategy is optimized, the time consumption of grid division is reduced, and the calculation power requirement is reduced.
The specific scheme is as follows, a complex model three-dimensional space curved grid division method comprises the following steps:
And 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 based on a wavefront advancing method and a three-dimensional constrained Delaunay (Delaunay) method in sequenceTAn initial linear tetrahedral unit, according toTThe mesh topological relation of the initial linear tetrahedral units constructs a hierarchical structure of tetra (body) → face) → edge → node, and stores mesh information.
Step 4, according to the initial linear four sides generated in the step 3The mesh information of the body unit and the calculation domain boundary information in step 2 are searched one by one, which is generated in step 3TAn 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 traversalTAnd obtaining a total of psi boundary elements from the initial linear tetrahedral units, wherein ѱ E is { 1,2,3, …, psi }.
And 5, retrieving the 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 tetrahedral unit comprises three boundary edges and three vertexes, and the edges on all the boundary surfaces in the psi boundary units are truncated to generate a boundary surface truncation point. And 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.
And 6, retrieving the surface, edge and point information of the non-boundary surfaces of the psi boundary units one by one according to the psi boundary unit information obtained in the step 5.
And 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 of the tetrahedral cells sharing a face are defined as adjacent cells, four of the adjacent cells are non-boundary cells, and fewer than four but at least one of the adjacent cells are boundary cells.
And 8, traversing all the boundary surface interception points and the corresponding projection points thereof, and calculating the cosine values of the included angles between the boundary surface interception points and the corresponding projection points thereof.
And 9, substituting all the truncation points obtained by truncation in the steps 5, 6 and 7 into the radial basis function calculation formula, wherein the value obtained by 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, the control point displacement is calculated by traversing the boundary surface interception points according to the sequence of the global serial number of the outer layer cycle of the boundary surface and the sequence of the local serial number of the inner layer cycle of the boundary surface, and the points are skipped repeatedly.
Step 10, substituting the deformation quantity calculated in the step 9 into the coordinate of the truncation point to calculate the coordinate (x) of the space curved gridc,yc,zc）。
And 11, inserting the space curved grid points obtained by calculation in the step 10 into grid data as newly added grid nodes, and updating 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.
Further, the specific process of step 3 is as follows: firstly, a wavefront advancing method is adopted to advance the calculation domain surface serving as an initial wavefront to the inside to generate a W-layer hierarchical grid, wherein W is larger than or equal to 1, and can be selected by a user according to the actual situation (the larger W is, the higher the precision of the space curve grid is, but the larger the calculation amount is, the larger W is, the value of W belongs to { 1,2,3, …,10 }).
And then generating initial linear tetrahedral units by adopting a three-dimensional constrained Delaunay (Delaunay) method by taking the current hierarchical grid surface as an initial constraint condition in the residual computation domain.
Further, when the position of the truncation point is calculated in the step 5, the proportional parameter λ is 1, which is the midpoint of the edge.
Further, if a curved grid with higher precision is to be constructed in step 5, a new truncation point is continuously generated on the basis of the first truncation of each boundary cell. And then projecting the respectively obtained new first-order truncation points to the boundary of the calculation domain to obtain corresponding projection points, 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.
Further, the basis function in the radial basis function in step 9 is a tight-branch radial basis function.
Further, the radial basis function in step 9 uses a quadratic polynomial.
Further, in step 11, the format output of the mesh data file adopts a CGNS format.
Further, in step 6, the non-boundary surface cutoff point is obtained by the same method as in step 5.
Further, in step 7, each edge is truncated by the same method in step 5 to obtain an adjacent cell truncation point.
Further, in step 7, for the multi-layer space curved grid, the neighboring cells of the neighboring cells are further retrieved in a loop.
The method for dividing the complex model three-dimensional space curved mesh is realized by cutting off the original linear tetrahedron to generate the curved surface tetrahedron unit, the mesh division requirement of the complex model can be well met under the condition of not improving the mesh magnitude, and the initial linear mesh is based on a front surface propelling method and a three-dimensional constrained Delaunay (Delaunay) method. In the process of generating the space curved grid point, only the boundary surface truncation point in the step 5 is selected as a control point, and the deformation point only selects the truncation points of the boundary unit and the adjacent units of the boundary unit, so that the number of the control point and the deformation point is greatly reduced, and the related calculation of the deformation point, namely the space curved grid point is greatly reduced. Compared with the prior art, step 8 innovatively provides a basis function selection mode based on the cosine value of the included angle; 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 deformation point calculation reduction and the simplification of the basis function form greatly improve the stability and the speed of grid division.
In summary, the invention provides a new control point and deformation point selection strategy, which reduces the selection of control points and deformation points while ensuring the grid quality, reduces a large amount of position calculation in the grid division process, optimizes the grid division strategy, reduces the time consumption of grid division and reduces the calculation power requirements. The mesh division requirement of the complex model is well met under the condition that the mesh magnitude is not improved.
Drawings
FIG. 1 is a schematic diagram of a three-dimensional curved mesh of an embodiment target model;
FIG. 2 is a plot of a blunt (10) cone linear/curved grid contrast;
FIG. 3 is a schematic cross-sectional (two-dimensional) view of an adjacent cell;
FIG. 4 is a schematic diagram of an X-51A linear grid;
FIG. 5 is a schematic diagram of an X-51A space curved grid;
FIG. 6 is a flow chart of the present invention.
Detailed Description
The technical solution of the present invention is described in detail below with reference to the accompanying drawings and examples.
Referring to the attached drawings, the method for dividing the three-dimensional curved grid of the complex model comprises the following steps:
And drawing a digital model based on the ACIS kernel according to the three-dimensional physical model (size parameters, formulas, drawings, real objects and the like), and reading the grid parameter settings such as grid size, boundary type and the like through a program.
And 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 establishing a Model class, importing the three-dimensional space geometric information in the digital Model into the Model class, and establishing a boundary marker (a member BCtype of the Model class in a program) according to the set boundary type.
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 toTAnd constructing a tetra body → face surface → edge → node point hierarchical structure by the mesh topological relation of the initial linear tetrahedron units to store mesh information.
Reading the grid parameters and the three-dimensional space geometric information in the Model class, firstly adopting a wavefront advancing method to advance the surface of the calculation region inwards to generate multi-layer (two layers in the spherical embodiment) grid units, and then utilizing a three-dimensional constrained Delaunay (Delaunay) method to generate grids of the residual region by taking the current wavefront as a constraint condition. And simultaneously constructing a tetra array class, a faceArray class, an edgeArray class and a nodeArray class to store the unit information.
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 traversalTAnd obtaining a total of psi boundary elements from the initial linear tetrahedral units, wherein ѱ E is { 1,2,3, …, psi }.
In this embodiment: and inserting the boundary mark BCtype into the tetraArray class, and establishing a corresponding relation. And then, the boundary surface, the boundary edge and the boundary point in the unit are marked in the insertion faceArray class, the edgeArray class and the nodeArray class of the BCtype (member BCtype of the tetraArray class in the program) in the tetraArray class.
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 tetrahedral unit comprises three boundary edges and a vertexP 1（x1,y1,z1）、P 2（x2,y2,z2) AndP 3（x3,y3,z3) X, y and z respectively represent coordinate values of three dimensions of points in a three-dimensional Cartesian rectangular coordinate system; the first edge corresponds to the vertexP 1AndP 2the second side corresponds toP 1AndP 3the third side corresponds toP 2AndP 3. And (3) cutting off the edges on all the boundary surfaces (the tetrahedral cell surface comprises three edges) in the psi boundary cells to generate the boundary surface cut-off points.
The position of the truncation point of the first boundary edge (the second and third boundary edges are the same) passes through the vertexP 1（x1,y1,z1) AndP 2（x2,y2,z2) And (3) calculating according to the calculation formula:
and the generated boundary surface truncation pointP H （xn,1,yn,1,zn,1）（HRepresenting boundary surface truncation point label) to the calculation domain boundary in the step 2 to obtain a projection pointP sub （xn,1,yn,1,zn,1）（subThe index represents the label of the projection point, the projection point and the cut-off point of the boundary surface have a one-to-one correspondence relationship), λ is a proportional parameter (in most cases, λ =1 can be directly taken as the midpoint of the edge for omitting the calculation step, the adjustment can be carried out based on the curvature of the physical model and the radial basis function in the special complex model with high curvature and the like, and the adjustment range is λ ∈ [0.5,2]。
Every time an interception point is generated, the interception point is cut offP H （xn,1,yn,1,zn,1) Storing cutoff point data classesP H （xN,1,yN,1,zN,1) And count (a)P H Where n represents the nth truncation point, and the maximum count value isN). All the projections are projected to the boundary surface at the same time to obtain projection pointsP sub （xn,1,yn,1,zn,1) Also storing projection point data classP sub （xN,1,yN,1,zN,1) And count (a)P sub Where n represents the nth projection point, and the maximum count value is alsoN）。
Further, if a higher precision curved grid is to be constructed, new truncation points continue to be generated on the basis of the first truncation of each border cell. By pointP H （xn,1,yn,1,zn,1) AndP 1（x1,y1,z1) And anP H （xn,1,yn,1,zn,1) AndP 2（x2,y2,z2) Respectively continue as end pointsAccording to
Generating an intercept point, at which (x)i,yi,zi) Is a point of one-off cutP H （xn,1,yn,1,zn,1），（xj,yj,zj) Are respectively a vertexP 1（x1,y1,z1) AndP 2（x2,y2,z2) The value of (c).
Then the cut-off points obtained respectivelyP H （xn,21,yn,21,zn,21) AndP H （xn,22,yn,22,zn,22) Projection onto the boundary of the computation domainP sub （xn,21,yn,21,zn,21) AndP sub （xn,22,yn,22,zn,22) And storing the corresponding interception point data class and the projection point data class. More high-order truncation points and high-order projection points are obtained in a circulating mode, and generally, each edge is truncated once to obtain one first-order truncation point or truncated twice to obtain three second-order truncation points, so that grid precision requirements can be met. Although more truncation points can be circularly intercepted to improve the grid precision, the more truncation times bring larger calculation amount, the slower calculation speed and the larger calculation force requirement. It is therefore not recommended to use a calculation scheme above the second order truncation point in the meshing.
The embodiment is as follows: and retrieving the boundary mark BCtype in the tetra array class, and retrieving the edge array class and the nodeArray class corresponding to the unit when the unit is the boundary unit. And further searching the boundary mark BCtype of the edge, and solving an interception point according to the point coordinate and the formula corresponding to the boundary edge when the boundary edge is searchedP H （xn,1,yn,1,zn,1) And a projection pointP sub （xn,1,yn,1,zn,1). Determining the point of interruptionP H （xn,1,yn,1,zn,1) And a projection pointP sub （xn,1,yn,1,zn,1) Post-entry truncation point data classP H （x N,1,y N,1,z N,1) AndP sub （x N,1,y N,1,z N,1) In (c), and count (n is the maximum count value ofN）。
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; the same method in the step 5 is adopted to obtain the non-boundary surface truncation pointP inr （xm,1,ym,1,zm,1) Storing cutoff point data classesP inr （xM,1,yM,1,zM,1) And then the number of the counting is counted,inrrepresents a non-boundary surface truncation point index,P inr the maximum count value of M is recorded as M.
And 7, retrieving adjacent units of the psi boundary units one by one according to the psi boundary unit information obtained in the step 5, extracting the information of surfaces, edges and points of non-adjacent surfaces in the adjacent units, and further circularly retrieving the adjacent units of the adjacent units if a multi-layer space curved grid is needed.
The adjacent unit: two of the tetrahedral cells sharing a face are defined as adjacent cells, four of the adjacent cells are non-boundary cells, and fewer than four but at least one of the adjacent cells are boundary cells.
The same method in the step 5 is adopted to truncate each edge to obtain the truncation point of the adjacent unitP adj （xk,1,yk,1,zk,1) Storing cutoff point data classesP adj （xK,1,yK,1,zK,1) In the step (2), counting,P adj the maximum count value of the medium K is K,adjrepresenting adjacent cell truncation point labels.
In this embodiment, class members tetra [ a, B ] of the faceArray class are retrieved, and a and B are two unit global numbers to which the faces belong. And searching adjacent units in the tetraArray class according to the global number A or B, and taking out the information of the edges and points of the non-adjacent surfaces.
Step 8, traversing boundary surface truncation point data classesP H （xN,1,yN,1,zN,1) And projection point data classP sub （xN,1,yN,1,zN,1) Calculating the cut-off point of the boundary surfaceP H （xn,1,yn,1,zn,1) Corresponding projection point thereofP sub （xn,1,yn,1,zn,1) The cosine of the angle between them is denoted as cos (θ) H-sub And storing the corresponding position of each point in the nodeArray class.
Wherein xn,H,yn,H,zn,HIs as followsnA three-dimensional coordinate value of a boundary surface truncation point representing x subn,,y subn,,z subn,Is as followsnThree-dimensional coordinate values of the respective projected points.
Step 9, all the truncation points obtained by truncation in the steps 5, 6 and 7 are substituted into the radial basis function calculation formula, and the calculation is obtained through the radial basis function calculation formulaf(x,y,z) The value of (1) is the deformation of the truncation point; store in temporary variable mesh-discrete current (total number)N+M+K) In (1). The radial basis function calculation formula is:
whereinαAs a function of the number of the coefficients,
in order to be the weight, the weight is,I N is the total number of control points andI N =N；
in order to be a radial basis function,
representing the Euclidean distances from all the truncation points to the corresponding control points
For the cut-off point of the boundary surface in step 5P H （xN,1,yN,1,zN,1) I.e. a subset of all boundary surface truncation points;
in the embodiment, based on the maximum included angle cosine value of the truncation point (i.e. the spatial curved grid point), the radial basis function calculation formula used in the present invention is:
the radial basis function only uses quadratic polynomial under most conditions, and the calculation required by the deformation quantity can be greatly reduced by using the low-order polynomial basis function with a simple form; and only in the complex situations of high curvature and the like of the complex model, the high-order polynomial is adopted to meet the grid quality requirement.
The weight is calculated by the control point (boundary surface interception point) and the projection point
When calculating the weight
InXRepresentative boundary surface truncation point (namely boundary surface truncation point in step 5)P H （xN,1,yN,1,zN,1) So as to total numberI N =N(ii) a When calculating the weightf (x,y,z) Is a known quantity, represents the displacement quantity of the control point,f (x,y,z) The minimum value of the Euclidean distance from the truncation point to the projection point on the same boundary surface is taken to avoid the condition of trivial solution(ii) occurs;
f (x,y,z)=min{║X b -X sub ║}
in the case of a first-order truncation pointX b B belongs to { 1,2,3 } and is three truncation points corresponding to three boundary edges on the same boundary surface,X sub the coordinate value of the projection point corresponding to the truncation 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;
finding weights
Post-substitution intercept point calculationf (x,y,z) At this time
In order to be of a known quantity,f (x,y,z) Is the amount to be calculated; when calculating the displacement of the grid point
InX=X j ，X j Representing all the truncation points obtained by the truncation points of the boundary surface in the step 5, the non-boundary surface in the step 6 and the adjacent unit in the step 7, and the total number of the truncation pointsI SUM =N+M+K；
Step 10, the deformation quantity calculated in step 9f (x,y,z) Substituted truncation point coordinates (x)J,1,yJ,1,zJ,1) To calculate the spatial curved grid coordinates (x)c,yc,zc) The calculation formula is as follows:
（xc,yc,zc）=（xJ,1,yJ,1,zJ,1）+f (x,y,z),J∈{N},{M},{K}
and traversing the temporary variable mesh-discrete to take out the deformation amount corresponding to each truncation point, and adding the deformation amount to the initial position of the truncation point to obtain the position of the transformed truncation point (the position of the added truncation point is the spatial curved grid point), namely the spatial curved grid point coordinate. And then, storing the obtained space curved grid point coordinates into the interception point class to replace the original interception point coordinates.
Step 11, calculating the high-order grid points (x) in step 10c,yc,zc) As new grid nodes, namely, the grid point coordinates of the space curve calculated in the step 10 are read from the interception point class, and are inserted into grid data as new grid nodes 5 to node10, and the grid topological relation is updated; and then outputting the new spatial curved grid data by using a corresponding grid data file format (such as CGNS format), and importing the output grid data file into visualization.
The present embodiment takes 3 target objects as an example to verify the validity of the present invention; the spatial curved meshing results of the final meshing of the embodiments are shown in the right diagrams of fig. 1 and 2, and fig. 5.
FIG. 1 is a global surface view of a left graph of a space curved grid with a sphere as a target object, a middle graph of a space curved grid sectional view, and a right graph of a local enlarged view.
FIG. 2 is a graph of the results of a blunt (10) cone spatial curved grid. The left graph is the initial linear tetrahedral mesh obtained in the third step, and the right graph is the schematic diagram of the finished space curved mesh.
Fig. 4 and 5 are graphs of results for an X-51A aircraft. Fig. 4 is the initial linear tetrahedral mesh obtained in step three, and fig. 5 is a schematic diagram of the completed spatial curved mesh.
From the comparison between the left diagram and the right diagram in fig. 2 and the comparison between fig. 4 and fig. 5, it can be seen that the space curved grid division method provided by the present invention realizes high-precision space curved grids under the condition of ensuring that the grid quality and the basic topology structure are not changed. The existing curved grid technology only basically realizes the object plane curved grid, although the object plane curved grid improves the calculation precision at the boundary, the numerical calculation precision is often limited by the precision of the internal grid. The invention adopts a brand-new strategy for selecting control points and deformation points (namely spatial curved grid points) to expand curved grids into spatial inner layer grids to realize high-precision spatial curved grid division of a complex model. In conclusion, the invention is practical and effective.
In conclusion, the invention provides a new control point and deformation point selection strategy, and well meets the mesh division requirement of the complex model under the condition of not improving the mesh magnitude. Compared with the prior art: 1. the invention realizes the space curved mesh division, and the prior art is basically an object plane curved mesh; 2. by adopting a new control point and deformation point selection strategy and a calculation method, the grid division requirement of the complex model is well met under the condition of not improving the grid magnitude; 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 grid quality, reduces a large amount of position calculation in the grid division process, optimizes the grid division strategy, reduces the time consumption of grid division and reduces the calculation force requirement. The deformation point calculation reduction and the simplification of the basis function form greatly improve the stability and the speed of grid division.
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.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111637505.5A CN114510775B (en) | 2021-12-30 | 2021-12-30 | Method for dividing three-dimensional space curved grids of complex model |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111637505.5A CN114510775B (en) | 2021-12-30 | 2021-12-30 | Method for dividing three-dimensional space curved grids of complex model |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114510775A true CN114510775A (en) | 2022-05-17 |
| CN114510775B CN114510775B (en) | 2023-06-27 |
Family
ID=81548505
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202111637505.5A Active CN114510775B (en) | 2021-12-30 | 2021-12-30 | Method for dividing three-dimensional space curved grids of complex model |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114510775B (en) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115186074A (en) * | 2022-06-24 | 2022-10-14 | 深圳市规划和自然资源数据管理中心 | Meta analysis-based method for simulating spatial distribution pattern of pH value of soil |
| CN116401916A (en) * | 2023-03-20 | 2023-07-07 | 北京云境智仿信息技术有限公司 | Method, device, medium and equipment for generating high-quality three-dimensional grid |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20130300735A1 (en) * | 2012-05-14 | 2013-11-14 | Autodesk, Inc. | Mesh boundary smoothing |
| CN109636913A (en) * | 2018-12-04 | 2019-04-16 | 山东理工大学 | Triangle gridding increment topology joining method based on Delaunay subdivision |
| CN112015735A (en) * | 2020-08-20 | 2020-12-01 | 西安数峰信息科技有限责任公司 | Data storage structure and data storage method of unstructured grid |
-
2021
- 2021-12-30 CN CN202111637505.5A patent/CN114510775B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20130300735A1 (en) * | 2012-05-14 | 2013-11-14 | Autodesk, Inc. | Mesh boundary smoothing |
| CN109636913A (en) * | 2018-12-04 | 2019-04-16 | 山东理工大学 | Triangle gridding increment topology joining method based on Delaunay subdivision |
| CN112015735A (en) * | 2020-08-20 | 2020-12-01 | 西安数峰信息科技有限责任公司 | Data storage structure and data storage method of unstructured grid |
Non-Patent Citations (1)
| Title |
|---|
| 骆冠勇 等: "用逐点插入法生成Delaunay四面体自适应网络", 计算力学学报, vol. 24, no. 6, pages 917 - 922 * |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115186074A (en) * | 2022-06-24 | 2022-10-14 | 深圳市规划和自然资源数据管理中心 | Meta analysis-based method for simulating spatial distribution pattern of pH value of soil |
| CN116401916A (en) * | 2023-03-20 | 2023-07-07 | 北京云境智仿信息技术有限公司 | Method, device, medium and equipment for generating high-quality three-dimensional grid |
| CN116401916B (en) * | 2023-03-20 | 2024-01-26 | 北京云境智仿信息技术有限公司 | Method, device, medium and equipment for generating high-quality three-dimensional grid |
Also Published As
| Publication number | Publication date |
|---|---|
| CN114510775B (en) | 2023-06-27 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN100468418C (en) | Method and program for generating volume data from boundary representation data | |
| De Floriani et al. | A survey on data structures for level-of-detail models | |
| CN102509339B (en) | Method for simplifying vertex clustering of three-dimensional models with texture constraint | |
| CN114510775A (en) | Complex model three-dimensional space curved grid division method | |
| WO2021203711A1 (en) | Isogeometric analysis method employing geometric reconstruction model | |
| Anderson et al. | Smooth, volume-accurate material interface reconstruction | |
| Cascón et al. | Comparison of the meccano method with standard mesh generation techniques | |
| Yuan et al. | Feature preserving multiresolution subdivision and simplification of point clouds: A conformal geometric algebra approach | |
| CN111695281B (en) | Particle rapid positioning method for tetrahedral mesh division finite element particle simulation | |
| Liu et al. | Memory-efficient modeling and slicing of large-scale adaptive lattice structures | |
| CN113610983A (en) | Automatic subdivision method for discrete point space curved surface triangular mesh | |
| CN114492250A (en) | Curved surface mesh generation method and system based on recursive decomposition and computer equipment | |
| Natarajan et al. | Simplification of three-dimensional density maps | |
| Marinić-Kragić et al. | Superimposed RBF and B-spline parametric surface for reverse engineering applications | |
| Qin et al. | Surface reconstruction from unorganized point clouds based on edge growing | |
| Langbein et al. | An efficient point location method for visualization in large unstructured grids. | |
| Makhlouf et al. | Reconstruction of the cad model using tps surface | |
| Kim et al. | Efficient encoding and decoding extended geocodes for massive point cloud data | |
| Khamayseh et al. | Use of the spatial kD-tree in computational physics applications | |
| CN113342999A (en) | Variable-resolution-ratio point cloud simplification method based on multi-layer skip sequence tree structure | |
| Hou et al. | Knot optimization for biharmonic b-splines on manifold triangle meshes | |
| Cascón et al. | The meccano method for automatic tetrahedral mesh generation of complex genus-zero solids | |
| Hao et al. | Surface reconstruction based on CAD model driven priori templates | |
| Brandts et al. | A geometric toolbox for tetrahedral finite element partitions | |
| CN111047687A (en) | Three-dimensional T-spline-based heterogeneous material solid modeling method |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |