Chapter 3. Interpolation of Functions

Vladimir Rovenski, Haifa

> restart:

3.1 Polynomial Interpolation of Functions

Example 1 . Plot over the range x=1..9 and look at its strange behaviour outside of the segment

containing the given nodes.

> xx:=[2,3,4,5,6,7,8]: yy:=[32,34,36,34,39,40,37]:

> f:=interp(xx, yy, x); INTERPOL:=plot(f, x=1..9):

> points1:=[seq([xx[j], yy[j]], j=1..nops(xx))]:

> POINTS2:=plot(points1, style=point, symbol=circle):

> plots[display]([POINTS2, INTERPOL]);

Then we plot the graph of on the segment = 1.8..8.2 and look at its oscillations.

Example 2 . [ K. Runge , 1901 ]

> f:=x -> 1/(1+25*x^2): x:=i -> -1+2*i/n: y:=i -> f(x(i)):

> n:=6: g:=interp([seq(x(i), i=0..n)], [seq(y(i),i=0..n)],x);

> plot([f(x), g], x=-1..1);

Example 3. [ S.Bernstein , 1912 ]

> f:=x -> abs(x): x:=i -> -1+2*i/n: y:=i -> f(x(i)):n:=6: g:=interp([seq(x(i), i=0..n)], [seq(y(i),i=0..n)],x);

> plot([f(x), g], x=-1..1);

> xx:=[2,3,4,5,6,7,8]: zz:=[0,0,0,0,0,0,0]: n:=nops(xx)-1:

> yy:=[32,34,36,34,39,40,37]: Points1:=plot([xx, yy], style=point,symbol=circle): delta:=(i,j)->if i=j then 1 else 0 fi:

> for k from 1 to n+1 do L[k]:=unapply(interp(xx,[seq(delta(k,j),j=1..n+1)],x),x); Lx[k]:=unapply(diff(L[k](x),x),x) od:

> HERM:=sum(yy[i]*(1-2*(x-xx[i])*Lx[i](xx[i]))*L[i](x)^2,i=1..n+1)+ sum(zz[i]*(x-xx[i])*L[i](x)^2, i=1..n+1):

> HERMITE:=plot(HERM, x=1.9..8.1):

> plots[display]([Points1, HERMITE],view=[1.9..8.1,25..41]);

3.2 Spline Interpolation of Functions

For example:

> readlib(spline): cubic1:=spline([0,1,2,3], [0,1,4,3], x);

Examples.

1 . Plot the " triangular waves " using a linear spline and then smooth them by a cubic spline.

> v:=spline([seq(i, i=0..9)], [seq((-1)^i, i=0..9)], x, linear):

> c:=spline([seq(i, i=0..9)], [seq((-1)^i, i=0..9)], x):

> plot([v, c], x=0..9, scaling=constrained);

3 . Compare the behavior of the Lagrange interpolation polynomial and polynomial splines derived by the "measurements

of air temperature over 10 days" given in Section 1.1.

> Days:=[12,13,14,15,16,17,18,19,20,21]; TT:=[15,17,17.5,19,20,19.5,18,17,17,19];

> f2:=spline(Days, TT, x, quadratic): f3:=spline(Days, TT, x):

> h:=unapply(interp(Days, TT, x), x):

> a:=plot([f2, f3], x=11.8..21.2): b:=plot(h, 11.8..21.2):

> c:=plot([seq([Days[i], TT[i]],i=1..10)], style=point, symbol=diamond):

> plots[display]([a, b, c], view=[11.7..21.2, 12..22]);

3.3 Construction of Curves Using Spline Functions

Example 1 . Try to tie the following knot , well known to sailors.

> tt:=[0,1/3,2/3,1,4/3,3/2,5/3,2,7/3,8/3,3.1]: xx:=[-6,-4,-2,0,4,4.8,4,0,-2,-4,-6]: yy:=[4,2,0,-2,-2,0,2,2,0,-2,-4]: zz:=[1,1,-1,1,-1,0,1,-1,1,-1,-1]:

> tt:=[0,1/3,2/3,1,4/3,3/2,5/3,2,7/3,8/3,3.1]: xx:=[-6,-4,-2,0,4,4.8,4,0,-2,-4,-6]: yy:=[4,2,0,-2,-2,0,2,2,0,-2,-4]: zz:=[1,1,-1,1,-1,0,1,-1,1,-1,-1]:

> readlib(spline): X:=spline(tt, xx, t): Y:=spline(tt, yy, t): Z:=spline(tt, zz, t):

> plots[tubeplot]({[X(t), Y(t), Z(t)], [-X(t), Y(t), -Z(t)]}, t=0..3, radius=.6, tubepoints=30, numpoints=190, orientation=[90, 26], color=[.8, .6, .4], projection=.5, style=patchnogrid, scaling=constrained, ambientlight=[.6, .6, .5], light=[75, 50, 1,.9,.7]);

To understand the key idea, plot the figure in the plane.

> XY:=plot([[X(t), Y(t), t=0..3], [-X(t), Y(t), t=0..3]], color=[red, blue]):

> Points:=plot([seq([op(i, xx), op(i, yy)], i=1..nops(xx)), seq([-op(i, xx), op(i, yy)], i=1..nops(xx))], style=point):

> plots[display](XY, Points, scaling=constrained);

>