A Triangle Mesh Standardization Method Based on Particle Swarm Optimization

To enhance the triangle quality of a reconstructed triangle mesh, a novel triangle mesh standardization method based on particle swarm optimization (PSO) is proposed. First, each vertex of the mesh and its first order vertices are fitted to a cubic curve surface by using least square method. Additionally, based on the condition that the local fitted surface is the searching region of PSO and the best average quality of the local triangles is the goal, the vertex position of the mesh is regulated. Finally, the threshold of the normal angle between the original vertex and regulated vertex is used to determine whether the vertex needs to be adjusted to preserve the detailed features of the mesh. Compared with existing methods, experimental results show that the proposed method can effectively improve the triangle quality of the mesh while preserving the geometric features and details of the original mesh.


Introduction
Triangle mesh reconstruction based on physical measurement data is one of important researches in reverse engineering. There are often some non-standardization mesh cells in a reconstructed triangle mesh because of the complexities of the physical measurement data and the shortages of a mesh reconstruction algorithm, such as long-narrow triangles, holes, redundancy triangles and T-vertex triangles. When the finite element analysis is carried out on the physical model, the non-standardization mesh directly influences the precision, efficiency and convergence of the calculation result. In addition, the efficiency and stability of rapid prototyping manufacturing are also reduced because of the existence of the non-standardization mesh. Therefore, to obtain a high-quality mesh, it is necessary to standardize the mesh reconstructed by industrial computed tomography (ICT) images.
Because it is relatively easy to standardize a redundancy triangle and a T-vertex triangle, existing research on mesh standardization is based on two aspects, hole repairing and long-narrow triangle processing.
Hole repairing standardizes the mesh by filling the holes of the reconstructed triangle mesh. There are two main categories for the hole repairing method, repairing based on the mesh and repairing based on the volume data. Repairing based on the mesh partially fills the holes of the mesh without changing the original mesh topology, such as the algorithms of Zhao [1], Attene [2], Jun [3] and Panchetti [4]. Repairing based on the volume data [5][6][7] is divided into the following 3 steps: (1) converting the mesh reconstructed to the volume mesh; (2) filling the holes; and (3) reconstructing the triangle mesh by the extracting iso-surface for the volume mesh. The central contributions of the above-mentioned methods are triangulation of the holes and making a new triangle coordinate with the triangle area around the hole. Hence, the quality of the repaired triangle is not a serious consideration.
The standardization of a long-narrow triangle aims to achieve a regular triangle as much as possible by relocating the position of the triangle vertices. The Laplace standardization algorithm [8,9] is popular and the most widely employed, and its key idea is to adjust the vertex position at a certain rate along the direction of the Laplace operator. The algorithm generally works quite well for improving the quality of the triangles. However, it can result in some shrinkage in the size of the mesh. Taubin [10] constructed a low-pass filter to effectively restrain the noise and mesh shrinkage, but it brought some disturbances. Vollmer [11] employed volume normalization for inhibiting the mesh shrinkage, yet the quality of the resulting mesh was general. Duguet [12] proposed an optimization algorithm combining the advantages of bilateral filtering and image filtering, and good mesh standardization was obtained. However, the error is bigger on the vertex curvature calculation. Chen [13] modified bilateral filtering and combined a quasi-Laplacian operator to relocate the vertex position. The algorithm can effectively make a uniform mesh while preserving the geometric detailed features of the original mesh. However, many parameters need to be set in the algorithm. There are many other methods for the long-narrow triangle standardization, such as those proposed by Zhu [14], Zheng [15], Wei [16], Gao [17], etc. These approaches of standardization improved the triangle quality of the mesh to some extent. However, the main considerations of the several approaches above preserved the geometric feature, reduced the volume shrinkage and improved the efficiency of the algorithm.
Each of the above-mentioned methods can only obtain good results in a specific aspect of mesh standardization. Some methods for mesh standardization are employed to obtain the uniform mesh, which can avoid the errors of finite element analysis. Some methods for mesh standardization are applied to get smaller error between standardized mesh and original mesh, which can reduce the reconstructed errors of CAD model. Some methods for mesh standardization are used to preserve the geometric features and details of the original mesh, which can improve the precision of finite element analysis. Some methods for mesh standardization are adopted to obtain the smoothed mesh, which can avoid stress concentration. However, methods described above cannot effectively improve the triangle quality of the mesh while preserving the geometric features and details of the original mesh, reducing the errors between standardized mesh and original mesh.
Because of the shortcoming of above-mentioned methods for mesh standardization, in this paper, we present a novel method for mesh standardization. First, a local curve surface for the vertex and its first-order adjacent vertices of the mesh is fitted using the least square method. Additionally, to improve the triangle quality of the mesh, we modify the PSO and apply it to achieve the optimum vertex on the fitted surface. Finally, threshold of the normal angle is used to determine whether the vertex needs to be regulated. The presented method can effectively improve the triangle quality of the mesh while preserving the geometric features and details of the original mesh.
The remainder of this paper is organized as follows. Section 2 describes a cubic curve surface fitting method for the vertex of the mesh and its adjacent vertices. Section 3 presents our PSO algorithm. Section 4 discusses the proposed method for mesh standardization. Section 5 shows some experimental results and analysis. The final section presents our conclusion.

