Polynomial Interpolation - Purdue University

[Pages:14]Chapter 4

Polynomial Interpolation

In this chapter, we consider the important problem of approximating a function f (x), whose values at a set of

distinct points x0, x1, x2, . . . , xn are known, by a polynomial P (x) such that P (xi) = f (xi), i = 0, 1, 2, . . . , n. Such a polynomial is known as an interpolating polynomial. An approximation of this function is desirable

if f (x) is difficult to evaluate or manipulate like the error function

erf(x)

=

2

x e-t2 dt,

0

for example. Note that polynomials are easily evaluated, differentiated or integrated. Another purpose of

such an approximation is that a very long table of function values f (xi) may be replaced by a short table and a compact interpolation subroutine.

4.1 Linear Interpolation

Let f (x) be given at two distinct points xi and xi+1. Now, given f (xi) = fi and f (xi+1) = fi+1, we wish to determine the interpolating polynomial P1(x) of degree 1 (i.e., a straight line) that approximates f (x) on the interval [xi, xi+1].

Let P1(x) = x + . From the conditions P1(xi) = fi and P1(xi+1) = fi+1 we obtain the two equations

xi + = fi xi+1 + = fi+1,

i.e.,

xi 1 xi+1 1

=

fi fi+1

.

This system is clearly nonsingular since xi = xi+1, and it has a unique solution:

Therefore,

= fi+1 - fi xi+1 - xi

= fixi+1 - fi+1xi . xi+1 - xi

P1(x) = fi +

fi+1 - fi xi+1 - xi

(x - xi).

Example 4.1 Using linear interpolation, approximate sin 6.5 from a table which gives (a) sin 6 = 0.10453; sin 7 = 0.12187. (b) sin 0 = 0; sin 10 = 0.17365.

Chap. 4. Polynomial Interpolation

CS414 Class Notes

57

Solution

(a) P1(6.5) = 0.10453 +

0.12187-0.10453 7-6

0.11320.

(b)

P1(6.5)

=

0

+

0.17365 10

(6.5)

0.11287.

The first answer is correct to 5 decimals whereas the second answer is correct only to 2 decimals!

2

Conclusion: Linear interpolation is suitable only over small intervals.

4.2 Polynomial Interpolation

Since linear interpolation is not adequate unless the given points are closely spaced, we consider higher order interpolating polynomials. Let f (x) be given at the selected sample of (n + 1) points: x0 < x1 < ? ? ? < xn, i.e., we have (n + 1) pairs (xi, fi), i = 0, 1, 2, . . . , n. The objective now is to find the lowest degree polynomial that passes through this selected sample of points.

f(x)

Pn(x)

f(x)

x0

x1

x2

x3 ... xn

Figure 4.1: Interpolating the function f (x) by a polynomial of degree n, Pn(x).

Consider the nth degree polynomial

Pn(x) = a0 + a1x + a2x2 + ? ? ? + anxn.

We wish to determine the coefficients aj, j = 0, 1, . . . , n, such that

Pn(xj ) = f (xj), j = 0, 1, 2, . . . , n.

These (n + 1) conditions yield the linear system

a0 + a1x0 + a2x20 + ? ? ? + anxn0 = f0

a0 + a1x1 + a2x21 + ? ? ? + anxn1 = f1

a0 ...

+

a1x2 ...

+

a2...x22

+

?

?

?

+

an...xn2

=

f2 ...

a0 + a1xn + a2x2n + ? ? ? + anxnn = fn

or V a = f , where aT = (a0, a1, . . . , an), f T = (f0, f1, . . . , fn), and V is an (n + 1) ? (n + 1) matrix given by

1 x0 x20 ? ? ? xn0

V

=

1 1 ...

x1 x2 ...

x21 x...22

??? ???

xn1 x...n2

.

1 xn x2n ? ? ? xnn

Chap. 4. Polynomial Interpolation

CS414 Class Notes

58

The matrix V is known as the Vandermonde matrix. It is nonsingular with nonzero determinant

For n = 3, for example,

n-1

det(V ) = (xi - xj ).

j=0

i>j

det(V ) = (x1 - x0)(x2 - x0)(x3 - x0)(x2 - x1)(x3 - x1)(x3 - x2)

Hence a can be uniquely determined as a = V -1f . Another proof of the uniqueness of the interpolating polynomial can be given as follows.

Theorem 4.1 Uniqueness of interpolating polynomial. Given a set of points x0 < x1 < ? ? ? < xn, there exists only one polynomial that interpolates a function at those points.

Proof Let P (x) and Q(x) be two interpolating polynomials of degree at most n, for the same set of points x0 < x1 < ? ? ? < xn. Also, let

R(x) = P (x) - Q(x).

Then R(x) is also a polynomial of degree at most n. Since

P (xi) = Q(xi) = fi,

we have, R(xi) = 0 for i = 0, 1, . . . , n. In other words R(x) has (n + 1) roots. From the Fundamental Theorem of Algebra however, R(x) cannot have more than n roots. A contradiction! Thus, R(x) 0, and

P (x) Q(x).

2

