# Smoothing Algorithm for Planar and Surface Mesh Based on Element Geometric Deformation

Authors: Shuli Sun; Minglei Zhang; Zhihong Gou
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

Smoothing is one of the basic procedures for improvement of mesh quality. In this paper, a novel and efficient smoothing approach for planar and surface mesh based on element geometric deformation is developed. The presented approach involves two main stages. The first stage is geometric deformation of all the individual elements through a specially designed two-step stretching-shrinking operation (SSO), which is performed by moving the vertices of each element according to a certain rule in order to get better shape of the element. The second stage is to determine the position of each node of the mesh by a weighted average strategy according to quality changes of its adjacent elements. The suggested SSO-based smoothing algorithm works efficiently for triangular mesh and can be naturally expanded to quadrilateral mesh, arbitrary polygonal mesh, and mixed mesh. Combined with quadratic error metric (QEM), this approach may be also applied to improve the quality of surface mesh. The proposed method is simple to program and inherently very suitable for parallelization, especially on graphic processing unit (GPU). Results of numerical experiments demonstrate the effectiveness and potential of this method.


---

## Body

## 1. Introduction

As an elementary technique, mesh quality improvement has been widely used in many applications, such as mesh generation [1–5], mesh simplification [6–8], dynamic grid or mesh deformation [9–14], and other mesh processing procedures. Due to its importance, mesh quality improvement technique has received a lot of attention and great progress has already been made in recent years [15–34]. Mesh improvement approach can be basically divided into two main categories, topological optimization and smoothing (also called geometrical optimization). Topological optimization changes the topology of a mesh, that is, the node-element connectivity relationship, while smoothing or geometrical optimization improves mesh quality by simply moving or adjusting node positions without changing the topology of mesh. This paper will focus on the latter, smoothing.Theoretically the optimization-based smoothing procedure [16, 21, 22, 29] would lead to better results. However, the capacity of improvements is limited when applying such kind of computationally expensive procedure to practical problem. In many works, for example, in [18], smoothing was commonly performed by Laplacian smoothing technique, that is, simply shifting each interior node (free vertex) to the centroid of the surrounding polygon or polyhedron. Although it is computationally inexpensive and easy to implement, however, it may produce occasionally invalid or illegal elements that are unacceptable in numerical analysis. Chen and Xu [25] presented the concept of optimal Delaunay triangulation and developed a quality smoothing scheme for triangular mesh [26] by minimizing the interpolation error among all triangulations with the same number of vertices in the local patch. This optimization-based smoothing method [26–28] could solve the corresponding optimization problem explicitly and thus the computational cost is as low as that of Laplacian smoothing. However, this method is only suitable for triangular mesh and cannot be expanded to other types of mesh.Recently Vartziotis et al. [33, 34] developed a novel smoothing method named GETMe (geometric element transformation method) for triangular and quadrilateral mesh based on element geometric transformation. This method first improves element shape by two-step geometric transformation and then updates positions of moved nodes to achieve the purpose of smoothing. It is not necessary to construct or solve optimization problem, so high efficiency of the method is expected. Test results show that the method achieves excellent performance in effect of smoothing and algorithm efficiency. Unfortunately, it is unsuccessful to expand this approach to polygonal elements with more than four edges (e.g., pentagon and hexagon).Inspired by the work of Vartziotis et al. [33, 34], this paper proposes a new and efficient smoothing approach for planar and surface mesh based on element geometric deformation. The approach includes two main stages. In the first stage, for each individual element, all vertices are moved by a purposely designed two-step stretching-shrinking operation (SSO) in order to get better shape of the element. The second stage is to determine the position of each node of the mesh by a weighted average strategy according to quality changes of its adjacent elements. The proposed SSO-based smoothing approach has the universal scheme and can be conveniently expanded from triangular mesh to quadrilateral mesh, arbitrary polygonal mesh, and mixed mesh. Combined with QEM, this approach may be also applied to improve the quality of surface mesh. Furthermore, the suggested method is simple to program and inherently very suitable for parallelization, especially on GPU.In the rest of this paper we first review the idea and basic procedures of GETMe [33, 34] and then the two-step SSO and the corresponding smoothing approach for planar mesh are presented. Next the smoothing approach for surface mesh is proposed by introducing QEM. Finally, some numerical experiments are given to demonstrate the effectiveness and feasibility of the method.


## 2. Related Work

