Plane parameterization method suitable for complex geometric region and having G1 smooth characteristic

Abstract

The invention provides a plane parameterization method which is suitable for a complex geometric region and has a G1 smooth characteristic. According to the technology, firstly, a subdivision structure is generated based on center axis information of a given area; then, by applying a G1 continuity constraint condition to adjacent surface patches, initial parameterization meeting G1 continuity on the whole is constructed; and finally, further optimizing the parameterization quality of each patch by using a quasi-conformal mapping technology. According to the method, rectangular and triangular Bessel patches are adopted to represent the geometric region at the same time, and the method can be effective for any topological complex domain. Compared with other methods in the prior art, not only is the design flexibility enhanced, but also the generation and parameterization quality of the multi-patch structure is remarkably improved, and the method is more efficient.

Description

Technical Field

[0001] The invention relates to the field of isogeometric analysis, in particular to a plane parameterization method suitable for complex geometric regions and having G1 smooth characteristics. Background Art

[0002] Isogeometric analysis (IgA) aims to bridge the gap between computer-aided design (CAD) and computer-aided engineering (CAE). Its core concept is to use a unified basis function with high smoothness, usually non-uniform rational B-splines (NURBS) to simultaneously represent the physical field in CAE and the CAD geometric model. Compared with the traditional finite element method, IgA shows many significant advantages, including but not limited to the accurate representation of geometric shapes, the high-order smoothness of the basis function, and the improvement of accuracy in each degree of freedom. These advantages make isogeometric analysis show its unique value and potential in dealing with engineering problems with complex geometry and high-order accuracy requirements.

[0003] A fundamental task in the IgA framework is to construct an analysis-aware parameterization from a given boundary representation of a physical domain, with basic requirements including injectivity, low distortion, and applicability to complex domains. Although this topic has been intensively studied in academia and engineering, most methods are limited to processing geometries with simple topological structures.

[0004] In theory, multi-piece parameterization technology can handle arbitrary topological structures, but in practical applications, most existing technologies still face challenges when dealing with complex domains. For example, in research based on multi-cubes, it is impossible to avoid the generation of singular points inside the region, which greatly affects the patch layout. In addition, skeleton-based methods also have limitations in obtaining high-quality parameterization results for complex shapes. Some methods that attempt to construct planes from complex boundaries through domain segmentation and global optimization, although feasible in theory, often lead to a significant increase in the number of patches in the parameterized structure in practice. The parameterization achieved by these methods can usually only reach global C 0 Continuity cannot meet the needs of high-precision analysis. Therefore, it is necessary to seek better new methods to overcome the above defects. Summary of the invention

[0005] In order to solve the problems of the prior art, the present invention provides a planar parameterization method suitable for complex geometric regions and having G1 smooth characteristics, which can be effective for any topologically complex domain. Compared with other existing methods, it not only enhances the flexibility of design, but also achieves significant improvements in the generation and parameterization quality of multi-faceted structures, and is more efficient.

[0006] In order to achieve the above object, the present invention provides a method having G 1The plane parameterization method of smooth characteristics, the specific steps are as follows:

[0007] 1) Automatic domain decomposition with gradual refinement, using triangular and rectangular Bezier patches to represent complex domains, including:

[0008] 1.1) Extracting the skeleton based on the region boundary, using a sparse optimization-based method; the skeleton extraction process is as follows:

[0009] A sparse optimization-based approach is used to extract the skeleton and compute a compact medial axis transform for any 2D object represented by a point cloud containing noise or missing data. The branch points, branches, the maximum inscribed circle centered at each branch point, and the corresponding tangent points on the skeleton are then determined.

[0010] Assumption c t (u) is the terminal branch spline curve connected to the branch point b, r b is the radius of the largest inscribed circle with b as the center. Then delete the terminal branches that meet the following conditions

[0011]

[0012] The left part represents the arc length of the terminal branch, r max It is the maximum value of the radius of the circle centered at each branch point, and α and β are positive parameters.

[0013] 1.2) Use the skeleton to guide the initial decomposition process, connect the corresponding tangent points on the domain boundary, and automatically generate triangles and rectangular patches to build the basic geometric structure;

[0014] 1.3) Based on the preliminary decomposition, the shapes of the remaining patches are further identified, and appropriate corner points are selected for each patch. The process of obtaining the corner points is as follows:

[0015] After completing the initial decomposition of the region Ω, further analysis and identification of the remaining patches are performed to determine whether there are patches that can be approximated to triangles. for A vertex of The boundary (edge )for An edge of , which is inside Ω, is also called the internal boundary. Representation patch The perimeter, L p represents the perimeter of Δpp′p″, and δ>0 is used as a threshold. satisfy

