Structural Optimization 5, 225-232 (1993)

StructuralOptimization © Springer-Verlag 1993

A m e t h o d for structural optimization which combines secondorder a p p r o x i m a t i o n s and dual techniques* T. L a r s s o n a n d M. R S n n q v i s t Department of Mathematics, LinkSping Institute of Technology, S-581 83 LinkSping, Sweden Abstract A structural optimization problem is usually solved iteratively as a sequence of approximate design problems. Traditionally, a variety of approximation concepts are used, but lately second-order approximation strategies have received most interest since high quality approximations can be obtained in this way. Furthermore, difficulties in choosing tuning parameters such as step-size restrictions may be avoided in these approaches. Methods that utilize second-order approximations can be divided into two groups; in the first, a Quadratic Programming (QP) subproblem including all available second-order information is stated, after which it is solved with a standard QP method, whereas the second approach uses only an approximate QP subproblem whose underlying structure can be efficiently exploited. In the latter case, only the diagonal terms of the second-order information are used, which makes it possible to adopt dual methods that require separability. An advantage of the first group of methods is that all available second-order information is used when stating the approximate problem, but a disadvantage is that a rather difficult QP subproblem must be solved in each iteration. The second group of methods benefits from the possibility of using efficient dual methods, but lacks in not using all available information. In this paper, we propose an efficient approach to solve the QP problems, based on the creation of a sequence of fully separable subproblems, each of which is efficiently solvable by dual methods. This procedure makes it possible to combine the advantages of each of the two former approaches. The numerical results show that the proposed solution procedure is a valid approach to solve the QP subproblems arising in second-order approximation schemes.

1

Introduction

Structural optimization is concerned with finding the best possible structure according to a given objective while simultaneously satisfying certain design restrictions. A model of a structural optimization problem consists of two related optimization problems (see, e.g. Kitsch 1981): the design problem that determines sizes, materials and shapes, and the analysis problem that determines the response of the structure when affected by loads. A standard solution concept is to use iterative methods where the design and analysis problems are approximated and solved alternately; this is known as the nested approach and works as follows. Given a design x(k), an approximate analysis problem is stated. Then, by using the well-known Finite Element Method (FEM) this problem is restated as an equation system known as the state (or equilibrium) equations. If linear elastic behaviour is assumed, it *Presented at NATO ASI "Optimization of Large Structural Systems", Berchtesgaden, Gerraany, Sept. 23 - Oct. 4, 1991

becomes a set of linear equations, g [ x ( k ) ] u = p[x(k)],

(1)

where K[x( k)] is the structural stiffness matrix, u is a vector of node deflections and p[x(k)] is a vector of node forces. (The vector p[x(k)] is usually assumed constant, p.) If the assumption of linear elasticity is not valid, a system of nonlinear equations is obtained. In order to find an improved design, an approximate design problem expressed explicitly in the design variables is formulated. To achieve this, the objective and constraint functions are approximated by means of Taylor series expansions. However, this is computationally very expensive because the derivatives needed are implicitly determined via the analysis variables (i.e. node deflections) through the state equations (1). In order to reduce the computational effort, it is therefore of the highest necessity to decrease the number of iterations. A key element to accomplish this is the construction of accurate approximations of the design problem. To state the approximate design problem, a first- or second-order approximation strategy is used. First-order structural optimization methods are essentiMly based on the generM sequential linear programming approach developed by Griffith and Stewart (1961), and known as the Method of Approximation Programming (MAP). In order to prevent the solution of the linear subproblem from being very remote from the current iteration point, which might easily happen since the solution to a linear problem is attained at a vertex of the feasible region, some step-size restrictions (or move limits) must be imposed. These restrictions are bounds on how much each variable value may alter in each iteration. They are crucial since a proper selection may give good convergence behaviour, whereas a bad selection may cause slow convergence. In structural optimization, when the design variables are, for instance, cross-sectional areas of bar-truss structures, it has been found very advantageous to use reciprocals of the design variables when stating the first-order approximate design problem, such reciprocal variables generally give superior linear approximations of stress and displacement responses (see, e.g. Fleury and Braibant 1985). A hybrid approach, where both direct and reciprocal variables are used, as proposed by Starnes and Haftka (1979), can also be used. Fleury (1979), and Fleury and Schmit (1980) introduced the combination of approximation by reciprocal variables and dual methods; this approach is referred to as the CONvex" LINearization (CONLIN) method. CONLIN has been fur-

