Marc Schoenauer

Leila Kallel

François Jouve

Fran Cois Jouve

Mechanical Inclusions Identification by Evolutionary Computation

Keywords: eme inverse, inclusions m echaniques, elasticit e lin eaire Genetic algorithms, stochastic optimization, inverse problem, mechanical inclusions, linear elasticity. 1

R ESUM E. L'identi cation d'inclusions est un probl eme th eoriquement mal pos e, pour lequel les m ethodes num eriques d eterministes se r ev elent instables. Les algorithmes evolutionaires sont des m ethodes d'optimisation stochastiques d'ordre 0 fournissant un cadre g en eral a la r esolution de probl emes inverses. Cependant, l'application des algorithmes evolutionaires au probl eme de l'identi cation d'inclusions doit faire face au probl eme crucial du choix de la repr esentation des solutions, c'est-a-dire de l'espace de recherche. Trois repr esentations sont propos ees. Les r esultats num eriques constituent une avanc ee importante dans le domaine de l'identi cation d'inclusions m ecaniques, et la robustesse des r esultats en pr esence de bruit est exp erimentalement etablie. Enn, les m erites respectifs des di erentes repr esentations sont discut es.

ABSTRACT. The problem of the identi cation of mechanical inclusion is theoretically ill-posed, and to-date numerical algorithms have demonstrated to be inaccurate and unstable. On the other hand, Evolutionary Algorithms provide a general approach to inverse problem solving. However, great care must be taken during the implementation: The choice of the representation, which determines the search space, is critical. Three representations are presented and discussed. Whereas the straightforward meshdependent representation su ers strong limitations, both mesh-independent representation provide outstanding results on simple instances of the identi cation problem, including experimental robustness in presence of noise. MOTS-CL ES : Algorithmes g en etiques, optimisation stochastique, probl

Introduction

In Structural Mechanics, the { non-destructive { identi cation of inclusions is a di cult problem, resisting to-date numerical methods: in its simplest instance, a structure is known to be made of two di erent materials of given mechanical characteristics, but their repartition in the structure is unknown. The available data consist of records of the mechanical behavior of the structure under known loadings. The goal is to nd the geometrical repartition of both materials from these experimental data. In steel manufacturing plants, for instance, it is of vital importance to check if coal scories are included in steel parts, and if their repartition does not dangerously weaken the whole part. For a given repartition of both materials, the computation of the simulated mechanical behavior of the structure is straightforward, using any Finite Element Analysis (FEM) package. The identi cation can then be viewed as an inverse problem.

Evolutionary Algorithms (EAs) are stochastic optimization methods that have been demonstrated useful to solve di cult, yet unsolved optimization problems. The most well-known algorithms referring to the Evolutionary Computation paradigm are the Genetic Algorithms (GAs), but other variations on the \survival of the ttest" theme also gave birth to powerful optimization algorithms. Requiring no regularity of the objective function (or on the constraints), EAs are able to tackle optimization problems on di erent kinds of search spaces, such as continuous, discrete or mixed spaces, as well as spaces of graphs or lists. The only prerequisite is the de nition of evolution operators such as crossover and mutation, satisfying as much as possible heuristically derived requirements. The two main drawbacks of EAs are rst the large number of evaluations of the objective function they usually require before eventually reaching a good, if not optimal, solution; and second, their stochastic aspect, weakening their robustness. Hence, EAs should be used with care, on problems beyond the reach of standard deterministic optimization methods.

This paper addresses the inverse problem of Mechanical Inclusions Identication using EAs. A possible objective function for such inverse problems is the di erence between the simulated mechanical behavior of a tentative repartition of both materials, and the actual (experimental) behavior of the real structure. However, the main di culty is to de ne the search space in which the EA will search. Considering past works on the Optimum Design problem, (a closely related problem from the Evolutionary Computation point of view, where the goal is to nd a partition of a design domain into material and void), a straightforward representation is de ned from a xed mesh of the structure, leading to a xed-length bitstring representation well-suited to standard Genetic Algorithms. This \bitarray" approach allowed signi cant breakthroughs in Topological Optimum Design 34,[START_REF] Kane | Genetic operators for two-dimensional shape optimization[END_REF]. However, this approach su ers from strong limitations due to the dependency of its complexity on the underlying mesh: the optimization problem rapidly becomes intractable when this mesh is Figure 1: A sample inclusion problem re ned. To cope with that di culty, two non-standard representations (termed the Vorono representation and the H-representation) were introduced 48, 49], independent of any a priori discretization, but leading to variable-length \individuals" for which speci c operators were designed and implemented.

This paper is concerned in applying Evolutionary Algorithms working on these representations to the inclusion identi cation problem in the context of linear elasticity, and is organized as follows:

The mechanical problem of inclusion identi cation is presented in details in Section 2. Section 3 gives a brief overview of the state-of-the-art in Evolutionary Algorithms (EAs), introducing the di erent components of Evolutionary Computation. In Section 4, a general framework to address inverse problems with EAs is introduced: The only prerequisite is a good (i.e. accurate and robust) numerical model of the direct problem, : : : and, of course, a good representation of the target search space. Three representations for that problem are presented in Section 5, together with their speci c evolution operators (crossover and mutations). The rst results, using the Vorono representation, are presented in Section 6, demonstrating outstanding performance on arti cial instances of the inclusion identi cation problem. Finally, comparative results of all three representations are discussed in Section 6.3, rst reinforcing the a priori arguments against the bitstring representation, then giving hints on how to a priori choose between the Vorono representation and the H-representation on a given instance of a problem. Further directions of research are sketched in the conclusive Section 7.

The Mechanical Problem

Background

Consider an open bounded domain IR N (N = 2; 3), with a su ciently smooth boundary @ , lled with a linear elastic material (as in Figure 1).

Under the hypothesis of small deformations (linear elasticity) the following equations hold:

"(x) := 1 2 (ru(x)ru T (x)) the strain tensor;

(1) and (x) := A(x)"(x) the stress tensor;

where u(x) is the displacement eld at point x and A(x) is the elasticity tensor (a fourth order tensor) involved in the Hooke's law (2): A is supposed to be inhomogeneous, which means that its value depends on the point x, and A(x) is a positive de nite tensor which satis es some symmetry conditions.

When A is given, one can state two kinds of boundary value problems, respectively of Dirichlet and Neumann type: div = 0 in ; u = u 0 on @ ;

(3) and div = 0 in ; :n = g 0 on @ ; [START_REF] Beyer | Toward a theory of evolution strategies: The ( ; )-theory[END_REF] where u 0 and g 0 are respectively a given displacement eld and a given external force eld on the boundary @ . It is well known (see e.g. 13] that each of these problem has a unique solution (for the Neumann's problem, one has to impose an integral condition on g 0 to ensure existence, and an integral condition on u to eliminate rigid displacements).

In the following, the inverse problem will be considered: nd A such that ; 8i 2 f1; :::; ng; 9u; 8 > > < > > : div = 0 in ; u = u i on @ ; :n = g i on @ ; = A";

(

where (g i ) i=1;::;n and (u i ) i=1;::;n are given. Problem ( 5) is a discrete version of the \ideal" inverse problem:

nd A, given the application: A : u j@ ! j@ .

(

) 6
The underlying physical problem is still lacking much more of known data than problem (5) since A is only known through a nite number of experimental measurements performed at a nite number of points. Hence, the \real" identi cation problem treated by mechanical engineers can be stated as: nd A such that ; 8i 2 f1; :::; ng; 9u; 8 > > < > > : div = 0 in ; u(x j ) = u j i 8j 2 f1; :::; pg; (x j ):n = g j i 8j 2 f1; :::; pg; = A"; [START_REF] Beyer | On the asymptotic behavior of multi-recombinant evolution strategies[END_REF] where information on the boundary is only known at a nite number of experimental points (x j ) j=1;:::;p for a nite number (n) of experiments. In addition, these data may be known with a certain amount of experimental error or noise.

State of the art

If the aim is the numerical treatment of problem (6) (and a fortiori problem (5) or ( 7)) by \classical" (i.e. non-stochastic) methods, two theoretical points are crucial: existence and uniqueness of A as a function of A , continuity of the dependency of A with respect to A .

Existence is of course essential to the pertinence of the identi cation problem, but uniqueness and continuity are only needed to ensure the reliability and the stability of deterministic numerical algorithms. On the other hand, EAs can deal with non-continuous functionals and non-unique solutions.

Problem [START_REF] Beyer | Toward a theory of evolution strategies: Self-adaptation[END_REF] is the elastic equivalent of the so-called tomography problem where the elliptic operator is the conductivity operator (div(Aru), u scalar eld) instead of the elasticity one (div(A"(u)), u vector eld).

The tomography problem has been widely studied. Under some hypothesis, existence and uniqueness have been proved. However, the continuity of the functional is only known in a weak sense, that cannot help numerical simulations.

The elasticity problem ( 6) is more di cult. Existence and uniqueness have been proved for isotropic Hooke's laws, but there is no continuity result (see 8, 9], or 15] for a comprehensive bibliographical discussion on this subject).

Numerical simulations by classical methods have shown that both tomography 37] and elastic identi cation problems 15] are ill-posed, and thus EAs are good tentative choice for a robust numerical method.

The simpli ed problem

In this preliminary paper, attention has been focused on representation and on speci c operators for EAs. To highlight speci c problems involved in these algorithms, the mechanical problem [START_REF] Beyer | On the asymptotic behavior of multi-recombinant evolution strategies[END_REF] was restricted to a two-dimensional simpler class of problems:

Let A 1 and A 2 be two isotropic elasticity tensors, fully de ned by Young's moduli E 1 and E 2 and Poisson ratios 1 and 2 . The aim is to solve Problems (5) and [START_REF] Beyer | On the asymptotic behavior of multi-recombinant evolution strategies[END_REF], restricting allowable values of A(x) to

A(x) = A 1 if (x) = 0 A 2 if (x) = 1 (8)
where is a characteristic function de ned on . These problems, although less general than ( 5) and ( 7) are still beyond the capabilities of deterministic algorithms (see 15]).

However, Problem (7) can be treated in the same way, as well as identi cation in non-linear elasticity (this is part of on-going work, see Section 7).

The direct elasticity problem is solved by a classical nite element method. All simulations (see Section 6) have been performed on a square domain, held xed on its left side, discretized by a regular mesh of size 24 24 elements. All details and mechanical constants are speci ed in Section 6.

Evolutionary Algorithms

This section presents a brief overview of Evolutionary Algorithms, both giving a modern generic and pragmatic point of view, and presenting a historical perspective on Evolutionary Computation.

The Generic Evolutionary Algorithm

Evolutionary Algorithms are stochastic optimization algorithms based on a crude mimic of the Darwinian principle of the Survival of the ttest.

Let the search space be a metric space E, and let F be a function from E on IR called the tness function. The problem is to nd the optimum of F on E.

A population of size P 2 IN is a set of P individuals (points of E) not necessarily distinct. This population is generally initialized randomly, uniformly on E. The tnesses of all individuals are computed. The population then undergoes a succession of generations, that can be described the following way:

Selection: Some P individuals are selected to give birth to o spring. Numerous selection processes can be used, either deterministic or stochastic. However, all implement in one way or another the Darwinian principle of the Survival of the ttest: the selection step is based on tness comparisons among individuals. Depending on the selection scheme used, some individuals can be selected more than once. At that point, the selected individuals give birth to copies of themselves (clones).

Application of evolution operators:

To each one of these copies is applied some evolution operator, giving birth to one or more o spring. The choice among possible operators is stochastic, according to user-supplied probabilities. These operators are always stochastic operators, and one usually distinguish between crossover (or recombination) operators, when two (or more) parents exchange part of their genetic material to build up one o spring, and mutation operators, in which one parent is randomly modi ed.

Evaluation: Stopping criterion: The process stops when: A target value for the best tness { when it is known { is reached, or some a priori maximal number of generations is run, or the best tness value in the population does not improve during a user-de ned number of generations, which ever comes rst. If the stopping criterion is not ful lled, another generation starts. Ideally, the evolution operators should be de ned on the same space than the tness function, called phenotype space, or behavioral space. Unfortunately, in most cases, one must introduce an intermediate space, called genotype space, or representation space. The mappings from the phenotype space in the genotype space is termed coding. The inverse mapping from the genotype space in the phenotype space is termed decoding. Genotypes undergo evolution operators, and their tness is evaluated on the corresponding phenotype. The properties of the coding mappings can greatly modify the dynamical behavior of the evolutionary algorithm.

Historical EAs

The most widely known Evolutionary Algorithms probably are Genetic Algorithms (GAs). First introduced by Holland 27], then popularized by Goldberg 24], historical GAs are based on a binary representation of the individuals (then often termed chromosomes) and heavily rely on the cross-over operator (i.e. stochastic operators generating one o spring from two Darwinian-selected parents). Parents are selected based on their tness, the operator considered crucial is the crossover operator, and the replacement step simply replaces all parents by all o spring.

On the other hand, historical Evolution Strategies 50, 44] and Evolutionary Programming 20] try to use the \natural" search space for the problem at hand (e.g. IR n when it comes to real parameter optimization). All parents are selected randomly, the main operator is considered to be the mutation operator, and the replacement step is a deterministic choice (based on tness values) among o spring and possibly parents in Evolution Strategies, and a stochastic tournament between parents and o spring in Evolutionary Programming. Moreover, for real-valued parameter optimization, all these algorithms use what is called Gaussian mutation:

x mutation ! x + N(0; ); [START_REF] Bui | An introduction to inverse problems in Structural Mechanics[END_REF] where N(0; ) is a normal random variable with user-de ned deviation > 0.

Nevertheless, modern Evolutionary Algorithms are de ned on a pragmatic basis: be they real-encoded GAs 40] or generalized Evolution Strategies using recombination 50]. Moreover, any possible evolution scheme (i.e. choice of selection-replacement methods) can be used, \as long as it works". More details about the recent advances in Evolutionary Computation can be found in 1, 18]. These references also describe examples of successes of Evolution Algorithms on di cult global optimization problems, on which standard classical methods (e.g. gradient-based methods for numerical optimization) fail.