Recently a novel smoothing approach named GETMe (geometric element transformation method) for triangular and quadrilateral mesh was proposed by Vartziotis et al. [33, 34]. This approach first improves element shape by two-step geometric transformation and then updates positions of moved nodes to achieve the purpose of smoothing. It is not necessary to construct or solve optimization problem, so high efficiency of the method is expected.The method first operates on a counterclockwise oriented triangle with pairwise disjoint verticeszi(0), i∈0,1,2. Taking the edge z0(0)z1(0) of the triangle, for example, rotating it by an angle θ clockwise with z0(0) as the center, then a new vertex z1(1/2) defined as the intersection of two subsidiary lines (the perpendicular of z0(0)z1(0) in z1(0) and the edge z0(0)z1(0) rotated clockwise by angle θ with z0(0) as the center) is obtained by(1)z11/2=z10-wz10-z00=wz00+1-wz10,where w=itan⁡θ and complex numbers are used to simplify notation.Applying the same procedure to the other two edges, a clockwise rotated new triangle with the verticesz0(1/2), z1(1/2), and z2(1/2) is obtained. To avoid the rotational effect, a similar counterclockwise transformation with angle θ by (2)zi1=zi1/2-wzi+1mod31/2-zi1/2=1+wzi1/2-wzi+1mod31/2,i=0,1,2,is applied to get a new triangle with the vertices z0(1), z1(1), and z2(1).It is proved that such geometric element transformations preserve the element centroid, and an arbitrary triangle will gradually converge to an equilateral triangle after enough successive transformations [33].The above transformations usually enlarge the size of the original triangle. In the event that the element does not have any boundary nodes, scaling is applied so that the new element preserves its original perimeter. If an edge of the element belongs to a boundary, the element is scaled according to the length of this edge and is moved so that the nodes belonging to this edge return to their original positions [33].This geometric transformation can be generalized leading to an additional degree of freedom in controlling the speed in which an element becomes regular [34]. For quadrilateral elements, similar transformation can be also conducted [34]. After successively performing such transformation, an arbitrary quadrilateral will become a regular square.In [33], the authors suggested a sequential GETMe approach based on successively improving elements with the lowest quality to complete smoothing of the mesh. Initially, the poor quality elements are rated. During the sequential process, the transformation algorithm is applied to the “worst” element, until the mesh quality reaches the desired level or termination criteria are met. After each transformation, all the affected elements are updated. Alternatively, a simultaneous GETMe approach [34] was presented for smoothing mixed planar meshes. It is based on transforming all elements simultaneously. Updated node positions are obtained as weighted means of the transformed element nodes.Testing shows that the method achieves excellent performance in the effect of smoothing and algorithm efficiency. Unfortunately, it is unsuccessful to expand this approach to polygonal elements with more than four edges (e.g., pentagon and hexagon).


## 3. Smoothing Algorithm for Planar Mesh

A new and universal smoothing approach for planar triangular mesh based on element geometric deformation is first proposed in this section. The presented approach includes two main stages. The first stage is geometric deformation of all the individual triangles. For each triangle, the vertices are moved by a specially designed two-step stretching-shrinking operation (SSO) in order to get better shape of the element. The second stage is to determine the position of each node of the mesh by a weighted average strategy according to the quality changes of its adjacent triangles. Next, the geometric deformation operation for triangle is first proposed and then is followed by its natural extensions to quadrilateral and arbitrary polygon.


### 3.1. Geometric Deformation of Single Triangles

