- Email: [email protected]
Optimum Design in Distributed Parameter Systems
Copyright © IFAC 3rd Symposium
SESSION 14 -
OPTIMAL DESIGN (2)
Contro l of Distributed Paramet e r Systcms
T ou louse. F ra nee . 1982
OPTIMUM DESIGN IN DISTRIBUTED PARAMETER SYSTEMS M. Koda Department of A eronautics, Faculty of Engineering, University of Tokyo, Hongo , Bunkyo-ku, Tokyo 113, japan
Abstract. A new numerical method for the computation of sensitivity of design performance functional with respect to boundary variations is used to study the optimum design problems in general distributed parameter systems. The method is based on a modification of schemes in functional derivative sensitivity analysis and includes a specialization of the optimal control algorithms in distributed parameter systems. A design sensitivity theorem is derived for the determination of the optimum boundary location that minimizes a design performance functional . Necessary conditions for minimum drag profiles in unsteady Navier-Stokes flows are obtained and a numerical desi gn algorithm based on the gradient method is formulated. Keywords. Distributed parameter systems; optimum design; sensitivity analysis; boundary control; computational method; optimisation.
I NTRODUCTI ON
is defined on a two-dimensional spatial domain. The extension of the results to higher dimension s is straightforward.
De s i gn of di s tributed parameter systems can be often characterized as free boundary problems in which part of the boundary, the free boundary, i s unknown and must be determined as a part of the problem. Pironneau (1973, 1974 ) worked on optimum (minimum drag) design problems in a st eady viscous fluid by utilizing the optimal control theory in distributed parameter systems (Lions, 1971). In Pironneau (1973, 1974), however, no comprehen s ive optimization theory wa s derived and applications of the theor y to the full (timedependent ) ~av i er - S tokes equations remain to be studied.
We consider the fluid mechanical system Ut = f(x,y,t,u, Vu, V2u)
(1)
where u(x,y,t) is the n-dimensional state vector defined on the spatial domain Q c R2 , and f(x,y,t,u, Vu,V 2u) denotes f(x,y,t,u,ux'u y , u
A numerical optimi zation theor y for the estimation of unkno~TI boundaries in general di s tributed parameter sys tems was developed by Chen and Seinfeld ( 19 75) . In this paper, we extend that ana l ysis to related problems of optimum desi gn in f luid mechanical distributed parameter sys t ems . Thu s , the purpose of this paper is to formulate in a comprehens ive manner and deve lop computational algorithms for the de rivation of the optimum shape or the optimum location of the boundar y in a cla ss of general distributed parameter systems. OPTnIU~1 DES IG\ PROBLHIS 1\ DISTRIBUTED PARAMETER SYSTEMS
In thi s paper, we confine our attention to optimum design problems in nonlinear systems of the first-order hyperbolic or parabolic partial differential equations and restrict our attention to the case in which the system 353
,u ) where the partial derivatives are xx yy indicated by a subscript, i.e., Ut' u ' and x u denote 3u/3 t, 3u/3x , and 32u/3x 2 , respecxx tively. The initial and boundary conditions to (1) are given in the general form u(x,y,O) = uo(x,y),
where u
n
x, )'
E Q
x, Y
E
3Q
( 2) (3)
denotes the derivative of u in the
direction normal to the boundary 3Q of Q. The system evolves from t=O to t=T, the period in which we are interested. In our problem, the domain Q E R2 is bounded by two surfaces rand S ( r , S E 3Q ); the outer surface r is known and the inner surface S is unknown and is to be determined optimally. Let V denote an area surrounded by S (see Fig. 1). Then the optimum design problem is to determine the coordinate location of each point of the boundary 3Q=S such that the design performance functional J defined as follows:
M. Koda
354
I:I
J
D(U,IJU,1J 2U)dXd Ydt Q
(4)
is minimized, and where D(u,lJu,1J 2u)=D(u,u x ,u y , Uxx,Uyy»O.
The problem may be viewed as an
optimal control problem of a class of distributed parameter systems, ~here the (open loop) control is the location of the boundary.
relate 6 ( oQ) to a variation in the design performance functional, 6J . A The variation 6 ( oQ ) transforms local coordinate system in the original domain Q to the global coordinate system in the new perturbed domain ~ (coordinates denoted by ~,~). The essence of the present variational approach is to find a new statement of the problem by suitable coordinate transformations. Suppose that ~ is mapped onto Q by transformations ~
(x,y,u, lJu,1J2 u ;0)
(6)
~ (x,y,u, lJu, 1J 2u; 0)
(7)
1 ~
2
where 0 is a parameter of the transformations and the new value of the state u = u(~,~,t) is mapped onto the original state function u(x,y,t) by the transformations (6) and (7) and the transformation u = ~ (x,y,u, lJu,1J 2u ; 0).
(8)
3
We assume that the transformations (6), (7), and (8) are one-to-one continuously differentiable and the value 0=0 corresponds to the identit y map:
n Fig. 1. Problem geometry.
x =
We transform the constrained optimal control problem presented here into an unconstrained problem through the introduction of adjoint functions. We define the augmented design performance functional by J
ITf ~T[f(x,y,t,u,lJu,1J2u)-u o
Q
]dxdydt (5) t
where ~(x,y,t) is the n-dimensional adjoint function and superscript T denotes transpose. Then if 50 is the boundary that minimizes J, it also clearly minimizes J . A In our approach to this optimization problem, we use a variational method on a variable region (Gelfand and Fomin, 1963). In particular, we wish to derive a design sensitivity of J to variations in the boundary location by A using the functional derivative sensitivity arguments (Koda, Dogru, and Seinfeld, 1979; Koda and Seinfeld, 1982). The constructive approach presented in this paper is based on the idea of design sensitivities to seek a design improvement by decreasing the design performance functional as much as possible. VARIATIO ~ AL
APPROACH
We consider a variation in the boundary location, 6(oQ)=6S, which transforms the boundary oQ to
Y
~
u =
~
(x,y,u,lJu, 1J 2u;0)
(9)
(x,y,u,lJu,1J 2u;0)
(10)
(x,y,u, lJu, IJLu;O) .
(11 )
1
2 3
The first variations of x, y, and u are defined as
A=J +
~
0 ~=o Q + 6 (o Q ) ~ith
a consequent mapping of
the domain Q into ~ = Q + 6Q , the domain corresponding to o~. The variation 6(o Q) introduces a variation in J . Our goal here is to A derive a design sensitivity 6J /6 (o Q), or to A
a
2 .
6x
~-x
u,O)
(12)
6y
~-y = o~o ~ 2 (x,y,u, lJu, 1J 2U;0)
(13)
6u
u(e~, t )
o~ ~ (x,y,u,lJu, 1J 00 1
u(x,y,t)
o ~o ~3 (x,y,u, lJu, 1J 2U;0)
(14)
where we have used the Taylor expansions of (6), (7), and (8), assuming that 0 is a small quantity. We may note that (14) can be expressed as 6u = 6u + u x 6x + UY6",
(15 )
where 6u=u(x,y,t)-u(x,y,t) is a perturbation of the state function in usual variational sense. Equation (15) implies that the total variation 6u is decomposed into a contribution due to the variation in the state, 6u, and into a contribution due to the coordinate transformations, u 6x+u 6y. The Jacobian of the transformationxis Yapproximated as
l o ( ~ ,~ )
l o (x,y)
I =1
+
~x OX
+
o
~y.
(16)
The first variation of the design performance functional J is defined as the principal A linear part (relative to 0) of the difference
Optimum Design
~n
Distributed Parameter Systems
J:J~{O(U'VU'V2U)
oJ A
+ tjJT (f;,ll,t) [f(f;,T),t,u,vu,V 2u)-u ) }df;dT)dt t -
J:J~{O(U'VU'V2U)
+ tjJT(x,y,t) [f(x,y,t,u,vu,V 2u)-u )}dxdyd P .(17) t
355
follow some computational aspects of this optimal design problem. The change of the design performance functional due to a shift of a boundary node tangential to a~=s is certainly smaller than the change of the design performance due to a normal shift. Thus, we shall restrict our arguments to the boundary variation (control) which is normal to the boundary surface S. We therefore specify
We specify that tjJ(x,y,t) satisfy the following adjoint equation:
ox (s)
on(s)cos8
(22)
ay (s)
on(s)sin8
(23)
= _fTtjJ+(fT tjJ) +(fT tjJ) _efT tjJ) _efT tjJ) where n(s) is the normal to a~ at s. u U x u y u xx Uyy yy xx x y For the purposes of numerical implementation, _OT+(OT ) +(OT ) _(OT) (OT) OU, ou , and ou in (20) should be approxiu u x u y u xx- u yy (18) x y x y xx yy mated by ou and ou. Expressing OU and ou n x y with the terminal condition in terms of OU and OU , the normal and tann s tjJ(x,y, T) = O. (19) gential derivatives of OU on S respectively, and noting that S is closed, (20) becomes tjJ
t
Making use of (12)-(16) to transform the domain of integration and employing the Green's theorem, (17) can be transformed to
=
+(Hu
ITJ o
TJ [-O+(Ecos8+Fsin8)u Jo a~ n xx
[O(oxcos8+oysin8)
cos 28+H u
+H
u
u
- (H x
u
) }oucos8+{H - (H ) }Ousin8 x u U y xx y yy
OU cos8+H xx
sin 28)u nn
a~
-H +{H
yy
x
U
OU sin8)dsdt yy
y
TJ
(20)
Jo as"2
where 8 is the angle between the positive xaxis and the outward normal at a point of the boundary a~, and s is the arc length along a~ (see Fig. 2). In (20), the Hamiltonian is defined as
u
)sin8cos8u yy
sn
)ondsdt
[Ecos8+Fsin8
~(H as u
-H xx
u
)sin8cos8})oudsdt yy
H(x,y,t,u,Vu,V 2u) = O(u,Vu,V 2u) +tjJ T (x,y,t)f(x,y,t,u,Vu,V 2 u)
+H
(21)
where
Hu = aH/au x ' Hu = aH/au xx ' etc. It is x xx important to note that (20) need not follow the boundary of the unknown perturbed domain ~.
u
sin 28)ou dsdt n
yy
(24 )
where
E=Hu -(H u )x and F=Hu -(H u )y' x xx y yy In deriving (24) we have used the relations similar to (15), i.e., ou=ou+unon, ou n =oun +unnon, etc.
n
~IAI:\
THEOREM - OESIGN SENSITIVITY
We now choose the appropriate boundary conditions for the adjoint system (18) and (19). From (3) we have x, y£3S"2
0,
an Fig. 2. Coordinate system at boundary. Since this paper intends to give mathematical background for a numerical solution, we now
When
~
gu
(25 )
0, the boundary condition is given
n
by Ecos8+Fsin8+
-H )sin8cos8) ~(H as u u xx
yy
M. Koda
356
cos 2 6+H
u
sin 2 6)g yy
-1
u
n
g = O. u
(26)
depending on the appropriate form of the boundary condition (3), where the adjoint function ~ is defined by (18) and (19) with the boundary condition given by (26) or (28), respectively.
For this case, (24) becomes =
I {fT[ an
-D+(Ecos6+Fsin6)u
0
- (Hu
n
-Hu xx
)sin6coS8usn]dt}on(s)ds. yy
(27) When the boundary condition (3) is not an explicit function of un' i.e., when g-1 does not exist, we have the boundary condiu
n
tion as follows: H u
cos 2 6 + H
u
xx
( 31)
-(H u -Hu )sin 6cos 6u sn ]dt xx yy
O.
(28)
yy
Remark 1. The functional derivative sensItIvities (30) and (31) can be viewed as the sensitivity of the design performance functional with respect to a variation in the normal coordinate of the boundary location at each position s per unit arc length. Remark 2. The functional derivative sensitivities (30) and (31) are actually sensitivity density functions since appropriate surface integrals (27) and (29) are needed to relate them to a total variation oJ . A Remark 3. In our optimum design problems, the functional derivatives (30) and (31) are essentially equal to variational derivatives.
For this case, (24) becomes
=f
an
{fT[
-D+(Ecos6+Fsin6)u
0
n
-H
u
)sin6cos6u yy
sn
]dt}on(s)ds.
(29) In (27) and (29), the variation in boundary coordinates, ox and oy, has been replaced by ones), the variation of the normal to the boundary surface S (where the direction angle of the normal is redifined by 6+TI, i.e. outward normal with respect to the domain V, see Fig. 1). The relations (27) and (29) naturally lead to the consideration of functional derivative sensitivities (Koda, Dogru, and Seinfeld, 1979; Koda and Seinfeld, 1982). The analysis can be summarized in the following design sensitivity theorem. Theorem. The functional derivative sensItIvIty oJA/on(s), the sensitivity of the design performance functional (4) with respect to the variation in boundary location o(an) = ones) to satisfy the state equations (1), (2), and (3), is given by
T[ -D+(Ecos6+Fsin6)u
fo
cos 2 6+H -(H
u
u
sin 2 6)u yy
n
In this section, the derivation of necessary conditions for minimum drag profiles in full (time-dependent) Navier-Stokes flow problems is presented. In an earlier work of Pironneau (1974), a particular minimum drag problem in a steady viscous fluid has been studied. Since there exist few results of a general nature for this class of optimum design problems, it is very instructive to reduce our general results to a familiar earlier case. We consider the two-dimensional Savier-Stokes equations Ut
u
(30)
or n
t
uu x - vu y - px + v (u xx +u yy )
(32)
= - uv x - vv y - py + v( v xx +v yy )
(33)
and the continuity equation
nn
To [ -D+(Ecos6+Fsin6)u
NECESSARY CONDITIONS FOR MINIMUM DRAG PROFILES
v
-Hu )sin6cos6u ]dt sn xx yy
f
It is clear that these design sensitivities inherently contain most useful information to solve the optimum design problems in distributed parameter systems. We should note that we do not consider any design constraints so far. Such design constraints, in addition to state constraints, arise in most applications and must be treated.
x
+
vy
=0
(34)
where u and v are the fluid velocity components in x- and y-coordinate directions, respectively, p denotes the pressure, and v is the viscosity of the fluid. The boundary conditions are (see Fig. 1) u
= v = 0,
x, YE: S
av
(35)
Optimum Design in Distributed Parameter Systems
v = 0,
x, ye: r
357
(36)
where Uoo denotes the fluid velocity at infinity
a(~,n)
Iy l a(x,y)
and the appropriate initial condition of the form (2) is assumed.
IdXdY -
Jydxd y
We choose the design performance functional for minimum drag profiles as J =
I ITI°
~
[2 (u 2+v2)+(u +v )2]dxdydt. x y y x
(38)
Thus, the problem is the design of a minimum drag body of given volume in two-dimensional unsteady flow fields. Applying the results obtained in the previous sections, we have the following set of linear partial differential equations for the adjoint functions wand w . 1
(oxcos9+oysin9)ds
a~
This is the total energy dissipated in an incompressible fluid in the time interval [O ,T]. Our problem here is to determine the coordinate location of the boundary of the domain Y such that (37) is minimized subject to the design constraint = IydXdY = const.
I
(37)
IsOn(S)dS where we have used (16), (22), and (23). From (31) we compute the functional derivative sensitivity density function as oJ / ones) V
-2 I
T
°[U n2+V n2+2{(W In)un +(w 2n)vn }]dt+A
where A is the Lagrange multiplier associated with (38). At the optimum, the condition oJA/on(s)= must hold; this implies
°
Corollary l,where
(w ) = u w +v w -u(w ) -v(w ) -uu -vu It xl X2 IX lY x Y -v[(w) 1
xx
+(w) 1
yy
]+q
(39)
x
In the same manner, we can obtain the necessary condition for a minimum drag body of given surface area subject to the design constraint const.
(w ) = u w +v I'i -u(w ) -v(w ) -uv -vv 2 t Y1 Y2 2 x 2 Y x Y -v[ (w) 2
(w ) +(w ) = 1
x
2 Y
xx
+(w) 2
yy
] +q
(40)
Y
°
(41)
w (x,y,T)=w (x,y,T)= 2
°
w (X,y,t)=w (X, y ,t)= 0, 1
2
x,ye:a ~ =
s,r
( 45)
Corollary 2. The necessary condition for the minimum drag profile given a unit surface area is T p(S)I [U 2+V 2+2{(W ) u +(w ) v }]dt n n 1 n n 2 n n
°
with the terminal and boundary conditions 1
~=-2A/V .
2
( 42)
const.
(43)
almost everywhere on S, where pes) is the radius of curvature of S at a point s, and the adjoint system is as in Corollary 1. Proof. The extension in n of the design constraint (45) can be written as
where we have used the boundary condition (28).
( 46)
Using the sensitivity results and the main theorem, we obtain the necessary condition for optimality for the problem (32)-(38) as follows . Corollary 1. The necessary condition for the minimum drag profile given a unit volume is that there exists a constant ~ such that T I
°[Un2+Vn2+2{(W)1 nun +(w 2) nvn }]dt = ~
almost everywhere on 5, where the adjoint system is defined by (39)-(43). Proof. The first variation of the design constraint (38) can be computed as
(44)
Thus we have o£ =
I
on(s)ds sP (s)
where the relation 1/p(s)=d9/ds is used. With the help of the main theorem, the rest of the proof proceeds as in Corollary 1. Remark 4. In steady-state case, the necessary conditions (44) and (46) reduced to those of
M. Koda
358
implemented for the necessary condition (46).
Pironneau (1974). NU~IERICAL
CONCLUSIONS
IMPLEMENTATION
The design sensitivities (44) and (46) are derived as functionals of the system state and adjoint functions. Thus, simultaneous solution of the system state and adjoint equations facilitates the computation of infinite dimensional gradients (functional derivative sensitivity density functions) of the design performance functional and a modification of the gradient method can be used as an iterative unconstrained optimization technique. An iterative numerical algorithm to compute
the optimum location of the boundary is derived directly from the necessary condition (44). The basis of the method is to choose ones) such that oJ is always negative. This is A accomplished by setting ones) as ones) = -a(s){
2 +v2 IT[U o n n
+2{(w ) u +(w ) v }]dt lnn 2nn
-~}
(47)
where a(s) is an arbitrary (small) positive weighting function. In practical numerical implementations, ~ in (47) is replaced by the mean value of
T
Io
[U 2+V 2+2{(W) u +(w ) v }]dt n n 1 n n 2 n n
( 48)
on S at each iteration. The algorithm proceeds as follows: Step 1.
Make an initial guess of the boundary location.
Step 2.
Solve the Navier-Stokes equations (32)-(36) from t=O to t=T. Evaluate J.
Step 3.
Solve the adjoint equations (39)-(43) from t=T to t=O.
Step 4.
Choose an appropriate weighting function a(s) in (47) and compute the normal shift ones) by using (47) and (48).
Step 5.
~lodify
Step 6.
Return to Step 2. Repeat until subsequent changes in J are less than a preset criterion.
the boundary location.
Detailed applications would, of course, require numerical integration of the relevent timedependent Navier-Stokes equations and reasonable numerical schemes developed by Fortin and Thomasset (1979) may be used in the present algorithm. Since the present method only involves the iteration of the normal shift in boundary location, the algorithm is considerably simpler and numerically advantageous than the other algorithms. We remark that a similar algorithm can be effectively
The optimum design problems in distributed parameter systems have been considered from an optimal control standpoint. Necessary conditions for minimum drag profiles in unsteady Kavier-Stokes flows have been derived and the results are of interest in their own right. In particular, it has been shown that the general results of the paper, i.e. design sensitivity theorem, can be reduced to the familiar purely steady cases of Pironneau (1974). The present method, based on the functional derivative sensitivity analysis, is clearly applicable to a broad class of optimum design problems that are described by free boundary problems. REFERENCES Chen, W. H., and J. H. Seinfeld (1975). Estimation of the location of the boundary of a petroleum reservoir. Society of Petroleum Engineers Journal, 15, 19-38. Fortin, M., and F. Thomasset (1979). ~lixed finite element methods for incompressible flow problems. J. of Computational Physics, 31, 114-145. Gelfand, I. M~ and F. D. Fomin (1963). Calculus of Variations. Prentice-Hall, Englewood Cliffs, New Jersey. Koda, M., A. H. Dogru, and J. H. Seinfeld (1979). Sensitivity analysis of partial differential equations with application to reaction and diffusion processes. J. of Computational Physics, 30, 259-282. Koda, M., and J. H. Seinfeld (1982). Sensitivity analysis of distributed parameter systems. IEEE Trans. on Autom. Control, in press. Lions, J. L. (1971). Optimal Control of Systems Governed by Partial Differential Equations. Springer-Verlag, New York. Pironneau, O. (1973). On optimum profiles in Stokes flow. J. of Fluid ~lechanics, ~, 117-128. Pironneau, O. (1974). On optimum design in fluid mechanics. J. of Fluid Mechanics, 64, 97-110.
SESSION 14 -
OPTIMAL DESIGN (2)
Contro l of Distributed Paramet e r Systcms
T ou louse. F ra nee . 1982
OPTIMUM DESIGN IN DISTRIBUTED PARAMETER SYSTEMS M. Koda Department of A eronautics, Faculty of Engineering, University of Tokyo, Hongo , Bunkyo-ku, Tokyo 113, japan
Abstract. A new numerical method for the computation of sensitivity of design performance functional with respect to boundary variations is used to study the optimum design problems in general distributed parameter systems. The method is based on a modification of schemes in functional derivative sensitivity analysis and includes a specialization of the optimal control algorithms in distributed parameter systems. A design sensitivity theorem is derived for the determination of the optimum boundary location that minimizes a design performance functional . Necessary conditions for minimum drag profiles in unsteady Navier-Stokes flows are obtained and a numerical desi gn algorithm based on the gradient method is formulated. Keywords. Distributed parameter systems; optimum design; sensitivity analysis; boundary control; computational method; optimisation.
I NTRODUCTI ON
is defined on a two-dimensional spatial domain. The extension of the results to higher dimension s is straightforward.
De s i gn of di s tributed parameter systems can be often characterized as free boundary problems in which part of the boundary, the free boundary, i s unknown and must be determined as a part of the problem. Pironneau (1973, 1974 ) worked on optimum (minimum drag) design problems in a st eady viscous fluid by utilizing the optimal control theory in distributed parameter systems (Lions, 1971). In Pironneau (1973, 1974), however, no comprehen s ive optimization theory wa s derived and applications of the theor y to the full (timedependent ) ~av i er - S tokes equations remain to be studied.
We consider the fluid mechanical system Ut = f(x,y,t,u, Vu, V2u)
(1)
where u(x,y,t) is the n-dimensional state vector defined on the spatial domain Q c R2 , and f(x,y,t,u, Vu,V 2u) denotes f(x,y,t,u,ux'u y , u
A numerical optimi zation theor y for the estimation of unkno~TI boundaries in general di s tributed parameter sys tems was developed by Chen and Seinfeld ( 19 75) . In this paper, we extend that ana l ysis to related problems of optimum desi gn in f luid mechanical distributed parameter sys t ems . Thu s , the purpose of this paper is to formulate in a comprehens ive manner and deve lop computational algorithms for the de rivation of the optimum shape or the optimum location of the boundar y in a cla ss of general distributed parameter systems. OPTnIU~1 DES IG\ PROBLHIS 1\ DISTRIBUTED PARAMETER SYSTEMS
In thi s paper, we confine our attention to optimum design problems in nonlinear systems of the first-order hyperbolic or parabolic partial differential equations and restrict our attention to the case in which the system 353
,u ) where the partial derivatives are xx yy indicated by a subscript, i.e., Ut' u ' and x u denote 3u/3 t, 3u/3x , and 32u/3x 2 , respecxx tively. The initial and boundary conditions to (1) are given in the general form u(x,y,O) = uo(x,y),
where u
n
x, )'
E Q
x, Y
E
3Q
( 2) (3)
denotes the derivative of u in the
direction normal to the boundary 3Q of Q. The system evolves from t=O to t=T, the period in which we are interested. In our problem, the domain Q E R2 is bounded by two surfaces rand S ( r , S E 3Q ); the outer surface r is known and the inner surface S is unknown and is to be determined optimally. Let V denote an area surrounded by S (see Fig. 1). Then the optimum design problem is to determine the coordinate location of each point of the boundary 3Q=S such that the design performance functional J defined as follows:
M. Koda
354
I:I
J
D(U,IJU,1J 2U)dXd Ydt Q
(4)
is minimized, and where D(u,lJu,1J 2u)=D(u,u x ,u y , Uxx,Uyy»O.
The problem may be viewed as an
optimal control problem of a class of distributed parameter systems, ~here the (open loop) control is the location of the boundary.
relate 6 ( oQ) to a variation in the design performance functional, 6J . A The variation 6 ( oQ ) transforms local coordinate system in the original domain Q to the global coordinate system in the new perturbed domain ~ (coordinates denoted by ~,~). The essence of the present variational approach is to find a new statement of the problem by suitable coordinate transformations. Suppose that ~ is mapped onto Q by transformations ~
(x,y,u, lJu,1J2 u ;0)
(6)
~ (x,y,u, lJu, 1J 2u; 0)
(7)
1 ~
2
where 0 is a parameter of the transformations and the new value of the state u = u(~,~,t) is mapped onto the original state function u(x,y,t) by the transformations (6) and (7) and the transformation u = ~ (x,y,u, lJu,1J 2u ; 0).
(8)
3
We assume that the transformations (6), (7), and (8) are one-to-one continuously differentiable and the value 0=0 corresponds to the identit y map:
n Fig. 1. Problem geometry.
x =
We transform the constrained optimal control problem presented here into an unconstrained problem through the introduction of adjoint functions. We define the augmented design performance functional by J
ITf ~T[f(x,y,t,u,lJu,1J2u)-u o
Q
]dxdydt (5) t
where ~(x,y,t) is the n-dimensional adjoint function and superscript T denotes transpose. Then if 50 is the boundary that minimizes J, it also clearly minimizes J . A In our approach to this optimization problem, we use a variational method on a variable region (Gelfand and Fomin, 1963). In particular, we wish to derive a design sensitivity of J to variations in the boundary location by A using the functional derivative sensitivity arguments (Koda, Dogru, and Seinfeld, 1979; Koda and Seinfeld, 1982). The constructive approach presented in this paper is based on the idea of design sensitivities to seek a design improvement by decreasing the design performance functional as much as possible. VARIATIO ~ AL
APPROACH
We consider a variation in the boundary location, 6(oQ)=6S, which transforms the boundary oQ to
Y
~
u =
~
(x,y,u,lJu, 1J 2u;0)
(9)
(x,y,u,lJu,1J 2u;0)
(10)
(x,y,u, lJu, IJLu;O) .
(11 )
1
2 3
The first variations of x, y, and u are defined as
A=J +
~
0 ~=o Q + 6 (o Q ) ~ith
a consequent mapping of
the domain Q into ~ = Q + 6Q , the domain corresponding to o~. The variation 6(o Q) introduces a variation in J . Our goal here is to A derive a design sensitivity 6J /6 (o Q), or to A
a
2 .
6x
~-x
u,O)
(12)
6y
~-y = o~o ~ 2 (x,y,u, lJu, 1J 2U;0)
(13)
6u
u(e~, t )
o~ ~ (x,y,u,lJu, 1J 00 1
u(x,y,t)
o ~o ~3 (x,y,u, lJu, 1J 2U;0)
(14)
where we have used the Taylor expansions of (6), (7), and (8), assuming that 0 is a small quantity. We may note that (14) can be expressed as 6u = 6u + u x 6x + UY6",
(15 )
where 6u=u(x,y,t)-u(x,y,t) is a perturbation of the state function in usual variational sense. Equation (15) implies that the total variation 6u is decomposed into a contribution due to the variation in the state, 6u, and into a contribution due to the coordinate transformations, u 6x+u 6y. The Jacobian of the transformationxis Yapproximated as
l o ( ~ ,~ )
l o (x,y)
I =1
+
~x OX
+
o
~y.
(16)
The first variation of the design performance functional J is defined as the principal A linear part (relative to 0) of the difference
Optimum Design
~n
Distributed Parameter Systems
J:J~{O(U'VU'V2U)
oJ A
+ tjJT (f;,ll,t) [f(f;,T),t,u,vu,V 2u)-u ) }df;dT)dt t -
J:J~{O(U'VU'V2U)
+ tjJT(x,y,t) [f(x,y,t,u,vu,V 2u)-u )}dxdyd P .(17) t
355
follow some computational aspects of this optimal design problem. The change of the design performance functional due to a shift of a boundary node tangential to a~=s is certainly smaller than the change of the design performance due to a normal shift. Thus, we shall restrict our arguments to the boundary variation (control) which is normal to the boundary surface S. We therefore specify
We specify that tjJ(x,y,t) satisfy the following adjoint equation:
ox (s)
on(s)cos8
(22)
ay (s)
on(s)sin8
(23)
= _fTtjJ+(fT tjJ) +(fT tjJ) _efT tjJ) _efT tjJ) where n(s) is the normal to a~ at s. u U x u y u xx Uyy yy xx x y For the purposes of numerical implementation, _OT+(OT ) +(OT ) _(OT) (OT) OU, ou , and ou in (20) should be approxiu u x u y u xx- u yy (18) x y x y xx yy mated by ou and ou. Expressing OU and ou n x y with the terminal condition in terms of OU and OU , the normal and tann s tjJ(x,y, T) = O. (19) gential derivatives of OU on S respectively, and noting that S is closed, (20) becomes tjJ
t
Making use of (12)-(16) to transform the domain of integration and employing the Green's theorem, (17) can be transformed to
=
+(Hu
ITJ o
TJ [-O+(Ecos8+Fsin8)u Jo a~ n xx
[O(oxcos8+oysin8)
cos 28+H u
+H
u
u
- (H x
u
) }oucos8+{H - (H ) }Ousin8 x u U y xx y yy
OU cos8+H xx
sin 28)u nn
a~
-H +{H
yy
x
U
OU sin8)dsdt yy
y
TJ
(20)
Jo as"2
where 8 is the angle between the positive xaxis and the outward normal at a point of the boundary a~, and s is the arc length along a~ (see Fig. 2). In (20), the Hamiltonian is defined as
u
)sin8cos8u yy
sn
)ondsdt
[Ecos8+Fsin8
~(H as u
-H xx
u
)sin8cos8})oudsdt yy
H(x,y,t,u,Vu,V 2u) = O(u,Vu,V 2u) +tjJ T (x,y,t)f(x,y,t,u,Vu,V 2 u)
+H
(21)
where
Hu = aH/au x ' Hu = aH/au xx ' etc. It is x xx important to note that (20) need not follow the boundary of the unknown perturbed domain ~.
u
sin 28)ou dsdt n
yy
(24 )
where
E=Hu -(H u )x and F=Hu -(H u )y' x xx y yy In deriving (24) we have used the relations similar to (15), i.e., ou=ou+unon, ou n =oun +unnon, etc.
n
~IAI:\
THEOREM - OESIGN SENSITIVITY
We now choose the appropriate boundary conditions for the adjoint system (18) and (19). From (3) we have x, y£3S"2
0,
an Fig. 2. Coordinate system at boundary. Since this paper intends to give mathematical background for a numerical solution, we now
When
~
gu
(25 )
0, the boundary condition is given
n
by Ecos8+Fsin8+
-H )sin8cos8) ~(H as u u xx
yy
M. Koda
356
cos 2 6+H
u
sin 2 6)g yy
-1
u
n
g = O. u
(26)
depending on the appropriate form of the boundary condition (3), where the adjoint function ~ is defined by (18) and (19) with the boundary condition given by (26) or (28), respectively.
For this case, (24) becomes =
I {fT[ an
-D+(Ecos6+Fsin6)u
0
- (Hu
n
-Hu xx
)sin6coS8usn]dt}on(s)ds. yy
(27) When the boundary condition (3) is not an explicit function of un' i.e., when g-1 does not exist, we have the boundary condiu
n
tion as follows: H u
cos 2 6 + H
u
xx
( 31)
-(H u -Hu )sin 6cos 6u sn ]dt xx yy
O.
(28)
yy
Remark 1. The functional derivative sensItIvities (30) and (31) can be viewed as the sensitivity of the design performance functional with respect to a variation in the normal coordinate of the boundary location at each position s per unit arc length. Remark 2. The functional derivative sensitivities (30) and (31) are actually sensitivity density functions since appropriate surface integrals (27) and (29) are needed to relate them to a total variation oJ . A Remark 3. In our optimum design problems, the functional derivatives (30) and (31) are essentially equal to variational derivatives.
For this case, (24) becomes
=f
an
{fT[
-D+(Ecos6+Fsin6)u
0
n
-H
u
)sin6cos6u yy
sn
]dt}on(s)ds.
(29) In (27) and (29), the variation in boundary coordinates, ox and oy, has been replaced by ones), the variation of the normal to the boundary surface S (where the direction angle of the normal is redifined by 6+TI, i.e. outward normal with respect to the domain V, see Fig. 1). The relations (27) and (29) naturally lead to the consideration of functional derivative sensitivities (Koda, Dogru, and Seinfeld, 1979; Koda and Seinfeld, 1982). The analysis can be summarized in the following design sensitivity theorem. Theorem. The functional derivative sensItIvIty oJA/on(s), the sensitivity of the design performance functional (4) with respect to the variation in boundary location o(an) = ones) to satisfy the state equations (1), (2), and (3), is given by
T[ -D+(Ecos6+Fsin6)u
fo
cos 2 6+H -(H
u
u
sin 2 6)u yy
n
In this section, the derivation of necessary conditions for minimum drag profiles in full (time-dependent) Navier-Stokes flow problems is presented. In an earlier work of Pironneau (1974), a particular minimum drag problem in a steady viscous fluid has been studied. Since there exist few results of a general nature for this class of optimum design problems, it is very instructive to reduce our general results to a familiar earlier case. We consider the two-dimensional Savier-Stokes equations Ut
u
(30)
or n
t
uu x - vu y - px + v (u xx +u yy )
(32)
= - uv x - vv y - py + v( v xx +v yy )
(33)
and the continuity equation
nn
To [ -D+(Ecos6+Fsin6)u
NECESSARY CONDITIONS FOR MINIMUM DRAG PROFILES
v
-Hu )sin6cos6u ]dt sn xx yy
f
It is clear that these design sensitivities inherently contain most useful information to solve the optimum design problems in distributed parameter systems. We should note that we do not consider any design constraints so far. Such design constraints, in addition to state constraints, arise in most applications and must be treated.
x
+
vy
=0
(34)
where u and v are the fluid velocity components in x- and y-coordinate directions, respectively, p denotes the pressure, and v is the viscosity of the fluid. The boundary conditions are (see Fig. 1) u
= v = 0,
x, YE: S
av
(35)
Optimum Design in Distributed Parameter Systems
v = 0,
x, ye: r
357
(36)
where Uoo denotes the fluid velocity at infinity
a(~,n)
Iy l a(x,y)
and the appropriate initial condition of the form (2) is assumed.
IdXdY -
Jydxd y
We choose the design performance functional for minimum drag profiles as J =
I ITI°
~
[2 (u 2+v2)+(u +v )2]dxdydt. x y y x
(38)
Thus, the problem is the design of a minimum drag body of given volume in two-dimensional unsteady flow fields. Applying the results obtained in the previous sections, we have the following set of linear partial differential equations for the adjoint functions wand w . 1
(oxcos9+oysin9)ds
a~
This is the total energy dissipated in an incompressible fluid in the time interval [O ,T]. Our problem here is to determine the coordinate location of the boundary of the domain Y such that (37) is minimized subject to the design constraint = IydXdY = const.
I
(37)
IsOn(S)dS where we have used (16), (22), and (23). From (31) we compute the functional derivative sensitivity density function as oJ / ones) V
-2 I
T
°[U n2+V n2+2{(W In)un +(w 2n)vn }]dt+A
where A is the Lagrange multiplier associated with (38). At the optimum, the condition oJA/on(s)= must hold; this implies
°
Corollary l,where
(w ) = u w +v w -u(w ) -v(w ) -uu -vu It xl X2 IX lY x Y -v[(w) 1
xx
+(w) 1
yy
]+q
(39)
x
In the same manner, we can obtain the necessary condition for a minimum drag body of given surface area subject to the design constraint const.
(w ) = u w +v I'i -u(w ) -v(w ) -uv -vv 2 t Y1 Y2 2 x 2 Y x Y -v[ (w) 2
(w ) +(w ) = 1
x
2 Y
xx
+(w) 2
yy
] +q
(40)
Y
°
(41)
w (x,y,T)=w (x,y,T)= 2
°
w (X,y,t)=w (X, y ,t)= 0, 1
2
x,ye:a ~ =
s,r
( 45)
Corollary 2. The necessary condition for the minimum drag profile given a unit surface area is T p(S)I [U 2+V 2+2{(W ) u +(w ) v }]dt n n 1 n n 2 n n
°
with the terminal and boundary conditions 1
~=-2A/V .
2
( 42)
const.
(43)
almost everywhere on S, where pes) is the radius of curvature of S at a point s, and the adjoint system is as in Corollary 1. Proof. The extension in n of the design constraint (45) can be written as
where we have used the boundary condition (28).
( 46)
Using the sensitivity results and the main theorem, we obtain the necessary condition for optimality for the problem (32)-(38) as follows . Corollary 1. The necessary condition for the minimum drag profile given a unit volume is that there exists a constant ~ such that T I
°[Un2+Vn2+2{(W)1 nun +(w 2) nvn }]dt = ~
almost everywhere on 5, where the adjoint system is defined by (39)-(43). Proof. The first variation of the design constraint (38) can be computed as
(44)
Thus we have o£ =
I
on(s)ds sP (s)
where the relation 1/p(s)=d9/ds is used. With the help of the main theorem, the rest of the proof proceeds as in Corollary 1. Remark 4. In steady-state case, the necessary conditions (44) and (46) reduced to those of
M. Koda
358
implemented for the necessary condition (46).
Pironneau (1974). NU~IERICAL
CONCLUSIONS
IMPLEMENTATION
The design sensitivities (44) and (46) are derived as functionals of the system state and adjoint functions. Thus, simultaneous solution of the system state and adjoint equations facilitates the computation of infinite dimensional gradients (functional derivative sensitivity density functions) of the design performance functional and a modification of the gradient method can be used as an iterative unconstrained optimization technique. An iterative numerical algorithm to compute
the optimum location of the boundary is derived directly from the necessary condition (44). The basis of the method is to choose ones) such that oJ is always negative. This is A accomplished by setting ones) as ones) = -a(s){
2 +v2 IT[U o n n
+2{(w ) u +(w ) v }]dt lnn 2nn
-~}
(47)
where a(s) is an arbitrary (small) positive weighting function. In practical numerical implementations, ~ in (47) is replaced by the mean value of
T
Io
[U 2+V 2+2{(W) u +(w ) v }]dt n n 1 n n 2 n n
( 48)
on S at each iteration. The algorithm proceeds as follows: Step 1.
Make an initial guess of the boundary location.
Step 2.
Solve the Navier-Stokes equations (32)-(36) from t=O to t=T. Evaluate J.
Step 3.
Solve the adjoint equations (39)-(43) from t=T to t=O.
Step 4.
Choose an appropriate weighting function a(s) in (47) and compute the normal shift ones) by using (47) and (48).
Step 5.
~lodify
Step 6.
Return to Step 2. Repeat until subsequent changes in J are less than a preset criterion.
the boundary location.
Detailed applications would, of course, require numerical integration of the relevent timedependent Navier-Stokes equations and reasonable numerical schemes developed by Fortin and Thomasset (1979) may be used in the present algorithm. Since the present method only involves the iteration of the normal shift in boundary location, the algorithm is considerably simpler and numerically advantageous than the other algorithms. We remark that a similar algorithm can be effectively
The optimum design problems in distributed parameter systems have been considered from an optimal control standpoint. Necessary conditions for minimum drag profiles in unsteady Kavier-Stokes flows have been derived and the results are of interest in their own right. In particular, it has been shown that the general results of the paper, i.e. design sensitivity theorem, can be reduced to the familiar purely steady cases of Pironneau (1974). The present method, based on the functional derivative sensitivity analysis, is clearly applicable to a broad class of optimum design problems that are described by free boundary problems. REFERENCES Chen, W. H., and J. H. Seinfeld (1975). Estimation of the location of the boundary of a petroleum reservoir. Society of Petroleum Engineers Journal, 15, 19-38. Fortin, M., and F. Thomasset (1979). ~lixed finite element methods for incompressible flow problems. J. of Computational Physics, 31, 114-145. Gelfand, I. M~ and F. D. Fomin (1963). Calculus of Variations. Prentice-Hall, Englewood Cliffs, New Jersey. Koda, M., A. H. Dogru, and J. H. Seinfeld (1979). Sensitivity analysis of partial differential equations with application to reaction and diffusion processes. J. of Computational Physics, 30, 259-282. Koda, M., and J. H. Seinfeld (1982). Sensitivity analysis of distributed parameter systems. IEEE Trans. on Autom. Control, in press. Lions, J. L. (1971). Optimal Control of Systems Governed by Partial Differential Equations. Springer-Verlag, New York. Pironneau, O. (1973). On optimum profiles in Stokes flow. J. of Fluid ~lechanics, ~, 117-128. Pironneau, O. (1974). On optimum design in fluid mechanics. J. of Fluid Mechanics, 64, 97-110.