Ashok V.Kumar,JOURNAL OF MECHANICAL DESIGN,1999,122(3).271.277

Ashok V.Kumar

Assistant Professor Department of Mechanical Engineering, University of Florida,Gainesville,FL32611

e-mail:akumar@https://www.360docs.net/doc/e7441169.html, A Sequential Optimization Algorithm Using Logarithmic Barriers:Applications to Structural Optimization

A sequential approximation algorithm is presented here that is particularly suited for problems in engineering design and structural optimization,where the number of vari-ables is very large and function and sensitivity evaluations are computationally expensive.

A sequence of sub-problems are generated using a linear approximation for the objective function and setting move limits on the variables using a barrier method.These sub-problems are strictly convex and computation per iteration is signi?cantly reduced by not solving the sub-problems exactly.Instead a few Newton-steps are taken for each sub-problem generated.A criterion,for setting the move limit,is described that reduces or eliminates step size reduction during line search.The method was found to perform well for unconstrained and linearly constrained optimization problems.It is particularly suit-able for application to design of optimal shape and topology of structures by minimizing their compliance since it requires very few function evaluations,does not require the hessian of the objective function and evaluates its gradient only once for every sub-problem generated.?S1050-0472?00?01603-2

?

1Introduction

Engineering design and structural optimization problems are often stated as non-linear programming problems.Typically these problems have a large number of variables and the evaluation of the objective function and the constraints involve computationally expensive analysis.Due to these reasons Schmit and Farshi?1?proposed sequential approximation methods for structural https://www.360docs.net/doc/e7441169.html,monly used approximations include linear approxima-tion leading to Sequential Linear Programming?SLP?and linear-ization in reciprocal variables as well as many hybrids of these two approximations?1–3?.

Sequential approximation techniques create a sequence of sub-problems that are easier to solve than the original problem.The approximate objective functions and constraints of the sub-problem are much cheaper to evaluate and therefore the sub-problems are relatively inexpensive to solve.The function and gradient evaluations of the original optimization problem are typi-cally performed only once per sub-problem.Unlike non-linear optimization algorithms derived by extending Newton’s method, methods based on sequential linearization do not require the hes-sian matrix of the objective function.The hessian matrix can be very large in applications involving large number of variables, therefore,constructing this matrix and its LU decomposition rep-resents a large overhead.As a result,sequential approximation methods often perform better for applications where the number of variables are very large or the second derivatives of the func-tion are very expensive to evaluate.

An important area for application of sequential programming has been structural optimization,where function evaluation in-volves computationally expensive structural analysis.The evalua-tion of gradient of the objective function and its hessian matrix each represent signi?cant computation.In particular,for shape and topology optimization of structures,the optimization problem involves very large number of variables.The optimal shape and topology of structures is obtained by computing the optimal ma-terial distribution within a given feasible region?see,Bends?e and Kikuchi?4?,Bends?e?5?and Kumar and Gossard?6,7??.Sequen-tial approximation algorithms such as CONLIN and MMA?8,9?have recently been used for this problem.

In Sequential Linear Programming,a solution to the non-linear programming problem is obtained by solving a sequence of linear programs created by linear approximations of the objective and constraints.Schmit and Farshi?1?suggested linearization in the reciprocal of the variables.Svanberg?10?further extended this idea,where he suggests a method of moving asymptotes?MMA?. One interpretation of this method can be that move limits are set on the variables using asymptotes for each sub-problem.These asymptotes are moved to make the method stable and to expedite convergence.Svanberg suggests a set of heuristic rules for mov-ing these asymptotes.The dual of the sub-problem is solved?rst and the solution for the primary problem is computed using the dual solution.The dual problems have concave objective function and hence can be solved easily using conventional gradient meth-ods such as the conjugate gradient method.

In this paper a novel sequential approximation method is de-scribed,in which a sequence of sub-problems are generated by linearizing the objective function and setting move limits on the variables using logarithmic barrier functions.This method can be interpreted as sequential linear programming using primal-dual technique?11,12?to generate a descent direction and step size. Move limits on the variables are?exible and are modi?ed each iteration in such a fashion that the upper and lower limits con-verge towards each other.The resulting algorithm is a non-linear programming technique that is particularly useful for applications involving large number of variables and linear constraints.The reason for using logarithmic barriers to set move limits on the variables is that the resulting sub-problems have a diagonal hes-sian matrix and therefore it is inexpensive to use Newton’s method for solving the sub-problems.The principal advantage of the method is that it does not require the exact solution of the sub-problems and hence requires very little computation per itera-tion.It is insensitive to the scaling of the variables and performs well even for ill-conditioned problems especially when the hes-sian matrix is nearly diagonal.

