Quasiparticle Swarm Optimization for Cross-Section Linear Profile Error Evaluation of Variation Elliptical Piston Skirt

Abstract

Variation elliptical piston skirt has better mechanical and thermodynamic properties and it is widely applied in internal combustion engine in recent years. Because of its complex form, its geometrical precision evaluation is a difficult problem. In this paper, quasi-particle swarm optimization (QPSO) is proposed to calculate the minimum zone error and ellipticity of cross-section linear profile, where initial positions and initial velocities of all particles are generated by using quasi-random Halton sequences which sample points have good distribution properties and the particles’ velocities are modified by constriction factor approach. Then, the design formula and mathematical model of the cross-section linear profile of variation elliptical piston skirt are set up and its objective function calculation approach using QPSO to solve the minimum zone cross-section linear profile error is developed which conforms to the ISO/1101 standard. Finally, the experimental results evaluated by QPSO, particle swarm optimization (PSO), improved genetic algorithm (IGA) and the least square method (LSM) confirm the effectiveness of the proposed QPSO and it improves the linear profile error evaluation accuracy and efficiency. This method can be extended to other complex curve form error evaluation such as cam curve profile.

1. Introduction

Piston skirt is the key parts of internal combustion engines (ICEs). Because internal combustion engines usually run under the circumstance of higher speed, larger pressure, and heavier load, it makes piston skirt work in overload conditions. Besides, piston skirt is also an important source of lubrication failures that will in turn lead to noise and power loss arisen from friction forces [1]. Therefore, it is necessary to improve its design, manufacture, measurement and evaluation method. With the rapid development of ICE industry, higher and higher design requirement for piston skirt is proposed for realizing high speed, high efficiency, low consumption, and low noise. The piston skirts of traditional formal cylinder and formal cone are seldom used, and they are mostly replaced with complex new-style profile piston skirts, especially for automobile engine. Among them convex variation elliptical piston skirt is widely applied that has the characteristics that the cross-section profiles at different heights are different ellipses and the axis-section profile is a convex curve. Compared with traditional one, the convex variation elliptical piston skirt has much better mechanical and thermodynamic properties. But the geometrical precision evaluation is a difficult problem because of its complex form. Recently many researchers have devoted themselves to develop different algorithms to compute the cross-section linear profile error of piston skirt. Huang et al. [2] proposed the algorithm based on single optimum seeking and moment method to evaluate cross-section linear profile error. In the hypothesis conditions of small deviation and small error with the measured data, least square cross-section linear profile error was calculated. Based on the combination of moment method and least square method, Huang and Wang [3] also put forward the evaluation method of piston skirt cross-section linear profile error. Liu et al. [4] set up the mathematical model of cross-section linear profile and deduced relevant designing formulae. The least square cross-section linear profile error was calculated. Nowadays most algorithms are based on LSM because of its ease of computation and the uniqueness of its solution. Because ISO 1101 (1996) recommends the form tolerance being evaluated based on the concept of minimum zone [5], the minimum zone method (MZM) has received much attention in recent years, and it has been applied to solve circularity (roundness), straightness, flatness, sphericity, and cylindricity error. Cheraghi et al. [6] formulated straightness and flatness errors by nonlinear optimization problems with linear objective function and nonlinear constraints. Samuel and Shunmugam [7] established the minimum circumscribed limacoid, maximum inscribed limacoid, and minimum zone limacoid in literature based on the computational geometry to evaluate sphericity error from coordinate measurement data. Weber et al. [8] proposed a unified linear approximation technique for use in evaluating the forms of straightness, flatness, circularity, and cylindricity. Non-linear equation for each form was linearized using Taylor expansion, and then it was solved as a linear program. Zhu and Ding [9] established the equivalence between the width of a point set and the inner radius of the convex hull of the self-difference of the set. An algorithm was proposed to calculate the “almost exact” minimum zone solution, which is implemented by solving a single linear programming problem. Li and Shi [10] applied the curvature technique for solving problems of roundness evaluation from coordinate data measured by CMM. Above methods are effective in solving simple form errors such as straightness, flatness, circularity, sphericity and cylindricity. Because the cross-section linear profile of variation elliptical piston skirt is more complex, it is difficult for traditional method to calculate the minimum zone error.

