A feature-preserving quadrilateral-dominant mesh generation method based on geometrically constrained frame field guidance

Abstract

A feature-preserving quadrilateral-dominant mesh generation method based on geometric constraint frame field guidance, comprising: obtaining a subdivision triangular mesh; initializing a frame field direction of each vertex on the subdivision triangular mesh to obtain an initial frame field direction; constructing a multi-resolution hierarchy based on the initial frame field direction; iteratively optimizing each layer frame field of the multi-resolution hierarchy by constructing a neighborhood consistency update rule and a direction constraint term to generate a final direction field; establishing a multi-level boundary fixed point constraint mechanism according to the multi-resolution hierarchy; randomly initializing a position field according to the multi-resolution hierarchy, and then iteratively optimizing to obtain a position field satisfying geometric constraints; adaptively partitioning the subdivision triangular mesh based on vertex principal curvature directions and feature boundary connectivity, so that the partition boundary is consistent with the geometric feature, and a structured region division result is generated; and generating a quadrilateral-dominant mesh by using a two-stage processing framework according to the final direction field, the position field and the structured region division result.

Description

Technical Field

[0001] This invention relates to the field of computer-aided engineering analysis and numerical simulation technology, specifically to a method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field. Background Technology

[0002] In numerical simulations such as finite element analysis, computational fluid dynamics, and structural simulation, it is often necessary to discretize continuous geometric models into mesh models. The element shape, size distribution, and alignment with geometric features of the mesh directly affect the accuracy, stability, and computational efficiency of the numerical calculations.

[0003] Triangular meshes are widely used due to their high degree of automation in generation and strong adaptability to complex geometries. However, in regions with significant curvature variations or geometric features such as sharp edges and polygonal lines, triangular elements are prone to shape distortion. Furthermore, to ensure computational accuracy, it is often necessary to significantly increase the mesh density, leading to increased computational costs. In contrast, quadrilateral elements are easier to align along the principal directions or boundary directions of the surface, exhibiting better numerical convergence and computational efficiency in structured or semi-structured mesh construction. Therefore, the demand for high-quality quadrilateral-dominated meshes is increasing in engineering simulations.

[0004] Among existing quadrilateral mesh generation methods, a common approach is to first construct a frame field (e.g., frame field / quadruple rotational symmetry field) on a triangular mesh, then generate a position field based on the frame field and extract the quadrilateral mesh. While this method can achieve mesh alignment with curvature or boundary directions to some extent, it still has the following shortcomings in practical applications:

[0005] (1) Near sharp edges, broken lines and feature lines, the frame field optimization process is easily affected by the global smoothing term and the direction drift occurs, which weakens or even disappears the geometric features;

[0006] (2) During the construction of the location field, when the boundary region lacks stable anchor point constraints, problems such as uneven distribution of boundary vertices and unstable corner positions are likely to occur, which will lead to a decrease in the geometric fidelity of the boundary.

[0007] (3) For complex topological regions such as holes, thin sheets or local flatness, existing methods are prone to local mesh collapse or extraction failure due to vertices being too close, which affects the quality and robustness of quadrilateral meshes. Summary of the Invention

[0008] To address the shortcomings of existing quadrilateral dominant mesh generation methods based on frame fields in preserving geometric features, this invention provides a feature-preserving quadrilateral dominant mesh generation method guided by geometrically constrained frame fields. By improving the boundary frame field initialization strategy, introducing boundary fixed points and a multi-resolution constraint propagation mechanism, and combining curvature-guided partitioning and two-stage mesh extraction, this method improves the regularity and geometric quality of the generated mesh while ensuring boundary stability and feature preservation.

[0009] To achieve the above objectives, the present invention adopts the following technical solution:

[0010] A method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field includes the following steps:

[0011] S1. Obtain and subdivide the preprocessed triangular surface mesh model to generate a subdivided triangular mesh;

[0012] S2. Initialize the frame field direction of each vertex on the subdivided triangular mesh to obtain the initial frame field direction;

[0013] S3. Based on the initial frame field orientation, construct a multi-resolution hierarchical structure; iteratively optimize each layer of the multi-resolution hierarchical structure by constructing neighborhood consistency update rules and orientation constraint terms to generate the final orientation field;

[0014] S4. Establish a multi-level boundary fixed-point constraint mechanism based on the multi-resolution hierarchical structure;

[0015] S5. Randomly initialize the position field according to the multi-resolution hierarchical structure, and then obtain the position field that satisfies the geometric constraints through iterative optimization;

[0016] S6. Adaptively partition the subdivided triangular mesh based on the connectivity between the principal curvature direction of the vertices and the feature boundary, so that the partition boundary is consistent with the geometric features, and generate structured region partitioning results.

[0017] S7. Based on the final orientation field, position field, and structured region division results, a two-stage processing framework is used to generate a quadrilateral dominant mesh.

[0018] As a preferred option, in S2:

[0019] For an internal vertex, a unit tangent vector is randomly generated in the tangent plane of the vertex as the initial direction;

[0020] For boundary vertices, dihedral angle-related directions are constructed based on the normal vectors of adjacent facets to reflect crease features; boundary polygonal direction is constructed based on adjacent boundary edges to reflect boundary orientation; the fusion weight is adaptively calculated according to the size of the dihedral angle, and the dihedral angle-related directions and boundary polygonal direction are weighted and fused, then projected onto the tangent plane and normalized to obtain the initial direction of the boundary vertex.

[0021] Preferably, S3 includes:

[0022] S31. A multi-resolution hierarchical structure is constructed using a bottom-up graph downsampling method. The product of the normal similarity of adjacent vertices and the dual area ratio is used as the merging priority. Vertices are merged layer by layer to generate multiple hierarchical meshes, and the mapping relationship between vertices between layers is established.

[0023] S32. Iteratively optimize the frame fields of each layer of the multi-resolution hierarchical structure in sequence. In each iteration, each vertex updates its own orientation according to the frame field orientation of its adjacent vertices, and applies orientation constraints to the boundary vertices.

[0024] Preferably, S32 includes:

[0025] S321. For each layer of mesh, each vertex is weighted and averaged according to the frame field direction of its adjacent vertices, and its own direction is updated smoothly in the neighborhood. The fourfold rotational symmetry of the frame field is handled by a compatibility function to eliminate the inconsistency of direction caused by rotational symmetry.

