https://ntrs.nasa.gov/search.jsp?

R=19850025225 2019-11-25T04:17:25+00:00Z

NASA Technical Memorandum 86382

. : NASA-TM-86382 19850025225
:
!
;
I
\
.'-~~---

Selecting Step Sizes _

in Sensitivity Analysis
by Finite Differences

Jocelyn Iott, Raphael T. Haftka,

and Howard M. Adelman

AUGUST 1985

NI\S/\ e"!!II!I!I!! ~ !!~I !!~ljl!!IIIIII!III!!I- j

NASA Technical Memorandum 86382

Selecting Step Sizes

in Sensitivity Analysis
by Finite Differences

Jocelyn Iott and Raphael T. Haftka

Virginia Polytechnic Institute and State University
Blacksburg, Virginia

Howard M. Adelman
Langley Research Center
Hampton, Virginia

NI\S/\
National Aeronautics
and Space Administration
Scientific and Technical
Information Branch

1985
Summary of subtracting large numbers which are very nearly
equal. Condition errors in finite difference derivatives
This paper deals with methods for obtaining near-
generally increase with decreasing step sizes.
optimum step sizes for finite difference approxima-
In the absence of methods for selecting a suitable
tions to first derivatives with particular application
step size, a trial-and-error approach, which requires
to sensitivity analysis. A technique denoted the fi-
many function evaluations for each derivative, must
nite difference (FD) algorithm, previously described
be used (ref. 1). Automatically calculating an op-
in the literature, is reviewed and applied in a tuto-
timum step size can decrease the time and cost of
rial manner to the derivative of a sine function. The
the calculations. One such method is described in
original FD algorithm is applicable to one derivative
reference 1. This method, referred to as the finite
at a time. This paper describes an extension of the
. difference (FD) algorithm, is formulated on the basis
FD method to the calculation of several derivatives
that the finite difference error is minimized when the
of matrix equations.
truncation error is equal to the condition error. The
Both the original and extended FD algorithms
FD algorithm has previously been demonstrated for
were applied to sensitivity analysis for a data-fitting
single derivatives and proved to be effective in pro-
problem in which a polynomial is passed through
ducing step sizes that gave good approximations to
data points, and the derivatives of the coefficients of
the exact derivatives (refs. 1 and 2). The need to
the polynomial with respect to uncertainities in the
apply the method to functions governed by matrix
data were calculated. The methods also were applied
equations has led to an extension of the FD algo-
to sensitivity analysis of the structural response of
rithm and is the main subject of this paper.
a finite-element-modeled swept wing. In a previous
In this paper, the FD algorithm is applied to
paper, this sensitivity analysis of the swept wing
three examples: (1) the sine function, (2) a matrix
required a time-consuming trial-and-error effort to
equation (ref. 3, pp. 12-13) involving data fitting by
obtain a suitable step size, but it proved to be a
a polynomial in which derivatives of coefficients with
routine application for the FD algorithm herein.
respect to uncertainties in the data are calculated,
and (3) a finite element structural analysis of a
Introduction swept wing. The polynomial and wing problems
One of the key areas in the development of meth- are used to test both the original and the extended
ods for optimization of engineering systems is sensi- versions of the FD algorithm. The finite element
tivity analysis-the calculation of derivatives of sys- model of the swept wing is included in this study
tem response with respect to design parameters. In because in a previous sensitivity study, described in
gradient-based optimization methods, derivatives are reference 4, the uncertainty of step size became a
used to determine design changes to move toward major concern, and a succession of trials was required
an optimum design. Sensitivity derivatives also are to determine an acceptable step size. This paper
used for rapid assessment of the effect of a change in contains a summary of the original FD algorithm,
a design parameter on the response of an engineer- the development of the extended FD algorithm, and
ing system. In these typical applications, especially the results of applying the algorithms to the above
when a large number of repetitive calculations are re- problems.
quired, the finite difference method is often a useful
technique for computing the derivatives. Symbols
One problem which arises in connection with the
finite difference method is the choice of step size. ai coefficient of ith term of Nth-order
The step size can contribute to two types of errors polynomial
in the finite difference approximation-truncation er- C (h) condition error in approximating first
ror and condition error. Truncation error is the dif- derivative by finite difference method
ference between the exact value and the computed
value of a perturbed function if it is assumed that c( cp) condition error in approximating
there is no loss of numerical precision in the calcu- second derivative by finite difference
lation. Truncation error is generally represented by method
neglected terms in the Taylor series representation
of a perturbed function. Truncation errors are in- e bound on total error in approximating
creased by overly large step sizes. Condition errors derivatives by finite difference method
are associated with numerical noise and are caused
F, F force vectors
by loss of numerical precision. These errors may re-
sult from computer round-off error or the operation f exact value of a function
 I' exact first derivative of a function One contribution to the condition error is round-off