Example 4.2 Given the following table for the function f (x), obtain the lowest degree interpolating polynomial.

xk ?1 2 4 fk ?1 ?4 4

Solution Since the number of nodes = 3 = n + 1, therefore n = 2. Let

P2(x) = a0 + a1x + a2x2,

that satisfies the following conditions at the points x0 = 1; x1 = 2, x2 = 4:

P (-1) = -1; P (2) = -4; P (4) = 4.

1 -1 1 a0 -1

1 2 4 a1 = -4 .

1 4 16

a2

4

Using Gaussian elimination with partial pivoting, we compute the solution

a0

-4

a = a1 = -2 ,

a2

1

Therefore, the interpolating polynomial is given by

P2(x) = -4 - 2x + x2.

We can use this approximation to estimate the root of f (x) that lies in the interval [2, 4]:

=

2

+

4 2

+

16

=

1

+

5.

Chap. 4. Polynomial Interpolation

CS414 Class Notes

59

2

An important remark is in order. One, in general, should not determine the interpolating polynomial by solving the Vandermonde linear system. These systems are surprisingly ill-conditioned for n no larger than 10. For example, for 0 < x0 < x1 < ? ? ? < xn = 1 uniformly distributed in [0,1], large n yields a Vandermonde matrix with almost linearly dependent columns, and the Vandermonde system becomes almost singular.

A more satisfactory form of the interpolating polynomial is due to Lagrange,

Pn(x) = f (x0)L0(x) + f (x1)L1(x) + ? ? ? + f (xn)Ln(x).

Here Lj(x), for 0 j n, are polynomials of degree n called fundamental polynomials or Lagrange polynomials. From the (n + 1) conditions

Pn(xi) = f (xi), 0 i n,

we see that

Lj(xi) =

0 i=j 1 i=j

The above property is certainly satisfied if we choose

Lj (x)

=

(x - x0)(x - x1) ? ? ? (x - xj-1)(x - xj+1) ? ? ? (x - (xj - x0)(xj - x1) ? ? ? (xj - xj-1)(xj - xj+1) ? ? ? (xj

xn) - xn)

i.e.,

Lj (x)

=

n

i=0

(x (xj

- xi) - xi)

.

i=j

The fact that Lj(x) is unique follows from the uniqueness of the interpolation polynomial. Recall that there is one and only one polynomial of degree n or less that interpolates f (x) at the (n + 1)

nodes x0 < x1 < ? ? ? < xn. Thus, for this given set of nodes, both forms yield the same polynomial!

Vandermonde Approach.

n

Pn(x) = aixi

Va=f

i=0

Operation count: It can be shown that the Vandermonde system V a = f can be solved in O(n2) arithmetic operations rather than O(n3) operations because of the special form of the Vandermonde matrix.

Lagrange Approach.

n

Pn(x)

=

n j=0

fj Lj (x)

=

n j=0

gj

(x - xi

i=0

(x - xj )

)

=

n j=0

gj

?

j

where

gj =

f (xj)

n

.

(xj - xi)

i=0

i=j

Operation count: The Lagrange form avoids solving the ill-conditioned Vandermonde linear, systems, and

evaluates the polynomials

j

=

(x

1 - xj)

n

(x

i=0

-

xi)

Chap. 4. Polynomial Interpolation

CS414 Class Notes

60

for a given x. The cost of these operations for each j = 0, 1, 2, . . . , n are:

n

u(x) = (x - xi)

j

=

ui=(0x) (x - xj )

gj = n fj

(xj - xi)

i=0

i=j

(2n + 1) ops (n + 1) ops 2n(n + 1) ops (evaluated once only)

Now, the cost for evaluating the polynomial

n

Pn(x) = gjj

j=0

(2n + 1) ops

Therefore, the total number of operations are 5n + 3 once gj's are evaluated.

4.3 Error in Polynomial Interpolation

Let en(x) be the error in polynomial interpolation given by

en(x) = f (x) - Pn(x).

Since Pn(xi) = f (xi), we have

en(xi) = 0 x0 < x1 < ? ? ? < xn.

In other words en(x) is a function with at least (n + 1) roots. Hence, we can write

f(x) Pn(x)

f(x)

x0

x1

x2

x3 ... xn

Figure 4.2: Computing the error in polynomial interpolation.

en(x) = (x - x0)(x - x1) ? ? ? (x - xn)h(x) = (x) ? h(x)

where h(x) is a function of x. Let us define

g(z) = en(z) - (z)h(x)

where x = xj , for 0 j n. Thus,

g(x) = 0, g(xj ) = 0,

0 j n,

Chap. 4. Polynomial Interpolation

CS414 Class Notes

61

i.e., g(z) has (n+2) distinct roots. Assuming that f (x) has at least n+1 continuous derivatives, the (n+1)th

derivative of g(z), i.e.,

d(n+1) dz(n+1)

g(z)

g(n+1)(z)

has at least one root x. Hence,

x0

x2

x4

g(z)

x1

x3

g!(z)

g!!(z)

g!!!(z)

giv(z) x*

Figure 4.3: The (n + 1)th derivative of a function g(x) with n + 2 roots (n = 3) must have at least one root.