As for the algorithm in the first stage, geometric deformation of each individual triangle can be performed by two operations, stretching and shrinking. The stretching operation will enlarge the size of the original triangle but change the element towards regular shape. In order to preserve its original size, a scaling or shrinking operation must be applied to the enlarged triangle with better shape. The detail of SSO is given as follows.For a counterclockwise oriented triangle with pairwise disjoint verticespi(0)=xi(0),yi(0)T, i∈0,1,2, the stretching operation is performed by pulling each vertex pi(0) out along the perpendicular direction ni of its opposite edge (see Figure 1).Figure 1
Geometric deformation of a triangle (stretching operation).Then we get an enlarged new triangle with the verticespi(1′) calculated by(3)pi1′=pi0+λni,i=0,1,2,where λ can be intuitively understood as a stretching factor to control the degree of geometric stretching. Large value of λ means more significant improvement for element shape or quality. There are many choices to define λ. Here, the following expression in terms of the average length l of the original element and the quality measurement q  (q∈[0,1], q takes the value of 1 for regular element and 0 for degenerated element) is adopted to reflect the fact that only small stretching is needed when the element approaches to regular shape (q gets close to 1):(4)λ=1-ql.In (4), the quality measurement q for triangle can be taken as 2r/R, where r and R represent the radius of incircle and circumcircle of the triangle, respectively.After preforming the stretching operation, the shape or quality of the original poor triangle can be significantly improved. This procedure, however, will usually enlarge the size of the original triangle. In order to preserve its original size, a scaling or shrinking operation shown in Figure2 must be applied to the enlarged triangle with better shape.Figure 2
Geometric deformation of a triangle (shrinking operation).The shrinking operation does not change the shape of the element. In general, there are two ways for shrinking, keeping perimeter or area of the element unchanged. In practice, test results show that there is no big difference between these two ways. In this paper, for the reason of computational cost, we keep the perimeter (the average length) of the element unchanged before and after deformation. In addition, the shrinking operation needs to keep the position of centroid unchanged for each element. Thus, we obtain the shrinking triangle with new vertices of(5)pi1=c0+kpi1′-c1′,i=0,1,2,where c(0)=∑i=02pi(0)/3 and c(1′)=∑i=02pi(1)/3 represent the centroid of the triangle before and after stretching, respectively. The shrinking factor k is defined as the ratio of the average length of the element before and after stretching.Stretching and shrinking operations make up a complete geometric element deformation process. The final shape of the triangle is only determined by the stretching factorλ and the shrinking factor k together, which can be easily calculated by coordinates of the vertices.Figure3 shows the geometric deformation process of a triangle. A poor shape triangle gradually becomes regular after 5 successive cycles of SSO.Figure 3
Geometric deformation process of a triangle.The proposed SSO is numerically efficient and is more simple to implement compared with the GETMe method. More important, the approach has the universal scheme and can be conveniently expanded to quadrilateral elements and arbitrary polygonal elements.


### 3.2. Geometric Deformation of Single Quadrilaterals

The SSO for quadrilaterals is a natural extension of the above procedure. For a counterclockwise oriented quadrilateral with four verticespi(0), i∈0,1,2,3, shown in Figure 4, stretching operation is performed by pulling each vertex pi(0) out along the perpendicular direction ni of the diagonal which links its two adjacent vertices. This operation will usually enlarge the size of the original element; in the meantime, the shape or quality of the original poor element can be significantly improved. Similarly, a shrinking operation is also applied to the enlarged element with better shape to preserve its original size.Figure 4
Geometric deformation of a quadrilateral (stretching operation).The complete SSO is formulated by(6)pi1′=pi0+λni,pi1=c0+kpi1′-c1′,i=0,1,2,3,where λ, k denote the stretching and shrinking factor and c(0) and c(1′) represent the centroid of the element before and after deformation, respectively. Also the stretching factor λ is used to control the degree of geometric stretching. For quadrilateral elements, the stretching factor λ takes the same expression in (4) but replacing q by a quality measurement for quadrilaterals.Figure5 shows the geometric deformation process of a quadrilateral. A poor shape quadrilateral gradually becomes regular after 5 successive cycles of SSO.Figure 5
Geometric deformation process of a quadrilateral.


### 3.3. Geometric Deformation of Arbitrary Polygonal Elements

In accordance with the deformation of quadrilateral elements, the extension of SSO to arbitrary polygonal elements is quite straightforward. For instance, the sketch of stretching deformation for a pentagon is illustrated in Figure6. One can easily derive the deformation formulation similar to (6).Figure 6
Geometric deformation of a pentagon (stretching operation).For such kind of mesh smoothing approach discussed in the following subsection, each element is deformed individually in the first stage, and the only thing we need to do is to perform SSO for each element independently without caring about the connection with its neighbor elements. Therefore, the proposed algorithm has the inherent advantage of dealing with mixed mesh.


### 3.4. SSO-Based Smoothing Algorithm