[0026] S322. Apply directional constraints to boundary vertices or feature line vertices; fuse the preset constraint directions with the neighborhood smoothing results to maintain the boundary and feature orientations.

[0027] S323. On each mesh layer, minimize the objective function that includes the neighborhood consistency term and the orientation constraint term, and iteratively optimize the frame field.

[0028] Preferably, S4 includes:

[0029] S41. Obtain the boundary edge set and corner point set and perform boundary path tracing to generate boundary path segments; perform segmented sampling on the boundary path segments based on the target mesh scale to obtain the fine-layer fixed point set.

[0030] S42. Implement a layer-by-layer sparse sampling strategy for the set of fixed points in the fine layer and establish inter-layer mapping relationships to realize the propagation of fixed point constraints in the multi-resolution hierarchical structure.

[0031] S43. For each fixed point, search for grid vertices that can serve as fixed points within its neighborhood and assign fixed point constraints to the matching candidate vertices. Fixed point constraints are divided into two categories: corner fixed points are treated as hard constraints, and their positions remain unchanged during optimization; non-corner fixed points are treated as soft constraints, and their deviation is limited by introducing a weighted penalty term into the position field optimization objective function.

[0032] As a preferred embodiment, S5 includes:

[0033] S51. Based on the multi-resolution hierarchical structure, establish a local grid coordinate system for each vertex and randomly initialize the position field in the tangent plane;

[0034] S52. Perform position field iterative optimization layer by layer, and use the optimization result of each layer as the initial value of the next fine layer to propagate layer by layer.

[0035] Preferably, S52 includes:

[0036] For each vertex, position field suggestions based on position compatibility are collected from its adjacent vertices, and the target position of the vertex is obtained by weighted fusion; the target position is corrected by combining the constraint position of the vertex's geometric constraints and weights, and then reprojected onto the tangent plane;

[0037] The continuous target positions are projected onto the local coordinate system, and the coordinates of discrete grid points are determined by rounding to update the position field.

[0038] For vertices with geometric constraints, a fixed-point penalty term and an orientation alignment constraint term are introduced into the objective function to maintain feature consistency.

[0039] As a preferred embodiment, S6 includes:

[0040] S61. For each vertex, construct a discrete curvature tensor matrix in its neighborhood; obtain the principal curvature direction and principal curvature magnitude through eigenvalue decomposition based on the discrete curvature tensor matrix.

[0041] S62. Determine the characteristic edge based on the dihedral angle of adjacent facets;

[0042] S63. Construct a curvature alignment cost function and calculate path quality indicators for candidate paths; select partitioned guidance paths based on path quality indicators.

[0043] S64. Starting from the seed face, perform connectivity expansion along adjacent faces that satisfy the unassigned and do not cross the feature edge constraints, assign partition identifiers to the faces, and establish mapping relationships from faces to partitions and from vertices to partitions.

[0044] S65. Post-process and optimize the initial partitioning results to obtain the structured region partitioning results.

[0045] As a preferred embodiment, the two-stage processing framework in S7 includes:

[0046] The first stage involves graph structure extraction and vertex merging. Compatibility analysis of frame field and position field is performed on adjacent vertices to identify vertex pairs that need to be merged, connected, or ignored. Reasonable vertex connection relationships are established, and quality optimization processing is performed on flat regions.

[0047] The second stage involves face extraction and mesh generation. A local coordinate system is established for the vertices and adjacent edges are sorted by angle. A ring tracing strategy is used to extract the target face. For cases where tracing fails, a backtracking search strategy is used to select the optimal tracing starting point.

[0048] As a preferred option, the first stage of processing includes:

[0049] For any pair of adjacent vertices in the original triangular mesh Perform frame field compatibility processing; calculate the discrete position compatibility relationship of adjacent vertex pairs based on the position field and target scale parameters, and output the corresponding index offset; manage vertex pairs that need to be merged using a disjoint-set data structure, and calculate the merging weight based on the position field error to determine the processing priority; calculate the geometric attributes of each representative vertex after merging; and perform quality constraint processing on the graph structure based on the structured region partitioning results.

[0050] The second phase includes:

[0051] A local coordinate system is established for each vertex, and its adjacent edges are arranged in order of local angles. A ring tracing strategy is used to generate a facet of the target size. When the local topology is complex and the tracing fails, a backtracking search is used to switch between candidate branches. For faces that are too large and have more vertices than the target size, a centralization or fan-shaped partitioning strategy is used to split them into several quadrilateral faces. The quadrilateral faces are checked for orientation consistency and the geometric quality is optimized using Laplacian smoothing. Topology consistency is then checked.

[0052] Compared with the prior art, the beneficial effects of the present invention are reflected in:

[0053] 1. Unlike traditional techniques that initialize the frame field direction for all vertices in the same way, this invention uses a technique that initializes the frame field direction for boundary vertices separately based on the size of the dihedral angle. This ensures that the direction field at the boundary and sharp edges is aligned with the geometric features from the beginning, reducing the possibility of excessive smoothing of feature directions during subsequent optimization.

[0054] 2. Unlike traditional techniques where boundary vertices lack fixed constraints during position field optimization, this invention employs a technique of sampling fixed points along the boundary path and propagating constraints step-by-step across multiple resolution layers. This results in more stable and uniformly distributed boundary vertices during the optimization process, making it less prone to boundary point shifts.

[0055] 3. Unlike traditional techniques that indiscriminately perform global mesh processing on complex regions such as holes and thin sheets, this invention adopts an adaptive partitioning technique based on the connectivity constraints of the principal curvature direction and feature boundaries. This ensures that the partition boundaries are consistent with the geometric features, thereby providing a structured basis for subsequent mesh extraction and reducing the difficulty of processing complex topological regions.

[0056] 4. Unlike traditional techniques that lack geometric partition constraints during mesh extraction, leading to improper merging across regions and a single tracking strategy that is prone to failure in topologically complex regions, this invention introduces partition constraints in the vertex merging stage to prohibit improper merging operations across partitions. In the face extraction stage, it adopts a strategy of combining ring tracking with backtracking search to select the optimal tracking starting point, thereby improving the tracking success rate and final mesh quality while ensuring topological correctness. Attached Figure Description

[0057] Figure 1 This is a flowchart of the method in Embodiment 1 of the present invention;