With the emergence of computational intelligence, the intelligence-oriented algorithms such as genetic algorithms (GAs) and particle swarm optimization (PSO) have been successfully employed to evaluate form error such as flatness, straightness, cylindricity, and so forth [11–13]. Wen attempted to calculate the minimum zone solution of piston skirt cross-section linear profile error by PSO [14]. Because there are several approximations in establishing profile error mathematical model, the solution is not very accurate.

In order to solve the minimum zone error of piston skirt cross-section linear profile accurately and efficiently, its mathematical model is founded, and QPSO is proposed. The paper is organized as follows: the design formula and mathematical model of the cross-section linear profile of variation elliptical piston skirt are set up. Then, QPSO for piston skirt cross-section linear profile error evaluation is developed. Finally, the examples are given, and conclusions are drawn.

2. Mathematical Model of Piston Skirt Cross-Section Linear Profile Error Evaluation

2.1. Design Formula

Set up the design coordinate system 𝑥^𝑂^𝑦^ of the cross-section linear profile of variation elliptical piston skirt, shown in Figure 1. 𝑄_𝑖 is the design point on the cross-section linear profile. The radial reduction Δ𝑙_𝑖 of the point 𝑄_𝑖 is usually formulated as [4]Δ𝑙_𝑖=𝐷−𝑑41−cos2𝛼^𝑖+2𝑏51−cos4𝛼^𝑖,(2.1) where 𝐷 is the diameter of long axis, 𝑑 is the diameter of short axis, 𝛼_𝑖 is the polar angle of the point 𝑄_𝑖, 𝑏 is the coefficient of plump degree, and 𝐺=𝐷−𝑑 is the ellipticity.

Figure 1 

The cross-section profile of variation elliptical piston skirt.

The design formula of the point 𝑄_𝑖 on cross-section linear profile is formulated as𝑙_𝑖=𝐷2−Δ𝑙^𝑖=𝐷2−𝐷4−𝑑1−cos2𝛼^𝑖+2𝑏51−cos4𝛼^𝑖,(2.2) where 𝑙_𝑖 is the radius of the point 𝑄_𝑖.

2.2. Mathematics Model of Cross-Section Linear Profile

The measurement model of cross-section linear profile of variation elliptical piston skirt is shown in Figure 2.

Figure 2 

The measurement model of cross-section profile.

𝑂 is the revolving centre of the measurement platform and 𝑂^ is the design center of piston skirt, 𝑒 is the setting eccentricity (𝑒=𝑂𝑂^) and 𝜃₀ is the eccentric angle, and 𝜙₀ is the angle between the measurement coordinate axis 𝑂𝑥 and the long axis 𝑂^𝑥^ of design profile (−10^∘≤𝜙₀≤10^∘). Assuming that 𝑃_𝑖(𝑟_𝑖,𝜃_𝑖) (𝑖=1,2,…,𝑛, 𝑛 is the number of measured point) is the measured point of the cross-section linear profile corresponding to the revolving centre 𝑂, and 𝑟_𝑖 and 𝜃_𝑖 are the radius and polar angle of point 𝑃_𝑖, respectively. 𝑃^𝑖(𝑟^𝑖,𝜃^_i) is the mapping point of 𝑃_𝑖, and 𝑟^𝑖 and 𝜃^𝑖 are the radius and polar angle of point 𝑃^𝑖 in the design coordinate system 𝑥^𝑂^𝑦^. 𝛽_𝑖 is the angle between 𝑃^𝑖𝑂^ and 𝑂^𝑥^, 𝛿_𝑖 is the angle between 𝑂𝑃^𝑖 and 𝑂^𝑃^𝑖. Because the setting eccentricity 𝑒 is very small, and 𝛿_𝑖 is also very small.