This subsection gives the algorithm to complete smoothing approach based on SSO. The weighted average strategy is selected in the second stage to determine the position of each node of the mesh according to the quality changes of its adjacent elements. This strategy is aimed at improving the average quality of the whole mesh, while the quality of bad elements can be also improved. Before updating the positions of nodes, each element needs to conduct SSO independently in the first stage to get the temporary vertices (the vertices of the deformed element with better shape). Then the position of nodepi can be updated by (7)pi=∑eΔqi,epi,eT∑eΔqi,e,where pi,eT and Δqi,e denote the corresponding temporary vertex and quality change of the adjacent element e, respectively. According to (7), the element with larger change in quality accounts for larger weight in determining the final position of the node, which leads to more emphasis on the transformed nodes of low quality elements.The main procedure of the SSO-based smoothing algorithm for planar mesh is described as follows.Algorithm 1 (SSO-based smoothing algorithm for planar mesh).
(1)
For each node of the mesh, find all of its adjacent elements.
(2)
For each element, conduct geometric deformation independently by SSO and save the coordinates of all temporary vertices.
(3)
For each node of the mesh, calculate the weighted average of the corresponding temporary vertices by (7) to get the position of this node.
(4)
If the mesh quality reaches the desired level or termination criteria are met, stop; otherwise, return to step(2) for next cycle.The proposed smoothing procedure has a fast convergence speed; usually 10~20 cycles will achieve quite good result. In addition, the main calculations, including SSO for elements and position updating for nodes with the time complexity ofO(n), have very good parallelizability inherently. Therefore, it is suitable for developing parallel scheme of the method, especially on GPU. Actually, a GPU implementation is completed by the authors. Some test results given in Section 5 show remarkable speedup gained by GPU parallel technology.


## 3.1. Geometric Deformation of Single Triangles

As for the algorithm in the first stage, geometric deformation of each individual triangle can be performed by two operations, stretching and shrinking. The stretching operation will enlarge the size of the original triangle but change the element towards regular shape. In order to preserve its original size, a scaling or shrinking operation must be applied to the enlarged triangle with better shape. The detail of SSO is given as follows.For a counterclockwise oriented triangle with pairwise disjoint verticespi(0)=xi(0),yi(0)T, i∈0,1,2, the stretching operation is performed by pulling each vertex pi(0) out along the perpendicular direction ni of its opposite edge (see Figure 1).Figure 1
Geometric deformation of a triangle (stretching operation).Then we get an enlarged new triangle with the verticespi(1′) calculated by(3)pi1′=pi0+λni,i=0,1,2,where λ can be intuitively understood as a stretching factor to control the degree of geometric stretching. Large value of λ means more significant improvement for element shape or quality. There are many choices to define λ. Here, the following expression in terms of the average length l of the original element and the quality measurement q  (q∈[0,1], q takes the value of 1 for regular element and 0 for degenerated element) is adopted to reflect the fact that only small stretching is needed when the element approaches to regular shape (q gets close to 1):(4)λ=1-ql.In (4), the quality measurement q for triangle can be taken as 2r/R, where r and R represent the radius of incircle and circumcircle of the triangle, respectively.After preforming the stretching operation, the shape or quality of the original poor triangle can be significantly improved. This procedure, however, will usually enlarge the size of the original triangle. In order to preserve its original size, a scaling or shrinking operation shown in Figure2 must be applied to the enlarged triangle with better shape.Figure 2
Geometric deformation of a triangle (shrinking operation).The shrinking operation does not change the shape of the element. In general, there are two ways for shrinking, keeping perimeter or area of the element unchanged. In practice, test results show that there is no big difference between these two ways. In this paper, for the reason of computational cost, we keep the perimeter (the average length) of the element unchanged before and after deformation. In addition, the shrinking operation needs to keep the position of centroid unchanged for each element. Thus, we obtain the shrinking triangle with new vertices of(5)pi1=c0+kpi1′-c1′,i=0,1,2,where c(0)=∑i=02pi(0)/3 and c(1′)=∑i=02pi(1)/3 represent the centroid of the triangle before and after stretching, respectively. The shrinking factor k is defined as the ratio of the average length of the element before and after stretching.Stretching and shrinking operations make up a complete geometric element deformation process. The final shape of the triangle is only determined by the stretching factorλ and the shrinking factor k together, which can be easily calculated by coordinates of the vertices.Figure3 shows the geometric deformation process of a triangle. A poor shape triangle gradually becomes regular after 5 successive cycles of SSO.Figure 3
Geometric deformation process of a triangle.The proposed SSO is numerically efficient and is more simple to implement compared with the GETMe method. More important, the approach has the universal scheme and can be conveniently expanded to quadrilateral elements and arbitrary polygonal elements.


## 3.2. Geometric Deformation of Single Quadrilaterals