Contributed by the Design Automation Committee for publication in the J OUR-NAL OF M ECHANICAL D ESIGN.Manuscript received Sept.1999.Associate Techni-cal Editor:A.Diaz.

In section 2,the sub-problem generated at each iteration is de-scribed along with the solution strategy.The criterion of selecting the move limits for each sub-problem is described in section 3.The algorithm requires a feasible starting point.A simple method for ?nding such a point using moving logarithmic barriers is de-scribed in section 4.A few examples are described in section 5,to illustrate the application of this sequential approximation method.

2The Sub-Problem and Solution Strategy

Consider an optimization problem with linear constraints and side constraints as described below,

?P ?:Min f ?x ?,

x ?R n ,f :R n →R subject to,

Ax ?b ,

A ?R m ?n ,

b ?R m

(1)

l i рx i рu i ,

i ?1,...,n

The objective function f is a nonlinear function.The vectors l

and u set the side constraints on the variables,written succinctly as,l рx рu .A sequence of sub-problems is generated for the above optimization problem.At the k th iteration,the sub-problem is de?ned as:

?P k ?:Min F k ?x ?,such that,x ?S k

S k ??x ?Ax ?b ;l ??k ?рx рu ??k ??

(2)

F k ?x ???f ?x ?k ??t

x ??k

?

i ?1

n

ln ?x i ?l ?i ?k ?

???k

?

i ?1

n

ln ?u ?i ?k ?

?x i ?

(3)

The objective function of the sub-problem,F k (x ),is strictly

convex.Assuming that S k is a non-empty and bounded set,there exists a unique global minimum for the above sub-problem.Note that the side constraints for the sub-problem (P k ),are not the same as that of the original problem ?P ?.Instead,l ?(k )and u ?(k )are ?exible move limits that satisfy the condition,

l рl ??k ??x ?u ??k ?рu

(4)

In this paper we refer to the side constraints of the sub-problem

as move limits.These move limits serve to limit the step size taken at each iteration k .Criteria for selecting their value is given later.If the gradient vector c ??f (x (k ))is taken to be a constant,the above sub-problem can be solved as a linear program,using barrier methods to enforce the move limits.The logarithmic bar-rier function has been found to perform well for the primal-dual method ?12?for linear programming.The primal-dual method is an interior point method for linear programming that uses loga-rithmic barriers to impose side constraints.The dual problem for the corresponding linear program,min c t x ,x ?S k can be derived as,

?D k ?:max q ???

(5)

q ????in?mum 1рx рu

?c t x ??t ?Ax ?b ???z t x

???t

b where,??R m and

z t ??c t ??t A ??R n and

x ?i ?k ?

?

?

l i ?k ?

if z i ?0u i

?k ?otherwise

?

(6)

The optimality criteria for the sub-problem (P k )can be derived

as,

?x i ?l ?i ?k ???u ?i ?k ??x i ?z i ??k ?u ?i ?k ??l ?i ?k ??2x i ??0,

i ?1,...,n

Ax ?b z ?A t ??c ,

(7)

where,z ?R n ,c ??f (x k )?R n and ??R m .

The optimality criteria are a set of nonlinear simultaneous equa-tions.Due to the convexity of the sub-problem,we know that the solution of these equations will yield the unique global minimum of the sub-problem.These equations can be therefore solved using Newton’s method.In doing so we treat c ??f (x (k ))as a constant vector,thereby making a linear approximation of the objective https://www.360docs.net/doc/e7441169.html,ing the standard ?rst order Taylor series expansion,we get the following linear simultaneous equations for evaluating the Newton descent direction.

?x i ?k ??l ?i ?k ???u ?i ?k ??x i ?k ???z i ??Z i ?k ???x i ?2?k ?x i ??v i

?k ?

A ?x ?0(8)

?z ?A t ???0

where,

v i ?k ???x i ?k ??l ?i ?k ???u ?i ?k ??x i ?k ??z i ?k ???k ?u ?i ?k ??l ?i ?k ??2x i ?k ?

?

(9)

