Micro-genetic approach for surface meshing on a set of unorganized points

Micro-genetic approach for surface meshing on a set

of unorganized points

Da-Pan Chen

, Zhen-Yu Hsu

, Fang-Pang Lin

a_{Department of Mechanical Engineering, National Chiao Tung University, 1001 Ta Hsueh Road, Hsinchu 300, Taiwan, ROC} b_{National Center for High-performance Computing, 7, R&D Rd. VI, Hsinchu Science Park, Hsinchu 300, Taiwan, ROC}

Received 18 November 2004; received in revised form 2 March 2006; accepted 1 November 2006

Abstract

Surface meshing plays a crucial role in mesh generation. Usually, surface meshing methods in three dimensions generate meshes rely-ing on prescribed patch interpolation. In this study, an approach of surface meshrely-ing directly on a set of unorganized points is developed, which consists of a mesh triangulation and a conversion scheme for primary triangular and quadrilateral surface meshes, a high order C1 surface function reconstruction, and a micro-genetic algorithm (MGA) to smooth the meshes. Practical cases are given to demonstrate its successful performance and its versatility.

 2006 Elsevier B.V. All rights reserved.

Keywords: Surface meshing; Unorganized points; C1surface function reconstruction; Micro-genetic algorithm

1. Introduction

Mesh generation is a critical step in a wide range of engi-neering practices, such as the study of scientific simulation by finite element analysis (FEA), of design automation with geometrical modeling, and of data presentation with graphics and visualization. It is well known that a quality mesh is imperative to the above study. In the engineering practices, meshes can be generated on the boundary and/ or in the interior of the object, which correspond to surface meshing and volume meshing, respectively. Usually, a well-defined geometrical model is requested before the mesh generation. In this paper, we present a different approach that surface meshing is on a set of unorganized points in which the surface function is yet to be defined.

For surface meshing, most existing methods[29] gener-ate the surface mesh based on pre-defined surface func-tions, either in parametric patches or in algebraic form,

by using existing mesh generation schemes such as Dela-unay Triangulation [20] or Advancing Front Methods

[21]. However, in many real cases, the given data may just be a set of unorganized points, in which its surface function cannot be devised in a usual fashion. It is often encoun-tered in the applications of biomechanics, which require a geometrical reconstruction from a set of sampling points that is extracted from a sequence of scanned images, e.g. histological sections in tomography. An immediate way to generate a surface mesh for this application is to triangu-late the given points[1–4]. However, some sets of such sam-pling points can be scattered and irregular, and results in locally ill-posed meshes, which may not acceptable for finite element analysis. In this study, a posteriori approach is adopted to tackle this problem: the given points are first triangulated to generate the triangular mesh and/or con-vert to the quadrilateral mesh and then an additional pro-cedure is introduced to enhance the quality of the mesh. For finite element meshes, the common used procedure for this enhancement is mesh smoothing.

The most popular mesh smoothing method is Laplacian smoothing[5], in which every internal grid node is

reposi-0045-7825/$ - see front matter  2006 Elsevier B.V. All rights reserved. doi:10.1016/j.cma.2006.11.003

*

Corresponding author. Tel.: +886 3 5776085x360; fax: +886 3 5773538.

E-mail address:[email protected](F.-P. Lin).

tioned at the geometrical center of the adjacent nodes. Generally, for surface mesh smoothing, Laplacian smooth-ing is first employed on a parametric plane and then maps the result onto the physical surface domain. It is well known that the mapping between a reconstructed surface and its parametric plane for Laplacian smoothing strongly affects the resultant mesh quality. A given unorganized point set is usually lack of geometrical regularity in distri-bution. It cannot guarantee to form a proper primary mesh for the mapping of the Laplacian smoothing.

Optimization-based mesh smoothing technique is another way to smooth the finite element mesh, in which the new location of node is found by using the optimization algorithms [6–9]. Freitag [6] presented several articles described about the finite element mesh smoothing tech-niques, which contained smart Laplacian smoothing, opti-mization-based smoothing and the combination of both, in two-dimensional plane mesh and three-dimensional tetrahedral mesh. For the combination method, the smart Laplacian smoothing is used to adjust every internal node and is followed by the optimization-based algorithm, the steepest decent method[10], in only the poorest-quality ele-ments[6]. The steepest decent method searches an optimal solution by a given initial step along a search direction, which is calculated from the gradients of the object func-tion, i.e. mesh quality measurement. However, the calculat-ing of gradients accordcalculat-ing to the real coordinates is inconvenient for this application, in which the surface func-tion is yet to be defined. Besides, the search direcfunc-tion may let the nodes to deviate from the original surface.

