Download The Generalized Lagrange Neville Interpolation | PHYS 6720 and more Study notes Physics in PDF only on Docsity!
Physics 6720 { Neville Interp olation { Novemb er 13, 2004
When it is exp ensive or dicult to evaluate a function f (x) at an arbitrary
value of x, we might consider, instead, interp olating from a table of values.
Consider the exp onential integral,
Ei(x) = P
Z
x
e
t
dt=t
tabulated b elow for small x.
x Ei(x)
This function diverges as log(x) at x = 0. To evaluate it requires doing
the integral, or summing a series, so it is somewhat exp ensive. Supp ose we
want the value of Ei(0:15). A very crude approximation cho oses either of the
two nearby values Ei(0:1) = 1 : 6228 or Ei(0:2) = 0 :8218. A common and
somewhat b etter approach makes a linear interp olation. Since 0.15 is midway
b etween 0.1 and 0.2 the linear interp olation just averages the two values, giving
1 Lagrange Interp olation
Let's generalize the linear interp olation by denoting the values in the table
by (x i ; y i ) for i = 0 ; : : : ; n. Then a linear function interp olating the rst two
values can b e written as
P
(x) = y 0
x x 1
x 0
x 1
x x 0
x 1
x 0
Since every term here is a constant except x, it is easy to see that the p olyno-
mial P 0 ; 1 has degree 1 at most. It also satis es
P
(x 0
) = y 0
P
(x 1
) = y 1
as is required for an interp olation of the two values.
The linear interp olation can also b e written in the form
P
(x) = y 0
L
(x) + y 1
L
(x)
where
L
(x) =
x x 1
x 0
x 1
L
(x) =
x x 0
x 1
x 0
These are examples of Lagrange interp olating p olynomials. They have the
prop erty
L
(x 0
) = 1 ; L
(x 1
) = 0 ; L
(x 0
) = 0 ; L
(x 1
that is, L
i
(x j
i;j
1.1 Generalized Lagrange Interp olation
Returning to our table, if we examine the change in the function b etween 0.
and 0.2 and b etween 0.2 and 0.3 we see that there is a dramatic change in
slop e. So clearly our result would b e improved if we selected an interp olating
p olynomial with curvature. It takes three p oints to x a parab ola. Why
stop there? We could include all n + 1 p oints, resulting in a p olynomial of
degree n. Such a p olynomial can b e constructed from higher degree Lagrange
interp olating p olynomials
L
(n)
k
(x) =
n Y
i=0;i 6 =k
x x i
x k
x i
They satisfy the sp ecial prop erty,
L
(n)
k
(x j
j;k
The construction gives
P
0 ;:::;n
(x) =
n X
k =
L
(n)
k
(x)y k
If we set x = x j
for any j , we see that all the p olynomials except L
(n)
j
vanish,
leaving P 0 ;:::;n
(x j
) = y j
. Thus the p olynomial passes through all the p oints, as
is the case for the Ei(x) function, as we can see in the gure. So we should
avoid high degree p olynomials when interp olating evenly spaced p oints.
2 Neville algorithm
Constructing the interp olating p olynomial is somewhat tedious. If our goal is
merely to get the interp olated value, and we don't care to know the co ecients
of the p olynomial, we may use the Neville algorithm. This algorithm starts
from the requested interp olation p oint x and generates a table of the form
x 0
P
(x)
P
(x)
x 1
P
(x) P 012
(x)
P
(x) P 0123
(x)
x 2
P
(x) P 123 (x) P 01234 (x)
P
(x) P 1234
(x)
x 3
P
(x) P 234 (x)
P
(x)
x 4
P
(x)
where P i;:::;k
(x) is the p olynomial interp olating the p oints x i
; x i+
: : : ; x j
, evalu-
ated at the p oint x. That is to say, the second column comes from interp olating
only one p oint each, meaning that it repro duces trivially the original values:
P
i
(x) = y i
. The third column comes from interp olating pairs of adjacent
p oints. Notice that P 0 ; 1
is the p olynomial we constructed ab ove. The fourth
column comes from interp olating three contiguous p oints. We illustrated P 123
ab ove. The last column has only one entry and comes from interp olating all
the p oints.
The identity
P
i;:::;`
(x) =
(x x j
)P
i;:::;j 1 ;j +1;:::`
(x) (x x k
)P
i;:::;k 1 ;k +1;:::`
(x)
x k
x j
is used to calculate the values in any column from pairs of values in one column
to the left. For example the rule gives
P
(x) =
(x x 3
)P
(x) (x x 1
)P
(x)
x 1
x 3
The student should verify that P 123
(x i
) = y i
for i = 1 ; 2 ; 3, as it should.
Here we apply the Neville scheme to the example of the exp onential inte-
gral, interp olated at x = 0 :15.
3 Algorithm
As we add one more p oint to the list of interp olated p oints, Eq. (1) tells us
how to compute the value of the new interp olating p olynomial. The matrix
structure of the interp olating table suggests intro ducing a matrix Q i;j
to hold
the ith row and j th column in the format
x 0
Q
Q
x 1
Q
Q
Q
Q
x 2
Q
Q
Q
Q
Q
x 3
Q
Q
Q
x 4
Q
so that Q i;j
= P
i j;i j +1;:::;i
. Converting Eq (1) to the Q's gives
Q
i;j
(x) =
(x x i j
)Q
i;j 1
(x) (x x i
)Q
i 1 ;j 1
(x)
x i
x i j
for i = j; : : : n 1 for each j = 1 ; 2 ; : : : and with the initial values Q i; 0
= y i
We leave the construction of the pseudo co de and the computer co de as an
exercise. With care the algorithm can b e designed so it keeps only one column
of the matrix, saving computer memory. The trick is to run the inner lo op
over i backwards so the fresh values for the j th column do not overwrite values
still needed from the j 1st column.