[0016]

[0017] It is believed that At this time, point p and the two endpoints p′ and p″ of the corresponding internal boundary are selected as the three corner points required for subsequent parameterization.

[0018] After the triangular patches are identified and distinguished, the remaining patches are classified as quadrilateral patches according to the default rule, and four corner points need to be selected on the domain boundary for each quadrilateral patch. Two of the tangent points have been determined to be used as corner points. For the remaining two corner points, an automated method based on the optimal transport principle is used for selection.

[0019] 1.4) Continue to finely decompose the region and construct a mathematical model of the obtained surface boundary segment using Bezier curves; the process of creating the mathematical representation of the Bezier curve is as follows:

[0020] For N discrete points The external boundary segment represented by the chord length parameterization method is used to determine the relevant parameters And use the control point c to fit a p-order Bezier curve C(u), which is done by minimizing

[0021]

[0022] Where λ is a positive weight and the regularization term

[0023]

[0024] It is used to measure the smoothness of the obtained curve. By verifying the mean square error defined by the following formula

[0025]

[0026] Is it greater than 4×10 -8 The fitting accuracy is evaluated with the specified tolerance. If the error value exceeds the specified tolerance, the corresponding patch needs to be further decomposed. The patches with large boundary fitting errors are usually quadrilaterals, so the decomposition process is as follows: for the longer of the two boundaries, a point k1 is selected, which has approximately the same geodesic distance to the two end points, and its closest point k2 on the opposite boundary is found. By connecting k1 and k2, the patch is decomposed into two parts. This decomposition and fitting process will be iterated until the boundary fitting errors of all patches are reduced to below the specified tolerance.

[0027] Then, we process the internal boundary segments. For the internal boundary segments (straight line segments) whose two endpoints are p1 and p4, we interpolate at p1 and p4 to construct a cubic Bezier curve to replace the original line segment. Then, we construct two other control points p2 and p3 according to the tangent directions of p1 and p4. Specifically, let is the tangent vector of two adjacent Bezier curves at p1, and is the tangent vector of two adjacent Bezier curves at p4, and the control points p2 and p3 are constructed as follows:

[0028] ∠p2p1p4≤σ,

[0029] and

[0030] ∠p3p4p1≤σ,

[0031] Among them, ρ and are the degrees of p1 and p4 respectively, and θ is given by and The included angle formed, Is and The included angle formed, σ is the threshold. When the above two equations contradict each other, the following conditions are adopted

[0032] ∠p2p1p4=σ,

[0033] ∠p3p4p1=σ,

[0034] 2) Build a global G 1 Smooth parameterized model:

[0035] 2.1) Using the discrete Coons interpolation method, the initial multi-slice parameterization is calculated; the process of the multi-slice parameterization is as follows:

[0036] For the first quadrilateral patches, Its boundary control points are Initial control point The structure is

[0037]

[0038] For Triangular Bezier surface represented by triangle patches The calculation formula for the internal control point that satisfies i+j+k=p is

[0039]

[0040] 2.2) Apply G on the patch interface 1 Continuity conditions, calculated to satisfy the global G 1 Smooth control points; the G is applied at the patch interface1 The process of continuity conditions is as follows:

[0041] Let a p-order Bezier curve

[0042]

[0043] is the common boundary of two adjacent rectangular Bezier patches, and its control points Then these two patches are G 1 Continuous, if there are scalar weight functions λ(u) and μ(u) such that

[0044] (1-λ(u))s(u)+λ(u)t(u)=(1-μ(u))q l (u)+μ(u)q r (u)

[0045] holds for all u where

[0046]

[0047] and are two rows of control points adjacent to a common boundary.

[0048] Here, the functions λ(u) and μ(u) are set as linear polynomials and quadratic polynomials as follows:

[0049]

[0050] Substituting λ(u) and μ(u) into the previous equation, we get the following G 1 Continuity Constraints

[0051]

[0052] They form a linear system.

[0053] When processing adjacent triangles and rectangular patches, first increase the order of the triangle patch from p to p+1, and then use the method described above to construct the above G 1 Restraint system.

[0054] In this linear system, the boundary control points are known and fixed s0,s p ,t0,t p , and then the remaining control points need to be calculated. The coefficients λ0, λ1 of λ(t) and the coefficients μ0, μ2 of μ(t) are determined by the first and last equations of the system respectively. Set the coefficient μ1 of μ(t) to (μ0+μ2) / 2. 1 The continuous system is divided into two groups, namely, the vertex G1 Continuous system (corresponding to k = 0, 1, p, p + 1) and edge G 1 Continuous system (corresponding to k=2,...,p-1).