Garimella et al. [8] presented a surface mesh quality optimization procedure that the nodes are repositioned based on element-based local parametric spaces. They employ the conjugate gradient method [10], whose search direction is estimated by computing the gradients of the mesh quality measurement with respect to the local para-metric spaces, to reposition the nodes to enhance the mesh quality. It is beneficial that mesh smoothing based on local parametric spaces can remain the nodes close to the origi-nal surface. The calculation of the gradients of object func-tion, which is based on the local parametric spaces and without the need of surface functions, is suitable for this application. However, two problems will be arisen: (1) the repositioned points are re-allocated on the planes of corresponding elements, that are not projected onto the original surface, and this will affect the geometrical accu-racy of model. (2) The gradient search methods as a local search method may encounter local optimum problem.

In this study, we propose an innovation surface mesh smoothing procedures that mesh nodes were repositioned by using the micro-genetic algorithm (MGA) based on local reconstructed surface. The MGA [13–16] is similar to the genetic algorithms (GA) [11,12], which is a global search method that searches optimal solution by employing natural evolution without calculating search direction and step size. The MGA works with small population size and reaches new optimal regions much earlier than the

con-ventional GA implementation[13]. It has been successfully applied to many fields [13–16]. Moreover, in order to ensure the geometrical accuracy of the analytical model, we projected the repositioned nodes onto the original sur-face based on an interpolation sursur-face function[18], which is reconstructed from the primary triangular elements. The one drawback of our approaches is that the computational cost of MGA is larger than the gradient search methods

[10]but it is feasible by applying parallel computation algo-rithms [19] to accelerate its computational efficiency. The computational efficiency will not be discussed in this study. The outline for the rest of the paper is as follows. Section 2 presents surface mesh generation and surface function reconstruction algorithm. Section3describes the details and procedures of our proposed MGA approach. Numerical results and discussions are given in Section 4

and the conclusions in Section5.

2. Mesh generation and surface function reconstruction 2.1. Mesh generation

The generation of surface mesh based on unorganized points set is necessary in some science and engineering fields, where geometrical data are often measured or gen-erated at isolated and unorganized positions, such as mentioned earlier of the biomedical research. In this application, the given data are just the points set and the surface function is yet to be defined. Therefore, the surface mesh cannot be generated based on its surface function directly. For this application, a common way to recon-struct the surface model is to triangulate the given points

[1–4]. The triangulation procedure aims at generating a pri-mary triangular surface mesh as well as creating back-ground triangular patches for the use in surface function reconstruction procedure. Furthermore, since the distribu-tion of the given points may be irregular over the surface of the model, the primary triangular meshes always contain some ill-posed triangles. Therefore, some mesh cleanup operators[23]were introduced to improve the topological connectivity of the triangular meshes, and then MGA approach was applied to enhance the mesh quality further. Once the primary triangular surface mesh was created, the quadrilateral surface mesh can be generated based on the triangular one. The conversion scheme [26–30] was employed to serve this purpose. It is a common and conve-nient way to generate an unstructured quadrilateral mesh. The quadrilateral mesh was created by a careful process to merge two adjoining triangles to form a quadrilateral ele-ment. However, the conversion scheme usually introduces plenty of ill-posed quadrilaterals. To improve the mesh quality, a two-stage procedure is required. First, mesh structure modification (topological improvement) opera-tors[24,25], such as cleanup, edge swapping, node elimina-tion and element collapse, were employed to refine the mesh connectivity. Then the mesh quality was further enhanced by applying the MGA approach.

2.2. Surface function reconstruction

In order to generate a proper surface mesh and ensure the geometrical accuracy of the analytical model, the sur-face function is necessary during the sursur-face mesh smooth-ing procedures. In some previous articles, the mesh smoothing procedures were applied to reposition the inter-nal nodes on the plane (or tangent plane) of the primary elements, which is C0-continue [8,26]. It is based on an assumption that the primary surface elements were well matched to the original surface. However, in this study, since the given data points may be chosen irregularly, the primary surface elements, which were triangulated directly from the given data points, may not match the original sur-face well. Therefore, in this study, a C1 continue surface function reconstruction algorithm was adopted to ensure the geometrical accuracy during the nodes repositioning.

There are many surface function reconstruction meth-ods developed [18,31–36]. For a finite element analytical model, a C1-continuous surface function is necessary for sufficient numerical accuracy. Here, a C1triangular patch interpolation method developed by Goodman and Said