Z i ?k ???u ?i ?k ??l ?i ?k ??2x i ?k ??z i

?k ?

The above equations can be solved using a linear equation solver.However,since the hessian matrix of the sub-problem is diagonal and we can use this special structure of the equations to solve them more ef?ciently.The solution can be derived in the following form by rearranging Eqs.?8?and ?9?.

?????ADA t ??1AV

?z ?A t ??

(10)

?x ??V ?D ?z ,V ?R n and D ?R n ?n

where,

D ii ?

?x i ?k ??l ?i ?k ???u ?i ?k ??x i ?k ?

???u ?i ?k ??l ?i ?k ??2x i ?k ??z i ?k ?

??2?k

,

D i j ?0,i j ,i ,j ?1,...,n

V i ?

v i

??u ?i ?k ??l ?i ?k ??2x i ?k ??z i ?k ?

??2?k

,i ?1,...,n (11)

This formulation does not require the hessian of the sub-problem to be assembled and hence enables us to use Newton’s

method for the sub-problem even when the number of variables are very large.In Eq.?10?,note that we need to solve only ‘m ’simultaneous equations to calculate ??.Therefore,in applications where the number of variables is very large as compared to the number of constraints this leads to signi?cant reduction in com-putation.Dual methods ?13?also have the same advantage when applied to such problems since the dual problem would have fewer variables than the primary problem.However,unlike in the dual methods we do not completely solve the sub-problem in this algorithm even though the Newton iterations derived above could be used to obtain the global minimum of the sub-problem (P k ).Since the minimum of the sub-problem is not the same as the minimum of the original problem (P ),we need not solve the sub-problem exactly.Fleury ?14?has suggested similar ideas in the context of dual methods.We restrict the number of Newton iterations to the minimum required to generate a descent direction.For unconstrained optimization such a direction is found at the very ?rst iteration.However,for constrained optimization,if our initial guess for the Lagrange multipliers are not close to the real values,more than one iteration may be required.In this case,we use the Eqs.?10?–?11?to compute ??and ?z to update ?,z and ?k but leave x unchanged.In fact,we found no advantage in updating x more than once per iteration ?or sub-problem de?ni-tion ?.Once a descent direction is found,the variables are updated as follows,

x?k?1??x?k???k?x

z?k?1??z?k???k?z(12)

??k?1????k???k??

The step size?k is selected such that the new values of the variables are feasible,i.e.,l?(k)?x?u?(k).The step size may there-fore be computed as follows,

?k?min

i?1,...,n

???i?(13) where,

??i??l?i?k??x i?k??x i?k?if?x i?0

u?i?k??x i?k?

?x i?k?

otherwise?,i?1,...,n(14)

After the variables are updated,the sub-problem is rede?ned by re-evaluating the functions and gradients and by resetting the move limits.Note that only one function and gradient evaluation is required per iteration to de?ne the sub-problem.However,once the descent direction and step size are found it is bene?cial to use Armijo’s rule?12?to check whether further step size reduction is necessary.When Armijo’s rule is used to reduce step size,one additional function evaluation is required per step size reduction. The move limits set by the barrier function usually minimizes the need for step size reduction so that often Armijo’s rule serves merely as a safety check.The Armijo rule used in our implemen-tation is given below.We use?m k?k?instead of?k in Eq.?12??as the step size,where we chose a?xed scalar?such that0???1and m k is the smallest positive integer for which the follow-ing relation holds.

f?x?k???f?x?k???m k?k?x?k??у???m k?k?f?x?k??t?x?k?

(15) In our implementation,we have used??0.5and??0.01.The magnitude of?x from Eq.?10?depends inversely on?k.As the variable approaches the optimal value and the upper and lower move limits converge towards each other,?k should tend towards zero so that step size does not become too small.We have used the following rule to set the value of?k.

?k?g k

2n

,where g k?z t?x?x??k??(16)

z and x?(k)are de?ned in Eq.?6?.Note that g k is the duality gap between the linear program min c t x,x?S k and its dual max q(?),

de?ned in Eq.?5?.At the optimal point,the duality gap vanishes as expected because z i?0if x i?x?i(k) 0.At the optimal solution of the original problem?P?,the duality gap of sub-problem must

be zero,g k?0.We make use of this to de?ne a convergence criterion as follows:

C?x??g?x?

f*nnр?,where,g?x??z

t?x?x??and