In any case, the price to pay for the wide application area is a fairly high computational cost. Moreover, many di erent parameters are involved (choice of the evolution scheme, population size, probability of application of di erent operators, stopping criterion). These parameters can have dreadful consequences on the overall convergence of the algorithm, they highly depend on the problem at hand, and no attempt to automatically adjust them has proved to be general enough [START_REF] Goldberg | Genetic algorithms in search, optimization and machine learning[END_REF][START_REF] Grefenstette | Optimization of control parameters for genetic algorithms[END_REF] but the usual trial-and-error systematic testing 47].

Theoretical results have been derived regarding the convergence of some speci c instances of EAs. The strongest of such results deal with the convergence of ESs on convex functions. Not only do the algorithms converge to the global optimum with probability 1 50], but also Beyer obtained very precise quantitative convergence rates for various schemes in the ES framework 2, 3, 4, 5, 6]. On the GA side, a nice global convergence result in probability has been proved by Cerf 10,[START_REF] Cerf | Une th eorie assymptotique des algorithmes g en etiques[END_REF] on a modi ed bitstring GA, while di erent partial results { including a non-convergence result on the non-elitist GA { have been shown by Rudolph 46,[START_REF] Rechenberg | Evolutionstrategie: Optimierung Technisher Systeme nach Prinzipien des Biologischen Evolution[END_REF]. However, these theoretical results are still of poor use when it comes to practical implementation (e.g. no convex function needs EA to be optimized, and the minimum population sizes obtained by Cerf are too large to be actually used). Hence, the user can only refer to past experimentations described in the literature.

Implementation hints

As said above, one strength of EAs is their ability to work on weird search spaces (i.e. with non-standard representations), as long as evolution operators (crossover and mutation) are provided. Moreover, though no theoretical results exist to show the way to design a successful evolutionary algorithm, some heuristically derived basic principles have been stated: the variance of the tness values of all individuals sharing some genetic materials (the schemata in GA terminology) should decrease when the amount of common material increases 43]; crossing-over two individual should respect their common genetic material 42]; mutation should be ergodic 42], i.e. a nite number of mutations should be able to join any two points of the search space; this point is also crucial for all theoretical convergence results based on Markov chain analysis 11]; mutation should respect the principle of strong causality 44], i.e. small changes of the genotype should provide small variations of the tness1 . Other important features to be considered are the the degeneracy of the representation (how many di erent genotypes map onto the same phenotype) and its redundancy (genetic material that could be omitted without modifying the phenotype). Whereas redundancy can be bene cial (see for instance the metaphor of natural introns 38]), it implies some degeneracy, which is usually considered harmful to evolutionary optimization. Finally, the degree of epistasis of the representation (the interaction between di erent part of the genotype to build the phenotype) must also be taken into account. Section 7 will consider these questions for the speci c representations used for the inclusion identi cation problem, in the light of the experimental results of the evolutionary optimization method described in the rest of the paper.

Evolutionary Inclusion Identi cation

This section presents the general numerical framework that will be used to solve the mechanical problem introduced in Section 2, using Evolutionary Algorithms.

An inverse problem