[18]was adopted to reconstruct the surface function. It is a simpler and efficient method for the surface function reconstruction. In this method, surface function is recon-structed by local cubic Bezier triangular patches. The gra-dients of vertices are necessary for this surface function reconstruction procedure. We adopted a local derivative estimation method, which is also developed by Goodman et al.[17], to calculate the gradients of vertices. Please refer to[17]and[18]for the details of the vertices gradients cal-culation and surface functions reconstruction, respectively. The comparison of geometrical accuracy for the original C0-continue and our enhanced C1-continue surface can be shown as follows: 36 points and triangulation in Whelan

[31]were chosen and two test functions were employed: F1ðx; yÞ ¼ 0:75 expððð9x  2Þ 2 þ ð9y  2Þ2Þ=4Þ þ 0:75 expðð9x þ 1Þ2=49 ð9y þ 1Þ=10Þ þ 0:5 expððð9x  7Þ2þ ð9y  3Þ2Þ=4Þ  0:2 expðð9x  4Þ2 ð9y  7Þ2Þ F2ðx; yÞ ¼ ð1:25 þ cosð5:4yÞÞ=ð6 þ 6ð3x  1Þ2Þ

The interpolated values of the test functions at a 25· 25 uniform mesh points in a unit square were computed and the maximum and mean errors were listed in Table 1. The surface error, shown in Table 1, indicates that the

C1-continuous surface is more accuracy than the C0 -con-tinue surface both on maximum and mean surface errors. Therefore, surface mesh smoothing based on the recon-structed C1surfaces is beneficial to enhance the geometrical accuracy of the analytical model.

3. Genetic algorithm—micro-genetic approach

The genetic algorithm is one of the recent popular evo-lutionary algorithms. The GA begins with a set of initial individuals called the population that represent the candi-date solutions, which were coded as a bit string to simulate the gene of nature, of the problem. At the beginning, sev-eral strings are created randomly to form an initial popula-tion. As per Darwin’s theory of evolution, the population is evolved generation after generation to search an optimal solution and historical information is then exploited to speculate on new search points with expected performance during the iteration [11]. The GA works with the popula-tion constituted by the coding parameter set and searches for a global optimal solution instead of a local one.

The GA yields good results via three major genetic oper-ators: selection, crossover, and mutation. First, the selec-tion operator selects the fittest strings of the previous generation to be the parents of the new generation. This is used to preserve the better historical information to survive at the new generation. After the selection, the crossover operator generates new strings by crossing pair parents of the old. The new strings succeed and exchange the best information of parents to be the new individuals. To avoid the loss of some important genes and increase the variation of the individuals, mutation operator is imported to add new information occasionally.

3.1. Micro-genetic algorithm (MGA)

The standard genetic algorithm is successfully applied to many different applications [11,12]. However, one major drawback is that the iterative global searching of the algo-rithm is time consuming. It will be deteriorating when addi-tional iterations are needed in the smoothing procedure. There are many approaches to tackle the problem. For the genetic algorithm per se, reduction of the population size is an effective one. For conventional GA, the usual choice of population size is based on the conception that larger population relates to better schema processing, lesser chance of premature convergence, and better optimal results [13]. This imposes a considerable loading on the computational time. To trade-off, the micro-genetic algo-rithm [13,16]is particularly adopted to accelerate the con-vergence of the conventional GA. The MGA is similar to the GA that proceeds with binary coded population and employs the selection and crossover operations to evolve population for generations, but with smaller population size than conventional GA. It had been reported that MGA reaches near optimal regions much earlier than the standard GA does[13].

Table 1

Surface errors of the test functions for C0and C1surface reconstruction Test function C0surface C1surface Maximum error Mean error Maximum error Mean error F1(x, y) 0.215184 0.027594 0.120067 0.022945 F2(x, y) 0.059769 0.008789 0.035420 0.005101

It is well known that the GA works poorly with small population size due to insufficient information processing, which results in premature convergence to local optimal solutions. For the MGA, the best individual is passed to the new generation to ensure that the good individual is held, and it requires multiple convergences. The best indi-vidual of the old is remained and the others are randomly generated after each convergence. This operation is used to add new information and avoid premature convergence, and the mutation rate is set to zero.

In this study, an unorganized points set were triangu-lated to form a primary surface triangular mesh and/or converted the primary triangular mesh into quadrilateral mesh. The primary surface mesh was first refined by the mesh structure modification operators [23–25], and then the MGA mesh smoothing procedures were applied to enhance the surface mesh quality further. The procedures of our surface mesh smoothing and the MGA adopted in our scheme were summarized as follows:

Step 1. Input data: input the data of primary mesh, which includes surface function, node positions and element connectivity. If the surface function is not given, the gradients of nodes will estimated by the gradients estimation procedures[17]. Step 2. MGA mesh smoothing begins:

Step 3. Search the optimal solutions within each adjacent element by using steps 4–7 and then choose the best one to be the new position of the node.