Using cosine theorem in the triangle Δ𝑃^𝑖𝑂𝑂^, we get the following:𝑟^𝑖2=𝑟^2𝑖+𝑒²−2𝑒𝑟_𝑖cos𝜃_𝑖−𝜃₀,(2.3) Equation (2.3) can be rewritten𝑟^𝑖=𝑟_𝑖2+𝑒²−2𝑒𝑟_𝑖cos𝜃_𝑖−𝜃₀.(2.4) From Figure 2, we can learn𝜃𝛽_𝑖𝑖=𝛿_𝑖+𝜂_𝑖,=𝜂_𝑖+𝜙₀.(2.5) So, we get the following:𝛽_𝑖=𝜃_𝑖−𝜙₀+𝛿_𝑖.(2.6) According to Taylor series expansion, we have the followingccooss42𝛽𝛽_𝑖𝑖=cos2𝜃_𝑖−𝜙₀+𝛿_𝑖≈cos2𝜃_𝑖−𝜙₀−2sin2𝜃_𝑖−𝜙₀⋅𝛿_𝑖,=cos4𝜃_𝑖−𝜙₀+𝛿_𝑖≈cos4𝜃_𝑖−𝜙₀−4sin4𝜃_𝑖−𝜙₀⋅𝛿_𝑖.(2.7) When 𝛽_𝑖=𝛼_𝑖, the radius design value 𝑙_𝑖 corresponding to the polar 𝛼_𝑖 can be rewritten𝑙_𝑖=𝐷2−𝐷4−𝑑1−cos2𝜃^𝑖−𝜙₀+2sin2𝜃_𝑖−𝜙₀⋅𝛿_𝑖+𝑏251−cos4𝜃^𝑖−𝜙₀+4sin4𝜃_𝑖−𝜙₀⋅𝛿_𝑖.(2.8) Using sine theorem in the triangle Δ𝑃^′_𝑖𝑂𝑂^, we get the following𝑒sin𝛿_𝑖=𝑟^𝑖sin𝜃_𝑖−𝜃₀.(2.9) Because 𝛿_𝑖 is very small, 𝛿_𝑖≈sin𝛿_𝑖, and it is substituted into (2.9), then (2.9) can be approximated as𝛿_𝑖=𝑒sin𝜃^𝑖−𝜃₀𝑟^𝑖.(2.10) Substituting (2.10) into (2.8), we have the following:𝑙_𝑖=𝐷2−𝐷4−𝑑1−cos2𝜃^𝑖−𝜙₀+2sin2𝜃^𝑖−𝜙₀⋅𝑒sin𝜃_i−𝜃₀𝑟^′_𝑖+𝑏251−cos4𝜃^𝑖−𝜙₀+4sin4𝜃^𝑖−𝜙₀⋅𝑒sin𝜃_𝑖−𝜃₀𝑟^𝑖.(2.11)

2.3. The Objective Function in Using QPSO to Calculate the Minimum Zone Error

The deviation 𝜀_𝑖 between the polar radius 𝑟^𝑖 of the mapping point 𝑃^𝑖 and the polar radius 𝑙_𝑖 of the design point 𝑄_𝑖 corresponding to the same polar angle is𝜀_𝑖=𝑟^𝑖−𝑙_𝑖=𝑟^𝑖−𝐷2−𝐷4−𝑑1−cos2𝜃^𝑖−𝜙₀+2sin2𝜃^𝑖−𝜙₀⋅𝑒sin𝜃_𝑖−𝜃₀𝑟^𝑖+𝑏251−cos4𝜃^𝑖−𝜙₀+4sin4𝜃^𝑖−𝜙₀⋅𝑒sin𝜃_𝑖−𝜃₀𝑟^𝑖,(2.12) where 𝑟^′_𝑖=𝑟^2𝑖+𝑒²−2𝑒𝑟_𝑖cos(𝜃_𝑖−𝜃₀).

According to the ISO/1101 standard, the minimum zone solution of linear profile error is the minimum width value of two ideal equidistance design profiles which encompass the measured real profile. Therefore, the objective function in using QPSO to calculate the minimum zone error of cross-section linear profile can be expressed as:𝑓𝜃₀,𝑒,𝜙₀,𝛽,𝐷,𝑑=minmax𝜀_𝑖−min𝜀_𝑖.(2.13) Equation (2.13) is a function of (𝜃₀,𝑒,𝜙₀,𝛽,𝐷,𝑑). Consequently, solving the minimum zone cross-section linear profile error of variation elliptical piston skirt is translated into searching the values of the parameters set (𝜃₀,𝑒,𝜙₀,𝛽,𝐷,𝑑), so that the objective function 𝑓(𝜃₀,𝑒,𝜙₀,𝛽,𝐷,𝑑) is the minimum.

3. Using QPSO to Evaluate Piston Skirt Cross-Section Linear Profile Error

3.1. Pseudorandom Numbers and Quasirandom Halton Sequences

3.1.1. Pseudorandom Numbers and Quasirandom Sequences