The problem to solve is the inverse problem corresponding to the direct problem described by Equation 8 of Section 2. Having a good simulation of a direct problem through Finite Elements Analysis, Evolutionary Computation can be used to address the inverse problem in a straightforward manner, as described schematically in Figure 2: the results of the FE Analyses performed Experimental conditions FEM on the structure at hand are compared to original experimental results, and the goal is to reach the smallest possible error.

Experimental results

Current structure

Numerical results

Fitness of current structure

The experimental data will consist in given displacements of the structure under given loading cases.

Experimental settings

All numerical results for the problem at hand (Section 2) will be obtained on the two-dimensional square structure presented in Figure 3: the structure is held xed on its left-side, and normal point-wise forces are applied at given points of the three other boundaries.

The aim is to identify the repartition of two materials into the square: a hard material (E = 1:5 and = 0:3) and a soft material (E = 1 and = 0:3). A xed mesh of size 24 24 will be used used throughout the experiments, unless otherwise mentioned.

The reference experimental loading cases used to compute tness values were actually computed from a known con guration of materials inside the structure. The optimal solution is thus known, which allows a better insight and understanding during the evolution of the population. Moreover, much exibility was required during the tuning of the overall process, that actual experimental results could not have bought. Finally, considerations about the noise in the experimental reference recording also favor simulated results: there still is a bias (due to numerical error in the Finite Element Analysis), but, this bias is hopefully the same for the experimental results and during the tness computation. Hence, it should not weaken the results of the evolutionary algorithm, as could unpredictable noise in actual measures. Of course, further work will have to consider results of actual experiments to fully validate the approach.

Figure 4 shows the four reference structures that will be used throughout this paper: the fairly easy corner problem and the checkerboard problem, dicult from the mechanical engineers point of view, are taken from 15], while the asymmetric and symmetric problems are purposely designed to test di erent behaviors of the representations introduced further in this section.

The tness functions

In real world situations, the design of the tness function should take into account all available loading cases experimented on the real structure, using as much information as possible about the behavior of the structure to be identi ed. However, the choice made in these preliminary experiments of using \simulated experiments" makes it possible to use as many loading cases as needed. In an attempt to have a sampling of the mechanical behavior of the structure as uniform as possible over the domain, 37 di erent loading cases were used, each loading case consisting of applying a given normal force at one point of the boundary of the structure (as shown on Figure 3).

A very important issue in inverse problem solving is the generalization capability of the solution: How good is the result of the optimization process when used under experimental conditions that are di erent from the ones used during the optimization itself? The usual answer is to use, during the identication, many di erent experimental conditions, also termed tness cases. The tness is then the average of the error over all tness cases. Needless to say, the total computational cost increases with the number of tness cases.

However, the use of simulated \experimental" values allows to test the overall algorithm with di erent points of view: In real world situations, some gauges are placed to the boundary of the structure, and only the displacements at those points are available. Figure 3 is an example of such a situation, where 9 gauges are represented. However, it seems clear that the more measure points, the easier is the identi cation task. Hence, three di erent tness functions have been used throughout the numerical experiments presented in Section 6, using the 37 loading cases described above, with di erent numbers of measure points:

The most informative tness function, referred hereafter to as the total tness, takes into account the displacements at all nodes of the mesh. An intermediate tness function uses only the displacements of all nodes lying at the boundary fo the structure, and is termed the boundary tness.

The real-world tness uses only the 9 measure points represented on Figure 3, incorporating much less information to the tness function.

The Evolutionary Algorithm

The Evolutionary Algorithm used for all experiments presented in this paper is a standard GA scheme: rank-proportional selection (only the rank in the population of all individuals is taken into account in the selection process, regardless of the actual tness values), crossover rate of 0.6, mutation rate per individual of 0.2, all o spring replace all parents. The population size is set to 100, and at most 300 generations of the algorithms are allowed { but it stops whenever 50 generations are run without any improvement of the overall best tness. Hence, between 10000 and 25000 Finite Element Analyses (each involving 37 loading cases, but with the same global matrix) were performed for each run, requiring around 2 hours of CPU time on a high range HP workstation (HP-PA 8000) for the 24 24 mesh used in most cases.

However, the stochastic nature of the algorithm imposes that all experiments are run many times with di erent independent initial populations: as clearly stated in 30], \You should never draw any conclusion of a single run of any Evolutionary Algorithm".

Hence, for each setting of experimental conditions (problem and representation), 21 independent runs have been performed, and the forthcoming results of Section 6.3 presenting the online performances ( tness along generations) or the o -line performances ( tness reached in the end of evolution) of these runs are median results: At each generation, the best individuals of all runs are sorted according to tness, and the median value (here, the 11 th out of 21) is shown. Note that EA results are generally presented through their averages. But, as discussed in 21], whereas using the average values (together with standard deviation and/or T-test results) make the implicit assumption of normal repartitions for the sampled values, choosing the median values seems more adapted to non-parametric statistics. Moreover, when dealing with tnesses that can take very small values (as the mean square errors of the displacements), averages hardy re ect actual behavior of the algorithm: Error values for di erent runs can be of fairly di erent orders of magnitude. Hence, a few large values might hide a majority of very small errors.

The genotype space

The most critical step in evolutionary algorithms is the choice of the representation, which de nes the search space. In any case, the algorithm will nd a good { if not the best { solution in the search space. Hence it seems that larger search spaces allow better solutions to be found. However, the larger is the search space, the more di cult is the optimization task, and a trade-o has to be found. Next section will present di erent tentatives to address this issue.

Representations for mechanical inclusions

This section introduces three representations for the problem described in Section 4.2. The rst one is the obvious mesh-based bitstring representation. The next two ones are variable-length representations, and use real numbers as main components. Hence, speci c operators (e.g. crossover and mutation) have to be designed for each of these representation.

Prerequisites

A solution to the inclusion identi cation problem is a partition of the domain of the structure into two subsets, each subset representing one of the materials involved. Moreover, all connected components of any subset should have a non-void interior and a regular boundary.

A theoretical framework has been developed by Ghaddar & al. 23] in the context of Structural Optimum Design: The search space is restricted to partitions with polygonal boundaries. Theoretical results are proven, approximation spaces are introduced and corresponding approximation results are obtained. Though the objective function considered in this paper is quite di erent from the one in 23], the same search space will be used here.

However, a signi cant di erence between the objective functions in 23] and the one to be used here is that the inclusion identi cation problem requires a Finite Element Analysis on the direct problem to compute the tness of a point of the search space (i.e. a given repartition of both materials), as detailed in Section 4.3. It is well-known that meshing is a source of numerical errors 14]. Hence, for any Evolutionary Algorithm, using a tness function based on the outputs of two Finite Element Analyses performed on di erent meshes is bound to failure, at least when the actual di erences of behavior will become smaller than the unavoidable numerical noise due to remeshing. The use of the same mesh for all Finite Element Analyses is thus mandatory in order to obtain signi cant results, at least for structures whose mechanical behaviors are to be compared, i.e. inside the same generation.

The bitarray representation

Once the decision to use a xed mesh has been taken, and with even very little knowledge of EAs, a straightforward representation for a partition of the given domain uses bitstrings. Each element of the xed mesh belongs to either one of the subsets of the partition, which can be symbolically labeled 0 or 1. The resulting representation can be viewed as a bitstring or, more precisely, as a bitarray. Indeed, considering such array of bits as a one-dimensional bitstring is bound to failure, and speci c two-dimensional crossover operators have to be designed (see [START_REF] Kane | Structural topology optimization in linear and nonlinear elasticity using genetic algorithms[END_REF]). Nevertheless, almost all previous works using Genetic Algorithms on the related Optimum Design problem did use that representation [START_REF] Holland | Adaptation in natural and arti cial systems[END_REF][START_REF] Cerf | An asymptotic theory of genetic algorithms[END_REF][START_REF] Kallel | Fitness distance correlation for variable length representations[END_REF][START_REF] Kane | Genetic operators for two-dimensional shape optimization[END_REF].

However, the limits of this bitstring representation clearly appear when it comes to re ne the mesh (in order to reach better precision of the resulting shape) or if one wants to apply evolutionary techniques to solve 3-dimensional problems: this would imply a huge bitstring, as the size of the bitstring is that of the underlying mesh. In order to reach convergence with large bitstrings, larger populations are required, increasing linearly with the size of the bitstring, according to both theoretical results of Cerf 11] and empirical studies 51].

Moreover, more generations are also needed to reach convergence, and the resulting algorithm rapidly becomes intractable. Some numerical experiments of section 6 con rm this phenomenon.

These considerations show the need for other representations, not relying on a given mesh { even if a xed mesh is used thereafter to compute the tness function. Two of such representations have been designed, and successfully used on the Optimum Design problem 48, 49].