[0058] Figure 2 This is a schematic diagram of the subdivision preprocessing in Embodiment 1 of the present invention;

[0059] Figure 3 This is a schematic diagram of boundary marker field initialization and marker field generation in Embodiment 1 of the present invention;

[0060] Figure 4 A schematic diagram of the position field generation after setting the boundary fixed point in Embodiment 1 of the present invention;

[0061] Figure 5 This is a schematic diagram of the position field generation without the addition of boundary fixed points in Embodiment 1 of the present invention;

[0062] Figure 6 This is a diagram illustrating the adaptive mesh partitioning effect based on principal curvature guidance in Embodiment 1 of the present invention.

[0063] Figure 7 This is a schematic diagram of the quadrilateral mesh extraction image in Embodiment 1 of the present invention;

[0064] Figure 8 This is a comparison diagram of the geometric model mesh generation effect in Embodiment 1 of the present invention;

[0065] Figure 9 A magnified view of the local details of the feature retention effect in Embodiment 1 of the present invention. Figure 1 ;

[0066] Figure 10 A magnified view of the local details of the feature retention effect in Embodiment 1 of the present invention. Figure 2 . Detailed Implementation

[0067] To make the technical means, inventive features, objectives, and effects of the invention readily understandable, the invention is further described below with reference to specific illustrations. However, the invention is not limited to the embodiments described below.

[0068] It should be noted that the structures, proportions, sizes, etc., illustrated in the accompanying drawings of this specification are only used to complement the content disclosed in the specification for those skilled in the art to understand and read, and are not intended to limit the conditions under which the present invention can be implemented. Therefore, they have no substantial technical significance. Any modifications to the structure, changes in the proportions, or adjustments to the size, without affecting the effects and objectives that the present invention can produce, should still fall within the scope of the technical content disclosed in the present invention.

[0069] Existing technologies typically suffer from the following problems: In sharp edges, feature lines, and boundary regions, the frame field is prone to deviating from the boundary direction and being weakened by the global smoothing term during iterative optimization, resulting in the blunting or even disappearance of sharp features; In the location field construction stage, there is a lack of effective control over the distribution of boundary points, which can easily lead to poor boundary geometric fidelity and unstable corner positions; In addition, in complex topological surfaces such as hole regions and flat regions, mesh collapse or cell quality degradation can easily occur due to vertices being too close together.

[0070] Example 1:

[0071] like Figure 1 The method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field is shown below. The specific steps are as follows:

[0072] S1. Obtain and subdivide the preprocessed triangular surface mesh model to generate a subdivided triangular mesh;

[0073] like Figure 2 As shown, the input triangular surface mesh model is subdivided and preprocessed to improve the mesh uniformity and element quality, and to provide a stable topological and geometric basis for the subsequent calculation of the frame field and position field.

[0074] Let the input triangular surface mesh be... ,in For the set of vertices, Let be the set of edges. This is a set of triangular facets. First, a directed edge or half-edge data structure is constructed to obtain vertex adjacency relationships and adjacent face information of the edges. Then, the mesh feature scale parameters are determined based on the mesh surface area and the number of target vertices. Then iterate through the edge set. Calculate the side length for each side. The dihedral angle is obtained by the angle between the normal vectors of the two adjacent faces and the dihedral angle. When the side length satisfies Or the dihedral angle satisfies When the edge is subdivided, it is marked as a subdivision edge and a subdivision operation is performed. This involves inserting a new vertex at the midpoint of the edge and splitting and updating the relevant triangular facets to obtain a subdivided triangular mesh. Preferably, the subdivision threshold can be set in association with the target scale parameter; the process stops when the maximum number of iterations is reached or no edge satisfies the subdivision condition.

[0075] During the subdivision process, the boundary edges and feature edges in the original mesh are inherited in terms of attributes. If the original edge is a boundary edge or a feature edge, the new edge generated by the subdivision will retain the same label to avoid the loss of boundary or sharp features due to subdivision. After the subdivision is completed, the vertex normal vectors and adjacency relationships are updated to provide the correct tangent plane information for subsequent steps.

[0076] S2. Initialize the frame field direction of each vertex on the subdivided triangular mesh to obtain the initial frame field direction;

[0077] The orientation of the frame field is the unit tangent vector located in the vertex tangent plane. For internal vertices, the orientation is initialized randomly or quasi-randomly to reduce the risk of getting trapped in local optima in subsequent optimizations. For boundary vertices, the initial orientation is constructed based on geometric features, combining the dihedral-related orientation with the boundary polygonal direction, and adaptively determining the fusion weight according to the size of the dihedral angle. This makes the crease feature orientation more emphasized in sharp-edge regions and the boundary direction continuity more emphasized in smooth-boundary regions, thereby improving the feature preservation stability of subsequent frame field optimizations.

[0078] For internal vertices Select a set of orthogonal unit bases in its tangent plane. , And generate the initial direction on this basis. In this embodiment, a random angle is generated. And order:

[0079]

[0080] For boundary vertices First, candidate directions related to dihedral angles and boundary polygons are constructed, and then normalized after being projected onto the vertex tangent plane. Let the unit normal vectors of the two faces adjacent to the vertex be... The relevant directions of the dihedral angle are constructed from it. Let the direction vectors of the two boundary edges adjacent to this vertex on the boundary polyline be... The direction of the broken line constructed from it. To ensure that the obtained direction lies within the tangent plane, tangent plane projection and normalization are performed separately:

[0081] ,

[0082] in For boundary vertices The unit normal vector.

[0083] In this embodiment, it is possible to take

[0084]

[0085] Preferably, when When the modulus is lower than the preset threshold, you can select or As an alternative boundary direction.

[0086] dihedral It can be obtained from the dot product of the normal vectors of adjacent faces:

[0087] Fusion weights The dihedral angle adapts to change, controlling the fusion ratio between the dihedral angle direction and the boundary line direction. In this embodiment, the following can be selected:

[0088]

[0089] in The preset acute angle threshold parameter is used; in this embodiment, For example, take The initial frame field direction at the final boundary vertex is defined as:

[0090]

[0091] S3. Based on the initial frame field orientation, construct a multi-resolution hierarchical structure; iteratively optimize each layer of the multi-resolution hierarchical structure by constructing neighborhood consistency update rules and orientation constraint terms to generate the final orientation field;

