Lab03．三次样条插值实验
【实验目的和要求】
1．使学生深入理解三次样条插值法，深入进行程序设计能力训练；
2．对第一与第二种边界条件，按三弯矩法，通过用Matlab语言设计计算三次样条插值的程序，以提高学生程序设计的能力。 【实验内容】
1．根据Matlab语言特点，描述三次样条插值法。
2．对第一与第二种边界条件，按三弯矩法，用Matlab语言设计计算三次样条插值的程序。 3．已知数据： xi 0.2 0.4 0.6 0.8 1.0 f(xi) 0.9798652 0.9177710 0.8080348 0.6386093 0.3843735 对于： (1) 自然边界条件S??(0.2)?S??(1.0)?0；
, S?(1.0)?1.55741. (2) 第一种边界条件S?(0.2)?0.20271输出用追赶法解出的弯矩向量(M1,M2,?M5)和S(0.2?0.1i) (i=0，1，…，8）的值，并画出y?S(x)的图形。
4．完成教材P45例8的计算，并将计算结果与Langrage插值法计算的结果
进行比较，由此说明三次样条插值的优越性。 【实验仪器与软件】
1．CPU主频在1GHz以上，内存在128Mb以上的PC；
2．Matlab 6.0及以上版本。 实验讲评：
实验成绩：
评阅教师： 2011 年 月 日
1
Lab03．三次样条插值实验
一、算法描述
1．根据Matlab语言特点，描述三次样条插值法. 答：
S(x) 在 [xj, xj+1]（j=1,2,?,n-1）上是三
次多项式，于是S\在[xj, xj+1] 上是一次多项式， 如果S\在 [xj,xj+1]（j=1,2,?,n-1）两端点上的值 已知，设S\j)=Mj，S\j+1)=Mj+1，则S\的表达 式为：?
??， 其中hj ＝
xj+1-xj，对S\进行两次积分，则得到1 个具有2
二、程序设计
2．对第一与第二种边界条件，按三弯矩法，用Matlab语言设计计算三次样条插值的程序。 答：
function s=selfspline(x0,y0,x) %x为待求值，按第二种边界条件求解 %x0,y0为插值节点，df为一介或二介导数 if (length(x0)==length(y0)) n=length(x0); else
disp('x0.y0维数不相等!'); return; end
%对x0进行由小到大排序
2
for j=1:n for i=1:n-j
if(x0(i)>x0(i+1)) t=x0(i);
x0(i)=x0(i+1); x0(i+1)=t; end end end
%对x0进行由小到大排序
%判断待求值x属于哪个区间 for m=1:n
if x0(m)<=x&x<=x0(m+1) m=m; break end end
%判断待求值x属于哪个区间
%求每个区间的长度 for t=1:n-1
h(t)=x0(t+1)-x0(t);
3
end
%求每个区间的长度
%求niu,lamda,d b=ones(1,n)*2; niu(1)=1; lamda(1)=1; for j=2:n-1
niu(j)=h(j-1)/(h(j-1)+h(j)); lamda(j)=h(j)/(h(j-1)+h(j)); end
d(1)=6/h(1)*((y0(2)-y0(1))/h(1)-y0(1)); d(n)=6/h(n-1)*(y0(n)-(y0(n)-y0(n-1))/h(n-1)); for i=2:n-1
d(i)=6*((y0(i+1)-y0(i))/h(i)-(y0(i)-y0(i-1))/h(i-1))/h(i-1)+h(i); end
%求niu,lamda,d
%消元
for k=1:n-1
niu(k)=-b(k)*(niu(k)/b(k))+niu(k);
4
b(k+1)=-lamda(k)*(niu(k)/b(k))+b(k+1); d(k+1)=-d(k)*niu(k)/b(k)+d(k); end %消元 %回代
M(n)=d(n)/b(n);
for a=n-1:-1:1
M(a)=(d(a)-lamda(a)*M(a+1))/b(a); end %回代 %求s(x)
s=M(m)*(x0(m+1)-x)^3/6*h(m)+M(m+1)*(x-x0(m))^3/6*h(m)...
+(y0(m)-M(m)*h(m)^2/6)*(x0(m+1)-x)/h(m)+(y0(m+1)-M(m+1)... *h(m)^2/6)*(x-x0(m)/h(j)) %求s(x)
三、实例求解
四、算法总结
5