Step 4. Initial population: the MGA requires multiple con-vergences. According to Ref.[13], a population size of five is chosen in each convergence. The best indi-vidual of the previous generation will be held. The others are generated randomly.

Step 5. Decode the strings and calculate their node posi-tions based on the reconstructed surface function (Section 2.2). Calculate their fitness values and then carry the best string to the next generation. Step 6. Select four strings (contains the best string) for

reproduction by employing the roulette wheel strategy[11]. Generate four individuals by employ-ing the crossover operator with probability of one

[13].

Step 7. Check the convergence criterion. If it is not conver-gence, go to step 5 or else go to step 3 or step 4. Step 8. Go to step 2 to smooth next node until the end of

the smoothing.

Step 9. Check the convergence of the whole mesh smooth-ing procedures. If it is not convergence, go to step 2 or else end off the smoothing procedures.

3.2. Design parameters

In the MGA mesh smoothing procedures, the design parameters were chosen to represent the nodes position. In order to avoid degenerate elements, the search space was restricted within a triangular area for each adjacent

element in both triangular and quadrilateral mesh smooth-ing. To represent a node lay on a triangle, we choose two parameters r1and r2, which relate to the barycentric

coor-dinates as the following: Consider a triangle T (Fig. 1) with vertices V1and V3in barycentric coordinates u, v and w,

such that any point on the triangle can be expressed as V ¼ uV1þ vV2þ wV3; uþ v þ w ¼ 1 ð1Þ

Now, with 0 < ri< 1, i = 1, 2, the barycentric coordinates

can be given as

u¼ r1ð1  r2Þ þ ð1  r1Þð1  r2Þ ¼ 1  r2

v¼ r2ð1  r1Þ

w¼ r1r2

ð2Þ

Substituting(2) into(1), we get

V ¼ ð1  r2ÞV1þ r2ð1  r1ÞV2þ r1r2V3 ð3Þ

As shown inFig. 1, the vertices V1,V2and V3are collinear

when r2= 1. To avoid that, we let 0 < r2< 0.5. After the

position of point V is obtained, the exact position is calcu-lated by mapping it to the original surface based on its cor-responding triangular patch and the local reconstructed surface function (Section2.2).

3.3. Fitness function

In this study, the fitness function is a criterion to judge the mesh quality and defined with the mesh quality mea-surement index. For the triangular mesh quality measure-ment, the common used quality index is

a¼ 2pffiffiffi3 kAB  ACk

kABk2þ kBCk2þ kACk2 ð4Þ where A, B, C are the vertices of a triangle. According to this mesh quality measurement, the fitness function for the triangular mesh smoothing Ftis defined as

Ft¼

Xn i¼1

ai ð5Þ

where n is the number of adjacent triangular elements.

For the quadrilateral mesh quality measurement, Knupp [37] presented an algebraic mesh quality metrics from the Jacobian matrix. The definition of quadrilateral shape quality metric is as follows: for a plane quadrilateral mesh, let the coordinates of the four nodes be (xk, yk),

k = 0, 1, 2, 3. The Jacobian matrices, Ak, one at each node

of the quadrilateral: Ak ¼ xkþ1 xk xkþ3 xk y_kþ1 yk ykþ3 yk   ð6Þ where the indices k + 1 and k + 3 are taken modulo four, for example, if k = 1 then k + 3 becomes 0. Four metric tensors are obtained by the combinations AT_kAk. Let

kk

ij; i; j¼ 1; 2, be the ijth component of the kth metric

ten-sor. Geometrically, at the kth node, kk₁₁ is the square of the length of the side connecting nodes k and k + 1, kk₂₂ is the square of the length of the side connecting nodes k and k + 3. Let hkbe the angle between the two sides joined

at the kth node, the quadrilateral shape quality metric can be expressed as b¼_P3 8 k¼0ð1 þ r2kÞ=ðrksin hkÞ ð7Þ where r¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffik22=k11 p

is the length ratio of kth node. This concept can be extended to measure quadrilateral surface mesh quality directly. Then the fitness function for the quadrilateral mesh smoothing Fqis defined as

Fq¼

Xn i¼1

bi ð8Þ

where n is the number of adjacent quadrilateral elements. 4. Results and discussion

The proposed approach is tested in several unorganized point datasets. A simple geometry of saddle shape with pre-scribed surface function is used to validate the procedure and to compare the performance of our MGA approach and the conjugate gradient method (Table 2). It is then applied to complicated geometries, such as wing-fuselage, which is often used in preliminary aircraft design, and biological dataset of shapes of a foot and a rat, which are constructed from contours by a three-dimensional laser scanning and from a sequence of segmented bio-images, respectively. The geometrical reconstruction of bio-images is a crucial pre-processing for physical modeling in biology and biomechanics. So, the applications are used to