error in evaluating equation (1), which is compara-
I" exact second derivative of a function
tively small for most mainframe computers and is,
I~' bound on second derivative therefore, considered negligible. However, if I (x) is
computed by a lengthy or ill-conditioned numerical
h step size process, the round-off contribution to the condition
hopt optimum step size error can be substantial. Additionally, condition er-
rors may result if I(x) is calculated by an iterative
hs step size used for estimating second process which is terminated early. The condition er-
derivative ror C(h) can be conservatively estimated as being
K stiffness matrix inversely proportional to the step size. The following
expression, used herein, is given in reference 2:
N degree of polynomial
T(h) truncation error in approximating first (3)
derivative by finite difference method
u,u exact and computed displacement
vectors, respectively Determination of fA
x independent variable In equation (3), fA is a bound on the absolute
Xi, Yi ith data point in Vandermonde error in the computed function I. Obtaining an ap-
problem propriate value of fA for use in equation (3) gen-
erally requires some knowledge of the problem be-
fA bound on absolute error in computed ing solved. A straightforward (but not always con-
function values venient) method of obtaining fA is to perform the
<I> finite difference approximation to calculation of I in single and double precision and
second derivative subtract the results. (See ref. 2 for examples of this
technique.) An alternative method for matrix prob-
Subscripts: lems is described in appendix A of this paper. A
max maximum method for estimating fA also is given in the ap-
pendix of reference 5.
mill minimum
Development of Optimum Step Size Formula
Summary of FD Algorithm
A bound e on the total error is the sum of the
Truncation and Condition Errors truncation and condition errors
The FD algorithm (ref. 2) is used to calculate an
optimum step size for a first derivative approximated (4)
by forward differences
where I/~'I is a bound on the second derivative in the
I'(x) = I(x + h) - I(x) (1) interval [x,x+hj. In many cases, such a bound on
h the second derivative is not available, and a central
Whenever finite difference formulas are used to ap- finite difference approximation
proximate derivatives, there are two sources of error:
truncation error and condition error. The trunca- ,,, '" ~ _ I(x + hs ) - 2 I(x) + I(x - hs )
Jb ",'j!- 2 (5a)
tion error T(h) is a result of the neglected terms in hs
the Taylor series expansion of the perturbed func-
tion and is represented in references 1 and 2 as being is used, so that
proportional to the step size as follows:
h 2
h
e= -
2
1<1>1 + -fA
h
(5b)
T(h) = "2 I" (x + ~h) (2)
It should be noted that equation (5a) is an approx-
The condition error is the difference between the nu- imation to" I"(x) and is a reasonable bound on I"
merical evaluation of a function and the"exact value. in the interval [x, x+h] only if I" varies slowly