Local Curve Surface Fitting of the Mesh
Although the reconstructed triangle mesh is discrete, the mesh can be considered as continuous from a local perspective. Therefore, any one of the local meshes can be regarded as a local curve surface, and the vertices of the mesh are distributed on the surface. Hence, the standardization of the triangle is a procedure in which the triangle vertices of the mesh are regulated on the local surface.
References [18,19] provided an approach in which a vertex's local surface is denoted by a cubic polynomial that fits the vertex of the mesh and its first-order adjacent vertices. Inspired by the method, in this paper, we employ the least squares to fit a cubic surface for the vertex and its adjacent vertices of the triangle mesh. To avoid geometry deformation and volume shrinkage of the reconstructed triangle mesh, the fitted surface is used for the constraint to the vertex adjustment. We regard vertex V i as the origin of the coordinates. Consider the normal n i of V i as the positive of the z axis, and create a local coordinate system to fit vertex V i and its first-order adjacent vertices N(V i ) into a cubic surface. The fitting curve surface is defined under the new coordinate system as: The specific steps for solving the coefficient of the local surface equation (Eq 1) are as follows: Step 1 Solve the partial derivative of the curve surface equation f(x,y) and obtain the normal equation, N(x,y) of f(x,y): Step 2 Set the normal vector of point (x i ,y i ,z i ) as (a i ,b i ,c i ). The normal (a i ,b i ,c i ) is rewritten as (−a i /c i ,−b i /c i ,−1) in to be consistent with the form of Eq 2.
Step 3 Write those coefficients of Eq 1 as a column vector, that is: Step 4 List the following equation for each point (x i ,y i ,z i ) under the local coordinate system: Step 5 Combine Eqs 4, 5 and 6, and form a linear system of equations, that is: is a column vector, and n represents the number of vertices N(V i ). If n is less than 4, the nearer second-order adjacent vertices to V i can be used.
Step 6 Solve Eq 7, and obtain a feasible solution of the least square method, that is, the coefficients of the fitting surface equation f(x,y).
In mesh standardization of this paper, the fitted surface is used for the region constraint of PSO.

Introduction of standard PSO
PSO [20,21] has been widely used in all types of engineering optimization problems because it is simple, easy to realize and has fewer adjustable parameters. The key idea of PSO is to find the optimal solution according to collaboration and information sharing between individuals in the group. The mathematical principle of PSO can be described as follows. In solution space, N particles present N possible solutions. The moving process of the particles is the searching process of the solution, rate of the particle stands for the searching direction, and each particle can determine the speed and position according to their own previous information and social information. The iterative formula for the speed and position of each particle can be defined as: The constraint of the particle speed can be presented as: The individual previous optimal position of a particle can be calculated according to Eq 11: where P t i and P t g are defined as: where V tþ1 i and X tþ1 i represent the speed and position of the i-th particle carried out t+1 times, respectively; V t i and X t i stand for the speed and position of the i-th particle carried out t times, respectively; P t i and P tþ1 i denote the individual previous optimal solution of the i-th and the (i+1)-th iteration, respectively; P t g is the group optimal solution after t times; c 1 and c 2 are learning factors; m is the population size; ω is a weight factor; r 1 2 (0,1) and r 2 2 (0,1) are the random distributions; V max indicates the maximum limit speed of the particle; and f(Á) is a fitness function.

