Description
We investigate interpolation of the function
f : [−1, 1] 7→ R, f(x) = 1
1 + 30x
2
,
very similar to a famous example, the “Runge” function, used to exhibit – yes you guessed it – the
“Runge phenomenon”.
f = @(x) 1./(1 + 30*x.^2 );
M1 Polynomial interpolation with equidistant points, xk = −1+2k/N, k = 0, 1, 2, . . . N. Weights
for the barycentric formula are wk = (−1)k
N
k
!
.
M2 Polynomial interpolation with Chebyshev points of the second kind,
xk = cos(kπ/N), k = 0, 1, 2, . . . N.
Weights for the barycentric formula are wk = (−1)k
, except w0 =
1
2
, wN = (−1)N 1
2
.
M3 Mix and match: Use barycentric formula with equidistant points, but with Chebyshev
weights. In this case the interpolating function is no longer a polynomial, but a rational
function.
Your coding tasks.
a. You will need to write code to compute the weights for the barycentric formula. Matlab has a
built-in command nchoosek; You may use that, or compute binomial coefficients recursively,
N
0
= 1,
N
k+1
=
N
k
N−k
k+1 .
b. You will have to write code implementing the barycentric formula,
p(x) =
PN
m=0
wmfm
P
x−xm
N
m=0
wm
x−xm
.
This can be done to varying degree of efficiency and flexibility. Remember, that Matlab likes
vectors; for loops not so much. Ideally, your function can take an array of evaluation points
as inputs, and return the values at all of those points simultaneously.
One detail you will have to pay attention to is the case x − xk = 0, i.e., an evaluation point
coinciding with an interpolation point. In that case your function should simply assign the
given value fk to p(xk). If your barycentric code allows vectors as input arguments, then the
Matlab function find is useful for this task. Otherwise, a simple if statement will do.
mr
