- Email: [email protected]
Optimum design of statically indeterminate beams under multiple loads
p
OPTIMUM
~~.?~9~8? 1987 Pcrgamon
53.00 + 0.00
Joumalr
Ltd.
DESIGN OF STATICALLY INDETERMINATE BEAMS UNDER MULTIPLE LOADS B. L.
KARIHALOO and
S. KANAGASUNDARAM
School of Civii and Mining Engineering, The University of Sydney, Sydney, N.S.W. 2006, Australia
(Received 7 ~~e~~er
1986)
paper presents minimum-weight designs of statically indeterminate beams under multiple loads. Under each load system consisting of external loads and self-weight the normal and shear stresses in the beams are restricted from above. The variation of the stiffness along the span of the beams is
Abstract-The
restricted to splines of order zero, one or two, thus reducing the optimization problem to a linear or non-linear programming problem depending on the cross-sectional shape. The latter is reduced to a sequence of linear programming problems which are solved by the highly-efficient active set method. Several examples are given to illustrate the potential of this method. A computer code in standard
FORTRAN 77 is also provided.
The minimum-weight designs of efastic frames under multiple loads are known [l, 21. The designs meet stability [I] or stiffness criteria [2] but would have contained sections with vanishing stiffness at the positions of vanishing bending moment had the stiffness not been restricted from below. In contrast to the above maximum-stiffness designs, the maximum-strength designs [3-51 do not exhibit sections with vanishing stiffness provided that both the normal and the shear stresses are restricted from above. It should however be stressed that the maximum-strength design of a statically indeterminate flexural member or a frame is the least-weight design only if both the stress constraints are active. This is unlikely to be always the case. In fact, such a design will generally be understressed except at a few sections. The ideal maximum-strength designs (as, indeed the ideal maximum-stiffness designs) are also difficult to fabricate economically because of their awkward geometries which reflect the variations in bending moment and/or shear force. To avoid awkward geometry the span of a flexural member may be divided into several segments, in each of which the stiffness (mass) is permitted only a limited variation as specified by spfines of order zero, one or two. This approach has been taken in the present work. The member to be optimized is subject to multiple loads. Under each load system consisting of external transverse loads and self-weight the normal and shear stresses are restricted to certain allowable values. The corresponding problem under a single load system was recently reported [6]. The minimax optimization problem is first reduced to a mathematical pro~amming problem by repIacing the analytical expressions for the stress constraints and objective function with splines of desired
order. The resulting mathemati~l programming problem is linear or non-linear depending on whether the mass-stiffness relationship of the member is linear or non-linear. In the latter case the sequential linear programming technique with move limits [7, S] was found to be the most efficient solution method. Moreover, the active set method (91 for solving the linear programming problems was found to be far superior to the traditional (reduced) Simplex algorithm. Several examples are given to illustrate the potential of the present technique. A c.omputer code in standard FORTRAN 77 is also provided.
Consider an isolated elastic member in flexure (Fig. 1) that forms part of a multispan beam or frame. The member is subjected to reactive forces Pi (i= I,... ,4; axial forces have been neglected) from the joints of the beam/frame to which it is connected, as well as external transverse forces f(x) and selfweight g(x). The external forces f(x) consist of several individual force systems {f;(x), I = 1, . . . , NL ) each acting on the member on separate occasions. The equilib~um state of the member of given stiffness under any of the load systems is completely described by the two force resultants, namely, bending moment M,(x) and shear force Q,(x), together with the appropriate boundary conditions. The force resultants are related as follows:
dQ/
= dx in which member. 521
x refers
-lxx> -g,(x),
to an arbitrary
section
along
the
522
B. L. KARIHALOOand
S.
KASAGASUNDARAM
I
I
Fig. I. An elastic flexural member subjected to external load system f,(x) and reactive forces (P,, i=l,... .4) from joints. (Gravitational force due to self-weight is not shown.)
Throughout this paper, it will be assumed that the mass of the member, which is proportional to the cross-sectional area A(x), is related to its flexural stiffness, which is proportional to the second moment of area I(x) through I(x) = CA”(X),
(3)
in which c and n are constants determined by crosssectional shape. For most member shapes, the constant n takes an integer value between one and three. Thus, for instance, n = 1 represents a rectangular section of constant depth, n = 2, a geometrically similar section, say circular or square, and n = 3, a rectangular section of constant width. Here we report results for beams of rectangular section only (n = 1 or 3); the solution of the optimization problem with n = 2 follows along similar lines. The strength (stress) constraints on the member of known cross-sectional shape may be written in the form of the following differential inequalities [6]: ui{x, s, f,, g,, M,, Q,, A,, (A,),} G 0 (i=l,2),
(I=1 ,...,
NL),
(4)
in which x and s are the co-ordinates of the crosssectional fibre. The two stresses refer to normal and shear stresses, respectively. In the sequel, we will omit reference to the load system unless otherwise required. The optimization problem for a flexural member isolated from a multi-span beam (or frame) may now be stated as follows. Determine the control A (x) and, thus, I(x) that meets the strength requirements (4) (the differential inequalities) under all loading conditions and minimizes the mass, W, of the member given by
in which y is the mass density of the material. The optimization problem for every member of the continuous beam or frame may be stated in a similar way. When the compatibility conditions at the member connections (joints) are invoked, the problem for the continuous beam or frame as a whole is formulated. The differential game problem described by eqns (4) and (5) can, however, be reduced to a minimax variational problem. MINIMAX PROBLEM
Assuming that the equilibrium equations, eqns (1) and (2), with the appropriate boundary conditions, can be solved in a closed form, the expressions for the force resultants M(x) and Q(x) are known for each load system. It should be mentioned that M(x) and Q(x) depend not only on x, f and g but also on the control A(x) itself through the compatibility conditions and self-weight. The control A(x), in turn, depends on all the load systems. Substituting the (implicit) analytical expressions for M and Q, i.e. M=M(x,j;g,A,A,), Q=Q(~,frg,A4), into the differential inequalities [eqn (4)] gives, for each load system, ai]x* s, .L g, M(x, J g, A, A,),
Q(x,f,g,A,A,),A,A,l~O,
(64
or, for brevity, S-&(x,s,f, A, A,) < 0
(i = 1,2).
(6b)
It should be mentioned that as g(x) = yA (x) it has not been repeated separately in eqn (6b). Next, the maxima of CIi with respect to s and f are determined for fixed values of x and A (x). Let these be attained at s = s,’ and f =f ,‘, i.e. Di(x, stq f:, A, A,)
w=y
A (~1
dx,
(5)
= max m/“” f&(x, s, J A, A,), I
(7)
Optimum design of statically indeterminate
d,=O
beams
$+1-L
4
d2
523
la1
I-*
lb)
d,=O
d2
(cl
l--x Fig. 2. Representation
dN+t=L
di
dg
of cross-sectional
and denote
area A (x, pi) using splines.
(9) become
Wx, A, A,) = @lx, st(x, A, A,),ft(x),
A, A& (8)
In practice, as the number of load systems {i.e. NL) is discrete, the maximum with respect to f is to be chosen among a finite set. Using this notation, the differential inequalities (4) reduce to
f-&,*(X,&~>,<0(i=l,2),(k=I
,*.., N).
(11)
Thus, the original game problem [eqns (4) and (5)] reduces to the following variational one: Minimize the functional (mass) L
i2y(x,
A,A,)
(9)
Next, each member is divided into a given number of segments (NSEG) and the order of variation (NOR) of cross-sectional area A (x) is prescribed using splines. For example, splines of order zero are used to represent constant cross-sectional area within a segment and those of order one and two to represent linear and parabolic variations, respectively (Fig. 2). Irrespective of the order of spline chosen to represent A (x) in each segment, it will be of the form
Y A(x,pk)dx 50
(k=I,...,N),
(‘12)
subject to the stress inequalities (11). It should tie emphasized that the latter include a maximisation operation with respect to load systems. App~i~tion of the above reduction procedure is iflustrated on a dexural member in the next section. APPLICATION TO FLEXURAL MEMBERS
To iflustrate the concepts intr~u~ in the preceding section, consider an arbitrary, horizontal, A(x)=A(x,p,) k=l,...,N, (10) flexural member subjected to external loading, as well as reactive forces from the joints of the beam (or where pr: denote the unknown (discrete) design frame) of which it forms a part (Fig. 1). Let the variables, and N is the total number of variables; reactive for,rces due to load system f;(x) be P,i 4) and let the external transverse loading N = NSEG, if NOR =Q and N= NSEC + 1, if (i=I,..., NOR = I or 2. The design variables could be the consist of several individual load systems f;(x), I=1 ,...,NL. It is assumed that a typical load variable width b&x)@ = 1) at the knots (k = 1, . . . , N), or the variable depth /~(x)(n = 3) at the system I;(x) will consist of a distributed load q&), knots. Using the representation (W), the inequalities concentrated load Wfjlocated at positions uj along the
524
B. L. KARMALOOand
member length and concentrated couples C,j acting at locations b,, as shown in the figure. If the reactive forces P,i (i = 1, . ,4) are assumed to be known, the bending moments ~,(~~ and Q,(x) due to external load system h(x) are given by
z w,,
M,(x)=P,,x-
-cl,>
S. KANAGASUNDARAM b(x), the strength constraints with eqns (16) give 12 _max mFlM,(x)s(x)l ” = b(x)h’ s 6
I
62
1-1
(IS) in conjunction
=
1
h*
b(x)h’
my
(17a)
-co G 0
7-$2(X)
[
x maxIQ,(x)l - 7. G 0.
(17b)
f
-
O~hhhN~~l~~
-v)dtl
(13)
I
Q,(x)=
Pir- f W,,0 j-t
in which m, and m, are respectively the total number of concentrated loads and couples acting on the members. The integral term involving the section area A(x) in each of the above two expressions is the contribution from the self-weight of the member. The strength constraints limiting the intensity of the normal (longitudinal) stress, cxxr and transverse (shear) stress, T__.at the section x are
The maxima of B, and a2 with respect to s are attained at s: = + h/2 and s: = 0, respectively, and so give 6 a: E-maxIM,(x)j-g,GO b(x)h2 / ny-_ - 3 Zb(x)h
(18a)
maxlQ,(x)l/
(18b)
r,G 0.
It is interesting to note that the differen~al inequalities (9) degenerate into algebraic inequalities (18) for the cross-sectional shape under consideration. Therefore, the lower bound on the optimum width of the member may be written, by inspection,
(19) where in which a, and to are known positive constants (allowable working stresses). From elastic beam theory it is known that
o(x) =max{lQ,(x)l,
Q2(xll,. . . * lQdxx>i>. CW
If it is assumed that the constant-depth beam consists of segments of linearly varying (NOR = 1) width, b(x) (Fig. 2b), it follows that for the ith segment, b(x)=(l-T)bi+
in which s(x) determines the position of the crosssectional fibre and K and a (x, s) are dependent on the cross-sectional shape. For a rectangular section of constant depth (h) and variable width b(x), i.e. for n = 1, K = l/~(~),u(x,~)=~(~)~~z/4-~z(x)]/2and -h/2 Q s G h/2, and for a rectangular section of constant width (b) and variable depth h(x), i.e. for n=3, K=l/b, a(x,s)=b[h2(x)/4-s2(x)]/2 and -h(x)/2Gs(x)dh(x)/2. In view of the dependence of a (x, s) and solution technique on the cross-sectional shape, i.e. on the value of n, it is convenient to illustrate the method separately for two typical values of n. Rectangular section of constant depth (n = 1) In the case of a flexural member of rectangular section with constant depth, h, and variable width
Z’bi+l,
dfGxGd,+,,
(21)
where T = T(x) = (x - di)/(di+ , - d,), and b, are the values of the optimum width at the knots. The stress constraints (19) for the segment d, < x g di+l can, therefore, be rewritten as , i=l
The functional
, . . . , NSEG.
(22)
to be minimized (12) becomes
W=~~(b,+b,+,)(d,+,-d,). 1-I
(23)
The problem of mi~m~ation of W (23) subject to the inequalities (22) is a linear pro~amming problem
525
Optimum design of statically indeterminate beams which is solved as follows:
where
1. Assume (bi)j, i = 1,. . . , N, in the first iteration 1. 2. Analyse the beam by standard stiffness matrix method (e.g. [3]) for all the load systems f;(x), and obtain n?(x) and e(x) [eqns (13, 14,20)]. 3. Divide each segment into m intervals, and solve the following linear programming problem using the active set algorithm [9]:
JI,= -
j=
Minimize
W=~‘~(~,+&+,)(CI~+, r-l
- 0
(24)
3
2bh2(x)
max Q,(x)h - M,(x); /
sf = 0,
(27)
lj12=3
sf2 = h’/4,
bh2(x)
G(x)= max(lM,(x)l,. The algebraic inequalities
,
(28)
. . , IMNL(x)l). (26a) reduce to
subject to [I - T(Q)]& +
T(Xk)Ei+
i=l,...,
h(x)>
I Z max(F,, F2), NSEG, k = 1,. . . , m + 1,
d,dxk
N.
(29)
The minimization of mass W [eqn (5)] subject to linear (29) and non-linear (26b) constraints results in a non-linear programming problem. It is however efficient to solve it as a sequence of linear programming problems as follows. First, the differential inequalities (26b) are rearranged such that if JI, > $2r we have
4. Repeat steps 2-3 until l(bi),+ 1-(b,),]g~(=0.5%),
i=l,...,
5. Calculate the minimum
mass (23).
LmaxlQ,h
N.
26P(x)
/
-M,$~-~o
or
The above steps were programmed on VAX 11/780 and several examples were solved. The results will be presented and discussed later in the paper. If the optimum beam is desired to have segments of constant (NOR = 0) or parabolically (NOR = 2) varying width, the optimization problem in step 3 above will still be a linear programming problem (see
fl = 3/(26r,). Similarly, if ti2 > I/,, we have
WI). (31)
Rectangular section of constant width (n = 3)
Substituting eqns (16) into the stress inequalities eqn (IS), we have 12 ul =
3
bh (x)
max m,a” IM,(x)s (x) ( - co G 0 s
(25a)
Assuming that the quantities within the square root sign in inequalities (30) and (31) are known, we may combine the inequalities (29)-(31) and write h(x)>max{F,,&F,},
(32)
where F,, F2 and F, are, respectively, the right hand sides of inequalities (29)-(31). Depending on the relative magnitudes of F,, F2 and F,, the inequalities (32) may be linear or non-linear. The minimization of mass The maxima of u, with respect to s are attained at s: = + h (x)/2, whereas that of c2 are attained at either sf2 = 0 or s:2 = h2(x)/4. The differential inequalities [eqn (7)] therefore reduce to *+
=6&l --
’Q: C.A.S I(,)_”
E
bh2(x)
00 Q cl
max($,. q2) - r. d 0,
(26d (26b)
L
W = yb
50
h(x)&
subject to the stress inequalities out in the following steps.
(33)
(32) is then carried
1. Assume h,(x) in the first iteration j = 1. 2. Analyse the structure by standard stiffness
8. L. KAIUHALCXJ and S. KANACASUNDARAM
526
BOkN
ZOkN/m Irrflrrtrlrr)
3m
6m 6m
I.-
3m
Load System 2 --
Load System 1 -PrlsmatlcBeam
0.5O.&030.2=---.,__ 0.1-
1 --____
00, 0.50Ab 0.3lm) 0.2\ 0.1O.OL 0
I 1
I 2
I 3
1 c
I 5
6
x Im)
Fig. 3. Minimum-weight design of pinned-clamped beam subjected to load systems as shown plus self-weight (member has constant depth h = 0.65 m, i.e. n = 1).
method [3] for all the load systemsfi(x), and obtain the force resultants M,(x) and Q,(x), (I = 1,. ..,NL), from eqns (13) and (14). 3. Solve the non-linear programming problem [eqns (32), (33)] using a sequential programming
approach with move-limits [7,8] to obtain h,, ,(x). 4. Repeat steps 2-3 until Jh,,, -h,)=26(=0.5%). 5. Calculate the minimum mass.
As with the case for n = 1, the optimum depth h(x) is approximated by splines of desired order [eqn (1O)]. It should be noted that if stiffness is constant within a segment (NOR = 0) then dh/dx = 0, and eqns (30) and (31) reduce to algebraic inequalities. The optimization probiem is linear in this case and step 3 in the above procedure reduces to the solution of a linear programming problem.
1.0PrismaticBeam
0.8-h
I
0.6-
______----
_,_---_
(ml 0.4__________________----------0.20.0.
1.0i h fm)
0.80.6-
-
O.L-u 0.20.0. 0
I 1
I 2
1 3
I Ir
I 5
x (rn)
Fig. 4. Minimum-weight
design of the beam shown in Fig. 3, but having rectangular cross-section of constant width b = 0.4m, i.e.n = 3.
6
521
Optimum design of statically indeterminate beams
30kN/m A GA
1
30kN
&IS
4m
6m
I
t---x 0
Load --
System
1
-Load
sPrismatic
0 .L-
Prismatic
Beam r
b o 3(ml 0 .2-
’
.*
System
2
Beam
.*
0 .l-‘, 0 .o0 Ab
o3 -
tm) 0.2 0 .I-\ 0 .o-
0
2
1
3
.4
5
6
7
9
8
x Im) Fig. 5. Minimum-weight design of two-span continuous beam subjected to load systems as shown plus self-weight (members are rectangular in cross-section with constant depth h = 0.65 m, i.e. n = 1). EXAMPLES AND DISCUSSION
The results of several numerical examples are shown in Figs 3-8. The following material properties were assumed in all the examples: u,, = 5 MPa, 7. = 0.5 MPa, E = 20 GPa and y = 2450 kg/m’. These values are typical for concrete. Data relating to the
O.E0.6-
Pritmatlc
Beam
I I _/__--- --
04 -
I A---__
external load systems acting on the beams are shown on the respective figures. Each member was designed to consist of l-m long segments, and the stress levels were calculated at sections spaced at 0.05-m intervals (i.e. m = 20). The integrals were evaluated using Simpson’s rule, and the iterations were terminated with c = 0.5%.
Prismatic
Beam
I 6
I 7
I --__ ---___
//* 0.2 I------L0.0, 0.6h lm)
0.4 0.2 = 0.0
0
I 1
I 2
I 3
I L
I 5
I 8
I 9
x Im)
Fig. 6. Minimum-weight design of the two-span continuous beam shown in Fig. 5, but having rectangular cross-section of constant width 6 s 0.4 m, i.e. n = 3.
10
528
B. L. KARIHALOO and S. KALAGASUNDAIUM
MkN
Load System
0.4
Prismattc
/’ -c,-’
oq-
Load System --
2
Beam
0.3 -1.
1
each
,’
/’
.’
0.4 b
0.3-
fm) 0.2 0.1 0.0
0
1
I 2
I 3
I 5
I I
I 6
I 7
I a
I 9
x (ml
Fig. 7. Minimum-weight design of two-span continuous beam subjected to load systems as shown plus self-weight. The concentrated forces in load system 2 represent a truck loading. The members are rectanguIar in cross-section with constant depth h = 0.65m, i.e. n = 1. Figures 3 and 4 show the minimum-weight design of propped cantilever beams of variable width (n = 1) and depth (n = 3), respectively. Figures 5-8 show the corresponding designs of two-span continuous beams. The load systems consist of dist~buted loads, concentrated Ioads and couples. The designs using
splines of different orders are compared. The designs using first and second order splines have been drawn separately since they would have been indistinguishable on the same diagram. A comparison of the volumes of op~mum beams is presented in Table 1. Also indicated are the volumes of the corresponding
0.6h
lm)
0.4 0.2 0.0
1 1
I
2
I
3
I
I
1
I
I
I
Ir
5
6
7
8
9
Fig. 8. Minimum-weight design of the two-span continuous beam shown in Fig. 7, but having a rectangular cross-section of constant width b = 0.4 m, i.e. n = 3.
Optimum design of statically indeterminate Table 1, Comparison
of volumes of optimum designs Volume (m’) Prismatic
Figure
3 (n = 1) 4 (R =3) 5 (n = 1) 6 (n = 3) 7 (n = 1) 8 (n = 3)
NOR=0
NOR-I
NOR=2
design
1.0695 1.1845
0.9733 1.1190
0.9996 1.1358
I .8639 1.8612
I .6943 1.8858 0.9449 1.2286
f.4411 I .6504 0.7864 0.8766
1.4149 I .6452 0.7913 0.8749
2.8037 2.7783 I so70
I .4933
prismatic designs. It should be mentioned that the total number of fully stressed sections in the minimum-weight beams is equal to the number of (discrete) design variables, the only exceptions being the beams of variable depth (n = 3) with first or second order splines in which fewer sections are fully stressed. All of the minimum-weight beams presented here satisfy the serviceability requirement in that the maximum deflection to span ratio is less than l/375. The non-linear programming problem associated with n = 3 was solved by a sequential linear programming technique with move-limits [7,8]. This technique was found to be the most efficient and reliable. In fact, the same optimum design was reached from different starting designs (both uniform and non-uniform), thereby pointing towards the global nature of the optimum. (The same statement is also valid for the optimum designs for n = 1.) The move-limits at a given iterative step were set equal to the 80% of the corresponding values at the previous step. The linear programs were solved by active set analysis [9j which was found to be very efficient and to use much less CPU time compared to the traditional approach using the (reduced) Simplex algorithm. A listing of the computer code in FORTRAN 77 is given in the Appendix. It should be mentioned that the sequential linear programming algorithm that was so efficient for the present problem may not be
529
beams
as efficient for solving other non-linear
programming problems. In fact, WChave found it rather unreliable for solving non-linear programming problems resulting from minimization problems under stress and deflection constraints. For these latter problems, the convex linearization and penalty function techniques seem to be more efficient and reIiabIe. This is currently under investigation and the results will be reported in a future communication.
Acknowledgement-The work reported in this paper was supported by grant F84/16078 from the Australian Research Grants Scheme to B. L. Karihaloo.
REFERENCES
I.
H.
K. Turner and R. H. Plaut, Optimal design for
stability under multiple loads. J. &ng AXE 2.
3.
4.
5.
6.
Mech: Die.,
106, 1365-1382 ~19801.
P. Pedersei and L. Jorgken:
Minimum mass design of elastic frames subjected to multiple load cases. Compur. Srrucr. 18, 147-157 (1984). S. Kanagasundaram, Optimum structural design (application of differential game theory). Ph.D. Thesis, The University of Newcastle, Australia (1984). S. Kanag~undaram and B. L. Karihaloo, Ootimal strength design of beam-columns. fnr. J. Sokk kr~. 19.937-953(1983h S. -Kanagasundaram and B. L. Karihaloo, Maximum strength design of structural frames. J. Strucr. Div., ASCE 111,1267-1287 (1985). B. L. Karihaloo and S. Kanagasundaram, Computeraided minimum-weight design of statically indeterminate beams. &grtg Opiimiz~f~on 10,139-156 (1986). P. Pedersen, The integrated approach of FEM-SLP for solving problems of optimal design. In Opfimizarion of Distribured Parameter Sysrems, Vol. I (Edited by E. 3. Haug and J. Cea), pp. 739-756. Sijthoff % Noordhoff, The Netherlands (I 98 1). P. Pedersen, On the optimal layout of multipurpose trusses. Comput. S~ruct. 2, 695-712 (1972). M. J. Best and K. Ritter, Linear Progrumming (Active Set Analysti and Computer Programs). Prentice-Hall, Englewood Cliffs, NJ (1985).
B. L. KARIHALOOand S. KANAGASUNDARAM
530
APPENDIX:FORTRAh 77 LISTINGOF SLP
0001
0002 0003 0004 0005 0006 0007 0008 0009 0010 0011 0012 0013 0014 0015 0014 0017 0018 0019 0020 0021 0022 0023 0024 0025 002L 0027 0028 0029 0030 0031 0032 0033 0034 0035 0036 0037 0038 0039 0040 0041 0042 0043 0044 0045 0046 0047 0048 0049 0050 0051 0052 0053 0054 0055 0056 0057 0058 0059 0060 OO6I 0062 0063 0064 0065 0066 0067 0068 0069 0070
SUBROUTINE
SLP(NV,tlINEQ,HEQ,DHIN,DMAX,TOL,IHAX,DELX, X,CC,IC,A,IROWA,6,C,DX,DINV,IROWD,S,V,WORK,JJ, AI,IROWA1,8I,Cl,DX1,DINVI,IROWDI,SI,VI,WORKI, JJ1,O8J,IOUT,JOlJT)
1 2 3
l
PURPOSE The program solves the following optimization problem: Hinimire
F
subject
to Ii (X ) j i
LE
G (X ) k i
EQ
* + 0 0 + + + t
non-linear
0,
0,
j-1
,..,MINEQ
k = i,..,HEQ
in which F is the objective function, H 's and G 's are the constraint functions k j and X , i=I,..,NV are the unknown design variables. i DESCRIPTION The program is based on Sequential Linear Programming (SLP) dpprOdCh using move-limits. In this approdeh the original problem is reduced to a sequence of linear programming problems with ohdngcs in the design variables (i.e., dX i i-1,. .,NV) ds the unknown vdridbles.The linear sub-problems are solved using active set strategy. The program is written in standard FORTRAN 77 for ease of portability. No derivdtivas are required. PARAMETERS
AND USER SUPPLIED
l
t + + t + t t + + + + t + + + + + + + + + + +
+ +
ROUTINES
4
I. Input
parameters
+ 4
NV -
INTEGER. NV specifies the number On entry, vdridbles. Unehdnged on exit.
+
of
design
4 + +
HINEQ
- INTEGER. HINEQ speoifies the number of On entry, inequality constraints. Unohdngcd on exit.
4 4
4 4
MEQ - INTEGER. On entry, MEQ specifies the number Unchanged on exit. constraints.
+
of
equality
4 4 +
DMIN
-
REAL array On entry,
of
DIMENSION GE NV. specifies DtlIN(I)(I-I ,..,NV)
lower move-limit design variable.
on the
Chdnge
the in the Ith
+ + l + l
DMAX -
REAL
array
of DIMENSION GE NV. On entry, DMAX(I)(I=I ,..,NV) specifies the upper move-limit on the ohdnge in the Ith design variable.
+ l 4 l +
TOL - REAL On entry, TOL must be set to specify the dccurdcy to which of change in the design variables the solution is desired. Unchdnged on exit.
l * + 4
+
Op(~~urn design of staucalfy jnde~errnIn~~e beams 0071 0072 0073 OQ74 0075 0076 0077 0076 0079 0080 0081 0082 0033 0004 0085 0006 0087 0080 0089 0090 0091 0092 0093 0094 0095 0096 0097 0098 0099 moo 0101 0102 0103 0104 0105 0106 0107 0108 0109
* *
IHAX
l
* * * *
-
l
* *
IBAX specifics the maximum to be performed. Unchanpcd
* +
number of on exit.
l
l
REAL On entry, DELX specifies an increment in design variables. Used in the calculation graditnts of the constraint and objtctive functions. For example,
the
l
af
* * l
*
l
dH(X
+ * *
-*-___-
* * * * * * * * + * * *
H
+DELX)
-
H(X
l
1
*
i i ________c_______I____
l
DELX
*
dX
l
i
l
d rtfcrs
whcrt
to
partial
l
derivative.
l
IOUT
-
*
INTEGER. on entry,
IOUT
associated
specifits
with
the
logical
output
l
unit numbtr Unchanged on
fife.
* * * *
exit. 2.
Input/Output
prramaters
l
REAL 8rray of l?fHENSiON Before entry, X(Il,i=i,..INV starting vrlues of the
X -
l
* * * * * * e * *
) i
l
On exit, X contains the dtsign VdridbieS
3.
Output 08J
*
GE NV.
must be set design
l
to
* *
vdridbles.
the optimum values
of
l
0 * * * *
l
parameter -
REAL. On exit 083 contains the objective function.
the
optimized
valut
of
4 * 4
l
*
+
oiii
*
0112 0113 0114 0115
* * * *
0116
l
4.
Workspaoe
(and
*
0119 cl120 012% 0122 0123 0124 at25 0126 0127 0123 0129 0130 0131 0132 0133
* * * * * * * * * * + * * * *
* * *
-
Dimension)
4
Parameters
INTEGER. On entry, IROWA must specify the first __ dimension of array A as declared in the CallanQ (sub) program. IROWA must bt GE (HINEQ + IlEQ +P*NVf.*
l 4 l
4
IAOWAl
-
INTEGER. On entry, IROWAI must dimension of array Al (sub) program. IROWAl + 2eNVl.
4
specify the as declared must
be
6E
first in the
4
callinq + FiEQ +
(XINEQ
*
1
l
* l
IRDWD
-
4
INTEGER. On entry,
lROWD must specify of array DINV as dtclared in
tht the CE
first calling CHINEQ
dimension
l l
+
l%EQ +2*NV).* l
IROWDi
-
*
INTEGER. On entry*
IROWDf must specify of array DINVl as declared in (sub) ‘program. XROWDl must be + 2*NVl.
the first dimension the calling 66
(HINEQ
+ UEQ
+
* l
1
l l 4
JJ,JJl
-
INTEGER PIWI
of
arrays
DIMENSION
* respectively,
#
GE NV dnd
+ 4
A,Al
-
l
* * * * * *
Associattd
l
IROWA
*
oiie
0134 0135 0136 0337 0138 0139 0140 0141 Oi42 0143
INTEGER On entry, iterations
* DELX
OlfO
0117
-
531
REAC arraym renpeotivtly,
of of
DIMENSION CIROWA,p) and where p must be GE NV.
fIRDWA,p+l)
4 4 4
B,Bl
-
REAL arrays and CHINEQ
DIHENSION + HEQ+l+?*NVl,
DINV, DINVl
-
REAI. arcrys of rtspeatfvely,whtre
GE (HfNEQ+MEQ+2*NV) respectively.
4 4 4
DItiENSION (1ROWD.p) and p must bt GE NV.
(IROWOi,p+l)
*
4
B. L. KARIHALOO and S. KANAGASLWDARAM
532 0144 0145 0146 0147 0148 0149 0150 0151 0152 0153 0154 0155 0154
0157 0158 0159 OibO 0161 0162 0163 0164. 0165 Oibb 0147 OIL0 0169 0170 0171 0172 0173 0174 0175 0176 0177 0170 0179 0100 0181 0182 0183 0184 Oi85 0186 0187 0188 0189 0190 0191 0192 0193 0194 0195 0196 0197 0198 0199 0200 0201 0202 0203 0204 0205 0206 0207 0208 0209 0210 02ii 0212 0213 0214 0215 0216 0217
a 4
a a + a a + a + a + * + + + + +
DX,C,S,V,WORK - REAL
arrays,
each
of
DItlENSION
GE NV.
DXl,C1,Sl,V1,WORK1 - REAL +rr+ys,
eroh
of
DIMENSION
GE NV+l.
INTEGER. entry, XC must specify the CC. XC must GE mrx(HINEQ,HEQ).
dimension
IC
-
cc 5.
on
-
REAL
User-supplied FUN1
-
a +
a + + + * + + + + a a a a + a
+rrry
of
Routines
(XC).
Example
et
end
of
proqr+m)
SUBROUTINE. This routine must be supplied to calculrte the values of all the “inequality” constraint functions for a given Set of design variables. The specifio+tion is :
NV -
MINEQ
-
X -
+
a + a + +
CI
FUNE
-
-
INTEGER. On entry, NV specifies design variables. Its be chrnqed in FUNI.
the value
number of must not
INTEGER. On entry, HINEQ specifies the number of inequality constraints. Its value must not 'be ohanqed in FUNI. REAL array of DItlENSION On entry, X contains the of the design variables to which the oonstraint FUNI must are required. values in X.
(NV). ourrent values corresponding function values not’ change the
REAL array of DIMENSION (HINEQ). FUN1 must set the elements of CI lqurl to the values of the “inequality” constraint functions corresponding $0 the design variables contained in X.
SUBROUTINE. This routine must be supplied to oaloulate the values of all the “equality” constraint functions for a given set of deriqn variables. The specification is I SUBROUTINE FUNE (NV,MEQ,X,CEl IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIHENSION X(NV),CEfflEQ) NV -
a +
a + + + + + + + + + a
(See
array
SUBROUTINE FUNI (NV,flINEQ,X,CIl IMPLICIT DOUBLE PRECISION (A-H,O-21 DIMENSION X(NV),CIfHINEQ)
4
a a + a a + + a + a a + + + + a + + +
DIflENSION
of
MEQ -
X -
INTEGER. NV specifies On entry, Its design variables. be oh+nqed in FUNE.
the value
number of must not
l
+ + +
+ + + + + + + + + + + + + 4 + 4 + + + + + + + + 4 + + + + a + l
I
I II . 1 t t 1 ,
4 4
I
XNTEGER. On entry, ME4 specifies the number of inequality constraints. Its value must not be changed in FUNE. REAL rrrry of DIMENSION X contains the On entry, of the design variables to whioh the constraint FUNE must l rc required. values in X.
a a + + + + + + + + + + + + +
(NV). current values oorrespondinq function values not change the
4
i :’
Optimum design of statically 0218 0219 0220 0221 0222 0223 0224 0225 0226 0227 0220 0229 0230 0231 0232 0233 0234 0235 0236 0237 0230 0233 0240 0241 0242 0243 0244 0245 0246 0247 0248 0249 0250 0251 0252 0253 0254 0255 0256 0257 0258 0259 0260 0261 0262 0263 0264 0265 0266
533
beams
*
l
CE -
l l
REAL FUNE
of set
array
must
DItlENSION
l
(HEQ).
the elements of CI equal to the values of the "equality" constraint functions corresponding to the design variables contained in X.
l l l
b 4l * l I)
l
OBJF
l l l l
-
SUBROUTINE.
l
Thie routine must be supplied to calculate the value of the objective function F(X) for a given set of design variables X. The specification is :
l
l
l l l
l
SUBROUTINE
l
IKPLICIT DOUBLE DIflENSION X(NV)
l
OBJF
(NV,X,FO) PRECISION
+
(A-H.O-2)
4 I)
l
*
l
NV - INTEGER. On entry, NV specifies the design variables. Its value changed in OBJF.
* l l
4t
number
of
must not be
l l l
*
l
x-
l
REAL
array
of DIIIENSION (NV).
On entry, X contains the of the design variables
l l l
to which the constraint rre required. OBJF must values in X.
l l
l
current
values correspondinq function values not change the
l * it * *
l
l
FO - REAL.
l
t
l
OBJF
l
the objective to the design
l
must
set
FO equal
to the value of function corresponding variables contained in X.
l
(I it I il
l l l
* l
IIOVE - SUBROUTINE. This routine calculates new move-limit9 on the changes in the design variables rcoording to user-specified criterion. It is not called in the firet itcrdtion.
* *
+ * * l l +
SUBROUTINE HOVE (NV,DHIN,DtlAX) ItlPLICIT DOUBLE PRECISION (A-H,O-2) DIMENSION DHIN(NV),DflAX(NV)
l l
l + *
l
*
NV - INTEGER. On entry NV specifies the number of design vrridbles. Its value must not be ohdnged in ROVE.
l l l l
4
DHIN - REAL array of DIMENSION (NV). On entry oontdins lower move-limits used in the previous iteration. On exit oontdins new lower move-limits calculdted dooordinq to user-specified LE 0, I-l,..NV. criterion. OKIN
l
0276 0277 0278 0279 0260 0201 0202 0283 0204 0285
l
0206 0207 0200 0289 0290
l
Best,
l
Aative
l l l l
DKAX - REAL rrrry of CIUENSION (NV). On entry oontdins upper move-limits used in the previous iterdtion. On exit oontdins new upper move-liritr OdIoUldted dooording to user-epecified GE 0, I=i,..~V. criterion. OKAX
l l l
l l
6. Other
routines
referenced
PIPZC,
l
* l l
l l l
* l l
l
The routines
PlPPC
dnd SOLN oan be obtdined
fron
* l
l
*
l
l
SOLN.
l
l
t
l
l
l
l
t 4
l
4
l
l
l
*
4 * 4 il
l
0267 0260 0269 0270 0271 0272 0273 0274 0275
CA.S x.l-4
indeterminate
tt.J., dnd Ritter, K., LINEAR PROGRAMINS: Set Analysis dnd Computer Programs, Englewood Cliff* New Jersey Prentice-Hdll InO., P-P- 254,283, 1985.
# l l
I) 4
534 0291 0292 0293 0294 0295 0296 0297 0298 0299 0300 0301 0302 0303 0304 0305 0306 0307 0308 0309 0310 0311 0312 0313 0314 031s 0316 0317 0318 0319 0320 0321 0322 0323 0324 0325 0326 0327 0328 0329 0330 0331 0332 0333 0334 0335 0336 0337 0338 0339 0340 0341 0342 0343 0344 0345 0366 0347 0348 0349 0350 0351 0352 0353 0354 0355 0356 0357 0358 0359 0360 0361 0362 0363
B. L. KARIHALOO Note that SOLN must parameters
4 l
it
and S.KANAGASUNDARAM
the be
argument list of routines PIPPC expanded td include the extra like IROWA,IROWAl,IROWD,IROUDl and
and
l
IOUT.
l
l
4t
l
*
7.
4 4 4 4 4 4 4 4
Error
Indicator JOUT
=
2
The from
linear subproblem is unboupded below.
JOUT
-
3
The linear subproblem has no feasible solution. Dr.
I
S.Kanagasundaram
C IllPLICIT
DOUBLE
PRECISION
(A-H,O-2)
C DIHENSION
DHIN~NV~,DflAX~NV~,CC~IC~IA(IROWA,NV~, 6(llINEQ+tlEQ+2+NV),C(NV),DX(NV),DINV(IRObJD,NV~, S(NV),V(NV),WORK(NV),JJ(NV),X(NV) Al(IROWAI,NV+1),8l(MINEQ+tlEQ+l+2*NV~,Cl(NV+l~, DXl(NV+1),DINVl(IROWDl,NV+l),Sl~NV+l~,Vl~NV+l~, WORK1 (NV+l) , JJl (NV+1 )
1 2 DItlENSXON 1 2 C
tiTOT=HINEQ+flEQ WRITE(IOUT,100)NV,klINEQ,tlEQ,TOL,ItlAX
ITE=l WRITE(IOUT,lIO)ITE PRINT 4, SLP_ITE-’ WRITE(IOUT,12O)(X(I),I=l,NV)
, ITE
Celculdtc move-limits IF(ITE.NE.I)CALL Calculdte functions
MOVE(NV,DMIN,DMAX)
values of 'inequrlity" constraint and their gradients
CALL FUNI(NV,HINEQ,X,CC) DO 2000 I=l,NV X(I)=X(I)+DELX CALL FUNI(NV,tlINEQ,X,Bl) X(1)=X(I)-DELX DO 2000 J=l,HINEQ
2000 3010
A(J,I)=(B?(J)-CC(J))/DELX DO 2010 i=l,HINEQ SUtl=O. DO DO 2020 J=l,NV
2020 2010 C C C
SUH-SUH+DHIN(J)tA(I,J) 8(L) --CC(L)-sull Impose
upper
limits on DX
HI =MINEQ+NV J-O DO 2030
I-HINEQ+l,Hl J-J+1 A(I,J)=l.DO 2030
B(I)=DflAX(J)-DHIN(J) tl2-tll +NV
C C C C
l
4
PROCRAWlER
4 4 LAST HODIFIED I 27 OCT. 1986 4 444444444444*444444444444444444444444444444444444+44444444444444444444444
C 1000
4 4 4 4 4 4
Include
non-negativity
J-0 DO
L=tll+l,fl2
2040
J-J+1
constraints
l
4 4 4
Optimum design of statically indeterminate 0364 0365 03cllr 0367 0368 0369 0370 0371 0372 0373 0374 0375 0376 a377 0378 0379 0380 0381 0382 0383 0384 0305 038& 0387 0380 0389 0390 0391 0392 0393 0394 039s 0396 0397 0398 0399 0400 0401 0402 0403 0404 0405 0406 0407 0408 0409 0410 0411 0412 0413 0414 0415 0416 0417 0418 0419 0420 0421 0422 0423 0424 0425 0426 0427 0428 0429 0430 0431 0432 0433 0434 0435 0436
2040 c c C c
CIlau1~t.L
values
funotionr
and
of their
“equality”
beams
ocnstraint
gradients
IF CMEQ.EQ.0) 60 TO 2050 CALL FUNE CNV,MEQ,X,CC) DO 2060 I-1,NV X(I)-XCI)*DELX CALL FUNE CNV,MEQ,X,BI> DO 2060 J-1 ,MEQ A~M2+J,I~-C0iCJ)-CCCJ~)/DELX DO 2070 I*l,MEQ SUM-O‘DO DO 2000 J-1,NV SUMmSUH+DMINCJ~~ACH2+I~J~ 6cN2+I)=-CCCI~-SUM
2060
20&O 2070 C C C 2050
Caloulete
of
coefficients
the
objective
function
CALL OB.YFCNV,X,FO) DO 2090 I*l,NV X(I)rXCf)+DELX CALL OBYF-DELX C(I)-CF-FO)/OELX
2090 c c c
Solution
of
tlT=H2+MEP Nl-NV+1 fill-NT+1 MTOTi -Ml
the
linear
progrrm
1
c CALL
PlP2C
CALL
SOLN
I 2
CA,IROWA,B,C,DX,DINV,IRO~D,S,V,WORK,JJ,O0JS,NV, tl2,HEQ,MT,Al,IROWA1,Bl,Cl,DXl,DINVl,IROWDl, Sl,V1,WORKI,JJI,O0JJI,Nl,Mll,MTOTl,JOUT,IOUT~ CJOUT,IOUT,DX,NV,OBJJ)
Check if solution of thr linear program is feasible. If not, double the value of the move-limitr and repeat the iteration.
3000
C c c 3020
3030 3050
c C :040
IFCJOUT.EQ.3)THEN DO 3000 I-1 .NV DMINCI~-2.O~*DHINCX~ DHAXCi~-2.DO~DMAXCI~ 60 TO 3010 END IF IFCJOUT.EQ.2)THEN WRITEC7,lSO) RETURN END IF
Check for
convergence
DO 3020 I*l,NV DXCII-DXCI)+DHINCI) DO 3030 I-1 ,NV XF~DABSCDXCI)).GT.TOt)GO CONTINUE DO 3050 I-1,NV XCil=DXCI)+XCI) CALL OBJF
iterations
ITE-ITE+l IFCITE.GT.IMAX>THEN
TO
3040
535
B.L.KARIWLOCI and S.KANAGASLJNDARAM
536
0437 0438 0439 0440 0441 0442 0443 0444 0445 0446 0447 0448 0449 0450 0451 0452 0453 0454 0455 0456
0001
0002 0003 0004 0005 0006 0007 0008 0009 0010 0011 0012 0013 0014 0015 0016 0017 0018 0019 0020 0021 0022 0023 0024 0025 0026 0027 0028 0029 0030 0031 0032 0033 0034 0035 0036 0037 0038 0039 0040 0041 0042 0043 0044 0045 0046 0047
WRITE(IOUT,180) STOP END IF DO 3060 I-l,NV X(I)-X(I)+DX(I) GO TO 1000 FORllAT(//,2X,‘NUHBER
3060 100 1 2 3 4 110 120 130 140 150 160 170 180
OF VARIABLES -‘,14,/2X, ‘NUHSER OF INEQUALITIES -’ ,14,/2X, ‘NUMBER OF EQUALITIES =‘,14,/2X, ‘CONVERGENCE CRITERION =‘,D12.6,/2X, ‘HAXIHUH NUMBER OF ITERATIONS -‘,X4) FORhAT(/SX,‘SLP_LTE -‘,X3) FORtlAT(/SX,‘VALUEB OF VARIABLEB’,/,/5fSX,DI2.6)) FORtlAT(/SX,‘LOWER HOVE-LIHITS ON VARIABLEB’,/,/5(5X,Dl2.6)) FORtlAT(/SX,‘UPPER HOVE-LIMITS ON VARIABLES’,/,/S(SX,Dl2.6)) FOR?lATf/2X,‘*** PROBLEPl IS UNBOUNDED FROM BELOW +**‘) FORtlATf/2X,‘OPTIHUH VALUE OF VARIABLES ..‘,/5(2X,Dl2.6)) FOR?‘lAT(//2X,‘OBJECTIVE FUNCTION VALUE AT OPTIliUt! =‘,D12.6) FORHAT(/PX,‘*** IlAXIHUfl NUHBER OF ITERATIONS PERFORtlED END
PROGRAM MAIN C C Calling program for routine C **************t********************************~*************************
SLP.
4
4 4
EXAHPLE
PROBLEM
4
4
4
Hinimirr
F
-
(X
-
+
1.0)
+
-X)
(X
1
4
3
2
2
l 4
(X
2
1
-
x
) 3
2
+
4
ix
4 4 4
-
x
3
4
4 4 4
4
4
l
4 4 4 4 4 4 4 4 4 4 4 4
)
+
-
cx 4
4
x
4 4
4
)
4
5
4 4
subjcot
to
4 4
2 X 2
4
+x
-x 3
+
2.0
-
2.0
*
SQRT(2.0)
GE
0.0
4
4 4 4
X
l
x
1
LE’2.0
4
5
4 4
X
LE
I.5
I
4 4 4
LE
X
1.5
3
4 4 4
LE
X
1.5
4
4 4 4
X
LE
4 5 4 4 3 2 4 l SQRT(2.0) +x - 2.0 - 3.0 X +x 4 3 1 2 4 4 This lx~mpla is taken from 4 Nummrioal Algorithms Group, NAG Library Manual, 4 1985. England, 4 *4444444*44444444*444444*44444*444444444*44444*444*444444444*4444*4444444
C
4
4
4
4 4
***‘I
1.5
4
4 4 4
0.0
4 4
4 4 4 4 4
Optimum design of statically indeterminate IRPLICIT
0048 0049
0050
0051 0052 0053 0054 0055 0056 0057 0050 0059 0060 0061 0062 0063 0064 0065 0066 0067 006.9 0069 0070 0071 0072 0073 0074 0075 0076 0077
c c C
Dlrfinr
problem
PARAHETER DIflENSION
DIMENSION 1 2 C C
Open
PRECISION(A-H,O-2) parameters
and
arrays
(NV=5,HINEQ=b,HEQ-1,IClb) DMXN(NV),DtlAX(NV),A~HINEQ+MEQ+2*NV,NV~, B~HINEQ+MEQ+2+NV~,C~NV~,DX~NV~,DINV~HINEQ+HEQ+2*NV,NV~ ,S(NV),V(NV),WORK~NV),JJ(NV),X(NV),CC(IC~ Al(flINEQ+HEQ+l+2+NV,NV+l~,Sl~HINEQ+nEQ+I+2*NV~, C1(NV+l),DXi(NV+1),DINVl(nINEQ+HEQ+l+2+NV,NV+l~, S1(NV+I),VI(NV+l),WORKI(NV+I),JJl~NV+I)
I 2
output
file in unit IOUT
10UT-7 OPEN~lJNIT=IOlJT,FILE=‘SLP.RES’,STATUSr’NEW’~ C C
Assign
initial move-limits to variables
C DO IO I=l,NV IO
OHIN(I1.00 DHAX(I)=l.DO
C
Define other
parameters
C TOL=l.D-04 IMAXDELX=I.D-05 IROWA=HINEQ+llEQ+2~NV IROWAl=HINEQ+HEQ+1+2*NV IROWD=HINEQ+tlEQ+2*NV IROWDl=HINEQ+HEQ+1+2+NV
0078 0079
ooao ooal 0082 0083 ooa4 0085 0086 0017 0088 0089
DOUBLE
beams
c C C
Starting values for the
variables
X(1)=0.500 X(2)=1 .ODO X(3)-0.6DO X(4)-1.4DO X(51-1 .ODO c 1 2
0090 0091 0092 0093 0094
3
CALL SLP(NV,HINEQ,HEQ,DHIN,DHAX,TOL,IllAX,DELX, X,CC,IC,A,IROWA,B,C,DX,DINV,IROWD,S,V,WORK,JJ,
A1,IROWAl,BI,C1,DXl,DINVI,IROWDI,Sl,Vl,WORKl, JJI,OSJ,IOUT,JOUT) STOP END
0001 0002 0003 0004 0005 0006 0007
SUBROUTINE ItlPLICIT DIHENSION C C
FUNI(N,li,X,B) DOUBLE PRECISION X(N),B(tl)
(A-H,O-Z)
Calculate "inequality" constraint functions
C
oooa
B(1) --X(2)+(X(3)**2)-X(4)-2.D0+2.DO+2.DO*DSQRTf2.DOy
0009 0010
BfP)=X(I)*X(S)-2.DO 8(3)-X(1)-1.5DO B(4)-X(3)-1.5DO
0011 0012
S(5)-X(4)-1.5DO B(b)=X(S)-1.500
0013 0014 0015
RETURN END
ooai 0002 0003 0004 0005 0006 0007
0008
SUBROUTINE IIlPLICIT DIMENSION C C C
FUNE(N,ll,X,B) DOUBLE PRECISION X(N),S(tl)
(A-H,O-Z)
Calculate the "equality" constraint function B(l)=X(l)+(X(2)r~2)+(X(3)~~3~-2.DO-3.DO~DSQRT~2.D0~ RETURN END
531
B. L. KARMALOO
538
and S. KANAGMUNDARAM
0001 0002 0003 0004 0005
0006 0007 0008 0009 0010 0011 0012 0013 0014 0015 0001 0002 0003 0004 0005 0006 0007 0008 0009 0010 0011 0012
SUBROUTINE IilPLICIT DItlENSION C C
C
Calculate
OBJF(N,X,F) DOUSLE PRECISION X(N)
objective
(A-H,O-2)
function
Al-X(I)-I.DO A29X(1)-X(2) A3-X(2)-X(3) A4-X(3)-%(4) A5-X(4)-X(5) F=(Alrr2)+(A2+r2)+(A3~~3)+(AS1+4) RETURN END SUBROUTINE ItlPLICIT DIMENSION C C
C 10
Assign
flOVE(N,DHIN,DHAX) DOUBLE PRECISION(A-H,O-2) DHIN(N),DllAX(N)
new move-limits
DO 10 I-I,N D~IN~I~-O.8D0~DHIN~I~ D~AX(1)-0.8DO~D~AX(I) RETURN END
OPTIMUM
~~.?~9~8? 1987 Pcrgamon
53.00 + 0.00
Joumalr
Ltd.
DESIGN OF STATICALLY INDETERMINATE BEAMS UNDER MULTIPLE LOADS B. L.
KARIHALOO and
S. KANAGASUNDARAM
School of Civii and Mining Engineering, The University of Sydney, Sydney, N.S.W. 2006, Australia
(Received 7 ~~e~~er
1986)
paper presents minimum-weight designs of statically indeterminate beams under multiple loads. Under each load system consisting of external loads and self-weight the normal and shear stresses in the beams are restricted from above. The variation of the stiffness along the span of the beams is
Abstract-The
restricted to splines of order zero, one or two, thus reducing the optimization problem to a linear or non-linear programming problem depending on the cross-sectional shape. The latter is reduced to a sequence of linear programming problems which are solved by the highly-efficient active set method. Several examples are given to illustrate the potential of this method. A computer code in standard
FORTRAN 77 is also provided.
The minimum-weight designs of efastic frames under multiple loads are known [l, 21. The designs meet stability [I] or stiffness criteria [2] but would have contained sections with vanishing stiffness at the positions of vanishing bending moment had the stiffness not been restricted from below. In contrast to the above maximum-stiffness designs, the maximum-strength designs [3-51 do not exhibit sections with vanishing stiffness provided that both the normal and the shear stresses are restricted from above. It should however be stressed that the maximum-strength design of a statically indeterminate flexural member or a frame is the least-weight design only if both the stress constraints are active. This is unlikely to be always the case. In fact, such a design will generally be understressed except at a few sections. The ideal maximum-strength designs (as, indeed the ideal maximum-stiffness designs) are also difficult to fabricate economically because of their awkward geometries which reflect the variations in bending moment and/or shear force. To avoid awkward geometry the span of a flexural member may be divided into several segments, in each of which the stiffness (mass) is permitted only a limited variation as specified by spfines of order zero, one or two. This approach has been taken in the present work. The member to be optimized is subject to multiple loads. Under each load system consisting of external transverse loads and self-weight the normal and shear stresses are restricted to certain allowable values. The corresponding problem under a single load system was recently reported [6]. The minimax optimization problem is first reduced to a mathematical pro~amming problem by repIacing the analytical expressions for the stress constraints and objective function with splines of desired
order. The resulting mathemati~l programming problem is linear or non-linear depending on whether the mass-stiffness relationship of the member is linear or non-linear. In the latter case the sequential linear programming technique with move limits [7, S] was found to be the most efficient solution method. Moreover, the active set method (91 for solving the linear programming problems was found to be far superior to the traditional (reduced) Simplex algorithm. Several examples are given to illustrate the potential of the present technique. A c.omputer code in standard FORTRAN 77 is also provided.
Consider an isolated elastic member in flexure (Fig. 1) that forms part of a multispan beam or frame. The member is subjected to reactive forces Pi (i= I,... ,4; axial forces have been neglected) from the joints of the beam/frame to which it is connected, as well as external transverse forces f(x) and selfweight g(x). The external forces f(x) consist of several individual force systems {f;(x), I = 1, . . . , NL ) each acting on the member on separate occasions. The equilib~um state of the member of given stiffness under any of the load systems is completely described by the two force resultants, namely, bending moment M,(x) and shear force Q,(x), together with the appropriate boundary conditions. The force resultants are related as follows:
dQ/
= dx in which member. 521
x refers
-lxx> -g,(x),
to an arbitrary
section
along
the
522
B. L. KARIHALOOand
S.
KASAGASUNDARAM
I
I
Fig. I. An elastic flexural member subjected to external load system f,(x) and reactive forces (P,, i=l,... .4) from joints. (Gravitational force due to self-weight is not shown.)
Throughout this paper, it will be assumed that the mass of the member, which is proportional to the cross-sectional area A(x), is related to its flexural stiffness, which is proportional to the second moment of area I(x) through I(x) = CA”(X),
(3)
in which c and n are constants determined by crosssectional shape. For most member shapes, the constant n takes an integer value between one and three. Thus, for instance, n = 1 represents a rectangular section of constant depth, n = 2, a geometrically similar section, say circular or square, and n = 3, a rectangular section of constant width. Here we report results for beams of rectangular section only (n = 1 or 3); the solution of the optimization problem with n = 2 follows along similar lines. The strength (stress) constraints on the member of known cross-sectional shape may be written in the form of the following differential inequalities [6]: ui{x, s, f,, g,, M,, Q,, A,, (A,),} G 0 (i=l,2),
(I=1 ,...,
NL),
(4)
in which x and s are the co-ordinates of the crosssectional fibre. The two stresses refer to normal and shear stresses, respectively. In the sequel, we will omit reference to the load system unless otherwise required. The optimization problem for a flexural member isolated from a multi-span beam (or frame) may now be stated as follows. Determine the control A (x) and, thus, I(x) that meets the strength requirements (4) (the differential inequalities) under all loading conditions and minimizes the mass, W, of the member given by
in which y is the mass density of the material. The optimization problem for every member of the continuous beam or frame may be stated in a similar way. When the compatibility conditions at the member connections (joints) are invoked, the problem for the continuous beam or frame as a whole is formulated. The differential game problem described by eqns (4) and (5) can, however, be reduced to a minimax variational problem. MINIMAX PROBLEM
Assuming that the equilibrium equations, eqns (1) and (2), with the appropriate boundary conditions, can be solved in a closed form, the expressions for the force resultants M(x) and Q(x) are known for each load system. It should be mentioned that M(x) and Q(x) depend not only on x, f and g but also on the control A(x) itself through the compatibility conditions and self-weight. The control A(x), in turn, depends on all the load systems. Substituting the (implicit) analytical expressions for M and Q, i.e. M=M(x,j;g,A,A,), Q=Q(~,frg,A4), into the differential inequalities [eqn (4)] gives, for each load system, ai]x* s, .L g, M(x, J g, A, A,),
Q(x,f,g,A,A,),A,A,l~O,
(64
or, for brevity, S-&(x,s,f, A, A,) < 0
(i = 1,2).
(6b)
It should be mentioned that as g(x) = yA (x) it has not been repeated separately in eqn (6b). Next, the maxima of CIi with respect to s and f are determined for fixed values of x and A (x). Let these be attained at s = s,’ and f =f ,‘, i.e. Di(x, stq f:, A, A,)
w=y
A (~1
dx,
(5)
= max m/“” f&(x, s, J A, A,), I
(7)
Optimum design of statically indeterminate
d,=O
beams
$+1-L
4
d2
523
la1
I-*
lb)
d,=O
d2
(cl
l--x Fig. 2. Representation
dN+t=L
di
dg
of cross-sectional
and denote
area A (x, pi) using splines.
(9) become
Wx, A, A,) = @lx, st(x, A, A,),ft(x),
A, A& (8)
In practice, as the number of load systems {i.e. NL) is discrete, the maximum with respect to f is to be chosen among a finite set. Using this notation, the differential inequalities (4) reduce to
f-&,*(X,&~>,<0(i=l,2),(k=I
,*.., N).
(11)
Thus, the original game problem [eqns (4) and (5)] reduces to the following variational one: Minimize the functional (mass) L
i2y(x,
A,A,)
(9)
Next, each member is divided into a given number of segments (NSEG) and the order of variation (NOR) of cross-sectional area A (x) is prescribed using splines. For example, splines of order zero are used to represent constant cross-sectional area within a segment and those of order one and two to represent linear and parabolic variations, respectively (Fig. 2). Irrespective of the order of spline chosen to represent A (x) in each segment, it will be of the form
Y A(x,pk)dx 50
(k=I,...,N),
(‘12)
subject to the stress inequalities (11). It should tie emphasized that the latter include a maximisation operation with respect to load systems. App~i~tion of the above reduction procedure is iflustrated on a dexural member in the next section. APPLICATION TO FLEXURAL MEMBERS
To iflustrate the concepts intr~u~ in the preceding section, consider an arbitrary, horizontal, A(x)=A(x,p,) k=l,...,N, (10) flexural member subjected to external loading, as well as reactive forces from the joints of the beam (or where pr: denote the unknown (discrete) design frame) of which it forms a part (Fig. 1). Let the variables, and N is the total number of variables; reactive for,rces due to load system f;(x) be P,i 4) and let the external transverse loading N = NSEG, if NOR =Q and N= NSEC + 1, if (i=I,..., NOR = I or 2. The design variables could be the consist of several individual load systems f;(x), I=1 ,...,NL. It is assumed that a typical load variable width b&x)@ = 1) at the knots (k = 1, . . . , N), or the variable depth /~(x)(n = 3) at the system I;(x) will consist of a distributed load q&), knots. Using the representation (W), the inequalities concentrated load Wfjlocated at positions uj along the
524
B. L. KARMALOOand
member length and concentrated couples C,j acting at locations b,, as shown in the figure. If the reactive forces P,i (i = 1, . ,4) are assumed to be known, the bending moments ~,(~~ and Q,(x) due to external load system h(x) are given by
z w,,
M,(x)=P,,x-
-cl,>
S. KANAGASUNDARAM b(x), the strength constraints with eqns (16) give 12 _max mFlM,(x)s(x)l ” = b(x)h’ s 6
I
62
1-1
(IS) in conjunction
=
1
h*
b(x)h’
my
(17a)
-co G 0
7-$2(X)
[
x maxIQ,(x)l - 7. G 0.
(17b)
f
-
O~hhhN~~l~~
-v)dtl
(13)
I
Q,(x)=
Pir- f W,,
in which m, and m, are respectively the total number of concentrated loads and couples acting on the members. The integral term involving the section area A(x) in each of the above two expressions is the contribution from the self-weight of the member. The strength constraints limiting the intensity of the normal (longitudinal) stress, cxxr and transverse (shear) stress, T__.at the section x are
The maxima of B, and a2 with respect to s are attained at s: = + h/2 and s: = 0, respectively, and so give 6 a: E-maxIM,(x)j-g,GO b(x)h2 / ny-_ - 3 Zb(x)h
(18a)
maxlQ,(x)l/
(18b)
r,G 0.
It is interesting to note that the differen~al inequalities (9) degenerate into algebraic inequalities (18) for the cross-sectional shape under consideration. Therefore, the lower bound on the optimum width of the member may be written, by inspection,
(19) where in which a, and to are known positive constants (allowable working stresses). From elastic beam theory it is known that
o(x) =max{lQ,(x)l,
Q2(xll,. . . * lQdxx>i>. CW
If it is assumed that the constant-depth beam consists of segments of linearly varying (NOR = 1) width, b(x) (Fig. 2b), it follows that for the ith segment, b(x)=(l-T)bi+
in which s(x) determines the position of the crosssectional fibre and K and a (x, s) are dependent on the cross-sectional shape. For a rectangular section of constant depth (h) and variable width b(x), i.e. for n = 1, K = l/~(~),u(x,~)=~(~)~~z/4-~z(x)]/2and -h/2 Q s G h/2, and for a rectangular section of constant width (b) and variable depth h(x), i.e. for n=3, K=l/b, a(x,s)=b[h2(x)/4-s2(x)]/2 and -h(x)/2Gs(x)dh(x)/2. In view of the dependence of a (x, s) and solution technique on the cross-sectional shape, i.e. on the value of n, it is convenient to illustrate the method separately for two typical values of n. Rectangular section of constant depth (n = 1) In the case of a flexural member of rectangular section with constant depth, h, and variable width
Z’bi+l,
dfGxGd,+,,
(21)
where T = T(x) = (x - di)/(di+ , - d,), and b, are the values of the optimum width at the knots. The stress constraints (19) for the segment d, < x g di+l can, therefore, be rewritten as , i=l
The functional
, . . . , NSEG.
(22)
to be minimized (12) becomes
W=~~(b,+b,+,)(d,+,-d,). 1-I
(23)
The problem of mi~m~ation of W (23) subject to the inequalities (22) is a linear pro~amming problem
525
Optimum design of statically indeterminate beams which is solved as follows:
where
1. Assume (bi)j, i = 1,. . . , N, in the first iteration 1. 2. Analyse the beam by standard stiffness matrix method (e.g. [3]) for all the load systems f;(x), and obtain n?(x) and e(x) [eqns (13, 14,20)]. 3. Divide each segment into m intervals, and solve the following linear programming problem using the active set algorithm [9]:
JI,= -
j=
Minimize
W=~‘~(~,+&+,)(CI~+, r-l
- 0
(24)
3
2bh2(x)
max Q,(x)h - M,(x); /
sf = 0,
(27)
lj12=3
sf2 = h’/4,
bh2(x)
G(x)= max(lM,(x)l,. The algebraic inequalities
,
(28)
. . , IMNL(x)l). (26a) reduce to
subject to [I - T(Q)]& +
T(Xk)Ei+
i=l,...,
h(x)>
I Z max(F,, F2), NSEG, k = 1,. . . , m + 1,
d,dxk
N.
(29)
The minimization of mass W [eqn (5)] subject to linear (29) and non-linear (26b) constraints results in a non-linear programming problem. It is however efficient to solve it as a sequence of linear programming problems as follows. First, the differential inequalities (26b) are rearranged such that if JI, > $2r we have
4. Repeat steps 2-3 until l(bi),+ 1-(b,),]g~(=0.5%),
i=l,...,
5. Calculate the minimum
mass (23).
LmaxlQ,h
N.
26P(x)
/
-M,$~-~o
or
The above steps were programmed on VAX 11/780 and several examples were solved. The results will be presented and discussed later in the paper. If the optimum beam is desired to have segments of constant (NOR = 0) or parabolically (NOR = 2) varying width, the optimization problem in step 3 above will still be a linear programming problem (see
fl = 3/(26r,). Similarly, if ti2 > I/,, we have
WI). (31)
Rectangular section of constant width (n = 3)
Substituting eqns (16) into the stress inequalities eqn (IS), we have 12 ul =
3
bh (x)
max m,a” IM,(x)s (x) ( - co G 0 s
(25a)
Assuming that the quantities within the square root sign in inequalities (30) and (31) are known, we may combine the inequalities (29)-(31) and write h(x)>max{F,,&F,},
(32)
where F,, F2 and F, are, respectively, the right hand sides of inequalities (29)-(31). Depending on the relative magnitudes of F,, F2 and F,, the inequalities (32) may be linear or non-linear. The minimization of mass The maxima of u, with respect to s are attained at s: = + h (x)/2, whereas that of c2 are attained at either sf2 = 0 or s:2 = h2(x)/4. The differential inequalities [eqn (7)] therefore reduce to *+
=6&l --
’Q: C.A.S I(,)_”
E
bh2(x)
00 Q cl
max($,. q2) - r. d 0,
(26d (26b)
L
W = yb
50
h(x)&
subject to the stress inequalities out in the following steps.
(33)
(32) is then carried
1. Assume h,(x) in the first iteration j = 1. 2. Analyse the structure by standard stiffness
8. L. KAIUHALCXJ and S. KANACASUNDARAM
526
BOkN
ZOkN/m Irrflrrtrlrr)
3m
6m 6m
I.-
3m
Load System 2 --
Load System 1 -PrlsmatlcBeam
0.5O.&030.2=---.,__ 0.1-
1 --____
00, 0.50Ab 0.3lm) 0.2\ 0.1O.OL 0
I 1
I 2
I 3
1 c
I 5
6
x Im)
Fig. 3. Minimum-weight design of pinned-clamped beam subjected to load systems as shown plus self-weight (member has constant depth h = 0.65 m, i.e. n = 1).
method [3] for all the load systemsfi(x), and obtain the force resultants M,(x) and Q,(x), (I = 1,. ..,NL), from eqns (13) and (14). 3. Solve the non-linear programming problem [eqns (32), (33)] using a sequential programming
approach with move-limits [7,8] to obtain h,, ,(x). 4. Repeat steps 2-3 until Jh,,, -h,)=26(=0.5%). 5. Calculate the minimum mass.
As with the case for n = 1, the optimum depth h(x) is approximated by splines of desired order [eqn (1O)]. It should be noted that if stiffness is constant within a segment (NOR = 0) then dh/dx = 0, and eqns (30) and (31) reduce to algebraic inequalities. The optimization probiem is linear in this case and step 3 in the above procedure reduces to the solution of a linear programming problem.
1.0PrismaticBeam
0.8-h
I
0.6-
______----
_,_---_
(ml 0.4__________________----------0.20.0.
1.0i h fm)
0.80.6-
-
O.L-u 0.20.0. 0
I 1
I 2
1 3
I Ir
I 5
x (rn)
Fig. 4. Minimum-weight
design of the beam shown in Fig. 3, but having rectangular cross-section of constant width b = 0.4m, i.e.n = 3.
6
521
Optimum design of statically indeterminate beams
30kN/m A GA
1
30kN
&IS
4m
6m
I
t---x 0
Load --
System
1
-Load
sPrismatic
0 .L-
Prismatic
Beam r
b o 3(ml 0 .2-
’
.*
System
2
Beam
.*
0 .l-‘, 0 .o0 Ab
o3 -
tm) 0.2 0 .I-\ 0 .o-
0
2
1
3
.4
5
6
7
9
8
x Im) Fig. 5. Minimum-weight design of two-span continuous beam subjected to load systems as shown plus self-weight (members are rectangular in cross-section with constant depth h = 0.65 m, i.e. n = 1). EXAMPLES AND DISCUSSION
The results of several numerical examples are shown in Figs 3-8. The following material properties were assumed in all the examples: u,, = 5 MPa, 7. = 0.5 MPa, E = 20 GPa and y = 2450 kg/m’. These values are typical for concrete. Data relating to the
O.E0.6-
Pritmatlc
Beam
I I _/__--- --
04 -
I A---__
external load systems acting on the beams are shown on the respective figures. Each member was designed to consist of l-m long segments, and the stress levels were calculated at sections spaced at 0.05-m intervals (i.e. m = 20). The integrals were evaluated using Simpson’s rule, and the iterations were terminated with c = 0.5%.
Prismatic
Beam
I 6
I 7
I --__ ---___
//* 0.2 I------L0.0, 0.6h lm)
0.4 0.2 = 0.0
0
I 1
I 2
I 3
I L
I 5
I 8
I 9
x Im)
Fig. 6. Minimum-weight design of the two-span continuous beam shown in Fig. 5, but having rectangular cross-section of constant width 6 s 0.4 m, i.e. n = 3.
10
528
B. L. KARIHALOO and S. KALAGASUNDAIUM
MkN
Load System
0.4
Prismattc
/’ -c,-’
oq-
Load System --
2
Beam
0.3 -1.
1
each
,’
/’
.’
0.4 b
0.3-
fm) 0.2 0.1 0.0
0
1
I 2
I 3
I 5
I I
I 6
I 7
I a
I 9
x (ml
Fig. 7. Minimum-weight design of two-span continuous beam subjected to load systems as shown plus self-weight. The concentrated forces in load system 2 represent a truck loading. The members are rectanguIar in cross-section with constant depth h = 0.65m, i.e. n = 1. Figures 3 and 4 show the minimum-weight design of propped cantilever beams of variable width (n = 1) and depth (n = 3), respectively. Figures 5-8 show the corresponding designs of two-span continuous beams. The load systems consist of dist~buted loads, concentrated Ioads and couples. The designs using
splines of different orders are compared. The designs using first and second order splines have been drawn separately since they would have been indistinguishable on the same diagram. A comparison of the volumes of op~mum beams is presented in Table 1. Also indicated are the volumes of the corresponding
0.6h
lm)
0.4 0.2 0.0
1 1
I
2
I
3
I
I
1
I
I
I
Ir
5
6
7
8
9
Fig. 8. Minimum-weight design of the two-span continuous beam shown in Fig. 7, but having a rectangular cross-section of constant width b = 0.4 m, i.e. n = 3.
Optimum design of statically indeterminate Table 1, Comparison
of volumes of optimum designs Volume (m’) Prismatic
Figure
3 (n = 1) 4 (R =3) 5 (n = 1) 6 (n = 3) 7 (n = 1) 8 (n = 3)
NOR=0
NOR-I
NOR=2
design
1.0695 1.1845
0.9733 1.1190
0.9996 1.1358
I .8639 1.8612
I .6943 1.8858 0.9449 1.2286
f.4411 I .6504 0.7864 0.8766
1.4149 I .6452 0.7913 0.8749
2.8037 2.7783 I so70
I .4933
prismatic designs. It should be mentioned that the total number of fully stressed sections in the minimum-weight beams is equal to the number of (discrete) design variables, the only exceptions being the beams of variable depth (n = 3) with first or second order splines in which fewer sections are fully stressed. All of the minimum-weight beams presented here satisfy the serviceability requirement in that the maximum deflection to span ratio is less than l/375. The non-linear programming problem associated with n = 3 was solved by a sequential linear programming technique with move-limits [7,8]. This technique was found to be the most efficient and reliable. In fact, the same optimum design was reached from different starting designs (both uniform and non-uniform), thereby pointing towards the global nature of the optimum. (The same statement is also valid for the optimum designs for n = 1.) The move-limits at a given iterative step were set equal to the 80% of the corresponding values at the previous step. The linear programs were solved by active set analysis [9j which was found to be very efficient and to use much less CPU time compared to the traditional approach using the (reduced) Simplex algorithm. A listing of the computer code in FORTRAN 77 is given in the Appendix. It should be mentioned that the sequential linear programming algorithm that was so efficient for the present problem may not be
529
beams
as efficient for solving other non-linear
programming problems. In fact, WChave found it rather unreliable for solving non-linear programming problems resulting from minimization problems under stress and deflection constraints. For these latter problems, the convex linearization and penalty function techniques seem to be more efficient and reIiabIe. This is currently under investigation and the results will be reported in a future communication.
Acknowledgement-The work reported in this paper was supported by grant F84/16078 from the Australian Research Grants Scheme to B. L. Karihaloo.
REFERENCES
I.
H.
K. Turner and R. H. Plaut, Optimal design for
stability under multiple loads. J. &ng AXE 2.
3.
4.
5.
6.
Mech: Die.,
106, 1365-1382 ~19801.
P. Pedersei and L. Jorgken:
Minimum mass design of elastic frames subjected to multiple load cases. Compur. Srrucr. 18, 147-157 (1984). S. Kanagasundaram, Optimum structural design (application of differential game theory). Ph.D. Thesis, The University of Newcastle, Australia (1984). S. Kanag~undaram and B. L. Karihaloo, Ootimal strength design of beam-columns. fnr. J. Sokk kr~. 19.937-953(1983h S. -Kanagasundaram and B. L. Karihaloo, Maximum strength design of structural frames. J. Strucr. Div., ASCE 111,1267-1287 (1985). B. L. Karihaloo and S. Kanagasundaram, Computeraided minimum-weight design of statically indeterminate beams. &grtg Opiimiz~f~on 10,139-156 (1986). P. Pedersen, The integrated approach of FEM-SLP for solving problems of optimal design. In Opfimizarion of Distribured Parameter Sysrems, Vol. I (Edited by E. 3. Haug and J. Cea), pp. 739-756. Sijthoff % Noordhoff, The Netherlands (I 98 1). P. Pedersen, On the optimal layout of multipurpose trusses. Comput. S~ruct. 2, 695-712 (1972). M. J. Best and K. Ritter, Linear Progrumming (Active Set Analysti and Computer Programs). Prentice-Hall, Englewood Cliffs, NJ (1985).
B. L. KARIHALOOand S. KANAGASUNDARAM
530
APPENDIX:FORTRAh 77 LISTINGOF SLP
0001
0002 0003 0004 0005 0006 0007 0008 0009 0010 0011 0012 0013 0014 0015 0014 0017 0018 0019 0020 0021 0022 0023 0024 0025 002L 0027 0028 0029 0030 0031 0032 0033 0034 0035 0036 0037 0038 0039 0040 0041 0042 0043 0044 0045 0046 0047 0048 0049 0050 0051 0052 0053 0054 0055 0056 0057 0058 0059 0060 OO6I 0062 0063 0064 0065 0066 0067 0068 0069 0070
SUBROUTINE
SLP(NV,tlINEQ,HEQ,DHIN,DMAX,TOL,IHAX,DELX, X,CC,IC,A,IROWA,6,C,DX,DINV,IROWD,S,V,WORK,JJ, AI,IROWA1,8I,Cl,DX1,DINVI,IROWDI,SI,VI,WORKI, JJ1,O8J,IOUT,JOlJT)
1 2 3
l
PURPOSE The program solves the following optimization problem: Hinimire
F
subject
to Ii (X ) j i
LE
G (X ) k i
EQ
* + 0 0 + + + t
non-linear
0,
0,
j-1
,..,MINEQ
k = i,..,HEQ
in which F is the objective function, H 's and G 's are the constraint functions k j and X , i=I,..,NV are the unknown design variables. i DESCRIPTION The program is based on Sequential Linear Programming (SLP) dpprOdCh using move-limits. In this approdeh the original problem is reduced to a sequence of linear programming problems with ohdngcs in the design variables (i.e., dX i i-1,. .,NV) ds the unknown vdridbles.The linear sub-problems are solved using active set strategy. The program is written in standard FORTRAN 77 for ease of portability. No derivdtivas are required. PARAMETERS
AND USER SUPPLIED
l
t + + t + t t + + + + t + + + + + + + + + + +
+ +
ROUTINES
4
I. Input
parameters
+ 4
NV -
INTEGER. NV specifies the number On entry, vdridbles. Unehdnged on exit.
+
of
design
4 + +
HINEQ
- INTEGER. HINEQ speoifies the number of On entry, inequality constraints. Unohdngcd on exit.
4 4
4 4
MEQ - INTEGER. On entry, MEQ specifies the number Unchanged on exit. constraints.
+
of
equality
4 4 +
DMIN
-
REAL array On entry,
of
DIMENSION GE NV. specifies DtlIN(I)(I-I ,..,NV)
lower move-limit design variable.
on the
Chdnge
the in the Ith
+ + l + l
DMAX -
REAL
array
of DIMENSION GE NV. On entry, DMAX(I)(I=I ,..,NV) specifies the upper move-limit on the ohdnge in the Ith design variable.
+ l 4 l +
TOL - REAL On entry, TOL must be set to specify the dccurdcy to which of change in the design variables the solution is desired. Unchdnged on exit.
l * + 4
+
Op(~~urn design of staucalfy jnde~errnIn~~e beams 0071 0072 0073 OQ74 0075 0076 0077 0076 0079 0080 0081 0082 0033 0004 0085 0006 0087 0080 0089 0090 0091 0092 0093 0094 0095 0096 0097 0098 0099 moo 0101 0102 0103 0104 0105 0106 0107 0108 0109
* *
IHAX
l
* * * *
-
l
* *
IBAX specifics the maximum to be performed. Unchanpcd
* +
number of on exit.
l
l
REAL On entry, DELX specifies an increment in design variables. Used in the calculation graditnts of the constraint and objtctive functions. For example,
the
l
af
* * l
*
l
dH(X
+ * *
-*-___-
* * * * * * * * + * * *
H
+DELX)
-
H(X
l
1
*
i i ________c_______I____
l
DELX
*
dX
l
i
l
d rtfcrs
whcrt
to
partial
l
derivative.
l
IOUT
-
*
INTEGER. on entry,
IOUT
associated
specifits
with
the
logical
output
l
unit numbtr Unchanged on
fife.
* * * *
exit. 2.
Input/Output
prramaters
l
REAL 8rray of l?fHENSiON Before entry, X(Il,i=i,..INV starting vrlues of the
X -
l
* * * * * * e * *
) i
l
On exit, X contains the dtsign VdridbieS
3.
Output 08J
*
GE NV.
must be set design
l
to
* *
vdridbles.
the optimum values
of
l
0 * * * *
l
parameter -
REAL. On exit 083 contains the objective function.
the
optimized
valut
of
4 * 4
l
*
+
oiii
*
0112 0113 0114 0115
* * * *
0116
l
4.
Workspaoe
(and
*
0119 cl120 012% 0122 0123 0124 at25 0126 0127 0123 0129 0130 0131 0132 0133
* * * * * * * * * * + * * * *
* * *
-
Dimension)
4
Parameters
INTEGER. On entry, IROWA must specify the first __ dimension of array A as declared in the CallanQ (sub) program. IROWA must bt GE (HINEQ + IlEQ +P*NVf.*
l 4 l
4
IAOWAl
-
INTEGER. On entry, IROWAI must dimension of array Al (sub) program. IROWAl + 2eNVl.
4
specify the as declared must
be
6E
first in the
4
callinq + FiEQ +
(XINEQ
*
1
l
* l
IRDWD
-
4
INTEGER. On entry,
lROWD must specify of array DINV as dtclared in
tht the CE
first calling CHINEQ
dimension
l l
+
l%EQ +2*NV).* l
IROWDi
-
*
INTEGER. On entry*
IROWDf must specify of array DINVl as declared in (sub) ‘program. XROWDl must be + 2*NVl.
the first dimension the calling 66
(HINEQ
+ UEQ
+
* l
1
l l 4
JJ,JJl
-
INTEGER PIWI
of
arrays
DIMENSION
* respectively,
#
GE NV dnd
+ 4
A,Al
-
l
* * * * * *
Associattd
l
IROWA
*
oiie
0134 0135 0136 0337 0138 0139 0140 0141 Oi42 0143
INTEGER On entry, iterations
* DELX
OlfO
0117
-
531
REAC arraym renpeotivtly,
of of
DIMENSION CIROWA,p) and where p must be GE NV.
fIRDWA,p+l)
4 4 4
B,Bl
-
REAL arrays and CHINEQ
DIHENSION + HEQ+l+?*NVl,
DINV, DINVl
-
REAI. arcrys of rtspeatfvely,whtre
GE (HfNEQ+MEQ+2*NV) respectively.
4 4 4
DItiENSION (1ROWD.p) and p must bt GE NV.
(IROWOi,p+l)
*
4
B. L. KARIHALOO and S. KANAGASLWDARAM
532 0144 0145 0146 0147 0148 0149 0150 0151 0152 0153 0154 0155 0154
0157 0158 0159 OibO 0161 0162 0163 0164. 0165 Oibb 0147 OIL0 0169 0170 0171 0172 0173 0174 0175 0176 0177 0170 0179 0100 0181 0182 0183 0184 Oi85 0186 0187 0188 0189 0190 0191 0192 0193 0194 0195 0196 0197 0198 0199 0200 0201 0202 0203 0204 0205 0206 0207 0208 0209 0210 02ii 0212 0213 0214 0215 0216 0217
a 4
a a + a a + a + a + * + + + + +
DX,C,S,V,WORK - REAL
arrays,
each
of
DItlENSION
GE NV.
DXl,C1,Sl,V1,WORK1 - REAL +rr+ys,
eroh
of
DIMENSION
GE NV+l.
INTEGER. entry, XC must specify the CC. XC must GE mrx(HINEQ,HEQ).
dimension
IC
-
cc 5.
on
-
REAL
User-supplied FUN1
-
a +
a + + + * + + + + a a a a + a
+rrry
of
Routines
(XC).
Example
et
end
of
proqr+m)
SUBROUTINE. This routine must be supplied to calculrte the values of all the “inequality” constraint functions for a given Set of design variables. The specifio+tion is :
NV -
MINEQ
-
X -
+
a + a + +
CI
FUNE
-
-
INTEGER. On entry, NV specifies design variables. Its be chrnqed in FUNI.
the value
number of must not
INTEGER. On entry, HINEQ specifies the number of inequality constraints. Its value must not 'be ohanqed in FUNI. REAL array of DItlENSION On entry, X contains the of the design variables to which the oonstraint FUNI must are required. values in X.
(NV). ourrent values corresponding function values not’ change the
REAL array of DIMENSION (HINEQ). FUN1 must set the elements of CI lqurl to the values of the “inequality” constraint functions corresponding $0 the design variables contained in X.
SUBROUTINE. This routine must be supplied to oaloulate the values of all the “equality” constraint functions for a given set of deriqn variables. The specification is I SUBROUTINE FUNE (NV,MEQ,X,CEl IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIHENSION X(NV),CEfflEQ) NV -
a +
a + + + + + + + + + a
(See
array
SUBROUTINE FUNI (NV,flINEQ,X,CIl IMPLICIT DOUBLE PRECISION (A-H,O-21 DIMENSION X(NV),CIfHINEQ)
4
a a + a a + + a + a a + + + + a + + +
DIflENSION
of
MEQ -
X -
INTEGER. NV specifies On entry, Its design variables. be oh+nqed in FUNE.
the value
number of must not
l
+ + +
+ + + + + + + + + + + + + 4 + 4 + + + + + + + + 4 + + + + a + l
I
I II . 1 t t 1 ,
4 4
I
XNTEGER. On entry, ME4 specifies the number of inequality constraints. Its value must not be changed in FUNE. REAL rrrry of DIMENSION X contains the On entry, of the design variables to whioh the constraint FUNE must l rc required. values in X.
a a + + + + + + + + + + + + +
(NV). current values oorrespondinq function values not change the
4
i :’
Optimum design of statically 0218 0219 0220 0221 0222 0223 0224 0225 0226 0227 0220 0229 0230 0231 0232 0233 0234 0235 0236 0237 0230 0233 0240 0241 0242 0243 0244 0245 0246 0247 0248 0249 0250 0251 0252 0253 0254 0255 0256 0257 0258 0259 0260 0261 0262 0263 0264 0265 0266
533
beams
*
l
CE -
l l
REAL FUNE
of set
array
must
DItlENSION
l
(HEQ).
the elements of CI equal to the values of the "equality" constraint functions corresponding to the design variables contained in X.
l l l
b 4l * l I)
l
OBJF
l l l l
-
SUBROUTINE.
l
Thie routine must be supplied to calculate the value of the objective function F(X) for a given set of design variables X. The specification is :
l
l
l l l
l
SUBROUTINE
l
IKPLICIT DOUBLE DIflENSION X(NV)
l
OBJF
(NV,X,FO) PRECISION
+
(A-H.O-2)
4 I)
l
*
l
NV - INTEGER. On entry, NV specifies the design variables. Its value changed in OBJF.
* l l
4t
number
of
must not be
l l l
*
l
x-
l
REAL
array
of DIIIENSION (NV).
On entry, X contains the of the design variables
l l l
to which the constraint rre required. OBJF must values in X.
l l
l
current
values correspondinq function values not change the
l * it * *
l
l
FO - REAL.
l
t
l
OBJF
l
the objective to the design
l
must
set
FO equal
to the value of function corresponding variables contained in X.
l
(I it I il
l l l
* l
IIOVE - SUBROUTINE. This routine calculates new move-limit9 on the changes in the design variables rcoording to user-specified criterion. It is not called in the firet itcrdtion.
* *
+ * * l l +
SUBROUTINE HOVE (NV,DHIN,DtlAX) ItlPLICIT DOUBLE PRECISION (A-H,O-2) DIMENSION DHIN(NV),DflAX(NV)
l l
l + *
l
*
NV - INTEGER. On entry NV specifies the number of design vrridbles. Its value must not be ohdnged in ROVE.
l l l l
4
DHIN - REAL array of DIMENSION (NV). On entry oontdins lower move-limits used in the previous iteration. On exit oontdins new lower move-limits calculdted dooordinq to user-specified LE 0, I-l,..NV. criterion. OKIN
l
0276 0277 0278 0279 0260 0201 0202 0283 0204 0285
l
0206 0207 0200 0289 0290
l
Best,
l
Aative
l l l l
DKAX - REAL rrrry of CIUENSION (NV). On entry oontdins upper move-limits used in the previous iterdtion. On exit oontdins new upper move-liritr OdIoUldted dooording to user-epecified GE 0, I=i,..~V. criterion. OKAX
l l l
l l
6. Other
routines
referenced
PIPZC,
l
* l l
l l l
* l l
l
The routines
PlPPC
dnd SOLN oan be obtdined
fron
* l
l
*
l
l
SOLN.
l
l
t
l
l
l
l
t 4
l
4
l
l
l
*
4 * 4 il
l
0267 0260 0269 0270 0271 0272 0273 0274 0275
CA.S x.l-4
indeterminate
tt.J., dnd Ritter, K., LINEAR PROGRAMINS: Set Analysis dnd Computer Programs, Englewood Cliff* New Jersey Prentice-Hdll InO., P-P- 254,283, 1985.
# l l
I) 4
534 0291 0292 0293 0294 0295 0296 0297 0298 0299 0300 0301 0302 0303 0304 0305 0306 0307 0308 0309 0310 0311 0312 0313 0314 031s 0316 0317 0318 0319 0320 0321 0322 0323 0324 0325 0326 0327 0328 0329 0330 0331 0332 0333 0334 0335 0336 0337 0338 0339 0340 0341 0342 0343 0344 0345 0366 0347 0348 0349 0350 0351 0352 0353 0354 0355 0356 0357 0358 0359 0360 0361 0362 0363
B. L. KARIHALOO Note that SOLN must parameters
4 l
it
and S.KANAGASUNDARAM
the be
argument list of routines PIPPC expanded td include the extra like IROWA,IROWAl,IROWD,IROUDl and
and
l
IOUT.
l
l
4t
l
*
7.
4 4 4 4 4 4 4 4
Error
Indicator JOUT
=
2
The from
linear subproblem is unboupded below.
JOUT
-
3
The linear subproblem has no feasible solution. Dr.
I
S.Kanagasundaram
C IllPLICIT
DOUBLE
PRECISION
(A-H,O-2)
C DIHENSION
DHIN~NV~,DflAX~NV~,CC~IC~IA(IROWA,NV~, 6(llINEQ+tlEQ+2+NV),C(NV),DX(NV),DINV(IRObJD,NV~, S(NV),V(NV),WORK(NV),JJ(NV),X(NV) Al(IROWAI,NV+1),8l(MINEQ+tlEQ+l+2*NV~,Cl(NV+l~, DXl(NV+1),DINVl(IROWDl,NV+l),Sl~NV+l~,Vl~NV+l~, WORK1 (NV+l) , JJl (NV+1 )
1 2 DItlENSXON 1 2 C
tiTOT=HINEQ+flEQ WRITE(IOUT,100)NV,klINEQ,tlEQ,TOL,ItlAX
ITE=l WRITE(IOUT,lIO)ITE PRINT 4, SLP_ITE-’ WRITE(IOUT,12O)(X(I),I=l,NV)
, ITE
Celculdtc move-limits IF(ITE.NE.I)CALL Calculdte functions
MOVE(NV,DMIN,DMAX)
values of 'inequrlity" constraint and their gradients
CALL FUNI(NV,HINEQ,X,CC) DO 2000 I=l,NV X(I)=X(I)+DELX CALL FUNI(NV,tlINEQ,X,Bl) X(1)=X(I)-DELX DO 2000 J=l,HINEQ
2000 3010
A(J,I)=(B?(J)-CC(J))/DELX DO 2010 i=l,HINEQ SUtl=O. DO DO 2020 J=l,NV
2020 2010 C C C
SUH-SUH+DHIN(J)tA(I,J) 8(L) --CC(L)-sull Impose
upper
limits on DX
HI =MINEQ+NV J-O DO 2030
I-HINEQ+l,Hl J-J+1 A(I,J)=l.DO 2030
B(I)=DflAX(J)-DHIN(J) tl2-tll +NV
C C C C
l
4
PROCRAWlER
4 4 LAST HODIFIED I 27 OCT. 1986 4 444444444444*444444444444444444444444444444444444+44444444444444444444444
C 1000
4 4 4 4 4 4
Include
non-negativity
J-0 DO
L=tll+l,fl2
2040
J-J+1
constraints
l
4 4 4
Optimum design of statically indeterminate 0364 0365 03cllr 0367 0368 0369 0370 0371 0372 0373 0374 0375 0376 a377 0378 0379 0380 0381 0382 0383 0384 0305 038& 0387 0380 0389 0390 0391 0392 0393 0394 039s 0396 0397 0398 0399 0400 0401 0402 0403 0404 0405 0406 0407 0408 0409 0410 0411 0412 0413 0414 0415 0416 0417 0418 0419 0420 0421 0422 0423 0424 0425 0426 0427 0428 0429 0430 0431 0432 0433 0434 0435 0436
2040 c c C c
CIlau1~t.L
values
funotionr
and
of their
“equality”
beams
ocnstraint
gradients
IF CMEQ.EQ.0) 60 TO 2050 CALL FUNE CNV,MEQ,X,CC) DO 2060 I-1,NV X(I)-XCI)*DELX CALL FUNE CNV,MEQ,X,BI> DO 2060 J-1 ,MEQ A~M2+J,I~-C0iCJ)-CCCJ~)/DELX DO 2070 I*l,MEQ SUM-O‘DO DO 2000 J-1,NV SUMmSUH+DMINCJ~~ACH2+I~J~ 6cN2+I)=-CCCI~-SUM
2060
20&O 2070 C C C 2050
Caloulete
of
coefficients
the
objective
function
CALL OB.YFCNV,X,FO) DO 2090 I*l,NV X(I)rXCf)+DELX CALL OBYF
2090 c c c
Solution
of
tlT=H2+MEP Nl-NV+1 fill-NT+1 MTOTi -Ml
the
linear
progrrm
1
c CALL
PlP2C
CALL
SOLN
I 2
CA,IROWA,B,C,DX,DINV,IRO~D,S,V,WORK,JJ,O0JS,NV, tl2,HEQ,MT,Al,IROWA1,Bl,Cl,DXl,DINVl,IROWDl, Sl,V1,WORKI,JJI,O0JJI,Nl,Mll,MTOTl,JOUT,IOUT~ CJOUT,IOUT,DX,NV,OBJJ)
Check if solution of thr linear program is feasible. If not, double the value of the move-limitr and repeat the iteration.
3000
C c c 3020
3030 3050
c C :040
IFCJOUT.EQ.3)THEN DO 3000 I-1 .NV DMINCI~-2.O~*DHINCX~ DHAXCi~-2.DO~DMAXCI~ 60 TO 3010 END IF IFCJOUT.EQ.2)THEN WRITEC7,lSO) RETURN END IF
Check for
convergence
DO 3020 I*l,NV DXCII-DXCI)+DHINCI) DO 3030 I-1 ,NV XF~DABSCDXCI)).GT.TOt)GO CONTINUE DO 3050 I-1,NV XCil=DXCI)+XCI) CALL OBJF
iterations
ITE-ITE+l IFCITE.GT.IMAX>THEN
TO
3040
535
B.L.KARIWLOCI and S.KANAGASLJNDARAM
536
0437 0438 0439 0440 0441 0442 0443 0444 0445 0446 0447 0448 0449 0450 0451 0452 0453 0454 0455 0456
0001
0002 0003 0004 0005 0006 0007 0008 0009 0010 0011 0012 0013 0014 0015 0016 0017 0018 0019 0020 0021 0022 0023 0024 0025 0026 0027 0028 0029 0030 0031 0032 0033 0034 0035 0036 0037 0038 0039 0040 0041 0042 0043 0044 0045 0046 0047
WRITE(IOUT,180) STOP END IF DO 3060 I-l,NV X(I)-X(I)+DX(I) GO TO 1000 FORllAT(//,2X,‘NUHBER
3060 100 1 2 3 4 110 120 130 140 150 160 170 180
OF VARIABLES -‘,14,/2X, ‘NUHSER OF INEQUALITIES -’ ,14,/2X, ‘NUMBER OF EQUALITIES =‘,14,/2X, ‘CONVERGENCE CRITERION =‘,D12.6,/2X, ‘HAXIHUH NUMBER OF ITERATIONS -‘,X4) FORhAT(/SX,‘SLP_LTE -‘,X3) FORtlAT(/SX,‘VALUEB OF VARIABLEB’,/,/5fSX,DI2.6)) FORtlAT(/SX,‘LOWER HOVE-LIHITS ON VARIABLEB’,/,/5(5X,Dl2.6)) FORtlAT(/SX,‘UPPER HOVE-LIMITS ON VARIABLES’,/,/S(SX,Dl2.6)) FOR?lATf/2X,‘*** PROBLEPl IS UNBOUNDED FROM BELOW +**‘) FORtlATf/2X,‘OPTIHUH VALUE OF VARIABLES ..‘,/5(2X,Dl2.6)) FOR?‘lAT(//2X,‘OBJECTIVE FUNCTION VALUE AT OPTIliUt! =‘,D12.6) FORHAT(/PX,‘*** IlAXIHUfl NUHBER OF ITERATIONS PERFORtlED END
PROGRAM MAIN C C Calling program for routine C **************t********************************~*************************
SLP.
4
4 4
EXAHPLE
PROBLEM
4
4
4
Hinimirr
F
-
(X
-
+
1.0)
+
-X)
(X
1
4
3
2
2
l 4
(X
2
1
-
x
) 3
2
+
4
ix
4 4 4
-
x
3
4
4 4 4
4
4
l
4 4 4 4 4 4 4 4 4 4 4 4
)
+
-
cx 4
4
x
4 4
4
)
4
5
4 4
subjcot
to
4 4
2 X 2
4
+x
-x 3
+
2.0
-
2.0
*
SQRT(2.0)
GE
0.0
4
4 4 4
X
l
x
1
LE’2.0
4
5
4 4
X
LE
I.5
I
4 4 4
LE
X
1.5
3
4 4 4
LE
X
1.5
4
4 4 4
X
LE
4 5 4 4 3 2 4 l SQRT(2.0) +x - 2.0 - 3.0 X +x 4 3 1 2 4 4 This lx~mpla is taken from 4 Nummrioal Algorithms Group, NAG Library Manual, 4 1985. England, 4 *4444444*44444444*444444*44444*444444444*44444*444*444444444*4444*4444444
C
4
4
4
4 4
***‘I
1.5
4
4 4 4
0.0
4 4
4 4 4 4 4
Optimum design of statically indeterminate IRPLICIT
0048 0049
0050
0051 0052 0053 0054 0055 0056 0057 0050 0059 0060 0061 0062 0063 0064 0065 0066 0067 006.9 0069 0070 0071 0072 0073 0074 0075 0076 0077
c c C
Dlrfinr
problem
PARAHETER DIflENSION
DIMENSION 1 2 C C
Open
PRECISION(A-H,O-2) parameters
and
arrays
(NV=5,HINEQ=b,HEQ-1,IClb) DMXN(NV),DtlAX(NV),A~HINEQ+MEQ+2*NV,NV~, B~HINEQ+MEQ+2+NV~,C~NV~,DX~NV~,DINV~HINEQ+HEQ+2*NV,NV~ ,S(NV),V(NV),WORK~NV),JJ(NV),X(NV),CC(IC~ Al(flINEQ+HEQ+l+2+NV,NV+l~,Sl~HINEQ+nEQ+I+2*NV~, C1(NV+l),DXi(NV+1),DINVl(nINEQ+HEQ+l+2+NV,NV+l~, S1(NV+I),VI(NV+l),WORKI(NV+I),JJl~NV+I)
I 2
output
file in unit IOUT
10UT-7 OPEN~lJNIT=IOlJT,FILE=‘SLP.RES’,STATUSr’NEW’~ C C
Assign
initial move-limits to variables
C DO IO I=l,NV IO
OHIN(I1.00 DHAX(I)=l.DO
C
Define other
parameters
C TOL=l.D-04 IMAXDELX=I.D-05 IROWA=HINEQ+llEQ+2~NV IROWAl=HINEQ+HEQ+1+2*NV IROWD=HINEQ+tlEQ+2*NV IROWDl=HINEQ+HEQ+1+2+NV
0078 0079
ooao ooal 0082 0083 ooa4 0085 0086 0017 0088 0089
DOUBLE
beams
c C C
Starting values for the
variables
X(1)=0.500 X(2)=1 .ODO X(3)-0.6DO X(4)-1.4DO X(51-1 .ODO c 1 2
0090 0091 0092 0093 0094
3
CALL SLP(NV,HINEQ,HEQ,DHIN,DHAX,TOL,IllAX,DELX, X,CC,IC,A,IROWA,B,C,DX,DINV,IROWD,S,V,WORK,JJ,
A1,IROWAl,BI,C1,DXl,DINVI,IROWDI,Sl,Vl,WORKl, JJI,OSJ,IOUT,JOUT) STOP END
0001 0002 0003 0004 0005 0006 0007
SUBROUTINE ItlPLICIT DIHENSION C C
FUNI(N,li,X,B) DOUBLE PRECISION X(N),B(tl)
(A-H,O-Z)
Calculate "inequality" constraint functions
C
oooa
B(1) --X(2)+(X(3)**2)-X(4)-2.D0+2.DO+2.DO*DSQRTf2.DOy
0009 0010
BfP)=X(I)*X(S)-2.DO 8(3)-X(1)-1.5DO B(4)-X(3)-1.5DO
0011 0012
S(5)-X(4)-1.5DO B(b)=X(S)-1.500
0013 0014 0015
RETURN END
ooai 0002 0003 0004 0005 0006 0007
0008
SUBROUTINE IIlPLICIT DIMENSION C C C
FUNE(N,ll,X,B) DOUBLE PRECISION X(N),S(tl)
(A-H,O-Z)
Calculate the "equality" constraint function B(l)=X(l)+(X(2)r~2)+(X(3)~~3~-2.DO-3.DO~DSQRT~2.D0~ RETURN END
531
B. L. KARMALOO
538
and S. KANAGMUNDARAM
0001 0002 0003 0004 0005
0006 0007 0008 0009 0010 0011 0012 0013 0014 0015 0001 0002 0003 0004 0005 0006 0007 0008 0009 0010 0011 0012
SUBROUTINE IilPLICIT DItlENSION C C
C
Calculate
OBJF(N,X,F) DOUSLE PRECISION X(N)
objective
(A-H,O-2)
function
Al-X(I)-I.DO A29X(1)-X(2) A3-X(2)-X(3) A4-X(3)-%(4) A5-X(4)-X(5) F=(Alrr2)+(A2+r2)+(A3~~3)+(AS1+4) RETURN END SUBROUTINE ItlPLICIT DIMENSION C C
C 10
Assign
flOVE(N,DHIN,DHAX) DOUBLE PRECISION(A-H,O-2) DHIN(N),DllAX(N)
new move-limits
DO 10 I-I,N D~IN~I~-O.8D0~DHIN~I~ D~AX(1)-0.8DO~D~AX(I) RETURN END