[0055] For each vertex (singular corner or regular corner) in the constructed multi-faceted structure, that is, the intersection of ρ Bezier facets, control vertices o1, o2, ..., o ρ The twist vector at the vertex should satisfy the G corresponding to k = 1 or p in the system 1 Continuity constraints. All these constraints form the linear system Vx v =y v , where the vector x v Represents control point The set of is obtained by solving the following problem

[0056]

[0057] VxD v =y v ,

[0058] x v Initial value of It is constructed by discrete Coons interpolation method.

[0059] For each common edge (boundary) shared by two Bezier patches, the two rows of control points adjacent to that boundary and Should satisfy G corresponding to k=2,...,p-1 1 Continuity constraints. These constraints form a linear system Ex e =y e , where the vector x e Represents control point and The set of is obtained by solving the following problem

[0060]

[0061] stEx e =y e

[0062] x e Initial value of It is also constructed by discrete Coons interpolation method.

[0063] 3) Using QC mapping technology to obtain high-quality parameterization of each patch; the process of calculating high-quality parameterization of the patch is as follows:

[0064] First The parameterization of a quadrilateral patch (i.e. a rectangular Bezier surface) is rewritten as a complex function

[0065]

[0066] in and here and It is a control point In this parameterization, the boundary control points and the control points adjacent to the common boundary are known, and the remaining interior control points are determined by solving the following optimization problem

[0067]

[0068] in yes The Beltrami coefficient of , the parameter γ is a positive weight, and Respectively and The Hessian matrix of .

[0069] When solving this nonlinear and non-convex problem, an auxiliary variable v is introduced and used to replace This leads to an interleaving problem

[0070]

[0071] ||v|| ∞ <1.

[0072] Here we further relax it into such a model

[0073]

[0074] st||v|| ∞ <1.

[0075] The weight η is relatively large. Here, an energy descent method is used to alternately solve and optimize the two sub-problems. The specific iterative process is as follows: 1 Parameterization Start by solving for v alternately. k and

[0076]

[0077] Among them, v k By fixed After solving the model, we get; By fixing v kThen the model is solved. This process continues until the solution of the model satisfies the condition ||v for a given threshold ε k+1 -v k Until ||<ε.

[0078] The beneficial effects of the present invention are:

[0079] 1. The present invention uses a combination of triangular Bezier surfaces and rectangular Bezier surfaces to perform a hybrid representation of complex domains, which is more flexible than traditional representation methods, can adapt to complex domains of arbitrary topology, and has a wider range of applications.

[0080] 2. The present invention combines a skeleton-guided coarse decomposition with a more refined decomposition procedure to achieve a good patch layout. Unlike the existing QMGLGO method and PeDP method, an advantage of the skeleton-guided domain decomposition process proposed in the present invention is that there are no internal singularities in the generated multi-faceted structure. This helps to generate better parameterizations for complex domains, thereby reducing distortion. In addition, by taking a series of operations to optimize the layout of the multi-faceted structure, such as automatically merging or further subdividing adjacent triangular patches as appropriate during the coarse decomposition process, the number of patches is reduced while maintaining geometric accuracy, which is particularly desirable in the case of IgA.

[0081] 3. The present invention uses QC mapping technology to ensure the injectivity and low distortion of the final parameterization. QC mapping is locally injective, and its angular distortion can be measured by the Beltrami coefficient. This is used as a framework to construct constraint conditions for parameterization, ensuring the injectivity of parameterization. In addition, the present invention uses an energy descent method to alternately solve difficult optimization problems in the parameterization process, so that when the same number of patches are generated, higher quality patches are constructed than other technologies. BRIEF DESCRIPTION OF THE DRAWINGS

[0082] In order to more clearly illustrate the technical solutions in the embodiments of the present invention, the drawings required for use in the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0083] Figure 1 It is a schematic diagram of the overall process of the technology of the present invention;

[0084] Figure 2 This is a diagram of the preliminary decomposition process based on the skeleton in the present invention;

[0085] Figure 3 It is the classification process of adjacent triangular facets in the preliminary decomposition process of the present invention;

[0086] Figure 4 This is a detailed decomposition process diagram in the present invention;

[0087] Figure 5 A schematic diagram of control points in the preliminary parameterization process of the present invention;

[0088] Figure 6 The present invention preliminarily parameterizes and uses QC mapping technology to obtain a high-quality facet schematic diagram. DETAILED DESCRIPTION

[0089] The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0090] The process of the present invention is as follows Figure 1 As shown in the figure, for complex domains with arbitrary topology, the region is first automatically decomposed level by level based on Bezier patches, and then the global G is constructed. 1 Smooth parameterization, and finally use QC mapping technology to achieve high-quality parameterization for each patch. It can be divided into three steps:

[0091] 1) Automatic region decomposition with gradual refinement

[0092] The skeleton of an object describes its complete shape and reflects its structure. Therefore, the present invention performs domain decomposition based on the skeleton information of the input shape.

[0093] First, the present invention extracts its skeleton from the boundary representation of the physical domain, which is a set of dense points. To this end, a sparse optimization-based method is used to extract the skeleton, and a compact medial axis transformation is calculated for any two-dimensional object represented by a point cloud containing noise or missing data. The key innovation that distinguishes the present invention from the prior art is that spline curves are used instead of traditional broken lines to represent skeleton branches. This representation method is more compact because it only needs to store the control points of the spline curves, thereby avoiding the need to store a large number of discrete points and optimizing the compactness and efficiency of the data. Then determine the branch points of the skeleton, the branches, the maximum inscribed circle centered at each branch point, and the corresponding tangent points on the skeleton, such as Figure 2 Shown in the left part.

[0094] It is worth noting that some small terminal branches (such as the one located at the bottom right of the backbone) may provide only minimal structural information about the domain. Therefore, they may not be essential elements in the subsequent decomposition process. t (u) is the terminal branch spline curve connected to the branch point b, r bis the radius of the largest inscribed circle with b as the center. To simplify the original skeleton of the region, delete the terminal branches that meet the following conditions

[0095]

[0096] The left part represents the arc length of the terminal branch, r max is the maximum value of the radius of the circle centered at each branch point. α and β are positive parameters. In the example, α is set to 1.0 and β is set to 0.35.

[0097] Use as Figure 2 The simplified skeleton in the second part on the left is then connected to the corresponding tangent points on the domain boundary to obtain the preliminary decomposition result, and automatically generate triangles and rectangular patches to construct the basic geometric structure, such as Figure 2 Shown in the second part from the right.

[0098] For a general domain, the degree of the branch point is usually 3, which is also the case in this example, resulting in the generation of several triangular facets in this step. In order to optimize the layout of the generated multi-facet structure and reduce the number of facets, the present invention adopts a strategy of pairing and merging adjacent triangular facets into quadrilateral facets, such as Figure 2 However, there are two cases that require special handling: Figure 2 In the right part, there is only one pair of adjacent triangle patches. In a coarse decomposition, a particular domain may be Figure 3 As shown on the left, there are more than two adjacent triangular patches. To this end, these patches should be grouped first, and then the adjacent patches are merged to form quadrilateral patches, such as Figure 3 As shown in the middle figure. Alternatively, you can choose to merge the middle adjacent faces to maintain symmetry, as shown in Figure 3 As shown on the right, this is more conducive to subsequent parameterization processing.

[0099] For domains with rotational symmetry, the degree of some branch points may exceed 3. Connecting the corresponding tangent points of each branch point in a clockwise direction along the domain boundary may result in a patch with more than four edges. If the number of edges of a patch is odd, it should be further subdivided into multiple quadrilateral patches and one triangular patch. On the contrary, if the number of edges is even, it is subdivided into several quadrilateral patches.

[0100] Then, based on the preliminary decomposition, the shapes of the remaining patches are further identified, and appropriate corner points are selected for each patch.

[0101] The coarse decomposition process automatically generates quadrilateral patches (denoted by "Q") and triangular patches (denoted by "T") from the input shape, such as Figure 4 The first row shows the remaining faces. Some of them are triangular, such as Figure 4The area defined by the red curve is shown. The present invention identifies these facets by analyzing the shape differences between these facets and triangles. for A vertex of The boundary (edge )for An edge of , which is inside Ω, is also called the internal boundary. Representation patch The perimeter, L p represents the perimeter of Δpp′p″, and δ>0 is used as a threshold. satisfy

[0102]

[0103] It is believed that At this time, point p and the two endpoints p′ and p″ of the corresponding internal boundary are selected as the three corner points required for subsequent parameterization.

[0104] After identifying and distinguishing the triangular patches, the remaining patches are classified as quadrilateral patches according to the default rule. For these quadrilateral patches, the prerequisite for subsequent parameterization is to select four corner points on the domain boundary. Two of the tangent points have been determined to be used as corner points. For the remaining two corner points, the present invention adopts an automated method based on the optimal transmission principle for selection. Figure 4 As shown in the first row, the corner points of all patches have been determined. Specifically, each triangular patch has three boundary segments, and each quadrilateral patch has four boundary segments.

[0105] Continue to refine the decomposition of the region and use Bezier curves to construct a mathematical model of the obtained surface boundary segment.

[0106] For N discrete points The external boundary segment represented by the chord length parameterization method is used to determine the relevant parameters And use the control point c to fit a p-order Bezier curve C(u), which is done by minimizing