Modified PSO
Since PSO was proposed, it gains significant popularity and improvement. Shi et al. [22] introduced the inertia weight factor to constrict the velocity and better regulate the capability of search. Clerc et al. [23] presented an improved PSO with a constraint factor by analyzing the speed iterative equation of classical PSO. Clerc et al. found that introduction of constraint factor can prevent oscillation of particles, ensure the convergence of the algorithm, and analyzed it from the algebraic point of view (discrete time) and the analytical view (continuous time) respectively. Nowadays, the PSO with constraint factor becomes the canonical PSO algorithm because of good performance. Mendes et al. [24] proposed a fully informed PSO (FI-PSO) to improve the optimization performance. Liu et al. [25] addressed a novel PSO with the scalefree topologies and gave rise to a better balance between the convergence speed and the optimum quality. Motivated by the shortcoming of FI-PSO, Du et al. [26] proposed a PSO with limited information (LI-PSO), which provides adequate information of each particle yet avoids the waste of information. Liang et al. [27] presented a comprehensive learning PSO (CLPSO), which employs all other particles' historical best information to update a particle's velocity. This strategy preserves the diversity of the swarm and discourages premature convergence. Besides, there have been many improved PSO which focus on the model coefficients and the population structure [28][29][30][31][32][33].
However, Most of the existing PSO algorithms use the inertia factor depending on iterations, which thus ignores the effect of the error between fitness and the optimal solution. Besides, in most of PSO algorithms, each particle is usually influenced by the best particle among its neighborhood, which thus neglects some useful information from other neighbors. Promoted by above works, especially FI-PSO, LI-PSO, CLPSO and canonical PSO, we modify the PSO by introducing the center of particles and modifying the inertia factor.
First, standard PSO often falls into the local optimum, especially in high order optimization. To avoid the shortcomings of standard PSO, inspired by references [24,26], we append a disturbance to the speed iteration formula (8) in this paper. The center of the current iterative particles is regarded as the disturbance, which indicates the influence of the most particles. The disturbance can effectively increase the information for particle searching optimization and disperse the appeal of the excellent particle to the whole particle swarm. The center of particles can be calculated by Eq 14.
where P t c is the center of the particle swarm in the t-th iteration, X t i is the position of i-th particle carrying out the t-th iteration, f(Á) is a fitness function, Med(Á) is a median filtering function, and m represents the population size.
Second, to adjust the searching region in real-time, we modify the inertia factor of PSO, which reflects the extent of inheritance to the speed of the particle itself. The search capability of PSO can be dynamically regulated by changing the value of ω. If ω is large, the PSO has strong global search capability; otherwise, the PSO has strong local search ability [22]. Linearly decreasing inertia weight (LDIW), as proposed by Shi [22,34], and its variants are widely employed at present. However, the value of ω depends only on the iterations under the strategy of LDIW, and the searching region of PSO cannot be adjusted according to the result of each calculation in real-time. To achieve a balance between global search and local search in PSO, in this paper, we apply the error between particle fitness and the optimal solution of the objective function to define an inertia factor, o t i which can be adaptively regulated, that is: where ω min and ω max indicate the minimum and maximum of the inertia factor, respectively. After many experiments and calculations, Shi [34] found that the convergence speed and precision of the PSO are higher when ω min = 0.4 and ω max = 0.9. In this paper, the values of these two parameters are the same as the verified values by Shi. E(x) max and E(x) min represent, respectively, the maximum and minimum of the errors between the fitness of each particle and the optimal solution of the objective function in the process of iteration. E t i ðxÞ denotes the error between the fitness of the current particle carrying out the t-th iteration and the optimal solution of the objective function. Eq 15 shows that the value of the inertia factor, o t i , is proportional to the distance between the particle's current position and the objective position. The global search ability and local search ability of the algorithm can be adjusted in real time in terms of the current position of the particle. Then, the solving precision of PSO is improved.
Finally, inspired by reference [23], to speed up the convergence rate, reduce the particle oscillation amplitude and avoid ineffective iteration, we introduce a constraint factor ξ. Because of the center of particles is added to the speed iterative equation described by Eq 17, the constraint factor is redefined as: where φ = c 1 + c 2 + c 3 and φ > 6, c 1 , c 2 and c 3 are learning factors.
To sum up, PSO is modified by introducing the disturbance of the particles' center, the constraint factor and the adaptive inertia factor. The speed formula and position formula of the modified PSO can be rewritten as: where P t c is the center of the particle swarm in the t-th iteration, c 3 is a center learning factor, and r 3 2 (0,1) is a random distribution.
The performance comparisons among the modified PSO, classical PSO [22], canonical PSO [23] and comprehensive learning PSO [27] are particularly described in S1 Appendix.
In mesh standardization of this paper, the modified PSO is used to regulate the vertices of triangles on the fitted local curve surface and more triangles with good quality can be obtained. Here, the particles of PSO stand for various possible positions of the target vertex. Movements of the particles represent optimization process of the target vertex. And the best particle obtained by the PSO is the final position of the target vertex.