demon-strate not only the effectiveness but also the practical use of our approach. In order to show the capability of our MGA mesh smoothing approach, all of the smoothed results are obtained by smoothing the primary meshes directly without any additional mesh treatment. The mesh quality is mea-sured according to Section 3.3 that the triangular mesh quality index is a and the quadrilateral mesh quality index is b. The results are summarized inTables 3 and 4for trian-gular and quadrilateral surface meshes, respectively. Signif-icant enhancement is found by our approach. Furthermore,

in Table 3, the worst quadrilateral mesh quality is set as

0.0001, which represents that one of the interior angle is greater than 179 and the codes are run in Linux PC with a dual core AMD 2 GHz CPU and 3 GHz RAM.

The first example is a basic mathematical function of saddle shape (Fig. 2), which allows us to scrutinize the per-formance of the approach. The surface function is denoted as the following:

z¼ x2

=10 y2

=4 ð9Þ

The saddle shape surface, unlike the usual well-behaved elliptic shape, has negative curvature, which is appealing to be used as a test case for the surface smoothing

[38]. The triangular surface mesh is generated based on

Table 2

The comparison between conjugate gradient method and MGA method using the saddle geometry

Method Object Number of node Number of element Worst quality Mean quality CPU time (s)

Conjugate gradient Triangular saddle mesh 601 1104 0.6834 0.9562 0.8000

MGA 0.4591 0.9509 3.8199

Conjugate gradient Quadrilateral saddle mesh 597 548 0.4099 0.8736 0.8400

MGA 0.5696 0.9001 2.3199

Table 3

Mesh quality improvement for the triangular surface meshes

Object Number of node Number of element Worst quality Mean quality CPU time (s) Wing-fuselage 1419 2432 0.2535 0.8483 3.0500 0.2535 0.9046 Foot 4039 8000 0.3803 0.8275 23.8200 0.3469 0.9192 Rat 25,670 51,354 0.2225 0.8963 150.0509 0.2460 0.9453 Table 4

Mesh quality improvement for the quadrilateral surface meshes

Object Number of node Number of element Worst quality Mean quality CPU time (s) Wing-fuselage 1405 1200 0.0001 0.6982 4.1000 0.4215 0.8580 Foot 4010 3971 0.0001 0.7049 19.6400 0.4111 0.8632 Rat 25,374 25,383 0.0001 0.7568 88.8099 0.2553 0.8899

Delaunay triangulation and the quadrilateral surface mesh converted from the triangular mesh (Section2.1). The per-formance of MGA can be observed by comparing with that of the conjugate gradient method in the saddle shape.

Ta-ble 2shows that the improvement of the mean quality of

triangular mesh reaches to 0.9562 and 0.9509 by both approaches, but the CPU time of our MGA approach is approximate five times of the conjugate gradient method. This is one drawback of GA, but can be easily tackled by parallelism and there are several parallel GA[19]developed with great success. The approach is similarly employed to the quadrilateral mesh. Our quadrilateral meshes are

gener-ated using a popular conversion scheme or fission scheme

[26–30]. The conversion scheme essentially merges two

neighboring triangles to form a new quadrilateral. This method may introduce ill-posed quadrilaterals, and needs further treatments[24,25]for practical use.Table 2shows that the improvement of mean quality of quadrilateral mesh reaches to 0.8736 and 0.9001 by employing the con-jugate gradient method and our MGA approach, res-pectively. It is clear that our performance as a global method achievement for mesh quality improvement is better than the conjugate gradient method. The enhance-ments of the mesh quality can be also visually observed

fromFigs. 2 and 3for triangular and quadrilateral surface

mesh, respectively.

Fig. 2. Triangular surface mesh of saddle: (a) original surface mesh; (b) surface mesh after MGA smoothing; and (c) surface mesh after conjugate gradient smoothing.

Fig. 3. Quadrilateral surface mesh of saddle: (a) original surface mesh; (b) surface mesh after MGA smoothing; and (c) surface mesh after conjugate gradient smoothing.

Fig. 4. Triangular surface mesh of wing-fuselage: (a) original surface mesh and (b) surface mesh after smoothing.

The following example is a complicated geometry of wing-fuselage configuration of NCKU-ILD-101[39]. Such kind of geometrical model is often used in CFD simulation for aircraft preliminary design. The usual practice for gen-erating a surface mesh for such a geometrical model is based on a pre-defined patch interpolation. The current data are prepared from a point datasets that is well gener-ated by usual mesh generation methods of Delaunay trian-gulation. The point datasets are deliberately reduced from 9962 points to 1314 points, but the feature of the geometry is carefully preserved by using AMIRA[22]. The reduction makes mesh so coarse that it is difficult to maintain quality