x?i??l i if z i?0

u i otherwise?(17) For the examples in this paper,we have used??1?10?6. Notice that the Lagrange Multipliers are also obtained as by-products during the solution process if the side constraints do not become active.

The algorithm presented above can also handle inequality con-straints of the form Cxрd,C?R rxn,d?R r by using slack vari-ables to convert them to equality constraints of the form Cx?s ?d,s?R r,s iу0.We illustrate this through an example in sec-tion5,where we convert an inequality constraint into an equality constraint by using an addition variable x n?1as the slack variable, where n is the number of variable in the problem.The side con-straints for this variable are set such that the lower bound is zero and the upper bound is an arbitrary large value.Notice that even though the objective function of the original problem is indepen-dent of this slack variable,the objective function of the sub-problem would have barriers to ensure that the slack variable re-mains positive.Hence using slack variables is identical to using barrier method to enforce inequality constraint.

3Criteria for Setting Move Limits

The move limits are set for each sub-problem in such a way that the variables stay within the feasible domain of the original opti-mization problem?P?.The upper and lower limits are moved such that they converge towards each other and the feasible region for the sub-problems shrinks progressively.The move limits help to stabilize the algorithm and prevents the step size from being ex-cessive.In addition,the curvature of the sub-problem objective function depends on these move limits.The value of the move limits are initialized to be the same as the side constraints on the original optimization problem?P?,that is,l?i0?l i,and u?i0?u i. The move limits are then updated every iteration using one of the following criteria.

Criteria I:At the k th iteration,the move limits are updated as follows.

if k?1,l i?1??l i,u i?1??u i

if k?1,a i?

z t?xàx??

?n?z i?k??z i?k?1???(18) l?i?k??max?x i?k??a i,l i,x i?k???u i?x i?k???

u?i?k??min?x i?k??a i,u i,x i?k???x i?k??l i??,i?1,...,n(19) Note that the move limits on each variable are set such that the upper and lower limits are equidistant from the current value of the variable.When the move limits are equidistant,the descent direction generated every iteration is not dependent on?k?even though its magnitude does depend on?k?as can be veri?ed from Eqs.?8?–?11?.This update criterion decreases the distance be-tween the bounds,2a i,as the duality gap z t(xàx?)decreases.At the optimal solution,when the duality gap vanishes,the upper and lower move limits converge on to the optimal values of the vari-ables.The inverse relation between a i and z i(k)?z i(k?1)in Eq.?18?ensures that the sub-problem is properly re-scaled during the Newton iterations.If the solution is oscillating during the itera-tions,the value of z i will also oscillate between positive and nega-tive values making the denominator in Eq.?18?larger and there-fore the value of a i smaller.

Another criterion for updating move limits is listed below, which is a modi?ed form of the criterion used by Svanberg?10?for the method of moving asymptotes.The main difference here is that we require the moving limits to be within the side constraints of the original problem.

Criteria II:At the k th iteration,the move limits are updated as follows.

If kр2,

l?i?k??l i

u?i?k??u i

(20)

If k?2,

l?i?k??max?x i?k??s?x i?k?1??1i?k?1??,

l i,x i?k???u i?x i?k??}

u?i?k??min{x i?k??s?u i?k?1??x i?k?1??,

u i,x i?k???x i?k??l i?}

(21)

where,0?s?1if(x i(k)?x i(k?1))(x i(k?1)?x i(k?2))р0,else s ?1.

In our implementation we have used s ?3/4or 4/3depending on the sign of (x i (k )?x i (k ?1))(x i (k ?1)?x i (k ?2)).

4Finding an Initial Feasible Point

The initial starting point should be inside the feasible region so that it satis?es the constraints.During subsequent iterations,the variable stays within the feasible region due to the move limits set by the barrier function.An initial feasible point can be found using a simple modi?cation to the moving barrier technique.An arbitrary starting point x 0can be projected on the hyper-plane Ax ?b by minimizing ?x ?x 0?2subject to Ax ?b .The pro-jected point x ?0is obtained as

x ?

0?x 0?A t ?AA t ??1?Ax 0?b ?(22)

The analytical center of the polyhedra S ??x ?Ax ?b ;1?(k )?x ?u ?(k )?is de?ned as the unique minimum over S of the convex function F b de?ned below.

P b :Min F b ?x ?,x ?S ,

(23)

F b ?x ???

