Numerical accuracy for the boundary and domain methods of the material derivative approach to shape design sensitivities is investigated through the use of mesh refinement. The results show that the domain method is generally more accurate than the boundary method, using the finite element technique. It is also shown that the domain method is equivalent, under certain assumptions, to the implicit differentiation approach not only theoretically but also numerically

Botkin, M. E.

Yang, R. J.

English

NASA Technical Reports Server

N87-18877
Accuracy of the Domain Method for the Material Derivative
Approach to Shape Design Sensitivities
R. J. Yang and M. E. Botkin
General Motors Research Laboratories
Warren, Michigan
Abstract
Numerical accuracy for the boundary and domain methods of the material derivative
approach to shape design sensitivities is investigated through the use of mesh refinement.
The results show that the domain method is generally more accurate than the boundary
method, using the finite element technique. It is also shown that the domain method is
equivalent, under certain assumptions, to the implicit differentiation approach not only
theoretically but also numerically.
Introduction
Haug and Choi et al. 1-4 developed a unified theory of structural design sensitivity
analysis for linear elastic structures, using a variational formulation of the structural state
equations. This theory allows one to take the total derivative, or material derivative, of
the variational state equation and to use an adjoint variable technique for design sensitiv-
ity analysis. The main attraction of this approach is that it allows one to compute the
derivatives of structural performances analytically. No discretization approximations are
involved during the derivation, and a step size need not be specified in the calculation.
However, the formulation requires evaluating accurate stress quantities on the boundaries
which are often difficult to obtain.
Accuracy of the shape design sensitivity theory was studied in Ref. 5 through the
equilibrium condition for different types of finite elements. However, a systematic study
through the refinement of the finite element mesh was still not found in the literature.
To improve the accuracy of shape design sensitivities. Choi et al. 6 proposed a new
domain method which transforms the boundary integrals into domain integrals and there-
fore is less influenced by the the inaccurate boundary stress evaluation. This method takes
advantage of the averaging nature of the finite element method, and is found to be more
accurate than the boundary approach which evaluates the derivatives using boundary in-
formation only 1-3 . However. a velocity field for the physical domain needs to be defined.
The necessity of defining the domain velocity may indicate that this method is closely
related to the implicit differentiation approach which also requires knowledge about the
domain change for differentiation of the elemental stiffness matrix 7
In this paper, the accuracy of the design sensitivity is studied through finite element
mesh refinement for a cantilever thin plate. Results of the domain and boundary methods
for the material derivative approach and the implicit differentiation approach are shown
and compared.
In a previous paper 8, the boundary integral formulation was shown to be equivalent
to the implicit differentiation approach, theoretically. In this report, the domain method
is shown to be equivalent to the implicit differentiation method, both in theoretical and
numerical aspects.
Shape Design Sensitivity Analysis
Two approaches for shape design sensitivity analysis are found in the literature. One is
the well known implicit differentiation approach and the other is the variational or material
FF£OE)iNG PAGE BtANK _ R_.MED
347
https://ntrs.nasa.gov/search.jsp?R=19870009444 2020-03-20T11:33:12+00:00Z
derivative approach. Detailed information for these two approaches is available in Refs. 4,
7, and 8. Only a brief background is provided in the following.
For the implicit differentiation approach, the displacement sensitivity is obtained by
differentiating the discretized structural system of equations
gz=F (1)
By assuming that the force vector F is independent of design, this leads to
Oz OK
Obi - K- 1-- = (2)
where K is the global stiffness matrix, z is the displacement vector, and bi is the design
variable.
The variational design sensitivity theory uses the material derivative concept of con-
tinuum mechanics and an adjoint variable method to obtain computable expressions for
the effect of shape variation on the functionals arising in the shape design problem. The
variation of the displacement functional _b with respect to shape change is derived by differ-
entiating the variational equilibrium equation and employing the adjoint variable method,
to obtain 1--4
i_ = [_ ai3(z)eiJ()_)Vknk dr (3)
./1"
where ¢ is defined by
= fa z$(x-_) d_ (4)
in which _ is the point of interest,/5 the Dirac-measure at zero, 12 the physical domain, a ij
and e i_ the stress and strain tensors, respectively. V the design perturbation and can be
thought of as velocity, and n k the unit normal vector of moving boundary r. The vectors
z and A are the displacement vectors for state and adjoint equations, respectively, which
can be expressed as follows:
/_1 aiJ(z)eq (-2)dfl = fr2 Ti_i dr (5)
f aiJ(A)eii(-A) dfl= fa b(x-_)Adfl (6)
where T, is a traction vector, 1"2 a loaded boundary, and - indicates the virtual displace-
ments that satisfy the kinematically admissible displacement field. The Einstein summa-
tion convention for a repeated index is used throughout this paper.
To obtain Eq. 3, the traction vector T,, the kinematically constrained boundary, and
the loaded boundary are assumed to be fixed during the design process, i.e., they are
independent of design. Note that the variation of the displacement functional of Eq. 3 is
only affected by the normal movement of the boundary of the physical domain.
Physically, the adjoint equation of Eq. 6 is interpreted by applying a unit load at the
poinl 5, where the displacement sensitivity is of interest. In Eq. 3, one sees that only
the boundary stress information is needed for evaluating the variation of the displacement
functional. Unfortunately, finite element analysis usually does not provide high quality
stress results, especially on the boundary. It has been shown that better finite element
results give better design sensitivity estimates, by examining the equilibrium equations for
different finite elements 5
348
Domain Method
The basic idea for the domain integral method is to take advantage of the averaging
nature of the finite element method, instead of evaluating the less accurate stresses on the
boundary. Since the finite element method is well known to provide better solutions inside
the finite element, the domain method has the advantage of predicting better sensitivities.
Applying the same procedure as in obtaining Eq. 3 with the domain method, the first
variation of the displacement constraint functional of Eq. 4, is obtained as 4,6
(7)
One should note that Eq. 7 is more general than Eq. 3, since only the loaded boundary
is assumed to be fixed, while both loaded and kinematically constrained boundaries are
assumed to be unchanged in Eq. 3. The kinematically constrained boundary and interface
boundary terms appear when the divergence theorem is used to transform the domain
integral to the boundary integral. It was shown in Ref. 6 that for an interface or built-up
structure problem this method simplifies the formulation and avoids specifying tedious
interface conditions and provides increased accuracy for shape design sensitivities.
To have a better understanding of Eq. 7, each term will be examined individually.
First, since the stress tensor a i3 is symmetric, the first term of Eq. 7 is divided into two
parts and then integrated by parts to obtain
l [z_,kVk,j + zj,kVk,i] df_
1 [z,,k  nj + =,,kVkn,]er= fr (8)
By assuming that only the free boundary is varied, the first term of Eq. 8 disappears,
since the traction vector is zero, i.e., oiJ(A)nj = ai_(A)ni = 0. The second term of Eq. 8
then can be further modified to
/o da (9)
where the velocity Vk can be parametrized as (c3xk/Obi)6bi, in which xk is the position
vector and bi is the design variable. Since all the integrals are linear in design, one can
eliminate 6b i or choose the value as a unit number. By interpreting the adjoint variable
A as the inverse of the reduced stiffness matrix K if all the displacement sensitivities are
desired, Eq. 9 is discretized, using the finite element formulation, as s
/f_ N_ /_ T ¢9xkaO(A)e'_(z)Vk, df_ = K-I __ B DB,k-_i z df_
1 e
Ne
= K-1 _ fa BTDB_z dfl
1 e
(10)
349
where the subscript i with a prime superscript indicates the derivative with respect to
the i th design variable. K, D, and B represent the stiffness, elasticity, and strain recovery
matrices, respectively. The stress-displacement and strain-displacement relationships are
employed in obtaining Eq. 10, which are defined in the following
a'J(A) = DBA (11)
 iJ(z) = Bz