[0107]

[0108] Where λ is a positive weight and the regularization term

[0109]

[0110] It is used to measure the smoothness of the obtained curve. By verifying the mean square error defined by the following formula

[0111]

[0112] Is it greater than 4×10 -8 The fitting accuracy is evaluated with the specified tolerance. If the error value exceeds the specified tolerance, the corresponding patch needs to be further decomposed. Figure 4 The second row, left image, shows the fitting results of the outer boundary segments, indicating that the four patches marked with "×" need to be decomposed. From the tested examples, it can be seen that the patches with large corresponding boundary fitting errors are usually quadrilaterals. This is because the boundary line segments of triangular patches are relatively simple in shape and can be accurately represented by Bezier curves.

[0113] Take one of the patches with a large fitting error as an example to illustrate the process of further decomposition. Figure 4 In the middle of the second row, for the longer of the two boundaries, select a point p1 that is roughly the same geodesic distance from both ends, and find its closest point p2 on the opposite boundary. This patch is decomposed into two parts by connecting p1 and p2. This decomposition and fitting process is iterated until the boundary fitting error of all patches falls below the specified tolerance. Figure 4 The second row, right figure, shows the fine decomposition results for the giraffe domain.

[0114] The shape of the internal boundary greatly affects the subsequent parameterization results. For an internal boundary segment (straight line segment) with two endpoints p1 and p4, interpolation is performed at p1 and p4 to construct a cubic Bezier curve to replace the original segment. Then, two other control points p2 and p3 are constructed based on the tangent directions of p1 and p4. Figure 4 As shown in the third row, is the tangent vector of two adjacent Bezier curves at p1, and is the tangent vector of two adjacent Bezier curves at p4, and the control points p2 and p3 are constructed as follows:

[0115] ∠p2p1p4≤σ,

[0116] and

[0117] ∠p3p4p1≤σ,

[0118] Among them, ρ and are the degrees of p1 and p4 respectively, and θ is given by and The included angle formed, Is and The included angle formed is σ, and σ is the threshold. The second inequality is to prevent the constructed cubic Bezier curve from being inconsistent with the original internal boundary. The deviation exceeds the predetermined range, because excessive deviation is likely to destroy the overall multi-faceted structure generated in the previous stage. The last few steps ensure that p2 and p3 are in The projection on is located at one third of the distance on this line. However, in some specific cases, the above two equations contradict each other, so the following condition is adopted

[0119] ∠p2p1p4=σ,

[0120] ∠p3p4p1=σ,

[0121] Once the inner boundary represented by a cubic Bezier curve is constructed, its degree can be further increased to p without changing its basic shape. Figure 4 The last figure shows the final breakdown of the Giraffe domain.

[0122] 2) Build a global G 1 Smooth parameterized model

[0123] First, the initial multi-slice parameterization is calculated using the discrete Coons interpolation method.

[0124] For the first quadrilateral patches, Its boundary control points are Initial control point The structure is

[0125]

[0126] For Triangular Bezier surface represented by triangle patches The calculation formula for the internal control point that satisfies i+j+k=p is

[0127]

[0128] Then apply G on the patch interface 1 Continuity conditions, calculated to satisfy the global G 1 Smooth control points.

[0129] Let a p-order Bezier curve

[0130]

[0131] is the common boundary of two adjacent rectangular Bezier patches, and its control points Then these two patches are G 1Continuous, if there are scalar weight functions λ(u) and μ(u) such that

[0132] (1-λ(u))s(u)+λ(u)t(u)=(1-μ(u))q l (u)+μ(u)q r (u)

[0133] holds for all u where

[0134]

[0135]

[0136] and are two rows of control points adjacent to a common boundary, such as Figure 5 shown.

[0137] Functions λ(u) and μ(u) can be expressed in a variety of mathematical forms. Here they are set as linear polynomials and quadratic polynomials, as shown below:

[0138]

[0139] Substituting λ(u) and μ(u) into the previous equation, we get the following G 1 Continuity Constraints

[0140]

[0141] They form a linear system.

[0142] Clearly, for some high-order rational functions λ and μ, almost all pairs of patches with a shared boundary satisfy the basic G 1 Continuity condition. The present invention only uses linear polynomials and quadratic polynomials for parameterization. Only when the patch pairs satisfy the conditions obtained here, they are considered to meet G 1 Continuity standards.

[0143] When processing adjacent triangles and rectangular patches, first increase the order of the triangle patch from p to p+1, and then use the method described above to construct the above G 1 Constraint system. There are no two adjacent triangle patches here because they have been merged into a quadrilateral patch.