?

i ?1

n

ln ?x i ?l ?i ?k ?

??

?

i ?1

n

ln ?u ?i ?k ?

?x i ?

The analytical center is an interior point of the polyhedra S .If

the point x ?0does not satisfy the feasibility requirement 1?x ?0?u of the original problem ?P ?,then move limits of the problem

(P b )are set such that 1

?(k )?x ?0?u ?(k )using Eq.?21?.The method used in section 2to ?nd a descent direction for the sub-problems (P k )can also be used for the optimization problem (P b ).Indeed,after setting c ?0and ?k ?1,Eqs.?10?and ?11?yield a descent direction for (P b ).In practice,larger values for ?k lead to faster convergence towards the analytical center.After updating the value of x ?0using the descent direction,the move limits may be reset every k th iteration as:

l ?i ?k ??x ?0i ?k ???;u ?i ?k ?

?u i if x i рl i u ?i ?k ??x ?0i ?k ???;

l ?i ?k ??l i if x i уu i

u ?i ?k ?

?u i ;

l ?i ?k ??l i ;

if l i рx i рu i ,

i ?1,...,n (24)

In the above equations,?is a small positive real number.At each iteration,the value of x ?0is changed so that it moves away from the move limits along the hyper-plane Ax ?b .By resetting the move limits using Eq.?24?,the move limits are again moved closer to the variable x ?0.As the iterations continue,the move limits approach the actual bounds on the variables ?l and u ?.The iterations are stopped when the feasibility conditions l ?x ?0?u are satis?ed.

5Examples and Application in Structural Optimiza-tion

In this section a few examples are given to illustrate the appli-cation of the algorithm presented in this paper,which we will refer to as the Moving Barrier Method ?MBM ?.We have com-pared the convergence rate of the Moving Barrier Method with the Method of Moving Asymptotes ?MMA ?for all the examples.In our implementation of MMA algorithm,two criteria suggested by Svanberg ?10?are available for setting the position of the moving asymptotes.These are listed below.Moving asymptote criterion I:

L i ?k ??x i ?k ???u i ?l i ?

U i ?k ??x i ?k ?

??u i ?l i ?

for k р2(25)L i ?k ??x i ?k ??s ?x i ?k ?1??l i

?k ?1??U i ?k ??x i ?k ??s ?u i ?k ?1??x i ?k ?1?

?

for k ?2,(26)

where,i ?1,...,n and s ?3/4if (x i (k )?x i (k ?1))(x i (k ?1)?x i (k ?2))у0,else s ?4/3.

Moving asymptote criterion II:

L i ?k ??tx i ?k ?

,

U i ?k ??x i ?k ?

/t ,where t ?1/3.

(27)

Example 1.Unconstrained Optimization.In this example

we consider a simple unconstrained nonlinear program with two variables.The optimization problem may be stated as:

min f ?x 1,x 2??

?x 2?g 2??x 1g 1

1000.

,

(28)

where,g 1(x 1,x 2)?50000.?5000.x 1?40.x 2?x 1x 2?0.002x 22and

g 2(x 1,x 2)?100000.?2g 1

The initial bounds on the variables were set as 0.?x 1?15.,.0.?x 2?10000.and the initial points were given as (x 1,x 2)?(1.,10.).The contours of the function f (x 1,x 2)are plot-ted in Fig.1.The function has a very ill-conditioned hessian ma-trix.This optimization problem was solved using moving barrier method ?MBM ?,Newton’s method with modi?ed Cholesky fac-torization ?12?and the method of moving asymptotes ?MMA ?.The graph shown in Fig.2presents a comparison of the conver-gence history of the MBM with the other algorithms.As expected,the fastest convergence was obtained using the Newton’s method,which converged within four iterations with an objective function value of ?796.07.However,Newton’s method requires the hes-sian matrix and step size reduction at each iteration leading to approximately 35function evaluations.The MBM algorithm with Armijo’s step size reduction and criteria II ?Eqs.20–21?for set-ting move limits gives very similar rate of convergence.It con-verged to the same value of objective function within eight itera-tions and required 29function evaluations.Criteria I ?equation ?lead to slower convergence and more function evaluations ?

ap-

Fig.1Contours of the function f …x 1,x 2

…

Fig.2Convergence history for example 1

