4a interpolation
pdf
keyboard_arrow_up
School
University of British Columbia *
*We aren’t endorsed by this school
Course
316
Subject
Mathematics
Date
Nov 24, 2024
Type
Pages
43
Uploaded by 1234sppencer
4a. Polynomial interpolation
MACM 316
1/43
Prelude: Interpolation versus Fitting
•
The human brain and vision system have a remarkable natural ability to recog-
nize the patterns underlying graphical data
•
The aim of Section 4 is to design algorithms that mimic this process . . .
Points from two different
Interpolate piecewise with
Interpolate with a single
Fit to a smooth function
data sources
lines (or polynomials)
smooth function
(no interpolation)
Interpolation
curve passes through all points
versus
Fitting
curve is “closest to” all points
r
4a. Polynomial interpolation
MACM 316
2/43
In Sections 4a/b/c we will study
interpolation
and
fitting
:
•
for a set of
n + 1
data points
(x
i
, y
i
)
with
i = 0, 1, 2, ... , n
(
=
n
1
!!
•
using polynomials of degree
p
>
1
(
= p
as small as possible!!
Our plan of attack:
4a:
piecewise linear interpolation (polynomial degree
p = 1
)
[ too simple since most real processes are smooth ]
4a:
smooth interpolation using a single polynomial (special case
p = n
)
[ imposing too much smoothness can generate “wiggles” ]
4b:
piecewise interpolation with polynomials of degree
1
<
p
⌧
n
(
=
splines
[ compromise on smoothness to get better results ]
4c:
fit with low-order polynomials (
1
6
p
⌧
n
) & other functions
(
=
least squares
[ when interpolating doesn’t make any sense – REAL DATA IS NOISY! ]
simple
4a. Polynomial interpolation
MACM 316
3/43
Taylor Polynomials are Unsuited for Interpolation
•
We’ve already used polynomials to approximate functions –
Taylor polynomials!
•
This approach is useful only if we know values of a function
f(x)
::::::::::::::::
and all its
::::::::::::::::::
derivatives at a single point
x = a
:
P
n
(x) = f(a) + f
0
(a)(x
-
a) +
f
00
(a)
2!
(x
-
a)
2
+
· · ·
+
f
(n)
(a)
n!
(x
-
a)
n
•
Disadvantages of the Taylor polynomial:
⇤
In practice, we might have measurements of some quantity
f(a)
, but seldom
its rate of change (first derivative) . . . not to mention higher derivatives!
⇤
Error in
P
n
(x)
grows rapidly away from
a
– error
/
(x
-
a)
n+1
⇤
The approximation is exact only in the limit
n
! 1
, and is often poor even
for moderate-to-large
n
.
Conclude:
Taylor polynomials are not useful for interpolating a list of data points
. . .
BUT
Taylor’s Theorem will still play an important role!
x
at
I
error
a
x
a
htt
FYI
for
error
estimates
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
4/43
Section 4a. Polynomial Interpolation
Consider this common scenario:
•
You have a list of measurements or “table of values” for 2 quantities,
x
and
y
•
Data derives from a
:::::::::::::::::::
continuous process
y = f(x)
BUT
f(x)
is unknown!
•
You want to find an approximation to
f(x)
that:
(a)
passes through all data points AND
(b)
fills in the gaps between the data
=
)
called
interpolation
Example 1:
Given a list of points
(x
i
, y
i
)
for
i = 0, 1, ... , 12
:
x
i
y
i
0.0
0.0000
0.5
1.5163
1.0
3.6788
-
at
x =
1.75?
2.0
5.4134
2.5
5.1303
3.0
4.4808
4.0
2.9305
5.0
1.6845
6.0
0.8924
7.0
0.4468
8.0
0.2147
9.0
0.1000
10.0
0.0454
Basic question:
How can we approximate the underlying function
f(x)
at another
point that isn’t in the list, say
x =
1.75?
4a. Polynomial interpolation
MACM 316
5/43
Basic question:
How can we approximate the underlying function
f(x)
at another
point that isn’t in the list (say
x =
1.75)?
Two-step solution:
1. Interpolate data with a function
P(x)
2. Use
P(1.75)
to approximate
f(1.75)
There are many ways to construct a continuous interpolating function
P(x)
:
•
piecewise straight line segments (
), continuous but not smooth
•
a single smooth function (
)
=
)
infinitely many choices for polynomial degree or functional form
Key point:
Interpolation comes in 2 main flavours:
piecewise
and
smooth
4a. Polynomial interpolation
MACM 316
6/43
Brief Look: Piecewise Linear Interpolation
Example 1 (again):
•
Find linear functions that interpolate
f(x)
on sub-intervals
[1, 2]
and
[2, 2.5]
•
Use these interpolants to approximate the two values
f(1.75)
and
f(2.25)
•
Also known as a
linear spline
x
i
y
i
0.0
0.0000
0.5
1.5163
1.0
3.6788
2.0
5.4134
2.5
5.1303
3.0
4.4808
4.0
2.9305
5.0
1.6845
6.0
0.8924
7.0
0.4468
8.0
0.2147
9.0
0.1000
10.0
0.0454
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
7/43
First “interpolant” on interval
[1, 2]
:
•
Straight line joining
(1.0,
3.6788)
to
(2.0, 5.4134)
is
[ check it! ]
P(x) =
1.9442 + 1.7346x
•
Interpolate – use
P(x)
to approximate
f(1.75)
f(1.75)
⇡
P(1.75) =
4.9798
Zoomed-in view
•
Assume we know the exact value
f(1.75) = 5.3218
. The interpolation error is
|
4.9798
-
5.3218
|
|
5.3218
|
= 0.0643
[ relative error ]
•
We could approximate
f(2.25)
by
extrapolating
P(x)
outside its interval of valid-
ity
[1, 2]
using
P(2.25) =
1.9422 + 1.7346(2.25)
⇡
5.8471
[ larger error – see plot ]
. . . extrapolation is always risky!
Second interpolant on interval
[2, 2.5]
:
(safer than extrapolating)
Q(x) =
6.5458
-
0.5662x
=
)
Q(2.25) = 5.2719
[ much smaller error ]
7
Qu
3.6788
5.4131318871
1
1.9442
1.7436
Pix
4.9798
relative
15.3218
P
1.7511
5.3218
0.0643
6
error
1.9442
1.734612251
5.8471
larger
654
smaller
Q
2.25
4a. Polynomial interpolation
MACM 316
8/43
Summary: Interpolation versus Extrapolation
•
With interpolation, error can (sometimes) be controlled and so there’s a chance
that the approximation is a good one
•
Extrapolation is always risky!
A
extrapolation
linear
fit
4a. Polynomial interpolation
MACM 316
9/43
Higher Order Interpolation (
p = n
)
Observe:
piecewise linear interpolation is not smooth (at each
x
i
there is a “cusp”).
Contrast with real processes that typically give rise to smooth functions.
Question:
How to find a single smooth function
f(x)
that passes through all points?
=
)
use a higher degree polynomial!
[ how high? ]
Suppose we have
n + 1
points labelled
(x
0
, y
0
), (x
1
, y
1
), ... , (x
n
, y
n
)
where
y
i
= f(x
i
)
for
i = 0, 1, ... , n
AND
x
0
<
x
1
<
· · ·
<
x
n
.
Theorem:
(Thms. 3.2 & 3.3, p. 108–109)
Given
n + 1
distinct points
x
i
,
there exists a
:::::::::::
unique
polynomial of
::::::::::::::
degree
n
P
n
(x) = a
0
+ a
1
x + a
2
x
2
+
· · ·
+ a
n
x
n
-
“interpolating polynomial”
that satisfies
P
n
(x
i
) = y
i
for
i = 0, 1, 2, ... , n
. Furthermore
f(x) = P
n
(x) + E
n
(x)
where the (absolute) interpolation error is given by
E
n
(x) =
f
(n+1)
(c)
(n + 1)!
(x
-
x
0
)(x
-
x
1
)
· · ·
(x
-
x
n
)
-
“error term”
for some
c
2
[x
0
, x
n
]
.
[
E
n
(x)
looks similar to Taylor remainder term! ]
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
10/43
Example 1 (Yet Again)
Example 1 (slide 7):
The quadratic interpolating polynomial for
f(x)
on the interval
[1, 2.5]
with three points (degree 2) is
P
2
(x) =
-
1.1235 + 6.3362x
-
1.5339x
2
[ we’ll soon have an algorithm to derive this ]
Check (by substitution):
P
2
(
1
) = 3.6788
,
P
2
(
2
) = 5.4134
,
P
2
(
2.5
) = 5.1303
•
Then we can approximate
f(1.75)
using
P
2
(1.75) =
5.2674
•
Relative error
⇡
5.3218
-
5.2674
5.3218
⇡
0.0102
=
)
smaller than for linear approx’n!
[ Recall:
P
1
(x) = 1.9442 + 1.7346x
,
P
1
(1.75) = 4.9798
,
relative error
⇡
0.0643
]
Graphical Comparison:
linear
and
quadratic
interpolants
Warning:
Increasing the polynomial degree does
NOT
mean the approximation improves!
(look ahead to slide 42)
4a. Polynomial interpolation
MACM 316
11/43
Linear Algebra of Interpolating Polynomials
•
Start with 3
ordered
points:
(x
0
, y
0
), (x
1
, y
1
), (x
2
, y
2
)
with
x
0
<
x
1
<
x
2
•
Fit a quadratic (degree
n = 2
):
y = c
0
+ c
1
x + c
2
x
2
•
Substitute the 3 points into the polynomial:
y
0
= c
0
+ c
1
x
0
+ c
2
x
2
0
,
y
1
= c
0
+ c
1
x
1
+ c
2
x
2
1
,
y
2
= c
0
+ c
1
x
2
+ c
2
x
2
2
•
Rewrite equations in matrix-vector form:
2
6
6
4
1
x
0
x
2
0
1
x
1
x
2
1
1
x
2
x
2
2
3
7
7
5
2
4
c
0
c
1
c
2
3
5
=
2
4
y
0
y
1
y
2
3
5
[ pattern is obvious! ]
•
Generalize to
n + 1
points:
(x
0
, y
0
), (x
1
, y
1
), ... , (x
n
, y
n
)
2
6
6
6
6
6
6
4
1
x
0
x
2
0
· · ·
x
n
0
1
x
1
x
2
1
· · ·
x
n
1
.
.
.
.
.
.
.
.
.
.
.
.
1
x
n
x
2
n
· · ·
x
n
n
3
7
7
7
7
7
7
5
|
{z
}
V
2
6
6
6
6
6
4
c
0
c
1
c
2
.
.
.
c
n
3
7
7
7
7
7
5
|
{z
}
c
=
2
6
6
6
6
6
4
y
0
y
1
y
2
.
.
.
y
n
3
7
7
7
7
7
5
|
{z
}
y
or
Vc = y
•
This is an
(n + 1)
⇥
(n + 1)
linear system for the polynomial coefficients
c
.
•
V
is called a
Vandermonde matrix
.
4a. Polynomial interpolation
MACM 316
12/43
Aside: Is
V
Invertible?
For any set of points
x
0
,
x
1
, . . . ,
x
n
, the Vandermonde matrix is
V =
2
6
6
6
6
6
4
1 x
0
x
2
0
· · ·
x
n
0
1 x
1
x
2
1
· · ·
x
n
1
.
.
.
.
.
.
.
.
.
.
.
.
1 x
n
x
2
n
· · ·
x
n
n
3
7
7
7
7
7
5
Fact:
If any two of the
x
i
are equal, then
V
has two identical rows and det
(V) = 0
!
Theorem:
V
is always invertible if the
x
i
are distinct.
Proof (by contradiction):
Suppose
V
is singular (non-invertible).
•
Then
Vc = 0
for some coefficient vector
c
6
= 0
.
•
It follows that the polynomial
q(x) = c
0
+ c
1
x + c
2
x
2
+
· · ·
+ c
n
x
n
is zero at
x = x
0
,
x = x
1
, . . . and
x = x
n
.
•
In other words,
q(x)
is a polynomial of degree
n
that has
(n + 1)
distinct roots.
•
Therefore,
q(x)
⌘
0
and
c
0
= c
1
= c
2
=
· · ·
= c
n
= 0
.
•
This is a contradiction because
c
6
= 0
.
Conclude:
V
must be invertible.
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
13/43
Example 2: Vandermonde Matrix
Consider these (out-of-date)
gasoline price data:
Year
1986
1988
1990
1992
1994
1996
Price ($/gal)
0.927
0.946
1.164
1.127
1.112
1.147
(a) Find a smooth polynomial that interpolates the data.
(b) Estimate the gas price in 1991.
Solution:
6 data points
=
)
use a polynomial of degree
n = 5
•
Set up Vandermonde matrix using some fancy M
ATLAB
’ing:
year = [1986 : 2 : 1996]’;
% column vector
V = repmat(year, 1, n);
% replicate / tile an array
V(:,n) = 1;
V = cumprod(V, 2, ’reverse’); % cumulative row product
•
OR an alternate method uses M
ATLAB
’s built-in
vander
function:
V = vander(year);
Warning:
⇤
vander
orders columns opposite to these notes – largest to smallest!
⇤
Entries in coefficient vector are likewise reversed:
(c
n
, ..., c
2
, c
1
, c
0
)
.
⇤
Advantage: This is the same order that
polyval
expects coefficients!
⇤
Code
vanderex.m
is posted on Canvas. Play with it!
4a. Polynomial interpolation
MACM 316
14/43
Example 2: M
ATLAB
Code
year
= [1986 : 2 : 1996]’;
vanderex.m
price = [0.927, 0.946, 1.164, 1.127, 1.112, 1.147]’;
V = vander(year)
% built-in function
c = V
\
price
y = linspace(min(year), max(year));
% dense points for plotting
p = polyval(c, y);
% evaluate polynomial with coefficients c
cond(V)
plot(year, price, ’o’, y, p, ’-’)
Part (b): Estimate 1991 price
>> polyval(c, 1991)
ans = 1.1797
Q:
Why the wiggles??
4a. Polynomial interpolation
MACM 316
15/43
Eliminating Wiggles
The “
::::::::::::
wiggles” result from build-up of floating point errors because . . .
•
Condition number of
V
is enormous
=
)
9.7738e+30
•
PLU decomposition in “
\
” introduces large errors in the coefficients
c
i
•
Coefficients
c
i
alternate in sign
=
)
causes big subtractive cancellation errors when evaluating
the polynomial!
Solution:
Large entries in
V
can be avoided
by “
::::::::::::::::::::::
shifting years”
[1986, 1996]
!
[0, 10]
>> year2 = year - min(year);
>> V = vander(year2);
>> cond(V)
ans = 5.6465e+05
% much smaller!
Wiggles are eliminated!
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
16/43
Exercise (Homework):
The code provides a more complete list of annual gas price
data from 1986 to 2017.
1. Interpolate the data for the entire period, again using only even-numbered years.
How do the results change?
2. Interpolate all price data, including every year. Describe how your interpolating
polynomial changes.
4a. Polynomial interpolation
MACM 316
17/43
Disadvantage of Standard Form
When using the “standard” expanded form of the polynomial
q(x) = c
0
+ c
1
x + c
2
x
2
+
· · ·
+ c
n
x
n
and solving for polynomial coefficients
c
i
directly:
•
The Vandermonde matrix is
:::::::::::::::::::::::
ill-conditioned, even for moderate
n
, which can lead
to large errors in the coefficients
c
i
.
•
Solving a
::::::::::
dense matrix system requires
O(n
3
)
floating point operations, which
is expensive.
Question:
Recognizing that
q(x)
is
UNIQUE
. . . Is there an alternate method to
determine
q(x)
that is cheaper and less prone to build-up of round-off errors?
Restate the Problem:
View
q(x)
as a sum of
monomial basis polynomials
M
j
(x)
:
q(x) =
n
X
j=0
c
j
M
j
(x)
with
M
j
(x) = x
j
Is there a
better basis
than the monomials
x
j
?
4a. Polynomial interpolation
MACM 316
18/43
Lagrange Form of Interpolating Polynomial
Text
§3.1
•
Given
n + 1
ordered
points:
(x
0
, y
0
), (x
1
, y
1
), ... , (x
n
, y
n
)
•
Idea:
Construct a
:::::::::::
simple
:::::::::
basis
:::::::::::::::::::
polynomial that
equals
1
at
x
0
and
0
at every other point:
L
0
(x) :
x
0
x
1
x
2
· · ·
x
n
#
#
#
#
1
0
0
· · ·
0
•
The following does the trick for
x
0
:
L
0
(x) =
(x
-
x
1
)
(x
0
-
x
1
)
(x
-
x
2
)
(x
0
-
x
2
)
· · ·
(x
-
x
n
)
(x
0
-
x
n
)
=
n
Y
i=1
x
-
x
i
x
0
-
x
i
[ degree
n
]
•
Construct similar basis polynomials for all
x
j
:
L
j
(x) =
n
Y
i=0
i
6
=j
x
-
x
i
x
j
-
x
i
=
)
L
j
(x
i
) =
⇢
1
if i = j
0
if i
6
= j
[ “
j
th
Lagrange polynomial” ]
•
The
n
th
degree interpolating polynomial for
f(x)
is then
P
n
(x) = y
0
L
0
(x) + y
1
L
1
(x) +
· · ·
+ y
n
L
n
(x) =
n
X
j=0
y
j
L
j
(x)
•
Check that
P
n
actually interpolates the function at
x = x
k
:
P
n
(x
k
) = y
0
L
0
(x
k
) + y
1
L
1
(x
k
) +
· · ·
+ y
n
L
n
(x
k
) = y
k
L
k
(x
k
) = y
k
Uniqueness:
When expanded, this
P
n
(x)
must be identical
to the “standard” form:
P
n
(x) =
c
0
+ c
1
x + c
2
x
2
+
· · ·
c
n
x
n
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
19/43
Accuracy / Cost of Lagrange Interpolation
The
::::::::::::::::::::::::
Lagrange form of the interpolating polynomial is a big improvement over the
:::::::::::::::::::::::
standard form (using the Vandermonde matrix method):
•
Accuracy:
⇤
No linear solves involved
⇤
No systematic problem with subtractive cancellation error
•
Cost:
⇤
Constructing each of the
(n + 1)
Lagrange polynomials requires
n
multipli-
cations and
n
divisions
⇤
Evaluating the interpolating polynomial requires
n + 1
multiplications.
⇤
Total cost:
2n(n + 1) + (n + 1) = 2n
2
+ 3n + 1 = O(n
2
)
=
)
For large
n
, this is a huge improvement over
O(n
3
)
!
•
Simplicity:
Lagrange code is very easy to implement.
4a. Polynomial interpolation
MACM 316
20/43
Example 3: Lagrange Interpolation
Return to the data from Example 1
and find the quadratic interpolating
polynomial
P
2
(x)
on
[1, 2.5]
x
i
y
i
.
.
.
.
.
.
1.0
3.6788
2.0
5.4134
2.5
5.1303
.
.
.
.
.
.
•
First set up the Lagrange polynomials:
At
x
0
= 1
:
L
0
(x) =
(x
-
2)
(1
-
2)
(x
-
2.5)
(1
-
2.5)
=
(x
-
2)(x
-
2.5)
1.5
[ leave unexpanded ]
At
x
1
= 2
:
L
1
(x) =
(x
-
1)
(2
-
1)
(x
-
2.5)
(2
-
2.5)
=
(x
-
1)(x
-
2.5)
-
0.5
At
x
2
= 2.5
:
L
2
(x) =
(x
-
1)
(2.5
-
1)
(x
-
2)
(2.5
-
2)
=
(x
-
1)(x
-
2)
0.75
•
Then construct the interpolating polynomial:
[ same result as slide 10 ]
P
2
(x) =
3.6788
(x
-
2)(x
-
2.5)
1.5
+ 5.4134
(x
-
1)(x
-
2.5)
-
0.5
+ 5.1303
(x
-
1)(x
-
2)
0.75
=
2.4525 (x
-
2)(x
-
2.5)
-
10.8268 (x
-
1)(x
-
2.5) + 6.8404(x
-
1)(x
-
2)
=
)
P
2
(1.75)
⇡
2.4525
·
0.25
·
0.75 + 10.8268
·
0.75
·
0.75
-
6.8404
·
0.75
·
0.25
⇡
5.2674
Check:
Never,
EVER
expand the Lagrange form (except to compare)
P
2
(x) =
-
1.1235 + 6.3362x
-
1.5339x
2
. . . which
MUST BE
the same polynomial from slide 10 –
UNIQUE!
a
i
3
terms
3.67886
2452.5
5.4134
115
2.5
5.1303
1
1
3
618,8
0.2511
0.751
54
0.7571
0.7s
5.637
0.75110.25
I
5.26735
4a. Polynomial interpolation
MACM 316
21/43
Example 3 (continued)
(a) Now, add one more point
(3.0, 4.4808)
and find the
Lagrange interpolating polynomial on the interval
[1.0, 3.0]
containing 4 data points.
(b) Use your result to estimate the original function at
x = 2.8
.
x
i
y
i
1.0
3.6788
2.0
5.4134
2.5
5.1303
3.0
4.4808
Solution:
(a) There are four points (
n = 3
), so the polynomial must be
:::::::::
cubic:
P
3
(x) =
3.6788(x
-
2)(x
-
2.5)(x
-
3)
(1
-
2)(1
-
2.5)(1
-
3)
+
5.4134(x
-
1)(x
-
2.5)(x
-
3)
(2
-
1)(2
-
2.5)(2
-
3)
+
5.1303(x
-
1)(x
-
2)(x
-
3)
(2.5
-
1)(2.5
-
2)(2.5
-
3)
+
4.4808(x
-
1)(x
-
2)(x
-
2.5)
(3
-
1)(3
-
2)(3
-
2.5)
(b) Leave in unexpanded form when evaluating – it’s much simpler!!
P
3
(2.8) =
3.6788(0.8)(0.3)(
-
0.2)
(
-
1)(
-
1.5)(
-
2)
+
5.4134(1.8)(0.3)(
-
0.2)
(1)(
-
0.5)(
-
1)
+
5.1303(1.8)(0.8)(
-
0.2)
(1.5)(0.5)(
-
0.5)
+
4.4808(1.8)(0.8)(0.3)
(2)(1)(0.5)
= 3.6788(0.016)
-
5.4134(0.216) + 5.1303(0.768) + 4.4808(0.432)
⇡
4.7653
367884,3
35
3
5.4134
Y
2
more
terms
3
6788
19
5
3
5.4134
is
t
2
more
3
6788
co
016145.413410.2161
t
2
move
4.7653
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
22/43
Example 3 (continued): Estimating Error
Recall Theorem on slide 9:
Error estimate is
f(x) = P
n
(x) + E
n
(x)
with
E
n
(x) =
f
(n+1)
(c)
(n + 1)!
(x
-
x
0
)(x
-
x
1
)
· · ·
(x
-
x
n
)
for some
c
2
[x
0
, x
n
]
In Example 3 (
n = 3
):
•
We approximated
f(2.8)
by
P
3
(2.8) = 4.7653
•
Suppose you know that
|
f
(4)
(x)
|
6
0.12
for all
x
=
)
this information is not usually available!!
•
How big can the absolute error
|
f(2.8)
-
4.7653
|
be?
x
i
y
i
1.0
3.6788
2.0
5.4134
2.5
5.1303
3.0
4.4808
E
3
(2.8)
=
1
4!
f
(4)
(c)(2.8
-
1)(2.8
-
2)(2.8
-
2.5)(2.8
-
3)
for some
c
2
[1, 3]
6
0.12
24
(1.8)(0.8)(0.3)(0.2)
= 4.32
⇥
10
-
4
=
)
PRETTY GOOD!
•
We’ll have a more reliable method for estimating error soon!
IE
312.871
4,1414741
2.8
1712.8
2312.8
2
5112.8311
for
some
ce
1.3
4
t
1
8110.81103110.21
I
4.32
10
4
20.0001
rel
error
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
23/43
Algorithm for Lagrange Interpolation
xlist = [0.0, 0.5, 1.0, 2.0, 2.5, 3.0, ...];
ylist = [0.0, 1.5163, 3.6788, 5.4134, 5.1303, 4.4808, ...];
np1 = length(xlist);
x = 1.75;
psum = 0.0;
for k = 1 : np1,
psum = psum + ylist(k)
*
eval_Lk(xlist, k, x);
end
psum
% - - - - - - - - - - - - - - - - - - - - - - - - - - - -
% EVAL_LK: Evaluate k’th Lagrange polynomial at x
function Lk = eval_Lk(xlist, k, x)
np1 = length(xlist);
Lk = 1.0;
for i = 1 : np1,
if i ˜= k,
Lk = Lk
*
(x-xlist(i)) / (xlist(k)-xlist(i));
end
end
end
A super-simple implementation that could be much more efficient!
lagsimple.m
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
24/43
Pros and Cons of the Lagrange Form
Advantages:
•
Lagrange polynomials can be written down easily for hand computation.
•
The algorithm is very simple to implement.
•
Cost is
O(n
2
)
– much better than
O(n
3
)
for the Vandermonde matrix method.
Disadvantages:
•
The error estimate from the Theorem (slide 9) is
not usable
:
E
n
(x) =
f
(n+1)
(c)
(n + 1)!
(x
-
x
0
)(x
-
x
1
)
· · ·
(x
-
x
n
)
. . . because
f
(n+1)
(c)
is unknown!
[ we only know a few values of
f(x)
]
•
If we add another point
(x
n+1
, y
n+1
)
and increase degree, we have to start over
again. No previously-calculated
L
k
(x)
for
k
6
n
can be reused.
Question:
Is there another method that avoids these two shortcomings?
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
25/43
Newton Divided Difference Form
Text
§3.3
•
Idea:
Build up
P
n
(x)
one point and one degree at a time
•
Start at degree
0
to interpolate the single point
(x
0
, y
0
)
:
P
0
(x) = a
0
[ constant function ]
•
Increase to degree
1
by adding a term
a
1
(x
-
x
0
)
P
1
(x) = y
0
+ a
1
(x
-
x
0
)
[ linear ]
Then solve for
a
1
so that
P
1
(x)
interpolates the point
(x
1
, y
1
)
•
This process extends easily to higher degrees:
P
2
(x) = a
0
+ a
1
(x
-
x
0
) +
a
2
(x
-
x
0
)(x
-
x
1
)
[ quadratic ]
.
.
.
P
n
(x) = a
0
+ a
1
(x
-
x
0
) + a
2
(x
-
x
0
)(x
-
x
1
)
+
· · ·
+
a
n
(x
-
x
0
)(x
-
x
1
)
· · ·
(x
-
x
n
-
1
)
=
n
X
k=0
a
k
"
k
-
1
Y
i=0
(x
-
x
i
)
#
Newton basis (degree
k
):
N
k
(x) =
k
-
1
Y
i=0
(x
-
x
i
)
(
⇤
)
•
In each step, determine
a
k
by substituting
P
n
(x
k
) = y
k
:
y
0
= P
0
(x
0
) = a
0
=
)
a
0
= y
0
y
1
= P
1
(x
1
) = a
0
+
a
1
(x
1
-
x
0
)
=
)
a
1
=
y
1
-
y
0
x
1
-
x
0
[ first divided difference ]
x2
ya
xo
yo
Ao
yo
a
1
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
26/43
•
The next divided difference comes from
y
2
= P
2
(x
2
) = a
0
+ a
1
(x
2
-
x
0
) +
a
2
(x
2
-
x
0
)(x
2
-
x
1
)
=
)
a
2
=
⇣
y
2
-
y
1
x
2
-
x
1
-
y
1
-
y
0
x
1
-
x
0
⌘
x
2
-
x
0
[ second divided difference ]
•
We now introduce some notation and give a
::::::::::::::::::::::::::::::::
recursive definition:
Step
0
:
f[x
i
] = y
i
Step
1
:
f[x
i
, x
i+1
] =
f[x
i+1
]
-
f[x
i
]
x
i+1
-
x
i
Step
2
:
f[x
i
, x
i+1
, x
i+2
] =
f[x
i+1
, x
i+2
]
-
f[x
i
, x
i+1
]
x
i+2
-
x
i
.
.
.
Step
s
:
f[x
i
, ... , x
i+s
] =
f[x
i+1
, ... , x
i+s
]
-
f[x
i
, ... , x
i+s
-
1
]
x
i+s
-
x
i
•
Then the coefficient
a
k
in equation
(
⇤
)
is just
a
k
= f[x
0
, x
1
, ... , x
k
]
(
i = 0
,
s = k
)
and is called the
k
th
divided difference
as
4
FY
XzNo
confusing
unclear
BUT
can
be
computed
EASILY
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
27/43
Newton Divided Difference Table
•
The notation on the previous slide is a bit confusing
•
These calculations are easily organized in tabular form:
Step 0
Step 1
Step 2
Step 3
Step
n
x
i
y
i
= a
0
a
1
a
2
a
3
...
a
n
x
0
f[x
0
]
f[x
0
, x
1
]
x
1
f[x
1
]
f[x
0
, x
1
, x
2
]
f[x
1
, x
2
]
f[x
0
, x
1
, x
2
, x
3
]
x
2
f[x
2
]
f[x
1
, x
2
, x
3
]
.
.
.
f[x
2
, x
3
]
f[x
1
, x
2
, x
3
, x
4
]
x
3
f[x
3
]
f[x
2
, x
3
, x
4
]
.
.
.
f[x
3
, x
4
]
f[x
2
, x
3
, x
4
, x
5
]
x
4
f[x
4
]
f[x
3
, x
4
, x
5
]
f[x
4
, x
5
]
x
5
f[x
5
]
.
.
.
.
.
.
fix
xD
this
I
fix
xs xx
this.gg
xz.xs
Tiggy
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
28/43
Example 4: Newton Divided Differences
Use data from Example 1 to construct the divided difference table, taking the subset
of points
x
i
2
{
0.0, 0.5, 1.0, 2.0, 2.5, 3.0
}
:
x
i
y
i
= a
0
a
1
a
2
a
3
a
4
0.0
0.0000
3.0326
0.5
1.5163
1.2924
4.3250
-
1.5097
1.0
3.6788
-
1.7269
0.6425
1.7346
0.0965
2.0
5.4134
-
1.5339
0.1216
-
0.5662
0.4005
2.5
5.1303
-
0.7328
-
1.2990
3.0
4.4808
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
29/43
Cubic polynomial starting at
x = 0.0
:
P
0
(x) =
0
P
1
(x) =
0 + 3.0326(x
-
0.0)
P
2
(x) =
3.0326(x
-
0.0) + 1.2924(x
-
0.0)(x
-
0.5)
P
3
(x) =
P
2
(x)
-
1.5097(x
-
0.0)(x
-
0.5)(x
-
1.0)
etc.
Read coefficients
a
k
directly off diagonals!
Quadratic polynomial starting at
x = 1.0
:
P
0
(x) =
3.6788
P
1
(x) =
3.6788 + 1.7346(x
-
1.0)
P
2
(x) =
P
1
(x)
-
1.5339(x
-
1.0)(x
-
2.0)
etc.
This is same polynomial as in earlier
examples! (check by expanding)
3 1
pts
x
0.0.5.112
3
6788
i
m
Read
coefficients
off
diagonals
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
30/43
Cost of Newton Interpolation
x
i
y
i
= a
0
a
1
a
2
a
3
a
4
a
5
x
0
= 0.0
0.0000
3.0326
x
1
= 0.5
1.5163
1.2924
4.3250
-
1.5097
x
2
= 1.0
3.6788
-
1.7269
0.6425
1.7346
0.0965
-
0.1736
x
3
= 2.0
5.4134
-
1.5339
0.1216
-
0.5662
0.4005
x
4
= 2.5
5.1303
-
0.7328
-
1.2990
x
5
= 3.0
4.4808
(n=5)
•
Cost of computing the Newton divided difference table:
Column 0:
no work (just fill in the data values)
Column 1:
n
divisions
Column 2:
n
-
1
divisions
.
.
.
Column
n
:
1 division
All columns:
1
2
n(n + 1)
•
Cost of evaluating
P
n
(x) =
n
X
k=0
a
k
k
-
1
Y
i=0
(x
-
x
i
)
=
)
also
1
2
n(n + 1)
•
Total cost:
n
2
+ n
operations
•
Compare with Lagrange interpolation:
2n
2
+ 3n + 1
•
Both are
O(n
2
)
BUT
Newton polynomial is about
half the cost!
no
work
a
p
division
In
Intl
n't
n
op's
O
n
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
31/43
Advantages of Newton Form
•
We can easily determine the interpolating polynomial for any desired interval
and degree once the divided difference table is computed.
•
Total cost is half that of the Lagrange polynomial.
•
Big Advantage:
The divided difference table provides an
easy error estimate
since the (absolute) error can be written as
E
n
(x) = f[x
0
, x
1
, ... , x
n+1
]
n
Y
i=0
(x
-
x
i
)
[ see argument from Exercise 22, page 133 ]
and then matching to the Theorem on slide 9 gives
f[x
0
, x
1
, ... , x
n+1
]
⇡
f
(n+1)
(c)
(n + 1)!
for some
c
2
[x
0
, x
n+1
]
•
In practice, error can estimated by adding one more data point and computing
one extra diagonal in the table – simple and cheap!
[ Extra O(n) cost ]
Exercise:
Derive the quadratic polynomial interpolant
P
2
(x)
in Newton form, for
Example 1 on the interval
[1, 2.5]
. Then use it to approximate
f(1.75)
and estimate
the error.
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
32/43
Solution:
We already found
P
2
(x)
in Example 4 (slide 28) but here it is again:
P
2
(x) =
3.6788 + 1.7346(x
-
1)
-
1.5339(x
-
1)(x
-
2)
=
)
P
2
(1.75) =
3.6788 + 1.7346(0.75)
-
1.5339(0.75)(
-
0.25)
⇡
5.2673
•
To estimate error, add to the table the extra di-
agonal corresponding to
x = 3.0
•
The error term has the form
E
2
(x) = f[1, 2, 2.5, 3](x
-
1)(x
-
2)(x
-
2.5)
where
a
3
= f[1, 2, 2.5, 3] =
0.4005
and so the
absolute error is
|
E
2
(1.75)
|
⇡
0.4005(0.75)(
-
0.25)(
-
0.75)
⇡
0.0563
x
i
y
i
= a
0
a
1
a
2
a
3
1.0
3.6788
1.7346
2.0
5.4134
-
1.5339
-
0.5662
0.4005
2.5
5.1303
-
0.7328
-
1.2990
3.0
4.4808
•
Relative error
⇡
0.0563
/
5.2673
⇡
0.0107
or about
1
% – pretty good!
•
The error estimate
|
E
2
|
6
0.0563
also gives us an
interval
containing
f(1.75)
:
P
2
(1.75)
-
0.0563
6
f(1.75)
6
P
2
(1.75) + 0.0563
5.2673
-
0.0563
6
f(1.75)
6
5.2673 + 0.0563
=
)
5.2110
6
f(1.75)
6
5.3236
•
Extra diagonal also provides a higher degree polynomial approximation:
P
3
(x) =
3.6788 + 1.7346(x
-
1)
-
1.5339(x
-
1)(x
-
2) + 0.4005(x
-
1)(x
-
2)(x
-
2.5)
3
678841.73461
1
1.53391 171
21
3.678841734610.751
1.533910
15,114.25
o
3280.4006
0.4006
ELI
7511
10400610.7511
0.2511
0.75
70.0563
0.05612
11.75
P
1.75
10.0523
5.2673
0.0563
Ef
1.751
5.2673
0.0563
Bonus
Paix
10.4006
x
1
x
2
x
2.5
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
33/43
Algorithm for Newton Interpolation
function a = newttable(x, y)
Canvas:
newttable.m
n = length(x);
a = zeros(n, n+1);
a(:,1) = x(:);
a(:,2) = y(:);
for j = 3 : n+1,
% build columns one at a time
for i = 1 : n-j+2,
a(i,j) = (a(i+1,j-1) - a(i,j-1)) / (a(i+j-2,1) - a(i,1));
end
end
>> xlist = [
1.0,
2.0,
2.5,
3.0];
>> ylist = [3.6788, 5.4134, 5.1303, 4.4808];
>> a = newttable(xlist, ylist)
a =
1.0000
3.6788
1.7346
-1.5339
0.4005
2.0000
5.4134
-0.5662
-0.7328
0
2.5000
5.1303
-1.2990
0
0
3.0000
4.4808
0
0
0
>> x0 = 1.75;
>> pval = a(1,2) + a(1,3)
*
(x0-a(1,1)) + a(1,4)
*
(x0-a(1,1))
*
(x0-a(2,1)) ...
+ a(1,5)
*
(x0-a(1,1))
*
(x0-a(2,1))
*
(x0-a(3,1))
pval = 5.3237
OR
>> pval = newtpoly(a, x0) % helper function!
Canvas:
newtpoly.m
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
34/43
A “Who’s Who” of Interpolation
Several “giants” of approximation and interpolation lived around the same time:
Isaac
Brook
Alexandre-Th´eophile
Joseph-Louis
Newton
Taylor
Vandermonde
Lagrange
(1643–1727)
(1685–1731)
(1735–1796)
(1736–1813)
There are many others who have contributed to polynomial interpolation:
Gauss, Lobatto, Chebyshev, Lebesgue, Weierstrass, Runge, Hermite, . . .
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
35/43
Example 5: Three Interpolants
Interpolate the four point values in the table below
x
-
3
-
2
0
1
5
f(x)
5
3
4
1
-
1
with a cubic polynomial that is written in:
(a) standard / monomial form
(b) Lagrange form
(c) Newton form
Then approximate
f(
-
1)
and plot the polynomial interpolant.
Note:
The last point is included only for error estimation (in Newton form)
Vandermonde
x
5
Work
out
the
details
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
36/43
Example 5: Standard form
x
-
3
-
2
0
1
f(x)
5
3
4
1
The Vandermonde system
Vc = y
is
2
6
6
6
4
1
-
3
9
-
27
1
-
2
4
-
8
1
0
0
0
1
1
1
1
3
7
7
7
5
·
2
6
6
6
4
c
0
c
1
c
2
c
3
3
7
7
7
5
=
2
6
6
6
4
5
3
4
1
3
7
7
7
5
=
)
c =
2
6
6
6
6
4
4
-
5
6
-
5
3
-
1
2
3
7
7
7
7
5
Then the cubic interpolant in standard form is
P
3
(x) = 4
-
5
6
x
-
5
3
x
2
-
1
2
x
3
=
)
P
3
(
-
1) = 4 +
5
6
-
5
3
+
1
2
=
11
3
⇡
3.6667
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
37/43
Example 5: Lagrange form
x
-
3
-
2
0
1
f(x)
5
3
4
1
The cubic interpolant in Lagrange form is:
P
3
(x) = 5
(x + 2)(x
-
0)(x
-
1)
(
-
3 + 2)(
-
3
-
0)(
-
3
-
1)
+ 3
(x + 3)(x
-
0)(x
-
1)
(
-
2 + 3)(
-
2
-
0)(
-
2
-
1)
+ 4
(x + 3)(x + 2)(x
-
1)
(0 + 3)(0 + 2)(0
-
1)
+ 1
(x + 3)(x + 2)(x
-
0)
(1 + 3)(1 + 2)(1
-
0)
=
-
5
12
(x + 2)x(x
-
1) +
1
2
(x + 3)x(x
-
1)
-
2
3
(x + 3)(x + 2)(x
-
1) +
1
12
(x + 3)(x + 2)x
Must have the same result
P
3
(
-
1) =
11
3
Exercise:
Expand and show that this is identical to the standard form.
1
in
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
38/43
Example 5: Newton form
x
-
3
-
2
0
1
5
f(x)
5
3
4
1
-
1
The Newton divided difference table is:
x
i
y
i
= a
0
a
1
a
2
a
3
a
4
-
3
5
-
2
-
2
3
5
6
1
2
-
1
2
0
4
-
7
6
31
336
-
3
5
21
1
1
1
2
-
1
2
5
-
1
Then the cubic interpolant in Newton form is
P
3
(x) = 5
-
2(x + 3) +
5
6
(x + 3)(x + 2)
-
1
2
(x + 3)(x + 2)x
P
3
(
-
1) =
11
3
⇡
3.6667
AND
the error term is
E
3
(x) =
31
336
(x + 3)(x + 2)x(x
-
1)
=
)
E
3
(
-
1) =
31
336
(
-
1 + 3)(
-
1 + 2)(
-
1)(
-
1
-
1)
⇡
0.3690
[ relative error
⇡
10% ]
abs
error
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
39/43
Summary of Polynomial Interpolation
The
n
th
degree interpolating polynomial
P
n
(x)
on
n + 1
points
(x
0
, y
0
), (x
1
, y
1
), ... (x
n
, y
n
)
with
x
0
<
x
1
<
· · ·
<
x
n
distinct &
ordered
satisfies the following properties:
•
it
::::::::::::::::::::
interpolates:
P
n
(x
i
) = y
i
for all
i = 0, 1, ... , n
•
it’s
::::::::::::::::::::::::::::::
maximally smooth: continuous and infinitely differentiable (polynomials)
•
it’s
:::::::::::
unique (Theorem on slide 9): standard, Lagrange and Newton forms are all
identical when expanded
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
40/43
Comparing Basis Polynomials
M
ATLAB
:
plotbases.m
Plots basis polynomials for
n = 2, 3, 4, 5, 6
equally-spaced points on
[0, 1]
:
Middix
Liam
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
41/43
Comparing Accuracy & Cost
1. “Standard” or monomial form:
uses the Vandermonde matrix
V
to obtain coeffi-
cients
c
i
for the expanded polynomial
P
n
(x) =
n
X
j=0
c
j
Mj(x)
with
M
j
(x) = x
j
“standard” or
monomial basis
•
matrix
V
is extremely ill-conditioned, so coefficients can be very inaccurate.
•
O(n
3
)
cost is expensive compared with other methods.
2. Lagrange form:
written in the form
P
n
(x) =
n
X
j=0
y
j
L
j
(x)
with
L
j
(x) =
n
Y
i=0
i
6
=j
x
-
x
i
x
j
-
x
i
Lagrange
basis
•
avoids subtractive cancellation errors.
•
O(n
2
)
total cost is a big improvement.
3. Newton form:
written in terms of divided difference table entries as
P
n
(x) =
n
X
j=0
f[x
0
, x
1
, ... , x
j
] N
j
(x)
with
N
j
(x) =
j
-
1
Y
i=0
(x
-
x
i
)
Newton
basis
•
also avoids subtractive cancellation errors.
•
cost is also
O(n
2
)
, but roughly half of that for the Lagrange form.
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
42/43
Afterword: Polynomial Wiggle
As the polynomial order increases, “wiggles” or oscillations tend to arise
[ These are very different from the smaller random-looking wiggles we saw on slide 14, arising from round-off error ]
Example 5:
Fit various order polynomials to data points that experience a “jump”:
•
Resulting interpolant doesn’t
necessarily approximate the
data well.
•
Wiggles typically increase in
number and amplitude as
n
increases.
•
High-order polynomials are
particularly terrible at
approximating functions with
jumps or rapid changes!
M
ATLAB
:
demoWiggle.m
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
4a. Polynomial interpolation
MACM 316
43/43
What’s Next?
•
Interpolating data with a single polynomial only works well when the data is
:::::::::::::::::::::
very smooth – that is, if it’s already close to polynomial!
•
Real data aren’t smooth – measurements are inherently random!
Conclude:
Smooth interpolation using a single high-order polynomial will almost
always generate oscillatory wiggles
=
)
smooth polynomial interpolation is not really a good idea
except in special cases
Question:
Can we relax the smoothness somewhat and still interpolate?
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help