A function, cosine ,which will use a quadratic spline to compute approximations for the cosine function at any real number.

A function, cosine,which will use a quadratic spline to compute approximations for the cosine function at any real number.

It should compute the quadratic spline S(x) using evenly-spaced nodes x1, x2, ..., xn , where x(k+1) - xk = d for k =1,2,...,n-1 , and whose derivative is specified at x1. The numbers y1 , y2 ,..., yn are the corresponding function values to be interpolated. On each interval [xk ,x(k+1)]

S(x) = Sk (x) = ak + bk (x - xk) + ck (x - xk)^2; k = 1, 2, ..., n-1

After applying the interpolation conditions, one obtains

ak = yk ; k = 1, 2, ..., n-1 (1)

and

ck = ((y(k+1)-yk)/d^2) - bk/d (2)

After applying the smoothness conditions, one obtains

bk + b(k+1) = (2/d) (y(k+1)-yk); k = 1, 2, ..., n-2 (3)

If S'(x) = v1 (derivative is specified to be v1 ), then one obtains

b1 = v1 (4)

a) Use Eqs. (3) and (4) to set up a linear system to solve for vector b that contains the coefficients b1, b2, ..., b(n-1).

Write down an example system for n = 5 in matrix form.

Also, find det(A) for the coefficient matrix, A, for any number of nodes, n, to show that there is a unique solution for any n. Recall that det(A) is the product of the diagonal elements for any lower triangular matrix.

(quadspline.m)

%%%%%%%%%%%%%%%%%%%%%  compute spline coefficients here  %%%%%%%%%%%%%%%

format long e

%%%%%%%%%%%%%%%%%% call sineval function here %%%%%%%%%%%%%%%%%

X = 250;

sineval(X,a,b,c)

X = -100;

sineval(X,a,b,c)

%%%%%%%%%%%%%%%%%%%%%  plot spline  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%number of points on which to plot. n = number of nodes

nplot = (n-1)*19+1;

xplot = zeros(nplot,1);

yplot = zeros(nplot,1);

%spacing between plot points

nspace = (x(n)-x(1))/(nplot-1);

k = 0;

for i = 1:n-1

  for j = 1:19

    k = k+1;

    xplot(k) = x(i) + (j-1)*nspace;

    yplot(k) = a(i) + b(i)*(xplot(k) - x(i)) + c(i)*(xplot(k) - x(i))^2;

  end

end

xplot(nplot) = x(n);

yplot(nplot) = a(i) + b(i)*(x(n) - x(n-1)) + c(i)*(x(n) - x(n-1))^2;

plot(xplot,yplot)

figure

abserr = abs(yplot - sin(xplot));

plot(xplot,abserr)

(forsub.m)

function [y] = forsub(L,b)

%Forward-substitution

%accepts an nX1 vector b, an nXn lower triangular matrix L

%generates an nX1 solution vector y

n = size(b,1);

y = zeros(n,1);

y(1) = b(1);

for i = 2:n 

  y(i) = b(i);

  for j = 1:i-1 

    y(i) = y(i) - L(i,j)*y(j); 

  end 

end