prox.42?.The method of moving asymptotes ?MMA ?using cri-teria I ?Eqs.25–26?took approximately 30iterations to converge to the same solution but since it uses only one function evaluation per iteration,the computational cost is similar to MBM.Notice that the objective function value does not decrease monotonically during the iterations for MMA https://www.360docs.net/doc/e7441169.html,ing criteria II ?Eq.27?,MMA does not converge within reasonable number of iterations.

Example 2.Cantilever Beam.This example is from Svan-berg ?10?but has been modi?ed by a change of variables to con-vert the non-linear inequality constraint into a linear inequality constraint.It involves weight minimization of a cantilever beam.As illustrated in the Fig.3,the beam consists of ?ve elements whose cross-section is square and hollow with constant thickness.The length of the side of the square cross-section for each element is treated as the design variable.The objective is to minimize the weight of the structure,subject to constraints on the vertical de-?ection at the end of the beam where the load is applied.The problem was stated in Svanberg ?10?as:

Minimize 0.0624?x 1?x 2?x 3?x 4?x 5?,

(29)

subject to:61/x 13?37/x 23?19/x 313?7/x 43?1/x 53

р1and x j ?0,

The non-linear inequality can be converted to linear inequality constraint by restating the problem in terms of a different variable that we de?ne as y i ?1/x i 3.In terms of y i the problem can be stated as:

Minimize 0.0624

??1y 1

?1/3

?

?1y 2

?1/3

?

?1y 3

?1/3

?

?1y 4

?1/3

?

?1y 5

?1/3

?

(30)

subject to 61y 1?37y 2?19y 3?7y 4?y 5р1

This problem was solved using both MBM and MMA algo-rithms.To use MBM,the inequality constraint was converted to an equality constraint by using a slack variable y 6increasing the total number of variables to six.In the Fig.4,the convergence history for MBM and MMA are presented.The problem was solved starting from the initial guess for the variables y i ??0.006,0.006,0.006,0.006,0.006?.The initial value for the slack variable was set as y 6?0.00001.The upper bound for all vari-ables were set as y i ?100,i ?1to 6.The convergence graph

shown in ?gure corresponds to criterion I for MBM ?Eqs.18,19?and criterion II ?Eq.27?for MMA.MBM converged to an objec-tive function value of 1.340within 6iterations requiring 13ob-jective function evaluations while MMA converged to the same value in 8iterations ?and 8function evaluations ?.In this example,MBM with criterion II ?Eqs.20–21?took more iterations to con-verge.MMA had very poor convergence rate using criterion I ?Eqs.25–26?for moving asymptotes.

Minimum Compliance Design.Structural optimization problems are nonlinear programs with a very large number of variables.Both function and gradient evaluations require ?nite element analysis of the structure and hence are computationally expensive.Here we consider an example of compliance minimi-zation where both the shape and topology of the structure are optimized.The structural optimization problem is stated as the minimization of compliance L (u )subject to a constraint on weight.The problem may be stated as,Minimize L (u )

L ?u ?????

?

?

f ?u ???d ??

?

?t

t ?u ???d ?(31)

subject to,

W ????

?

?

?d ?рW O

(32)

?

?0

?t ?u ˉ?D ?????u ?d ?0?L ?u ˉ?,

0р?р1

In the above equations,?(x )is the shape density function and

u ˉ

is the displacement ?eld.We de?ne the shape of the structure to be the region where shape density function has a value above a threshold value ?15,7?.Equation ?32?describes the constraint that the weight ?volume ?W of the component should be less than or equal to W O .The externally applied forces and traction are de-noted as f and t respectively.D is the elasticity matrix whose elements are functions of the material properties.Since our shape is represented using the density function ?,therefore,the material properties and the D matrix depend on the density function.We assume relations between Young’s Modulus and density of the form E ?E 0?n .This type of relation has been used before ?5,16,6?and the corresponding material has been referred to as Solid Isotropic Material with Penalization ?or SIMP ?.

A feasible region is de?ned within which the optimal shape must ?t and this region is divided into a quadrilateral mesh.The density function is represented using piece-wise bilinear interpo-lation over the quadrilateral elements of the mesh.The variables of the design are the density function values at the nodes and therefore the number of variables is equal to the number of nodes in the mesh.

Example 3:Shape and Topology Optimization.Figure 5illustrates a planar shape and topology optimization example.The feasible region,the applied loads,the support boundary conditions and the mesh used to represent the shape are displayed in Fig.5?a ?.The structure being designed is supported at the two lower corners and has to carry the three sets of loads shown in the ?gure.The initial geometry is assumed to be the rectangular