The SSO for quadrilaterals is a natural extension of the above procedure. For a counterclockwise oriented quadrilateral with four verticespi(0), i∈0,1,2,3, shown in Figure 4, stretching operation is performed by pulling each vertex pi(0) out along the perpendicular direction ni of the diagonal which links its two adjacent vertices. This operation will usually enlarge the size of the original element; in the meantime, the shape or quality of the original poor element can be significantly improved. Similarly, a shrinking operation is also applied to the enlarged element with better shape to preserve its original size.Figure 4
Geometric deformation of a quadrilateral (stretching operation).The complete SSO is formulated by(6)pi1′=pi0+λni,pi1=c0+kpi1′-c1′,i=0,1,2,3,where λ, k denote the stretching and shrinking factor and c(0) and c(1′) represent the centroid of the element before and after deformation, respectively. Also the stretching factor λ is used to control the degree of geometric stretching. For quadrilateral elements, the stretching factor λ takes the same expression in (4) but replacing q by a quality measurement for quadrilaterals.Figure5 shows the geometric deformation process of a quadrilateral. A poor shape quadrilateral gradually becomes regular after 5 successive cycles of SSO.Figure 5
Geometric deformation process of a quadrilateral.


## 3.3. Geometric Deformation of Arbitrary Polygonal Elements

In accordance with the deformation of quadrilateral elements, the extension of SSO to arbitrary polygonal elements is quite straightforward. For instance, the sketch of stretching deformation for a pentagon is illustrated in Figure6. One can easily derive the deformation formulation similar to (6).Figure 6
Geometric deformation of a pentagon (stretching operation).For such kind of mesh smoothing approach discussed in the following subsection, each element is deformed individually in the first stage, and the only thing we need to do is to perform SSO for each element independently without caring about the connection with its neighbor elements. Therefore, the proposed algorithm has the inherent advantage of dealing with mixed mesh.


## 3.4. SSO-Based Smoothing Algorithm

This subsection gives the algorithm to complete smoothing approach based on SSO. The weighted average strategy is selected in the second stage to determine the position of each node of the mesh according to the quality changes of its adjacent elements. This strategy is aimed at improving the average quality of the whole mesh, while the quality of bad elements can be also improved. Before updating the positions of nodes, each element needs to conduct SSO independently in the first stage to get the temporary vertices (the vertices of the deformed element with better shape). Then the position of nodepi can be updated by (7)pi=∑eΔqi,epi,eT∑eΔqi,e,where pi,eT and Δqi,e denote the corresponding temporary vertex and quality change of the adjacent element e, respectively. According to (7), the element with larger change in quality accounts for larger weight in determining the final position of the node, which leads to more emphasis on the transformed nodes of low quality elements.The main procedure of the SSO-based smoothing algorithm for planar mesh is described as follows.Algorithm 1 (SSO-based smoothing algorithm for planar mesh).
(1)
For each node of the mesh, find all of its adjacent elements.
(2)
For each element, conduct geometric deformation independently by SSO and save the coordinates of all temporary vertices.
(3)
For each node of the mesh, calculate the weighted average of the corresponding temporary vertices by (7) to get the position of this node.
(4)
If the mesh quality reaches the desired level or termination criteria are met, stop; otherwise, return to step(2) for next cycle.The proposed smoothing procedure has a fast convergence speed; usually 10~20 cycles will achieve quite good result. In addition, the main calculations, including SSO for elements and position updating for nodes with the time complexity ofO(n), have very good parallelizability inherently. Therefore, it is suitable for developing parallel scheme of the method, especially on GPU. Actually, a GPU implementation is completed by the authors. Some test results given in Section 5 show remarkable speedup gained by GPU parallel technology.


## 4. Smoothing Algorithm for Surface Mesh

The SSO-based smoothing algorithm can be also applied to surface mesh. Different from planar mesh, we must consider an important issue; that is, the movement of the node is not allowed to deviate from the original surface too much. In this section, we introduce the quadratic error metric (QEM) [6] which has been widely used in surface simplification to control the moving distance of the nodes away from the original surface in the process of geometric deformation.


### 4.1. The Basic Idea of QEM [6]

It is easy to observe that, in original surface mesh, each nodepi is the solution of the intersection of a set of planes, S(pi). Then the QEM of node pi, Δ(pi), is defined as the sum of squared distances to its planes:(8)Kp=a2abacadabb2bcbdacbcc2cdadbdcdd2,Qi=∑p∈SpiKp,Δpi=piTQipi.In (8)  a, b, c, d represent the coefficients of the plane defined by the equation ax+by+cz+d=0 with a2+b2+c2=1; pi=xi,yi,zi,1T denotes the augmented coordinate vector of node pi.The QEM can be used to describe the deviation degree of a node away from the original surface. Note that the initial value of QEM for each node is 0, since each node in original mesh lies in the planes of all its incident elements. In smoothing of surface mesh, we need to control the value of QEM to ensure that the element geometric deformation is acceptable.


### 4.2. SSO-QEM Based Smoothing for Surface Mesh

