\magnification=\magstep1
\centerline{Mathematics 171 -- Numerical Analysis 1}
\centerline{Lab Day 4 -- Cubic Spline Interpolation}
\centerline{November 5, 1998}
\bigskip
\noindent
{\bf Goals and Background}
\bigskip
The goals of today's lab are to introduce the Maple {\tt spline}
procedure for {\it free cubic splines} and another ``homegrown'' Maple
procedure for {\it clamped cubic splines}.  We will use
them to construct cubic spline interpolating functions.
Recall that on Tuesday we saw that given $x_0, \ldots, x_n$ and 
function values $f(x_0), \ldots, f(x_n)$ to be interpolated,  
the coefficients $a_j,b_j,c_j,d_j$
in the free cubic spline interpolating function 
$$S(x) = 
\cases{S_0(x) = a_0 + b_0(x-x_0) + c_0(x-x_0)^2 + d_0(x-x_0)^3,\ x_0\le x\le x_1&\cr
S_1(x) = a_1 + b_1(x-x_1) + c_1(x-x_1)^2 + d_1(x-x_1)^3,\ x_1\le x\le x_2&\cr
\qquad \vdots & \cr
S_{n-1}(x) = a_{n-1} + b_{n-1}(x-x_{n-1}) + c_{n-1}(x-x_{n-1})^2 + 
d_{n-1}(x-x_{n-1})^3,&\cr
\qquad \qquad \qquad x_{n-1}\le x\le x_n&\cr}$$
satisfy $a_j = f(x_j)$ for $0 \le j \le n - 1$ and
$$\eqalign{b_j &= {a_{j+1} - a_j\over h_j} - {h_j\over 3}(2c_j + c_{j+1})\cr
           d_j &= {c_{j+1} - c_j\over 3 h_j}\cr}$$
for all $j = 0,\ldots, n$.
The $c_j$ are found by solving the $(n-1)\times (n-1)$ 
system of linear equations:
$$\eqalign{
&\pmatrix{
2(h_0+h_1) & h_1 & 0 &\cdots &0&0\cr
h_1 & 2(h_1+h_2) & h_2  & \cdots &0 &0\cr
\vdots&  &  && \vdots & \vdots \cr
0 & \cdots & 0 & 0&h_{n-2} & 2(h_{n-2} + h_{n-1}) \cr}\cdot 
\pmatrix{c_1 \cr \vdots \cr c_{n-1} \cr}\cr
& = 
\pmatrix{ {3(a_2-a_1)\over h_1} - {3(a_1-a_0)\over h_0}\cr
\vdots\cr
{3(a_n-a_{n-1})\over h_{n-1}} - {3(a_{n-2}-a_{n-1})\over h_{n-1}} }\cr}
$$
where $h_j = x_{j+1} - x_j$ and $c_0 = c_n = 0$.  There
is a very similiar system of linear equations whose solutions
give the $c_j$ coefficients for clamped cubic splines.
\bigskip
\noindent
{\bf A New Maple Command}
\bigskip
Maple has a spline interpolation function that can be used to construct 
free cubic spline interpolating
functions.  It is invoked by first entering the command:
\medskip
\centerline{\tt readlib(spline);}
\medskip
\noindent
(This does something very close to what the {\tt with} command for reading
in external packages does.)
Then the format of the spline command is similar to the {\tt interp} command
from the last lab.  
\medskip
\centerline{\tt spline(xlist,ylist,var);}
\medskip
\noindent
The output from the spline command is the piecewise cubic definition of
the {\it natural cubic spline\/} function $S(x)$ -- in Maple terms it is an 
{\it expression, not a function\/}, though(!) You can assign this to a symbolic name, 
plot, use it to define a function, or whatever you like.
\bigskip
\noindent
{\bf Laboratory Questions}
\bigskip
\item{A)}  In class on Tuesday, we worked out an example where 
$x_0 = 1, x_1 = 1, x_2 = 3, x_3 = 4$ and $f(x_0) = a_0 = 5, f(x_1) = a_1 = 4, 
f(x_2) = a_2 = 0, f(x_3) = a_3 = 1$.  We found that the cubic spline was
given by 
$$\cases{
S_0(x) = 5 - {5\over 8}x - {3\over 8} x^3, & $0 \le x \le 1$\cr
S_1(x) = 4 - {7\over 4}(x-1) - {9\over 8}(x-1)^2 + {1\over 2} (x-1)^3, & $1\le x \le 3$\cr
S_2(x) = -{1\over 4}(x-3) + {15\over 8}(x-3)^2 - {5\over 8} (x-3)^3, & $3\le x \le 4$\cr}$$
Verify that $S$ satisfies the continuity, differentiability, and second-order 
differentiability conditions for cubic spline functions.  Plot $S(x)$ and the points 
$(0,5), (1,4), (3, 0), (4,1)$
together (see the Maple online help page for the {\tt plot} command for some
ideas about plotting sets of points and use the option {\tt style = point}).
\bigskip
\item{B)}  Use Maple's {\tt spline} command to do problems 3d and 4d in section 3.4
of Burden and Faires.  Also, plot the spline together with $f(x) = x\cos x - 2 x^2  + 3 x - 1$,
and discuss the goodness of the ``fit'' between the spline function and $f(x)$.
\bigskip
\item{C)}  Use Maple's {\tt spline} command to do problem 27 in section 3.4 of Burden
and Faires.
\bigskip
\item{D)}  Maple does not have a built-in procedure or a library
package for computing {\it clamped\/} cubic
splines.  However, you can get the source code for one by entering the command:
\medskip
\noindent
\centerline{\tt cp /home/fac/little/Classes/Numer-98/clampedspine.map clampedspline.map}
\medskip
\item{}
(The source listing is also available on the course homepage.)
Do this and
read the code into your Maple session.  The format for using the procedure is similar
to the one for {\tt spline}, except you also need to supply the values for 
the derivative of the spline function at the endpoints ({\tt df0, dfn} are
$f'(x_0)$ and $f'(x_n)$ respectively):
\medskip
\centerline{\tt clampedspline(xlist,ylist,df0,dfn,var);}
\medskip
\item{}
Use this procedure to solve problem 28 in section 3.4 of Burden and Faires.
\bigskip
\noindent
{\bf Assignment}
\bigskip
Individual lab write-ups due in class Thursday, November 12.
\bye
\bye