M S Berkani

S G I U R G E A

C Espanet

andC J L Coulomb

Kieffer

Study on Optimal Design Based on Direct Coupling Between a FEM Simulation Model and L-BFGS-B Algorithm

Keywords:

This paper investigates the capability of a Quasi-Newton optimization algorithm, the L-BFGS-B one, to reduce the drawbacks of the direct coupling with FEM models. After a short description of the L-BFGS-B algorithm, the authors have tested it to the TEAM Workshop Problem 25. L-BFGS-B algorithm issues linked to FEM estimation of the objective function have been analyzed. Finally the authors present the results of an electromagnetic device optimal design, using a large number of parameters (47).

I. INTRODUCTION

I N the context of a growing part of electrical power in the in- dustry, there is a strong need to develop fast and efficient design methodologies. In this situation, a numerical model based on the finite element method (FEM), used in the optimization process, provides both the benefits of a fast development and of the accuracy. On the other hand, several drawbacks have been identified in the past few decades; two of them seem to be the most important. First, the FEM model could require an important computation time to estimate the objective function, which doesn't enable to control the global computation cost of the optimization. A second drawback is linked to the objective function that may be noisy due to the effects of mesh variation [START_REF] Ratnajeevan | Fictitious minima of object functions, finiteelement meshes, and edge elements in electromagnetic device synthesis[END_REF]. For example, Fig. 1 shows a zoom on the variation of the objective function estimated by the FEM model of the TEAM Workshop Problem 25 (TEAM 25) [START_REF] Y O S H I D A | Investigation of simulated annealing method and its application to optimal design of die mold for orientation of magnetic powder[END_REF], for the parameter variation. In this example the others parameters are fixed at their initial value (see Table I). This noise could perturb the gradient-based algorithms. One way to reduce computation time is to substitute the objective function, estimated by the FEM, with an approximation model. This strategy requires less computation time, since it reduces the number of FEM simulations. Different approximation methods have been proposed such as polynomial regression [START_REF] Kim1 | Optimal design of slotless-type PMLSM considering multiple responses by response surface methodology[END_REF], radial basis function [START_REF] Ishikawa | A combined method for the global optimization using radial basis function and deterministic approach[END_REF], multiquadratic basis functions [START_REF] Alotto | A multiquadrics based algorithm for the acceleration of simulated annealing optimization procedures[END_REF], and kriging [START_REF] L E B E N S Z T A J N | Kriging: A useful tool for electromagnetic device optimization[END_REF]. These solutions provide a continuous and smooth approximation and are compatible with deterministic and stochastic optimization algorithms. Nevertheless, to calibrate an approximation model, several points are needed on each parameter that implies a number of FEM estimations, which depends exponentially on the parameters number. This negative feature can be avoided by choosing in the approximation model only the most significant parameters of the FEM model [START_REF] Costa | An adaptive method applied to the diffuse element approximation in optimization process[END_REF]. Approximation errorsmayalsoinfluence the precision of the final result.

The purpose of this contribution is to investigate the capability of a Quasi-Newton based algorithm, the L-BFGS-B one, to produce a fine result by using direct coupling with an interfaced FEM model. To achieve this study, we consider the TEAM 25 problem. In the Section II, the authors describe the L-FBGS-B algorithm. Then, in the Section III, they present the optimisation problem formulation of TEAM 25. In the Section IV, they analyse the particularities of L-BFGS-B algorithm coupled to a FEM model through the use case TEAM 25. Finally, in Section V, the authors present the results of an electromagnetic device optimal design; the interest of this example is that the optimization uses 47 parameters, which is quite large for a FEM model.

II. L-BFGS-B OVERVIEW

As described in [START_REF] Byr D | Al i mi t edmemo ryal gorithm for bound constrained optimization[END_REF] the L-BFGS-B algorithm is a Quasi-Newton based algorithm that allows finding the minimum of an n dimensional function:

(1)