226

ther developed in a number of papers (e.g. Fleury 1989b). The same concept of combining approximations and dual methods has also been used by Svanberg (1982, 1987). Svanberg (1987) introduced the Method of Moving Asymptotes 1 (MMA). Here, the intermediate variables ~ or where Uj and Lj are given values called asymptotes, are used instead of the reciprocal and direct variables. These intermediate variables can be viewed as a way to avoid trust regions as step-size restrictions. It should be noted that although non-linearities are introduced, both CONLIN and MMA are first-order approximations. Moreover, both have the difficulty how to choose proper move limits (in MMA the asymptotes Uj and L j). There are several fundamental reasons underlying the efficiency obtained by combining approximation concepts and dual methods. Two of them are that dual methods can exploit the special algebraic structure of the approximate problem generated at each iteration, and that only simple bounds restrict the dual variables. We will now turn to second-order approximations; to do this we consider the following general problem: [P]

mxin f ( x ) ,

s.t.(g(x) < 0 ~x E X. Here, g(x) = [gl(x),...,gm(x)] T, and the set X consists of simple bounds on the variables. Most optimization algorithms set forth to solve [P] generate a sequence of points, x(k), tending to an optimal solution, starting from an initial solution, x(0). First, an approximate problem, where /(x) and g(x) are replaced with approximations ](x) and ~(x) derived at x(k), is stated and solved, giving ~(k) and a search direction d(k) = ~(k) _ x(k). Thereafter, to obtain the new iteration point, a line-search on a merit (or descent) function ~ is performed, i.e. let

a (k) e arg min ff[x (k) + a d (k)] ~>_0

(2)

where a is a step-size, arid x(k+l) = x(k) + a(k)d(k). The general form of the merit function is ~(x, r) = f(x) + V(x, r), where V(x, r) is some penalty function that measures the infeasibility of the point x, and r is a vector of penalty parameters which controls the degree by which the objective function is penalized at infeasible points. The purpose of the descent function is to enforce steady progress to a locally optimal solution by balancing the conflicting aims of reducing the objective function value and to satisfy the nonlinear constraints. It is crucial to choose an approximation such that a descent direction is obtained for the merit function, i.e. so that

ff[x (k+l)] < ff[x(k)],

k = 0, 1,2 .... ,

(3)

can always be fulfilled. Unfortunately, it is highly impractical in reality to perform line-searches on general structural optimization problems, because of the implicit dependence between the design and analysis variables. A simplified approach is therefore to select the point from each approximate design problem as the new iteration point. In this procedure, however, it is not possible to ensure convergence. Straightforward approximations based on second-order Taylor series expansions give nonlinearities both in the objective and in the constraints. This will create a rather difficult subproblem for which efficient solution procedures are hard to

develop. It is, however, possible to avoid introducing nonlinearities in the constraints. The most promising and efficient approaches along this line are based on solving a sequence of QP subproblems; these methods are known as Sequential Quadratic Programming (SQP) or Recursive Quadratic Programming (RQP) methods. They use approximate or exact Hessian matrices to state a QP subproblem at each iteration, and can be interpreted as an application of Newton's method to the problem of finding a stationary point to the Lagrangian function m L(x, )~) = f(x) + ~ Aigi(x). (4) i=1 Historically, SQP methods start with Wilson (1963), who suggested a method called SOLVER, see also van de Panne (1975). This method involves a sequence of quadratic subproblems of the following type: [QP]

rn~n f[x(k)] + Vf[x(k)]T[x - x(k)] +

+~[x -

x(k)]TH[x(k) , ,k(k) ] [x - x(k)],