Pseudorandom numbers are deterministic, but they try to imitate an independent sequence of genuine random numbers. Common pseudorandom number generators include linear congruential, quadratic congruential, inversive congruential, parallel linear congruential, et.al. In contrast to pseudorandom numbers, the points in a quasirandom sequence do not imitate genuine random points. But they try to cover the feasible region in an optimal way. Quasirandom generators do not generate numbers, but sequences of points in the desired dimension. Common quasirandom sequence generators include Halton, Hammersley, Faure, Sobol, and Niederreiter generators [15]. In this paper we focus our attention on Halton sequence since it is conceptually very appealing, and it can be produced easily and fast with simple algorithms.

3.1.2. Halton Sequences

Halton sequences are not unique, and they depend on the set of prime numbers taken as bases to construct their vector components. Typically and most efficiently, the lowest possible primes are used.

Let 𝑏 be a prime number. Then any integer 𝑘,𝑘≥0, can be written in base-𝑏 representation as𝑘=𝑑_𝑗𝑏^𝑗+𝑑_𝑗−1𝑏^𝑗−1+⋯+𝑑₁𝑏+𝑑₀,(3.1) where 𝑑_𝑖∈{0,1,…,𝑏−1},𝑖=0,1,…,𝑗. Define the base-𝑏 radical inverse function, 𝜙_𝑏(𝑘), as𝜙_𝑏(𝑘)=𝑑⁰𝑏¹+𝑑¹𝑏²+⋯+𝑑^𝑗𝑏^𝑗+1.(3.2) Notice that for every integer, 𝑘≥0,𝜙_𝑏(𝑘)∈[0,1].

The 𝑘th element of the Halton sequence is obtained via the radical inverse function evaluated at 𝑘. Specifically, if 𝑏₁,…,𝑏_𝑑 are 𝑑 different prime numbers, then a 𝑑-dimensional Halton sequence of length 𝑚 is given by {𝑥₁,…,𝑥_𝑚}, where the 𝑘th element of the sequence is𝑋_𝑘=𝜙_𝑏1(𝑘−1),…,𝜙_𝑏𝑑(𝑘−1)^𝑇,𝑘=1,…,𝑚.(3.3)

3.2. QPSO for Evaluating Cross-Section Linear Profile Error of Piston Skirt

Particle swarm optimization (PSO) method is one of the most powerful methods for solving unconstrained and constrained global optimization problems. The method was originally proposed by Kennedy and Eberhart as an optimization method in 1995 [16], which was inspired by the social behavior of bird flocking and fish schooling. It utilizes a “population” of particles that fly through the problem hyperspace with given velocities [17]. In PSO initial position and initial velocity of particles are often randomly generated by using pseudorandom numbers [13, 18]. Because the positions of initial particles have influence on the optimization performance, Richard and Ventura [19] proposed initializing the particles in a way that they are distributed as evenly as possible throughout the problem space. This ensures a broad coverage of the search space. They concluded that applying a starting configuration based on the centroidal Voronoi tessellations (CVTs) improves the performance of the PSO compared with the original random initialization. As an alternative method, Campana et al. [20] proposed reformulating the standard iteration of PSO into a linear dynamic system. The system can then be investigated to determine the initial particles’ positions such that the trajectories over the problem hyperspace are orthogonal, improving the exploration mode and convergence of the swarm.

Quasirandom sequences have been successfully applied in numerical integration and in random search optimization methods [21]. The idea of a good initial population has also been used in genetic programming and genetic algorithm [22]. In this work, quasirandom Halton sequences are applied to generate the initial positions and velocities of particles in PSO for solving the minimum zone profile error of variation elliptical piston skirts. For short, the proposed PSO is called quasiparticle swarm optimization (QPSO).

QPSO algorithm begins by using quasirandom Halton sequences to initialize a swarm of 𝑁 particles (𝑁 is referred as particle size), each having 𝑠 unknown parameters (𝑠 is referred as the dimensionality of optimized variables) to be optimized at each iteration. The ideal cross-section linear profile can be decided by the set of six parameters (𝜃₀,𝑒,𝜙₀,𝛽,𝐷,𝑑) and the method takes (𝜃₀,𝑒,𝜙₀,𝛽,𝐷,𝑑) as a particle. Therefore, the dimension 𝑠 of the particle is six. The best particle with the minimum objective function value 𝑓(𝜃₀,𝑒,𝜙₀,𝛽,𝐷,𝑑) is considered as the minimum zone solution to the cross-section linear profile error. The flow of QPSO for evaluating cross-section linear profile error is as follows.