The Vorono representation

A possible way of representing partitions of a given domain comes from computational geometry, more precisely from the Vorono diagram theory. The ideas of Vorono diagrams are already well-known in the FEM community, as a powerful tool to generate good meshes 22]. However, the representation of partitions by Vorono diagrams to describe their evolutionary optimization seems to be original.

Vorono diagrams: Consider a nite number of points V 0 ; : : : ; V N (the Vorono sites) of a given subset of IR n (the design domain). To each site V i is associated the set of all points of the design domain for which the closest Vorono site is V i , termed Vorono cell. The Vorono diagram is the partition of the design domain de ned by the Vorono cells. Each cell is a polyhedral subset of the design domain, and any partition of a domain of IR n into polyhedral subsets is the Vorono diagram of at least one set of Vorono sites (see [START_REF] Michalewicz | Genetic Algorithms+Data Structures=Evolution Programs[END_REF][START_REF] Beyer | On the asymptotic behavior of multi-recombinant evolution strategies[END_REF] for a detailed introduction to Vorono diagrams, and a general presentation of algorithmic geometry).

The genotype: Consider now a (variable length) list of Vorono sites, each site being labeled 0 or 1. The corresponding Vorono diagram represents a partition of the design domain into two subsets, if each Vorono cell is labeled as the associated site (here the Vorono diagram is supposed regular, i.e. to each cell corresponds exactly one site). Example of Vorono representations can be seen in Figure 5. The Vorono sites are the dots in the center of the cells. Note that this representation does not depend in any way on the mesh that will be used to compute the mechanical behavior of the structure. Furthermore, Vorono diagrams being de ned in any dimension, the extension of this representation to IR 3 and IR n is straightforward.

An important remark is that this representation presents a high degree of epistasis (the in uence of one site on the physical shape is modulated by all neighbor sites). This will be discussed in more details in Section 7.

Decoding: Practically, and for the reasons stated in Section 5.1 the tness of all structures will be evaluated using the same xed mesh. A partition described by Vorono sites is thus mapped on this xed mesh: the subset an The crossover operators exchange Vorono sites on the basis of geometricallybased choice. In this respect it is similar to the speci c bitarray crossover described in 35]; moreover, this mechanism easily extends to any dimension 31]. Figure 5 demonstrates an application of this crossover operator. a rst mutation operator performs a Gaussian mutation on the coordinates of the sites, as in Evolution Strategies (see Equation 9) or randomly ips the boolean attribute of some sites. In all cases, the strength of the mutation (the standard deviation of the Gaussian mutation, or the number of sites whose labels are ipped) is adjusted proportionally to the relative tness of the individual at hand, as in Evolutionary Programming 17]: t individuals undergo weak mutations while un t individuals are more deeply modi ed.

Figure 6: The H-representation crossover operator. A random line is drawn across both structures, and the holes on one side are exchanged.

"standard" mutation operators for variable-length representations must include addition and destruction of some sites on the list.

H-representation

Another representation for partitions is based on an old-time heuristic method in Topological Optimum Design (TOD): from the initial design domain, considered as plain material, one removes material at locations where the mechanical stress is minimal, until the constraints are violated. However, the lack of backtracking makes this method useless in most TOD problems. Nevertheless, this idea gave birth to the \holes" representation 16], later termed H-representation.

The representation: The design domain is by default made of one material, and a (variable length) list of \holes" describes the repartition of the other material. These holes are elementary shapes taken from a library of possible simple shapes. Only homogeneous representations involving either rectangular or triangular holes are considered at the moment. Moreover, two di erent representations can be designed from each elementary form, depending on which material the \holes" actually stand for. Results of Section 6.3 will show that this choice does in uence the results. Three representations will hence be considered in the rest of the paper: the 1-rectangle representation, in which holes are rectangles representing material 1, the 0-rectangle representation where the rectangular holes stand for material 0, and the 1-triangle representation (with self-explained meaning).

Example of structures described in the 1-rectangle representation are presented in Figure 6. The rectangles are taken in a domain larger than the design domain, in order not to bias the boundary parts of the design domain toward the default value.

All holes representations, as the Vorono representation, is independent from any mesh, and hence its complexity does not depend on any required accuracy for the simulation of the mechanical behavior of the structure. Its merits and limitations will be discussed in the light of the experimental results presented in section 6.

Decoding: As for the Vorono representation, the simulated behavior of the shapes is computed on a given xed mesh, to limit the numerical noise due to re-meshing. The criterion to decide which subset an element does belong to, is based on whether its center of gravity belongs to a hole (in which case the whole element is considered made of the hole-material) or not.

Evolution operators: The evolution operators are quite similar to those of the Vorono representation: crossover by geometrical (2D or 3D) exchange of holes (see Figure 6 for an example); mutation by Gaussian modi cation of the characteristics (coordinates of the center, width and length) of some holes, as de ned by Equation 9; mutation by addition or destruction of some holes;

Initialization procedure: whereas the initialization procedure for the Vorono representation was straightforward (uniform choice of the number of Vorono sites between 1 and a user-supplied maximum number, and uniform choice of the sites in the structure, and of their label), that of the holes representation raises some interesting issues (discussed in details in 32]), as no obvious procedure exists that would be \uniform" over the space of holes. Based on results involving Fitness Distance Correlation analyses [START_REF] Jensen | Topological Structural Design using Genetic Algorithms[END_REF][START_REF] Kahng | Toward more powerful recombinations[END_REF] that are out of the scope of this paper, the following were used, assuming a rectangular structure = 0; L x ] 0; L y ] : In all cases, the number of holes is uniformly chosen in 1; N Max ] for some user-supplied parameter N Max (50 in all experiments presented here).

For the rectangle representations, uniform random choice of the center of the rectangle in , and uniform random choice of the width and height in 0; Lx 4 ] and 0; Ly 4 ] respectively. For the triangle representations, uniform random choice of a rst point in , then determination of the two other points by gaussian perturbations of mean deviation x Max =2 of the coordinates of the rst point.

Numerical Results

This section presents some numerical experiments on the inclusion identication problem using Evolutionary Algorithms. The rst results address the feasibility issue. In subsection 6.2, EAs demonstrate their robustness in the presence of noise in the experimental data. Finally, comparative experiments of subsection 6.3 bring a partial answer to the question of the a priori choice of a representation for a given instance of problem. In particular, the bitarray representation de nitely demonstrates its poor performances when compared to the mesh-independent representations (Vorono , 1-and 0-rectangles).

Feasibility Results

The very rst experiments were performed using the Vorono representation, described in Section 5.3. They were obtained on the reference partitions represented in Figure 4 of section 2.

The rst results on the corner problem were astonishingly good: the exact solution was found 3 times out of 10 independent runs, in less than 100 generations when using the \total" tness, and even twice (in about 250 generations) using the boundary tness. Figure 7 shows an example of a successful run, where the best individual at di erent stages of evolution is plot together with the error. The Vorono sites are represented on the gure by grey little dots.

On the symmetrical and asymmetric problems of Figure 4, some di erence begin to appear between the di erent tnesses, but these examples are still easy problems. And when it came to the checker problem, the total tness gave much better results than the real-world tness, as can be seen in Figure 8. However, the real-world tness gave some interesting results, as can be seen on Figure 8-b: the actual values are clearly identi ed along the free boundary, and badly along the xed part of the boundary, where very little information is available.

Results on noisy data

