function s=spline2(x0,y0,y21,y2n,x)
%三次样条插值，应用II型边界条件
%s=spline2(x0,y0,y21,y2n,x)
%x0,y0 are existed points,x are insert points,y21,y2n are the second
%dirivitive numbers given.
n=length(x0);
km=length(x);
a(1)=-0.5;
b(1)=3*(y0(2)-y0(1))/(2*(x0(2)-x0(1)));
for j=1:(n-1)
h(j)=x0(j+1)-x0(j);
end
for j=2:(n-1)
alpha(j)=h(j-1)/(h(j-1)+h(j));
beta(j)=3*((1-alpha(j))*(y0(j)-y0(j-1))/h(j-1)+alpha(j)*(y0(j+1)-y0(j))/h(j));
a(j)=-alpha(j)/(2+(1-alpha(j))*a(j-1));
b(j)=(beta(j)-(1-alpha(j))*b(j-1))/(2+(1-alpha(j))*a(j-1));
end
m(n)=(3*(y0(n)-y0(n-1))/h(n-1)+y2n*h(n-1)/2-b(n-1))/(2+a(n-1));
for j=(n-1):-1:1
m(j)=a(j)*m(j+1)+b(j);
end
for k=1:km
for j=1:(n-1)
if((x(k)>=x0(j))&(x(k)<x0(j+1)))
l(k)=j;
end
end
end
for k=1:km
sum=(3*(x0(l(k)+1)-x(k))^2/h(l(k))^2-2*(x0(l(k)+1)-x(k))^3/h(l(k))^3)*y0(l(k));
sum=sum+(3*(x(k)-x0(l(k)))^2/h(l(k))^2-2*(x(k)-x0(l(k)))^3/h(l(k))^3)*y0(l(k)+1);
sum=sum+h(l(k))*((x0(l(k)+1)-x(k))^2/h(l(k))^2-(x0(l(k)+1)-x(k))^3/h(l(k))^3)*m(l(k));
s(k)=sum-h(l(k))*((x(k)-x0(l(k)))^2/h(l(k))^2-(x(k)-x0(l(k)))^3/h(l(k))^3)*m(l(k)+1);
end