[0092] Based on the initial frame field orientation obtained in step S2, a multi-resolution hierarchical structure is constructed, and the frame field orientation of adjacent vertices is gradually made consistent through iterative optimization on each layer of the mesh. At the same time, directional constraints are applied to the boundary vertices to maintain geometric features, and finally a globally continuous frame field aligned with the boundary orientation is generated.

[0093] After constructing the multi-resolution hierarchical structure, the corresponding frame field representation is generated and optimized on each layer of the mesh. This step, through a multi-resolution iterative optimization strategy, ensures that the frame field remains continuous and smooth throughout the entire region while maintaining boundary and feature direction constraints.

[0094] S31. The multi-resolution hierarchical structure is constructed using a bottom-up graph downsampling method. The original subdivided triangular mesh is used as the finest layer. The product of normal similarity and dual area ratio is used as the merging priority for adjacent vertex pairs. Vertex pairs are merged sequentially from high to low priority, shrinking two vertices into a new vertex each time, thereby gradually reducing the number of vertices and generating a coarser mesh representation. The vertex merging priority is defined as:

[0095]

[0096] Where similarity is the dot product of the normal vectors and ratio is the area ratio. Merging is only considered when similarity is greater than a preset threshold; in this embodiment, the threshold is set to 0.5.

[0097] The normal vector and area of ​​the new coarse-layer vertices generated by merging are obtained by weighted averaging of the areas of the merged vertices. The frame orientation of the new vertex needs to be projected onto the tangent plane of the new vertex and normalized after weighted averaging to ensure that the result is a valid tangent vector. Merging stops when the number of vertices in a certain layer drops to about 1 / 3 of the previous layer or falls below a preset lower threshold.

[0098] In this invention, layers with a higher number of vertices and higher mesh resolution are called finer layers, while layers with a lower number of vertices and lower mesh resolution are called coarser layers. The original subdivided triangular mesh is the finest layer, and the mesh obtained when merging stops is the coarsest layer. Approximately [amount missing] are generated between the two. The layers are represented by a grid from fine to coarse, where This represents the original number of grid vertices.

[0099] The standard field optimization starts from the coarsest layer. Coarse layers have fewer vertices and a lower solution space dimension, making iterative optimization more likely to converge to a globally consistent solution. The convergence result of the coarse layer serves as the initial value for the optimization of adjacent fine layers, propagating layer by layer to the finer layers. This ensures that each fine layer only needs local refinement based on a good initial value to converge. Simultaneously, a complete inter-layer mapping relationship is established, recording the set of fine layer vertices corresponding to each coarse layer vertex for use in the coarse-to-fine propagation.

[0100] S32. Iteratively optimize the frame field of each layer in the multi-resolution hierarchical structure. In each iteration, each vertex updates its own orientation based on the frame field orientation of its adjacent vertices. Simultaneously, orientation constraints are applied to boundary vertices to ensure the frame field becomes smooth and consistent overall and remains aligned with geometric features at the boundaries. S321. During the frame field optimization process, iterative updates are performed on each mesh layer from coarse to fine. For vertices… According to its set of adjacent vertices The orientation of the frame field in the image is taken as the new orientation of the frame field at that vertex after weighted averaging.

[0101]

[0102]

[0103] in As vertex The current direction of the standard frame field, Its normal vector, For adjacency weight, This is the frame field compatibility function. Since the frame field has four-fold rotational symmetry, there are multiple equivalent rotational configurations for the frame field directions of adjacent vertices when comparing them. The compatibility function selects the one with the best alignment and then participates in the weighted average, thereby eliminating the inconsistency in direction caused by rotational symmetry. In this embodiment, four rotational configurations can be enumerated and the number of rotations required to rotate adjacent frame fields can be selected. S322: To maintain the consistency of direction at geometric boundaries and feature lines, a direction constraint term is introduced during the frame field optimization process. For boundary vertices or feature line vertices, their constraint directions are preset. And constraints are imposed by fusing the results with the neighborhood average:

[0104]

[0105] in The constraint weights are used to balance the effects of neighborhood smoothing and orientation constraints. In this embodiment, the boundary vertices and feature line vertices are selected as... To strictly maintain the boundary direction; smaller weights are applied to internal vertices to enhance overall smoothness.

[0106] The frame field optimization process can be equivalently understood as minimizing an objective function that includes neighborhood consistency terms and orientation constraint terms at each mesh layer.

[0107]

[0108] in This represents a neighborhood consistency term, used to reduce the difference in field orientation between adjacent vertices. This represents the boundary direction constraint term, used to keep the boundary and feature directions stable.

[0109] S323. Iterative optimization ends when a preset termination condition is met. In this embodiment, the maximum number of iterations is reached. Or the change in the field direction of the standard frame in two consecutive iterations is less than the threshold. When the target field of that layer is considered to have converged, the direction change can be defined as the maximum or average of the angle between the directions of the vertex target field in two adjacent iterations; threshold The settings can be adaptively adjusted based on grid resolution or target precision, for example, setting it to 10. -3 ~10 -2In radians. After the current layer converges, the optimization result is used as the initial value for the next finer layer, and the process continues to advance to finer layers until the finest layer is optimized.

[0110] like Figure 3 As shown, the green lines represent boundary edges, and the red crosses represent the four-directional frame field of each vertex. In the boundary region of the triangular mesh, the frame field direction is initialized along the geometric boundary to ensure that the frame field is consistent with the boundary orientation; in the internal region, the frame field is optimized and continuously distributed between adjacent triangles, thus providing a stable directional reference for subsequent position field generation and mesh extraction.

[0111] S4. Based on the multi-resolution hierarchical structure constructed in step S3, establish a multi-level boundary fixed-point constraint mechanism. In the position field optimization process, this step maintains the boundary features and corner positions by setting fixed points at the boundary and using multi-resolution hierarchical propagation of fixed points, thus avoiding boundary drift or incorrect alignment across features.

[0112] S41. Obtain the boundary edge set and corner point set and perform boundary path tracing to generate boundary path segments; perform segmented sampling on the boundary path segments based on the target mesh scale to obtain the fine-layer fixed point set.