Before performing element geometric deformation to surface mesh, the quadratic error matrixQi in (8) for each node must be calculated in advance. For each element, construct the plane according to the coordinates of its vertices and find the coefficients of a, b, c, d, and then calculate the matrix Kp for each element according to (8). By using the relationship between nodes and faces of the mesh, for each node, find all of its adjacent elements, and then add the matrix Kp of its adjacent elements up to get the error matrix Qi. With the aid of the quadratic error matrix Qi, we can easily calculate the error of the node pi by (8) when the node is moved during smoothing process. Please note that the quadratic error matrix Qi only needs to be calculated once.The SSO for element geometric deformation needs to be revised to meet the requirement of surface mesh. For simplicity, we take triangular surface mesh as an example. Figure7 is a local region of a triangular surface mesh, and the element p0(0)p1(0)p2(0) will be deformed by the geometric transformation. The deformation is restrained in the original plane during the stretching and shrinking operations; namely, the position of each vertex only changes in the original plane of the element.Figure 7
Element stretching for triangular surface mesh.The stretching operation is exactly the same as that in planar case. The stretching factorλ is taken as the same expression in (4). The three vertices of the element, respectively, move along the direction normal to their opposite edge. The new positions of pi(1′) can be calculated by (3).In the scaling or shrinking process, the basic principle remains to keep the centroid of the element unchanged as that in planar case, but the shrinking factork is no longer the same in (5). Representing the three vertices after shrinking with p0(1), p1(1), and p2(1), the total quadratic error due to element deformation can be formulated as (9)Δe=∑i=02pi1TQipi1.Recalling the SSO formulation in (3) and (5), one can see that the total quadratic error Δe is only decided by the stretching factor λ and shrinking factor k; that is, (10)Δe=Eλ,k.Since the stretching factor λ is given in terms of the average length l of the original element and the quality measurement q, the total quadratic error Δe of the element is only the function of shrinking factor k. The remaining thing we need to do is to find the optimal value of k by minimizing the total quadratic error Δe, which can be easily solved by analytical or numerical manner. With the optimal value of k obtained, the new positions of pi(1) can be calculated by (5).The main procedure of so-called SSO-QEM based smoothing algorithm for surface mesh is summarized as follows.Algorithm 2 (SSO-QEM based smoothing algorithm for surface mesh).
(1)
For each element, formulate the equation of the plane according to the coordinates of its vertices, and calculate the matrixKp according to (8).
(2)
For each node of the mesh, find all of its adjacent elements; then add the matrixKp of its adjacent elements up to get the error matrix Qi.
(3)
For each element, project the coordinates of its vertices to the plane, and then conduct geometric deformation independently by stretching operation; find the optimal value of the shrinking factork by minimizing the total quadratic error Δe; conduct shrinking operation independently and save the coordinates of all temporary vertices.
(4)
For each node of the mesh, calculate the weighted average of the corresponding temporary vertices by (7) to get the position of this node.
(5)
If the mesh quality reaches the desired level or termination criteria are met, stop; otherwise, return to step(3) for next cycle.With the aid of QEM to control the moving distance of the nodes away from the original surface, the proposed SSO-QEM based smoothing algorithm is able to improve the quality of surface mesh while maintaining the shape of surface as much as possible.


## 4.1. The Basic Idea of QEM [6]

It is easy to observe that, in original surface mesh, each nodepi is the solution of the intersection of a set of planes, S(pi). Then the QEM of node pi, Δ(pi), is defined as the sum of squared distances to its planes:(8)Kp=a2abacadabb2bcbdacbcc2cdadbdcdd2,Qi=∑p∈SpiKp,Δpi=piTQipi.In (8)  a, b, c, d represent the coefficients of the plane defined by the equation ax+by+cz+d=0 with a2+b2+c2=1; pi=xi,yi,zi,1T denotes the augmented coordinate vector of node pi.The QEM can be used to describe the deviation degree of a node away from the original surface. Note that the initial value of QEM for each node is 0, since each node in original mesh lies in the planes of all its incident elements. In smoothing of surface mesh, we need to control the value of QEM to ensure that the element geometric deformation is acceptable.


## 4.2. SSO-QEM Based Smoothing for Surface Mesh

