Atkinson K. An Introduction to Numerical Analysis
Подождите немного. Документ загружается.
228 APPROXIMATION OF FUNCTIONS
Thus
{4.7.11)
which has
2n
+ 4 relative maxim3;
and
minima,
of
alternating sign.
By
Theorem
4.9
and
the preceding inequalities,
2(.793)a
2
"+
3
2{1.207)a
2
"+
3
2
n +
3
~
P1n+l{f) =
P2n+2(/)
~
2
n +
3
{4.7.12)
For
the practical evaluation
of
the coefficients ci in (4.7.1) use the formula
21"'
ci
= - cos
(jO)f(cos
0)
dO
7T
0
(4.7.13)
and
the
midp0int
or the trapezoidal numerical integration rule. This was il-
lustrated earlier in (4.5.21).
Interpolation
at
the Chebyshev zeros
If
the error in the minimax approxima-
tions is nearly
c,+
1
T,+
1
(x),
as
derived in (4.7.5), then the error should be nearly
zero
at
the roots of T,+
1
(x)
on
[-1,1].
These are the points
[
2}
+ 1 ]
xj
= cos
TT
2n + 2
j=0,1,
...
,n
(4.7.14)
Let
I,(x)
be
the
polynomial of
degree~
n that interpolates
f(x)
at
these nodes
{xi}.
Since it has zero error at the Chebyshev zeros
xi
of
T,+
1
(x),
the continuity
of
an
interpolation polynomial with respect to the function values defining it
suggests
that
I,(x)
should approximately equal Cn(x), and hence also
q:(x).
More
precisely, write
I,(x)
in its Lagrange form and manipulate it as follows:
n
I,(x)
= L
f(x)i/x)
j-0
n n
= L
C,(x)l/x)
+ L
[!(xi)-
Cn(xj)}li(x)
j-0 j-0
= C,(x)
I,(x)
= C,(x)
(4.7.15)
because
f(x)-
C,(x
1
)
= c,+
1
T,+
1
(xi)
=
0.
The
term
I,(x)
can be calculated using the algorithms Divdif
and
Interp
of
Section 3.2. A note of caution:
If
c,+
1
= 0, then the error
f(x)-
q!(x)
is
likely
to have as approximate zeros those
of
a higher degree Chebyshev polynomial,
NEAR-MINIMAX APPROXIMATIONS 229
usually Tn+
2
(x). This is likely
to
happen if
f(x)
is either odd or even
on
[
-1,
1]
(see Problem
29).
[A
function
f(x)
is
even {odd} on an interval
[-a,
a]
if
/(
-x)
=
f(x)
{/(
-x)
=
-/(x)}
for all x in
[-a,
a].] The preceding example
with
f(x)
=
tan-
1
xis
an apt illustration.
A further motivation for considering
/n(x)
as
a near-minimax approximation
is
based on the following important theorem about Chebyshev polynomials.
Theorem 4.12 For a
fixed
integer n >
0,
consider the minimization problem:
Tn
= Infimum [ Max lxn +
Q(x)
1]
deg(Q):Sn-1
-1
:Sx
:s;
1
(4.7.16)
with
Q(x)
a polynomial. The minimum
Tn
is
attained uniquely
by
letting
(4.7.17)
defining
Q(x)
implicitly. The minimum
is
1
T=--
n
2n-!
(4.7.18)
Proof
We begin by considering some facts about the Chebyshev polynomials.
From the definition,
T
0
(x)
=
1,
T
1
(x)
= x. The triple recursion relation
is
the basis of an induction proof that
T,(x)
= 2n-lxn + lower degree terms
n
~
1
(4:7.19)
Thus
1
2n-l
Tn
(X)
= xn + lower degree terms
n
~
1 (4.7.20)
Since Tn(x) = cos(nO),
X=
cos
0,
0.::;;
e.::;;
'IT, the polynomial Tn(x)
attains relative maxima and minima at n + 1 points in [
-1,
1]:
j=
0,1,
...
,n
( 4.7 .21)
Additional values of
j do not lead to new values of
xj
because of the
periodicity of the cosine function. For these points,
Tn(x)
=
(-l)j
j=
0,1,
...
,n
(4.7.22)
and
-1
=
Xn
<
Xn-l
< · · · < x
1
<
Xo
= 1
230 APPROXIMATION OF FUNCTIONS
The
polynomial
Tn(x)j2n-t
has leading coefficient 1
and
Max
--T
x =
--
I
1 I 1
-1:s:x~1
2n-l
n(
)
2n-l
Thus
Tn
::::;;
lj2n-l.
Suppose that
1
T
<--
n
2n-1
(4.7 .23)
(4.7.24)
We
show that this leads to a contradiction.
The
assumption (4.7.24)
and
the definition (4.7.16) imply the existence
of
a polynomial
M(x)
= xn +
Q(x)
degree ( Q )
:;;;
n - 1
with
(4.7 .25)
Define
which
has
degree::::;;
n -
1.
We examine the sign
of
R(xi)
at
the points
of
(4.7.21). Using (4.7.22) and (4.7.24),
1
R(x
0
)
= R(1) =
2
n-l
-
M(l)
> 0
R(xr)
= -
2n1-l
-
M(xr)
= - [
2n1-l
+
M(xl)]
< 0
and
the
sign
of
R(xi)
is (
-1)i.
Since R has n + 1 changes of sign, R
must
have
at
least n zeros. But then degree (
R)
< n implies that R = 0;
thus M = (1/2n-l)Tn.
To
prove that no polynomial other than
(1/2n-l)Tn(x)
will minimize
(4.7.16), a variation of the preceding
proof
is used. We omit it.
ll
Consider
now
the problem
of
determining n + 1 nodes
xi
in [
-1,
1],
to be
used
in
constructing an interpolating polynomial Pn(x) that is to approximate
the
given function
f(x)
on [
-1,
1].
The
error in Pn(x)
is
!(
)
_ ( )
=
(x-
x
0
) • • •
(x-
xn)J<n+l)(t
)
X
Pn
X
(n
+ 1)!
'Ox
(4.7 .26)
The
value
J<n+
1
>(
~x)
depends on {xi},
but
the dependence
is
not
one that can
be
dealt
with explicitly. Thus to try to make
II/-
Pnlloo
as small as possible, we
NEAR-MINIMAX
APPROXIMATIONS
231
Table 4.7 Interpolation data for
J(x)
= ex
0
1
2
3
X;
.923880
.382683
-.382683
.923880
consider only the quantity
f(x;)
2.5190442
1.4662138
.6820288
.3969760
Max l{x- x
0
)
..
·
(x-
xJ
I
-l:sxsl
2.5190442
1.9453769
.7047420
.1751757
( 4.7 .27)
We
choose
{xj}
to minimize this quantity.
The
polynomial in (4.7.27)
is
of
degree
n + 1 and has leading coefficient
1.
From the preceding theorem, (4.7.27)
is
minimized by taking this polynomial
to
be Tn+
1
(x)j2n,
and
the minimum
value
of
(4.7.27) is
lj2n.
The nodes {xj} are the zeros
of
Tn+I(x),
and
these are
given in (4.7.14). With this choice of nodes,
Pn
=
In
and
(4.7 .28)
Example Let
f(x)
=ex
and x = 3.
We
use the Newton divided difference
form
of
the interpolating polynomial. The nodes, function values, and needed
divided differences are given in Table
4.7.
By
direct computation,
Max
lex-
1
3
(x)
I,;,
.00666
-ls;xs;l
(4.7.29)
whereas the bound m (4.7.28)
is
.0142. A graph
of
ex-
I
3
(x)
is shown in
Figure 4.7.
y
232 APPROXIMATION OF FUNCTIONS
The
error
II/-
/nlloo
is
generally not much worse than
Pn(/).
A precise result
is
that
n;::O
(4.7.30)
For
a proof, see Rivlin (1974, p.
13).
Actual numerical results are generally better
than that predicted by the bound in (4.7.29), as is illustrated
by
(4.7.29) in the
preceding examples
of
f(x)
= ex.
Forced
oscillation
of
the error Let
j(x)
E C[
-1,
1], and define
Let
II
F
11
(X) =
'2:'
C
11
• ,T,(x)
•=n
(4.7.31)
be nodes,
th~
choice of which
is
discussed below. Determine the coefficients
en.
k
by
forcing the error
f(x)
- Fn(x) to oscillate in the manner specified
as
neces-
sary by Theorem
4.10:
i=0,1,
...
,n+1
(4.7.32)
We have introduced another unknown,
En,
which
we
hope will be nonzero. There
are
n + 2 unknowns, the coefficients en,O•
•••
,
en,
n and
En,
and there are n + 2
equations in (4.7.32). Thus there
is
a reasonable chance that a solution exists.
If
a
solution exists, then by Theorem
4.9,
(4.7.33)
To
pick the nodes, note that if
en+
1
T
11
+
1
(x)
in (4.7.5)
is
nearly the minimax
error, then the relative maxima and minima in the minimax error,
f(x)
-
q!(x),
occur
at
the relative maxima and minima of
Tn+
1
(x).
These maxima and minima
occur when
Tn+
1
(x)
= ±
1,
and they are given by
x.
=cos(~)
' n + 1
i=O,l,
...
,n+1
(4.7.34)
These seem
an
excellent choice to use in (4.7.32) if F"(x)
is
to look like the
minimax approximation
q:(x).
The
system (4.7.32) becomes
n
}:'
cn.kTk(xi) + (
-l)iEn
=
J(xi)
k-0
i=0,1,
...
,n+1
(4.7.35)
NEAR-MINIMAX APPROXIMATIONS 233
Note
that
(
ki7T
)
Tk
(X;)
= cos ( k . cos
-l
(X;)) = cos
--
n + 1
Introduce
En=
cn,n+
1
/2.
The system (4.7.35) becomes
since
n+l
(
ki7T
)
2:"
cn,kcos
--1
=
f(x;)
k=O
n +
(
ki7T
) .
cos
--
= (
-1)'
n + 1
i=0,1,
...
,n+1
for k = n + 1
(4.7.36)
The
notation
I:"
means that the first
and
last terms are to be halved before
summation
begins.
To
solve (4.7.36), we need the following relations:
j = k = 0
or
n + 1
n
+,~
(
ij7T
) (
ik7T
) n + 1
I:
cos
--
cos
--
=
--
{
n
+ 1
i=O
n + 1 n + 1 2
O<j=k<n+1
(4.7.37)
0
j=t=k
k:s;n+1
The
proof
of
these relations is closely related
to
the proof of the relations in
Problem
42
of
Chapter 3. .
Multiply equation
i
in
(4.7.36) by
cos[(ij1r)j(n
+ 1)] for some 0 ::;j::::; n +
1.
Then
sum
on
i, halving the first and last terms.
This
gives
·n+l
n+l
(
ij7T
) (
ik7T
)
n+l
(
ij7T
)
2:"
en
k
I:"
cos
--
cos
--
=
2:"
f(x;)
cos
--
k=O
, i=O n + 1 n + 1
i=O
n + 1
Using the relations (4.7.37), all
but
one
of
the terms in the summation
on
k will
be
zero. By checking the two cases j =
0,
n + 1,
and
0 < j < n +
1,
the
same
formula is obtained:
2
n+
1 (
ij7T
)
cn,j
=
--1
I:"f(x;)cos
--
n + i=O n + 1
O::;j::;n+1
(4.7.38)
The
formula for
En
is
1
n+l
.
En=
n +
1
~"
(
-1)'f(x;)
t=O
(4.7.39)
There
are
a number
of
connections between this approximation
Fn(x)
and
the
Chebyshev expansion Cn(x). Most importantly, the coefficients
cn,j
are ap-
proximations
to
the coefficients
cj
in Cn(x). Evaluate formula (4.7.13) for
cj
by
234 APPROXIMATION OF FUNCTIONS
using the trapezoidal
ruk
(5.1.5) with n + 1 subdivisions:
21'"
c
1
= -
/(cos
8)
cos
(JO)
dO
'IT
0
2 n + 1 ( ( i
'fT
) ) (
ji
'IT
)
'IT
,;,
-I:"!
cos
--
· cos
--
·
--
= c .
'IT
;-o
n + 1 n + 1 n + 1
n,;
( 4.7 .40)
It
is
well-known that for periodic integrands, the trapezoidal numerical integra-
tion rule
is
especially accurate-(see Theorem
5.5
in Section 5.4). In addition, it
can be shown that
cnj
=
cj
+
c2(n+l)-j
+
c2(n+1)+j
+
c4(n+1)-j
+
...
(4.7.41)
If
the Chebyshev coefficients
en
in (4.7.2) decrease rapidly, then the approxima-
tion
Fn(x) nearly equals Cn(x), and it is easier to calculate than Cn(x).
Example Use
f(x)
=ex,
-1
~
x
~
1,
as
before. For n = 1, the nodes are
{x;}
=
{-1,0,1}
E
1
= .272, and
F
1
(
x)
= 1.2715 + 1.1752x (4.7.42)
For
the error, the relatiye maxima for
ex-
F
1
(x)
are given in Table 4.8. For
n = 3, £
3
= .00547 and
F
3
(x)
= .994526 + .995682x + .543081x
2
+ .179519x
3
(4.7.43)
The points of maximum error are given in Table 4.9. From Theorem 4.9
.00547
~
p3(/)
~
.00558
and this says that F
3
(x)
is
an excellent approximation to
qj(x).
To see that F
3
(x)
is
an approximation
to
C
3
(x),
compare the coefficients c
3
•
1
with the coefficients c
1
given in Table 4.3 in example (4.5.22) at the end of
Section 4.5.
The
results are
given
in
Table 4.10.
Table 4.8 Relative
maxima
of
lex-
F
1
(x)i
x
e"-
F
1
(x)
- 1.0
.272
.1614
-.286
1.0
.272
NEAR-MINIMAX
APPROXIMATIONS
23S
Table 4.9 Relative maxima of
lex-
Fix)l
X
ex-
FJ(x)
- 1.0 .00547
-.6832
-.00552
.0493 .00558
. 7324 - .00554
1.0 .00547
Table 4.10 Expansion coefficients for C
3
(x)
and F
3
(x)
to
e.r
j
cj
cn.j
0
2.53213176
2.53213215
1 1.13031821
1.13032142
2 .27149534
.27154032
3 .04433685
.04487978
4 .00547424
E
4
=
.00547424
As with the interpolatory approximation In(x ), care must be taken
il
j(x)
is
odd
or
even in [
-1,
1].
In such a case, choose n as follows:
If
f
is
{
~~~},
then chooses n
to
be
{
~v~~}
(4.7.44)
This ensures en+ I
=I=
0 in (4.7.5), and thus the
nodes
chosen
will
be
th~
correcl
ones.
An
analysis
of
the convergence of Fn(x) to
f(x)
is given in Shampin1(1970),
resulting in a bound similar to the bound (4.7.30) for In(x):
(4.7.45)
with
w(n)
empirically nearly equal to the bounding coefficient in
(4.7.3l).
Botll
In(x)
and Fn(x) are practical near-minimax approximations.
We
now give an algorithm for computing Fn(x), which can then be
ealuated
using the algorithm Chebeval
of
Section
4.5.
Algorithm Approx (c, E,
J,
n)
1.
Remark: This algorithm calculates the coefficients
cj
in
n
Fn(x)
=
f..'
cjTj(x)
j=O
-l::s;x:s;l
according to the formula (4.7.38), and E
is
calculat~
from
(4.7.39). The term c
0
should be halved before using
apri
thm
Chebeval.
236
APPROXIMATION
OF
FUNCTIONS
2.
Create
X;:=
cos(i1rj(n
+
1))
/;
:=
f(x;)
i=0,1,
...
,n+1
3. Do through step 8 for j =
0,
1,
...
, n +
1.
4.
sum
:=
[/
0
+ (
-1)
1
/n
+
1
]/2
5. Do through step 6 for i =
1,
...
, n.
6. sum
:=
sum
+/;cos
(ij7T
/(n
+
1))
7.
End loop on i.
8. c
1
= 2
sumj(n
+
1)
9.
End loop on
j.
10.
E
:=
cn+
1
/2
and exit.
The cosines in step 6 can be evaluated more efficiently
by
using the trigono-
metric addition formulas for the sine and cosine functions, but
we
have opted for
pedagogical simplicity, since the computer running time
will
still be quite small
with our algorithm. For the same reason,
we
have not used the
FFT
techniques of
Section 3.8.
Discussion
of
the Literature
Approximation theory
is
a classically important area of mathematics, and it
is
also an increasingly important tool in studying a wide variety of modern
problems in applied IlJathematics, for example, in mathematical physics and
combinatorics. The variety of problems and approaches
in
approximation theory
can be seen in the books of Achieser (1956), Askey (1975b), Davis (1963),
Lorentz (1966), Meinardus (1967), Powell (1981), and Rice (1964), (1968). The
classic work
on
orthogonal polynomials
is
Szego (1967), and a survey of more
recent work is given in Askey (1975a). For Chebyshev polynomials and their
many uses throughout applied mathematics and numerical analysis, see Fox and
Parker (1968) and Rivlin (1974). The related subject of Fourier series and
approximation
by
trigonometric polynomials was only alluded to in the text, but
it
is
of central importance in a large number of applications. The classical
reference is Zygmund
(1959).
There are many other areas of approximation
theory that we have not even defined. For an excellent survey of these areas,
including
an
excellent bibliography, see Gautschi (1975). A major area omitted in
the present text is approximation by rational functions. For the general area, see
Meinardus (1967, chap.
9)
and Rice (1968, chap.
9).
The generalization
to
BIBLIOGRAPHY
237
rational functions of the Taylor polynomial
is
called the Pade approximation; for
introductions, see Baker
(1975)
and Brezinski (1980). For the related area of
continued fraction expansions of functions,
see
Wall (1948). Many of the functions
that are of practical interest are examples of what are called the
special functions
of mathematical physics. These include the basic transcendental functions (sine,
log, exp, square root), and in addition, orthogonal polynomials, the Bessel
functions, gamma function, and hypergeometric function. There
is
an extensive
literature on special functions, and special approximations have been devised for
most of them. The most important references for special functions are in
Abramowitz and
Stegun (1964), a handbook produced under the auspices of the
U.S. National Bureau of Standards, and Erdelyi et
al.
(1953), a three volume set
that
is
often referred to
as
the "Bateman project."
For
a general overview and
survey of the techniques for approximating special functions, see Gautschi
(1975). An extensive compendium of theoretical results for special functions and
of methods for their numerical evaluation
is
given in Luke (1969), (1975), (1977).
For
a somewhat more current sampling of trends in the study of special
functions, see the symposium proceedings in Askey (1975b).
From the advent of large-scale
use
of computers
in
the 1950s, there has been a
need for high-quality polynomial or rational function approximations of the basic
mathematical functions and other special functions. As pointed out previously,
the approximation of these functions requires a knowledge of their properties.
But it also requires an intimate knowledge of the arithmetic of digital computers,
as surveyed in Chapter
1.
A general survey of numerical methods for producing
polynomial approximations
is
given in Fraser (1965), which has influenced the
organization of this chapter. For a very complete discussion of approximation of
the elementary functions, together with detailed algorithms, see Cody and Waite
(1980); a discussion of the associated programming project
is
discussed in Cody
(1984). For a similar presentation of approximations,
but
one that also includes
some of the more common special functions, see Hart et al. (1968). For an
extensive set of approximations for special functions, see Luke (1975), (1977).
For
general functions, a successful and widely used program for generating
minimax approximations
is
given in Cody et al. (1968). General programs for
computing minimax approximations are available in the
IMSL and NAG libraries.
Bibliography
Abramowitz, M., and I. Stegun (eds.) (1964). Handbook
of
Mathematical Func-
tions.
National Bureau of Standards, U.S. Government Printing Office,
Washington, D.C. (It
is
now published by Dover, New
York)
Achieser, N. (1956). Theory
of
Approximation (trans!. C. Hyman). Ungar, New
York.
Askey,
R.
(1975a). Orthogonal Polynomials and Special Functions. Society for
Industrial and Applied Mathematics,
Philadelphia.