[0113] First, boundary path tracing is performed based on the input set of boundary edges and corner points. Starting from each corner point, the path is continuously traced along the boundary edge to the next corner point, resulting in several boundary path segments. For each boundary path segment, the density of fixed points is determined based on the path length and the target mesh scale, and a set of fixed points is generated on the path segment. In this embodiment, the number of fixed points can be adaptively determined by the path length and the target scale, for example:

[0114]

[0115] in The path segment length, The target grid scale; in this embodiment, it can be... The number is adjusted to the closest even number to improve the uniformity of the segments, but this invention is not limited to this adjustment.

[0116] S42. A coarse-layer fixed point set is constructed by adopting a layer-by-layer sparse sampling strategy, and an inter-layer mapping relationship is established to realize the propagation of fixed point constraints in a multi-resolution hierarchical structure.

[0117] In the finest layer, samples are taken along the boundary path segment at preset intervals, and fixed points are marked, with corner points being marked as fixed points first to ensure the stability of critical boundary positions. For coarser layers, a sparser sampling strategy is used to generate the set of fixed points for that layer, in order to reduce constraint density and improve optimization efficiency.

[0118] In this embodiment, the first The sampling interval of a layer can be set in an exponential increment, for example:

[0119]

[0120] in The layer numbering is used, but the specific interval function form is not limited in this invention. For each layer, the set of three-dimensional coordinates of fixed points is recorded, and inter-layer mapping relationships are established for subsequent cross-layer propagation and consistency maintenance.

[0121] S43. Next, constraint point identification and propagation are performed. For each fixed point, a search is conducted within its neighborhood for mesh vertices that can serve as "actual fixed points", and the fixed point constraint is assigned to the best-matching candidate vertex.

[0122] In this embodiment, a spatial neighborhood-based search method can be used, such as constructing a local search region and selecting the nearest candidate vertex by distance;

[0123] Furthermore, the search range can be appropriately expanded at different levels to enhance cross-level stability; for example, the search radius can be increased with each level.

[0124]

[0125] However, this invention does not limit the specific search structure (such as bounding box, sphere neighborhood, etc.) or the radius setting form. For candidate vertices with repeated competition, the conflict can be resolved by rules such as "nearest priority" or "unoccupied priority", thereby ensuring the uniqueness and uniformity of the fixed point set.

[0126] Fixed-point constraints in position field optimization can be categorized into hard and soft constraints: corner fixed points are used for hard constraints, ensuring that the position of their corresponding vertices remains unchanged during optimization; non-corner fixed points are used for soft constraints, introducing a weighted penalty term to limit the degree to which the position field deviates from the fixed point. Through this multi-level propagation mechanism combining hard and soft constraints, coarse-level constraints remain consistent across fine-level constraints, improving the convergence and robustness of the overall optimization while ensuring the stability of boundary features. In this embodiment, soft constraints can be implemented by introducing a fixed-point penalty term into the position field optimization objective function, for example, by adding the following constraint term:

[0127] ,

[0128] in For a fixed set of points, The target location is a fixed point. The corresponding constraint weight parameters are used to balance the maintenance of fixed points with the smoothness of the overall position field; this invention does not limit the specific form of the constraint terms or the way the weight values ​​are taken.

[0129] S5. Randomly initialize the position field and obtain a position field that satisfies the geometric constraints through iterative optimization;

[0130] This step generates the position field using a multi-resolution iterative optimization strategy. Starting with the random initialization of the position field, a multi-level mesh structure is constructed, and finally, a discrete position field satisfying geometric constraints is obtained through iterative optimization.

[0131] S51. Based on the multi-resolution hierarchical structure constructed in step S3, establish a local grid coordinate system for each vertex and randomly initialize the position field in the tangent plane.

[0132] First, establish a local grid coordinate system for each vertex. For vertex i, its normal vector is known. and the direction of the main stadium Calculate the intersection of the tangent plane and the plane. orthogonal unit vectors Construct a local orthogonal coordinate system. Generate random numbers within this coordinate system. And calculate the random offset vector. :

[0133]

[0134] Set initial position field And verify that it satisfies .

[0135] Subsequently, in each layer of the multi-resolution hierarchical structure in step S3, position field iterative optimization is performed sequentially from the coarsest layer to the finest layer, and the optimization result of each layer is used as the initial value for the next finest layer and propagated layer by layer.

[0136] S52. In each iteration, the optimal grid alignment between adjacent vertices is obtained through position compatibility calculation, and continuous positions are rounded to discrete grid positions, thereby generating a discrete and regular position field. For the boundary and its neighboring vertices, fixed point constraints or position constraints are introduced during the optimization process to improve the stability of the boundary point distribution and the geometric fidelity of the boundary.

[0137] S521. Collect position field suggestions based on position compatibility for each vertex and obtain the target position of the vertex through weighted fusion; correct the target position by combining the constraint position and weight of the geometric constraints of the vertex, and reproject it onto the tangent plane.

[0138] Location field optimization is performed layer by layer downwards, starting from the coarsest layer. For each mesh layer, location field suggestions from adjacent vertices are collected for each vertex, and the target location is obtained through weighted fusion.

[0139]

[0140] in, The positional suggestion for vertex i given by vertex j. For weighted weights.

[0141] If the vertex has geometric constraints, the target position is corrected by combining the constraint position and weight, and then reprojected onto the tangent plane.

[0142] Position compatibility calculations are used to determine the optimal position field configuration between adjacent vertices. For a pair of adjacent vertices, the optimal intersection point is first calculated, which lies simultaneously in the tangent plane of both vertices and minimizes the sum of the squared distances to both vertices. In this embodiment, the optimal intersection point can be solved analytically or iteratively using numerical methods; the invention is not limited to a specific solution method. The optimal intersection point is attached to the grid coordinate system of the two vertices to obtain the reference grid point position. A local coordinate system is constructed, and the relative position is calculated to obtain the grid point coordinates.

[0143]

[0144] in `floor` is the floor function. The base grid point position is calculated as follows:

[0145]

[0146] Starting from each reference grid point, considering the four possible offsets of 4-RoSy, a total of 16 location combinations are formed. All combinations are enumerated, and the distance cost for each combination is calculated. ,in and They respectively represent the combinations of candidates. and The calculated 3D reference positions of the candidate grid points at both ends are used. The combination that minimizes the distance between the two grid point positions is selected as the final compatible position configuration.