for simulation. However, it is well suited to test the perfor-mance of the current MGA procedure. The points at the edges and joints have to be marked and fixed to avoid sin-gularity during the smoothing. The results are shown in

Figs. 4 and 5 for triangular and quadrilateral meshes,

respectively. The improvement is significant, which can be seen from the quality indices measured in Tables 3

and 4. The mean quality index of triangular mesh is

enhanced from 0.8483 to 0.9046. It is well known that quadrilateral mesh is less stiff than triangular mesh and allows more degrees of freedom to move the mesh to change the shape of elements. As a result, it needs more

Fig. 6. Triangular surface mesh of foot: (a) original surface mesh and (b) surface mesh after smoothing.

smoothing treatments after generating the surface mesh. Our results show significant improvement in terms of worst and average property measured within elements. The mean and worst qualities are enhanced from 0.6982 to 0.8580 and 0.0001 to 0.4215, respectively.

Due to the advance of image processing, geometry reconstruction and mesh generation, nowadays larger and more complicated geometry models can be reconstructed from a sequence of bio-medical images data set. The need for such reconstruction is immense and it has become a common practice for bio-medical and bioengineering study. By scrutinizing the reconstruction procedure for such geometry models, it can be found that the surface geometry is usually defined by a set of unorganized points, or at least by a sequence of un-associated contours, that is

identified by pattern recognition methods on each image of interest. The surface triangulation is not straightforward. It is even more challenging for surface mesh quality enhance-ment in that the re-arrangeenhance-ment of the point distribution needs higher order interpolation methods for accuracy. Two test cases are given to demonstrate the capability of the proposed methods for tackling these challenging issues. The first example is the foot shape model of Polhemus

[40]. The original model was created by using FastSCAN with the FastRBF Extensions. The surface data points are measured by the laser scanner FastSCAN, and then the data points are reconstructed to form the geometrical model by the software FastRBF Extensions. It is a very convenient and efficient way to reconstruct a geometrical model, especially when the surface function is not

Fig. 8. Triangular surface mesh of rat: (a) original surface mesh; (b) surface mesh after smoothing; (c) original surface mesh (enlarged); and (d) surface mesh after smoothing (enlarged).

pre-defined, by using a laser scanner system. In this study, the original foot model contains 25,845 nodes. Similarly, the original data points are reduced to 4039 and it contains 8000 triangles on the surface. Our MGA approach con-structs a high order interpolation and iteratively optimizes the point distribution in terms of local evolution. The trian-gular mesh quality was enhanced as expected from 0.8275 to 0.9192 (Fig. 6). The original surface mesh obtained is actually well defined. However, the current approach can still give further enhancement of the mesh quality. For sys-tematical comparison, the quadrilateral mesh is generated and smoothed in a similar fashion. The results of obvious improvement are as expected that the surface quadrilater-als are more regular after the MGA approach (Fig. 7). This

also can be readily observed from the mean quality, which is improved from 0.7049 to 0.8632 (Table 4).

The last case is a rat shape that is reconstructed from a sequence of image slices of histological sections from Ryut-aro Himeno of RIKEN [41]. The original reconstructed model contains skin, bone, organs and so on, and it can help people to observe the whole model of rat virtually without dissect a real rat. This approach can be used to reconstruct more models of biology. In this study, we extract the skin model and reduce it by using AMIRA. As shown in Fig. 8, this model contains 25,670 nodes and 51,354 elements. The reconstructed shape is rough and irregular and the primary surface mesh contains many poor triangles elements. Some of the nodes connect to three

Fig. 9. Quadrilateral surface mesh of rat: (a) original surface mesh; (b) surface mesh after smoothing; (c) original surface mesh (enlarged); and (d) surface mesh after smoothing (enlarged).

elements and others connect to more than nine elements. The poor connectivity is treated first by the mesh structure modification operators [23–25]. To avoid interpolation error, the singular points, which locate at the boundary or contain large curvature variation, need be identified and fixed before the surface mesh smoothing. The MGA procedure is then used to re-construct surface interpolation function and to smooth the mesh by optimizing the point distribution. The result is shown inFig. 8b and the mean mesh quality improved from 0.8963 to 0.9453. The similar improvement for quadrilateral mesh can be found inFig. 9

andTable 4.

5. Conclusions