The feasibility domain accepts simple bounds on variables and can be defined by a hyper cube with minimal and maximal limit of each parameter as shown in the following: [START_REF] Y O S H I D A | Investigation of simulated annealing method and its application to optimal design of die mold for orientation of magnetic powder[END_REF] At each iteration, the last values of the function, and of its gradient, are computed; the hessian matrix approximation is updated. Quadratic approximation of the objective function is then used at different steps of the algorithm:

(3)
The updating mechanism of is realized by using the last p correction pairs defined by: (

representing the evolutions of the current point and of the objective function gradient between two successive iterations. This is one aspect that can explain the behavior of the algorithm coupled with FEM model. Three stages compose an iteration of the L-BFGS-B algorithm (Fig. 2): First Step-Cauchy Point Computation: Theobjectiveofthis step is to find the local minimum of quadratic approximation of the objective function, starting from the current point ,on the path defined by the projection of the steepest descent direction on the feasible domain. On this projection path, a point x can be expressed by its travel distance from the starting point . Traveling on the opposite gradient direction, one by one, the n variables of objective function will reach those limits at values of For the corresponding point x can be expressed: [START_REF] Alotto | A multiquadrics based algorithm for the acceleration of simulated annealing optimization procedures[END_REF] with ,a n d , the steepest descent direction, projected on the feasibility domain.

On the segment , the quadratic approximation can be expressed as a one-dimensional quadratic expression in : [START_REF] L E B E N S Z T A J N | Kriging: A useful tool for electromagnetic device optimization[END_REF] where the scalar values are obtained using the values of and . The optimal value can be estimated:

. Its corresponding point represents the Cauchy point: [START_REF] Costa | An adaptive method applied to the diffuse element approximation in optimization process[END_REF] The indexes of the variables fixed on their limits, during the Cauchy point computation, are grouped in a set . Second Step-Subspace Minimization: It represents a Quasi Newton stage starting from the Cauchy Point obtained previously. A new descent direction is obtained by searching the minimum of (3) on the free variable subspace, by using the following formulation: [START_REF] Byr D | Al i mi t edmemo ryal gorithm for bound constrained optimization[END_REF] In this step, the feasibility domain is reduced at a subspace of the feasibility domain , by considering as free variables, the variables that are not fixed on limits during the Cauchy point step: [START_REF] K I E F F E R | D e s i g n of an electrical machine for an integrated magnetocaloric coupling system[END_REF] The others variables are fixed on their boundary value obtained during the Cauchy point step. [START_REF] Centner | U. Schafer Comparison of High-Speed Induction Motors Employing Cobalt-Iron and Silicon Electrical Steel ICEM-2008 Conf[END_REF] TheThir dStep-LineSear ch: It achieves a linear search on to find a new point that satisfies a sufficient decrease condition and a curvature condition (known also as Wolfe's conditions). The line search implies successive one-dimensional interpolations (quadratic or cubic) to compute an acceptable step from point to point. This stage requires additional estimations of the function and of its gradient, , if doesn't satisfy the Wolfe conditions.

III. APPLICATION The algorithm has been tested with using the classical optimization problem TEAM 25 [START_REF] Y O S H I D A | Investigation of simulated annealing method and its application to optimal design of die mold for orientation of magnetic powder[END_REF]. The goal is to optimize the geometry of a die mold used in the production of permanent magnets. The die mold geometry is described by 4 parameters: R1,L2,L3andL4asshowninFig.3 [START_REF] Y O S H I D A | Investigation of simulated annealing method and its application to optimal design of die mold for orientation of magnetic powder[END_REF].Wehavetofind the optimal value of the four geometrical parameters, to obtain a radial magnetic fieldonthearcef. Reference value of the flux density presents a sinusoidal distribution of and components:

(11)

This objective function is obtained by minimizing the residual sum of squares between the reference values and FEM simulated values, on ten different points on the arc ef:

(12)
The objective function is strictly positive. Sufficiently small value of the objective function can be chosen as STOP criteria for this problem. The definition domain of each parameter has been normalized into the interval [0,1]. In the gradient evaluation, partial derivatives are obtained using the Newton's difference quotient. The same normalized variation has been considered for all normalized parameters. Five FEM computations

are used to compute

, for a point . The optimization parameters, their initial values and feasibility domains are given in Table I.

IV. RESULTS AND ALGORITHM ANALYSIS

In this section, some important issues of coupling L-BFGS-B algorithm with a FEM model are analyzed, such as the compatibility with the FEM model and the influence of and parameters. We analyze also two observed difficulties of this methodology due to the Line Search stage and STOP criteria of the algorithm.

A. Compatibility of the Algorithm With a FEM Model

Some previous work shows the capability of response surfaces to deal with the objective function presenting a mesh noise. On the same optimization problem a lower number of iteration has been obtained in [START_REF] Costa | An adaptive method applied to the diffuse element approximation in optimization process[END_REF] after a screening stage, by an adapting response surface applied on 3 parameters. The results obtainedin [START_REF] Costa | An adaptive method applied to the diffuse element approximation in optimization process[END_REF]areresumedinT ableII,(rs1andrs2results).

Using Quasi-Newton research directions, L-BFGS-B implicitly uses quadratic approximations that may have the same behaviour of a response surface, in the case of a noisy objective function. For each iteration, the algorithm operates two implicit local approximations. When using a direct connection with the FEM solver, we can anticipate that this mechanism blur the discontinuities due to the mesh variation. The method can take into consideration an important number of parameters. This enables a strong reduction of the objective function, when an important number of parameters are influent. As it is shown in the Table II, the minimum value of the objective function obtained for and , with 815 FEM estimations was which divides by ten the previous results obtained in [START_REF] Costa | An adaptive method applied to the diffuse element approximation in optimization process[END_REF] (rs2 in Table II). The results of this optimization are given in Table I. Computing time dedicated to the algorithm computation was (s) and the one dedicated to FEM estimation was 8532(s).

B. Influence of (The Gradient Computation Normalized Va r i a t i o n )

The influence of the normalized variation used in the gradient computation, , has been tested by considering 4 values. Fig. 4 shows the evolution of the objective function depending on the number of estimations of for the same number of correction pairs . It can be observed that is influent on the gradient accuracy and therefore on the speed of convergence. This dependence is not linear and the better performance has been obtained for .

C. Influence of the Number of the Last p Correction Pairs Considered in Hessian Approximation

The parameter influences the convergence too. The considered example shows a better convergence for a greater value of . In Fig. 4 we can observe also the evolution of the objective function when the algorithm use for and values. Expression [START_REF] Ishikawa | A combined method for the global optimization using radial basis function and deterministic approach[END_REF] shows that a number p of previous results, estimated on the FEM model, are used to actualize the Hessian approximation. Acting on the number of the last experiments will provide a local behavior of the quadratic approximation (for a small value of ) or global one (for a higher value of ). This parameter has also an influence on the approximation of the Hessian matrix, B. As shown in [START_REF] Byr D | Al i mi t edmemo ryal gorithm for bound constrained optimization[END_REF], a lower value of limits the memory and the time computation necessary to estimate the approximation of the Hessian matrix, used in the Quasi-Newton process. However, the algorithm computation time is negligible in comparison with the FEM computation time. Table II synthesizes the results obtained for different configurations and of the algorithm. In order to have a comparison between the different configurations, and also to compare the results obtained from L-BFGS-B with other results from literature, columns and from Table II shows the necessary number of iterations to obtain FEM computed value of the objective function .L a s t columns "Final " and "Final " represent the result obtained when the optimization process is stopped because the line search process doesn't succeed.

D. Analysis on Line Search Stage

During the line search step, the objective functions and its gradient are calculated successively to find a point that satisfies the Wolfe conditions. As in the line search, and are directly estimated on the Objective Function (FEM estimated), this step can correct the approximation errors linked to the quadratic approximation. Fig. 5(a) presents the L-BFGS-B evolution for and . We observe also the number of line search per iteration. The minimum point (presented in Table II) has been reached after 27 estimations and 19 iterations. A large number of iterations execute only one estimation for the line search step. This happens when satisfies both Wolfe conditions. For a positive definite , [START_REF] Byr D | Al i mi t edmemo ryal gorithm for bound constrained optimization[END_REF] has been proven that that the approximation solution of the quadratic problem defined in the subspace minimization step defines a descent direction for the objective function. One problem in coupling LBFGS algorithm with a FEM model corresponds to the situation where the line search doesn't find satisfying Wolfe conditions. This situation occurs when the objective function is close to the theoretical optimum point. In this situation, the algorithm stops before reaching a convergence criterion. Fig. 5(b) shows the objective function estimation on the last line search segment before the algorithm is stopped, when the line search process doesn't succeed, for and . The real objective function does not descend in this last interval.

One explanation could be the noise of the gradient which is more significant near the optimal point. Even if the quadratic approximation function ( 3) is descending, the real function may increase within the small interval of the line search interval. In Table II, "Final " and "Final " results were obtained by this criterion and it represents the last iteration that corresponds to a decrease of the objective function.

V. L ARGE NUMBER OF PARAMETERS OPTIMIZATION: MAGNETO CALORIC DEVICE OPTIMIZATION

To illustrate the behavior of L-BFGS-B coupled to a practical FEM model having a large number of parameters, we present here the case of a thermo magnetic device design presented in [START_REF] K I E F F E R | D e s i g n of an electrical machine for an integrated magnetocaloric coupling system[END_REF]. The sizing of the rotor magnets (Fig. 6(a)) is carried out in order to obtain a maximum flux density variation in the regenerator air gap, to ensure the performance of the magneto caloric device. To limit the saturation, the extremity of the rotor pole is realized on a cobalt-iron alloy "Vacoflux" with a high saturation flux density [START_REF] Centner | U. Schafer Comparison of High-Speed Induction Motors Employing Cobalt-Iron and Silicon Electrical Steel ICEM-2008 Conf[END_REF]. An initial optimal design has been obtained with using only two parameters. The objective function to minimize is the quadratic error sum between the reference value of the flux density module and the value of the FEM estimated value , for two test zones in )a n d in . The solution has been obtained by a direct coupling of the FEM model and L-BFGS-B algorithm. This optimal geometry and the flux distribution are shown in Fig. 6(b). The reference value and the FEM computed value are presented in Fig. 6(d). To improve the result, the problem has been modified, by increasing the number of parameters of the FEM model from 2 to 47 parameters characterizing the magnets geometry and their orientation. Then, the L-BFGS-B was coupled to this model. The results were obtained after 72 h of computation and about 7000 FEM estimations. The resulting geometry and the magnets orientation

VI. CONCLUSION

The analysis of the L-BFGS-B algorithm and the numerical experimentation proved the capability of this algorithm to produce an accurate result when is directly coupled with an interfaced FEM model. The influence of two parameters has been observed on the algorithm convergence: the gradient computation normalized variation and the last p correction pairs considered in Hessian approximation. A study of the line search stage shows a different behavior when the objective function is FEM estimated. Finally the algorithm capability to solve a practical OP with a large number of parameters (47) has been confirmed on the FEM model of special device for magnetocaloric refrigeration.

Fig. 1 .

 1 Fig. 1. Objective function noise observed on FEM simulation of TEAM Workshop Problem 25.

Fig. 2 .

 2 Fig. 2. L-BFGS-B algorithm overview.

Fig. 3 .

 3 Fig. 3. TEAM Workshop Problem 25 geometry.

Fig. 4 .

 4 Fig. 4. Evolution of the objective function for different p, .

Fig. 5 .

 5 Fig. 5. (a) Number of line search per iteration, (b) Objective function evolution on the last line search.

Fig. 6 .

 6 Fig. 6. (a) Two parameters design of the initial configuration, (b) Magnetic flux distribution, (c) Optimal 47 parameters design, (d) Flux density in the regenerator: reference and computed values.

TABLE I CONSIDERED

 I

GEOMETRIC PARAMETERS

TABLE II NUMBER

 II OF FEM ESTIMATIONS AND OPTIMIZATION RESULTS