[0147] S522. After collecting adjacent vertex suggestions and calculating the weighted average, continuous target positions need to be rounded to discrete grid positions. For vertex i and target position target, the relative displacement is first calculated in the local coordinate system. Then, the projection components in the q and t directions are calculated respectively, and the grid point coordinates are determined by rounding.

[0148]

[0149] The final rounding position field value is calculated as follows:

[0150]

[0151] in and This forms an orthogonal basis for the local grid coordinate system. The update result must satisfy... In this embodiment, it is allowed As a tolerance range.

[0152] S523. After completing the optimization at a certain level, the position field results of that level are propagated to a finer level and used as initial values ​​for further optimization. For vertices with geometric constraints, a weighted constraint penalty term is introduced during the optimization process to maintain feature consistency. The objective function can be expressed as:

[0153] ,

[0154] in For position-fixed constraints, For orientation alignment constraints, , These are the weight parameters.

[0155] Figure 4 (The green line represents the boundary edge, and the red dots represent the position field coordinates) and Figure 5 This is a comparison diagram of the position field generation effects, in which... Figure 4 The result of position field generation after introducing boundary fixed point constraints is shown. Green lines represent boundary edges and red dots represent discrete position field coordinates. Figure 5 The results of position field generation without the introduction of boundary fixed point constraints are shown to illustrate the situation where the position field shifts or is unevenly distributed near the boundary under the condition of lack of boundary constraints.

[0156] S6. Calculate the principal curvature information of each vertex of the subdivided triangular mesh obtained in step S1, and establish a curvature-guided path search mechanism based on feature boundary connectivity constraints to generate curvature-aligned partition boundaries. Establish a mapping relationship from vertices to partitions based on the partitioning results, and perform partition compatibility checks to ensure connectivity between adjacent regions, obtaining structured region partitioning results and providing partition constraints for mesh extraction in step S7. This step improves the stability of mesh generation in complex topological regions such as porous and flat regions, and reduces the risk of mesh collapse.

[0157] Mesh adaptive partitioning based on principal curvature. For example... Figure 6 As shown, this step adaptively partitions the mesh based on the connectivity between the principal curvature direction of the vertices and the feature boundary, so that the partition boundary is consistent with the geometric features, and provides a structured region division result for subsequent mesh extraction.

[0158] S61. First, obtain the principal curvature information of the mesh vertices. For vertex i, construct a discrete curvature tensor matrix in its neighborhood. Then, eigenvalue decomposition is performed to obtain the principal curvature direction and principal curvature magnitude:

[0159] ,

[0160] in Represents vertices The set of neighborhood edges, For the normalized area term, To be with the edge The change in normal of adjacent facets Let be the side length. Let be the edge direction vector. This is derived from the eigenvalue matrix. The magnitude of the principal curvature is obtained from the eigenvalues ​​of the eigenvector matrix. The principal curvature direction is obtained from the eigenvectors. In this embodiment, the principal curvature information can be obtained through the curvature tensor method; however, this invention does not limit the specific curvature calculation method.

[0161] S62. Perform feature boundary detection based on curvature information. This involves detecting shared edges between adjacent facets. Calculate the dihedral angle And based on this, determine the feature edges:

[0162]

[0163] in , Let the unit normal vectors be the adjacent faces on both sides of the shared edge. The preset threshold is used. Feature edges serve as partition boundary constraints and are prohibited from being crossed during partition expansion.

[0164] S63. Next, perform the curvature-guided path search. (For mesh edges) Constructing the curvature alignment cost function The cost is used to measure the consistency between the edge direction and the reference direction, as well as the principal curvature direction of the endpoints, while also considering the impact of curvature magnitude changes on path selection. Based on the cost, candidate paths are... Calculate path quality metrics:

[0165]

[0166] in The weights are the edge lengths of the edges along the path. To evaluate the consistency between the path and the principal curvature guidance criterion, priority is given to selecting the path that... The optimal path is used as the partition boot path.

[0167] In this embodiment, the curvature alignment cost It can be composed of a direction alignment term, a curvature sensitivity term, and an angle penalty term, for example, it can take:

[0168]

[0169] in This is the orientation alignment coefficient. The curvature combination value of the two vertices at the two ends of the edge. Angle penalty item, The weighting parameters are used to adjust curvature sensitivity; the specific calculation methods and parameter values ​​of the above items are not intended to limit the invention.

[0170] S64. Based on the above curvature guidance criterion, perform connectivity expansion starting from the seed facet, assigning partition identifiers to the facets. The partition expansion process must satisfy the following constraints:

[0171] ,

[0172]

[0173] in For the current face, For the candidate extended patches adjacent to it, the mesh is divided into several connected partitions according to the above extension rules, and a mapping relationship between patches and partitions is established. , This indicates that the patch has not yet been assigned a partition identifier, and the expansion process is not allowed to cross feature boundaries. For the partition mapping of vertex v.

[0174] S65. After completing the initial partitioning, the partitioning results are optimized and constraints are applied, including merging excessively small partitions, adjusting thin or anomalously shaped regions, and smoothing partition boundaries. Where geometric features need to be preserved, further partition boundary constraints can be applied to improve the stability and feature preservation capabilities of the partitioning results. The final structured region partitioning result will be used as the partitioning constraints for vertex merging and connection relationship establishment in step S7.

[0175] S7. Based on the structured region division results, a two-stage processing framework is used to generate a quadrilateral dominant mesh:

[0176] The first stage involves graph structure extraction and vertex merging. Compatibility analysis of frame field and position field is performed on adjacent vertices to identify vertex pairs that need to be merged, connected, or ignored. Reasonable vertex connection relationships are established, and quality optimization processing is performed on flat regions.

[0177] The second stage involves face extraction and mesh generation. A local coordinate system is established for the vertices and adjacent edges are sorted by angle. A ring tracing strategy is used to extract the target face. In the case of tracing failure, a backtracking search strategy is used to select the optimal tracing starting point, thereby generating a topologically correct and geometrically good quadrilateral dominant mesh.

[0178] Quadrilateral mesh extraction. This step uses a two-stage process to generate a quadrilateral-dominant mesh. The first stage is used to construct the graph structure and perform vertex merging and connection relationships, such as... Figure 7 As shown, the orange lines represent the extracted map. The second stage is used for patch tracking to generate patches and output a mesh, as shown below. Figure 8 As shown. This process is based on the standard frame field obtained in steps S3 and S5. With position field By analyzing the compatibility of adjacent elements and processing topological constraints, a graph structure that can be used for patch generation is obtained.