[0144] In this linear system, the boundary control points are known and fixed s0,s p ,t0,t p, and then the remaining control points need to be calculated. The coefficients λ0, λ1 of λ(t) and the coefficients μ0, μ2 of μ(t) are determined by the first and last equations of the system respectively. Set the coefficient μ1 of μ(t) to (μ0+μ2) / 2. 1 The continuous system is divided into two groups, namely, the vertex G 1 Continuous system (corresponding to k = 0, 1, p, p + 1) and edge G 1 Continuous system (corresponding to k=2,...,p-1).

[0145] For each vertex (singular corner or regular corner) in the constructed multi-faceted structure, that is, the intersection of ρ Bezier facets, control vertices o1, o2, ..., o r The twist vector at the vertex should satisfy the G corresponding to k = 1 or p in the system 1 Continuity constraints. All these constraints form the linear system Vx v =y v , where the vector x v Represents control point The set of can be obtained by solving the following problem

[0146]

[0147] VxD v =y v ,

[0148] x v Initial value of It is constructed by the discrete Coons interpolation method mentioned above.

[0149] For each common edge (boundary) shared by two Bezier patches, the two rows of control points adjacent to that boundary and Should satisfy G corresponding to k=2,...,p-1 1 Continuity constraints. These constraints form a linear system Ex e =y e , where the vector x e Represents control point and The set of can be obtained by solving the following problem

[0150]

[0151] stEx e =y e

[0152] x e Initial value of It is also constructed by discrete Coons interpolation. The constraints of the above two problems always have solutions, so both problems have global optimal solutions.

[0153] Except for the control points adjacent to the common boundary, the remaining control points of each patch remain the control points obtained by Coons interpolation or by degree raising method. Figure 6 The two figures on the left show the calculated G of the giraffe domain. 1 Smooth parameterization reveals that the resulting mapping is not globally injective (||μ g || ∞ =1.48>1).

[0154] 3) Use QC mapping technology to obtain high-quality parameterization of each patch

[0155] The QC mapping is locally injective, and its angular distortion can be measured by the Beltrami coefficient. The present invention uses this as a framework to generate high-quality parameterization for each face. Here, a quadrilateral face is taken as an example, and the case of a triangle should be treated similarly.

[0156] Since the QC mapping is defined in complex form, first The parameterization of a quadrilateral patch (i.e. a rectangular Bezier surface) is rewritten as a complex function

[0157]

[0158] in and here and It is a control point In this parameterization, the boundary control points and the control points adjacent to the common boundary are known, and the remaining interior control points are determined by solving the following optimization problem

[0159]

[0160] in yes The Beltrami coefficient of , the parameter γ is a positive weight, and Respectively and The Hessian matrix of . In the above model, the first term of the objective function measures the conformal deformation of the mapping, while the second term is related to the uniformity of the parameterization. The constraints ensure the injectivity of the parameterization.

[0161] because and are coupled together. Solving this nonlinear and non-convex problem is not an easy task. Therefore, an auxiliary variable v is introduced and used to replace This leads to an interleaving problem

[0162]

[0163] ||n|| ∞ <1.

[0164] Here we further relax it into such a model

[0165]

[0166] st||n|| ∞ <1.

[0167] The weight η is relatively large. In order to process the model efficiently, an energy descent method is used here to alternately solve and optimize the two sub-problems. The specific iterative process is as follows: 1 Parameterization Start by solving for v alternately. k and

[0168]

[0169] Among them, v k By fixed After solving the model, we get; By fixing v k Then the model is solved. This process continues until the solution of the model satisfies the condition ||v for a given threshold ε k+1 -v k Until ||<ε.

[0170] In each iteration of the entire alternating solution process, there are several types of integration, such as

[0171]

[0172]

[0173] Need to calculate, where f1 is a rational function, given by The partial derivative of The value of

[0174]

[0175] Where f2 is v k The real part of f3 is v kThe imaginary part of . In each iteration, the Bernstein basis functions appearing in the integral remain constant. Inspired by the validity evaluation of the integrals produced in the isogeometric discretization, the values ​​and derivatives of these basis functions are precomputed. Therefore, the integral is evaluated using a lookup strategy. Figure 6 The two figures on the right show the parameterization results of the giraffe domain based on QC mapping, which is much better than Figure 6 Initial results depicted in the two figures on the left.