g(n+1)(z) = e(nn+1)(z) - h(x)(n+1)(z), where (n+1)(z) = (n + 1)! Moreover, since en(z) = f (z) - Pn(z), therefore,

e(nn+1)(z) = f (n+1)(z) - Pn(n+1)(z) in which the last term Pn(n+1)(z) is zero. Using the above equations, we have

g(n+1)(z) = f (n+1)(z) - h(x)(n + 1)!

We had shown earlier that g(n+1)(z) has at least one root x in the interval containing x and xj , j = 0, . . . , n. This implies that g(n+1)(x) = 0, and therefore,

h(x)

=

f (n+1)(x) (n + 1)!

Consequently, we have the following representation of the error in polynomial interpolation

en(x)

=

f (x)

-

Pn(x)

=

(x

-

x0)(x

-

x1)

?

?

?

(x

-

xn)

f (n+1)(x) (n + 1)!

where x is in the interval containing x and xj, j = 0, . . . , n.

Chap. 4. Polynomial Interpolation

CS414 Class Notes

62

Example 4.3 Obtain the Lagrange interpolating polynomial from the data

sin 0 = 0,

sin

6

=

1 2,

sin

3

=

3 2,

sin

2

=

1,

and

use

it

to

evaluate

sin

12

and

sin

4

Solution The interpolating polynomial is given by

P3(x) = f0L0(x) + f1L1(x) + f2L2(x) + f3L3(x)

=

1 2

L1(x)

+

3 2

L2(x)

+

L3(x)

At

x

=

12

,

L1

12

=

12

-

0

6

-

0

12

-

3

6

-

3

12

-

2

6

-

2

=

15 16

L2

12

=

12

-

0

3

-

0

12

-

6

3

-

6

12

-

2

3

-

2

=

-

5 16

L3

12

=

12

-

0

2

-

0

12

-

6

2

-

6

12

-

3

2

-

3

=

1 16

Therefore,

P3

12

= 0.26062.

SincAetsxin=124=,

0.25882, we L1

have |error| =

4

= 0.5625,

0.0018. L2

4

= 0.5625,

L3

4

= 0.0625.

Therefore

P3

4

= 0.70589.

SincLeetsinus4e=sti0m.7a0te71t1h,eweerrhoarveen|(exr)roart|x==0.400: 122.

e3

4

=

4

-

0

4

-

6

4

-

3

4

-

2

f (4)(x) 4!

Since f (x) = sin x; f (4)(x) = sin x, and thus,

|f (4)(x)| 1.

Therefore,

e3 4

4

?

12

?

12

?

4

?

1 4!

=

1.76

?

10-3.

Similarly,

en(x)

at

x

=

12

is

bounded

by

e3 12

12

?

12

?

4

?

5 12

?

1 4!

=

2.94

?

10-3.

2

Chap. 4. Polynomial Interpolation

CS414 Class Notes

63

4.4 Runge's Function and Equidistant Interpolation

In this section, we discuss a problem that arises when one constructs a global interpolation polynomial on a

fairly large number of equally spaced nodes. Consider the problem of approximating the following function

on the interval [-1, 1]

f (x)

=

1 1 + 25x2

using the interpolation polynomial Pn(x) on the equidistant nodes

xj

=

-1

+

2j , n

j = 0, 1, . . . , n.

Figure ?? shows how the polynomials P4(x), P8(x), P12(x), and P16(x) approximate f (x) on [-1, 1]. We see

n= 4 1

n= 8 1

0.8 0.8

0.6

0.6

0.4

0.2 0.4

0

0.2 -0.2

0 -0.2

-0.4 -0.6 -0.8

-0.4

-1

-1

-0.8 -0.6 -0.4 -0.2

0

0.2

0.4

0.6

0.8

1

-1

-0.8 -0.6 -0.4 -0.2

0

0.2

0.4

0.6

0.8

1

n = 12 1

n = 16 2

0.5

0

0 -2

-0.5

-4 -1

-1.5

-6

-2 -8

-2.5 -10

-3

-12 -3.5

-4

-14

-1

-0.8 -0.6 -0.4 -0.2

0

0.2

0.4

0.6

0.8

1 -1

-0.8 -0.6 -0.4 -0.2

0

0.2

0.4

0.6

0.8

1

Figure 4.4: Equal spaced interpolation of Runge's function.

that as n increases, the interpolation error in the central portion of the interval diminishes, while increasing rapidly near the ends of the interval. Such behavior is typical in equidistant interpolation with polynomials of high degree. This is known as Runge's phenomenon, and the function f (x) is known as Runge's function.

If the (n + 1) interpolation nodes are not chosen equally spaced, but rather placed near the ends of the interval at the roots of the so called Chebyshev polynomial of degree (n + 1), the problem with Runge's function disappears (see Fig. ??). The roots of the Chebyshev polynomial of degree (n + 1) are given by

xk = - cos

k n

+ +

1

2

1

,

k = 0, 1, . . . , n.

 ................
................

In order to avoid copyright disputes, this page is only a partial summary.

Google Online Preview   Download