[0179] Phase 1: Graph Structure Extraction and Vertex Merging

[0180] S71. For any pair of adjacent vertices in the original triangular mesh Perform frame field compatibility processing to eliminate orientation inconsistencies caused by the fourfold rotational symmetry at adjacent vertices. Orientation compatibility can be expressed as:

[0181]

[0182] in , For the compatible direction vector, , This represents the unit normal vector of the corresponding vertex. Compatibility is achieved by enumerating four rotation configurations and selecting the configuration that optimizes the alignment. In this embodiment, the invention selects the configuration that maximizes the dot product of the frame field vectors, ensuring good continuity of the frame field between adjacent vertices.

[0183] S72. Based on the completion of the target field compatibility processing, and based on the position field and target scale parameters... (or Calculate the discrete position compatibility relationship of adjacent vertex pairs and output the corresponding index offset (or equivalent discrete displacement).

[0184]

[0185] in , This is the position offset. , The vertex position, , For position field value, The target mesh scale parameter. This compatibility relation is used to determine the relationship type of adjacent vertices: merge, connect, or ignore. In this embodiment, the determination can be based on the norm condition of the offset difference:

[0186]

[0187] When the infinite norm of the offset difference is 0, it indicates that the two vertices should be merged into one vertex; when the infinite norm is 1 and the L1 norm does not exceed 1, it indicates that the two vertices should be connected by an edge; otherwise, the relationship between the vertex pairs is ignored.

[0188] S73. For vertex pairs that need to be merged, a disjoint-set data structure is used for efficient management, and the merging weight is calculated based on the position field error to determine the processing priority:

[0189]

[0190] in and ′ represents the corresponding position field (or equivalent reference position) after compatibility processing. According to Merge from smallest to largest, and perform adjacency conflict detection before merging to avoid topological errors such as self-loops or multiple edges.

[0191] After merging, the geometric properties (position, normal, and partition information, etc.) of each representative vertex are recalculated. In this embodiment, a weighted average method can be used to calculate the merged properties, for example, for the vertex set within the merged cluster. :

[0192]

[0193] In this embodiment, weight The weights can be Gaussian weights or other monotonic weight functions that are related to the position field deviation; the present invention does not limit the form of the weights.

[0194] S74. After establishing the merging and connection relationships, perform quality constraint processing on the graph structure, such as identifying and handling anomalies like flattened triangles and thin regions, to improve the stability of subsequent patch generation. In this embodiment, partitioning information can be introduced to restrict cross-partition operations. The specific constraint rules are as follows:

[0195] Obtain the latest partition information of the vertices to be processed, using Represents vertices The number of vertices in the partition, when satisfying

[0196]

[0197] or

[0198]

[0199] Merge operations are disabled to prevent improper merges or joins across geometric partitions.

[0200] Phase 2: Patch Extraction and Mesh Generation

[0201] In the second stage, patch tracing is performed based on the graph structure obtained in the first stage to generate patches and output the mesh.

[0202] S75. Establish a local coordinate system for each vertex and arrange its adjacent edges in order of local angles to ensure consistent direction during the tracing process. For vertices... adjacent vertices It can calculate the azimuth angle of the side in the local coordinate system:

[0203]

[0204] in As vertices Orthogonal basis within the local tangent plane. The adjacent edges are ordered to support consistent ring tracing.

[0205] S76. Generate a face of the target size using a circular tracing strategy. Tracing starts from the initial vertex and initial edge, and proceeds step by step according to the edge usage status and adjacency order, recording the boundary sequence. When returning to the initial vertex and the sequence length meets the target face size, tracing is considered successful and the face is output. To avoid topological errors, used edges are marked during the tracing process, and self-intersections or non-closed boundaries caused by reuse are prohibited.

[0206] S77. When local topology complexity leads to tracking failure, this embodiment employs backtracking search (such as DFS backtracking) to switch between candidate branches to improve the success rate; and can calculate quality scores for candidate starting points or candidate branches, selecting the scheme with the best score. For example, a scoring function can be defined based on the deviation of the angle from the ideal angle (taking a quadrilateral as an example, the ideal angle is...). ):

[0207]

[0208] And select The invention does not limit the form of the scoring function to tracking starting from a smaller point.

[0209] S78. For excessively large faces with more vertices than the target size, in this embodiment, a centralization or sector partitioning strategy is used to split them into several quadrilateral faces in order to improve geometric quality and reduce the instability caused by large faces.

[0210] S79. The generated quadrilateral facets... It can calculate the surface normal vector and perform orientation consistency checks, for example:

[0211]

[0212] By comparison The consistency with vertex normals determines whether the face order needs to be flipped, thus avoiding normal flipping.

[0213] Finally, the mesh is geometrically optimized. In this embodiment, Laplacian smoothing can be used to improve the geometric quality, where boundary vertices remain unchanged and internal vertices are updated by averaging their neighborhood values; and a topological consistency check is performed on the generated mesh (e.g., each internal edge is shared by exactly two faces) to ensure that the output mesh topology is correct.

[0214] like Figure 8 The quadrilateral-dominant mesh generation results and such Figure 9 The mesh details shown in Figure 10, through the above two-stage processing, enable the stable generation of topologically correct and feature-preserving quadrilateral dominant meshes even in the presence of holes, thin sheets, and areas where feature lines intersect.

Claims