[0176] Each embodiment in this specification is described in a progressive manner, and the same or similar parts between the embodiments can be referred to each other, and each embodiment focuses on the differences from other embodiments. In particular, for the device embodiment, the above is only a preferred embodiment of the present invention. Since it is basically similar to the method embodiment, the description is relatively simple, and the relevant parts can be referred to the partial description of the method embodiment. The above is only a specific embodiment of the present invention, but the protection scope of the present invention is not limited to this. Any technician familiar with the technical field is within the technical scope disclosed by the present invention. For ordinary technicians in the technical field, changes or replacements that can be easily thought of should be covered within the protection scope of the present invention without departing from the principle of the present invention. Therefore, the protection scope of the present invention should be based on the protection scope of the claims.

Claims

1. A method for plane parameterization suitable for complex geometric areas and having G1 smooth characteristics, characterized by comprising the steps of: 1) Automatic domain decomposition with gradual refinement, using triangular and rectangular Bezier patches to represent complex domains, including: 1.1) Extract the skeleton based on the region boundary, using a sparse optimization-based method; 1.2) Use the skeleton to guide the initial decomposition process, connect the corresponding tangent points on the domain boundary, and automatically generate triangles and rectangular patches to build the basic geometric structure; 1.3) Based on the initial decomposition, further identify the shapes of the remaining patches and select appropriate corner points for each patch; 1.4) Continue to refine the region and use Bezier curves to construct a mathematical model of the resulting patch boundary segment, as follows: 1.4.1) For N discrete points The external boundary segment represented by the chord length parameterization method is used to determine the relevant parameters And use the control point c to fit a p-order Bezier curve C(u), which is done by minimizing get; among them, λ is a positive weight, regularization term Used to measure the smoothness of the resulting curve; 1.4.2) By verifying the mean square error defined by the following formula: Is it greater than 4×10 -8 The fitting accuracy is evaluated with the specified tolerance; if the error value exceeds the specified tolerance, the corresponding patch needs to be further decomposed; 1.4.3) Processing internal boundary segments: For internal boundary straight line segments with two endpoints p1 and p4, interpolation is performed at p1 and p4 to construct a cubic Bezier curve to replace the original line segment; then, two other control points p2 and p3 are constructed according to the tangent directions of p1 and p4; 2) Build a global G 1 Smooth parameterized model: 2.1) Using the discrete Coons interpolation method, the initial multi-slice parameterization is calculated; 2.2) Apply G on the patch interface 1 Continuity conditions, calculated to satisfy the global G 1 Smooth control points; 3) Use QC mapping technology to obtain high-quality parameterization of each patch, as follows: 3.1) Rewrite the parameterization of the lth quadrilateral patch, i.e. a rectangular Bezier surface, into a complex function in and here and It is a control point In this parameterization, the boundary control points and the control points adjacent to the common boundary are known, and the remaining interior control points are determined by solving the following optimization problem: in Yes l The Beltrami coefficient of , the parameter γ is a positive weight, and Respectively represent x l and l The Hessian matrix of 3.2) Introduce an auxiliary variable v and use it to replace This leads to an interleaving problem ||v|| ∞ <1. Here we further relax it into the following model: s.t.||v|| ∞ <1. Among them, the weight η is relatively large; 3.3) Using an energy descent method, the two sub-problems are optimized alternately; 3.4) This process continues until the solution of the model satisfies the condition ||v for a pre-given threshold ε k+1 -v k Until ||<ε.

2. The method for plane parameterization applicable to complex geometric regions and having G1 smoothness according to claim 1, characterized in that: The process of extracting the skeleton in step 1.1) is as follows: A sparse optimization-based approach is used to extract the skeleton and compute a compact medial axis transform of any 2D object represented by a point cloud containing noise or missing data. The skeleton's branch points, branches, the largest inscribed circle centered at each branch point, and the corresponding tangent points on the skeleton are then determined. Assumption c t (u) is the terminal branch spline curve connected to the branch point b, r b is the radius of the largest inscribed circle with b as the center, and delete the terminal branches that meet the following conditions The left part represents the arc length of the terminal branch, r max It is the maximum value of the radius of the circle centered at each branch point, and α and β are positive parameters.

3. The method for plane parameterization applicable to complex geometric regions and having G1 smoothness according to claim 1, characterized in that: The process of obtaining corner points in step 1.3) is as follows: 1.3.1) After completing the preliminary decomposition of the region Ω, further analyze and identify the remaining patches to determine whether there are patches that can be approximated to triangles; 1.3.2) After the triangular patches are identified and distinguished, the remaining patches are classified as quadrilateral patches according to the default rule. Four corner points need to be selected on the domain boundary for each quadrilateral patch, of which two tangent points have been determined to be used as corner points; for the remaining two corner points, an automated method based on the optimal transmission principle is used for selection.