s.t. [ g[x(k)] + Vg[x(k)]T[x- x(k)]

(

x

<

0,

E

X

Here, H[x(k), ),(k)] = vTVL[x(k),),(k)] is the Hessian matrix of the Lagrangian function. The solution of each quadratic programming subproblem is used for the approximation in the following iteration, and the Lagrangian multipliers associated with the linearized constraints are used to update the multipliers associated with the nonlinear constraints. Later, Han (1976, 1977) developed a quasi-Newton SQP method where the second-order derivatives were approximated using well-known quasi-Newton updating formulae. An improvement has been suggested by Schittkowski (1983), who suggested the use of an active (constraint) set strategy to avoid unnecessary calculations of gradients. An active set strategy will be preferred because it avoids the expensive calculation of the gradients of all the implicit constraint functions. In Han's algorithm no active set strategy is used; to prove global convergence of the algorithm, it is necessary to include approximations of all constraints when defining [QP]. An active set strategy is also used in Pshenichnyj's algorithm (Pshenichnyj 1970, 1987), for which global convergence has been proven. This algorithm has been successfully used for structural optimization problems (see Belegundu and Arora 1984; Choi et al. 1983). Pshenichnyj, however, used H = I, the identity matrix, thus ignoring the possibility of using second-order information. This can also been seen as a way of restricting the step-size by adding the Euclidean norm as a penalty term to the objective function; it can therefore be viewed as a first-order approximation scheme. Merits of Pshenichnyj's approach are that it preserves separabilitY and that it makes the objective function strictly convex; hence, efficient dual methods can be used. Renwei and Peng (1987) suggested a dual method where only the diagonal terms of the second-order information are used. Lim and Arora (1986) improved Pshenichnyj's algorithm by using second-order information in a quasi-Newton manner. The resulting QP subproblems are solved with a standard QP subroutine.

227 One difficulty of applying SQP methods to structural optimization problems is that straightforward calculations of the second-order derivatives are computationally expensive. However, Fleury and Smaoui (1988) and Ringertz (1988) have demonstrated that the Hessian of the Lagrangian function can be calculated using only one additional load case in the sensitivity calculation, thereby supplying an important tool for the use of SQP schemes in structural optimization. Fleury (1988) approximated the Hessian matrix with diagonal elements, which results in the separability of the objective function being kept; hence it is advantageous to use dual methods. This scheme is similar to that suggested by Renwei and Pent (1987). The use of diagonal elements only can be viewed as a way to impose properly chosen move limits, and Fleury (1989a) proposed an extension to MMA where the asymptotes are calculated utilizing the second-order diagonal terms. Ringertz (1988) used the full Hessian and solved the resulting nonseparable problem by a standard quadratic programming routine. A dual approach is not useful in this case since the Lagrangian subproblem is not separable. The contents of this paper are as follows. First, in Section 2, we state the design problem and describe the overall solution procedure, which is essentially equivalent to SOLVER. Further, we show how the second-order approximation of the design problem can be obtained; this QP subproblem is in the general case nonseparable. In Section 3, we describe our proposed solution procedure for the QP subproblem, based on a partial linearization of the objective function. This procedure gives rise to a series of fully separable subproblems, which can be solved with dual methods. Thereafter, in Section 4, we summarize the entire solution scheme. In the following section we give numerical results, and finally, in Section 6 conclusions will be drawn concerning the efficiency of the proposed algorithm. 2

Problem formulation and SQP scheme

To formulate the overall design problem we will use either direct or reciprocal variables. The notations used are given in Table 1. The design problem can then be stated as mxinw(x),

