- Email: [email protected]
Error evaluation of free-form surface based on distance function of measured point to surface
Computer-Aided Design 65 (2015) 11–17
Contents lists available at ScienceDirect
Computer-Aided Design journal homepage: www.elsevier.com/locate/cad
Error evaluation of free-form surface based on distance function of measured point to surface✩ Gaiyun He a , Mei Zhang b , Zhanjie Song b,c,∗ a
Key Laboratory of Mechanism Theory and Equipment Design of Ministry of Education, Tianjin University, Tianjin 300072, China
b
School of Science, Tianjin University, Tianjin 300072, China
c
Center for Applied Mathematics, Tianjin University, Tianjin 300072, China
article
info
Article history: Received 23 June 2013 Accepted 22 February 2015 Keywords: Free-form surface Distance function Profile error Error evaluation algorithm
abstract As free-form surface is widely used in engineering, it is urgently needed to develop advanced methodology of detecting and evaluating the profile error. To this end, the semantic of profile tolerance in ASME Y14.5.1M is reviewed and the mathematical definition of profile tolerance is discussed. Subsequently, a mathematical model for error evaluation is built. This mathematical model is augmented based on distance function by considering the second-order terms in the computation of the distance from point to surface. Then, a profile error evaluation algorithm, which combines Differential Evolution (DE) algorithm and Nelder–Mead (NM) algorithm, is developed to solve this model. The proposed model and optimization algorithm are validated with simulation results from a case study. Additionally, the model is superior to Least-squares (LS) model in simplicity, efficiency and robustness. © 2015 Elsevier Ltd. All rights reserved.
1. Introduction Functional or aesthetic requirements lead to the design of more complex free-form surfaces on product. Due to the complicated manufacturing process of free-form surface and the induced high cost, the false rejection rate in quality control procedures should be reduced. To this end, that the detection and evaluation methods must be improved to make the evaluated error sufficiently close to the real value. For free-form surface detection, contact measurement and noncontact measurement are two major methods to obtain the 3D coordinates of points on the surface. With the points, the measured surface can be represented so that the profile error is achieved by comparing the measured surface and the nominal one with the evaluation models and algorithms. The most commonly used approach in commercial Coordinate Measuring Machine (CMM) is the Least-Squares (LS) method. The objective function of the LS model is the sum of squared distances from the measured points to the designed surface, and the variables to be optimized are the six-degree-of-freedom motions of the measured surface (from a
✩ This paper has been recommended for acceptance by Dr. H. Suzuki.
∗ Corresponding author at: School of Science, Tianjin University, Tianjin, China. Tel.: +86 22 27402324; fax: + 86 22 27402324. E-mail address: [email protected] (Z. Song). http://dx.doi.org/10.1016/j.cad.2015.02.004 0010-4485/© 2015 Elsevier Ltd. All rights reserved.
nominal location). It is an optimization problem defined in the product space R3 × SO (3), where SO (3) and R3 denote the group of rotation transformation and the translation transformation in 3D space, respectively. To solve this optimization problem, Iterative Closest Point (ICP) algorithm is the most used method. ICP is a category of algorithms for establishing correspondences between two data sets, solving for the optimal transformation that causes alignment [1]. It is presented by Besl and McKay [2] firstly in 1992, that the ICP is sensitivity to the initial value and inefficiency. Later on, various types of algorithms have been proposed. An algebraic algorithm for workpiece location was proposed by Li and Yeung [3] in 1996 based on the singular value decomposition algorithm. Menq and Yau [4] presented an optimal match scheme applying the steepest descent algorithm that aligns the measurement data with the design data in CAD-directed dimensional inspection in 1992. Huang and Gu [5] adopted the same method to optimize the LS function in 1996. According to the steepest descent algorithm theory, the efficiency of the algorithm is low. Using freeform surface registration to inspect engineering components constructed by NURBS surfaces [6] was proposed by Ristic and Brujic [7] in 1997. In their work, the Multidimensional Downhill Simplex method is used to find the initial value of the ICP algorithm to improve its computational efficiency. Although the LS method is easy for implementation, efficient in computation and less sensitive to asperities (such as that caused by dirt or scratches on the surface), it does not conform to the ASME
12
G. He et al. / Computer-Aided Design 65 (2015) 11–17
standards [8] which recommend that the profile error should be evaluated based on the concept of minimum zone. It formulates an optimization problem called a minimax optimization problem corresponding to minimize the maximum error. In the optimization process, the maximum distance of the measured point to the nominal surface is calculated and the maximum value should be reduced after each iteration. The traditional optimization algorithm is not suitable for solving such problems. Therefore, many heuristic optimization methods are applied in solving this problem. For example, the genetic algorithm (GA) [9–11] is used in evaluating the circularity, and the particle swarm algorithm (PSA) [12,13] is employed to evaluate conicity and cylindricity. Since the minimum zone is determined by the distance from the sampling points to the measured profile, calculating the distance is crucial for the error evaluation using the above algorithm. This issue attracted attention from many scholars. The point-to-surface distance function expressed in l2 -norm is transformed into l1 -norm by some scholars [14,15]. Flory and Hofer [14] derived an approximation of the unsigned distance function from point clouds to B-spline surfaces, and employed the approximation to solve the registration of point clouds and fitting of point sets based on the l1 -norm. Recently, a signed point-to-surface distance function is defined by Zhu [16,17] who differentiated the distance function and achieved second-order Taylor series expansion. On this basis, the profile error evaluation model was presented. Though the second-order Taylor series expansion has higher precision, too much differential information is needed. Due to the difficulty of computation, the first-order Taylor expansion is used in practice. On the basis of analyzing the semantics of profile tolerance in ASME Y14.5.1 [8], and the work from Zhu [16,17], a new model of free-form surface profile error evaluation is presented in the paper. The rest of this paper is structured as follows. In Section 2, mathematical definition of profile tolerance based on the distance function is given, and the error evaluation model is built with the rigid-body motion variables added into the distance function. In this part, Zhu’s [16,17] idea on establishing the error evaluation model is taken and his method is improved. In Section 3, an advanced algorithm of free-form surface profile error is proposed. Section 4 deals with simulation and Section 5 with conclusions. 2. Free-form surface profile error evaluation model 2.1. Mathematical definition of free-form surface profile tolerance A profile tolerance zone is an area (profile of a line) or a volume (profile of a surface) generated by offsetting each point on the nominal surface in the direction normal to the nominal surface at that point. With this definition in ASME Y 14.5.1 [8], it is important to find the normal distance of a measured point to the nominal surface. Any measured point p xp , yp , zp has a unique corresponding point q xq , yq , zq in the nominal surface which is called theoretical point. The unit normal vector of q can be denoted as nq (xnq , ynq , znq ). With the ICP algorithm, p and q can compose a point pair. The signed distance from a measurement point to the nominal surface can be expressed as:
Fig. 1. Distance from point to surface.
2.2. Free-form surface profile error evaluation model The error evaluation model usually formulates an optimization problem. In this section, the objective function, the variables and the constraints of this optimization problem are proposed. Two coordinate systems CW and CM , which denote the reference frame of the nominal surface and the reference frame of the coordinate measuring device respectively, are built. The location of the CW relative to the CM is represented by the vector–matrix pair g (T , R), where T ∈ R3 is the position vector and R ∈ SO(3) is the rotation matrix, representing the position and orientation of CM relative to CW . If the measured surface has no manufacturing error, it will totally coincide with the transformed nominal surface. However, the manufacturing error usually exists. Profile error is usually used as a criterion to evaluate how much the measured points coincide with the surface. During the process of evaluation the profile error, the position and orientation of the nominal surface should be adjusted continuously to compare with the measured surface. The target position of the surface, where it minimizes the distance between the farthest measured point and the surface, can be represented by a small rotation and a small translation from the nominal one. In this paper, the rotation angle around X , Y and Z coordinate axes is expressed as α , β and γ respectively. Besides, the translation along X , Y and Z axes is expressed as δx , δy and δz respectively. Here we define two vectors
T v and ω where v = δx , δy , δz and ω= [α, β, γ ]T . The position of point q, which is the corresponding point of point p, is changed with the position of the nominal surface. It is turned into q′ , and q′ = g (T , R) = Rq + T .
d (p, S ) = (p − q) · nq
(2.1)
where S represents the surface. The signed distance from any measured point to the free-form surface and the relevant variables are shown in Fig. 1. If the measured point is located in the profile tolerance zone, the following equation should be satisfied:
(p − q) · nq 6 t
(2.2)
2 where t denotes the tolerance value usually given in a drawing.
(2.3)
Since the α , β and γ are small enough, the high order terms of Taylors expansions in R can be ignored. R* , the approximation of R, can be written as: 1
γ
−γ β
−α
∗
R =
−β α .
1
(2.4)
1
Matrix ω ˆ is introduced to simplify the calculation.
⌢
0
ω = R − I = −γ β ∗
γ 0
−β α
−α
0
(2.5)
where I is the identity matrix. For any vector x ∈ R3 , we have:
ωˆ x = ω × x.
(2.6) ′
q′
According to Eqs. (2.4) and (2.5) q and n can be expressed as q = I +ω ˆ q+v ′
′
nq = I + ω ˆ nq .
(2.7)
G. He et al. / Computer-Aided Design 65 (2015) 11–17
13
Combining Eq. (2.7) with Eq. (2.1), the following equation can be got:
∧
d (p, S ) = p − I + ω q − v ·
∧
I + ω nq .
(2.8)
With the Difference Operation, the differential of Eq. (2.8) can be given as follows:
1d (p, S ) ∧ ∧ = p − I + ω q − v · I + ω nq − (p − q) · nq ∧ ∧ = p − q − ω q − v · I + ω nq − (p − q) · nq ∧ ∧ = (p − q) · nq + (p − q) · ω nq − ω q + v · nq ∧ ∧ − ω q + v · ω nq − (p − q) · nq . ∧
(2.9)
Zhu et al. [17] deduces that theline (p − q) parallels nq , and the
Fig. 2. The signed distance function [17].
∧
ω is skew-symmetric, so (p − q)· ω nq in Eq. (2.9) is zero. Also he ∧ ∧ supposed that ω q + v · ω nq is a term higher than the second order which can be neglected. With these above hypothesis, the following equation is the approximation of Eq. (2.9):
/mm
∧ 1d (p, S ) = − ω q + v · nq = −nq · v − q × nq · ω.
Distance
1
(2.10)
Theoretical Value
0.95
DDF Approximation
0.9
FOTE Approximation
0.85 0.8 0.75 0.7
Their model for computing the distance from the point to the transformed surface was deduced as follows:
0.65
0
0.1
0.2
0.3
0.5
0.6
0.7
0.8
0.9
1
k
d (p, S ) = (p − q) · nq + 1d (p, S )
= (p − q) · nq − nq · v − q × nq · ω.
d (p, S ) = (p − q) · nq + 1d (p, S )
= (p − q) · nq − nq · v − q × nq · ω ∧ ∧ − ω q + v · ω nq .
Fig. 3. 3 Methods for calculation of distance from point to surface.
(2.11)
But it should be noted that the distance approximation will be more accurate if the second order term is accounted for. In our work, the distance function is illustrated as follows:
(2.12)
This function is called Differential Distance Function (DDF). Eq. (2.12) expressed in the form of matrix which is used as the mathematical model for computing the point to surface distance is deduced as follows:
δ
0.4
method (based on Eq. (2.13)) and First Order Taylor Expansion (FOTE) (based on Eq. (2.11) [18]) method, is shown in Fig. 3 to compare with the theoretical value. Fig. 3 shows that the distance from a point to surface described by DDF is more accurate. As presented in Section 2.2, the evaluation error will be found when the surface moves to the position where the distance is minimized between the measured point, whose distance to the surface is the farthest among all the measured points, and the surface [19] (shown in Fig. 4). So the evaluation model is built as follows: min
T
δ δ x x x δ δ y y δ y δz δz T δz d (p, S ) = α A α + b α + c β β β γ γ γ
d = max d pi , S α, β, γ , δx , δy , δz
α (2.13)
where c = xp xnq − xq xnq + yp ynq − yq ynq + zp znq − zq znq , b and A are given in Box I. Actually, the last term in Eq. (2.13) is a second-order term which is different from the work of Zhu [17]. To validate the effectiveness of the distance function, here a simulation is implemented. For a given sphere O, its radius is 1 and its center O is the origin of coordinate system. A given point Cp ∈ R3 has the coordinate (0, 0, 2). Cq is the point on the surface of sphere O closest to the point Cp . The initial distance from Cp to Cq is 1 (shown in Fig. 2 [17]). If the sphere has rigid body motion V = k α, β, γ , δx , δy , δz
T
=
k [0.07, 0.08, 0.09, 0.3, −0.2, 0.3]T , k is from zero to one, the variation of signed distance of two methods, that is the DDF
s.t.
β γ ∈ R6 . δx δy δz
(2.14)
According to the minimum zone principle, the profile error is
ε = 2d.
2.3. Convergence analysis of the model The convergence of the newly proposed model (Eq. (2.9) in Section 2.2) is given in this section that a optimum can be found. The norm of A in Eq. (2.13) can be written as a partitioned matrix as follows:
A |A| = 11 −A12
A12 A22
(2.15)
14
G. He et al. / Computer-Aided Design 65 (2015) 11–17
−x n q q −y n −znq b= yp znq − zp ynq zp xnq − xp znq xp ynq − yp xnq
0 0 0 A= 0 znq −ynq
0 0 0 −znq 0 xnq
0 0 0 ynq −xnq 0
0 −znq ynq −2 zq znq + yq ynq xq ynq + yq xnq xq znq + zq xnq
znq 0 −x n q xq ynq + yq xnq −2 xq xnq + zq znq yq znq + zq ynq
−ynq xnq 0 xq znq + zq xnq yq znq + zq ynq −2 yq ynq + xq xnq
Box I. Table 1 The parameter used in DE algorithm. Parameter
Value
NP F CR Maximum generation Termination threshold
30 0.5 0.5 300 10−6
Table 2 The parameter used in NM algorithm.
o
Fig. 4. Free-form surface profile error.
0 where A11 =
A22
0 0
0 0 0
0 0 0
, A12 =
−2 zq znq + yq ynq = xq ynq + yq xnq xq znq + zq xnq
0 −znq ynq
znq 0 −xnq
xq ynq + yq xnq −2 xq xnq + zq znq yq znq + zq ynq
−ynq xnq 0
Parameter
Value
Initial population Contraction coefficient Expansion coefficient
7 1.6 0.5
,
xq znq + zq xnq yq znq + zq ynq −2 xq xnq + yq ynq .
Then we have:
A A12 |A| = 11 −A12 A22 −A12 A22 = − A11 A12 = − |−A12 | A12 − A11 (−A12 )−1 A22 = (−1)2 |A12 |2 = |A12 |2 .
(2.16)
It can be easily proved that the sequential principal minor of A is zero. Matrix A is a positive semi-definite matrix, and the distance function is a convex function. So this objective function can reach the minimum value that complies with the definition of profile error. 3. Free-form surface profile error evaluation algorithm In the process of optimizing, the maximum distance should be found and should be reduced during each iteration. With a large number of measured points, direct search algorithm is timeconsuming and not suitable for this model. The combination of Differential Evolution algorithm (DE) [20,21] and Nelder–Mead (NM) simplex search algorithm is applied in this paper to ensure rapid and global convergence. DE algorithm is a kind of evolutionary algorithm first proposed by R. Storn and K. Price [20], which has the advantage of rapid
optimization. DE has three control parameters: the initial population size-NP, crossover constant-CR and the amplification factor of the difference vector-F. Also termination threshold and maximum generation should be set according to the optimum precision needed. The NM simplex algorithm, first published in 1965 [22], is an enormously popular direct search method for multidimensional unconstrained minimization. The NM algorithm has the control parameters including initial population, contraction coefficient and expansion coefficient. More details about these two algorithm can be found in [20,22]. As proposed in the two papers, simplex algorithm like NM algorithm usually depends on a good initial points, otherwise it takes long time to converge and easily fall into a local minimum. The minimum can be found globally with the DE algorithm. However, it is time-consuming in a small scale because it keeps a lot population in every generation. So the two algorithm are combined to solve the error evaluation model in this paper. The overall framework of the error evaluation algorithm is the DE algorithm. The NM algorithm is inserted in each generation of DE algorithm through the way that the best individual generated by DE algorithm is selected to be the initial population in NM algorithm and the result outputted by NM is used for the next generation of DE algorithm. Through analyzing the theory of DE and NM algorithm and considering the model in this paper, the parameters of the DE and the NM algorithm used are illustrated in Tables 1 and 2, respectively. The flowchart of the error evaluation algorithm is shown in Fig. 5, where ν = [α, β, γ , δx , δy , δz ]T . 4. Simulation of error evaluation algorithm In order to verify the error evaluation model and algorithm, the simulation is designed in this paper. Firstly, a mold surface, a blade surface, a cycloidal gear surface, a Gleason gear surface, an involute gear surface and an aircraft inlet surface are designed
G. He et al. / Computer-Aided Design 65 (2015) 11–17
15
Fig. 5. Flowchart of the evaluation algorithm.
a
b
c
d
e
f
Fig. 6. Simulation model diagram.
by 3D CAD software Siemens NX 8.0. Their maximum projection areas are 90 × 120, 33 × 80, 20 × 40, 5 × 26, 10 × 3 and 40 × 30 respectively (the maximum projection area is used to characterize the size of surface). They are marked from a to f and Fig. 6 is the diagram of the simulation model. Then, system error Err = (0.1, 0.198, 0.2, 1, 2, 3) and random error ε ∈ N (0, 0.04) which means a normal distribution with mean 0 and variance 0.04, are added into the measured information. Finally, enough points are sampled on the surfaces to reflect the quality information. Coordinates of these points and their corresponding unit normal vectors can be given by the UG/OPEN Grip (a API Tool for the
secondary development that provides some sentences to calculate the geometrical information in batch). Profile errors of surfaces in Fig. 6 are evaluated using the DDF model and algorithm we proposed. In contrast experiment these surfaces are evaluated by LS model and the same algorithm we proposed. The flowchart of the simulation process is shown in Fig. 7. The results are shown in Table 3 and Fig. 8. The advanced method raises the average evaluation precision by 35.6%. Especially, when the measured surface has large curvature or large change of curvature, such as surface (a), (b) and (f), the new method has better performance that the evaluation precision goes
16
G. He et al. / Computer-Aided Design 65 (2015) 11–17 Table 3 Comparison of the profile error evaluated by two methods. Profile error/mm
New method LS method
Surf.(a)
Surf.(b)
Surf.(c)
Surf.(d)
Surf.(e)
Surf.(f)
1.2065 2.3186
1.0073 2.1854
1.0579 2.3836
1.3280 1.8236
1.3462 1.8085
1.1203 2.3980
Surf. = Surface.
Table 4 Comparison of variables after optimization by two methods.
ν α/rad New method LS method
β/rad
−0.0075
0.0047
0.0077
−0.0102
γ /rad
δx /mm
δy /mm
δz /mm
0.0107 0.0088
−0.2341
−0.1306
0.4952
0.0641
−0.0709 −0.2131
Fig. 8. Comparison of the evaluated error of two methods used in the 6 surfaces.
Fig. 7. Flowchart of the simulation process.
up nearly 54%. Take surface (b) as an example, the initial error is 2.0682 mm. The error after optimization by the method is 1.0073 mm. The initial position and the position after optimization of points and measuring points are shown in Figs. 9 and 10.
The parameter values after optimization are shown in Table 4. The position transformation matrixes calculated by two methods are shown in Table 5 where we can find that the matrix obtained by new method is closer to the theoretical value. In order to verify the stability of the algorithm, the simulation experiments are implemented 50 times. Fig. 11 shows that the error evaluated by new method is stable and more close to the true tolerance. The new method is also time efficient. 50 times simulation experiments are implemented on a computer with i7 CPU and 32G memory, Windows 7 operation system and Matlab 7.12 platform. Applied on the surface (b), when 80 points are sampled, it costs averagely 25.4 s. According to the analysis, most of the time is consumed in the calculation of the distance from points to surface. We adopt a strategy to avoid this repetitive calculation. As shown in Eq. (2.13), matrix A, vector b, and c are constant during the optimization. If these constants are calculated and stored in the memory, they will be called in the program without computing. This strategy is time saving but memory consuming.
Fig. 9. Initial position of theoretical points and measured points.
G. He et al. / Computer-Aided Design 65 (2015) 11–17
17
Fig. 10. Theoretical points and measured points position after optimization.
Deviation /mm
3
Acknowledgments
the new method LS method
2.5
The research is supported by NSFC under Grant Nos. 51275348 and 51105271. The authors would like to express their sincere gratitude to the anonymous referees, the editor, Prof. Ying Li and Dr. Longzhen Guo for many valuable suggestions and comments that helped to improve the article.
2 1.5 1 0.5
References 0
5
10
15
20
25
30
35
40
45
50
Experiment Times Fig. 11. Analysis of algorithm stability. Table 5 Comparison of position transformation matrixes. Optimization method
Position transformation matrixes 0.9605 −0.1782 0.2136
Theoretical value
0 0.9628 −0.1620 0.2164 0 0.9593 −0.1626 0.2310 0
New method
LS method
0.1947 0.9791 −0.0586 0 0.1784 0.9822 −0.0584 0 0.1771 0.9832 −0.0432 0
−0.1987 0.0978 0.9752 0
−0.2030 0.0948 0.9745 0
−0.2201 0.0823 0.9720 0
1 2 3 1 1.4056 2.2403 2.8781 1 0.6763 2.0456 2.5940 1
5. Conclusions In this paper, firstly, profile tolerance was defined mathematically and the distance function presented by Zhu et al. [16,17] was improved. On this basis, a mathematical model for profile error evaluation was proposed. This model was developed based on the distance function by considering the second-order terms in computation. Less differential information of surfaces was needed in this model such as theoretical points, measured points and unit normal vectors of the surface in theoretical points. Compared with the distance function in the model proposed by Zhu et al. [16,17], the distance function in this new model has been proved to be convex, which provides a theoretical proof for the convergence of this model. The developed model was optimized with the algorithms combined the advantages of DE and NM algorithms. This optimization algorithm was validated with some simulations. At last, case study results showed that the average evaluation precision calculated with the advanced model was improved by 35.6% compared with LS method.
[1] Bouaziz S, Tagliasacchi A, Pauly M. Sparse iterative closest point. In: Computer graphics forum, Vol. 32. Wiley Online Library; 2013. p. 113–23. [2] Besl PJ, McKay ND. Method for registration of 3-d shapes. In: Robotics-DL tentative. International Society for Optics and Photonics; 1992. p. 586–606. [3] Li X, Yeung M, Li Z. An algebraic algorithm for workpiece localization. In: Robotics and automation, 1996. Proceedings., 1996 IEEE international conference on, Vol. 1. IEEE; 1996. p. 152–8. [4] Menq CH, Yau HT, Lai GY. Automated precision measurement of surface profile in cad-directed inspection. IEEE Trans Robot Autom 1992;8(2):268–78. [5] Huang X, Gu P, Zernicke R. Localization and comparison of two free-form surfaces. Comput Aided Des 1996;28(12):1017–22. [6] Fuhr RD, Hsieh L, Kallay M. Object-oriented paradigm for nurbs curve and surface design. Comput-Aided Des 1995;27(2):95–100. [7] Ristic M, Brujic D. Efficient registration of nurbs geometry. Image Vis Comput 1997;15(12):925–35. [8] Principles T. ASME/ANSI standard y14. 5.1 m-94. New York (NY) 1994. [9] Shopova EG, Vaklieva-Bancheva NG. Basic—a genetic algorithm for engineering problems solution. Comput Chem Eng 2006;30(8):1293–309. [10] Wen X, Song A. An improved genetic algorithm for planar and spatial straightness error evaluation. Int J Mach Tools Manuf 2003;43(11):1157–62. [11] Wen X, Xia Q, Zhao Y. An effective genetic algorithm for circularity error unified evaluation. Int J Mach Tools Manuf 2006;46(14):1770–7. [12] del Valle Y, Venayagamoorthy GK, Mohagheghi S, Hernandez JC, Harley RG. Particle swarm optimization: basic concepts, variants and applications in power systems. IEEE Trans Evol Comput 2008;12(2):171–95. [13] Wen X, Zhao Y, Xu Y, Sheng D. Quasiparticle swarm optimization for crosssection linear profile error evaluation of variation elliptical piston skirt. Math Probl Eng 2012;2012. [14] Flöry S, Hofer M. Surface fitting and registration of point clouds using approximations of the unsigned distance function. Comput Aided Geom Design 2010;27(1):60–77. [15] Pottmann H, Leopoldseder S, Hofer M. Registration without ICP. Comput Vis Image Underst 2004;95(1):54–71. [16] Zhu L, Xiong Z, Ding H, Xiong Y. A distance function based approach for localization and profile error evaluation of complex surface. J Manuf Sci Eng 2004;126(3):542–54. [17] Zhu LM, Zhang XM, Ding H, Xiong YL. Geometry of signed point-to-surface distance function and its application to surface approximation. J Comput Inf Sci Eng 2010;10:041003. [18] Ding H, Zhu L. Geometric theories and methods for digital manufacturing of complex surfaces. Beijing, China: Science Press; 2011. (in Chinese). [19] Gou J, Chu Y, Li Z. A geometric theory of form, profile, and orientation tolerances. Precis Eng 1999;23(2):79–93. [20] Storn R, Price K. Differential evolution—a simple and efficient heuristic for global optimization over continuous spaces. J Global Optim 1997;11(4): 341–59. [21] Okdem S. A simple and global optimization algorithm for engineering problems: differential evolution algorithm. Turk J Electr Eng 2004;12(1). [22] Nelder JA, Mead R. A simplex method for function minimization. Comput J 1965;7(4):308–13.
Contents lists available at ScienceDirect
Computer-Aided Design journal homepage: www.elsevier.com/locate/cad
Error evaluation of free-form surface based on distance function of measured point to surface✩ Gaiyun He a , Mei Zhang b , Zhanjie Song b,c,∗ a
Key Laboratory of Mechanism Theory and Equipment Design of Ministry of Education, Tianjin University, Tianjin 300072, China
b
School of Science, Tianjin University, Tianjin 300072, China
c
Center for Applied Mathematics, Tianjin University, Tianjin 300072, China
article
info
Article history: Received 23 June 2013 Accepted 22 February 2015 Keywords: Free-form surface Distance function Profile error Error evaluation algorithm
abstract As free-form surface is widely used in engineering, it is urgently needed to develop advanced methodology of detecting and evaluating the profile error. To this end, the semantic of profile tolerance in ASME Y14.5.1M is reviewed and the mathematical definition of profile tolerance is discussed. Subsequently, a mathematical model for error evaluation is built. This mathematical model is augmented based on distance function by considering the second-order terms in the computation of the distance from point to surface. Then, a profile error evaluation algorithm, which combines Differential Evolution (DE) algorithm and Nelder–Mead (NM) algorithm, is developed to solve this model. The proposed model and optimization algorithm are validated with simulation results from a case study. Additionally, the model is superior to Least-squares (LS) model in simplicity, efficiency and robustness. © 2015 Elsevier Ltd. All rights reserved.
1. Introduction Functional or aesthetic requirements lead to the design of more complex free-form surfaces on product. Due to the complicated manufacturing process of free-form surface and the induced high cost, the false rejection rate in quality control procedures should be reduced. To this end, that the detection and evaluation methods must be improved to make the evaluated error sufficiently close to the real value. For free-form surface detection, contact measurement and noncontact measurement are two major methods to obtain the 3D coordinates of points on the surface. With the points, the measured surface can be represented so that the profile error is achieved by comparing the measured surface and the nominal one with the evaluation models and algorithms. The most commonly used approach in commercial Coordinate Measuring Machine (CMM) is the Least-Squares (LS) method. The objective function of the LS model is the sum of squared distances from the measured points to the designed surface, and the variables to be optimized are the six-degree-of-freedom motions of the measured surface (from a
✩ This paper has been recommended for acceptance by Dr. H. Suzuki.
∗ Corresponding author at: School of Science, Tianjin University, Tianjin, China. Tel.: +86 22 27402324; fax: + 86 22 27402324. E-mail address: [email protected] (Z. Song). http://dx.doi.org/10.1016/j.cad.2015.02.004 0010-4485/© 2015 Elsevier Ltd. All rights reserved.
nominal location). It is an optimization problem defined in the product space R3 × SO (3), where SO (3) and R3 denote the group of rotation transformation and the translation transformation in 3D space, respectively. To solve this optimization problem, Iterative Closest Point (ICP) algorithm is the most used method. ICP is a category of algorithms for establishing correspondences between two data sets, solving for the optimal transformation that causes alignment [1]. It is presented by Besl and McKay [2] firstly in 1992, that the ICP is sensitivity to the initial value and inefficiency. Later on, various types of algorithms have been proposed. An algebraic algorithm for workpiece location was proposed by Li and Yeung [3] in 1996 based on the singular value decomposition algorithm. Menq and Yau [4] presented an optimal match scheme applying the steepest descent algorithm that aligns the measurement data with the design data in CAD-directed dimensional inspection in 1992. Huang and Gu [5] adopted the same method to optimize the LS function in 1996. According to the steepest descent algorithm theory, the efficiency of the algorithm is low. Using freeform surface registration to inspect engineering components constructed by NURBS surfaces [6] was proposed by Ristic and Brujic [7] in 1997. In their work, the Multidimensional Downhill Simplex method is used to find the initial value of the ICP algorithm to improve its computational efficiency. Although the LS method is easy for implementation, efficient in computation and less sensitive to asperities (such as that caused by dirt or scratches on the surface), it does not conform to the ASME
12
G. He et al. / Computer-Aided Design 65 (2015) 11–17
standards [8] which recommend that the profile error should be evaluated based on the concept of minimum zone. It formulates an optimization problem called a minimax optimization problem corresponding to minimize the maximum error. In the optimization process, the maximum distance of the measured point to the nominal surface is calculated and the maximum value should be reduced after each iteration. The traditional optimization algorithm is not suitable for solving such problems. Therefore, many heuristic optimization methods are applied in solving this problem. For example, the genetic algorithm (GA) [9–11] is used in evaluating the circularity, and the particle swarm algorithm (PSA) [12,13] is employed to evaluate conicity and cylindricity. Since the minimum zone is determined by the distance from the sampling points to the measured profile, calculating the distance is crucial for the error evaluation using the above algorithm. This issue attracted attention from many scholars. The point-to-surface distance function expressed in l2 -norm is transformed into l1 -norm by some scholars [14,15]. Flory and Hofer [14] derived an approximation of the unsigned distance function from point clouds to B-spline surfaces, and employed the approximation to solve the registration of point clouds and fitting of point sets based on the l1 -norm. Recently, a signed point-to-surface distance function is defined by Zhu [16,17] who differentiated the distance function and achieved second-order Taylor series expansion. On this basis, the profile error evaluation model was presented. Though the second-order Taylor series expansion has higher precision, too much differential information is needed. Due to the difficulty of computation, the first-order Taylor expansion is used in practice. On the basis of analyzing the semantics of profile tolerance in ASME Y14.5.1 [8], and the work from Zhu [16,17], a new model of free-form surface profile error evaluation is presented in the paper. The rest of this paper is structured as follows. In Section 2, mathematical definition of profile tolerance based on the distance function is given, and the error evaluation model is built with the rigid-body motion variables added into the distance function. In this part, Zhu’s [16,17] idea on establishing the error evaluation model is taken and his method is improved. In Section 3, an advanced algorithm of free-form surface profile error is proposed. Section 4 deals with simulation and Section 5 with conclusions. 2. Free-form surface profile error evaluation model 2.1. Mathematical definition of free-form surface profile tolerance A profile tolerance zone is an area (profile of a line) or a volume (profile of a surface) generated by offsetting each point on the nominal surface in the direction normal to the nominal surface at that point. With this definition in ASME Y 14.5.1 [8], it is important to find the normal distance of a measured point to the nominal surface. Any measured point p xp , yp , zp has a unique corresponding point q xq , yq , zq in the nominal surface which is called theoretical point. The unit normal vector of q can be denoted as nq (xnq , ynq , znq ). With the ICP algorithm, p and q can compose a point pair. The signed distance from a measurement point to the nominal surface can be expressed as:
Fig. 1. Distance from point to surface.
2.2. Free-form surface profile error evaluation model The error evaluation model usually formulates an optimization problem. In this section, the objective function, the variables and the constraints of this optimization problem are proposed. Two coordinate systems CW and CM , which denote the reference frame of the nominal surface and the reference frame of the coordinate measuring device respectively, are built. The location of the CW relative to the CM is represented by the vector–matrix pair g (T , R), where T ∈ R3 is the position vector and R ∈ SO(3) is the rotation matrix, representing the position and orientation of CM relative to CW . If the measured surface has no manufacturing error, it will totally coincide with the transformed nominal surface. However, the manufacturing error usually exists. Profile error is usually used as a criterion to evaluate how much the measured points coincide with the surface. During the process of evaluation the profile error, the position and orientation of the nominal surface should be adjusted continuously to compare with the measured surface. The target position of the surface, where it minimizes the distance between the farthest measured point and the surface, can be represented by a small rotation and a small translation from the nominal one. In this paper, the rotation angle around X , Y and Z coordinate axes is expressed as α , β and γ respectively. Besides, the translation along X , Y and Z axes is expressed as δx , δy and δz respectively. Here we define two vectors
T v and ω where v = δx , δy , δz and ω= [α, β, γ ]T . The position of point q, which is the corresponding point of point p, is changed with the position of the nominal surface. It is turned into q′ , and q′ = g (T , R) = Rq + T .
d (p, S ) = (p − q) · nq
(2.1)
where S represents the surface. The signed distance from any measured point to the free-form surface and the relevant variables are shown in Fig. 1. If the measured point is located in the profile tolerance zone, the following equation should be satisfied:
(p − q) · nq 6 t
(2.2)
2 where t denotes the tolerance value usually given in a drawing.
(2.3)
Since the α , β and γ are small enough, the high order terms of Taylors expansions in R can be ignored. R* , the approximation of R, can be written as: 1
γ
−γ β
−α
∗
R =
−β α .
1
(2.4)
1
Matrix ω ˆ is introduced to simplify the calculation.
⌢
0
ω = R − I = −γ β ∗
γ 0
−β α
−α
0
(2.5)
where I is the identity matrix. For any vector x ∈ R3 , we have:
ωˆ x = ω × x.
(2.6) ′
q′
According to Eqs. (2.4) and (2.5) q and n can be expressed as q = I +ω ˆ q+v ′
′
nq = I + ω ˆ nq .
(2.7)
G. He et al. / Computer-Aided Design 65 (2015) 11–17
13
Combining Eq. (2.7) with Eq. (2.1), the following equation can be got:
∧
d (p, S ) = p − I + ω q − v ·
∧
I + ω nq .
(2.8)
With the Difference Operation, the differential of Eq. (2.8) can be given as follows:
1d (p, S ) ∧ ∧ = p − I + ω q − v · I + ω nq − (p − q) · nq ∧ ∧ = p − q − ω q − v · I + ω nq − (p − q) · nq ∧ ∧ = (p − q) · nq + (p − q) · ω nq − ω q + v · nq ∧ ∧ − ω q + v · ω nq − (p − q) · nq . ∧
(2.9)
Zhu et al. [17] deduces that theline (p − q) parallels nq , and the
Fig. 2. The signed distance function [17].
∧
ω is skew-symmetric, so (p − q)· ω nq in Eq. (2.9) is zero. Also he ∧ ∧ supposed that ω q + v · ω nq is a term higher than the second order which can be neglected. With these above hypothesis, the following equation is the approximation of Eq. (2.9):
/mm
∧ 1d (p, S ) = − ω q + v · nq = −nq · v − q × nq · ω.
Distance
1
(2.10)
Theoretical Value
0.95
DDF Approximation
0.9
FOTE Approximation
0.85 0.8 0.75 0.7
Their model for computing the distance from the point to the transformed surface was deduced as follows:
0.65
0
0.1
0.2
0.3
0.5
0.6
0.7
0.8
0.9
1
k
d (p, S ) = (p − q) · nq + 1d (p, S )
= (p − q) · nq − nq · v − q × nq · ω.
d (p, S ) = (p − q) · nq + 1d (p, S )
= (p − q) · nq − nq · v − q × nq · ω ∧ ∧ − ω q + v · ω nq .
Fig. 3. 3 Methods for calculation of distance from point to surface.
(2.11)
But it should be noted that the distance approximation will be more accurate if the second order term is accounted for. In our work, the distance function is illustrated as follows:
(2.12)
This function is called Differential Distance Function (DDF). Eq. (2.12) expressed in the form of matrix which is used as the mathematical model for computing the point to surface distance is deduced as follows:
δ
0.4
method (based on Eq. (2.13)) and First Order Taylor Expansion (FOTE) (based on Eq. (2.11) [18]) method, is shown in Fig. 3 to compare with the theoretical value. Fig. 3 shows that the distance from a point to surface described by DDF is more accurate. As presented in Section 2.2, the evaluation error will be found when the surface moves to the position where the distance is minimized between the measured point, whose distance to the surface is the farthest among all the measured points, and the surface [19] (shown in Fig. 4). So the evaluation model is built as follows: min
T
δ δ x x x δ δ y y δ y δz δz T δz d (p, S ) = α A α + b α + c β β β γ γ γ
d = max d pi , S α, β, γ , δx , δy , δz
α (2.13)
where c = xp xnq − xq xnq + yp ynq − yq ynq + zp znq − zq znq , b and A are given in Box I. Actually, the last term in Eq. (2.13) is a second-order term which is different from the work of Zhu [17]. To validate the effectiveness of the distance function, here a simulation is implemented. For a given sphere O, its radius is 1 and its center O is the origin of coordinate system. A given point Cp ∈ R3 has the coordinate (0, 0, 2). Cq is the point on the surface of sphere O closest to the point Cp . The initial distance from Cp to Cq is 1 (shown in Fig. 2 [17]). If the sphere has rigid body motion V = k α, β, γ , δx , δy , δz
T
=
k [0.07, 0.08, 0.09, 0.3, −0.2, 0.3]T , k is from zero to one, the variation of signed distance of two methods, that is the DDF
s.t.
β γ ∈ R6 . δx δy δz
(2.14)
According to the minimum zone principle, the profile error is
ε = 2d.
2.3. Convergence analysis of the model The convergence of the newly proposed model (Eq. (2.9) in Section 2.2) is given in this section that a optimum can be found. The norm of A in Eq. (2.13) can be written as a partitioned matrix as follows:
A |A| = 11 −A12
A12 A22
(2.15)
14
G. He et al. / Computer-Aided Design 65 (2015) 11–17
−x n q q −y n −znq b= yp znq − zp ynq zp xnq − xp znq xp ynq − yp xnq
0 0 0 A= 0 znq −ynq
0 0 0 −znq 0 xnq
0 0 0 ynq −xnq 0
0 −znq ynq −2 zq znq + yq ynq xq ynq + yq xnq xq znq + zq xnq
znq 0 −x n q xq ynq + yq xnq −2 xq xnq + zq znq yq znq + zq ynq
−ynq xnq 0 xq znq + zq xnq yq znq + zq ynq −2 yq ynq + xq xnq
Box I. Table 1 The parameter used in DE algorithm. Parameter
Value
NP F CR Maximum generation Termination threshold
30 0.5 0.5 300 10−6
Table 2 The parameter used in NM algorithm.
o
Fig. 4. Free-form surface profile error.
0 where A11 =
A22
0 0
0 0 0
0 0 0
, A12 =
−2 zq znq + yq ynq = xq ynq + yq xnq xq znq + zq xnq
0 −znq ynq
znq 0 −xnq
xq ynq + yq xnq −2 xq xnq + zq znq yq znq + zq ynq
−ynq xnq 0
Parameter
Value
Initial population Contraction coefficient Expansion coefficient
7 1.6 0.5
,
xq znq + zq xnq yq znq + zq ynq −2 xq xnq + yq ynq .
Then we have:
A A12 |A| = 11 −A12 A22 −A12 A22 = − A11 A12 = − |−A12 | A12 − A11 (−A12 )−1 A22 = (−1)2 |A12 |2 = |A12 |2 .
(2.16)
It can be easily proved that the sequential principal minor of A is zero. Matrix A is a positive semi-definite matrix, and the distance function is a convex function. So this objective function can reach the minimum value that complies with the definition of profile error. 3. Free-form surface profile error evaluation algorithm In the process of optimizing, the maximum distance should be found and should be reduced during each iteration. With a large number of measured points, direct search algorithm is timeconsuming and not suitable for this model. The combination of Differential Evolution algorithm (DE) [20,21] and Nelder–Mead (NM) simplex search algorithm is applied in this paper to ensure rapid and global convergence. DE algorithm is a kind of evolutionary algorithm first proposed by R. Storn and K. Price [20], which has the advantage of rapid
optimization. DE has three control parameters: the initial population size-NP, crossover constant-CR and the amplification factor of the difference vector-F. Also termination threshold and maximum generation should be set according to the optimum precision needed. The NM simplex algorithm, first published in 1965 [22], is an enormously popular direct search method for multidimensional unconstrained minimization. The NM algorithm has the control parameters including initial population, contraction coefficient and expansion coefficient. More details about these two algorithm can be found in [20,22]. As proposed in the two papers, simplex algorithm like NM algorithm usually depends on a good initial points, otherwise it takes long time to converge and easily fall into a local minimum. The minimum can be found globally with the DE algorithm. However, it is time-consuming in a small scale because it keeps a lot population in every generation. So the two algorithm are combined to solve the error evaluation model in this paper. The overall framework of the error evaluation algorithm is the DE algorithm. The NM algorithm is inserted in each generation of DE algorithm through the way that the best individual generated by DE algorithm is selected to be the initial population in NM algorithm and the result outputted by NM is used for the next generation of DE algorithm. Through analyzing the theory of DE and NM algorithm and considering the model in this paper, the parameters of the DE and the NM algorithm used are illustrated in Tables 1 and 2, respectively. The flowchart of the error evaluation algorithm is shown in Fig. 5, where ν = [α, β, γ , δx , δy , δz ]T . 4. Simulation of error evaluation algorithm In order to verify the error evaluation model and algorithm, the simulation is designed in this paper. Firstly, a mold surface, a blade surface, a cycloidal gear surface, a Gleason gear surface, an involute gear surface and an aircraft inlet surface are designed
G. He et al. / Computer-Aided Design 65 (2015) 11–17
15
Fig. 5. Flowchart of the evaluation algorithm.
a
b
c
d
e
f
Fig. 6. Simulation model diagram.
by 3D CAD software Siemens NX 8.0. Their maximum projection areas are 90 × 120, 33 × 80, 20 × 40, 5 × 26, 10 × 3 and 40 × 30 respectively (the maximum projection area is used to characterize the size of surface). They are marked from a to f and Fig. 6 is the diagram of the simulation model. Then, system error Err = (0.1, 0.198, 0.2, 1, 2, 3) and random error ε ∈ N (0, 0.04) which means a normal distribution with mean 0 and variance 0.04, are added into the measured information. Finally, enough points are sampled on the surfaces to reflect the quality information. Coordinates of these points and their corresponding unit normal vectors can be given by the UG/OPEN Grip (a API Tool for the
secondary development that provides some sentences to calculate the geometrical information in batch). Profile errors of surfaces in Fig. 6 are evaluated using the DDF model and algorithm we proposed. In contrast experiment these surfaces are evaluated by LS model and the same algorithm we proposed. The flowchart of the simulation process is shown in Fig. 7. The results are shown in Table 3 and Fig. 8. The advanced method raises the average evaluation precision by 35.6%. Especially, when the measured surface has large curvature or large change of curvature, such as surface (a), (b) and (f), the new method has better performance that the evaluation precision goes
16
G. He et al. / Computer-Aided Design 65 (2015) 11–17 Table 3 Comparison of the profile error evaluated by two methods. Profile error/mm
New method LS method
Surf.(a)
Surf.(b)
Surf.(c)
Surf.(d)
Surf.(e)
Surf.(f)
1.2065 2.3186
1.0073 2.1854
1.0579 2.3836
1.3280 1.8236
1.3462 1.8085
1.1203 2.3980
Surf. = Surface.
Table 4 Comparison of variables after optimization by two methods.
ν α/rad New method LS method
β/rad
−0.0075
0.0047
0.0077
−0.0102
γ /rad
δx /mm
δy /mm
δz /mm
0.0107 0.0088
−0.2341
−0.1306
0.4952
0.0641
−0.0709 −0.2131
Fig. 8. Comparison of the evaluated error of two methods used in the 6 surfaces.
Fig. 7. Flowchart of the simulation process.
up nearly 54%. Take surface (b) as an example, the initial error is 2.0682 mm. The error after optimization by the method is 1.0073 mm. The initial position and the position after optimization of points and measuring points are shown in Figs. 9 and 10.
The parameter values after optimization are shown in Table 4. The position transformation matrixes calculated by two methods are shown in Table 5 where we can find that the matrix obtained by new method is closer to the theoretical value. In order to verify the stability of the algorithm, the simulation experiments are implemented 50 times. Fig. 11 shows that the error evaluated by new method is stable and more close to the true tolerance. The new method is also time efficient. 50 times simulation experiments are implemented on a computer with i7 CPU and 32G memory, Windows 7 operation system and Matlab 7.12 platform. Applied on the surface (b), when 80 points are sampled, it costs averagely 25.4 s. According to the analysis, most of the time is consumed in the calculation of the distance from points to surface. We adopt a strategy to avoid this repetitive calculation. As shown in Eq. (2.13), matrix A, vector b, and c are constant during the optimization. If these constants are calculated and stored in the memory, they will be called in the program without computing. This strategy is time saving but memory consuming.
Fig. 9. Initial position of theoretical points and measured points.
G. He et al. / Computer-Aided Design 65 (2015) 11–17
17
Fig. 10. Theoretical points and measured points position after optimization.
Deviation /mm
3
Acknowledgments
the new method LS method
2.5
The research is supported by NSFC under Grant Nos. 51275348 and 51105271. The authors would like to express their sincere gratitude to the anonymous referees, the editor, Prof. Ying Li and Dr. Longzhen Guo for many valuable suggestions and comments that helped to improve the article.
2 1.5 1 0.5
References 0
5
10
15
20
25
30
35
40
45
50
Experiment Times Fig. 11. Analysis of algorithm stability. Table 5 Comparison of position transformation matrixes. Optimization method
Position transformation matrixes 0.9605 −0.1782 0.2136
Theoretical value
0 0.9628 −0.1620 0.2164 0 0.9593 −0.1626 0.2310 0
New method
LS method
0.1947 0.9791 −0.0586 0 0.1784 0.9822 −0.0584 0 0.1771 0.9832 −0.0432 0
−0.1987 0.0978 0.9752 0
−0.2030 0.0948 0.9745 0
−0.2201 0.0823 0.9720 0
1 2 3 1 1.4056 2.2403 2.8781 1 0.6763 2.0456 2.5940 1
5. Conclusions In this paper, firstly, profile tolerance was defined mathematically and the distance function presented by Zhu et al. [16,17] was improved. On this basis, a mathematical model for profile error evaluation was proposed. This model was developed based on the distance function by considering the second-order terms in computation. Less differential information of surfaces was needed in this model such as theoretical points, measured points and unit normal vectors of the surface in theoretical points. Compared with the distance function in the model proposed by Zhu et al. [16,17], the distance function in this new model has been proved to be convex, which provides a theoretical proof for the convergence of this model. The developed model was optimized with the algorithms combined the advantages of DE and NM algorithms. This optimization algorithm was validated with some simulations. At last, case study results showed that the average evaluation precision calculated with the advanced model was improved by 35.6% compared with LS method.
[1] Bouaziz S, Tagliasacchi A, Pauly M. Sparse iterative closest point. In: Computer graphics forum, Vol. 32. Wiley Online Library; 2013. p. 113–23. [2] Besl PJ, McKay ND. Method for registration of 3-d shapes. In: Robotics-DL tentative. International Society for Optics and Photonics; 1992. p. 586–606. [3] Li X, Yeung M, Li Z. An algebraic algorithm for workpiece localization. In: Robotics and automation, 1996. Proceedings., 1996 IEEE international conference on, Vol. 1. IEEE; 1996. p. 152–8. [4] Menq CH, Yau HT, Lai GY. Automated precision measurement of surface profile in cad-directed inspection. IEEE Trans Robot Autom 1992;8(2):268–78. [5] Huang X, Gu P, Zernicke R. Localization and comparison of two free-form surfaces. Comput Aided Des 1996;28(12):1017–22. [6] Fuhr RD, Hsieh L, Kallay M. Object-oriented paradigm for nurbs curve and surface design. Comput-Aided Des 1995;27(2):95–100. [7] Ristic M, Brujic D. Efficient registration of nurbs geometry. Image Vis Comput 1997;15(12):925–35. [8] Principles T. ASME/ANSI standard y14. 5.1 m-94. New York (NY) 1994. [9] Shopova EG, Vaklieva-Bancheva NG. Basic—a genetic algorithm for engineering problems solution. Comput Chem Eng 2006;30(8):1293–309. [10] Wen X, Song A. An improved genetic algorithm for planar and spatial straightness error evaluation. Int J Mach Tools Manuf 2003;43(11):1157–62. [11] Wen X, Xia Q, Zhao Y. An effective genetic algorithm for circularity error unified evaluation. Int J Mach Tools Manuf 2006;46(14):1770–7. [12] del Valle Y, Venayagamoorthy GK, Mohagheghi S, Hernandez JC, Harley RG. Particle swarm optimization: basic concepts, variants and applications in power systems. IEEE Trans Evol Comput 2008;12(2):171–95. [13] Wen X, Zhao Y, Xu Y, Sheng D. Quasiparticle swarm optimization for crosssection linear profile error evaluation of variation elliptical piston skirt. Math Probl Eng 2012;2012. [14] Flöry S, Hofer M. Surface fitting and registration of point clouds using approximations of the unsigned distance function. Comput Aided Geom Design 2010;27(1):60–77. [15] Pottmann H, Leopoldseder S, Hofer M. Registration without ICP. Comput Vis Image Underst 2004;95(1):54–71. [16] Zhu L, Xiong Z, Ding H, Xiong Y. A distance function based approach for localization and profile error evaluation of complex surface. J Manuf Sci Eng 2004;126(3):542–54. [17] Zhu LM, Zhang XM, Ding H, Xiong YL. Geometry of signed point-to-surface distance function and its application to surface approximation. J Comput Inf Sci Eng 2010;10:041003. [18] Ding H, Zhu L. Geometric theories and methods for digital manufacturing of complex surfaces. Beijing, China: Science Press; 2011. (in Chinese). [19] Gou J, Chu Y, Li Z. A geometric theory of form, profile, and orientation tolerances. Precis Eng 1999;23(2):79–93. [20] Storn R, Price K. Differential evolution—a simple and efficient heuristic for global optimization over continuous spaces. J Global Optim 1997;11(4): 341–59. [21] Okdem S. A simple and global optimization algorithm for engineering problems: differential evolution algorithm. Turk J Electr Eng 2004;12(1). [22] Nelder JA, Mead R. A simplex method for function minimization. Comput J 1965;7(4):308–13.