2
near x. The step size hs used for estimating the sec- equation (10) for <I> in equation (7). This results in
ond derivative is not necessarily the same as that
used for estimating the first derivative. The optimum
step size is found by minimizing the error defined in
C!L
hs = 2(1 + IxDy 1 + III (11)
equation (4). Taking the derivative of e with respect
to h, substituting <I> for It',
and setting the resulting The interval hs is used to compute <I> from equa-
expression equal to zero gives tion (5a)j <1>, in turn, is used to calculate the er-
ror 2(<1» from equation (9). If the error is too large
8e 1<1>1 2 (greater than 0.1), hs is increased. If it is too small
- = - - -fA =0 (6)
8h 2 h2 (less than 0.001), hs is decreased. Once a value of hs
is found which produces an acceptable value for <1>,
Solving for h yields the value of <I> is used in equation (7) to produce an
optimum step size, hopt. The limit of 0.1 is chosen to
hopt = 2y~
.f!X (7) prevent the condition error from being excessive. The
limit of 0.001 is chosen because when the condition
error is very small, there is a risk of high truncation
Clearly equation (7) is not applicable to deriva- errors.
tives of functions at points where the second deriva-
tive is zero or close to zero-for example, when the Applications of FD Algorithm
function is nearly constant, linear, or an odd func-
tion. For these special cases, the FD algorithm fails Sine Function
to provide a satisfactory finite difference step size. As A simple example to demonstrate the FD algo-
discussed in reference 2, a computerized implementa- rithm is sin x. This application is described in a
tion of the algorithm terminates with an appropriate
tutorial fashion in appendix B to illustrate the steps
error message when 1<1>1 is less than a certain thresh- in the procedure.
old value. Equation (6) leads to the observation that
Polynomial-Fit (Vandermonde) Problem
hopt 1""1_
--
2
'J." - - - f A (8)
2 hopt The polynomial-fit example is based on the prob-
lem of passing an Nth-order polynomial through
Therefore, the optimum step size is the one which N + 1 data points, Xi = i, Yi = I(xi),
balances the truncation and condition errors. i = 0, 1, 2, ... , N. The form of the polynomial is

N
Estimation of Second Derivative
I(x) = Eaixi
Before evaluating equation (7) for h opt , the sec- i=O
ond derivative is evaluated from equation (5a). This
requires a value for hs, which is obtained by the fol- To solve for the coefficients ai, a system of N + 1
lowing trial-and-error procedure. The relative condi- equations in N + 1 unknowns is generated. It is
tion error in !" resulting from use of equation (5a) required to calculate the derivatives of the coefficients
is shown in references 1 and 2 to be bounded by ai with respect to the data point locations Xi. The
equation for computing the coefficients is
(9) 1 x 02 x 30 x 0N
Xo ao Yo
1 Xl x21 x31 xN 1 al Yl
The choice of an initial guess for hs is based on the
conditions that the function is well scaled and that 1 X2 x22 x32 xN 2 a2 Y2
x2 xN
=
the function and its second derivative are comparable 1 x3 3 x33 3
a3 Y3
in magnitude, such that

I" = 1 + III (10) 1 xN xJv x1v x~ aN YN

(1 + IxD 2 (12)
The matrix in equation (12) is referred to as the
where Ixl is of the order unity. The trial value of hs to Vandermonde matrix and is useful in the present
be used in equation (5a) is calculated by substituting work because the range of values (from unity to NN)

3
 TABLE 1. FINITE DIFFERENCE EVALUATION OF·
causes it to become ill-conditioned even for relatively
8all8x3 FOR POLYNOMIAL-FIT PROBLEM
small values of N. For example, N = 10 generates
values from 1 to lO lD .
In this example, the value for N was selected to be Derivative
12 because it produced condition errors which were Step size, Finite difference Derivative error,
significant but not catastrophic. The values of Xo h derivative ·errora percent
and Yo were set to zero, so that ao = O. The FD 10 1 -5.59003 X 106 1.23 X 106 28.2108
algorithm was used to find an optimum step size to 10- 2 -4.46303 1.03 X 10 5 2.3624
use in approximating the derivatives of each coeffi- 10- 3 -4.37003 1.00 X 104 0.2294
cient of the polynomial with respect to the location 1.71 X 10- 4 (hopt) -4.36028 0.25 X 10 3 0.0057
of the third data point, i.e., aai! aX3, i = 1, 2, ... , 12. 10- 4 -4.35313 0.69 X 104 0.1583
The data points Yi were generated from the equation 10- 5 -4.36095 0.92 X 103 0.0211
Yi = i + i 5 + i 9 at values of i ranging from 0 through 10- 6 -4.34963 1.04 X 104 0.2385
12. This method of generating data gives a means for 10- 7 -4.32533 3.47 X 104 0.7959
estimating fA because the polynomial which exactly 10- 8 -4.05803 3.02 X 10 5 6.9266
fits these data has coefficients al = a5 = a9 = 1.0, 10- 9 -5.31903 9.59 X 105 .21.9953
and the remaining coefficients should be 0.0. Any
deviation in the values of the coefficients is the result aRelative to exact value (-4.36003 X 106).
of numerical imprecision. A standard matrix subrou-
tine denoted MATOPS was used to solve the system
of equations. MATOPS offers three options (full piv- It was of interest to evaluate the error estimate
oting, partial pivoting, and no pivoting) correspond- of equation (5b) for this problem because the exact
ing to three levels of accuracy. The least accurate error was available. Using equations (5b) and (8)
option was used for evaluating derivatives, and a first gives the bound for the minimum error emin as
estimate for fA was obtained by subtracting the least
accurate solutions from the most accurate. The sec- (13)
ond method used for estimating fA was the iterative
improvement technique described in appendix A.
The two estimates of fA were very close, but both
So for fA = 0.153 and hopt = 1.71 X 10- 4 , emin =
were sensitive to the values of Xi. For the nominal
3580. This value is fourteen times larger than the ac-
point where Xi = i, fA was evaluated as 8.3 X 10- 4 .
tual error. To resolve this large difference between
However, at the nearby point of Xi = 0.01 +i, fA was
the error predicted by the FD algorithm and the ac-
estimated as 0.153. The low value of the error for
tual error, the derivative was calculated in a nar-
Xi = i may be explained by a reduction in round-off
row range of step sizes in the neighborhood of hopt.
errors when all the terms in the Vandermonde matrix
A sample of the results is given in table II. It was
are exact integers.
clear that near hopt, the value of the finite differ-
Based on an fA value of 0.153, equation (11) gave ence derivative exhibits a random oscillation, and the
hs = 2.2, and c(~) = 6.8 X 1011 from equation (9).
small error obtained at hopt is somewhat fortuitous.
Because c(~) was outside the interval [0.001,0.11, hs
was reduced, successively, by factors of 10 till a value The FD algorithm also was used for derivatives
of 2.2 X 10- 3 resulted in c( ~) = 0.06. The optimum of the remaining polynomial coefficients with respect
step size, hopt, was found from equation (7) to be to X3. Table III lists the second derivative step size
1.71 X 10- 4 • The derivative aaI/ax3 was computed hs, the optimum step size hopt, and the value of fA
by repeating the solution of equation (12) with X3 for all 12 derivatives. It is seen from table III that
replaced by X3 + hopt, while the remaining Xi values hopt varies from 1.642 X 10- 4 to 3.054 X 10- 4 •
as well as the Yi values were left unchanged from
their nominal values. The value of the derivative was Swept-Wing Problem
obtained as -4.36028 X 106. The exact value of the A finite element structural model of a swept wing
derivative was calculated analytically by Gaylen A. is shown in figure 1. The Engineering Analysis Lan-
Thurston of NASA Langley Research Center to be guage (EAL) finite element program (ref. 6) was used
-4.36003 X 106 , so that the error in the derivative to model the wing. The model contains 194 elements
at hopt is 250, or 0.0057 percent. The derivatives as consisting of rods, triangular membranes, and shear
a function of the step size are listed in table I for a panels connected at 88 nodes and constrained at the
wide range of h values bounding hopt. It is evident root. Additional details of the wing model are avail-
that the choice of hopt is indeed a good one. able in ref~rences 4 and 7.

4
 TABLE II. FINITE DIFFERENCE EVALUATION membrane elements 1 through 6 located at the root
OF 8a1l8x3 IN THE NEIGHBORHOOD OF
hopt = 1.71 X 10- 4 FOR POLYNOMIAL- of the wing model as shown in figure 1.
FIT PROBLEM The first step is to calculate an initial value for hs
from equation (11). This calculation requires values
of the function f and the independent variable x and
Derivative an estimate of the error bound fA for this problem.
Step size, Finite difference Derivative error, In this example, x is the nominal value of the design
h derivative error a percent variable (0.2 in.), f is the normal displacement of
0.4356 X 10 4 -4.35986 X 106 170 0.0039 node 44 (40.8 in.), and fA (obtained by the method
0.5445 -4.35980 230 0.0053 described in appendix A) is 1.0 X 10- 8 in. These
0.6806 -4.35911 920 0.0211 values result in an initial hs value of 3.775 X 10- 5 in.
0.8507 -4.36278 2750 0.0631 used in equation (5a) to approximate cI>, which was
1.063 -4.35696 3070 0.0704 then used to calculate c(cI» from equation (9). The
1.329 -4.35715 2880 0.0661 value for c( cI» was 0.06 and fell within the acceptable
1.662 -4.36103 1000 0.0229 range [0.001,0.1]. The value of the optimum step
1.714 (hopt) -4.36028 250 0.0057 size hopt, 9.4 X 10- 6 in., produced a 0.0059 percent
2.077 -4.36072 690 0.0158 error in the finite difference derivative relative to
2.596 -4.36110 1070 0.0245 the exact value obtained in reference 4. This step
3.245 -4.36372 3690 0.0846 is better than the one (2 X 10- 5 in.) obtained by
4.056 -4.36583 5800 0.1330 trial and error in reference 4 and was achieved with
5.071 -4.36493 4900 0.1124 far fewer function evaluations. The variation of the
6.338 -4.36316 3130 0.0718 error with step size is summarized in table IV and
indicates that indeed the value of 9.4 X 10- 6 in. is
aRelative to exact value (-4.36003 X 10 6). close to optimal. Unlike the case of the polynomial-
fit problem, for this problem equation (13) predicts a
minimum error which is very close to the actual error
TABLE III. VALUES OF fA, hs, AND hopt FOR
DERIVATIVES 8ai/8x3 FOR POLYNOMIAL- (0.004 in. compared with 0.003 in.).
. FIT PROBLEM
TABLE IV. VARIATION OF DERIVATIVE OF TIP
DISPLACEMENT WITH STEP SIZE FOR
SWEPT-WING PROBLEM
Coefficient fA hs hopt
0.1526 2.209 10 3 1.714X 10 4
al X

a2 0.4560 5.395 1.794 Error b in

a3 0.5643 0.600 1.642 Step size,a finite difference
a4 0.3886 4.982 2.004 h, in. derivative, percent
a5 0.1674 2.314 2.093 2 X 10- 3 0.9400
a6 4.777 X 10- 2 1.748 2.449 2 X 10- 4 0.0960
a7 9.273 X 10- 2 0.770 2.112 c2 X 10- 5 0.0098
as 1.231 X 10- 2 2.807 2.478 d9.4 X 10- 6 0.0059
ag 1.099 X 10- 3 5.931 2.884 2 X 10- 6 0.0200
alO 6.311 X 10- 5 2.010 2.633
all 2.102 X 10- 6 3.668 2.913 aNominal thickness = 0.2 in.
a12 3.086 X lO- S 4.444 3.054 bRelative to analytically calculated value of -50.907 in.
(ref. 4).
CBest value obtained by trial and error in reference 4.
dStep size calculated by FD algorithm.
The problem, as described in reference 4, con-
cerns derivatives of the wing displacements with re-
spect to element cross-sectional dimensions. In order
to achieve accurate results, in reference 4 a time- Extension of FD Algorithm to Vector
consuming trial-and-error process was used to deter- Derivatives
mine an appropriate step size. The derivative of the When the derivative of a vector with respect to
normal displacement at node 44 with respect to one a design parameter is required, the FD algorithm
of the design variables is chosen to investigate the is not always useful in its original form. As seen
application of the FD algorithm. The design vari- in the polynomial-fit example (table III), the FD
able selected controls the thickness of the triangular algorithm gives a different optimum step size for the

5
 I~.~-----------------------------600------------------------------~
r

---r-
50
t•
X

Figure 1. Geometry, element numbers, and node numbers for swept wing. Dimensions in inches .

.6
derivative of each component of the vector. From to the polynomial-fit problem discussed earlier pro-
the standpoint of efficiency, it is desirable in the case duced 12 optimum step sizes, one for each coefficient.
of a derivative of a vector to find a single step size (See table III.) To calculate all 12 finite difference
for all components which is the optimum compromise sensitivity derivatives with the same step size, the
among them. It is possible, of course, that no such extended FD algorithm was used. To apply the ex-
compromise exists if the components have widely tended FD algorithm, the second derivative step size
different condition and truncation errors. However, hs obtained for a! (2.209 X 10- 3 , see table III) was
in most practical cases it is expected that such a used for all 12 coefficients. All 12 c( <l?) values were
compromise will be possible, and an extension of the in the acceptable range of [0.001,0.11. The values
FD algorithm described below is appropriate. for hopt changed slightly from the values given in ta-
In seeking a compromise among the step size re- ble III and ranged from 1. 714 x 10-4 to 3.364 X 10-4 •
quirements of the different components, it is natu- Subdividing this range into 20 intervals resulted in
ral to seek to minimize some norm of the error vec- 21 vectors, each containing 12 error components. The
tor. Given an error vector E with components ei, IIEII! and IIElb norms (eqs. (14) and (15)) produced
i = 1, 2, ... , N, there are three norms which may be the same step size of 1.879 x 10- 4 , whereas the IIElioo
used: produced a step size of 1. 797 x 10-4. A check of the
N
IIEII! = L leil (14)
actual errors based on these step sizes showed that
all 12 errors are acceptably small.
i=!

~ v'N,~ ~ ~t.'1
Swept-Wing Problem
IIEII, (15) To demonstrate the application of the extended
FD algorithm to the swept wing, a set of 12 displace-
IIElioo = m~ leil (16) ment derivatives were chosen. These are the deriva-
t tives of the normal displacements at nodes 5 through
where erms denotes the root-mean-square of the 10 and 39 through 44. (See fig. 1.) The FD algo-
error. rithm was first used to find 12 optimum step sizes.
The extended FD algorithm proceeds through the For this purpose, hs = 3.775 X 10- 5 in. was used for
following steps: all 12 components, and one step size was extracted
1. Select an initial hs (e.g., by using eq. (11)), and for calculating all 12 derivatives. Unlike the previous
calculate <l? from equation (5a) for each component case, the three norm criteria produced different step
of the vector. sizes: 9.05 x 10-6 in. from IIEII!, 9.44 x 10- 6 in. from
2. Estimate the condition error in <l? (c( <l?) from IIEII2' and 9.6 X 10- 6 in. from IIElioo. Derivatives
eq. (9)) for each component, and vary hs to bring using these step sizes were compared with analytical
that error into the required interval of [0.001,0.1]. results from reference 4. All three gave sufficiently
Failure to find a common hs is an indication that it accurate results. For example, table V shows that a
is futile to use a single step size for all components. step size of 9.44 x 10- 6 in. produced excellent results
In that case, the components may be divided into compared with the analytical solution.
two or more groups and subsequent steps repeated
for each group separately. Effect of fA
3. Estimate the optimum step size hopt for each Because the error bound fA may be uncertain in
component from equation (7). The largest and some problems, the effect of fA on the step size cal-
smallest values of hopt are denoted hmax and hmilll culations and derivative results is of interest. To test
respectively. . this effect, the value of the nominal fA was varied by
4. Divide the interval [h milll hmaxl into m equal a factor of 10 for each of the single-derivative prob-
subintervals resulting in m + 1 step size values h! = lems. The results are shown in tables VI, VII, and
hmilll h2, h3, ... , hm+! = h max . Equation (4) is used VIII for the sine function, the polynomial-fit prob-
to estimate the error in each component for each hi lem, and the swept-wing problem, respectively. In
value and the associated values of <l? and fA. each case, the modified values of fA had a noticeable
5. Select as hopt that value of hi which minimizes but small effect on the results, and the errors were
the desired error norm. still acceptable.
Applications of Extended FD Algorithm Concluding Remarks
Polynomial-Fit Problem This paper deals with methods for obtaining near-
The application of the standard FD algorithm optimum step sizes for finite difference approxima-

7
 TABLE V. FINITE DIFFERENCE DERIVATIVES TABLE VIII. EFFECT OF €A ON DERIVATIVE IN
OF SWEPT-WING TIP DEFLECTION FOR SWEPT-WING PROBLEM
hopt = 9.44 X 10- 6 in.

Finite Derivative
Derivative difference error,a
Node Finite difference Analytical method Error, €A, in. hopt, in. derivative percent
number method of reference 4 percent Nominal 9.4 X 10- 6 -50.904 0.0059
5 -2.9372 -2.9374 0.00681 (0.10348 x 10- 7)
6 -3.7578 -3.7580 0.00532 0.1 x Nominal 2.525 x 10- 6 -50.914 0.0138
7 -4.5608 -4.561 0.00438 (0.10348 x 10- 8)
8 -5.5679 -5.5682 0.00538 10 x Nominal 2.929 X 10 5 -50.900 0.0138
9 -10.645 -10.646 0.00939 (0.10348 x 10- 6)
10 -11.597 -11.598 0.00862
39 -47.737 -47.739 0.00419
aRelative to exact value (-50.907).
40 -48.219 -48.222 0.00622
41 -49.773 -49.776 0.00603 nite difference (FD) algorithm, previously described
42 -50.112 -50.115 0.00599 in the literature, is reviewed and applied in a tuto-
43 -50.496 -50.499 0.00594 rial manner to the derivative of a sine function. The
44 -50.904 -50.907 0.00589 original FD algorithm is applicable to one derivative
at a time. This paper describes an extension of the
FD algorithm to the calculation of several derivatives
TABLE VI. EFFECT OF €A ON DERIVATIVE OF sin x of matrix equations.
The FD algorithm requires an estimate of the er-
ror in the computed value of the functions being dif-
Finite Derivative ferentiated. This paper describes a simple method
difference error,a for calculating the error estimate for matrix prob-
€A hopt derivative percent lems. Denoted the iterative improvement technique,
Nominal 1.108 x 10- 3 0.706679 0.06 the method involves a second solution of the matrix
(2.188 x 10- 7) equation. Since it uses the same decomposed matrix
0.1 x Nominal 1.20 x 10 4 0.708333 0.17 that was used in the original solution, it imposes a
(2.188 x 10- 8) minimal computational burden on the analysis.
10 x Nominal 3.464 x 10 3 0.705831 0.18 Both the original and extended FD algorithms
(2.188 x 10- 6 ) were applied to sensitivity analysis for a data-fitting
problem in which a polynomial is passed through
aRelative to exact value (0.707107). data points, and the derivatives of the coefficients
of the polynomial with respect to uncertainties in
the data were calculated. The algorithms also were
TABLE VII. EFFECT OF €A ON DERIVATIVE
8at/8x3 IN POLYNOMIAL-FIT PROBLEM applied to sensitivity analysis of the structural re-
sponse of a finite-element-modeled swept wing. In a
previous paper, this sensitivity analysis of the swept
Finite Derivative wing required a time-consuming trial-and-error effort
difference error,a to obtain a suitable step size, but it proved to be a
€A hopt derivative percent routine application for the FD algorithm herein.
Nominal 1.714 X 10- 4 -4.3603 X 106 0.006 One limitation of the FD algorithm is that it re-
(0.153) quires a finite lower bound on the second derivative
0.1 x Nominal 2.656 x 10- 6 -4.3596 X 106 0.009 of the function being differentiated. This bound is
(0.0153) generally estimated by a central difference approxi-
10 x Nominal 4.885 x 10- 5 4.3636 X 106 0.08 mation. In certain 'special cases, such as when the
(1.53) function is nearly constant, linear, or an odd func-
tion, the second derivative is zero or close to zero,
a Relative to exact value (-4.3600 x 106). and in such a case, the FD algorithm is not a reli-
able provider of an optimum step size for the first
derivative. In the absence of this special condition,
tions to first derivatives with particular application however, the FD algorithm has been found to be use-
to sensitivity analysis. A technique denoted the fi- ful, accurate, and reliable.

8
Appendix A and F is a force vector. If U is the solution obtained
numerically, F can be defined as
Iterative Improvement Technique for
Estimating € A (A2)
A method for estimating fA is based on the so-
called "iterative improvement" of a solution to a Subtracting equation (A2) from equation (AI) gives
system of linear equations. This process is explained
for the system K(U-U) =F-F (A3)

KU=F (AI) so that an estimate of U - U may be obtained by

solving equation (A3). The solution gives an order
where K is an n X n matrix, and F and U are of magnitude estimate for the error fA, and since it
n x I vectors. In finite element structural analysis, uses the same decomposed matrix that was used in
K is the stiffness matrix, U is a displacement vector, the original solution (eq. (AI)), it imposes a minimal
computational burden on the analysis.

9
 Appendix B Since c( cI» is outside the acceptable interval
[O.OOl,O.lJ, hs must be increased, and new values for
Application of FD Algorithm to Sine cI> and c( cI» must be calculated. For convenience, hs
Function . was increased by a factor of 10. The new value for
1cI>1 is 0.711656 for c( cI» = 0.008621. Since this value
A simple example to demonstrate the FD algo- falls within the acceptable interval, cI> may be used
rithm is sin x. The derivative of this function was to calculate hopt from equation (7).
evaluated with a hand calculator which carries up to
10 significant digits. The function sin x is convenient
because the exact derivatives are known and can be
compared with the approximations. The bound on
hopt =2/M
the absolute error in the computed function values fA
is controlled by the number of digits used on the cal-
culator. In this example, the number of digits used
=2 2.188 X 10-7
0.711656
was six. The function sin x and its first derivative = 1.108 X 10- 3
were evaluated at x = 0.785398 rad, or 45°. There-
fore, f(x) = sin x = 0.7071067812 to 10 significant Substituting hopt into the forward difference for-
figures, and sin x = 0.707107 to 6 significant fig- mula (eq. (1)) gives an approximation to the first
ures. The difference between these two numbers is derivative
€A = 2.188 X 10-7.
The first step is the calculation of the initial step f'(x) ~ f(x + hopt) - f(x)
size hs to be used in estimating the second derivative. hopt
From equation (11)
0.707890 - 0.707107
1.108 X 10- 3
hs = 2(1 + Ix!) Vr2r.
1 + If I = 0.706675

2.188 X 10-7 Comparing this result with the exact derivative in-
= 2(1 + 0.785398) 1 + 0.707107 dicates that the use of hopt gives a relative error of
0.06 percent. Table BI shows derivatives calculated
= 0.001278 with step sizes larger and smaller than the optimum
along with the error relative to the exact derivative.
Using hs in equation (5a) gives

f(x + hs ) - 2 f(x) + f(x - hs ) TABLE BI. FINITE DIFFERENCE DERIVATIVES OF

cI>= 2 sin x AT 45° FOR VARIOUS STEP SIZES
hs
0.708010 - 2(0.707107) + 0.706202
=
0.001278 2 Step size, Finite difference Derivative
h derivative error,a percent
= -1.225
1.108 X 10- 1 0.666525 5.74
1.108 X 10- 2 0.703159 0.56
The relative condition error is calculated from equa-
1.108 X 10- 3 (hopt) 0.706675 0.06
tion (9)
1.108 X 10- 4 0.703971 0.40
~ 4fA
1.108 X 10- 5
c(cI» = h; 1cI>1 0.722022 2.11

4(2.188 X 10-7) aRelative to exact value (0.707107).

= (0.001278)2(1.225)
= 0.4374

10
References 4. Camarda, Charles J.; and Adelman, Howard M.: Static
and Dynamic Structural-Sensitivity Derivative Calcula-
1. Gill, Philip E.; Murray, Walter; Saunders, Michael A.; tions in the Finite-Element-Based Engineering Analysis
and Wright, Margaret H.: Computing the Finite- Language (EAL) System. NASA TM-85743, 1984.
Difference Approximations to Derivatives for Numeri- 5. Lyness, J. N.: Has Numerical Differentiation a Future?
cal Optimization. SOL 80-6 (Contract DAAG29-79-C- Proceedings of the Seventh Manitoba Conference on Nu-
0110), Dep. Operations Research-SOL, Stanford Univ., merical Mathematics and Computing, Univ. of Manitoba,
May 1980. Winnipeg, 1977, pp. 107-129.
2. Gill, Philip E.; Murray, Walter; Saunders, Michael A.; 6. Whetstone, W. D.: EISI-EAL Engineering Analysis Lan-
and Wright, Margaret H.: Computing Forward- guage Reference Manual-EISI-EAL System Level 2091.
Difference Intervals for Numerical Optimization. SIAM Engineering Information Systems, Inc., July 1983.
J. Sci. fj Stat. Comput., vol. 4, no. 2, June 1983, 7. Schmit, Lucien A., Jr.; and Miura, Hirokazu: Ap-
pp. 310-321. proximation Concepts for Efficient Structural Synthesis.
3. Franklin, Joel N.: Matrix Theory. Prentice-Hall, Inc., NASA CR-2552, 1976.
c.1968.

11
 Standard Bibliographic Page
1. Report No. 12. Government Accession No. 3. Recipient's Catalog No.
NASA TM-86382
4. Title and Subtitle 5. Report Date
Selecting Step Sizes in Sensitivity Analysis by Finite Differences August 1985
6. Performing Organization Code

7. Author(s}
506-53-53-07
Jocelyn Iott, Raphael T. Haftka, and Howard M. Adelman 8. Performing Organization Report No.
L-15938
9. Performing Organization Name and Address 10. Work Unit No.
NASA Langley Research Center
Hampton, VA 23665 11. Contract or Grant No.

13. Type of Report and Period Covered

12. Sponsoring Agency Name and Address
National Aeronautics and Space Administration Technical Memorandum
Washington, DC 20546 14. Sponsoring Agency Code

15. Supplementary Notes

Jocelyn Iott and Raphael T. Haftka: Virginia Polytechnic Institute and State University, Blacksburg,
~ Virginia.
Howard M. Adelman: Langley Research Center, Hampton, Virginia.
16. Abstract
This paper deals with methods for obtaining near-optimum step sizes for finite difference approximations to
first derivatives with particular application to sensitivity analysis. A technique denoted the finite difference
(FD) algorithm, previously described in the literature and applicable to one derivative at a time, is extended
to the calculation of several derivatives simultaneously. Both the original and extended FD algorithms
are applied to sensitivity analysis for a data-fitting problem in which derivatives of the coefficients of an
interpolation polynomial are calculated with respect to uncertainties in the data. The methods are also
applied to sensitivity analysis of the structural response of a finite-element-modeled swept wing. In a
previous study, this sensitivity analysis of the swept wing required a time-consuming trial-and-error effort
to obtain a suitable step size, but it proved to be a routine application for the extended FD algorithm
herein.

17. Key Words (Suggested by Authors(s» 18. Distribution Statement

Sensitivity analysis Unclassified-Unlimited
Finite element analysis
Finite difference method
Optimization methods
Structural analysis
Subject Category 39
19. Security Classif.(of this report} 120. Security Classif.(of this page} 21. No. of Pages 122. Price
Unclassified Unclassified 12 A02

For sale by the National Technical Information Service, Springfield, Virginia 22161
NASA-Langley, 1985
End of Document