4. The method for plane parameterization applicable to complex geometric regions and having G1 smoothness according to claim 3, characterized in that: Step 1.31) Whether there is a face that can be approximated to a triangle. The specific process is as follows: for A vertex of The boundary (edge )for An edge of , which is inside Ω, is also called the internal boundary; Representation patch The perimeter, L p represents the perimeter of Δpp′p″, and δ>0 is used as a threshold; if the patch satisfy: It is believed that is triangular, in which case point p and the two endpoints p′ and p″ of the corresponding internal boundary are selected as the three corner points required for subsequent parameterization.

5. The method for plane parameterization applicable to complex geometric regions and having G1 smoothness according to claim 1, characterized in that: The specific process of decomposing the corresponding facets described in step 1.42) is as follows: for the longer of the two boundary segments, select a point k1, which has approximately the same geodesic distance to the two end points, and find its nearest point k2 on the opposite boundary; by connecting k1 and k2, decompose this facet into two parts; this decomposition and fitting process will be iterated until the boundary fitting errors of all facets are reduced to below the specified tolerance.

6. The method for plane parameterization applicable to complex geometric regions and having G1 smooth characteristics according to claim 1 or 5, characterized in that: Step 1.43) The specific process is as follows: set up is the tangent vector of two adjacent Bezier curves at p1, and is the tangent vector of two adjacent Bezier curves at p4, and the control points p2 and p3 are constructed as follows: and Among them, ρ and are the degrees of p1 and p4 respectively, and θ is given by and The included angle formed, Is and The included angle formed, σ is the threshold; When the above two equations contradict each other, the following condition is adopted 7. The method for plane parameterization applicable to complex geometric regions and having G1 smoothness according to claim 1, characterized in that: The process of multi-chip parameterization in step 2.1) is as follows: For the l-th quadrilateral patch represented by a rectangular Bezier surface, l = 1, ..., M q , whose boundary control points are Initial control point The construction is: For the triangular Bezier surface S represented by the lth triangular patch l (u,v,w),l=1,…,M t , the calculation formula of the internal control point that satisfies i+j+k=p is:

8. The method for plane parameterization applicable to complex geometric regions and having G1 smoothness according to claim 1, characterized in that: Step 2.2) Apply G to the patch interface 1 The process of continuity conditions is as follows: Let a p-order Bezier curve is the common boundary of two adjacent rectangular Bezier patches, and its control points Then these two patches are G 1 Continuous; if there exist scalar weight functions λ(u) and μ(u) such that: (1-λ(u))s(u)+λ(u)t(u)=(1-μ(u))q l (u)+μ(u)q r (u) holds for all u where and are two rows of control points adjacent to a common boundary; Here, the functions λ(u) and μ(u) are set as linear polynomials and quadratic polynomials as follows: Substituting λ(u) and μ(u) into the previous equation, we get the following G 1 Continuity constraints: They form a linear system; When processing adjacent triangles and rectangular patches, first increase the order of the triangle patch from p to p+1, and then construct the above G 1 Restraint system.

9. The method for plane parameterization applicable to complex geometric regions and having G1 smooth characteristics according to claim 8, characterized in that: In the above G 1 In the constraint system, the boundary control points are known and fixed s0,s p ,t0,t p , and then the remaining control points need to be calculated; The coefficients λ0, λ1 of λ(t) and the coefficients μ0, μ2 of μ(t) are determined by the first and last equations of the system respectively; the coefficient μ1 of μ(t) is set to (μ0+μ2) / 2; G 1 The continuous system is divided into two groups, namely, the vertex G 1 Continuous system (corresponding to k = 0, 1, p, p + 1) and edge G 1 Continuity system (corresponding to k = 2, ..., p-1); for each vertex, singular corner or regular corner in the constructed multi-faceted structure, that is, the intersection of ρ Bezier facets, control vertices o1, o2, ..., o ρ The twist vector at the vertex should satisfy G corresponding to k = 1 or p in the system 1 Continuity constraints; all these constraints form the linear system Vx v =y v , where the vector x v Represents control point The set of is obtained by solving the following problem s.t.Vx v =y v , x v Initial value of It is constructed by discrete Coons interpolation method; For each common edge shared by two Bezier patches, the two rows of control points adjacent to that edge and Should satisfy G corresponding to k=2,...,p-1 1 Continuity constraints; these constraints form a linear system Ex e =y e , where the vector x e Represents control point and The set of is obtained by solving the following problem s.t.Ex e =y e x e Initial value of Constructed via discrete Coons interpolation.

10. The plane parameterization method applicable to complex geometric areas and having G1 smooth characteristics according to claim 1, characterized in that: Step 3.3) The specific iterative process of the energy descent method is as follows: 1 Parameterization Start by solving for v alternately. k and Among them, v k By fixed After solving the model, we get; By fixing v k After solving the model, we get .

Images