[gi(x)=c~i(x)-~i s.t.~gi(x)=ui(x )-gi

xj

t

< -<

0 ,i=l,...,n, 0 ,i=n+l,...,m,

e

Xj,j= l,...,n .

The objective function is stated as n

w(x) =

ej

j=l if reciprocal variables are used, and as n

w(x) =

i ,

zj

the cross-sectional area (or reciprocal area) of element j of the structure; j = 1, ..., n lower bound of the area (or reciprocal area) of _xj element j; j = 1, . : . , n upper b o u n d o f the area (or reciprocal area) of element j; j = 1, ..., n w (x) weight of the structure cj constant for the j-th element's contribution to the total weight, i.e. length times density; j = 1, ..., n (ri(x) functions describing stresses; i = 1, ..., n ui(x ) functions describing deflections; i = n + 1, ..., m limit on in for constraint i; i = 1, ..., n ui limit on deflection in constraint i; i = n + 1, ..., rn xj {z.ilz_ i < x i _<~,:} ; j = 1, . . . I n

X

(6)

j=l in the case of direct variables. The inequalities ai(x ) < ~i,i = 1 , . . . , n and ui(x) < ui, i - n + 1 , . . . , m are behavioural constraints that impose limitations on structural

X=N~=IXj

response quantities. The functions c~i(x) and ui(x ) are usually non-linear and implicit with respect to the design variables. In order to evaluate these functions, it is necessary to solve an auxiliary analysis problem that determines the behaviour of the structure by minimizing its potential energy. In our formulation, the restrictions on structural response consist of stress and deflection limits, but it is of course also possible to introduce restrictions on eigenfrequencies, buckling, etc. Furthermore, when the deflections (or displacements) u(x) are known, it is possible to determine the stresses. Assuming linear elastic behaviour, the stresses can be written as gi(x) -- qTu(x), where qi are constant vectors. Since we consider only deflection and stress constraints, the functions gi(x) can all be written as gi(x) = q T u ( x ) - Ki, i = 1 , . . . , n and gi(x) = q T u ( x ) - u i , i = n + 1 , . . . , m . (In the latter case, the vectors qi contain only one non-zero element, whose value is one.)

2.2

2.1 Problemformulation

[DP]

Table 1. Notation

The SQP scheme

In order to solve problem [DP], an SQP scheme, which iteratively solves a sequence of quadratic approximations, is applied. Given a design x(k) and a set of Lagrangian multipliers ~(k), an analysis is performed in order to obtain the response of this design. Thereafter, a quadratic approximation is made, which results in a subproblem such as the problem [QP]. In order to obtain an improved design, this QP subproblem is solved giving the primal solution g and the corresponding dual solution A. In the next iteration, these solutions are used to state the next QP subproblem, i.e. x(k+l) = K and ~(k+l) = ~. This iterative scheme continues until some convergence criteria are fulfilled, e.g. conditions on weight reduction and feasibility. To state the problem [QP], it is necessary to evaluate first- and second-order derivatives of the functions w(x) and g(x) with respect to the design variables at the current point x(k). The gradient Vw(x) is easily calculated since w(x) is explicitly given, whereas Vg(x) is implicitly determined through the analysis variables. For one individual constraint gi(x), Vgi(x ) can be stated as

Vgi(x)T = vhT"(x)]r = qTdU' dx'

(7)

228 The Jacobian of u, du, can be obtained by first differentiatax ing the state equations (1) with respect to any x j, giving

introduction of the term 7 I can also be interpreted as a way to introduce move limits, since

OK K0U Oxj u q- Oxj = 0,

~[xl _ x(k)]T(H + ' y I ) [ x - x (k)] =

(8)

if p is asssumed independent of zj. Then - ~ is determined by solving (8) for each independent design variable. This is known as the pseudo-load technique (see, e.g. Adelman and Haftka 1986), and involves solving several systems of equations of the type K v = s i. Performing a matrix inversion of K is computationally inefficient; instead, a Cholesky decomposition of the stiffness matrix is usually done, which decomposes K as K = LTL. Then only forward and backward substitutions must be performed each time a system of equations must be solved. The tensor dK is in general easily calculated, and if, for instance, bar-truss structures are considered it is computed analytically. The Hessian matrix can be stated as H(x, X) = VTVw(x) + x T V T V g ( x ) . (9) Each individual element of H can be computed (see Ringertz 1988) via

O2w

~ i q T 02" m

Hjk - OxjOx~k +

i=1

OzjOxk --

02w T 02u - OxjOx~k -t- qtot OxjOXk • However, the second derivatives ~

(10)

02u

are unknown. Dif-

ferentiating (8) yields 02K 0K 0u OK 0u + K 02u OxjTz-£xku + OxjOXk + Ox~Oz~ O-~jOxk - 0 .

(11)

Introducing a virtual displacement vector Vtot, (10) may be written as 0u) 02w _ v~t (\ ~ 02K u + 0 Z -0u- + -OK Hjk OzjOx k Oxj Ox k Ox k ~ ' -

-

(12) where Vtot is obtained by solving

K v = qtot. In the following we will use the notation [DQP] s.t. { A5

m~n

(13)

f(5)=~sTRS-rTS,

_
for the approximated problem. Here, A = V g [x(k)], 1% = H [x(k), )~(k)], r = V w[x(k)], b = - - g [ x ( k ) ] + V g [x(k)]Tx (k), 5 = x - x(k) , and the set A is defined as A = {5 Ix_j - x ~.k) = ~_j < 5j ~_ ~j = "~j x ~.k), j = 1, .. .,n}. If positive definiteness of the matrix R is required, it is possible to ensure this by adding a positive definite matrix, e.g. by redefining i = n + ~I, (14) where I is the identity matrix and the positive scalar 7 is selected large enough. (It should be noted that, for a general problem, a proper selection of ~, is not an easy task.) The

i = ~[x-x(k)]TH[x-x(k)]-t-

i 12 , ~7llx-x(k)l

(15)

where the second term introduces an implicit step-size restriction. This move limit technique has earlier been used by Pshenichnyj (1970, 1987). A study of different realizations and interpretations of move limits can be found in the paper by Jonsson and Larsson (1990). 3

3.1

Solution p r o c e d u r e for t h e Q P s u b p r o b l e m

The partial linearization procedure

It is obvious that the objective function f(6) of problem [DQP] is not separable, and therefore, the merit of applying dual methods to solve it is lost. Many general non-linear programming algorithms can, of course, be used for obtaining a local optimal solution to the problem [DQP]. However, the principal drawback of these methods is that often they do not capture the inherent structure of the problem. One iterative method that takes the structure into account is given by Larsson and Migdalas (1990). It was developed for a nonseparable, non-linear programming problem over a Cartesian product set, where the objective function is assumed to be partially separable. The main concept of the method is to induce separability by performing an approximation of the non-separable part of the objective function using a firstorder Taylor series expansion, i.e. by a partial linearization of the objective function. The method can be viewed as a generalization of the Frank-Wolfe method (Frank and Wolfe 1956), since only a part of the objective is linearized. The non-separable part of the objective function f(6) of problem [DQP] is 1 5 T R & However, it is possible to rewrite I t according to 1%= R + D where D = diag(R), so that the truly non-separable part can be identified as h(5) = 15TRS. Suppose that the current feasible iteration point is 5(0, then h(5) is approximated by a first-order Taylor series expansion, which yields h(5)~25(0TRS(0- + [RS(/)]T[5 - 5(/)].

(16)

Using this approximation of the non-separable part of f(5) and dropping the constant term, we will obtain a separable subproblem [DQPS] s.t.

AS ~

m~n f(5) = 15TD5 + [RS(e)
r]r~,

.

Assume that the matrix t t is positive definite and that there is a feasible solution to the problem [DQP], and let ~(~) be the unique solution to the problem [DQPS] and d(~) = 5(~) -5(~). Then the two following propositions can be shown (see Larsson and Migdalas 1990). P r o p o s i t i o n 1. If d(£) = 0 then 5(£) is an optimal solution to the problem [DQP], otherwise d(l) is a feasible direction of descent.

229 Proposition 1 states that unless the optimal solution has been found, that d ( 0 is a feasible direction of descent, which makes it possible to perform a line-search from the current iteration point in the direction d(~) with respect to the original objective function in problem [DQP]. Moreover, this line-search is analytically solvable. A new point 6(t+l) is constructed according to ~f(/+1) = 8(g) + a(g)d(/), where the step-length c~(g) = arg mi>~{f(6(g) + ad(g))[6(g) + ad(g) is feasible in [DQP]}. Following this procedure, the second proposition, concerning the convergence behaviour, can be stated as follows. P r o p o s i t i o n 2. The algorithm either terminates with an optimal solution in a finite number of iterations or it generates an infinite sequence {6(~)} such that any accumulation point is an optimal solution to [DQP].

If R is not positive definite, although it holds for D, then only a local optimum can be guaranteed to be found. Since D is a diagonal matrix it is easy to ensure positive definiteness as described in the previous section. However, a selection of a too small parameter 7 in (14) will give rise to an illconditioned dual objective function. 3.2

A dual solution procedure

Instead of solving the problem [DQPS] with a primal method, the Lagrangian dual problem can be solved. This dual problem is defined as [LD]

~n>_a~001),

where @(N) is the dual objective function defined through the Lagrangian subproblem as [LDS]

O(n)= ~i~L(~,,)

= 16TDS-[r_gS(0]rs+

It is well-known that the dual function is concave. As L(6, ~1) is strictly convex and A is a convex set, 6(r/) is uniquely determined for any vector r/ > 0, and the dual function is differentiable so that it is possible to apply a gradient related search method for the solution of [LD]. The gradient of O(Tt) is

VO(~t) = A6(~t) - b . (19) The dual problem is now to maximize a 'differentiable and concave function subject to only non-negativity restrictions on the dual variables. Here a slightly modified steepest ascent method to take into account the non-negativity restrictions can be used. The optimal primal solution is obtained by inserting the optimal ~ values into (17). 4

The proposed algorithm

In this section, we summarize the proposed algorithm. The algorithm starts by selecting an initial design x(0) and initializing the Lagrangian multipliers to )~(0). In an arbitrary iteration, given the solution x(k), a structural analysis is first performed by solving the state equations (1). Then a second-order approximation of the design problem is made. We solve the resulting QP subproblem by applying the partial linearization concept, which creates a sequence of separable subproblems, each solved by utilizing the dual method. Thereafter, the new values of the design variables and Lagrangian multipliers are taken as the primal and dual optimal solutions, respectively. The process restarts according to the SQP scheme, which continues until the convergence criteria are fulfilled. The entire procedure can be formalized in the following way. • Initialization -

+~/T(A6 - b ) .

-

Here, L(6,~t) is the Lagrangian function to the problem [DQPS] and ~? are the Lagrangian multipliers associated with the.constraints A6 _< b; the set A is kept as explicit constraints. Under certair~ convexity assumptions, which hold in our case, it can be shown that solving the dual problem is equivalent to solving the primal problem (see, e.g. Bazaraa and Shetty 1979). The problem [LDS] separates into n onedimensional minimization problems of the type 1 2 [LDS/] min )~'(hj, ~?) --- -~Djj6j - -rj(~l)~j,

-

Initiafize iteration counter, k := 0. Let x(0) constitute an initial design. Initialize the Lagrangian multipliers, A(0) := 0.

• Major iteration k

Perform a structural analysis at x(k). - Formulate an approximate design problem [DQP]. -

-

Initialize subiteration counter, £ := 0, and set 8(0) := 0.

- Parlial linearizalion iteration g

* Approximate problem [DQP] around 6(0 giving [DQPS]. , Solve the dual problem [LD], which gives

[~(~1,~(~1].

where r(~t) = r - R6(g) - N T A . These are analytically solvable giving 6(~t), where each 6j (,) is defined as

* Make a line-search in [DQP] from 6(i) in the

5j(rl) = l S j if ~

direction 6(g) - 6(0, giving 6(g+l). * Check convergence criterion for the problem [DQP], if not satisfied update £ := £ + 1 and repeat.

<0,

(17)

- End partial linearization

Furthermore, the value of the dual objective function can now be stated explicitly according to 0(~/) = 2~07)TD607) + [~6(e) _ r]T6(rl) + + , T [ A ~ ( ~ ) - b].

(18)

, Let the new design x(k+l) := x(k) + 5(g+l). , Let )~(k+l) := ~(g). - Check convergence criteria for the problem [DP], if not satisfied update k := k + 1 and repeat. • Choose x(k+l) as the optimal design. • End

230

In each partial linearization iteration we choose the optimal dual solution from the previous iteration as the starting point. This selection has shown to be very efficient and only a few dual line-searches are needed to find the new optimal solution. It is easily shown that if the partial linearization procedure has converged, e.g. d(l) = 0, then ~(t) equals the optimal Lagrangian multipliers associated with the problem [DQP]. This follows from the fact that the first-order optimality conditions for [DQP] and [DQPS] coincide when ~(/) = ,(~). Therefore, the choice A(k+l) := ~(*) is valid. 5

Computational

0

results

0

In order to prove the feasibility of the suggested procedure and to draw preliminary conclusions about its efficiency, numerical results are presented. The implementation was done in FORTRAN 77 and executed on a SUN 4/390. Numerical results are presented for the well-known ten-bar truss problem given in Fig. 1, and a 56-bar truss structure given in Fig. 5. Material properties, stress limits, displacement limits, and loading conditions are given in Table 2.

L

.~¢

Iterations 15

20

25

30

Fig. 2. Weight histories using direct and reciprocal variables, respectively 1.18

1.16A 1.14!i 1.12 1.08

107 0.1 25000 2 100 360

1.06 1.04 1.02

l.oo

L

",~ :~ - 7":

0.98

We have tested formulations based on direct and reciprocal design variables, respectively. The overall convergence of the structural weight for the ten-bar truss is presented in Fig. 2 (dashed line: direct variables, solid line: reciprocal variables); the corresponding convergence of the maximal normalized constraint violation is presented in Fig. 3. The feasibility behaviour using reciprocal variables formulation is stable, which is due to the fact that these variables lead to accurate approximations of the feasible set. When using direct variables, the feasible set is not approximated as accurately as with reciprocal variables. Therefore, the intermediate designs may become quite infeasible with respect to the original constraint set. J

10

1.10

Table 2. Problem data Young's modulus (psi) Density (lbm/in 3) Stress limit (psi) Displacement limit (in) Load, P (kip) Length, L (in)

5

.¢

/ / / / / / / /

Fig. 1. Ten-bar truss The performance of the proposed partial linearization procedure for solving the QP subproblems is illustrated in Fig. 4, where the convergence behaviour from three different main

.

5

........ .

.

l0

~------- . . . . . .

.

Iterations 15

20

25

30

Fig. 3. Constraint feasibility histories iterations is illustrated. To be able to compare the convergence, the optimal [DQP] objective values have been normalized. The solid curve corresponds to an early main iteration and the dashed and dotted curves to later ones. (Here, the direct variables were used; however, the principal bellaviour is the same for reciprocal variables.) Clearly, only a few partial linearization iterations are needed to obtain a near-optimM solution. To study the performance of the proposed procedure for a larger problem, we considered the 56-bar truss problem. As above, the behaviour of the objective value convergence for three different main iterations is presented. Also here, few partial linearizations are needed to obtain a near-optimal solution. 6

Conclusions

The proposed solution procedure constitutes an efficient approach to the solution of the QP subproblems arising in SQP schemes applied to the optimization of bar-truss structures. Instead of solving the QP problem directly, a sequence of separable subproblems is solved utilizing efficient dual methods. In our application only a few iterations are needed to obtain convergence to an optimal QP solution. An explanation to the few iterations is that the Hessian matrix is rather diagonally dominant for the bar-truss structures of our test examples.

231 throughout the solution process. An important feature is that there are no difficulties in choosing proper tuning parameters such as move limits; this might be advantageous in real-world applications where there is no a priori information about the characteristics of the structural optimization problem under consideration.

1.12 o

1.10

~, 1.08 ~

1.06

~

1.02

Acknowledgements Z

1 . ~

0.98 0

.....

. 1

The research leading to this paper was partially supported by the Center for Industrial Information Technology (CENIIT). We thank Prof. Anders Klarbring for valuable suggestions to an earlier version of the paper.

°"

. . . . . . . 2 3 4 5 6 7 Iterations

10

8

References

Fig. "4. Number of partial linearization iterations

Adelman, H.M.; Haftka, R.T. 1986: Sensitivity analysis of discrete structural systems. A I A A J. 24, 823-832 Bazaraa, M. S.; Shetty, C. M. 1979: Nonlinear programming. New York: John Wiley ,~ Sons

\

Belegundu, A.D.; Arora, J.S. 1984: A recursive quadratic programming method with active set strategy for optimal design. Int. Y. Num. Meth. Eng. 20, 803-816 Choi, K.K.; Haug, E.J.; Hou, J.W.; Sohoni, V.N. 1983: Pshenichnyj's linearization method for mechanical system optimization. J. Mech. Trans. Auto. Des. 105, 97-103 Fleury, C. 1979: Structural weight optimization by dual methods of convex programming. Int. J. Num. Meth. Eng. 14, 1761-1783 P

P

Fleury, C. 1989a: First and second order approximation strategies in structural optimization. Struct. Optim. 1, 3-10

Fig. 5.56-bar truss

Fleury, C. 1989b: CONLIN: an efficient dual optimizer based on convex approximation concepts. Struct. Optim. 1, 81-89

1.12,

'• •

1"1°t

Fleury, C.; Bralbant, V. 1985: Structural optimization: a new dual method using mixed variables. Int. 3". Num. Meth. Eng. 23, 409-428

1.08

•

Fleury, C.; Schmit, L.A. 1980: Primal and dual methods in structural optimization. J. Struct. Div. 106, 1117-1133

~"~ 1.06 "~ 1.04

Fleury, C.; Smaoui, H. 1988: Convex approximation strategies in structural optimization. In: Eschenauer, H.A.; Thierauf, G. (eds.) Discretization methods and structural optimization - procedures and applications, pp. 118-126. Berlin, Heidelberg, New York: Springer

. ....q

1.02 Z

1.00 0.98 1

i

i

2

3

i

i

i

4 5 6 Iterations

i

i

i

7

8

9

10

Fig. 6. Number of partial linearization iterations We have- compared the performance of the SQP approach with that of first-order approximation schemes as well as Pshenichnyj's method. For these methods, it is possible to obtain faster weight convergence if suitable move limits are imposed, but steady behaviour with respect to feasibility cannot be achieved at the same time. The fact that firstorder approximations might reduce the number of iterations compared to second-order approximations has indeed been demonstrated by Woo (1986). The overall SQP scheme is conservative compared to first-order procedures with respect to constraint feasibility. Therefore, it would be advantageous for problems where it is important to maintain feasibility

Frank, M.; Wolfe, P. 1956: An algorithm for quadratic programming. Naval Research Logistics Quarterly 3, 95-110 Griffith, R.E.; Stewart, R.A. 1961: A non-linear programming technique for the optimization of continuous processing systems. Management Science 7, 379-392 Han, S.P. 1976: Superlinearly convergent variable metric algorithms for general nonlinear programming problems. Mathematical Programming 11,263-282 Han~ S.P. 1977: A globally convergent method for nonlinear programming. J. Optimiz. Theory Appl. 22, 297-309 Kirseh, U. 1981: Optimum structural design. Graw-ttfll

New York: Mc

Jonsson, ().; Larsson, T. 1990: A note on step-size restrictions in approximation procedures for strt~ctural optimization. Comp. ~A Struct. 37, 259-263

232 Larsson, T.; Migdalas, A. 1990: An algorithm for nonlinear programs over Cartesian product sets. Optimization 21, 535-542

approximation concepts and dual methods. AIAA J. 18, 12521260

Lim, O.K.; Arora, J.S. 1986: An active set RQP algorithm for engineering design optimization. Comp. Meth. Appl. Mech. Eng. 57, 51-65

Starnes, J.H.; Haftka, R.T. 1979: Preliminary design of composite wings for bucklings, strength and displacement constraints. J. Aircraft 16, 564-570

Pshenichnyj, B.N. 1970: Algorithms for the general problem of mathematical programming. Kibernetica 5, 120-125 Pshenichnyj, B.N. 1987: The linearization method. Optimization 18, 179-196 Renwei, X.; Peng, L. 1987: Structural optimization based on second-order approximations of functions and dual theory. Comp. Meth. Appl. Mech. Eng. 65, 101-114 Ringertz, U.T. 1988: A mathematical programming approach to structural optimization. Doctoral Thesis, Report No. 88-24, Royal Institute of Technology, Stockholm Schittkowski, K. 1983: On the convergence of a sequential quadratic programming method. Optimization 14, 197-216 Schmit, L.A.; Fleury, C. 1980: Structural synthesis by combining

Received Dec. 2, 1991 Revised manuscript received Jan. 22, 1993

Svanberg, K. 1982: An algorithm for optimum structural design using duality. Mathematical Programming Study 20, 161-177 Svanberg, K. 1987: The method of moving asymptotes - a new method for structural optimization. Int. J. Num. Meth. Eng. 24, 359-373 van de Panne, C. 1975: Methods for linear and quadratic programming. Amsterdam: North-Holland Wilson, R.B. 1963: A simplicial method]or concave programming. Ph.D. Dissertation, Harvard University, Cambridge, Mass. Woo, T.H. 1986: Space frame optimization subject to frequency constraints. Proc. AIAA/ASME/ASCE/AHS 27th Structures, Structural Dynamics, and Materials Conf., pp. 103-115