Our Approach for Mesh Standardization
To quantitatively analyze the triangular standardization of the mesh, we have to determine an evaluation standard. In this paper, the triangle average quality [35] of the local mesh is used for the standard to evaluate the standardization of the mesh. The triangle average quality of the local mesh is defined as: T sum represents the number of triangles in the mesh; Q k stands for the k-th triangle quality of the mesh; S k denotes the k-th triangle area of the mesh; l k1 , l k2 , and l k3 indicate the edge lengths of the k-th triangle; and Q k 2 (0,1) describes the quality of a triangle. When Q k = 1, the triangle is an equilateral triangle, whereas for smaller values of Q k , the triangle is longer-narrower in shape.
Therefore, the objective function of the proposed algorithm can be defined as: where N V i is a first-order neighborhood triangle set of vertex V i . The constraints of vertex V i (x i ,y i ,z i ) are: Dy ¼ maxfjy i j À jy V i 1 j; jy i j À jy V i 2 j; Á Á Á ; jy i j À jy V i V sum jg ð23Þ Eq 26 is a cubic surface equation, which is fitted by vertex V i and its first-order adjacent vertices. The fitted surface is set as the constraint of V i to ensure that the geometry size of the mesh does not shrink. In addition, to better preserve the detailed features of the mesh, the normal angle between the original vertex and the regulated vertex is used to determine whether the vertex needs to be adjusted when we have obtained the new position of the vertex using the improved PSO.
According to the above analyses, the presented approach about triangle mesh standardization can be implemented by the following steps.
Step 1 Use Eq 1 to fit a local cubic surface for vertex V i and its first-order adjacent vertices at the first.
Step 2 Initialize the particles, and set the current position of the particles as individual previous optimal values. Determine the group's optimal value using Eq 13, and calculate the center of particles according to Eq 14.
Step 3 If the constraints of Eqs 24, 25 and 26 are satisfied, update the speed and position of each particle according to Eqs 17 and 18. Otherwise, x or y coordinate of the particle is assigned the boundary value according to Eqs 24 or 25, the opposite direction of the speed is set, and then calculate z coordinate of the particle according to Eq 26, and then update speed and position of the particle according to Eqs 17 and 18.
Step 4 Calculate the fitness of each particle using Eq 21. Update the individual previous optimal positions and the group's previous optimal position according to Eqs 12 and 13, and recalculate the center of particles using Eq 14.
Step 5 If the end of the iteration or the error threshold meets the condition, go to the next step; otherwise, jump to step 3.
Step 6 Store the new coordinate value of the vertex, and calculate and store the normal angle / n i between the original vertex and the regulated vertex.
Step 7 If / n i is within the threshold range, update the vertex position and vertex normal vector information; otherwise, maintain the original position of the vertex.
Step 8 If all vertices of mesh have been traveled, finish the mesh standardization; otherwise, go to step 1.