After these rst satisfactory results on exact data, further validation of the proposed approach had to take into account possible errors and noise in the data, before using actual { noisy { experimental data. In order to test the sensitivity of the algorithm to noise, arti cial noise was purposely introduced in the reference \experimental" displacements.

Figure 9 shows the results obtained on the (easy) \corner" example, when the relative amount of noise is set to 2% and 5%: all reference displacements were multiplied by a term (1 + "), where " was a random variable uniformly distributed in 0:02; 0:02] and 0:05; 0:05] respectively. The results are (of course !) degraded, but they demonstrate a fairly good robustness, at least on this example.

Comparative results

This section compares the four representations de ned in section 5, namely the bitarray, the Vorono , the 1-and 0-rectangles and the triangle representations. As discussed in Section 4.4, all results are median over 21 independent runs.

Figure 10 presents the on-line results, i.e. the plots of error on the displacements along generations, while Table 1 presents the o -line results after 300 generations. Note that the 0-rectangle was not tested on the two symmetrical problems (the symmetric and checker-board problems): though the problems are not mechanically symmetric (the hard and soft materials are exchanged), the asymmetry is hardly distinguishable from the EA point of view, and the results of the 0-rectangle representation are broadly of the same nature than those of the 1-rectangle on these problems. Some clear tendencies emerge from these results:

The bitarray representation is clearly outperformed by all other representations, except on the checker-board problem, where the di erence is not signi cant. Moreover, the best run of the bitarray trials was evolved up to 1000 generations, with resulting errors of respectively 0:014, 0:0004, 0:0003 and 0:0005 for the 4 problems, and were worse than the results of most other representations in 300 generations. In particular they were never better than the Vorono representation. The Vorono representation performs globally better than all Holes rep-(a) (b) The results on the asymmetric problem demonstrates : : : the asymmetry of the holes representations: the 0-rectangle is widely outperformed by the 1-rectangle on that { purposely de ned { problem. The shape of the target subset to be matched by the holes does in uence the relative results of the Holes-representations: the triangle representation outperforms the rectangle representation on the corner problem, while the reverse is clearly true on the checker-board problem.

In uence of the mesh

This section presents preliminary investigations to a posteriori justify the a priori arguments against the bitarray representation related to its dependency on the underlying mesh (see Section 5.2). Three di erent meshes are considered: the 24 24 mesh used in the preceding sections, together with a two times coarser 12 12 mesh and a re ned 48 48 mesh are used. The same 37 loading cases are considered in the tness (see section 4.3), and the reference displacements are computed on the same mesh than the following runs.

Tables 2 and3 give the respective o -line results of the bitarray, the Vorono and the 1-rectangle representations, for the asymmetric and the checkerboard problems of Figure 4. Further, Figures 11 and12 show the result obtained by the bitarray-, the Vorono -and the rectangle-representation on the asymmetric and the checker-board problems for one run of 300 generations on the 48 48 mesh: whereas both variable-length representations obtain a solution resembling the actual structure, the bitarray representation ends up in a Asymmetric Checker best median best median Voronoi 0 0.00012 0.00013 0.00026 rectangle 0 0 0.00032 0.0017 bitarray 0.003 0.0054 0.00055 0.00074 totally irrelevant result. Although, for technical reasons (a single run for the 48 48 mesh take about 10 hours of CPU time on a powerful HP-PA8000 processor), the results for the 48 48 mesh are still incomplete, it is obvious that the performance of the bitarray representation does indeed decrease rapidly when the complexity of the mesh increases.

Discussion and further work

Due to the monopolistic use of the bitarray representation in all works of the EA literature dealing with similar problems, it is necessary to emphasize here its inadequation for problems involving the search for a partition of a domain in IR n . A priori considerations on the complexity of the representation for ne (or 3D) meshes, \Fitness Distance Correlation"-based heuristical results 32] and preliminary experimental results on the Topological Optimum Design problem 48, 49] gave a broad picture that is now clearly con rmed by the experimental results of Section 6. Note that the inclusion identi cation problems treated in this paper allow a ner control of the numerical experiments than the problems of Topological Optimum Design (e.g. the global optimum is known, the problem is unconstrained, : : :).

The validation results in the presence of noise of Section 6.2 are of course only a rst step toward using results of actual mechanical experiments as reference data. An intermediate step, part of on-going work, is to consider for reference data the results of simulations on a very ne mesh, and to test the reliability of the approach with respect to the accuracy of the numerical approximation.

Regarding the a priori choice of a representation for a given problem, the speci c characteristics of the problem at hand should be taken into account: if there is strong evidence that one material is present only in small quantities, then one hole representation should be chosen, with the holes representing the minority material. Further, if the general shape of the small inclusions can be guessed, choosing elementary shapes for holes that easily ll up the supposed inclusion shapes can help the identi cation process. Though this is generally beyond what can be asked to mechanical experts, detailed Fitness Distance Correlation analyses 32] can be of some help here again.

On the other hand, if both materials are approximately present in the same proportions, the Vorono representation clearly is to be preferred. Moreover, the Vorono representation seems to be more robust to wrong guesses of actual proportions, and hence should probably be used in any doubtful case.

Another important direction of research is to address the more general problem (6), i.e. the case where the mechanical characteristics of the materials are unknown. The proposed representations can { and will { be extended easily to that framework: the label of Vorono sites can easily be replaced by a set of real values (e.g. the E and coe cient of the material); in the same line, realvalued coe cients can be attached to each hole in the Hole representation. The decoded values for a point in the structure will then be the average of all values that cover it. First experiments on the problem of identi cation of spatial geological models in petroleum prospection 39] gave encouraging results.

Summary and Conclusion

This paper has presented the rst results using Evolutionary Computation on an important engineering problem: non destructive inclusion identi cation in Structural Mechanics. In spite of theoretical results on the identi ability in the linear elastic case 9], the standard deterministic methods had to face a ill-conditioned problem, and has proved inaccurate and unstable. On the opposite, the evolutionary method demonstrates powerful on the simpli ed problem of linear elasticity involving two materials of known characteristics, widely outperforming the to-date results of specialized numerical methods.

However, the { almost { standard representation of structures by arrays of bits (based on a given mesh) proved to be totally inaccurate to that type of problem, and the breakthrough results were obtained using non-standard representations together with speci cally designed genetic operators. The Vorono representation uses a variable number of labeled Vorono sites to de ne a partition, while the Holes representations uses a list of elementary shapes over the structure: these elementary shapes indicate which parts of the structure are not made of the default value.

First, the complexity of these representations are independent of any discretization of the structure. Hence the complexity of the algorithm itself, in terms of number of tness evaluations, only depends on the problem at hand.

Second, the Vorono representation seems more robust with respect to a priori knowledge about the quality of the solution: unless the solution is known to be made of very small parts of one material spread over in the structure, in which case a hole representation with default value to the majority material is the best choice, the Vorono representation performs as good or even better than the rectangle and triangle representations tested here. Moreover, the Vorono representation is independent of any guess regarding the shape of the inclusion itself, whereas the holes representation perform poorly if the wrong guess is made.

Figure 2 :

 2 Figure 2: Evolutionary approach for the inverse problem of mechanical inclusion identi cation: The tness of the individual at hand is achieved by comparing the actual experimental results with the numerical results obtained by Finite Element Analysis on the structure that the individual represents.

Figure 3 :

 3 Figure 3: The experimental structure is xed on its left boundary, and normal forces are applied on di erent points of the free boundary. The large dots are the 9 real-world-like measure points.

Figure 4 :

 4 Figure 4: The experimental structure, and the four test cases, respectively termed corner, checkerboard, asymmetric and symmetric.

Figure 5 :

 5 Figure 5: The Vorono representation crossover operator. A random line is drawn across both diagrams, and the sites on one side are exchanged

Figure 7 :Figure 8 :

 78 Figure 7: A successful run of the Vorono representation on the corner problem using boundary tness. Plots of the best individual at di erent generations of the evolution

Figure 9 :

 9 Figure 9: Robustness to noise. (a): with 2% noise. (b): with 5% noise.

Figure 11 :Figure 12 :

 1112 Figure 11: Best structure after 300 generations on the asymmetric problem, on the 48 48 mesh (colors are inverted on the rectangle plot).

Table 1 :

 1 O -line performances of GAs for all four 24 24 test problems and di erent representations. The best and median result from 21 runs of 300 generations each are presented. The gures are the L 2 error on the displacement at the measure points: A null best error means that the actual global solution was found exactly at least once. resentations, with signi cant di erence for the symmetric and checkerboard problems.

	Triangle
	Bitarray
	1-Rectangle
	Voronoi

Table 2 :

 2 O -line performances of GAs for two test problems on di erent representations using a 12 12 mesh after 300 generations.

	Asymmetric best median best median Checker Voronoi 2.6e-05 0.0003 0.0005 0.0005 Rectangle 1.6e-05 4e-05 0.0015 0.0015 Bitarray 0.57 0.57 0.01 0.01

Table 3 :

 3 Partial results of o -line performances for two test problems on di erent representations using a 48 48 mesh. From one or three runs are aggregated here.

This section gives a detailed presentation of the mechanical problem of inclusion identi cation, and states the simpli cation hypotheses made throughout this paper.

another approach to the relation between operators and tness can be found in 19].


Marc Schoenauer
email: marc.schoenauer@polytechnique.fr

Shape Representations and Evolution Schemes

come   L'archive ouverte pluridisciplinaire

Introduction

Representation is long acknowledged a central issue for Evolutionary Computation. From early discussions on the now out-of-date Michalewicz 1992] problem of encoding real parameters in GAs (e.g. binary vs Gray coding Caruna and Scha er 1988]) to recent works on comparing representations for the TSP problem Radcli e and Surry 1994], EC researchers have tried to characterize the desirable properties of the mapping from the genotype space, where evolution takes place, into the phenotype space, where environmental pressure acts Fogel 1995a, Fogel 1995b]. However, most of these works either present general recommendations and heuristics (e.g. the degree of degeneracy in the representation should be as small as possible to avoid loosing information) or focus on xed-length representations (e.g. binary encoding for the TSP problem).