feasible

Fig.3Weight minimization of cantilever

beam

Fig.4Convergence history for example

2

Fig.5Optimal shape and topology design

region.Figure 5?b ?displays the optimal geometry computed under the constraint that the ?nal geometry must have 30percent less mass than the initial geometry.A cubic relation was used between Young’s modulus and density.In this example,the mesh con-sisted of 3200elements and 3333nodes.The number of variables of the optimization problem is equal to the number of nodes in the ?nite element mesh.

The convergence history is shown in Fig.6for MBM using criteria I and II ?labeled as MBM-I and MBM-II respectively ?and MMA.At the ?rst iteration,the value of the compliance is very low because the weight constraint is not yet satis?ed.Even though MBM-I is slow to converge during the early iterations at the end of 25iterations all three curves converge to a objective function value of approximately 4100.MBM algorithms took 15.3minutes to complete 25iterations while MMA took 20.7minutes to com-plete the same number of iterations on the same computer.The number of objective function evaluations required was equal to the number of iterations for all three methods since no step size reduction was required for MBM algorithm at any iteration in this example.For MMA only criterion II ?Eq.27?was used for mov-ing asymptotes because the algorithm did not converge fast using criterion I.

Example 4:Michell’s Frame.In Fig.7,we consider an ex-ample where 70percent of the mass is removed from a rectangu-lar feasible region loaded like a cantilever beam.When large per-centage of mass is removed from the initial feasible geometry,the resultant optimal geometry tends to be frame-like.Therefore,a relatively ?ne mesh is required to represent the geometry.This implies that the optimization problem has a very large number of variables.The optimal geometry obtained is shown superimposed on the mesh used to represent the geometry in Fig.7.The mesh used here contains 9801nodes and 9600elements.Here a qua-dratic Young’s modulus—density relation was used.The optimal solution for this problem was obtained analytically by Michell

?17?.The similarity between the solution obtained here and the analytical solution con?rms that the optimization algorithm did converge to the correct solution.

The convergence history for MBM and MMA for this example is presented in Fig.8.Again MBM-I and MBM-II stands for MBM using criterion I and II respectively,for setting move limits.The MBM algorithm required only one function and gradient evaluation per sub-problem for this example since no step size reductions were necessary.The convergence criterion ?Eq.17?was satis?ed after 35iterations for MBM-II and the objective function value was reduced to 4165.At the end of 35iterations,the objective function value was 4776for MBM-I and 4216for MMA.The MBM algorithm took 78.5minutes to complete 35iterations while MMA took 95.5minutes for the same number of iterations.Even though,in Fig.8,the objective function value appears to have almost converged by the sixth or seventh iteration for MMA and MBM-II,shape and topology change drastically over the next 30iterations.The shape and topology converge and become well de?ned ?with little or no gray areas ?only after about 35iterations.Even though MBM-II converges in 35iterations,the iterations where continued up to 50iterations as shown in Fig.8,to show that eventually all three methods converge to almost the same value for the objective function.

6Conclusions

A sequential approximation technique for nonlinear program-ming is described that locally approximates the objective function linearly and sets move limits on the variables using logarithmic barrier functions.This technique,like other sequential lineariza-tion techniques,does not require the evaluation of the hessian of the objective function.The algorithm can handle inequality con-straints also by converting them into equality constraints by using slack variables.The Lagrange Multipliers are obtained as a by-product of constrained minimization if the side constraints are not active.

Overall the performance of MBM is similar to that of MMA for problems involving only linear constraints.The modest improve-ment in speed over MMA observed for the minimum compliance design examples can be attributed to the fact that computation per iteration is reduced for MBM by not solving the sub-problems completely.Even though there was only one constraint for these examples,evaluation of the dual objective function and its gradi-ent for MMA becomes expensive as the number of variables in the primary problem increases.Furthermore,many iterations are required to completely solve the dual problem.For MBM only one Newton step is taken for the sub-problem if the very ?rst step yields a descent direction.Additional iterations for the sub-problem did not improve the rate of convergence.

The two criteria presented here for setting the move limits for MBM appear to work well for all examples in this paper.Criterion II requires storing the solution from the last two iterations as well as values of the move limits from the previous iteration,whereas criterion I requires storage of only the Kuhn-Tucker vector from the previous iteration.The rate of convergence of these two crite-ria varied from example to example making it dif?cult to