Step 1. Input the measurement values (𝑟_𝑖,𝜃_𝑖) (𝑖=1,2,…,𝑛) of the cross-section linear profile. If the point is measured in the Cartesian coordinates, it needs to be transformed into the polar coordinates.

Step 2. Generate the initial positions and initial velocities of all particles by using quasirandom Halton sequences.

Step 3. Calculate the objective functions of all particles according to (2.13). The less the objective function value is, the better the particle is.

Step 4. Update velocity. Because constriction factor approach (CFA) ensures the convergence of the search procedures based on the mathematical theory and can generate higher-quality solutions [23], CFA is employed to modify the velocity. The velocity and position parameters of each particle (𝑝_𝑖) in the swarm are updated at iteration (𝑡) according to CFA:

𝑣^𝑖𝑡+1=𝐾𝑣^𝑡𝑖+𝐾𝑐₁𝑟=_1𝑗𝑝𝑏𝑒𝑠𝑡^𝑡𝑖−2𝑝^𝑡𝑖+𝑐₂𝑟_2𝑗𝑔𝑏𝑒𝑠𝑡^𝑡−𝑝^𝑡𝑖,|||2−𝜑−√𝜑²−4𝜑|||,(3.4) where 𝑣^𝑡𝑖 and 𝑝^𝑡𝑖 are the velocity and position of 𝑖th particle at iteration t, respectively. 𝑟_1𝑗 and 𝑟_2𝑗(𝑗=1,2,…,𝑠) are uniform random numbers between 0 and 1. 𝑐₁, and 𝑐₂ are acceleration factors that determine the relative pull for each particle toward its previous best position (𝑝best) and the group’s best position (𝑔best), respectively. 𝑐₁ and 𝑐₂ meet the conditions 𝜑=𝑐₁+𝑐₂.

Step 5. Update position. The position of each particle is modified by 𝑝^𝑖𝑡+1=𝑝^𝑡𝑖+𝑣^𝑡𝑖Δ𝑡.(3.5)

Step 6. Update 𝑝best. Calculate the objective function of all particles. If the current objective function value of a particle is less than the old 𝑝best value, the 𝑝best is replaced with the current position.

Step 7. Update 𝑔best. If the current objective function value of a particle is less than the old 𝑔best value, the 𝑔best is replaced with the current position.

Step 8. Go to Step 4 until the maximum iteration is satisfied.

Step 9. Output the computation results.

4. Results and Discussion

4.1. Optimizing Classical Testing Functions

In order to verify the optimization efficiency of QPSO, numerical experiments on some classical testing functions are carried out [24]. Two examples are given as follows.

Function 1
𝑓₁(𝑥)=100×𝑥₂−𝑥²¹²+1−𝑥₁²,−2.048≤𝑥₁,𝑥₂≤2.048.(4.1) It is hard to be minimized. The global minimum point is at (1.0, l.0), and the global minimum is zero.

Function 2
It is the Schaffer test function defined as 𝑓₂(𝑥)=0.5+sin²𝑥²¹+𝑥²²−0.51.0+0.001×𝑥²¹+𝑥²²²,−100≤𝑥_𝑖≤100.(4.2) This function has many circle ridges nearby the global minimum (1.0, l.0), and the function value of the nearest circle ridge (𝑥²¹+𝑥²²=3.138²) is 0.009716. It is very easy to trap in this value.

The proposed algorithms were written in MATLAB, and the experiments were run in Windows XP on an IBM ThinkPad X200-7457 A46 with 2.26 GHz main frequency and 1 GB memory. QPSO is also a stochastic optimization method and it is important to evaluate the average performance. For comparison, two stochastic optimization methods including PSO [13] and IGA [11] are employed. The popsize size was set 50, and 20 trials were performed in prescribed maximum iteration 200. In specified initial ranges, initial populations and initial positions were randomly generated by using pseudorandom numbers for PSO and IGA. Initial populations were generated by using Halton random sequences for QPSO. The mean values and the standard deviations are tabulated in Table 1.

Table 1 

Mean and standard deviation of functions.