1. A method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field, characterized in that, Includes the following steps: S1. Obtain and subdivide the preprocessed engineering simulation triangular surface mesh model to generate a subdivided triangular mesh; S2. Initialize the frame field direction of each vertex on the subdivided triangular mesh to obtain the initial frame field direction, including: for internal vertices, randomly generate a unit tangent vector in the tangent plane of the vertex as the initial direction; for boundary vertices, construct dihedral angle related directions based on the normal vectors of adjacent facets to reflect crease features; construct boundary polyline direction based on adjacent boundary edges to reflect the boundary direction; adaptively calculate the fusion weight according to the size of the dihedral angle, weight and fuse the dihedral angle related directions and the boundary polyline direction, and normalize them after projecting them onto the tangent plane to obtain the initial direction of the boundary vertices; S3. Based on the initial frame field orientation, construct a multi-resolution hierarchical structure; iteratively optimize each layer of the multi-resolution hierarchical structure by constructing neighborhood consistency update rules and orientation constraint terms to generate the final orientation field; S4. Establish a multi-level boundary fixed-point constraint mechanism based on the multi-resolution hierarchical structure; S5. Randomly initialize the position field according to the multi-resolution hierarchical structure, and then obtain the position field that satisfies the geometric constraints through iterative optimization; S6. Adaptively partition the subdivided triangular mesh based on the connectivity between the principal curvature direction of the vertices and the feature boundary, so that the partition boundary is consistent with the geometric features, and generate structured region partitioning results. S7. Based on the final orientation field, position field, and structured region division results, a two-stage processing framework is used to generate a quadrilateral dominant mesh.

2. The method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field according to claim 1, characterized in that, S3 include: S31. A multi-resolution hierarchical structure is constructed using a bottom-up graph downsampling method. The product of the normal similarity of adjacent vertices and the dual area ratio is used as the merging priority. Vertices are merged layer by layer to generate multiple hierarchical meshes, and the mapping relationship between vertices in the layers is established. S32. Iteratively optimize the frame fields of each layer of the multi-resolution hierarchical structure in sequence. In each iteration, each vertex updates its own orientation according to the frame field orientation of its adjacent vertices, and applies orientation constraints to the boundary vertices.

3. The method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field according to claim 2, characterized in that, S32 includes: S321. For each layer of mesh, each vertex is weighted and averaged according to the frame field direction of its adjacent vertices, and its own direction is updated smoothly in the neighborhood. The fourfold rotational symmetry of the frame field is handled by a compatibility function to eliminate the inconsistency of direction caused by rotational symmetry. S322. Apply directional constraints to boundary vertices or feature line vertices; fuse the preset constraint directions with the neighborhood smoothing results to maintain the boundary and feature orientations. S323. On each mesh layer, minimize the objective function that includes the neighborhood consistency term and the orientation constraint term, and iteratively optimize the frame field.

4. The method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field according to claim 1, characterized in that, S4 include: S41. Obtain the boundary edge set and corner point set and perform boundary path tracing to generate boundary path segments; perform segmented sampling on the boundary path segments based on the target mesh scale to obtain the fine-layer fixed point set. S42. Implement a layer-by-layer sparse sampling strategy for the set of fixed points in the fine layer and establish inter-layer mapping relationships to realize the propagation of fixed point constraints in the multi-resolution hierarchical structure. S43. For each fixed point, search for grid vertices that can serve as fixed points within its neighborhood and assign fixed point constraints to the matching candidate vertices. Fixed point constraints are divided into two categories: corner fixed points are treated as hard constraints, and their positions remain unchanged during optimization; non-corner fixed points are treated as soft constraints, and their deviation is limited by introducing a weighted penalty term into the position field optimization objective function.

5. The method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field according to claim 1, characterized in that, S5 include: S51. Based on the multi-resolution hierarchical structure, establish a local grid coordinate system for each vertex and randomly initialize the position field in the tangent plane; S52. Perform position field iterative optimization layer by layer, and use the optimization result of each layer as the initial value of the next fine layer to propagate layer by layer.

6. The method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field according to claim 5, characterized in that, S52 includes: For each vertex, position field suggestions based on position compatibility are collected from its adjacent vertices, and the target position of the vertex is obtained by weighted fusion; the target position is corrected by combining the constraint position of the vertex's geometric constraints and weights, and then reprojected onto the tangent plane; The continuous target positions are projected onto the local coordinate system, and the coordinates of discrete grid points are determined by rounding to update the position field. For vertices with geometric constraints, a fixed-point penalty term and an orientation alignment constraint term are introduced into the objective function to maintain feature consistency.

7. The method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field according to claim 1, characterized in that, S6 include: S61. For each vertex, construct a discrete curvature tensor matrix in its neighborhood; The principal curvature directions and magnitudes are obtained from the discrete curvature tensor matrix through eigenvalue decomposition. S62. Determine the characteristic edge based on the dihedral angle of adjacent facets; S63. Construct a curvature alignment cost function and calculate path quality indicators for candidate paths; select partitioned guidance paths based on path quality indicators. S64. Starting from the seed face, perform connectivity expansion along adjacent faces that satisfy the unassigned and do not cross the feature edge constraints, assign partition identifiers to the faces, and establish mapping relationships from faces to partitions and from vertices to partitions. S65. Post-process and optimize the initial partitioning results to obtain the structured region partitioning results.

8. The method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field according to claim 1, characterized in that, In S7, the two-stage processing framework includes: The first stage involves graph structure extraction and vertex merging. Compatibility analysis of frame field and position field is performed on adjacent vertices to identify vertex pairs that need to be merged, connected, or ignored. Reasonable vertex connection relationships are established, and quality optimization processing is performed on flat regions. The second stage involves face extraction and mesh generation. A local coordinate system is established for the vertices and adjacent edges are sorted by angle. A ring tracing strategy is used to extract the target face. For cases where tracing fails, a backtracking search strategy is used to select the optimal tracing starting point.

9. The method for generating a feature-preserving quadrilateral dominant mesh guided by a geometrically constrained frame field according to claim 8, characterized in that, The first phase of processing includes: Perform frame field compatibility processing on any adjacent vertex pair (i,j) in the original triangular mesh; calculate the discrete position compatibility relationship of adjacent vertex pairs based on the position field and target scale parameters, and output the corresponding index offset; manage vertex pairs that need to be merged using a disjoint-set data structure, and calculate the merging weight based on the position field error to determine the processing priority; calculate the geometric attributes of each representative vertex after merging; and perform quality constraint processing on the graph structure based on the structured region partitioning results. The second phase includes: A local coordinate system is established for each vertex, and its adjacent edges are arranged in order of local angles. A ring tracing strategy is used to generate a facet of the target size. When the local topology is complex and the tracing fails, a backtracking search is used to switch between candidate branches. For faces that are too large and have more vertices than the target size, a centralization or fan-shaped partitioning strategy is used to split them into several quadrilateral faces. The quadrilateral faces are checked for orientation consistency and the geometric quality is optimized using Laplacian smoothing. Topology consistency is then checked.

Images