recom-Fig.6Convergence history for example

3

Fig.7Michell’s

frame

Fig.8Convergence history for example 4

mend one over the other.The performance of both MBM and MMA were found to be signi?cantly dependent on the strategy used for setting move limits or asymptotes.

The main limitation of the algorithm presented here is that it is not directly applicable to solve problems involving nonlinear con-straint.Many different approaches for extending the algorithm to handle nonlinear constraints are possible and need to be explored. We are currently studying an approach that involves linearizing the nonlinear constraints using Taylor series expansion at every iteration so that the sub-problem only involves linear constraints. Encouraging results were obtained for some examples but further study is required before the results can be published. References

?1?Schmit,L.A.,and Farshi,B.,1974,‘‘Some Approximation Concepts for Structural Synthesis,’’AIAA J.,12,No.5,pp.692–699.

?2?Fleury,C.,1979,‘‘Structural Weight Optimization by Dual Methods of Con-vex Programming,’’Int.J.Numer.Methods Eng.,14,pp.1761–1783.

?3?Fleury,C.,and Braibant,V.,1986,‘‘Structural Optimization:A New Dual Method Using Mixed Variables,’’Int.J.Numer.Methods Eng.,23,pp.409–428.

?4?Bendso”e,M.P.,and Kikuchi,N.,1988,‘‘Generating Optimal Topologies in Structural Design Using a Homogenization Method,’’Comput.Methods Appl.

Mech.Eng.,71,pp.197–224.

?5?Bendso”e,M.P.,1995,Optimization of Structural Topology,Shape and Mate-rial,Springer Verlag,Berlin.

?6?Kumar,A.V.,and Gossard,D.C.,1993,‘‘Geometric Modeling for Shape and

Topology Optimization,’’Fourth IFIP WG5.2,Geometric and Product Mod-eling,Wilson,P.R.,Wozny,M.J.,Pratt,M.J.,eds.

?7?Kumar,A.V.,and Gossard,D.C.,1996,‘‘Synthesis Of Optimal Shape And Topology of Structures,’’ASME J.Mech.Des.,118,No.1,pp.68–74.

?8?Zhang,W.H.,and Fleury,C.,1997,‘‘A Modi?cation of Convex Approxima-tion Methods for Structural Optimization,’’Comput.Struct.,64,Nos.1–4,pp.

89–95.

?9?Sigmund O.,1997,‘‘On the Design of Compliant Mechanisms Using Topol-ogy Optimization,’’Mech.Struct.Mach.,25,No.4,pp.493–524.

?10?Svanberg,K.,1987,‘‘The Method of Moving Asymptotes—A New Method for Structural Optimization,’’Int.J.Numer.Methods Eng.,24,pp.359–373.?11?Monteiro,R.D.C.,and Adler,I.,1989,‘‘Interior Path Following Primal-Dual Algorithms.Part I:Linear Programming,’’Math.Program.,44,pp.27–41.?12?Bertsekas,D.P.,1999,Nonlinear Programming,2nd ed.,Athena Scienti?c, MA.

?13?Fleury,C.,1993,‘‘Mathematical Programming Methods for Constrained Op-timization:Dual Methods,’’Structure Optimization:Status and Promise,Ka-mat,M.P.,ed.,series on Progress in Astronautics and Aeronautics,AIAA, Chap.7,pp.123–150.

?14?Fleury,C.,1982,‘‘Reconciliation of Mathematical Programming and Optimal-ity Criteria methods,’’Foundations of Structural Optimization,Morris,A.,ed., Chap.10,pp.363–404,Wiley,New York.

?15?Kumar,A.V.,1993,‘‘Shape and Topology Synthesis of Structures Using a Sequential Optimization Algorithm,’’Ph.D.thesis,Massachusetts Institute of Technology,Cambridge,MA.

?16?Zhou,M.,and Rozvany,G.I.N.,1991,‘‘The COC Algorithm,Part II:Topo-logical,Geometrical and Generalized Shape Optimization,’’Comput.Methods Appl.Mech.Eng.,89,pp.309–336.

?17?Michell,A.G.M.,1904,‘‘The limit of Economy of Material in Frame Struc-tures,’’Philos.Mag.,8,No.4.