Experimental Results and Discussions
The proposed algorithm is implemented in Visual C++6.0 and tested on a PC with an Intel Core2 2.2G Processor and 3G RAM. The parameter settings of the algorithm are shown in Table 1. We compare the standardization performance of the proposed method with several traditional methods, such as Laplace, Taubin, Vollmer and the state-of-the-art method of Chen, by making use of the triangle mesh reconstructed by ICT serial images. The values of the parameter of the compared methods are the same as the original references.
When the error of adjacent iterations is set to 1.0e-4, we can find the optimal position of the vertex after about 40 iterations according to the experimental results. Hence, the iterations of the proposed method can be set to 50. The population size, m, of PSO is generally set as 20~40, but it is worth noting that a larger value of m gives a lower efficiency of PSO. Furthermore, m is set as 11 because of the smaller region of the fitted surface. Here, m is set as an odd number to conveniently calculate the center of the particles. c 1 is an individual learning factor denoting the cognitive ability of the particle itself; c 2 is a global learning factor, which is used to indicate the information sharing cooperation ability of all particles; and c 3 is a central learning factor standing for the balance ability of the center of the particles. Here, c 1 = 2.05 and c 2 = 2.05 are set empirically. c 3 is also set as 2.05 to balance the influence among the center of the particles,  the individual optimum and the global optimum. Reference [36] showed that the smaller the error threshold of the normal angle is, the better the preserving feature of mesh is and the worse the triangle quality of the mesh is. The error of the angle is commonly set as 5°*15°, so we set the error threshold of the normal angle as 8°to strike a balance between the features and the triangle quality of mesh.
The experimental results are shown in Figs 1-3. The comparisons of the standardization algorithms on the wheel hub model are shown in Fig 1. The quality of the mesh is better after Laplace standardization. However, some features of the mesh are lost, and the mesh is distorted (see Fig 1B). To some extent, the standardization results by Taubin and Vollmer preserve the feature of the mesh, whereas the triangle quality of the mesh is worse (see Fig 1C and 1D). The standardized mesh by Chen can nicely preserve the feature of the mesh, but the triangle quality of the mesh is general (see Fig 1E). The result of the proposed method shows that not only detailed features (e.g., wheel edges) of the mesh can be well preserved but that the triangle quality of the mesh is also better (see Fig 1F).  Fig 2B and 2D). The results of the local zoomed region are as shown in Fig 3. The triangle is uniform, and the quality of the triangle is good after that the reconstructed mesh is standardized by the proposed method and the method of Laplace (see Fig 3B and 3F). However, many detailed features of the mesh standardized by Laplace's method are lost (see Fig 3B), and the triangle quality of the mesh standardized by other methods are worse (see Fig 3C, 3D and 3E). Therefore, the proposed method can well standardize the mesh reconstructed by ICT serial images compared with existing methods.
To precisely compare the performances of the proposed and existing methods, the information of the mesh standardized by different methods is counted. Table 2 shows the related information of the original mesh. The bounding box size, surface area, volume, maximum error, average error of the meshes standardized by different methods and the consuming time of the  Table 3. Table 4 shows the triangle quality distribution of the different mesh standardization methods.
According to the data in Tables 2 and 3, for the wheel model, the error of bounding box size, surface area error and volume error of mesh standardized by Laplace are maximal. Error of bounding box size, surface area error and volume error of mesh standardized by Chen are much closer to original mesh. And our method nearly keeps the bounding box size, surface area and volume consistent with the original mesh's. The maximum error and mean error of mesh standardized by Laplace are maximal. For the maximum error of the standardized mesh, Chen's method is slightly less than the proposed method. However, the mean error of mesh standardized by the proposed method is minimal. For the cylinder head model, the errors of mesh standardized by Laplace are maximal, which contains the error of bounding box size, surface area error, volume error, maximum error and mean error. The errors of mesh standardized by Chen and Taubin are smaller than other methods except our method. According to the data in Table 4, compared with existing methods, the proposed method can obtain the largest number of triangles with good quality (0.6 < Q 1.0) and the least number of triangles with poor quality (Q 0.3). Where, the triangles with good quality on the wheel model is 84.41%, the triangles with good quality on the cylinder head model is 88.23%, the triangles with poor  quality on the wheel model is 2.32%, the triangles with poor quality on the cylinder head model is 1.01%.
In conclusion, our method nearly keeps the bounding box size, the volume and the surface area of the mesh consistent with the original mesh's, and better preserves the geometric features of the mesh. More importantly, the proposed method can sharply decrease the number of long-narrow triangles and greatly improve the quality of the triangles. However, we have to fit the curve surface according to the vertex and its first-order adjacent vertices in advance, and then search the optimal vertex. As a result, the consuming time of our method is relatively long.

Conclusions
In this paper, we present a novel method for triangle mesh standardization based on the PSO. The proposed method can effectively improve the triangle quality of the mesh while preserving the geometric features and details of the original mesh. The PSO is modified by introducing the center of the particles, the constraint factor and the inertia factor, which can effectively avoid the local optimal, accelerate the convergence speed and adjust the searching region in real time. The fitted surface is used for the searching region of the PSO, which solves the volume shrinkage of the mesh standardized by most of existing algorithms. The threshold of the normal angle between the original vertex and the regulated vertex is used to determine whether the vertex needs to be regulated, to ensure that the detailed features of the triangle mesh are not lost. However, how to set the threshold of the normal angle in terms of the geometric feature of the triangle mesh and how to further shorten the consuming time of the proposed method are still problems worth considering.
Supporting Information S1 Appendix. Performance comparisons of between the proposed PSO and the other PSO. (DOCX)