As seen in Table 1, the QPSO method for two examples could provide more accurate and stable solution. Figures 3 and 4 show the optimizing processes of these methods at one trial for two examples, respectively. As seen in the figures, it is evident that the optimization performance of QPSO is better than those of PSO and IGA.

Figure 3 

The evolution process of function 𝑓1 by three different methods.

Figure 4 

The evolution process of function 𝑓2 by three different methods.

4.2. Examples

Simulation Example
According to the design formula of cross-section linear profile of variation elliptical piston skirt, the simulation data with random noise are generated. The setting eccentricity and the angle between the measurement coordinate axis and the long axis of design profile are set by the coordinate translation and rotation transform. In the experiment, the design data of 𝜃₀,𝑒,𝜙₀,𝛽,𝐷,𝑑 are listed in Table 2, and the transformed simulation data are shown in Table 3. For comparison, IGA and PSO were employed. Considering the values of optimized parameters, 𝜃₀,𝑒,𝜙₀,and𝛽 are very small, and 𝐷 and 𝑑 are usually larger, in order to save searching time, the initial populations and initial positions were randomly generated by using pseudorandom numbers in: ([−𝜀,𝜀],[−𝜀,𝜀],[−𝜀,𝜀],[−𝜀,𝜀],[−𝜀+max(𝑟_𝑖),𝜀+max(𝑟_𝑖)],[−𝜀+min(𝑟_𝑖),𝜀+min(𝑟_𝑖)]) for IGA and PSO. For QPSO the initial positions were generated by using quasrandom Halton sequence in ([−𝜀,𝜀],[−𝜀,𝜀],[−𝜀,𝜀],[−𝜀,𝜀],[−𝜀+max(𝑟_𝑖),𝜀+max(𝑟_𝑖)],[−𝜀+min(𝑟_𝑖),𝜀+min(𝑟_𝑖)]) and the initial velocities were generated by using quasirandom Halton sequence in:([−𝜀,𝜀],[−𝜀,𝜀],[−𝜀,𝜀],[−𝜀,𝜀],[−𝜀,𝜀],[−𝜀,𝜀]). In our experiments, 𝜀 is all set 0.5. The searching process and optimization results of the minimum zone error of cross-section linear profile by different methods are shown in Figure 5 and Table 2. As seen in Figure 5 and Table 2, the minimum zone error by QPSO is 0.0775 mm and it is smaller than that by PSO and IGA. It takes about 40 iterations for the proposed QPSO to find the optimal solution and it is faster than PSO and IGA. The cross-section linear profile error by LSM is 0.0932 mm and it is larger than the minimum zone error.

Table 2 

Results of simulate example.

Table 3 

Simulated measurement data.

Figure 5 

The evolution process by different methods.

Practical Example
The cross-section profiles of piston skirt of a SL 105 diesel engine are variation ellipses and its main cross-section profiles (upper end, convexity and lower end) are inspected by Coordinate Measurement Machine (CMM). And the minimum zone error of every cross-section linear profile is calculated by the proposed QPSO and the results are listed in Table 4. For comparison, the ellipticities of three cross-sections are calibrated by 19JPC microcomputer-type all-purpose tool microscope and the values are also listed in Table 4. From the table, we can learn the cross-section linear profile error of MZM is less than that of LSM. And the ellipticity calculated by QPSO is almost the same as the calibration value.

Table 4 

Results of practical example.

5. Conclusions

In this paper, QPSO is proposed to calculate the minimum zone error and ellipticity of cross-section linear profile of variation elliptical piston skirt, which initial positions and initial velocities of all particles are generated by using quasirandom Halton sequences and the particles’ velocities are modified by constriction factor approach. The design formula and mathematical model of the cross-section linear profile are set up and its objective function calculation approach using QPSO to solve the minimum zone error of cross-section linear profile is developed. The simulation and practical examples confirm the optimization efficiency of QPSO is better than that of PSO and IGA for complex optimal problems. Compared with conventional evaluation methods, the proposed method not only has the advantages of simple algorithm and good flexibility, but also improves cross-section linear profile error evaluation accuracy. The proposed method can be extended to other complex curve profile error evaluation.

Acknowledgments

The research was supported by Natural Science Research Project of the People’s Republic of China (no: 51075198), Natural Science Research Project of Jiangsu Province (BK2010479) Jiangsu Provincial Project of 333 Talents Engineering of China, Jiangsu Provincial Project of Six Talented Peak of China.