This paper studies variable-length representations for two-or three-dimensional shapes, in the eld of Topological Optimum Design (TOD): the goal is to nd a structure (a shape in a given design domain) having prescribed mechanical properties and minimal weight. The straightforward representation for shapes that has been exclusively used in past evolutionary attempts on Optimum Design problems Jensen 1992, Chapman, Saitou and Jakiela 1994, Chapman and Jakiela 1995, Kane, Jouve and Schoenauer 1995] is a bitstring representation based on a mesh of the design domain. However, this representation hardly scales up, neither when considering 3-D shapes, nor when the accuracy of the mechanical behavior of the structure (depending on the size of the mesh) becomes a central issue.

Therefore two other representations are designed; both allow to dissociate the representation of the shape and the accuracy of the tness computation. These representations are variable-length representation with a high degree of degeneracy; nevertheless preliminary results show they signi cantly outperform the (non-degenerate) bitarray representation.

Forthcoming section 2 presents the Mechanical background of the Optimum Design problem and brie y recalls the results obtained by standard deterministic approaches.

Section 3 describes how Evolutionary Computation addresses some limitations of these standard methods, by handling 2-D shapes as bitstrings (rather standing for bit-arrays). Speci c evolution operators have been designed to overcome the geometrical bias of standard operators. Nevertheless this representation is intrinsically limited in the sense that the accuracy of the tness is linked to the size of the chromosome.

Section 4 therefore introduces two new representations for 2-D and 3-D shapes, overcoming the above limitation. These are experimentally studied on a benchmark problem of TOD, and compared using the three main evolutionary schemes, namely GAs, EP and ES. Our results suggest the performance depends on two main features: the degree of epistasis of a representation, i.e. the degree the expression of one gene depends on the other genes, and the symmetry of the representation with respect to that of the optimal solutions.

Last, some avenues for further research are presented.

2 Topological Optimum Design

Background

Optimum Design in Structural Mechanics consists in nding the best design for a structure in an initial domain of IR 2 or IR 3 .

More precisely, identifying a shape amounts to nd a partition of the design domain (a subset of IR n , n = 2 or n = 3), into two subsets; one subset is the structure while the other subset represents void 1 .

The optimality criterion is determined by the mechanical properties of the structure, depending on the constitutive law (behavior) of the material; only hyper-elastic materials [START_REF] Ciarlet | Mathematical Elasticity, Vol I : Three-Dimensional Elasticity[END_REF] will be considered in the following.

In most cases of Optimum Design, the goal is to minimize the weight of the structure (i.e. the amount of material), while meeting some engineering requirements for given loading cases (e.g. forces, pressure, prescribed displacements, : : :). This problem is of utter importance for part manufacturers, as a small decrease in the weight of some widely used part results in large cuts in the manufacturing cost.

Deterministic state of the art

Two contexts are distinguished in Structural Optimum Design:

When the solution is sought as the continuous deformation of a given initial shape, iterated small modications of the shape are a deterministic way toward the solution. The methods of domain variation, or sensitivity analysis Cea 1981], are based on gradientlike optimization techniques in that context. But continuous deformations do not a ect the topology of the structure (i.e. its number of holes). When the topology of the solution is unknown, the problem amounts to Topological Optimum Design (TOD); the only deterministic method in that context, to the best of the author`s knowledge, is the homogenization method [START_REF] Bendsoe | [END_REF]Kikushi 1988, Allaire and[START_REF] Allaire | Optimal design for minimum weight and compliance in plane stress using extremal microstructures[END_REF]. Homogenization proceeds by rst relaxing the problem and considering probabilistic shapes in the design domain: the density of material ranges in 0; 1] instead of being either 0 (for void regions) or 1 (for the structure itself). Theoretical results ensure that the optimal solution lies in this superset of probabilistic shapes, and can be approximated by deterministic gradient-based methods.

1 A formally equivalent problem is that of inclusion identi cation Constantinescu 1994], where the goal consists of identifying the repartition of two materials with di erent mechanical properties from the global behavior of the material.

M Ω

Figure 1: A TOD benchmark problem: The 2 1 cantilever plate. The design domain is xed on its left boundary and a force is applied at point M.

The optimal probabilistic shape is then mapped into a real shape. Spectacular results in two and three dimensions have been obtained by the homogenization method [START_REF] Allaire | Shape optimization by the homogenization method[END_REF]]. The greatest limitation of this approach is that dealing with probabilistic shapes is possible so far in the frame of linear elasticity only. Moreover, even in this restricted frame, homogenization method can handle neither multiple loading cases, nor loading applied on the unknown boundary (e.g. a uniform pressure on one side of the structure).

A TOD benchmark

The classical benchmark problem of topological optimization in two dimensions is the cantilever plate, described in Figure 1: a plate is xed on one of its boundary, and one (or more) punctual forces are applied on the other end of the design domain. The goal is to minimize the weight of the structure while staying below some upper bounds on the maximal displacement of the points where the forces are applied.

A shape is therefore evaluated from two criterions: its weight (or equivalently its area in the two-dimensional case), and its mechanical behavior under the prescribed loading. The numerical simulation of the mechanical behavior of the structure is achieved using the Finite Element Method (FEM) Zienkiewicz 1977[START_REF] Ciarlet | The Finite Element Method for Elliptic Problems[END_REF]], which implies a discretization of the shape into small elements, termed mesh.

Formally, the TOD problem is a constrained problem: one aims at minimizing the weight of the shape, while satisfying the prescribed mechanical requirements.

\Standard" Evolutionary TOD

This section brie y describes to-date results obtained on TOD problems in the frame of Evolutionary Computation.

Bitarray Representation

As stated above, the computation of the tness of a shape involves meshing that shape. In that context, the most natural representation for shapes is the bit-array representation, directly based on a xed regular mesh of the design domain into rectangular elements. Each element of the mesh is given a boolean value (0 or 1) indicating which subset it belongs to (material of void). This representation is straightforward from the point of view of both FEM (this is the simplest possible mesh of the domain) and GAs (the representation can be viewed as a bitstring). All previous stochastic works addressing the topological optimum design problem have adopted that representation, be they based on simulated annealing Ghaddar, Y. Maday, and A. T. Patera 1996] or GAs Jensen 1992, Chapman, Saitou and Jakiela 1994, Chapman and Jakiela 1995].

Evolution Operators

However, standard genetic operators appear ill-suited for the above representation of shapes, in the sense they suffer from geometrical bias. Bit-arrays are not bit-string: one-dimensional operators can only exchange horizontal parts of the design domain; if the good schemata involve vertical bands of the structure, these tend to be disrupted by any one-dimensional crossover, and the building block hypothesis therefore does not apply.

Two-dimensional crossover operators, exchanging twodimensional regions of the shapes (e.g. rectangular blocks or regions separated by random lines), have therefore been purposely designed. It has been demonstrated on the cantilever plate problem that these speci c twodimensional crossover operators clearly outperform the standard one-and two-points bitstring crossover operators Kane and Schoenauer 1995].

Constrained Optimization

Many methods have been designed to constrained evolutionary optimization (see Michalewicz 1995] for a survey of such methods). In particular, the method described in [START_REF] Schoenauer | Constrained GA optimization[END_REF] has been successfully applied to another problem of structural mechanics, the optimization of truss structures. But as the focus of this work is representation, the standard method of penalization was chosen, being a non speci c robust way of handling constraints.

The rst draft for the tness function thus has the following expression:

F = Area + (D Max D Lim ) + (1)
where D Max is the maximal displacement of the structure when the prescribed force is applied (computed using the FEM), D Lim the imposed limit value for the dis-placement and is a positive user-supplied penalty parameter (a + denotes the positive part of a). Though some mechanical di culties in fact lead to slightly more complex expression of the tness (e.g. structures not connecting the xed boundary and the loading are not valid solutions, see Kane andSchoenauer 1995, Schoenauer 1995]), equation (1) will be considered in the following for the sake of simplicity.

Successes

Using speci c two-dimensional evolution operators on bitarray representation, EC has successfully tackled TOD problems that could not be addressed by other techniques Kane, Jouve andSchoenauer 1995, Kane 1996] | at the expense of large computational time. The most signi cant results are the following:

Overall, EC can accommodate any mechanical model for which there exists a numerical simulation algorithm. This is con rmed as the results obtained for the large displacement model appear to be the rst results in Optimum Design of structures in nonlinear elasticity.

The optimization can take into account more than one loading case (as in the design of a bicycle), as well as loading applied on the unknown boundary of the structure (as in the case of the underwater dome).

In some situations, many optimal, or near-optimal solutions exist. Evolutionary algorithms, using for instance the sharing scheme [START_REF] Goldberg | Genetic algorithms with sharing for multi-modal function optimization[END_REF], are able to provide the engineer with a range of such solutions, allowing him to take into account inarticulate criterions.

Drawbacks

The limitations of these results come from the following fact.

The accuracy of Evolutionary TOD is dictated by the size of the mesh underlying the FEM analysis and the tness computation: the above mentioned results were obtained on rather coarse meshes (e.g. 10 20) whereas real-world problems and accurate analyses of the mechanical behavior of the shapes require much ner meshes (e.g. 100 200).

Increasing the size of the mesh would not only increase the cost of the tness computation (which is roughly quadratic in the size of the mesh) | but also the size of the chromosomes. And increasing the size of the individuals would require to increase in turn the size of the population and the number of generations to reach the same level of convergence [START_REF] Cerf | Une th eorie assymptotique des algorithmes g en etiques[END_REF], Cerf 1996] proved that the minimal size of the population for a GA to converge increases linearly in term of the size of the bitstring).

Finally the bit-array representation poorly scales up when re ning the mesh | not to speak of handling 3-D shapes... Therefore other representations have been designed in order to overcome this limitation, and dissociate the complexity of the representation and the accuracy of the evaluation.

4 Variable-length representations for shapes Consider now a (variable length) list of Vorono sites, each site being labeled 0 or 1. The corresponding Vorono diagram represents a shape (a partition of the design domain into two subsets), if each Vorono cell is labeled as the associated site (here the Vorono diagram is supposed regular, i.e. to each cell corresponds exactly one site). Example of Vorono representations can be seen in Figure 2. The Vorono sites are the dots in the center of the cells. Note that Vorono representation of shapes does not depend in any way on the mesh that will be used to compute the behavior of the shapes. Furthermore, Vorono diagrams being de ned in any dimension, the extension of this representation to IR 3 and IR n is straightforward.

Evolution operators

The evolution operators on the Vorono representation are inspired from both the two-dimensional crossover operators designed for the bit-array representation and the usual operators for variable length representations:

The crossover operators exchange Vorono sites on the basis of geometrically-based choice. In this respect it is similar to the speci c bitarray crossover described in Kane and Schoenauer 1995]; moreover, this mechanism easily extends to any dimension [START_REF] Kahng | Toward more powerful recombinations[END_REF]. Figure 2 demonstrates an application of this crossover operator. a rst mutation operator performs a Gaussian mutation on the coordinates of the sites, or randomly ips the boolean attribute of some sites; a "standard" mutation for variable-length representations adds or deletes some sites from the list. An important remark is that this representation presents a high degree of epistasis (the in uence of one site on the physical shape is modulated by all neighbor sites). This will be discussed in more details in section 5.1.

Practically, the tness of all shapes is evaluated using the same xed mesh 2 . A shape described by Vorono sites is thus mapped on this xed mesh: the subset (material or void) an element belongs to is determined from the label of the Vorono cell in which the center of gravity of that element lies.

H-representation

Another representation for shapes is based on an oldtime heuristic method in TOD: from the initial design domain, remove material where the mechanical stress is minimal, until the constraints are violated. However, the lack of backtracking makes this method useless in most TOD problems. Nevertheless, this idea gave birth to the \holes" representation Dejonghe 1993], later termed Hrepresentation.

The representation

The design domain is by default made of material, and a (variable length) list of \holes" describes the topology of the structure. These holes are elementary shapes taken from a library of possible simples shapes. Only rectangular holes are considered at the moment. On-going work Seguin 1995] is concerned with other elementary holes (e.g. triangles, circles).

Example of structures described in the H-representation are presented in Figure 3. The rectangles are taken in a domain larger than the design domain, in order not to bias the boundary parts of the design domain.

Evolution operators

The evolution operators are quite similar to those of the Vorono representation: crossover by geometrical (2D or 3D) exchange of holes (see Figure 3 for an example); mutation by Gaussian modi cation of the characteristics (coordinates of the center, width and length) of some holes; mutation by addition or deletion of a hole;

The H-representation, as the Vorono representation, is independent from any mesh, and hence its complexity does not depend on any required accuracy for the simulation of the underlying physical phenomenon. Its merits and limitations will be discussed in the light of the experimental results presented in next section.

As for the Vorono representation, the simulated behavior of the shapes is computed on a given xed mesh, to limit the numerical noise due to re-meshing. The criterion to decide which subset an element does belong to, the higher the computational cost). Hence, the tnesses of di erent structures that are to be compared should be performed with the same mesh. Otherwise, numerical noise due to re meshing might hide the actual di erences in the mechanical behavior of di erent structures.

Figure 3: The H-representation crossover operator. A random line is drawn across both structures, and the holes on one side are exchanged. is based on whether its center of gravity belongs to a hole (in which case the whole element is void) or not.

Preliminary comparative results

The same benchmark problem of Optimal Design has been used for the rst comparative results for the three representations of shapes presented in preceding section (see Schoenauer 1995]). As expected, and even for a rather coarse mesh, the bit-array representation is outperformed by all other representations in terms of computational cost as well as in terms of quality of the solution. Hence, it will not be considered any more in the rest of the paper.

Both the Vorono representation and the Hrepresentation derive a quasi-optimal solution. But these solutions are more rapidly found as H-based than Vorono -based shapes.

Tentative explanations for that are proposed in next section.

Epistasis and Evolutionary Schemes

Epistasis and Symmetry of Representations

In the biological context, the epistasis is a measure of how the expression in the phenotype of one single gene is in uenced by the other genes of the genotype. In Evolutionary computation, epistasis has a strong in uence on the way genetic material is transmitted from parents to o spring, and how it can be modi ed by evolution operators thereafter applied.

In the bit-array representation, the contribution of a given bit to the phenotype is clear, and it cannot modi ed by any other bit. The epistasis is here minimal, all genetic material transmitted to the o spring is expressed, i.e. is dominant. Such a situation can be termed strong transmission.

In the Vorono representation, the in uence of a site on the nal phenotype can be greatly modi ed by other neighbor sites. If one site labeled void becomes surrounded by sites labeled material, its in uence on the nal phenotype almost vanishes. Some part of genetic information transmitted to o spring is recessive, and, by contrast to the bit-array situation, this case will be termed weak transmission.

In the H-representation, the situation is even more complex: Only holes are strongly transmitted. When a hole is transmitted from parent to o spring, the whole area will stay a hole. On the opposite, no such strong transmission is achieved for \non-holes". The transmission is asymmetric, strong for the hole value, and weak for the default value (i.e. any hole covering the same area changes its value). One can say the hole-value represent dominant genetic material, while the default value is recessive. The bottom-up approach of GA, trying to recombine building blocks to reach the optimal solution should be much more hindered by epistasis than the top-down approach only relying on the competition between phenotypes to reach an optimal point of the search space. In order to check these points, experiments were conducted on the benchmark cantilever problem (section 2.3) using the two di erent approaches of EC.

Evolutionary Schemes

Evolutionary Algorithms (EA) crudely mimic the evolution of a population of points of the search space. This population is usually initialized randomly, and undergoes a succession of generations. The general outline of a generation can be viewed as SELECT parents from the population APPLY evolutionary operators to the selected parents to generate offspring REPLACE some parents by some offspring

The most widely used EAs, namely Evolutionary Programming Fogel, Owens and Walsh 1966], Evolutionary Strategies Schwefel 1981] Of course, numerous variations of these canonical algorithms exist and are being used in practical applications on a pragmatic basis Michalewicz 1992, Fogel 1995a[START_REF]Evolutionary Algorithms in theory and practice[END_REF]. However, the a priori adequation of an evolutionary scheme to a given application is still an open question.

Experimental results

Experiments were conducted using two evolutionary schemes:

The standard generational GA, with population size of 100, using ranked-based selection, crossover rate of 0.6, mutation rate per individual of 0.2, with Gaussian mutations of xed variances. A (15+100)-ES (15 parents generate 100 o spring, the best 15 among parents + o spring become the parents of the next generation) with strength of Gaussian mutations depending on the tness of the individual at hand as in EP Fogel 1992]. In both cases, the maximum number of generations allowed is set to 100, and the algorithm stops if no improvement is observed during 10 generations. These stopping criterions are severe, and most runs did not reach convergence. But, in contrast with Kane, Jouve andSchoenauer 1995, Kane 1996], the goal here is to observe the behavior of the algorithms with the perspective of fast convergence, rather than to reach the optimal shape by all means.

As one of the goals is to study the impact of the symmetry of the representations on the evolutionary algorithm, three representations are considered: the Vorono representation (section 4.1), the \Hole" representation and the \Plate" representation, which are both H-representation (section4.2), where the default value is material in the case of the Hole-representation, and void in the case of the Plate-representation (the structure is here an assembly of small plates).

In order to test the in uence of both the epistasis and the symmetry of the representation interacting with the .Comparative results between combinations of representations and evolutionary schemes: Averaged (over 30 independent runs) best tness along generations. The constraint on the displacement is very strong (limit value 0.2) and the penalty parameter is small. evolution scheme, di erent instances of the cantilever plate problem described in section 2.3 are considered. For all problems, the design domain is the 2 1 rectangle, the mesh is a 20 10 regular mesh (i.e. all elements are here squares) and the same force is applied at point (2; 0:5) (point M in Figure 1).

The maximum value of the displacement is assigned di erent values, increasing from 0.2, the actual displacements of the full horizontal beam joining point M and the xed boundary of thickness 0:8 to 37.10, the actual displacement of the full beam of thickness 0:1. Note that only in this latter case is the optimum known: due to the coarse mesh used here, the minimal structure connecting point M and the xed boundary is the full beam of thickness 0:1, which does indeed respect the constraint.

Finally, the penalty parameter is varied too: intuitively, the higher the penalty factor, the more dicult the problem (feasible regions are environed by steep slopes).

It is fair to say that only tendencies could be observed, and no absolute conclusion can be drawn.

Both H-representations globally outperform the

Vorono representation. So the high degree of epistasis seems to penalize the representation here.

2. ES and GA are hardly distinguishable, except { when the limit value for the displacement is small (strong constraint), in which case GA signi cantly outperforms ES for all representations (see Figure 4). However, this phenomenon decreases when the penalty parameter increases. A possible explanation could be that the elitist ES gets stuck more easily in the rst feasible local optimum for small values of the penalty parameter, whereas both algorithms encounter this same di culty for large values of the penalty parameter. Figure 5: Comparative results between combinations of representations and evolutionary schemes: Averaged (over 30 independent runs) best tness along generations. The constraint on the displacement is weak (limit value 6.40) and the penalty parameter is large.

{ for the Vorono representation with large limit for the displacement (weak constraint) in which case ES performs better than GA (see Figure 5). This is the only hint meeting the a priori expectation that epistasis should favor the top-down against the bottom-up approach. 3. The \Hole" representation demonstrates slightly better performances than the \Plate" representation in almost all cases (e.g. in Figure 5) except when the optimal solution contains about as much material as void, as can be seen in Figure 4: a value of 0:5 for the tness of feasible structures corresponds to exactly the same amount of void and material.

6 Further directions for the comparison of representations

Overall, these experiments show that far too many parameters are involved in the TOD problem to make it a good test-bed for shape representation ! A next step will be to design such an adequate test-bed, as an unconstrained problem the solution of which is known, with tunable amount of void and material, and scalable degree of di culty. On-going work is concerned with the problem of (non-destructive) identi cation of the repartition of two materials in a structure (e.g. scories in a steel piece) using mechanical experiments Constantinescu 1994]. Some criterions investigated in the literature will guide systematic experiments:

The tness variance theory of Radcli e Radcli e and Surry 1994] studies the variance of the tness as a function of the order of an extension of schemas called formae Radcli e 1991], and, simply put, shows that the complexity and di culties of evolution increases with the average variance of the tness. But if the formae and their order (or their precision) are well-de ned on any binary representation, including the bit-array representation of section 3, it is not straightforward to extend these de nitions to variable length representations discussed above. Moreover, Radcli e's tness variance does not take into account the possible evolution operators. Further step in that direction would be to study the variance of the change of tness with respect to a given evolution operator (e.g. the Gaussian mutation of Vorono sites for di erent standard deviations), as in Fogel 1995b] The tness distance correlation of [START_REF] Jones | Fitness distance correlation as a measure of problem difculty for genetic algorithms[END_REF] studies the correlation between the distance to the optimal point and the tness. Simply put again, the idea is: the stronger this correlation, the narrower the peak the optimum belongs to, and the more di cult the problem. Conjectures based on this remark are experimentally conrmed in the GA-frame. Nevertheless, the di culty in shape representation is to de ne a distance which is meaningful for both the representation and the problem at hand. The rst important issue is whether the considered distance should be purely genotype-based (i.e. de ned on the coded individual only) or partially or totally phenotype-based (i.e. de ned on the same space than the tness function).

Conclusion and further work

Di erent evolutionary approaches for shape optimization have been presented. The emphasis has been put on the representation: the simple bit-array representation allowed signi cant advances in the domain of topological shape optimization (e.g. the rst results in nonlinear elasticity), but hardly scales up, in the sense that the accuracy of the tness computation is commanded by the size of the individuals.

Two other, variable-length, representations, the Vorono and the H-representations have been designed to overcome this issue; they demonstrate good results on the Topological Optimum Design problem, signi cantly and consistently outperforming the bitarray representation.

The main di erence between those latter representations a priori appears their respective degree of epistasis, i.e. the way one gene (Vorono site or hole) may modify the phenotypic traits due to other genes. Systematic experiments have been conducted to see how the degree of epistasis interacts with the evolutionary scheme (bottom-up as in GAs or top-down as in ES) and the symmetry of the representation ... and mainly demonstrate that the classical benchmark in Topological Optimum Design is a problem far too complex to serve as a test-bed for evaluating shape representations.

A number of powerful representations of shapes remain to be investigated in the frame of evolutionary optimization:

The CAD community uses splines de ned from control points to describe and manipulate shapes. The L-system paradigm Prusinkiewicz and Lindenmeyer 1990] uses grammar rules to simulate plant growth. The resulting \plants" can be viewed as shapes. Moreover, it has been demonstrated that such grammar rules can be optimized using Evolutionary Computation algorithms Simms 1994].

The fractal theory of Iterated Function Systems Barnsley 1988] o ers another possible representation, for which the inverse problem has been successfully tackled by means of Evolutionary Computation [START_REF] Garigliano | [END_REF]al. 1993, Lutton andMartinez 1994].

It is emphasized however, that a key problem of evolutionary computation yet is the choice of an adequate representation. Widening this choice would only ask more loudly for judicious choice criterions; and it is our conviction that the rst, and maybe the more important step on this way, would be to elaborate a convenient test-bed, of tunable di culty.

On-going work considers another shape identi cation problem to this end, namely the identi cation of inclusions (e.g. scories in a piece of steel) Constantinescu 1994]. This latter problem presents two main advantages: it is unconstrained and fully tunable (the optimal solution can be xed a priori).

4. 1

 1 The Vorono representation A possible way of representing shapes comes from computational geometry, more precisely from the Vorono diagram theory. The ideas of Vorono diagrams are already well-known in the FEM community, as a powerful tool to generate good meshesGeorge 1991]. However, the representation of shapes by Vorono diagrams and their evolutionary optimization seems to be original.4.1.1 Vorono diagramsConsider a nite number of points V 0 ; : : : ; V N (the Vorono sites) of a given subset of IR n (the design domain). To each site V i is associated the set of all points of the design domain for which the closest Vorono site is V i , termed Vorono cell. The Vorono diagram is the partition of the design domain de ned by the Vorono cells. Each cell is a polyhedral subset of the design domain, and any partition of a domain of IR n into polyhedral subsets is the Vorono diagram of at least one set of Vorono sites (see[START_REF] Preparata | [END_REF] Shamos 1985, Boissonnat and M. Yvinec 1995] for a detailed introduction to Vorono diagrams, and a general presentation of algorithmic geometry).

Figure 2 :

 2 Figure 2: The Vorono representation crossover operator. A random line is drawn across both diagrams, and the sites on one side are exchanged

  Figure4: .Comparative results between combinations of representations and evolutionary schemes: Averaged (over 30 independent runs) best tness along generations. The constraint on the displacement is very strong (limit value 0.2) and the penalty parameter is small.

  and Genetic Algorithms[START_REF] Holland | Adaptation in natural and arti cial systems[END_REF], are instances of this general scheme:EP and ES do not use initial selection process while GA uses tness-based stochastic selection (e.g. roulette wheel).

	EP uses mutation only, each parent generating one o spring, ES uses mostly mutation, each parent gen-erating usually more than one o spring and GA use mostly recombination (crossover), two parents gener-ating two o spring.
	EP replaces the parents using a stochastic tourna-ment among parents and o spring, ES selects the best individuals among parents and o spring in ( + ) ES or among o spring only in ( ; ) ES as the new parents, and GA globally replaces all parents with all o spring.

This intends to limit the bias due to the unavoidable numerical noise of FEM (the ner the mesh, the lower the numerical error,

Acknowledgements

Grateful thoughts to Mich ele Sebag for useful discussions and helpful advices.