- Email: [email protected]
Shape design sensitivity analysis of thermoelasticity problems
Computer Methods in Applied Mechanics and Engineering 102 (1993) 41-60 North-Holland CMA 273
Shape design sensitivity analysis of thermoelasticity problems R.J. Yang Ford Motor Company, Village Plaza, 23400 Michigan Avenue, Dearborn, MI 48124, USA Received 8 April 1991 Revised manuscript received 31 October 1991 Analytical expressions for shape design sensitivities of three-dimensional, linear thermoelasticity problems are derived, using the material derivative approach. The surface integrals in the sensitivity expressions are simplified by using the natural coordinates for numerical integration. Three different strategies for calculating the sensitivity of the structural responses are described and discussed. The sensitivities of the displacement and stress computed from all three strategies are computed and match the analytical solutions.
1. Introduction
Thermal stress is one of the important concerns in designing powertrain components. High thermal stress may initiate a fatigue crack and subsequently cause failure of the components. The sensitivity of thermal stress provides useful information. It can be used for trade-off study, what-if study, and also for design optimization. The shape design sensitivity analysis predicts the effect of a shape change on the performance of the structure. Meric [1] used the boundary element method for design ~ensitivity analysis of thermoelastic solids. Dems [2] derived boundary sensitivity expressions for thermoelasticity problems. Hou et al. [3] used a domain method for design sensitivity analysis and optimization of thermoelastic solids. They analyzed a two-dimensional thermal fin and a fixed-fixed beam. TortoreUi et al. [4] derived analytical sensitivities for nonlinear 'ynamic thermoelastic systems. All the aforementioned research used the augmented functional approach for design sensitivity analysis which was first proposed by Dems and Mroz [5]. This paper proposes a shape design sensitivity method for general, three-dimensional, linear, thermoelastic solids. The derivation follows the material derivative concept developed in [6]. The analytical formulations are obtained in a variational form. Three different sensitivity strategies are described and their advantages and disadvantages are discussed. Numerical integration considerations for sensitivity equations are discussed. Finally, stresses, displacements, temperatures, and their sensitivities for four examples are computed and compared with the analytical solutions. Correspondence to: Dr. R.J. Yang, Ford Motor Company, Village Plaza, Suite 1100, 23400 Michigan Avenue, Dearborn, MI 48124, USA. 0045-7825/93/$05.00 © 1993 Elsevier Science Publishers B.V. All rights reserved
R.J. Yang, Shape design sensitivity analysis
42
2. Variational form of thermoelastic equations Consider a thermoelastic, isotropic and homogeneous solid in three dimensions. The steady state heat conduction equation and different boundary conditions are written as [7, 8]
-kO.-g
in n ,
0 = oo
o n r 0° ,
(1)
kOini = q
on F0l ,
koch' + h(O- 0®)= 0 on F~, where 0 = T - T 0, T is the temperature, To the reference temperature of the natural (stress-free) state of the solid body, 0o the prescribed temperature, 0= the ambient temperature, n ~ the unit normal vector on the boundary, k the heat conductivity of the body, h the convective heat transfer coefficier.', q the heat flux vector, g the internal heat source, F ° the temperature boundary, F0t the heat flux boundary, and F0~ the heat convection boundary ~. The weak form of the heat transfer problem is obtained by multiplying both sides of the differential equation (1) by a virtual temperature § and then integrating over the physical domain/2, to obtain
f kO,#,dn + fr:. hO~dr = f gOdl2 + frl.q#dr + fr~.hO~#dr
for all #~ I:,
(2)
where Y is the space of kinematically admissible temperature; i.e.,
v ~- (0 ¢ [H'(a)]: 0 =0, x ¢ r ° ) .
(3)
Define energy bilinear form and the load flmctional for the temperatures as follows:
A(O,O)~,f kO~d/2+ fr:.hO0dr L(O)E f g~ dK2+ fr~ q~d r +
f:.
(4) hO®OdF,
Equation (2) is then rewritten as the simple form A(e, ~)= L(~)
for all ~ E Y,
(s)
' The superscript and subscript indicate component and differentiation with respect to field variable, respectively. The Einstein summation convention for a repeated index is used throughout this paper.
R.J. Yang, Shape design sensitivity analysis
43
For three-dimensional elasticity problems, the differential equation for equilibrium and the boundary conditions are written as
o'~ =
in ~ ,
fi
zi= z°
on r ° ,
t'= tr'in i
on
(6)
r 2 , . .
where the stress tensor ,7'j is defined as (7) and f i is the body force, t i the traction on the boundary, z ~the displacement component in the ith direction, F ° the displacement boundary, r 2 the traction boundary, D ~jm" the elasticity tensor, a the coefficient of linear thermal expansion, fl the thermal modulus and/3 = aE/ ( 1 - 2v), E is Young's modulus and v is Poisson's ratio. The weak form of the elasticity problem is obtained by multiplying both sides of the differential equation (6) by virtual displacement i ' and then integrating over the physical domain fl to obtain
f tr~i(z, O)eiJ([)dr/=
f f'i' d~/+
fr 2 ,i~i dr
for all ~ E Z .
(8)
Using the stress-strain relationship, (8) is rewritten as
f
f
fr t'i'dF forall£~Z ,
(9)
where Z is the space of kinematicaUy admissible displacements; i.e.,
zffi (z~ [H'(o)]': z' =o, x ~ r ° ) .
(lo)
Note that this is the variational form of the elasticity problem and not the variational problem
governing thermoelasticity. Define energy bilinear form and load functional for the displacement vector as r
a(z, [) -
~ ~'J(z)~'J(~)da, (11)
eft)- f f'~' + 130z.idn +
ti~. i dF.
Equation (9) is then written as the form
a(z, i)
- 1(£) for all £ • Z .
(12)
44
R.J. Yang, Shape design sensitivity analysis
3. Design sensitivity analysis The following design sensitivity analysis uses the material derivative concept from continuum mechanics. Several material derivative formulas are given for deriving the analytical expressions for sensitivities of performance measures with respect to shape change. The derivation for these formulas and more information can be found in [5, 6, 9]. The pointwise material derivative is defined as 8 z Ox k = z ' + z k V k , = z ' + Ox k Ot
(13)
where t is considered as the design in this case and the velocity V is considered as the design change. Assuming reasonable smoothness, the partial derivative and derivative with respect to x' commute as follows:
az), az' -~
-
(14)
8X i •
Let ~1 be a domain functional defined as 0t = f F d,f2.
(15)
The material derivative of Ot is
~, = f [F' + (Fvk)k] dR =
f
(F' + FkV' + FV:)dO.
(16)
Let 02 be a boundary functional,
~ -- fr F dF
.
(17)
The material derivative of 02 is
(b2-- fr [F' +
(F~n' + HF)Vkn k]dF,
(18)
where H is the curvature of F in R 2 and twice the mean curvature in R 3, Two methods that are commonly used for design sensitivity analysis of elasticity problems are employed, the direct method and the adjoint variable method. The direct method, i.e., the implicit differentiation approach, differentiates the state equation and obtains sensitivity for different structural responses. The adjoint variable method first differentiates the constraint equation and then eliminates the unknowns by the adjoint variable method. Three different sensitivity strategies are proposed and discussed in the following.
45
R.J. Yang, Shape design sensitivity analysis
3.1. Direct method Taking the material derivative of (2), using (13), (14), (16) and (18), and rearranging terms to obtain
A(O, O) -- f kO~ da
+
r~ h6# d r
= frJ q[-O~V~ + I
+
(Okn k +
HO)Vmnm] dr
fr~o h[(O®Ok- 00~ - O~O)Vk + {[(0 - O=)O]~n~ + H(O f k[(O~
+
-
O=)o}Vmnm] dF
k + g ~ V ~ + gOV~- ~ d a , O~O~)V~-O~O~V~]-
(19)
where g'= q'= 0'® = 0 is assumed and the path for the virtual temperature 0 - 0 is chosen. Comparing (2) and (19), it is noted that the bilinear forms on the left side are identical and that the right sides are different; namely, they can be considered as two different loadings with the same stiffness matrix. Field equation (2) is solved first to obtain 0 and to keep the decomposed stiffness matrix. It is reused for the sensitivity equation (19) to obtain 0. For the elasticity part, one takes the material derivative of (9), assumes the fact that f~' - t ~' - ~ = O, and rearranges terms to obtain t"
a(~, £) -- J O°(~)e.°(~.)d~2 ^ ij( O ~-tv , k - ~°(O~°(~)v~l + f ~ ' v ~ + ,f, i -zi . ,~,~ k - f [~°(Oz~V~ + ,~
+ ~Io~',-oi~v~ + or,v~l da + fr,~-t~V ~ + [(t~J)kn~+ Ht'~]Vmn m dr, f
^ ii -
~
k
-i
k
•,i - i . r k +fkz r + y~.i-ilrk z v kdJ2
+ fr~-t~z.~ Vk + [(ti~)knk + Hti~.i]Vmnm dF.
(2o)
Similar to its thermal counterpart, the elasticity field equation (9) is solved first for displacement z and subsequently the sensitivity equation (20) is solved for i, with the identical stiffness matrix. Notice that the displacement sensitivity ~ depends on displacement z, temperature 0 and temperature sensitivir 0.
R.J. Yang, Shape design sensitivity analysis
46
3.2. A d j o i n t variable m e t h o d
Define a performance measure of structural responses @,
!1#= f F( zm, o'iJ(z, 0)) d n .
(21)
Taking the material derivative of (21), the variation of ~, is obtained as f
~- J F'(Z m, O'Q(Z, 0)) "~"(FVk)k d~']
_ f OF - - 07 =
OF
[~'(z, 0)]' + o 7
f O~OFq [eij(~ .
.
zkVk. )
z ~' + (~'~)~ d a
~iJ~(O .
OkVk)] j¢ ~Z m (~m .
zm k V k) + (FVk)k dO
Oo.o + OF
-'-oct' [-8'/30
+ O'(-zkV k) + 8"/30kVk]
OFm ZkmV"+ (FVk)~ d a OZ
(22)
where the displacement and temperature sensitivities ~ and 0 are unknowns. To obtain an explicit expression for @ in terms of the velocity field V, ~ has to be eliminated. An adjoint equation is introduced by replacing ~ E Z in (22) by a virtual displacement ,~ E Z, and equating terms involving A to the energy bilinear form (11), yielding the adjoint equation for the elasticity problem,
a(A, A) =
~
?r'~(,~) + ~
,~" dO
for all ,~ E Z.
(23)
Notice that if A is replaced by ~, the right-hand side of (23) is identical with the first two terms on the right-hand side of (22). Thus, m k = .(,, +) + f ~OF [-8'~0 + +'J(-=+v ~) + 8'+/~o+v ~] - --~OF zk V + (FVk)k dO Zm (24)
Using the symmetry condition of the bilinear form (20) and substituting A for £, the variation of ~ becomes = a(~, ~) + f o~,~ oF - " - [-~'00
i
+ o'(-=~v ~) + ~'Oo,,V*l
oF Oz"
zk,. V ,, + (FVk), dO
OF ]dO
i k +,¢'(z, "" i k" o),~'(,)v**] + f [o.,,i](,)z~v~ e)~,,V~ o.~J(z, + f '',,a V ~ + f "x
V, d~
R.J. "fang, Shape design sensitivity analysis
47
OF ztv* + (FV*), da + f °-LF [6-'~(-z,y ') + ,~'J/3o,,V"l oz" Oo.iJ
+ fF ~ - t .i,,^ ki rvrk + [ ( t ~ ) k n k + H t ~ ] V ~ n '~ d r .
(25)
Notice that (25) is explicit and computable except for terms including 0. Use the same approach for £ and define an adjoint equation for the thermal part as follows:
A(A, .4) = ~:/3fi, AI - 0or,---~. ~ dD
for all ft. E Y.
(26)
Replacing 0 for A in (26) and substituting it into (25), the variation of ~ can be rewritten as ~, = A(A, 0)
+ f [6.0( A)ZkV ~ ik + tr'itz, 0) AkV k _ trOtz, O)~ij(A)V~] + f kAiV k + :~'qiVk ,, --k d/'] + f -.--OFi] [~iJ(__zkVk ) dp ~iJ~Okgk ] ~Gr
~F m k Oz m z k V + (Fvk)k d n
+ fr 2 -t'A~V k + [(t'A')kn k + HtiA']Vmn m dF.
(27)
Using the symmetry condition of the bilinear form (19) and substituting A for 0, the variation of ~ becomes
-- f k[(OkA, + O~Ak)V~ - O,A,V~] + gkAV k + gaV~ d12
+ f [~"(a)z'~v~ + ,,'J(,, O)A,Vj'* -crU(z, +
~
Ie,j(_z,y,, ) + 8,jOo~v~]
O),'~(,~)V~,I + f,a' 'V k +f'A'V~ dO
oF 7.k v~ + (FVk), d O t~zm m
+ fr~o q[--Ak Vk + (Ak nk + HA)Vmn m] dF
fr~ h[(O.Ak - OAk - O,A)V k + + fr*-t'A'~Vk
{[(0
-
O.)A]kn k + H(O - O.)A} Vmn m] d r
+ [(t'X')knk + Ht'A']vmnm dr.
(28)
Notice that (28) is a computable expression for the sensitivity of a general performance functional. The solution sequence is the temperature 0, the displacement vector z, the adjoint displacement vector A, and finally the adjoint temperature A.
48
R.J. Yang, Shape design sensitivity analysis
3.3. Implementation considerations
In the direct method, the sensitivities are obtained by solving, the heat transfer problem 0 of (2), the elasticity problem z of (9), the temperature sensitivity 0 of (19), and the displacement sensitivity ~ of (20). Notice that the sequence for design sensitivity analysis is identical with the structural analysis, namely, temperature and then displacement, and that the decomposed stiffness matrix can be reused for both the thermal and elasticity analyses. If a single load case is considered, the computational overhead for design sensitivity analysis is NV heat transfer and NV elasticity analyses, where NV is the number of the design variables. As the decomposed stiffness matrices can be reused from the analyses, only the forward and backward substitution is required. The adjoint variable method is different from the direct method in obtaining the sensitivities. The heat transfer and thermal stress analyses are first solved, the constraint functional is differentiated, and subsequently the adjoint displacement of (23) is solved before the adjoint temperature of (26). Similar to the direct method, the decomposed stiffness matrices for temperature and displacement can be reused for the sensitivity analysis. The sensitivity overhead in this case is N heat transfer and N elasticity analyses, where N is the number of the active constraints. Neither the adjoint variable nor the direct method can provide the most efficient computation if used alone. The most efficient way for computing the sensitivities is the hybrid approach [10, 11]. The hybrid approach compares the active constraints N and the design variable number NV. If the active constraint number N is smaller, as in the overdesigned case, the adjoint variable method is employed. On the other hand, if N is greater than NV, indicating the underdesigned case or the design is infeasible, the direct method is used. Another way. of obtaining the sensitivities is to use (25) for sensitivity calculation, in which the unknown 0 is solved by the temperature sensitivity equation of (19) and the A by the adjoint displacement equation of (26). The advantage of this strategy is that the programming efforts are minimal compared to the adjoint displacement-adjoint temporature approach. This stratogy can be considered as the sensitivity analysis for the elasticity part only, provided that 0 and 0 are known. The hybrid approach can also be applied at this level, i,e., after # and 0 are solved. Table I summarizes these three different strategies. The first level hybrid approach which is most efficient combines strategy l and strategy II, while the second level combines strategy l and strategy III. The overhead for strategy lII is that NV heat transfer sensitivity analyses have to be solved for each iteration. That requires more computation resources than the direct-direct or Table 1 Implementation strategies for sensitivityanalysis Strategy i StrategyII Strategy ill Heat transfer Direct Adjoint Direct sensitivity analysis (NV)* (N) (NV)
Elasticity Direct Adjoint Adjoint sensitivity analysis (NV) (N) (N) *Number in parentheses is the analysis count. N, number of active constraints; NV, number of design variables.
R.J. Yang, Shape design sensitivity analysis
49
adjoint-adjoint approaches. However, because the degrees-of-freedom for the temperature problem are fewer than for the stress problem, it may still be viable. 4. Numerical integration
The variational forms of the structural and sensitivity analyses derived in the previous section are not limited to any specific number method, although the finite element method is commonly used. A numerical integration scheme must be used to form the matrix equation, when the finite element method is employed. In general, Gauss quadrature is employed in a numerical integration scheme for the finite element analysis. It first transforms the integral to the natural coordinates and integrates over the coordinate system. For example, the functional in (15) can be transformed to the form
~, = f FJ d ~ ,
(29)
where l is the determinant of the Jacobian matrix which transforms the undeformed configuration into the natural coordinate system, and d ~ the volume increment in the natural coordinates. Haber [_12] found that (29) could be exploited for sensitivity analysis, as the volume increment d n is constant for all time or for all designs in this case. The variation of (29) is written as
(30)
4,, = f (F' + F,v*)J + FJ d / ) ,
where J is the derivative with respect to the design variable. Comparing (30) and (16), it is found that J d~] -- V~ dD. That relationship was pointed out and proved in [13]. Arora et al. [14] used the same 'reference volume' approach for design sensitivity analysis with nonlinear structural response. Using the same approach, the variation of a general functional ~2 on the boundary is
4,2= fr (F' + f,V*)J + FJ d r ,
(31)
where dr' is the surface increment in the natural coordinates. Comparing (31) and (18), it is found that (31) is simpler than (18), as (18) contains the curvature term H. Using the relationships in (30) and (31) and applying the same approach described in the previous section, the temperature sensitivity equation of (19) is rewritten as
A(O, O)- f kO~J dD + fr 2 hOOJdI" -
(q,v j + qj) dP + f.o h(O o) J dP -
(32)
50
R.J. Yang, Shape design sensitivity analysis
The displacement sensitivity equation of (20) is rewritten as
a(~,[)- f ~J(~)~J(~)J d.0 = f ~[or ' J (Z)zkV - ' i~ + crii(z, O)zkV -~ j~ + /3oz:i+ fikvki, ilJ + [f'E'-- ¢riJ(Z, O)eiS(f.)]J dfJ +
z (tkV J
-k t~j)£. ' d.P.
(33)
The variation of a general structurai response ~, is in the following form:
= f k[(OkA ~+ O~Ak)VkJ- O~A,J] + (gkVkJ + gJ)a d~2 +f[~j( +
, k + ~r~J(Z, A)zkVj ~ jk + f ~kV kA~ ]J + [fiA' O)AkV
"-'---oF 8o.iJ
[Oo(_z,V,) + 8 , q J o ~ v ~1
oz m Zk V
+ frio (q, VkJ+ q J ) A d [ ' + fr~o h(O®-O)AJd['+
-
o'q(Z, o)EiJ(A)]J d~
"b ( F V k ) k
<,:v'J
J
d~ dP
<34)
It is observed that the surface integral can be simplified significantly, however, the domain integrals are almost unchanged. 5. Numerical examples
Four examples are used to demonstrate and verify the analytical sensitivity expressions by using all three strategies described in Section 3, direct-direct, adjoint-adjoint and directadjoint. The first example verifies the thermal stress sensitivities subjected to a temperature change. The second and the third examples verify the thermal displacement sensitivities considering the heat conduction and convection effects. The final example is a threedimensional tube which is used to show the applicability of the technique to realistic problems.
5.1. One.dimot~ional thermal stress problem A one-dimensional tapered beam shown in Fig, 1 is the first example, The beam is heated from 0 °C to T = 1000 °C. The analytical stress field is obtained by the compatibility condition between the thermal expansion and the deformation due to the mechanical load. The expression is shown as follows: e:a r(d~ - dl)
or = - [d, + (d 2 - d,)xtl] Iog(dJd~) "
(35)
Selecting the heights of two end points as the design variables d I and d 2, the stress
R.J. Yang, Shape design sensitivity analysis
51
w Fig. 1. Clamped beam.
sensitivities are do" = - E a T
(d2/d~ - 1)[d, + (d 2 - d , ) x / l ] - d 2 l o g ( d 2 / d t ) [d, + (d 2 - d , ) x / I ] l o g ( d 2 / d , ) 2
do" = - E a T
(d'/d2-
dd 2
1)[d, + (d 2 - d , ) x / l ] + d, l o g ( d 2 / d , )
(36)
[d~ + (d 2 - d~)x/l] log(d2/d~) 2
The one-dimensional problem is solved by the finite element method and is modeled as a three-dimensional beam problem. Insulation boundary conditions are imposed on all surfaces except the two ends. The finite element model of the beam contains 68 grid points and 5 20-noded HEXA elements. The left and right ends of the beam are clamped. The length I and the width w of the beam are chosen as 50 mm and I mm, respectively. Young's modulus E, Poisson's ratio v, and the coefficient of linear thermal expansion o. for the beam are 2.07 × 10s MPa, 0.3 and 1.2 × 10 -7 °C -1, respectively. Heat transfer analysis is not required for this problem, as the temperature distribution is known. However, to check the heat transfer equation in (2), temperatures of 1000°C are specified at both ends. As expected, the temperature results in the beam are identical with the boundary temperatures. This indicates that the heat transfer finite element analysis is capable of describing a uniform temperature profile. The stress and stress sensitivities, from the exact solutions and the finite element analysis, are listed in Table 2. No stress error greater than 0.05% from finite element analysis is observed in this case. All three senslt "iv"i~y strategies: direct-direct, adjoint-adjoint and direct-adjoint show identical accuracy in which the error is below 4%, compared with the analytical solutions. 5.2. Heat c o n d u c t i o n
The second example is a one-dimensional heat conduction problem. The beam has an identical configuration to the first example, as shown in Fig. 1. In this example, the design variable is the length of the beam, while the heights at either end remain unchanged. The one-dimensional heat conduction problem is governed by
-kOxx = g ,
x~[O,/l, (37)
0(0)=0,
0(1)=T .
t~ tO
Table 2 Stress sensitivity c o m ~
for damped beam
(1) x (mm)
(2) Exact stress
(3) FEM stress
(4)" Error (%)
(5) b DV val~e
(6) DSAf exact
(7) DSAs I
(8) c Error (%)
(9) DSAh 11
Error (%)
(11) DSAi III
(12)" Error (%)
0
-27.25
-27.24
0.0
1(1.0) 2(1.2)
14.0~ -11.7
13.95 -11.63
0.6 0.6
13.95 -11.63
0.6 0.6
13.95 -11.63
0.6 0.6
10
-26.2
-26.19
0.0
1(1.0) 2(1.2)
20
-25.23
-25.22
0.0
1(1.0) 2(1.2)
(10) d
7.452 -6.21
7.37 -6.142
1.1 1.1
7.37 -6.142
1.1 1.1
7.37 -6.142
1.1 1.1
1.785 -1.487
!.716 -1.43
3.9 3.9
1.716 -1.43
3.9 3.9
1.716 -1.43
3.9 3.9 ¢b
30
-24.33
-24.32
0.0
2(1.2)
-3.106 2.588
-3.165 2.637
-1.9 -1.9
-3.165 2.637
-1.9 -1.9
-3.165 2.637
-1.9 -1.9
l(1.0)
r~
40
-23.49
-23.49
0.0
1(1.0) 2(I.2)
-7.338 6.115
-7.388 6.157
-0.7 -0.7
-7.388 6.157
-0.7 -0.7
-7.388 6.157
-0.7 -0.7
50
-22.71
-22.7
0.0
1(1.0) 2(1.2)
-11.01 9.174
-11.05 9.210
-0.4 -0.4
-11.05 9.210
-0.4 -0.4
-11.05 9.210
-0.4 -0.4
" ( 4 ) = (1 - 1(3)/(2)1) × i 0 0 ~ .
b Design variables (d~ and d2) and values. ( 8 ) - (1 - 1(7)1(6)1) x 1 0 0 % . ( 1 0 ) = (1 - l ( 9 ) / ( 6 ) l ) x 1 0 0 % . ( 1 2 ) = (1 - 1 ( 1 1 ) / ( 6 ) 1 ) x 100%.
f Exact solution from (36). s Strategy [: dire.ct-dire.ct. a Strategy If: adjoint-adjoint. ~Strategy IH: direct-adjoint.
zs z~
R.J. Yang, Shape design sensitivity analysis
53
The exact solution for the temperature is obtained by integrating (37) as
t +
7
•
-
(38)
The differential equation for the elasticity field and boundary conditions for the onedimensional case are On
E Ox
/30
Constant
x~[0,1] (39)
u(O)-O,
u(O-o.
The exact solution for displacements and stresses are integrated to obtain u
- ot [ gx3
or = - E a
T
gl 2
(40)
(r~" + 1"i2"k/"
The temperature, displacement and stress sensitivities with respect to the length I are dO dl du dl d~ dl
O0 OI -
au OI Oor O!
O0 Ox #x Ol b
g x ( I - x) kl T
Ou Ox Ox OI Oor Ox Ox O!
3g
gl
(41)
Eagl 6k
The conduction problem is also modeled as a three-dimensional beam problem, as shown in Fig. 1. Insulation boundary conditions are imposed on all edges, 0 °C and 1000 °C temperatures are applied at the left and the right ends, respectively. Both the left and the right ends of the beam are clamped. Young's modulus E, Poisson's ratio v, and the coefficient of linear thermal expansion a for the beam are 7.07 x 105 MPa, 0.3 and 1.2 x 10 -~ °C -!, respectively. The ambient temperature T® and fin base temperature are 0°C and 1000°C. The thermal conductivity k is assumed to be unity and the internal heat source g = 0. The temperature and its sensitivity, computed from the finite element analysis and from the analytical expressions of (38) and (41), are listed in Table 3. The finite elemen* results show excellent agreement with the analytical solutions. The displacements and the sensitivities, computed from finite element analysis and from the analytical expressions of (40) and (41) are listed in Table 4. The finite element results again show excellent agreement with the analytical solutions. The maximum error in displacement sensitivity is below 1%. The stress is 12.42 MPa both from the finite element and from the analytical solution. The stress sensitivity diminishes, as no internal heat source is applied.
R.J. Yang, Shape design sensitivity analysis
54
Table 3 Temperature and sensitivity for heat conduction ,,
T
dT/dd
X
(ram)
Exact
FEM
Exact
FEM
0 10 20 30 40 50
0 200 400 600 800 1000
0 200 400 600 800 1000
0.0 0.0 0.0 0.0 0.0 0.0
0.3333E-11 - 0.1584E- 10 0.7804E- 12 0.4822E- 11 0.1986E- 10 -0.3333E- 11
5.3. Heat convection
A straight thermal fin of rectangular profile is chosen as the third example. Figure 2 illustrates the configuration and the notations. The symbol I denotes the length of the fin, k the thermal conductivity, h the heat transfer coefficient at the exposed surface, To the fixed temperature at the fin base, and T~, the temperature of the bulk of the ambient fluid. Without losing any generality, the heat losses from the side edges of the fin are neglected and the heat flow is expressed per unit width. If it is assumed that the fin is sufficiently slender, the fin can be treated as a one-dimensional problem with no heat convected out of the free end. With the assumptions, the analytical temperature of the fin can is written as [15] O _- T - T~. _ cosh m ( l - x) 0o T o - T. cosh ml '
(42)
where m - - V ~ A ", C is the perimeter and A is the cross sectional area. For a wide fin, assuming no heat losses from the side edges is reasonable, and m -V'2"-h/-~, where d is the
height of the fin. For a freely expanded fin, the differential equation for the elasticity part is as follows: E OU
ox - t3o = o ,
x~[O,l],
(43)
x(O) = o .
,
a
i|r
T./~,,
k mm
i
....
x
III W
h : heat transfer coefficient k : thermal conductivity Fig. 2, Thermal fin of rectangular profile.
R.J. Yang, Shape design sensitivity analysis
55
The displacement is obtained by integrating (43) as u(x)
a[sinh m l
0o
-
sinh m ( l
-
x)]
(44)
m cosh m l
Selecting the height of the fin d as the design variable, the sensitivities for the temperature and the displacement are d(O/Oo)
O(O/Oo) a m am Od
dd
0m (l - x) sinh m ( l - x) cosh m l - I sinh m l cosh m(1 - x) = ad cosh 2 m l ' d(u/eo)
O(u/O0) Om
dd
Om
Od
= a { l m [ 1 - sinh m ( l - x) sinh ml]
- cosh m l [ m ( l - x) cosh m ( l - x) + sinh m l - sinh m ( l - x)]} I ( m cosh m l ) 2 0 m ad '
Om
1
2k~~
(45)
Od = - 2
The thermal fin problem is also modeled as a three-dimensional beam problem, as shown in Fig. 1. Insulation boundary conditions are imposed on side edges, and only the top and bottom surface are allowed heat convection. The left end of the fin is clamped. Young's modulus E, Poisson's ratio v, and the coefficient of linear thermal expansion a for the fin at 2.07 x l0 s MPa, 0.3 and 1.2 x 10 -7 oC=~, respectively. The ambient temperature T® and fin base temperature is 0 oC and 1000 oC. The thermal conductivity k and heat transfer coefficient h are assumed to be unity. The temperature and its sensitivity, computed from the finite element analysis and from the analytical expressions of (42) and (45), are listed in Table 5. The finite element results show excellent agreement with the analytical solutions. The assumption of no heat convected at the right end is proved to be reasonable, as the temperature is diminishing at the fin end. Table 5 Temperature and sensitivity for heat convection
T
dT/dd
X
(mm)
Exact
FEM
Exact
FE M
0 10 20 30 40 50
0.1000E + 4 0.2431E+3 0.5912E+2 0.1442E+2 0.3700E + 1 0.1699E + 1
0.1000E + 4 0.2448E+3 0.6032E+2 0.1481E+2 0.3856E + 1 0.1779E + 1
0.0000E +0 0.1719E+3 0.8366E+2 0.3073E+2 0.1076E + 2 0.6006E + 1
0.5106E-09 0.1671E+3 0.8257E+2 0.3058E+2 0.1085E+2 0.6097E + 1
0.0 -0.48E-3 -0.72E-3 -0.72E-3 -0.48E-3 0.0
0 10 20 30 40 50
0.0 -0.481b-E-3 -0.']7.03E-3 -0.7203E-3 -0.4805E-3 0.0
(3) FEM disp. 0.0 -0.96E-5 -0.144E-4 -0.144E-4 -0.96E-5 0.0
(5) DSA ~ exact
a ( 4 ) = ( 1 - 1(3)/(2)1) x b ( 7 ) = (1 - ] ( 6 ) / ( 5 ) ] ) x " ( 9 ) = (1 - 1(8)/(5)1) x d (11) =(1 -J(10)/(5)D 100%. 100%. x 100%.
100%.
0.OO 0.6441E-3 0.8025E-3 0.84~9E-3 0.8540E-3 0.8645E-3
0.0 0.6422E-3 0.7984E-3 0.8364E-3 0.8457E-3 0.8485E-3
10 20 30 40 50
0
FEM disp.
Exact disp.
x (mm)
(3)
(2)
(t)
0.0 -0.9591E-5 -0.1439E-4 -0.1439E-4 -0.9591E-5 0.0
0.0 O.1753E-3 0.3283E-3 0.3925E-3 0.4152E-3 0.4243E-3
DSA c exact
(5)
0.0 0.1816E-3 0.3341E-3 0.3998E-3 0.4229E- 3 0.4374E-3
DSA t 1
(6)
" E x a c t solution from ( 4 5 ) . f Strategy I: direct-direct. * Strateb~ II: adjoint-adjoint. b Strategy HI: direct-adjoint.
-0.3 -0.5 -0.9 - 1.0 -1.9
m
Error (%)
(4)"
-0.0 -3.6 -1.7 -1.9 -1.8 -3.1
(%)
(7) b Error
0.1 0.0 0.0 0.1
(%)
(7) b
Error
(6) DSA f I
Exact solution from (41). f Strategy i: direct-direct. * Strategy II: adjoint-adjoint. b Strategy HI: direct-adjoint.
m -0.1 0.0 0.0 -0.1 ~
(4) ~ Error (%)
Table 6 Displacement sensitivity comparison for heat convection
" ( 9 ) = (I - I(8)/(5)I)x i o o % . d (II) = (I - ] ( 1 0 ) I ( 5 ) D x 100%.
*'(7) = (I - I(6)/(5)l) x I o o % .
~'( 4 ) = (I - [(3)/(2)I) x I O O % .
(2) Exact disp.
(1) x (mm)
Table 4 Displacement scnsitivity comparison for heat c o m l ~
(8)
-0.0 0.1816E-3 0.3341E-3 0.3998E-3 0.4229E-3 0.4374E-3
DSA g II
(8)
0.0 -0.9591E-5 -0.1439E-4 -0.1439E-4 -0.9591E-5 0.0
DSA g II
-3.1
-1.8
-3.6 -1.7 -1.9
(%)
(9)" Error
0.1 0.0 0.0 0.1
(%)
(9) ~ Error
0.0 0.1816E-3 0.3341E-3 0.3998E-3 0.4229E-3 0.4374E-3
(10) DSA h III
0.0 -0.9591E-5 -0.1439E-4 -0.1439E-4 -0.9591E-5 0.0
(10) DSA" III
-3.6 -1.7 -1.9 -1.8 -3.1
(%)
(11) e Error
0.1 0.0 0.0 0.1
t%)
(11) e Error
e~
g~
¢,q
,o
¢Jt
R.J. Yang, Shape design sensitivity analysis
57
The displacements and the sensitivities, computed from finite element analysis and from the analytical expressions of (44) and (45), are listed in Table 6. The finite element results also show excellent agreement with the analytical solutions. The maximum error in displacement is below 4%. As the fin is freely expanded, no stress occurs in this case. The finite element results show very small values in the stresses and their sensitivities, which are omitted.
5.4. Three-dimensional tube problem A three-dimensional tube problem is shown in Fig. 3. The tube is idealized to simulate an engine exhaust manifold in which high thermal stress may initiate a fatigue crack and subsequently cause failure of the components. A finite element model is generated for design sensitivity analysis, as shown in Fig. 4. This model consists of 720 grids, 2110 degrees-of-freedom and 448 linear H E X A elements. Zero displacement conditions are imposed in the .x-direction of the left end (section A). One point on the surface is fixed in all three directions to eliminate rigid body motion. The rest of the points are allowed to move in the yz plane to simulate the bolt condition between the runner and the engine head. A fixed temperature condition of 780 °C is imposed at the left end. The thermal conductivity k and heat transfer coefficient h for the inner and outer surfaces of the tube are assumed to be unity. Young's modulus E, Poisson's ratio v and the coefficient of linear thermal expansion a for the tube are 2.07 x 105 MPa, 0.3 and 1.2 x 10 -7 °C -I, A
(01t ..... (ill
V
1 I
.
B rI r
design variables b(I) = r b(2) = r
0
b(3) = r I b(4) = d I b($):
db
Fig. 3. Three-dimensional tube.
o
R.J. Yang, Shape design sensitivity analysis
58
720 grids 448 HEXAelements 2110 D. O. F.
870
Fig. 4. Finite element model of three-dimensional tube.
respectively. The gas temperature in the tube is 2000 °C. The ambient temperature outside the tube is 0 °C. The x-displacement at grid 1202 (Fig. 4), Von Mises stress at grid 870 in element 340 and at grid 1202 in element 277 (Fig. 4), and their sensitivities are listed in Table 7. Note that the central finite difference results with 0.1% design perturbation is used to calculate the Table 7 Sensitivity comparison for a three-dimensional tube
( 1) ~oltem
(2) DV'
(3) FD
(8) DSA~ (Ill)
(9) d Error (%)
. . 2.5613E-4 3.7556E~4 9.7297Em5
. . 10 0.81 0,29
2.5614E-4 3.7556E~4 9.7297E-5
10 0.81 0,29
2.0343E~-2 -3,8482 6,2727 1,3381 -1,4674E-2
2.4 0.16 0,23 0,88 1,52
1(3.0) 1.9868E-2 2(1.5) m 3 . 8 5 4 4 3(1.0) 6.2585 4(10) 1.353 5(~12) -1.4901E-2
2.0344E-2 -3.8482 6.2726 1,3381 -1,4674E-2
2.4 2,0343E-2 0,16 -3,8482 0.23 6,2726 0,88 1,3381 1.52 -1,4672E-2
2,4 0.16 0,23 0,88 1.54
1(3.0) 2(1.5) 3(1.0) 4(10) 5(12)
O,1732 7.9142 6,2726 -0,65371
8,97 O,1732 1,62 7,9141 0,23 6, 2726 1.55 -0,65371
8.97 O.1732 1.62 7,9142 0,23 6, 2727 1,55 -0,65370
i
Design variable and value, b (5)~ (t -1(4)/(~)1) × 100%.
(7) = (1 -1(6)/(3)1) × 100%. d (9)--(i - 1(8)/(3)1) × 100%. Strategy I: direct-direct, f Strategy lh adjoint,adjoint, g Strategy IIl: direct-adjoint,
. .
(7) ~ Error (%)
Von Mises stress
O,15895 7.78S3 6, 2585 -0.65373
. .
(6) DSAf 11
1(3.0) 2(1,5) 3(1,0) 4(10) 5(12)
Von Mises stress at point 1202 (4.2809 MPa)
. . 2,5613E-4 3,7556E-4 9,7297E-5
(5) b Error (%)
x-di~t ,acement at point 1202 (1.5806E ~ 3 mm)
at point 870 (0.65788 MPa)
. . . . 2,3283E-4 3,7253Em4 9,7013E-5
(4) DSA~ I
10 0,81 0,29
8.97 1.62 0,23 1,55
R.J. Yang, Shape design sensitivity analysis
59
percentage errors in sensitivities, as no exact solution is available in this case. The results show that the sensitivities calculated from all three strategies have excellent agreement with the finite difference solutions. The sensitivities with relatively small magnitudes compared to others are neglected.
6. Summary Analytical expressions for shape design sensitivities of three-dimensional, linear thermoelasticity problems were derived, using the material derivative approach. Numerical integration for obtaining the sensitivities using the finite element method was discussed and the advantages are identified. It was found that the surface integral can be simplified significantly by using the natural coordinates. Three different strategies, direct-direct, adjoint-adjoint and direct-adjoint, were discussed and verified by four example problems. It was found that the sensitivities computed from all three strategies are identical and have quite good agreement with the analytical solutions. The error in displacement and stress is below 4%. It is also shown that the temperature and its sensitivity have excellent agreement with the analytical solutions. The technique is further verified by a three-dimensional tube example. The most efficient way for computing the sensitivities is the hybrid approach. Two levels of the hybrid approach were suggested. The first level hybrid approach, which is the most efficient, combines the direct-direct and adjoint-adjoint method, while the second level combines the direct-direct and direct-adjoint methods.
Acknowledgment The author wishes to express appreciation to the reviewers, Professor J.W. Hou of Old Dominion University, S. Poe of PTPE, D. Tao of PMES, D. Wagner, K. Yeung, A. Mahajan and J. Liu of CAE department for valuable suggestions and discussion.
References [1] R.A. Meric, Boundary element in shape design sensitivity analysis of thermoelastic solids, in: C.A. Mota Soares, ed., Computer Aided Optimal Design: Structural and Mechanical System, NATO ASI Series (Springer, Berlin, 1987). [2] K. Dems, Sensitivity analysis in thermoelasticity problems, in: C.A. Mota Soares, ed., Computer Aided Optimal Design: Structural and Mechanical System, NATO ASl Series (Springer, Berlin, 1987). [3] G. Hou, J.S. Sheen and C.H. Chuang, Shape sensitivity analysis and design optimization of linear, thermoelastic solids, in: Proc. 31st AIAA SDM Conf., Long Beach, CA, 1990. [4] D.A. Tortorelli, R.B. Haber and S.C-Y. Lu, Adjoint sensitivity analysis for nonlinear dynamic thermoelastic systems, AIAA J. 29 (1991) 253-263. [5] K. Dems and Z. Mroz, Variational approach by means of adjoint systems to structural optimization and .~'nsitivity analysis II, Structure shape variation, Internat. 3. Solids and Structures 20 (1984) 527-552. [6] E. 2. Haug, K.K. Choi and V. Komkov, Design Sensitivity Analysis of Structural Systems (Academic Press, New York, 1986).
60
R.J. Yang, Shape design sensitivity analysis
[7l N. Kikuchi, Finite Element Methods in Mechanics (Cambridge Univ. Press, New York, 1986). [81 T.J.R. Hughes, The FirSte Element Method (Prentice Hall, Englewood Cliffs, NJ, 1987). [91 J.A. Bennett and M.E. Botkin, eds., The Optimum Shape: Automated Structural Design (Plenum, New York, 1986).
[101 RJ. Yang, A three-dimensional shape optimization system-SHOP3D, Comput. & Structures 31 (1989) 881-890.
[lll R.J. Yang, A hybrid approach for shape optimization, in: V.A. Tipnis and E.M. Patton, eds., Computers in Engineering (ASME, New York, 1988) 107-112.
[12] R.B. Haber, A new variational approach to structural shape design sensitivity analysis, in: C.A. Mota Soares, ed., Computer Aided Optimal Design: Structural and Mechanical System, NATO ASI Series (Springer, Berlin, 1987). [13] R.J. Yang hnd M.E. Botkin, Accuracy of the domain material derivative approach to shape design sensitivities, AIAA J. 25 (1987) 1606-1610. [14] J.B. Cardoso and J.S. Arora, Variational method for design sensitivity analysis in nonlinear structural mechanics, AIAA J. 26 (1988) 595-603. [lSl A.J. Chapman, Fundamentals of Heat Transfer (Macmillan, New York, 1987).
Shape design sensitivity analysis of thermoelasticity problems R.J. Yang Ford Motor Company, Village Plaza, 23400 Michigan Avenue, Dearborn, MI 48124, USA Received 8 April 1991 Revised manuscript received 31 October 1991 Analytical expressions for shape design sensitivities of three-dimensional, linear thermoelasticity problems are derived, using the material derivative approach. The surface integrals in the sensitivity expressions are simplified by using the natural coordinates for numerical integration. Three different strategies for calculating the sensitivity of the structural responses are described and discussed. The sensitivities of the displacement and stress computed from all three strategies are computed and match the analytical solutions.
1. Introduction
Thermal stress is one of the important concerns in designing powertrain components. High thermal stress may initiate a fatigue crack and subsequently cause failure of the components. The sensitivity of thermal stress provides useful information. It can be used for trade-off study, what-if study, and also for design optimization. The shape design sensitivity analysis predicts the effect of a shape change on the performance of the structure. Meric [1] used the boundary element method for design ~ensitivity analysis of thermoelastic solids. Dems [2] derived boundary sensitivity expressions for thermoelasticity problems. Hou et al. [3] used a domain method for design sensitivity analysis and optimization of thermoelastic solids. They analyzed a two-dimensional thermal fin and a fixed-fixed beam. TortoreUi et al. [4] derived analytical sensitivities for nonlinear 'ynamic thermoelastic systems. All the aforementioned research used the augmented functional approach for design sensitivity analysis which was first proposed by Dems and Mroz [5]. This paper proposes a shape design sensitivity method for general, three-dimensional, linear, thermoelastic solids. The derivation follows the material derivative concept developed in [6]. The analytical formulations are obtained in a variational form. Three different sensitivity strategies are described and their advantages and disadvantages are discussed. Numerical integration considerations for sensitivity equations are discussed. Finally, stresses, displacements, temperatures, and their sensitivities for four examples are computed and compared with the analytical solutions. Correspondence to: Dr. R.J. Yang, Ford Motor Company, Village Plaza, Suite 1100, 23400 Michigan Avenue, Dearborn, MI 48124, USA. 0045-7825/93/$05.00 © 1993 Elsevier Science Publishers B.V. All rights reserved
R.J. Yang, Shape design sensitivity analysis
42
2. Variational form of thermoelastic equations Consider a thermoelastic, isotropic and homogeneous solid in three dimensions. The steady state heat conduction equation and different boundary conditions are written as [7, 8]
-kO.-g
in n ,
0 = oo
o n r 0° ,
(1)
kOini = q
on F0l ,
koch' + h(O- 0®)= 0 on F~, where 0 = T - T 0, T is the temperature, To the reference temperature of the natural (stress-free) state of the solid body, 0o the prescribed temperature, 0= the ambient temperature, n ~ the unit normal vector on the boundary, k the heat conductivity of the body, h the convective heat transfer coefficier.', q the heat flux vector, g the internal heat source, F ° the temperature boundary, F0t the heat flux boundary, and F0~ the heat convection boundary ~. The weak form of the heat transfer problem is obtained by multiplying both sides of the differential equation (1) by a virtual temperature § and then integrating over the physical domain/2, to obtain
f kO,#,dn + fr:. hO~dr = f gOdl2 + frl.q#dr + fr~.hO~#dr
for all #~ I:,
(2)
where Y is the space of kinematically admissible temperature; i.e.,
v ~- (0 ¢ [H'(a)]: 0 =0, x ¢ r ° ) .
(3)
Define energy bilinear form and the load flmctional for the temperatures as follows:
A(O,O)~,f kO~d/2+ fr:.hO0dr L(O)E f g~ dK2+ fr~ q~d r +
f:.
(4) hO®OdF,
Equation (2) is then rewritten as the simple form A(e, ~)= L(~)
for all ~ E Y,
(s)
' The superscript and subscript indicate component and differentiation with respect to field variable, respectively. The Einstein summation convention for a repeated index is used throughout this paper.
R.J. Yang, Shape design sensitivity analysis
43
For three-dimensional elasticity problems, the differential equation for equilibrium and the boundary conditions are written as
o'~ =
in ~ ,
fi
zi= z°
on r ° ,
t'= tr'in i
on
(6)
r 2 , . .
where the stress tensor ,7'j is defined as (7) and f i is the body force, t i the traction on the boundary, z ~the displacement component in the ith direction, F ° the displacement boundary, r 2 the traction boundary, D ~jm" the elasticity tensor, a the coefficient of linear thermal expansion, fl the thermal modulus and/3 = aE/ ( 1 - 2v), E is Young's modulus and v is Poisson's ratio. The weak form of the elasticity problem is obtained by multiplying both sides of the differential equation (6) by virtual displacement i ' and then integrating over the physical domain fl to obtain
f tr~i(z, O)eiJ([)dr/=
f f'i' d~/+
fr 2 ,i~i dr
for all ~ E Z .
(8)
Using the stress-strain relationship, (8) is rewritten as
f
f
fr t'i'dF forall£~Z ,
(9)
where Z is the space of kinematicaUy admissible displacements; i.e.,
zffi (z~ [H'(o)]': z' =o, x ~ r ° ) .
(lo)
Note that this is the variational form of the elasticity problem and not the variational problem
governing thermoelasticity. Define energy bilinear form and load functional for the displacement vector as r
a(z, [) -
~ ~'J(z)~'J(~)da, (11)
eft)- f f'~' + 130z.idn +
ti~. i dF.
Equation (9) is then written as the form
a(z, i)
- 1(£) for all £ • Z .
(12)
44
R.J. Yang, Shape design sensitivity analysis
3. Design sensitivity analysis The following design sensitivity analysis uses the material derivative concept from continuum mechanics. Several material derivative formulas are given for deriving the analytical expressions for sensitivities of performance measures with respect to shape change. The derivation for these formulas and more information can be found in [5, 6, 9]. The pointwise material derivative is defined as 8 z Ox k = z ' + z k V k , = z ' + Ox k Ot
(13)
where t is considered as the design in this case and the velocity V is considered as the design change. Assuming reasonable smoothness, the partial derivative and derivative with respect to x' commute as follows:
az), az' -~
-
(14)
8X i •
Let ~1 be a domain functional defined as 0t = f F d,f2.
(15)
The material derivative of Ot is
~, = f [F' + (Fvk)k] dR =
f
(F' + FkV' + FV:)dO.
(16)
Let 02 be a boundary functional,
~ -- fr F dF
.
(17)
The material derivative of 02 is
(b2-- fr [F' +
(F~n' + HF)Vkn k]dF,
(18)
where H is the curvature of F in R 2 and twice the mean curvature in R 3, Two methods that are commonly used for design sensitivity analysis of elasticity problems are employed, the direct method and the adjoint variable method. The direct method, i.e., the implicit differentiation approach, differentiates the state equation and obtains sensitivity for different structural responses. The adjoint variable method first differentiates the constraint equation and then eliminates the unknowns by the adjoint variable method. Three different sensitivity strategies are proposed and discussed in the following.
45
R.J. Yang, Shape design sensitivity analysis
3.1. Direct method Taking the material derivative of (2), using (13), (14), (16) and (18), and rearranging terms to obtain
A(O, O) -- f kO~ da
+
r~ h6# d r
= frJ q[-O~V~ + I
+
(Okn k +
HO)Vmnm] dr
fr~o h[(O®Ok- 00~ - O~O)Vk + {[(0 - O=)O]~n~ + H(O f k[(O~
+
-
O=)o}Vmnm] dF
k + g ~ V ~ + gOV~- ~ d a , O~O~)V~-O~O~V~]-
(19)
where g'= q'= 0'® = 0 is assumed and the path for the virtual temperature 0 - 0 is chosen. Comparing (2) and (19), it is noted that the bilinear forms on the left side are identical and that the right sides are different; namely, they can be considered as two different loadings with the same stiffness matrix. Field equation (2) is solved first to obtain 0 and to keep the decomposed stiffness matrix. It is reused for the sensitivity equation (19) to obtain 0. For the elasticity part, one takes the material derivative of (9), assumes the fact that f~' - t ~' - ~ = O, and rearranges terms to obtain t"
a(~, £) -- J O°(~)e.°(~.)d~2 ^ ij( O ~-tv , k - ~°(O~°(~)v~l + f ~ ' v ~ + ,f, i -zi . ,~,~ k - f [~°(Oz~V~ + ,~
+ ~Io~',-oi~v~ + or,v~l da + fr,~-t~V ~ + [(t~J)kn~+ Ht'~]Vmn m dr, f
^ ii -
~
k
-i
k
•,i - i . r k +fkz r + y~.i-ilrk z v kdJ2
+ fr~-t~z.~ Vk + [(ti~)knk + Hti~.i]Vmnm dF.
(2o)
Similar to its thermal counterpart, the elasticity field equation (9) is solved first for displacement z and subsequently the sensitivity equation (20) is solved for i, with the identical stiffness matrix. Notice that the displacement sensitivity ~ depends on displacement z, temperature 0 and temperature sensitivir 0.
R.J. Yang, Shape design sensitivity analysis
46
3.2. A d j o i n t variable m e t h o d
Define a performance measure of structural responses @,
!1#= f F( zm, o'iJ(z, 0)) d n .
(21)
Taking the material derivative of (21), the variation of ~, is obtained as f
~- J F'(Z m, O'Q(Z, 0)) "~"(FVk)k d~']
_ f OF - - 07 =
OF
[~'(z, 0)]' + o 7
f O~OFq [eij(~ .
.
zkVk. )
z ~' + (~'~)~ d a
~iJ~(O .
OkVk)] j¢ ~Z m (~m .
zm k V k) + (FVk)k dO
Oo.o + OF
-'-oct' [-8'/30
+ O'(-zkV k) + 8"/30kVk]
OFm ZkmV"+ (FVk)~ d a OZ
(22)
where the displacement and temperature sensitivities ~ and 0 are unknowns. To obtain an explicit expression for @ in terms of the velocity field V, ~ has to be eliminated. An adjoint equation is introduced by replacing ~ E Z in (22) by a virtual displacement ,~ E Z, and equating terms involving A to the energy bilinear form (11), yielding the adjoint equation for the elasticity problem,
a(A, A) =
~
?r'~(,~) + ~
,~" dO
for all ,~ E Z.
(23)
Notice that if A is replaced by ~, the right-hand side of (23) is identical with the first two terms on the right-hand side of (22). Thus, m k = .(,, +) + f ~OF [-8'~0 + +'J(-=+v ~) + 8'+/~o+v ~] - --~OF zk V + (FVk)k dO Zm (24)
Using the symmetry condition of the bilinear form (20) and substituting A for £, the variation of ~ becomes = a(~, ~) + f o~,~ oF - " - [-~'00
i
+ o'(-=~v ~) + ~'Oo,,V*l
oF Oz"
zk,. V ,, + (FVk), dO
OF ]dO
i k +,¢'(z, "" i k" o),~'(,)v**] + f [o.,,i](,)z~v~ e)~,,V~ o.~J(z, + f '',,a V ~ + f "x
V, d~
R.J. "fang, Shape design sensitivity analysis
47
OF ztv* + (FV*), da + f °-LF [6-'~(-z,y ') + ,~'J/3o,,V"l oz" Oo.iJ
+ fF ~ - t .i,,^ ki rvrk + [ ( t ~ ) k n k + H t ~ ] V ~ n '~ d r .
(25)
Notice that (25) is explicit and computable except for terms including 0. Use the same approach for £ and define an adjoint equation for the thermal part as follows:
A(A, .4) = ~:/3fi, AI - 0or,---~. ~ dD
for all ft. E Y.
(26)
Replacing 0 for A in (26) and substituting it into (25), the variation of ~ can be rewritten as ~, = A(A, 0)
+ f [6.0( A)ZkV ~ ik + tr'itz, 0) AkV k _ trOtz, O)~ij(A)V~] + f kAiV k + :~'qiVk ,, --k d/'] + f -.--OFi] [~iJ(__zkVk ) dp ~iJ~Okgk ] ~Gr
~F m k Oz m z k V + (Fvk)k d n
+ fr 2 -t'A~V k + [(t'A')kn k + HtiA']Vmn m dF.
(27)
Using the symmetry condition of the bilinear form (19) and substituting A for 0, the variation of ~ becomes
-- f k[(OkA, + O~Ak)V~ - O,A,V~] + gkAV k + gaV~ d12
+ f [~"(a)z'~v~ + ,,'J(,, O)A,Vj'* -crU(z, +
~
Ie,j(_z,y,, ) + 8,jOo~v~]
O),'~(,~)V~,I + f,a' 'V k +f'A'V~ dO
oF 7.k v~ + (FVk), d O t~zm m
+ fr~o q[--Ak Vk + (Ak nk + HA)Vmn m] dF
fr~ h[(O.Ak - OAk - O,A)V k + + fr*-t'A'~Vk
{[(0
-
O.)A]kn k + H(O - O.)A} Vmn m] d r
+ [(t'X')knk + Ht'A']vmnm dr.
(28)
Notice that (28) is a computable expression for the sensitivity of a general performance functional. The solution sequence is the temperature 0, the displacement vector z, the adjoint displacement vector A, and finally the adjoint temperature A.
48
R.J. Yang, Shape design sensitivity analysis
3.3. Implementation considerations
In the direct method, the sensitivities are obtained by solving, the heat transfer problem 0 of (2), the elasticity problem z of (9), the temperature sensitivity 0 of (19), and the displacement sensitivity ~ of (20). Notice that the sequence for design sensitivity analysis is identical with the structural analysis, namely, temperature and then displacement, and that the decomposed stiffness matrix can be reused for both the thermal and elasticity analyses. If a single load case is considered, the computational overhead for design sensitivity analysis is NV heat transfer and NV elasticity analyses, where NV is the number of the design variables. As the decomposed stiffness matrices can be reused from the analyses, only the forward and backward substitution is required. The adjoint variable method is different from the direct method in obtaining the sensitivities. The heat transfer and thermal stress analyses are first solved, the constraint functional is differentiated, and subsequently the adjoint displacement of (23) is solved before the adjoint temperature of (26). Similar to the direct method, the decomposed stiffness matrices for temperature and displacement can be reused for the sensitivity analysis. The sensitivity overhead in this case is N heat transfer and N elasticity analyses, where N is the number of the active constraints. Neither the adjoint variable nor the direct method can provide the most efficient computation if used alone. The most efficient way for computing the sensitivities is the hybrid approach [10, 11]. The hybrid approach compares the active constraints N and the design variable number NV. If the active constraint number N is smaller, as in the overdesigned case, the adjoint variable method is employed. On the other hand, if N is greater than NV, indicating the underdesigned case or the design is infeasible, the direct method is used. Another way. of obtaining the sensitivities is to use (25) for sensitivity calculation, in which the unknown 0 is solved by the temperature sensitivity equation of (19) and the A by the adjoint displacement equation of (26). The advantage of this strategy is that the programming efforts are minimal compared to the adjoint displacement-adjoint temporature approach. This stratogy can be considered as the sensitivity analysis for the elasticity part only, provided that 0 and 0 are known. The hybrid approach can also be applied at this level, i,e., after # and 0 are solved. Table I summarizes these three different strategies. The first level hybrid approach which is most efficient combines strategy l and strategy II, while the second level combines strategy l and strategy III. The overhead for strategy lII is that NV heat transfer sensitivity analyses have to be solved for each iteration. That requires more computation resources than the direct-direct or Table 1 Implementation strategies for sensitivityanalysis Strategy i StrategyII Strategy ill Heat transfer Direct Adjoint Direct sensitivity analysis (NV)* (N) (NV)
Elasticity Direct Adjoint Adjoint sensitivity analysis (NV) (N) (N) *Number in parentheses is the analysis count. N, number of active constraints; NV, number of design variables.
R.J. Yang, Shape design sensitivity analysis
49
adjoint-adjoint approaches. However, because the degrees-of-freedom for the temperature problem are fewer than for the stress problem, it may still be viable. 4. Numerical integration
The variational forms of the structural and sensitivity analyses derived in the previous section are not limited to any specific number method, although the finite element method is commonly used. A numerical integration scheme must be used to form the matrix equation, when the finite element method is employed. In general, Gauss quadrature is employed in a numerical integration scheme for the finite element analysis. It first transforms the integral to the natural coordinates and integrates over the coordinate system. For example, the functional in (15) can be transformed to the form
~, = f FJ d ~ ,
(29)
where l is the determinant of the Jacobian matrix which transforms the undeformed configuration into the natural coordinate system, and d ~ the volume increment in the natural coordinates. Haber [_12] found that (29) could be exploited for sensitivity analysis, as the volume increment d n is constant for all time or for all designs in this case. The variation of (29) is written as
(30)
4,, = f (F' + F,v*)J + FJ d / ) ,
where J is the derivative with respect to the design variable. Comparing (30) and (16), it is found that J d~] -- V~ dD. That relationship was pointed out and proved in [13]. Arora et al. [14] used the same 'reference volume' approach for design sensitivity analysis with nonlinear structural response. Using the same approach, the variation of a general functional ~2 on the boundary is
4,2= fr (F' + f,V*)J + FJ d r ,
(31)
where dr' is the surface increment in the natural coordinates. Comparing (31) and (18), it is found that (31) is simpler than (18), as (18) contains the curvature term H. Using the relationships in (30) and (31) and applying the same approach described in the previous section, the temperature sensitivity equation of (19) is rewritten as
A(O, O)- f kO~J dD + fr 2 hOOJdI" -
(q,v j + qj) dP + f.o h(O o) J dP -
(32)
50
R.J. Yang, Shape design sensitivity analysis
The displacement sensitivity equation of (20) is rewritten as
a(~,[)- f ~J(~)~J(~)J d.0 = f ~[or ' J (Z)zkV - ' i~ + crii(z, O)zkV -~ j~ + /3oz:i+ fikvki, ilJ + [f'E'-- ¢riJ(Z, O)eiS(f.)]J dfJ +
z (tkV J
-k t~j)£. ' d.P.
(33)
The variation of a general structurai response ~, is in the following form:
= f k[(OkA ~+ O~Ak)VkJ- O~A,J] + (gkVkJ + gJ)a d~2 +f[~j( +
, k + ~r~J(Z, A)zkVj ~ jk + f ~kV kA~ ]J + [fiA' O)AkV
"-'---oF 8o.iJ
[Oo(_z,V,) + 8 , q J o ~ v ~1
oz m Zk V
+ frio (q, VkJ+ q J ) A d [ ' + fr~o h(O®-O)AJd['+
-
o'q(Z, o)EiJ(A)]J d~
"b ( F V k ) k
<,:v'J
J
d~ dP
<34)
It is observed that the surface integral can be simplified significantly, however, the domain integrals are almost unchanged. 5. Numerical examples
Four examples are used to demonstrate and verify the analytical sensitivity expressions by using all three strategies described in Section 3, direct-direct, adjoint-adjoint and directadjoint. The first example verifies the thermal stress sensitivities subjected to a temperature change. The second and the third examples verify the thermal displacement sensitivities considering the heat conduction and convection effects. The final example is a threedimensional tube which is used to show the applicability of the technique to realistic problems.
5.1. One.dimot~ional thermal stress problem A one-dimensional tapered beam shown in Fig, 1 is the first example, The beam is heated from 0 °C to T = 1000 °C. The analytical stress field is obtained by the compatibility condition between the thermal expansion and the deformation due to the mechanical load. The expression is shown as follows: e:a r(d~ - dl)
or = - [d, + (d 2 - d,)xtl] Iog(dJd~) "
(35)
Selecting the heights of two end points as the design variables d I and d 2, the stress
R.J. Yang, Shape design sensitivity analysis
51
w Fig. 1. Clamped beam.
sensitivities are do" = - E a T
(d2/d~ - 1)[d, + (d 2 - d , ) x / l ] - d 2 l o g ( d 2 / d t ) [d, + (d 2 - d , ) x / I ] l o g ( d 2 / d , ) 2
do" = - E a T
(d'/d2-
dd 2
1)[d, + (d 2 - d , ) x / l ] + d, l o g ( d 2 / d , )
(36)
[d~ + (d 2 - d~)x/l] log(d2/d~) 2
The one-dimensional problem is solved by the finite element method and is modeled as a three-dimensional beam problem. Insulation boundary conditions are imposed on all surfaces except the two ends. The finite element model of the beam contains 68 grid points and 5 20-noded HEXA elements. The left and right ends of the beam are clamped. The length I and the width w of the beam are chosen as 50 mm and I mm, respectively. Young's modulus E, Poisson's ratio v, and the coefficient of linear thermal expansion o. for the beam are 2.07 × 10s MPa, 0.3 and 1.2 × 10 -7 °C -1, respectively. Heat transfer analysis is not required for this problem, as the temperature distribution is known. However, to check the heat transfer equation in (2), temperatures of 1000°C are specified at both ends. As expected, the temperature results in the beam are identical with the boundary temperatures. This indicates that the heat transfer finite element analysis is capable of describing a uniform temperature profile. The stress and stress sensitivities, from the exact solutions and the finite element analysis, are listed in Table 2. No stress error greater than 0.05% from finite element analysis is observed in this case. All three senslt "iv"i~y strategies: direct-direct, adjoint-adjoint and direct-adjoint show identical accuracy in which the error is below 4%, compared with the analytical solutions. 5.2. Heat c o n d u c t i o n
The second example is a one-dimensional heat conduction problem. The beam has an identical configuration to the first example, as shown in Fig. 1. In this example, the design variable is the length of the beam, while the heights at either end remain unchanged. The one-dimensional heat conduction problem is governed by
-kOxx = g ,
x~[O,/l, (37)
0(0)=0,
0(1)=T .
t~ tO
Table 2 Stress sensitivity c o m ~
for damped beam
(1) x (mm)
(2) Exact stress
(3) FEM stress
(4)" Error (%)
(5) b DV val~e
(6) DSAf exact
(7) DSAs I
(8) c Error (%)
(9) DSAh 11
Error (%)
(11) DSAi III
(12)" Error (%)
0
-27.25
-27.24
0.0
1(1.0) 2(1.2)
14.0~ -11.7
13.95 -11.63
0.6 0.6
13.95 -11.63
0.6 0.6
13.95 -11.63
0.6 0.6
10
-26.2
-26.19
0.0
1(1.0) 2(1.2)
20
-25.23
-25.22
0.0
1(1.0) 2(1.2)
(10) d
7.452 -6.21
7.37 -6.142
1.1 1.1
7.37 -6.142
1.1 1.1
7.37 -6.142
1.1 1.1
1.785 -1.487
!.716 -1.43
3.9 3.9
1.716 -1.43
3.9 3.9
1.716 -1.43
3.9 3.9 ¢b
30
-24.33
-24.32
0.0
2(1.2)
-3.106 2.588
-3.165 2.637
-1.9 -1.9
-3.165 2.637
-1.9 -1.9
-3.165 2.637
-1.9 -1.9
l(1.0)
r~
40
-23.49
-23.49
0.0
1(1.0) 2(I.2)
-7.338 6.115
-7.388 6.157
-0.7 -0.7
-7.388 6.157
-0.7 -0.7
-7.388 6.157
-0.7 -0.7
50
-22.71
-22.7
0.0
1(1.0) 2(1.2)
-11.01 9.174
-11.05 9.210
-0.4 -0.4
-11.05 9.210
-0.4 -0.4
-11.05 9.210
-0.4 -0.4
" ( 4 ) = (1 - 1(3)/(2)1) × i 0 0 ~ .
b Design variables (d~ and d2) and values. ( 8 ) - (1 - 1(7)1(6)1) x 1 0 0 % . ( 1 0 ) = (1 - l ( 9 ) / ( 6 ) l ) x 1 0 0 % . ( 1 2 ) = (1 - 1 ( 1 1 ) / ( 6 ) 1 ) x 100%.
f Exact solution from (36). s Strategy [: dire.ct-dire.ct. a Strategy If: adjoint-adjoint. ~Strategy IH: direct-adjoint.
zs z~
R.J. Yang, Shape design sensitivity analysis
53
The exact solution for the temperature is obtained by integrating (37) as
t +
7
•
-
(38)
The differential equation for the elasticity field and boundary conditions for the onedimensional case are On
E Ox
/30
Constant
x~[0,1] (39)
u(O)-O,
u(O-o.
The exact solution for displacements and stresses are integrated to obtain u
- ot [ gx3
or = - E a
T
gl 2
(40)
(r~" + 1"i2"k/"
The temperature, displacement and stress sensitivities with respect to the length I are dO dl du dl d~ dl
O0 OI -
au OI Oor O!
O0 Ox #x Ol b
g x ( I - x) kl T
Ou Ox Ox OI Oor Ox Ox O!
3g
gl
(41)
Eagl 6k
The conduction problem is also modeled as a three-dimensional beam problem, as shown in Fig. 1. Insulation boundary conditions are imposed on all edges, 0 °C and 1000 °C temperatures are applied at the left and the right ends, respectively. Both the left and the right ends of the beam are clamped. Young's modulus E, Poisson's ratio v, and the coefficient of linear thermal expansion a for the beam are 7.07 x 105 MPa, 0.3 and 1.2 x 10 -~ °C -!, respectively. The ambient temperature T® and fin base temperature are 0°C and 1000°C. The thermal conductivity k is assumed to be unity and the internal heat source g = 0. The temperature and its sensitivity, computed from the finite element analysis and from the analytical expressions of (38) and (41), are listed in Table 3. The finite elemen* results show excellent agreement with the analytical solutions. The displacements and the sensitivities, computed from finite element analysis and from the analytical expressions of (40) and (41) are listed in Table 4. The finite element results again show excellent agreement with the analytical solutions. The maximum error in displacement sensitivity is below 1%. The stress is 12.42 MPa both from the finite element and from the analytical solution. The stress sensitivity diminishes, as no internal heat source is applied.
R.J. Yang, Shape design sensitivity analysis
54
Table 3 Temperature and sensitivity for heat conduction ,,
T
dT/dd
X
(ram)
Exact
FEM
Exact
FEM
0 10 20 30 40 50
0 200 400 600 800 1000
0 200 400 600 800 1000
0.0 0.0 0.0 0.0 0.0 0.0
0.3333E-11 - 0.1584E- 10 0.7804E- 12 0.4822E- 11 0.1986E- 10 -0.3333E- 11
5.3. Heat convection
A straight thermal fin of rectangular profile is chosen as the third example. Figure 2 illustrates the configuration and the notations. The symbol I denotes the length of the fin, k the thermal conductivity, h the heat transfer coefficient at the exposed surface, To the fixed temperature at the fin base, and T~, the temperature of the bulk of the ambient fluid. Without losing any generality, the heat losses from the side edges of the fin are neglected and the heat flow is expressed per unit width. If it is assumed that the fin is sufficiently slender, the fin can be treated as a one-dimensional problem with no heat convected out of the free end. With the assumptions, the analytical temperature of the fin can is written as [15] O _- T - T~. _ cosh m ( l - x) 0o T o - T. cosh ml '
(42)
where m - - V ~ A ", C is the perimeter and A is the cross sectional area. For a wide fin, assuming no heat losses from the side edges is reasonable, and m -V'2"-h/-~, where d is the
height of the fin. For a freely expanded fin, the differential equation for the elasticity part is as follows: E OU
ox - t3o = o ,
x~[O,l],
(43)
x(O) = o .
,
a
i|r
T./~,,
k mm
i
....
x
III W
h : heat transfer coefficient k : thermal conductivity Fig. 2, Thermal fin of rectangular profile.
R.J. Yang, Shape design sensitivity analysis
55
The displacement is obtained by integrating (43) as u(x)
a[sinh m l
0o
-
sinh m ( l
-
x)]
(44)
m cosh m l
Selecting the height of the fin d as the design variable, the sensitivities for the temperature and the displacement are d(O/Oo)
O(O/Oo) a m am Od
dd
0m (l - x) sinh m ( l - x) cosh m l - I sinh m l cosh m(1 - x) = ad cosh 2 m l ' d(u/eo)
O(u/O0) Om
dd
Om
Od
= a { l m [ 1 - sinh m ( l - x) sinh ml]
- cosh m l [ m ( l - x) cosh m ( l - x) + sinh m l - sinh m ( l - x)]} I ( m cosh m l ) 2 0 m ad '
Om
1
2k~~
(45)
Od = - 2
The thermal fin problem is also modeled as a three-dimensional beam problem, as shown in Fig. 1. Insulation boundary conditions are imposed on side edges, and only the top and bottom surface are allowed heat convection. The left end of the fin is clamped. Young's modulus E, Poisson's ratio v, and the coefficient of linear thermal expansion a for the fin at 2.07 x l0 s MPa, 0.3 and 1.2 x 10 -7 oC=~, respectively. The ambient temperature T® and fin base temperature is 0 oC and 1000 oC. The thermal conductivity k and heat transfer coefficient h are assumed to be unity. The temperature and its sensitivity, computed from the finite element analysis and from the analytical expressions of (42) and (45), are listed in Table 5. The finite element results show excellent agreement with the analytical solutions. The assumption of no heat convected at the right end is proved to be reasonable, as the temperature is diminishing at the fin end. Table 5 Temperature and sensitivity for heat convection
T
dT/dd
X
(mm)
Exact
FEM
Exact
FE M
0 10 20 30 40 50
0.1000E + 4 0.2431E+3 0.5912E+2 0.1442E+2 0.3700E + 1 0.1699E + 1
0.1000E + 4 0.2448E+3 0.6032E+2 0.1481E+2 0.3856E + 1 0.1779E + 1
0.0000E +0 0.1719E+3 0.8366E+2 0.3073E+2 0.1076E + 2 0.6006E + 1
0.5106E-09 0.1671E+3 0.8257E+2 0.3058E+2 0.1085E+2 0.6097E + 1
0.0 -0.48E-3 -0.72E-3 -0.72E-3 -0.48E-3 0.0
0 10 20 30 40 50
0.0 -0.481b-E-3 -0.']7.03E-3 -0.7203E-3 -0.4805E-3 0.0
(3) FEM disp. 0.0 -0.96E-5 -0.144E-4 -0.144E-4 -0.96E-5 0.0
(5) DSA ~ exact
a ( 4 ) = ( 1 - 1(3)/(2)1) x b ( 7 ) = (1 - ] ( 6 ) / ( 5 ) ] ) x " ( 9 ) = (1 - 1(8)/(5)1) x d (11) =(1 -J(10)/(5)D 100%. 100%. x 100%.
100%.
0.OO 0.6441E-3 0.8025E-3 0.84~9E-3 0.8540E-3 0.8645E-3
0.0 0.6422E-3 0.7984E-3 0.8364E-3 0.8457E-3 0.8485E-3
10 20 30 40 50
0
FEM disp.
Exact disp.
x (mm)
(3)
(2)
(t)
0.0 -0.9591E-5 -0.1439E-4 -0.1439E-4 -0.9591E-5 0.0
0.0 O.1753E-3 0.3283E-3 0.3925E-3 0.4152E-3 0.4243E-3
DSA c exact
(5)
0.0 0.1816E-3 0.3341E-3 0.3998E-3 0.4229E- 3 0.4374E-3
DSA t 1
(6)
" E x a c t solution from ( 4 5 ) . f Strategy I: direct-direct. * Strateb~ II: adjoint-adjoint. b Strategy HI: direct-adjoint.
-0.3 -0.5 -0.9 - 1.0 -1.9
m
Error (%)
(4)"
-0.0 -3.6 -1.7 -1.9 -1.8 -3.1
(%)
(7) b Error
0.1 0.0 0.0 0.1
(%)
(7) b
Error
(6) DSA f I
Exact solution from (41). f Strategy i: direct-direct. * Strategy II: adjoint-adjoint. b Strategy HI: direct-adjoint.
m -0.1 0.0 0.0 -0.1 ~
(4) ~ Error (%)
Table 6 Displacement sensitivity comparison for heat convection
" ( 9 ) = (I - I(8)/(5)I)x i o o % . d (II) = (I - ] ( 1 0 ) I ( 5 ) D x 100%.
*'(7) = (I - I(6)/(5)l) x I o o % .
~'( 4 ) = (I - [(3)/(2)I) x I O O % .
(2) Exact disp.
(1) x (mm)
Table 4 Displacement scnsitivity comparison for heat c o m l ~
(8)
-0.0 0.1816E-3 0.3341E-3 0.3998E-3 0.4229E-3 0.4374E-3
DSA g II
(8)
0.0 -0.9591E-5 -0.1439E-4 -0.1439E-4 -0.9591E-5 0.0
DSA g II
-3.1
-1.8
-3.6 -1.7 -1.9
(%)
(9)" Error
0.1 0.0 0.0 0.1
(%)
(9) ~ Error
0.0 0.1816E-3 0.3341E-3 0.3998E-3 0.4229E-3 0.4374E-3
(10) DSA h III
0.0 -0.9591E-5 -0.1439E-4 -0.1439E-4 -0.9591E-5 0.0
(10) DSA" III
-3.6 -1.7 -1.9 -1.8 -3.1
(%)
(11) e Error
0.1 0.0 0.0 0.1
t%)
(11) e Error
e~
g~
¢,q
,o
¢Jt
R.J. Yang, Shape design sensitivity analysis
57
The displacements and the sensitivities, computed from finite element analysis and from the analytical expressions of (44) and (45), are listed in Table 6. The finite element results also show excellent agreement with the analytical solutions. The maximum error in displacement is below 4%. As the fin is freely expanded, no stress occurs in this case. The finite element results show very small values in the stresses and their sensitivities, which are omitted.
5.4. Three-dimensional tube problem A three-dimensional tube problem is shown in Fig. 3. The tube is idealized to simulate an engine exhaust manifold in which high thermal stress may initiate a fatigue crack and subsequently cause failure of the components. A finite element model is generated for design sensitivity analysis, as shown in Fig. 4. This model consists of 720 grids, 2110 degrees-of-freedom and 448 linear H E X A elements. Zero displacement conditions are imposed in the .x-direction of the left end (section A). One point on the surface is fixed in all three directions to eliminate rigid body motion. The rest of the points are allowed to move in the yz plane to simulate the bolt condition between the runner and the engine head. A fixed temperature condition of 780 °C is imposed at the left end. The thermal conductivity k and heat transfer coefficient h for the inner and outer surfaces of the tube are assumed to be unity. Young's modulus E, Poisson's ratio v and the coefficient of linear thermal expansion a for the tube are 2.07 x 105 MPa, 0.3 and 1.2 x 10 -7 °C -I, A
(01t ..... (ill
V
1 I
.
B rI r
design variables b(I) = r b(2) = r
0
b(3) = r I b(4) = d I b($):
db
Fig. 3. Three-dimensional tube.
o
R.J. Yang, Shape design sensitivity analysis
58
720 grids 448 HEXAelements 2110 D. O. F.
870
Fig. 4. Finite element model of three-dimensional tube.
respectively. The gas temperature in the tube is 2000 °C. The ambient temperature outside the tube is 0 °C. The x-displacement at grid 1202 (Fig. 4), Von Mises stress at grid 870 in element 340 and at grid 1202 in element 277 (Fig. 4), and their sensitivities are listed in Table 7. Note that the central finite difference results with 0.1% design perturbation is used to calculate the Table 7 Sensitivity comparison for a three-dimensional tube
( 1) ~oltem
(2) DV'
(3) FD
(8) DSA~ (Ill)
(9) d Error (%)
. . 2.5613E-4 3.7556E~4 9.7297Em5
. . 10 0.81 0,29
2.5614E-4 3.7556E~4 9.7297E-5
10 0.81 0,29
2.0343E~-2 -3,8482 6,2727 1,3381 -1,4674E-2
2.4 0.16 0,23 0,88 1,52
1(3.0) 1.9868E-2 2(1.5) m 3 . 8 5 4 4 3(1.0) 6.2585 4(10) 1.353 5(~12) -1.4901E-2
2.0344E-2 -3.8482 6.2726 1,3381 -1,4674E-2
2.4 2,0343E-2 0,16 -3,8482 0.23 6,2726 0,88 1,3381 1.52 -1,4672E-2
2,4 0.16 0,23 0,88 1.54
1(3.0) 2(1.5) 3(1.0) 4(10) 5(12)
O,1732 7.9142 6,2726 -0,65371
8,97 O,1732 1,62 7,9141 0,23 6, 2726 1.55 -0,65371
8.97 O.1732 1.62 7,9142 0,23 6, 2727 1,55 -0,65370
i
Design variable and value, b (5)~ (t -1(4)/(~)1) × 100%.
(7) = (1 -1(6)/(3)1) × 100%. d (9)--(i - 1(8)/(3)1) × 100%. Strategy I: direct-direct, f Strategy lh adjoint,adjoint, g Strategy IIl: direct-adjoint,
. .
(7) ~ Error (%)
Von Mises stress
O,15895 7.78S3 6, 2585 -0.65373
. .
(6) DSAf 11
1(3.0) 2(1,5) 3(1,0) 4(10) 5(12)
Von Mises stress at point 1202 (4.2809 MPa)
. . 2,5613E-4 3,7556E-4 9,7297E-5
(5) b Error (%)
x-di~t ,acement at point 1202 (1.5806E ~ 3 mm)
at point 870 (0.65788 MPa)
. . . . 2,3283E-4 3,7253Em4 9,7013E-5
(4) DSA~ I
10 0,81 0,29
8.97 1.62 0,23 1,55
R.J. Yang, Shape design sensitivity analysis
59
percentage errors in sensitivities, as no exact solution is available in this case. The results show that the sensitivities calculated from all three strategies have excellent agreement with the finite difference solutions. The sensitivities with relatively small magnitudes compared to others are neglected.
6. Summary Analytical expressions for shape design sensitivities of three-dimensional, linear thermoelasticity problems were derived, using the material derivative approach. Numerical integration for obtaining the sensitivities using the finite element method was discussed and the advantages are identified. It was found that the surface integral can be simplified significantly by using the natural coordinates. Three different strategies, direct-direct, adjoint-adjoint and direct-adjoint, were discussed and verified by four example problems. It was found that the sensitivities computed from all three strategies are identical and have quite good agreement with the analytical solutions. The error in displacement and stress is below 4%. It is also shown that the temperature and its sensitivity have excellent agreement with the analytical solutions. The technique is further verified by a three-dimensional tube example. The most efficient way for computing the sensitivities is the hybrid approach. Two levels of the hybrid approach were suggested. The first level hybrid approach, which is the most efficient, combines the direct-direct and adjoint-adjoint method, while the second level combines the direct-direct and direct-adjoint methods.
Acknowledgment The author wishes to express appreciation to the reviewers, Professor J.W. Hou of Old Dominion University, S. Poe of PTPE, D. Tao of PMES, D. Wagner, K. Yeung, A. Mahajan and J. Liu of CAE department for valuable suggestions and discussion.
References [1] R.A. Meric, Boundary element in shape design sensitivity analysis of thermoelastic solids, in: C.A. Mota Soares, ed., Computer Aided Optimal Design: Structural and Mechanical System, NATO ASI Series (Springer, Berlin, 1987). [2] K. Dems, Sensitivity analysis in thermoelasticity problems, in: C.A. Mota Soares, ed., Computer Aided Optimal Design: Structural and Mechanical System, NATO ASl Series (Springer, Berlin, 1987). [3] G. Hou, J.S. Sheen and C.H. Chuang, Shape sensitivity analysis and design optimization of linear, thermoelastic solids, in: Proc. 31st AIAA SDM Conf., Long Beach, CA, 1990. [4] D.A. Tortorelli, R.B. Haber and S.C-Y. Lu, Adjoint sensitivity analysis for nonlinear dynamic thermoelastic systems, AIAA J. 29 (1991) 253-263. [5] K. Dems and Z. Mroz, Variational approach by means of adjoint systems to structural optimization and .~'nsitivity analysis II, Structure shape variation, Internat. 3. Solids and Structures 20 (1984) 527-552. [6] E. 2. Haug, K.K. Choi and V. Komkov, Design Sensitivity Analysis of Structural Systems (Academic Press, New York, 1986).
60
R.J. Yang, Shape design sensitivity analysis
[7l N. Kikuchi, Finite Element Methods in Mechanics (Cambridge Univ. Press, New York, 1986). [81 T.J.R. Hughes, The FirSte Element Method (Prentice Hall, Englewood Cliffs, NJ, 1987). [91 J.A. Bennett and M.E. Botkin, eds., The Optimum Shape: Automated Structural Design (Plenum, New York, 1986).
[101 RJ. Yang, A three-dimensional shape optimization system-SHOP3D, Comput. & Structures 31 (1989) 881-890.
[lll R.J. Yang, A hybrid approach for shape optimization, in: V.A. Tipnis and E.M. Patton, eds., Computers in Engineering (ASME, New York, 1988) 107-112.
[12] R.B. Haber, A new variational approach to structural shape design sensitivity analysis, in: C.A. Mota Soares, ed., Computer Aided Optimal Design: Structural and Mechanical System, NATO ASI Series (Springer, Berlin, 1987). [13] R.J. Yang hnd M.E. Botkin, Accuracy of the domain material derivative approach to shape design sensitivities, AIAA J. 25 (1987) 1606-1610. [14] J.B. Cardoso and J.S. Arora, Variational method for design sensitivity analysis in nonlinear structural mechanics, AIAA J. 26 (1988) 595-603. [lSl A.J. Chapman, Fundamentals of Heat Transfer (Macmillan, New York, 1987).