In this study, we present an evolutionary procedure, which integrates the micro genetic algorithm and a C1 sur-face function interpolation method, applied to sursur-face meshing on a set of unorganized points. Original surface meshes are generated based on Delaunay triangulation or reduced from the original models. Quadrilateral meshes are generated via the conversion of the triangular meshes and elementarily modified by the mesh structure modifica-tion operators. A C1interpolation function is constructed using the primary local elements to ensure the surface geo-metrical feature is preserved. The meshes are then smoothed using the micro-genetic algorithm on the recon-structed high order surface elements. The procedure is then tested in complicated and large-scale point datasets. The results show that the procedure successfully achieves better surface meshing with mesh quality significantly enhanced for different kinds of practical applications.

Generally, the generation of surface mesh is based on a pre-defined surface function, either in parametric patches or in algebraic form. In this study, the given data of model is just a set of unorganized points. In this situation, the major problem of surface mesh generation is of the mesh smoothing. As mention above, the given points are first tri-angulated to form a triangular mesh and/or converted them into quadrilateral mesh, and then surface mesh smoothing procedure is applied to improve the mesh qual-ity. A MGA mesh smoothing procedure is adopted in this study, which allows us to avoid the calculation of search direction and step size and to enable a global search for optimum. Furthermore, a C1surface function interpolation method is integrated into the MGA mesh smoothing proce-dure to ensure the geometrical accuracy of models during the surface mesh smoothing.

The extensions of the proposed method are to enhance the efficiency of MGA and the accuracy of the recon-structed surface functions. The MGA as one of the GA methods is by nature well suited to parallelization. There-fore, our further enhancement will be in performance issues of parallelism. The accuracy of the reconstructed local interpolation functions on surfaces determines the baseline of the accuracy of the entire model. Higher order

interpo-lation methods can be further developed to serve the purpose.

Acknowledgement

This work is supported by the National Center for High-performance Computing, Taiwan, ROC.

References

[1] H. Weimer, J. Warren, Fast approximating triangulation of large scattered datasets, Adv. Engrg. Software 30 (6) (1999) 389–400. [2] H.K. Zhao, S. Osher, R. Fedkiw, Fast surface reconstruction using

the level set method, in: Proceedings of the First IEEE Workshop on Variational and Level Set Methods, in conjunction with the Eighth international Conference on Computer Vision (ICCV), Vancouver, Canada, 2001, pp. 194–202.

[3] H. Park, K. Kim, Applications: an adaptive method for smooth surface approximation to scattered 3d points, Comp. Aided Des. 27 (1995) 929–939.

[4] H. Edelsbrunner, E.P. Mu¨cke, Three-dimensional alpha shapes, ACM Trans. Graphic (TOG) 13 (1) (1994) 43–72.

[5] D.A. Field, Laplacian smoothing and Delaunay triangulations, Commun. Appl. Numer. Meth. 4 (1988) 709–712.

[6] L. Freitag, On combining Laplacian and optimization-based mesh smoothing techniques, AMD Trends Unstruct. Mesh Gener. ASME 220 (1997) 37–43.

[7] P.M. Knupp, Achieving finite element mesh quality via optimization of the Jacobian matrix norm and associated quantities. Part I—a framework for surface mesh optimization, Int. J. Numer. Meth. Engrg. 48 (2000) 401–420.

[8] R.V. Garimella, M.J. Shashkov, P.M. Knupp, Triangular and quadrilateral surface mesh quality optimization using local parame-trization, Comput. Methods. Appl. Mech. Engrg. 193 (2004) 913–928. [9] S.A. Canann, J.R. Tristano, M.L. Staten, An approach to combined Laplacian and optimization-based smoothing for triangular, quadri-lateral, and quad-dominant meshes, in: Proceedings of the Seventh International Meshing Roundtable, 1998, pp. 479–494.

[10] J.S. Arora, Introduction to Optimum Design, McGraw-Hill Book Company, 1989.

[11] D.E. Goldberg, Genetic Algorithms in Search, Optimization and Machine Learning, Addison-Wesley Longman Publishing Co., Inc., Boston, MA, 1989.

[12] M. Gen, R. Cheng, Genetic Algorithms and Engineering Optimiza-tion, John Wiley & Sons, New York, 2000.

[13] K. Krishnakumar, Micro-genetic algorithms for stationary and non-stationary function optimization, Proc. SPIE Intell. Cont. Adapt. Syst. 1196 (1989) 289–296.

[14] W.A. Bedwani, O.M. Ismail, Genetic optimization of variable structure PID control systems, in: Proceedings of the ACS/IEEE International Conference on Computer Systems and Applications (AICCSA’01), Beirut, Lebanon, June, 2001, pp. 27–30.

[15] S. Chen, Y. Wu, S. McLaughlin, Genetic algorithm optimization for blind channel identification with higher order cumulant fitting, IEEE Trans. Evol. Comput. 1 (4) (1997) 259–265.