Before performing element geometric deformation to surface mesh, the quadratic error matrixQi in (8) for each node must be calculated in advance. For each element, construct the plane according to the coordinates of its vertices and find the coefficients of a, b, c, d, and then calculate the matrix Kp for each element according to (8). By using the relationship between nodes and faces of the mesh, for each node, find all of its adjacent elements, and then add the matrix Kp of its adjacent elements up to get the error matrix Qi. With the aid of the quadratic error matrix Qi, we can easily calculate the error of the node pi by (8) when the node is moved during smoothing process. Please note that the quadratic error matrix Qi only needs to be calculated once.The SSO for element geometric deformation needs to be revised to meet the requirement of surface mesh. For simplicity, we take triangular surface mesh as an example. Figure7 is a local region of a triangular surface mesh, and the element p0(0)p1(0)p2(0) will be deformed by the geometric transformation. The deformation is restrained in the original plane during the stretching and shrinking operations; namely, the position of each vertex only changes in the original plane of the element.Figure 7
Element stretching for triangular surface mesh.The stretching operation is exactly the same as that in planar case. The stretching factorλ is taken as the same expression in (4). The three vertices of the element, respectively, move along the direction normal to their opposite edge. The new positions of pi(1′) can be calculated by (3).In the scaling or shrinking process, the basic principle remains to keep the centroid of the element unchanged as that in planar case, but the shrinking factork is no longer the same in (5). Representing the three vertices after shrinking with p0(1), p1(1), and p2(1), the total quadratic error due to element deformation can be formulated as (9)Δe=∑i=02pi1TQipi1.Recalling the SSO formulation in (3) and (5), one can see that the total quadratic error Δe is only decided by the stretching factor λ and shrinking factor k; that is, (10)Δe=Eλ,k.Since the stretching factor λ is given in terms of the average length l of the original element and the quality measurement q, the total quadratic error Δe of the element is only the function of shrinking factor k. The remaining thing we need to do is to find the optimal value of k by minimizing the total quadratic error Δe, which can be easily solved by analytical or numerical manner. With the optimal value of k obtained, the new positions of pi(1) can be calculated by (5).The main procedure of so-called SSO-QEM based smoothing algorithm for surface mesh is summarized as follows.Algorithm 2 (SSO-QEM based smoothing algorithm for surface mesh).
(1)
For each element, formulate the equation of the plane according to the coordinates of its vertices, and calculate the matrixKp according to (8).
(2)
For each node of the mesh, find all of its adjacent elements; then add the matrixKp of its adjacent elements up to get the error matrix Qi.
(3)
For each element, project the coordinates of its vertices to the plane, and then conduct geometric deformation independently by stretching operation; find the optimal value of the shrinking factork by minimizing the total quadratic error Δe; conduct shrinking operation independently and save the coordinates of all temporary vertices.
(4)
For each node of the mesh, calculate the weighted average of the corresponding temporary vertices by (7) to get the position of this node.
(5)
If the mesh quality reaches the desired level or termination criteria are met, stop; otherwise, return to step(3) for next cycle.With the aid of QEM to control the moving distance of the nodes away from the original surface, the proposed SSO-QEM based smoothing algorithm is able to improve the quality of surface mesh while maintaining the shape of surface as much as possible.


## 5. Experimental Results