Finally, the first term of Eq. 7 is written as
The second term of Eq. 7 can also be derived in the same way to obtain a similar
expression as in Eq. 12 as
Ne
1 e
Substituting Eqs. 12 and 13 into Eq. 7, the following expressions are obtained
_z
= ], [eJ(_)z,,_V_,, + o'J(z)_,,_V_,,.- e"(z)_'"(A)V_,_] da
N_
= -K-' Z f_ BTDB:z + BT'DBz + B TDBzVk,k df_
1 •
N_
:-.-'z(io
I e
(14)
where [ J [ is the determinant of the Jacobian matrix J which transforms the undeformed
configuration into the natural coordinate system. The constraint functional change $_ is
equal to Oz/Ob_, since all displacement sensitivities are calculated and the design change
_b_ is chosen as a unit number.
For the implicit differentiation approach of Eq. 2, the derivative of the global stiffness
matrix can be evaluated at the elemental level, i.e.,
Oz _K_IOK
z (i5)
350
where K_' is calculated numerically in natural coordinates as 7
K e'-- f__l/__1 f__l [ BTIDB [ g ] At-B TDB' [ J I .+B TDB J i_ ] d_dTId_ (16)
Comparing Eqs. 14 and 15, one sees that both are equivalent if the following expression
is valid.
I g I_=l J [ Vk,k (17)
To prove Eq. 17 is true, one should notice that the determinant of the Jacobian can
be separated into two parts. The first contribution is from the deformed to undeformed
configuration, denoted by I J Id, which depends on design. The other is the contribution
of transformation from the undeformed global to local or natural configurations, denoted
by I J 1, which is independent of design. The relationship is expressed in the following
form
[J[r=l Jid I gt (18)
where r denotes the deformed configuration. Differentiating Eq. 18 with respect to design,
one obtains
I J I'=i J I_ I JI (19)
It was shown in Ref. 4 that [ J [_= Vk,k at r = 0, if the design change is assumed to be
equal to a unit vector. Thus, Eq. 19 is identical to the form of Eq. 17, and the equivalence
of Eqs. 14 and 15 is proved.
Another way to prove the validity of Eq. 17 is to carry out the differentiation directly
by the definition of the Jacobian. Consider a two-dimensional case as an example, the
right side of Eq. 17 is obained as
I g l Vk'k = c3ec3 7 c%?_ _--Ox- + _-] (20)
where the velocity _ and Vy are defined as
_X
by
(21)
Substituting Eq. 21 into Eq. 20, the following expression is obtained
O (Oxi)y OyOx)[ J IV k'k -- i)bi O( Or] 0_-_
=1J I;
(22)
This simple calculation also verifies that the relationship of Eq. 17 is valid. Note that the
design change 6b, is assumed to be unity in Eq. 17.
351
Numerical Verification and Comparison
In this section, the equivalence of the domain method and the implicit differentiation
is verified through a simple example. The accuracy of design sensitivities will be examined
and compared through the refinement of the finite element mesh.
A simple two-dimensional thin plate is considered as an example. The finite element
configuration, dimensions, material properties, loading condition, and design variable are
shown in Fig. 1. Design variable b is chosen to move the upper traction free boundary.
The load of 100 lb is applied parabolically at the right of the plate.
An 8-noded two-dimensional plane stress isoparametrie element is employed for analy-
sis. The boundary stresses and strains that appear in Eq. 3 are computed by extrapolat-
ing linearly from the stresses at Gauss points, where the optimal or the best approximate
stresses are located. Numerical results for design sensitivity of point A in the Y-direction
for lxl, 2x2, 3x3.4x4, 5x5 and 6x6 meshes are shown in Table 1.
In Table 1, column 1 represents different finite element meshes and column 2 the
displacement of point A in the Y-direction for the initial design, b. Columns 3 and 4
represent the displacement sensitivities at point A of Fig. 1, using the boundary method
(BM) of Eq. 3 and the domain integral (DM) of Eq. 7, respectively, for different meshes.
Column 5 has the results using the implicit differentiation approach (IDA) of Eq. 2. The
derivative of the global stiffness matrix is carried out by differentiating the element stiffness
matrix, analytically.
Fig. 2 shows the same results as in Table 1. From Fig. 2 and Table 1, one observes that
the displacements and the sensitivities for the implicit approach (IDA) do not change much
after 3x3 finite element mesh. However, the design sensitivity for the boundary method
of the variational approach (BM) is still increasing at the limit of mesh refinement. This
implies that the boundary method (BM) is more sensitive to the finite element results,
although it provides the analytical formulation for sensitivities. And one concludes that
the boundary method of the variational approach tends to yield better gradient estimates,
once a more accurate analysis is used and better boundary stresses are obtained. The
same conclusion is also found in Refs. 5 and 8.
Comparing column 4 with 5, one sees that the domain method results (DM) are very
close to those obtained from the analytical implicit differentiation approach (IDA). Clearer
interpretation can be observed from Fig. 2, which plots the displacement and displace-
ment sensitivity versus finite element mesh size. This numerical agreement verifies the
equivalence of the two approaches.
In Ref. 8. the boundary the method was shownto be theoretically equivalent to the
implicit approach, however, they yield slightly different results numerically as also shown
in Table 1 and Fig. 2. The difference results from different numerical schemes for these two
approaches, i.e., one uses the boundary, and the other the domain information. If consistent
numerical schemes are used for the domain method and the implicit approach as in this
report, they are shown to be equivalent, not only theoretically but also numerically.
The disadvantage of the domain method is in computational aspects. Numerical eval-
uation of Eq. 7 is more complicated than evaluation of Eq. 3, because Eq. 7 requires
integration over the entire domain, whereas Eq. 3 requires integration only over the vari-
able boundary. In addition, a velocity field which satisfies regularity properties must be
defined in the domain. The choice of velocity for the physical domain is more difficult than
that for the variable boundary. Although, a boundary layer scheme 9 and a displacement-
like velocity field 10 were proposed to alleviate these problems, the domain method still
requires more analyst and computational efforts.
352
Summary
It is shown that accurate finite element analysis results in accurate design sensitivities.
For the boundary method of the material derivative approach to shape design sensitivities,
the accuracy of the finite element is more crucial, since the finite element method usually
does not give accurate stresses on the boundary.
The domain method is generally more accurate than the boundary method in the
material derivative approach for evaluating the design sensitivities; however, a velocity
field for the physical domain needs to be defined. The necessity of defining a domain
velocity field and integrating the domain integral instead of the boundary integral, as in
boundary method, requires both more analyst time and computational time.
It is also shown that the domain method is equivalent, under certain assumptions,
to the implicit differentiation approach not only theoretically but also numerically. The
numerical equivalence is valid only if the numerical methods used for both approaches are
consistent.
353
References
1. Haug, E. J., Choi, K. K., Hou, J. W., and Yoo, Y. M., "A Variational Method for
Shape Optimal Design of Elastic Structures," New Directions in Optimum Structural
Design, Ed. E. Atrek et al., Wiley, New York, 1984.
2. Choi, K. K. and Haug, E. J., "Shape Design Sensitivity Analysis of Elastic Structures,"
Journal of Structural Mechanics, 11(2), pp. 231-269, 1983.
3. Choi. K. K., "Shape Design Sensitivity Analysis of Displacement and Stress Con-
straints," Journal of Structural Mechanics, 13(1), pp. 27-41, 1985.
4. Haug, E. J.. Choi, K. K., and Komkov, V., Design Sensitivity Analysis of Structural
Systems. Academic Press, 1986.
5. Yang, R. J. and Choi, K. K., "Accuracy of Finite Element Based Shape Design Sensi-
tivity Analysis." Journal of Structural Mechanics, 13(2), pp. 223-239, 1985.
6. Choi, K. K. and Seong, H. G., "A Domain Method for Shape Design Sensitivity Analy-
sis of Built-up Structures," Computer Methods in Applied Mechanics and Engineering,
Vol. 57. No. 1, pp. 1-16, 1986.
7. Ramakrishnam C. V. and Francavilla, A., "Structural Shape OptimizationUsing Penal-
ty Functions," Journal of Structural Mechanics, 3(4), pp. 403-422, 1974-1975.
8. Yang, R. J. and Botkin, M. E., "Comparison between the Variational and Implicit
Differentiation Approaches to Shape Design Sensitivities," AIAA, Vol. 24, No. 6, pp.
1027-1032, 1986.
9. Haug, E. J. and Choi, K. K., "Material Derivative Methods for Shape Design Sensitivity
Analysis," The Optimum Shape: Automated Structural Design, Ed. J. A. Bennett and
M. E. Botkin, 1986.
10. Choi, K. K. and Seong, H. G., "A Numerical Method for Shape Design Sensitivity
Analysis and Optimization of Built-up Structures," The Optimum Shape: Automated
Structural Design, Ed. J. A. Bennett and M. E. Botkin, 1986.
354
Table 1. Comparison of Design Sensitivity Accuracy
BM DM IDA
mesh displacement dv/db dv/db dv/db
lxl 2.495E-5 -4.196E-6 -4.843E-6 -5.248E-6
2x2 2.760E-5 -4.518E-6 -5.167E-6 -5.235E-6
3x3 2.841E-5 -4.856E-6 -5.369E-6 -5.375E-6
4x4 2.845E-5 -4.995E-6 -5.391E-6 -5.394E-6
5x5 2.854E-5 -5.093E-6 -5.412E-6 -5.413E-6
6x6 2.856E-5 -5.158E-6 -5.425E-6 -5.426E-6
b
/
/
/
/
/
i
E • 1.0 x 10? psi
V-0.3
_A
10"
I00 Ib
I0" Ixl 2x2
3x3 4x4
finite element meshes
(8-nodeplane stress element)
Y
Fig. 1 Square Thin Plate
355
,4--w
e-.
o
r-
E
om
_' JI •
=
5.5-
4.5
4
3.5
3
2.5
_°---.D ..........
/
Z_ INITIAL DISPLACEMENT
x 0iS,ONsr_smvn_A_
o 0_S,0NsrNsmvrrY
0 DESIGN SENSmVITY (0M)
A A A
io i5 20 2'5 i0
Finite Element Number
Fig. 2 Accuracy of Design Sensitivity
3'5 40
356


Accuracy of the domain method for the material derivative approach to shape design sensitivities

Abstract

Similar works

Full text

Available Versions

NASA Technical Reports Server