Numerical Methods in Software and Analysis

717 Pages·1993·35.731 MB·English

by  John R. Rice (Auth.)

N u m e r i c al M e t h o d s,
S o f t w a r e,
a nd A n a l y s is
Second Edition
John R. Rice
Department of Computer Science
Purdue University
West Lafayette, Indiana
ACADEMIC PRESS, INC.
Harcourt Brace Jovanovich, Publishers
Boston San Diego New York
London Sydney Tokyo Toronto
This book is printed on acid-free paper. ©
Copyright © 1993, 1983 by Academic Press, Inc.
All rights reserved.
No part of this publication may be reproduced or
transmitted in any form or by any means, electronic
or mechanical, including photocopy, recording, or
any information storage and retrieval system, without
permission in writing from the publisher.
IMSL and PROTRAN are the trademarks of IMSL Inc. for its proprietary computer
software. No materials describing such software may be produced or distributed
without the written permission of the owners of the trademark and license rights
in the software and the copyright in the published materials.
ACADEMIC PRESS, INC.
1250 Sixth Avenue, San Diego, CA 92101-4311
United Kingdom edition published by
ACADEMIC PRESS LIMITED
24-28 Oval Road, London NW1 7DX
ISBN 0-12-587755-2
Printed in the United States of America
92 93 94 95 EB 9 8 7 6 5 4 3 2 1
PREFACE
Objectives
Modem science depends on our capability to carry out substantial computations to explore the
consequence of the laws of nature. These laws are expressed as mathematical models, and numeri-
cal methods are used throughout the computation. The objective of this text is to introduce science
and engineering students to the methods, tools, and ideas of numerical computation. These stu-
dents are assumed to have taken the standard introductory calculus sequence, plus an elementary
introduction to matrix theory, ordinary differential equations, and Fortran programming.
An introductory course in numerical methods faces a fundamental problem: students who want
to do realistic computations need to know more than they can learn in one course. This text has
about twice as much material as can be covered in a typical one-term course, and yet most topics
are not covered in depth. The use of high-quality mathematical software is the leading way out
of the dilemma of too little time to learn too much knowledge and know-how. The mathematical
software approach is based on the premise that experts in numerical computation can incorporate
enough of their experience and know-how into a computer program so that it performs better than
the program a typical scientist can create with a reasonable effort. The validity of this approach has
been established in many areas of numerical computation; the objective of this text may be restated
as to present scientific problem solving using standard mathematical software.
Organization and Content
This book is organized into three parts, as shown in the diagram below. Part I is background and
preliminary material. Few instructors will cover this thoroughly; in a one-term course it is feasible
to summarize the main ideas in a few lectures. Part II is the core of the course, presenting the
principal methods and ideas of numerical computation. Part III contains very relevant, but somewhat
independent, material about software engineering and performance evaluation. There is rarely time
in a one-term course to cover this material in detail, yet it is important that young scientists gain
some appreciation for this material.
xi
xii PREFACE
PART I: Preliminaries
1. MATH AND CS BACKGROUND
2. NUMERICAL SOFTWARE
3. ERRORS
PART II: Core Material PART III: Software Engineering
12. SOFTWARE PRACTICE
4. MODELS AND FORMULAS
13. PERFORMANCE EVALUATION
5. INTERPOLATION
14. SOFTWARE VALIDATION
6. MATRIX COMPUTATION
15. PROTRAN
7. DERIVATIVES/INTEGRALS
8. NONLINEAR EQUATIONS
9. DIFFERENTIAL EQUATIONS
10. PARTIAL DIFFERENTIAL EQUATIONS
11. APPROXIMATION
A fairly uniform approach is used in each area of numerical computation. First, an intuitive
development is made of the problems and the basic methods for their solution. Then, relevant math-
ematical software is reviewed and its use outlined. In many areas rather extensive examples and case
studies are given. Finally, a deeper analysis of the methods is presented as in traditional numerical
analysis texts.
The use of mathematical software tools may resolve the dilemma facing the practicing scientist,
but it still leaves a dilemma for the instructor using this text: should one emphasize how to use the
software or understanding how it works? This is not an easy choice to make; ample material is pre-
sented so that an instructor may choose whatever balance seems appropriate for the students' needs
and the instructor's taste. Some sections and chapters are "starred" to help organize the material. The
double-starred items contain material that is either more advanced (for example, osculatory inter-
polation and adaptive quadrature) or somewhat peripheral to the main topics (for example, sources
of matrix problems and reporting computational experiments). Since the remaining material is still
more than can be covered in a one-term course, the single-starred items are suggested as candidates
to be skipped. Some single-starred chapters also have single-starred sections; this provides suggested
topics for an instructor who wants to cover these chapters briefly. The dependence of the chapters is
shown in the following diagram.
PREFACE xiii
Software
This text discusses numerous programs and software packages, but it is not intended to be a user's
guide for this software. The IMSL library and ACM Algorithms are discussed systematically as
they are the most widely available general collections of high-quality mathematical software. It is
not essential to have this particular software; any reasonably complete program library will have
programs similar to those used in most examples and problems. We note (see Chapter 12) that the
cost of developing just one good program is more than the cost of acquiring several hundred good
programs from IMSL, ACM, and similar sources. These sources should be viewed as providing the
initial step toward a software parts technology which might dramatically reduce the present heavy
burden of software development costs.
This text also uses the PROTRAN system introduced by IMSL. PROTRAN provides mathe-
matical problem-solving statements which substantially ease the use of library routines and simulta-
neously improve programming reliability. PROTRAN has two important additional attractions for
use in this text. First, example programs in PROTRAN are much shorter than in Fortran, and even
if PROTRAN is not available to the students, they can easily visualize how an equivalent Fortran
program would be. Second, PROTRAN is an example of the very high level languages that most
students will be using in a few years as this methodology spreads throughout scientific computation.
Acknowledgments
The most important acknowledgment is to the students at Purdue University who use preliminary
versions of this text. Their suggestions and patience with typos and rough spots are greatly appre-
ciated. Second, I thank Thomas Aird especially, for carefully reviewing the discussion of IMSL
software and for making the PROTRAN system available. I thank Robert Lynch for his invaluable
xiv PREFACE
assistance in developing the library applications presented here; I also thank Dinesh Manocha and
Thomas Stiller for their help in implementing these applications. I thank Carl de Boor for permis-
sion to include the text of two programs (INTERV and Example 5.6) from PPPACK and for many
valuable insights into the nature of scientific computation. I thank John Reid for material used in the
cover design and Figure 6.3. Many valuable comments were received from Neil Berman, Ronald
Boisvert, Henry Darilek, Wayne Dyksen, Don Kainer, Iftekhar Karimi, Robert Lynch, and Granville
Sewell. I thank them and the others who have contributed to the development of this text. I thank
Doris Moore for typing the text and several revisions, Richard Gaudette for preparing the art work,
and Gary Lelvis for coordinating the preparation, typing, and typesetting of the text.
John R. Rice
June, 1992
**1
MATHEMATICS AND COMPUTER
SCIENCE BACKGROUND
The background assumed in this text is a basic sequence in calculus, an introduction to differential
equations and linear algebra (or matrices), and a working knowledge of Fortran programming. The
purpose of this chapter is to give a highly condensed review of this background. The material pro-
vides references for the rest of the text, and the problems allow instructors or students to check their
knowledge against the assumed background in calculus and linear algebra. Many instructors will
choose to spend little, if any, time presenting this material in lectures.
Some material in this book requires more background and experience than presented here, for
example, the chapter on partial differential equations. An instructor must judge whether the students
have the necessary experience and maturity for such material; some instructors might choose to use
the material here as a first introduction to the topic. The sections marked by ** contain material
which is more advanced or peripheral and which can be skipped; all the material requiring more
background is in such sections.
1.1 CALCULUS
The principal topics in calculus are the real and complex number systems, the concept of limits and
convergence, and the properties of functions. Convergence of a sequence of numbers x, is defined
as follows: The sequence X{ converges to the limit x* if, given any tolerance c > 0, there is an index
Ν = N(c) so that for all i> Ν we have |a?- - χ* | < e. The notation for this is
t
lim X{ = x*
Ϊ—+ΌΟ
Convergence is also a principal topic of numerical computation, but with a different emphasis. In
calculus one studies limits and convergence with analytic tools; one tries to obtain the limit or to show
that convergence takes place. In computations, one has the same problem but little or no theoretical
1
2 MATHEMATICS AND COMPUTER SCIENCE BACKGROUND
knowledge about the sequence. One is frequently reduced to using empirical or intuitive tests for
convergence; often the principal task is to actually estimate the value of the tolerance c for a given
The study of functions in calculus revolves about continuity, derivatives, and integrals. A func-
tion f(x) is continuous if
lim./(*,•) = /(**)
holds for all x* and all ways for the Xi to converge to x*. We list six theorems from calculus which
are useful for estimating values that appear in numerical computation.
Theorem 1.1 (Mean value theorem for continuous functions) Let f(x) be continuous on the in-
terval [a, b]. Consider points XHI and XLOW in [a, b] and a value y so that f(XLOW) < y <
f(XHI). Then there is a point in [a, 6] so that
Theorem 1.2 (Mean value theorem for sums) Let f(x) be continuous on the interval [a, b], let
x\,X2^..,xbe points in [a, 6], and let w\,W2,.. - ,w be positive numbers. Then there is a point
n n
£ in [a, b] so that
η η
£>*/(**) =/(Ο 5>i
Theorem 1.3 (Mean value theorem for integrals) Let f(x) be continuous on the interval [a, b] and
let w(x) be a nonnegative function [w(x) > 0] on [a, 6]. Then there is a point ξ in [a, 6] so that
/ w(x)f(x)dx = /(£) / w(x)dx
Ja J a
Theorems 1.2 and 1.3 show the analogy that exists between sums and integrals. This fact derives
from the definition of the integral as
/ f(x)dx= lim y2f(xi)(xM ~Xi)
Ja max\xi-xi\-+0
i+
where the points a?», with X( < x \, are a partition of [a, 6]. This analogy shows up for many nu-
i+
merical methods where one variation applies to sums and another applies to integrals. Theorem 1.2
is proved from Theorem 1.1, and then Theorem 1.3 is proved by a similar method. The assumption
that w(x) > 0(wi > 0) may be replaced by w(x) < 0(w{ < 0) in these theorems; it is essential that
w(x) be of one sign as shown by the example w(x) = f(x) = χ and [a, b] = [-1,1].
Theorem 1.4 (Continuous functions assume max/min values) Let f(x) be continuous on [a, b]
with |a|, |6| < oo. Then there are points XHI and XLOW in [a, b] so that for all χ in [a, 6]
f(XHI) < f(x) < f(XLOW)
The derivative of f(x) is defined by
£ , / ( «+ *>-/(«>
=/ ( x ) = U m
dx Λ -Ο h
CALCULUS 3
As an illustration of the difference between theory and practice, the quantity [/(x + ft) - f(x)]/h
can be replaced by [f(x + ft) - f(x - h)]/(2h) with no change in the theory but with a dramatic
improvement in the rate of convergence; that is, much more accurate estimates of /'(x) are obtained
for a given value of ft. The kth derivative is the derivative of the (k - l)th derivative; they are de-
k k (4) (5)
noted by df/dx or /"(*), /'"(x), / (x), / (x),...
Theorem 1.5 (Mean value theorem for derivatives) Let f(x) be continuous and differentiate in
[a, 6] with \a\, |6| < oo. Then there is a point ζ in [a, b] so that
ο — a
Equivalently, given c and χ in [a, 6] there is a point ξ between χ and c so that
f(x) = f(c) + f'(0(x-c)
The special case of Theorem 1.5 with /(a) = f(b) = 0 is known as Rolle's theorem. It states that
a
if f() = fQ>) = 0> then there is a point ξ between a and b so that /'(£) = 0. This is derived from
Theorem 1.5 by multiplying through by 6 - a, renaming a, b as x, c, and then applying the first form
to the smaller interval [x, c] or [c, x], depending on the relation between χ and c.
A very important tool in numerical analysis is the extension of the second part of Theorem 1.5
to use higher derivatives.
Theorem 1.6 (Taylor's series with remainder) Let f(x) have η +1 continuous derivatives in [a, 6].
Given points χ and c in [a, b] we have
Kx) = /(c) + /'(cX* - c) + / " ( O ^^ + / ' " ( c ) ^^ +... + / < " ) ( ) ^^ +
C ß n + l) 0 r
ζ! 3! ni
where R \ has either one of the following forms (ξ is a point between χ and c):
n+
(n + 1)!
R i(x) = -^7 [\x-trfn+\lt)dt
n+
If a function / depends on several variables, one can differentiate it with respect to one variable,
say x, while keeping all the rest fixed. This is the partial derivative of / and is denoted by 6f/6x
or f. Higher-order and mixed derivatives are defined by successive differentiation. Taylor's series
x
for functions of several variables is a direct extension of the formula in Theorem 1.6, although the
number of terms in it grows rapidly. For two variables it is
2
/(*, V) = /(c, d) + f(x -c) + f(y-d)+ ^[f(x - cf + 2f(x - c)(y - d) + f (y - d)] + ...
x y xx xy yy
where all the partial derivatives are evaluated at the point (c, d).
4 MATHEMATICS AND COMPUTER SCIENCE BACKGROUND
Theorem 1.7 (Chain rule for derivatives) Let f(x, y,..., z) have continuous first partial derivatives
with respect to all its variables. Let χ = x(t), y = y(t),..., ζ = z(t) be continuously differentiate
functions oft. Then g{t) = /(x(t), 2/(0, · · ·, z{t)) is continuously differentiable and
g'(t) = fxx'(t) + fy'(t) + ... + fz'(t)
y 2
Finally we state
Theorem 1.8 (Fundamental theorem of algebra) Letpix) be a polynomial of degree η > 1, that is,
2 n
p(x) = ao + a\x + a,2X +... + a x
n
where the ai are real or complex numbers and a ^ 0. Then there is a complex number ξ so that
n
p(O = 0.
1.1 PROBLEMS
2 2
1.1.1 Determine the following limits: (C) fix) - x = \
(D) x = cosß, fix) = sin θ
_10_i
(A) lim^i Zp^f
1.L5 Give the first four terms of the Taylor's series at c
(B) lim ^o χ—s^in χ for the following functions:
x
(C) \imx^/2 ^z^f (A) fix) = i\+x)k,c = 0
(D) lim _
n n+y/n (B) fix) = log 27, c = 1
(E) lim^o ψ- (C) /Or) = e*,c = 0
(F) lim _oo[log(n + 1) - log(n)] (D) f(x)= l/x,c= -1
n
1.1.2 Find the mean value ξ of Theorem 1.3 for the fol- (E) fix) = sin χ, c = π/6
lowing integrals:
(F) fix) = \ / l - x\ c = -2
2 2 1.1.6 Obtain the Taylor's series for sin a; and cos χ with
(A) jlxdx, with/Or) = s
c = 0. Square each of them, add the results, and verify that
2 2
2 sin x + cos χ = 1.
(B) Jq xdx, with/Or) = s
1.1.7 Find f and f for each of the following:
x y
2 2
(C) J y/x~ix + \)dx with f(x) = Jx~ 2 2
0 (A) fix,y)=x -2xy+4y
2 2 2 (B) fix, y) = ix + y)/ix - y)
(D) J^(x + \)dx with f(x) = x + 1
2
(C) fix,y)=x\0£iy/x)
1.1.3 Find the mean value ξ of Theorem 1.5 for the fol-
(D) fix, y) = ln(y + z) where ^(x, 3/) = ye*'*
lowing:
1.1.8 Find ci//cft for each of the following:
2
(A) fix) = x - Ax + 3, a = 1, b = 4
2 3
(A) fix, y) = x + xy, x = 3t, j/ = *
(B)/(ar) = a?i, α = 0, 6 = 9
2 2 2
(B) fix,y)= y/x + y, x = t, y = e* - e~*
iQ fix) = x + l/x, a = 1,6 = 9
1.1.9 Give the first six terms of the Taylor's series for the
1.1.4 Find /' and /" for each of the following: following:
3
(A) fix) = y/2x~+3 iA)fix,y) = ix + y), = 0,d = 0
c
a;+1
(B) /(x) = sin2a: (B)/(x, /) = e /y, c = 0,d= 1
2
2 2
(C)/(^2/)=log0r + 2/), c=l,d = l