451
Numerical CMPSC/MATH Computations
Lecture 36
Nov 14, 2011
Prof. Kamesh Madduri
Class Overview
Numerical Differentiation
Finite difference schemes
Richardson Extrapolation
Slides from textbook follow.
2
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Numerical Differentiation
Finite Difference Approximations
Automatic Differentiation
Numerical Differentiation
Differentiation is inherently sensitive, as small
perturbations in data can cause large changes in result
Differentiation is inverse of integration, which is inherently
stable because of its smoothing effect
For example, two functions shown below have very similar
denite integrals but very different derivatives
Michael T. Heath
Scientic Computing
47 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Numerical Differentiation
Finite Difference Approximations
Automatic Differentiation
Numerical Differentiation, continued
To approximate derivative of function whose values are
known only at discrete set of points, good approach is to t
some smooth function to given data and then differentiate
approximating function
If given data are sufciently smooth, then interpolation may
be appropriate, but if data are noisy, then smoothing
approximating function, such as least squares spline, is
more appropriate
< interactive example >
Michael T. Heath
Scientic Computing
48 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Numerical Differentiation
Finite Difference Approximations
Automatic Differentiation
Finite Difference Approximations
Given smooth function f : R R, we wish to approximate
its rst and second derivatives at point x
Consider Taylor series expansions
f (x) 2 f (x) 3
h+
h +
2
6
f (x) 2 f (x) 3
h
h +
f (x h) = f (x) f (x)h +
2
6
Solving for f (x) in rst series, obtain forward difference
approximation
f (x + h) = f (x) + f (x)h +
f (x + h) f (x) f (x)
f (x + h) f (x)
h +
h
2
h
which is rst-order accurate since dominant term in
remainder of series is O(h)
f (x) =
Michael T. Heath
Scientic Computing
49 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Numerical Differentiation
Finite Difference Approximations
Automatic Differentiation
Finite Difference Approximations, continued
Similarly, from second series derive backward difference
approximation
f (x) f (x h) f (x)
+
h +
h
2
f (x) f (x h)
h
which is also rst-order accurate
f (x) =
Subtracting second series from rst series gives centered
difference approximation
f (x + h) f (x h) f (x) 2
h +
2h
6
f (x + h) f (x h)
2h
which is second-order accurate
f (x) =
Michael T. Heath
Scientic Computing
50 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Numerical Differentiation
Finite Difference Approximations
Automatic Differentiation
Finite Difference Approximations, continued
Adding both series together gives centered difference
approximation for second derivative
f (x + h) 2f (x) + f (x h) f (4) (x) 2
h +
h2
12
f (x + h) 2f (x) + f (x h)
h2
which is also second-order accurate
f (x) =
Finite difference approximations can also be derived by
polynomial interpolation, which is less cumbersome than
Taylor series for higher-order accuracy or higher-order
derivatives, and is more easily generalized to unequally
spaced points
< interactive example >
Michael T. Heath
Scientic Computing
51 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Numerical Differentiation
Finite Difference Approximations
Automatic Differentiation
Automatic Differentiation
Computer program expressing function is composed of
basic arithmetic operations and elementary functions, each
of whose derivatives is easily computed
Derivatives can be propagated through program by
repeated use of chain rule, computing derivative of function
step by step along with function itself
Result is true derivative of original function, subject only to
rounding error but suffering no discretization error
Software packages are available implementing this
automatic differentiation (AD) approach
Michael T. Heath
Scientic Computing
52 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Richardson Extrapolation
Romberg Integration
Richardson Extrapolation
In many problems, such as numerical integration or
differentiation, approximate value for some quantity is
computed based on some step size
Ideally, we would like to obtain value limiting as step size
approaches zero, but we cannot take step size arbitrarily
small because of excessive cost or rounding error
Based on values for nonzero step sizes, however, we may
be able to estimate value for step size of zero
One way to do this is called Richardson extrapolation
Michael T. Heath
Scientic Computing
53 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Richardson Extrapolation
Romberg Integration
Richardson Extrapolation, continued
Let F (h) denote value obtained with step size h
If we compute value of F for some nonzero step sizes, and
if we know theoretical behavior of F (h) as h 0, then we
can extrapolate from known values to obtain approximate
value for F (0)
Suppose that
F (h) = a0 + a1 hp + O(hr )
as h 0 for some p and r, with r > p
Assume we know values of p and r, but not a0 or a1
(indeed, F (0) = a0 is what we seek)
Michael T. Heath
Scientic Computing
54 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Richardson Extrapolation
Romberg Integration
Richardson Extrapolation, continued
Suppose we have computed F for two step sizes, say h
and h/q for some positive integer q
Then we have
F (h) = a0 + a1 hp + O(hr )
F (h/q ) = a0 + a1 (h/q )p + O(hr ) = a0 + a1 q p hp + O(hr )
This system of two linear equations in two unknowns a0
and a1 is easily solved to obtain
a0 = F (h) +
F (h) F (h/q )
+ O(hr )
q p 1
Accuracy of improved value, a0 , is O(hr )
Michael T. Heath
Scientic Computing
55 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Richardson Extrapolation
Romberg Integration
Richardson Extrapolation, continued
Extrapolated value, though improved, is still only
approximate, not exact, and its accuracy is still limited by
step size and arithmetic precision used
If F (h) is known for several values of h, then extrapolation
process can be repeated to produce still more accurate
approximations, up to limitations imposed by
nite-precision arithmetic
Michael T. Heath
Scientic Computing
56 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Richardson Extrapolation
Romberg Integration
Example: Richardson Extrapolation
Use Richardson extrapolation to improve accuracy of nite
difference approximation to derivative of function sin(x) at
x=1
Using rst-order accurate forward difference
approximation, we have
F (h) = a0 + a1 h + O(h2 )
so p = 1 and r = 2 in this instance
Using step sizes of h = 0.5 and h/2 = 0.25 (i.e., q = 2), we
obtain
sin(1.5) sin(1)
F (h) =
= 0.312048
0.5
sin(1.25) sin(1)
F (h/2) =
= 0.430055
0.25
Michael T. Heath
Scientic Computing
57 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Richardson Extrapolation
Romberg Integration
Example, continued
Extrapolated value is then given by
F (0) = a0 = F (h)+
F (h) F (h/2)
= 2F (h/2)F (h) = 0.548061
(1/2) 1
For comparison, correctly rounded result is
cos(1) = 0.540302
< interactive example >
Michael T. Heath
Scientic Computing
58 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Richardson Extrapolation
Romberg Integration
Example: Romberg Integration
As another example, evaluate
/2
sin(x) dx
0
Using composite trapezoid rule, we have
F (h) = a0 + a1 h2 + O(h4 )
so p = 2 and r = 4 in this instance
With h = /2, F (h) = F (/2) = 0.785398
With q = 2, F (h/2) = F (/4) = 0.948059
Michael T. Heath
Scientic Computing
59 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Richardson Extrapolation
Romberg Integration
Example, continued
Extrapolated value is then given by
F (h) F (h/2)
4F (h/2) F (h)
=
= 1.002280
22 1
3
which is substantially more accurate than values previously
computed (exact answer is 1)
F (0) = a0 = F (h)+
< interactive example >
Michael T. Heath
Scientic Computing
60 / 61
Numerical Integration
Numerical Differentiation
Richardson Extrapolation
Richardson Extrapolation
Romberg Integration
Romberg Integration
Continued Richardson extrapolations using composite
trapezoid rule with successively halved step sizes is called
Romberg integration
It is capable of producing very high accuracy (up to limit
imposed by arithmetic precision) for very smooth
integrands
It is often implemented in automatic (though nonadaptive)
fashion, with extrapolations continuing until change in
successive values falls below specied error tolerance
Michael T. Heath
Scientic Computing
61 / 61