Several examples are selected to demonstrate the effectiveness and efficiency of the proposed approach. Both examples are tested on Windows 7 system with Intel (R) Xeon (R) CPU (2.67 GHz) and 24 GB memory.First a triangular mesh with 7229 nodes and 13942 elements is tested (see Figure8). The quality of initial mesh is very poor. Several smoothing methods, including Laplacian method, GETMe method, and the proposed SSO-based method, are conducted to improve the quality of mesh. The experimental results are listed in Table 1. From the results, it can be seen that the proposed SSO-based method has slight advantages over GETMe method in average quality and gives much better result in quality of the worst element. In this test, Laplacian method even gives negative value of the quality, which means invalid or illegal elements are produced.Table 1
Quality statistics of triangular mesh by different methods.
Methods
Average quality
Quality of theworst element
Initial mesh
0.77
0.002
Laplacian
0.95
−0.860
GETMe
0.94
0.004
SSO
0.96
0.023Figure 8
A triangular mesh before and after smoothing.
(a)
Initial mesh
(b)
Smoothed mesh by suggested methodThe second example is a planar quadrilateral mesh shown in Figure9, with the boundary nodes fixed. The mesh contains 1089 nodes and 1024 elements, and the initial quality is also very poor. The experimental results are listed in Table 2. In this test the performances of the proposed method and GETMe method are almost the same. Both of them give better results in quality of the worst element compared with Laplacian method.Table 2
Quality statistics of quadrilateral mesh by different methods.
Methods
Average quality
Quality of theworst element
Initial mesh
0.37
0.03
Laplacian
0.74
0.08
GETMe
0.73
0.11
SSO
0.74
0.10Figure 9
A quadrilateral mesh before and after smoothing.
(a)
Initial mesh
(b)
Smoothed mesh by suggested methodIn the following two surface meshes of dragon and bunny [35] (Figure 10) are tested by the proposed SSO-QEM based smoothing approach. The dragon mesh consists of 50000 nodes and 19994 triangular surface elements, while the bunny mesh contains 34835 nodes and 7996 triangular surface elements. Table 3 lists the smoothing results. It can be seen that the proposed SSO-QEM algorithm has obviously improved the quality of each mesh while the surface is well maintained. Here, the geometric error in Table 3, defined as the ratio of Hausdorff distance [7] and the average size of the mesh, is adopted to measure the difference between the new and original surface meshes.Table 3
Smoothing results of surface meshes by the proposed method.
Mesh
Initial quality
Qualityafter smoothing
Running time (s)
Geometric error
Dragon
0.58
0.74
3.66
0.8696
Bunny
0.63
0.77
1.41
0.4022Figure 10
Surface mesh before and after smoothing.
(a)
Initial mesh
(b)
Smoothed mesh by suggested methodFinally, an example for testing the performance of GPU implementation of the proposed method is given. The code is realized on the platform CUDA 4.0. The type of GPU used for this test is Tesla C2050 (1.15 GHz) with 3034 MB memory. The testing planar mesh consists of 445182 nodes and 887622 triangles. In order to compare the efficiency, 50 cycles or iterations are run on CPU and GPU, respectively (actually it does not need that much of iterations to reach convergence). Table4 gives the comparison of smoothing result and running time. More than 300 times speedup is gained by GPU implementation which displays great superiority and potential of the proposed method in the aspect of computational efficiency.Table 4
Comparison of smoothing result and efficiency on CPU and GPU.
Methods
Average quality
Quality of theworst element
Running time (s)
Initial
0.876
0.246
—
Smoothing on CPU
0.952
0.504
276.05
Smoothing on GPU
0.955
0.504
0.8150


## 6. Conclusions

This paper presents a novel and efficient smoothing algorithm for 2D planar and surface mesh based on element geometric deformation. The main features of the proposed method are summarized as follows.(1)
The smoothing of mesh involves two main stages, namely, geometric deformation of all the elements and updating positions of all the nodes. These two stages can be performed independently.(2)
The element geometric deformation consists of two operations: stretching and shrinking. The stretching operation changes the shape and improves the quality of the element, while the shrinking operation preserves the centroid and size of the element. These operations mainly involve two parameters, the stretching factor and shrinking factor. As for the algorithm in this paper, the stretching factor is decided by the size and quality of the element, while the shrinking factor for planar mesh is defined as the ratio of the average length of the element before and after stretching. Single elements will converge to regular shape after successively geometric deformations.(3)
In the stage of node updating, the final position of each node is calculated by taking weighted average of the corresponding temporary nodes obtained from element stretching and shrinking operations. The element with larger change in quality accounts for larger weight in determining the final position of the node, which leads to more emphasis on the transformed nodes of low quality elements.(4)
The proposed SSO-based smoothing approach has the universal scheme and can be conveniently expanded from triangular mesh to quadrilateral mesh, arbitrary polygonal mesh, and mixed mesh. Moreover, each element is deformed individually in the geometric deformation stage without caring about the connection with its neighbor elements; therefore, the proposed algorithm has the inherent advantage of dealing with mixed mesh.(5)
The main computation efforts, including SSO for elements and position updating for nodes, have very good parallelizability inherently. It is very suitable for developing parallel scheme of the proposed method, especially on GPU. Some test results show remarkable speedup gained by GPU parallel technology, which is very important for dealing with large scale meshes.(6)
For smoothing surface mesh, quadratic error metric (QEM) is introduced to element geometric deformation to control the moving distance of the nodes away from the original surface. The SSO for element geometric deformation only needs slight revision by shifting the shrinking factor, which can be solved by minimizing the total quadratic error due to element deformation.To conclude, the proposed method has the universal scheme and very simple formulation. It works efficiently for smoothing triangular mesh, quadrilateral mesh, arbitrary polygonal mesh, and mixed mesh. It can be also applied to improve quality of surface mesh by incorporating with QEM. Moreover, it is very suitable for developing parallel implementation. Results of numerical experiments demonstrate the effectiveness and potential of this method. Further study will focus on extensions of the method to mixed volume meshes.



---

*Source: 435648-2015-05-24.xml*