[16] S. Chakravarty, R. Mittra, L. Fellow, N.R. Williams, Application of a microgenetic algorithm (MGA) to the design of broad-band microwave absorbers using multiple frequency selective surface screens buried in dielectrics, IEEE Trans. Antennas Propag. 50 (3) (2002) 284–296..

[17] T.N.T. Goodman, H.B. Said, L.H.T. Chang, Local derivative estimation for scattered data interpolation, Appl. Math. Comput. 68 (1995) 41–50.

[18] T.N.T. Goodman, H.B. Said, A C1_{triangular interpolant suitable for}

scattered data interpolation, Comm. Appl. Numer. Meth. 17 (1991) 479–485.

[19] V.S. Gordon, D. Whitley, Serial and parallel genetic algorithms as function optimizers, in: Proceedings of the Fifth International Conference on Genetic Algorithms (IGGA5), 1993, pp. 177–183. [20] A. Okabe, B. boots, K. Suqihara, Spatial Tessellations: Concepts and

Applications of Voronoi Diagrams, John Wiley & Sons, New York, 1992.

[21] P.L. George, Automatic Mesh Generation, Application to Finite Element Methods, J. Wiley & Sons, New York, 1991.

[22] Template Graphics Inc. (TGS).

[23] S.A. Canann, S.N. Muthukrishnan, R.K. Phillips, Topological refinement procedures for triangular finite element meshes, Engrg. Comp. 12 (1996) 243–255.

[24] S.A. Canann, S.N. Muthukrishnan, R.K. Phillips, Topological improvement procedures for quadrilateral finite element meshes, Engrg. Comp. 14 (1998) 168–177.

[25] P. Kinney, CleanUp: improving quadrilateral finite element meshes, in: Proceedings of the Sixth International Meshing Roundtable, 1997, pp. 437–447.

[26] T.S. Lau, S.H. Lo, C.K. Lee, Generation of quadrilateral mesh over analytical curved surfaces, Finite Elem. Anal. Des. 27 (1997) 251–272. [27] B. Cheng, B.H.V. Topping, Improved adaptive quadrilateral mesh generation using fission elements, Adv. Engrg. Software 29 (7–9) (1998) 733–744.

[28] S. Maza, F. Noel, J.C. Leon, Generation of quadrilateral meshes on free-form surfaces, Comp. Struct. 71 (1999) 505–524.

[29] S.A. Canann, Y.C. Liu, A.V. Mobley, Automatic 3D surface meshing to address today’s industrial needs, Finite Elem. Anal. Des. 25 (1997) 185–198.

[30] B.P. Johnston, J.M. Sullivan Jr., A. Kwasnik, Automatic conversion of triangular finite element meshes to quadrilateral elements, Int. J. Numer. Meth. Eng. 31 (1991) 67–84.

[31] T. Whelan, A representation of a C2 interpolant over triangles, Comp. Aided Geom. Des. 3 (1986) 53–66.

[32] A. Sheffer, E. de Sturler, Surface parametrization for meshing by triangulation flattening, in: Proceedings of the Ninth International Meshing Roundtable, 2000, pp. 161–172.

[33] C. Gotsman, X. Gu, A. Sheffer, Fundamentals of spherical para-meterization for 3D meshes, ACM Trans. Graphic. (Proc. SIG-GRAPH 2003) 22 (3) (2003) 358–363.

[34] S. Lee, G. Wolberg, S.Y. Shin, Scattered data interpolation with multilevel B-splines, IEEE Trans. Vis. Comp. Graphic. 3 (3) (1997) 228–244.

[35] P. Costantini, C. Manni, A local shape-preserving interpolation scheme for scattered data, Comp. Aided Geom. Des. 16 (1999) 385– 405.

[36] M.J. Lai, Convex preserving scattered data interpolation using bivariate C1cubic splines, J. Comput. Appl. Math. 119 (2000) 249– 258.

[37] P.M. Knupp, Algebraic mesh quality metrics for unstructured initial meshes, Finite Elem. Anal. Des. 39 (2003) 217–241.

[38] G. Guy, G. Medioni, Inference of surfaces, 3D curves, and junctions from sparse, noisy, 3D data, IEEE Trans. Pattern Anal. Mach. Intell. 19 (11) (1997) 1265–1277.

[39] D.L. Sheu et al., National Cheng Kung University Islanders 101, Technical Report, Department of Aeronautics and Astronautics, National Cheng Kung University, Taiwan, ROC.

[40] Polhemus FastSCAN Cobra, <http://fastscan3d.com/download/sam-ples/>.