发布时间 : 星期日 文章基于 - MATLAB - 的函数的插值方法更新完毕开始阅读c60b607f27284b73f2425021
运行后输出f(3)?s2(3)和f(8)?s4(8)的值如下
x2 = s2 = x 4= s4 = 3 24.6221 8 68.5581
例6.7.3和例6.7.2用的方法不同，所以计算的结果不同，但是两种方法计算的结果相近.
6.7.6 用MATLAB计算三次样条
例6.7.6 给出节点的数据如下： -2.46 -5.26 -1.87 0.05 1.65 2.69 4.56 7.89 10.31 分别求在下列条件下在插值点x1??0.02，x2?2.56处的压紧三次样条插值，并显示该样条函数的有关信息：
?(3.50)?29.16； ?(?1)?5，Sn（1）端点约束条件为Snx -1.00 -0.54 0.13 1.12 1.89 2.06 2.54 2.82 3.50 f(x) ?(3.50)?0. ?(?1)?0，Sn（2）端点约束条件为Sn解 （1）输入MATLAB程序
>> X=[-1.00 -0.54 0.13 1.12 1.89 2.06 2.54 2.82
3.50];
Y=[-2.46 -5.26 -1.87 0.05 1.65 2.69 4.56 7.89 10.31];
XI=[-0.02 2.56]; YI= spline (X, [5,Y,29.16],XI), PP = spline (X, [5,Y,29.16])
运行后屏幕显示压紧样条分别在x1??0.02，x2?2.56处的插值和该样条函数的有关信息如下
YI =
-3.1058 4.7834 PP = form: 'pp'
breaks: [-1 -0.5400 0.1300 1.1200 1.8900 2.0600 2.5400 2.8200 3.5000] coefs: [8x4 double] pieces: 8 order: 4 dim: 1
'（2）因为端点约束条件为Sn(?1)?0，S'n(3.50)?0，所以输入MATLAB程序
>> YI= spline (X, [0,Y,0],XI), PP= spline (X, [0,Y,0])
运行后屏幕显示压紧三次样条分别在x1??0.02，x2?2.56的插值和该样条函数的有关信息如下
YI =
-3.0192 4.7501 PP =
form: 'pp'
breaks: [-1 -0.5400 0.1300 1.1200 1.8900 2.0600 2.5400 2.8200 3.5000] coefs: [8x4 double] pieces: 8 order: 4 dim: 1
96．
例6.7.7 在默认端点约束条件下，用两种方法计算在下列插值点处的三次样条插值.
（1）给出节点的数据与例6.7.6相同,插值点XI=2.56； （2）节点(X(i),Y(i))的横坐标向量X从?5到5的整数，纵坐标向量Y=（-2.36，0.85，1.12，2.36，2.36，3，4，1.46，0.49，0.06, 0），插值点向量XI是从?4到4的11个等分点.
解 （1）输入MATLAB程序
>>X=[-1.00 -0.54 0.13 1.12 1.89 2.06 2.54 2.82
3.50];
Y=[-2.46 -5.26 -1.87 0.05 1.65 2.69 4.56 7.89
10.31];
XI=2.56;
Y1= spline (X, Y,XI),Y2=interp1(X,Y,XI,'spline')
运行后屏幕显示三次样条插值的两种结果如下
Y1 = 4.730 2，Y2 =4.730 2. （2）输入MATLAB程序
>> X= -5:5; Y= [-2.36 .85 1.12 2.36 2.36 3 4 1.46 .49 .06 0]; XI = linspace(-4,4,11); Y1= spline (X, Y,XI),
Y2=interp1(X,Y,XI,'spline') 运行后屏幕显示
Y1 =
0.8500 1.0203 1.8857 2.4779 2.3683 3.0000
4.0656 2.5935 0.8247 0.4074 0.0600 Y2 =
0.8500 1.0203 1.8857 2.4779 2.3683 3.0000
4.0656 2.5935 0.8247 0.4074 0.0600
例6.7.8 设函数f(x)?1定义在区间[?1,1]上，取n?10，按等距节点
1?25x2构造分段三次样条函数Sn(x).试用MATLAB程序计算在各小区间中点处分段三次样条插值Sn(xi?1/2)及其相对误差.
解 .在MATLAB工作窗口输入程序
>> h=0.2;x0=-1:h:1;y0=1./(1+25.*x0.^2);xi=-0.9:h:0.9; fi=1./(1+25.*xi.^2); yi=spline (x0,y0,xi);
Ri=abs((fi-yi)./fi); xi,fi,yi,Ri,i=[xi',fi',yi',Ri']
运行后屏幕显示各小区间中点xi处的函数值fi，插值si，相对误差值Ri(略).
6.7.7 几种作三次样条有关图像的MATLAB程序
求有关分段三次样条图形的MATLAB主程序 （一）限定端点约束条件的作图程序
function S=splinetx(x0,y0,xｊ,x,y,dy1,dyn) S = spline(x0,[dy1,y0,dyn],xｊ); Sn = spline(x0,[dy1,y0,dyn],x);
plot(x0,y0,'o',x,Sn,'-',xｊ,S,'*',x,y,'-.')
legend('节点(xi,yi)', '分段三次样条函数','插值点(x,S)','被插值
函数y')
（二）不限定端点约束条件的作图程序
function S=splinetx1(x0,y0,xi,x,y)
S= interp1(x0,y0,xi, 'spline'); Sn= interp1(x0,y0,x,
'spline');
plot(x0,y0,'o',x,Sn,'-',xi,S,'*',x,y,'-.')
legend('节点(xi,yi)', '分段三次样条函数','插值点(x,S)','被插值
函数y')
97．
（三）自由作图程序
直接在MATLAB工作窗口编程序，例如，
>>subplot(2,2,1),x1=-8:4/3:-4,c1=sin(x1);xx1 = -8:0.1:-4; pp1 = interp1 (x1,c1,xx1,'spline ');
cc1 =sin(xx1);%pp1 = spline (x1,c1,xx1); plot(x1,c1,'bo',xx1,pp1,'k-',xx1,cc1,'r-.') subplot(2,2,2)
x2=-4:4/3:-0;c2=sin(x2); xx2 = -4:0.1:0; pp2 = interp1 (x2,c2,xx2,'spline ');
cc2=sin(xx2);plot(x2,c2,'bo',xx2,pp2,'k-',xx2,cc2,'r-.') title('y=sinx及其三次样条插值函数,节点(xi,yi)的图形') subplot(2,1,2)
x=-8:4/3:8;c=sin(x);xx = -8:0.1:8; pp = spline(x,c,xx);
cc=sin(xx); plot(x,c,'bo',xx,pp,'k-',xx,cc,'r-.') legend('节点(xi,yi)','三次样条插值函数','y=sinx 的函数')
或 >> subplot(2,2,1),x1=-8:4/3:-4,c1=cos(x1);xx1 = -8:0.1:-4;
Y1=interp1(x1,c1,xx1, 'pchip');
pp1 = interp1 (x1,c1,xx1,'spline '); cc1 =cos(xx1);
plot(x1,c1,'bo',xx1,pp1,'k-',xx1,Y1,'r:',xx1,cc1,'g-.') subplot(2,2,2)
x2=-4:4/3:-0;c2= cos(x2); xx2 = -4:0.1:0; pp2 = interp1 (x2,c2,xx2,'spline ');
Y2=interp1(x2,c2,xx2, 'pchip');cc2=cos(xx2);
plot(x2,c2,'bo',xx2,pp2,'k-',xx2,Y2,'r:',xx2,cc2,'g-.') title('y=cosx及其三次样条插值函数,分段三次埃尔米特插值, 节点(xi,yi)的图形') subplot(2,1,2)
x=-8:4/3:8;c= cos(x);xx = -8:0.1:8;
pp = spline(x,c,xx);Y=interp1(x,c,xx, 'pchip'); cc= cos(xx);
plot(x,c,'bo',xx,pp,'k-',xx,Y,'r:',xx,cc,'g-.')
legend('节点(xi,yi)','三次样条插值函数','分段三次埃尔米特插值','y=cosx 的函数')
或 >>n=7,h=2*pi/n; x=(-2*pi+h)/2:h:(2*pi-h)/2;
y= tan(cos((3^(1/2)+sin(2*x))./(3+4*x.^2))); x0=-pi:h:pi;X=-pi:h/12:pi;
y0= tan(cos((3^(1/2)+sin(2*x0))./(3+4*x0.^2))); Y= tan(cos((3^(1/2)+sin(2*X))./(3+4*X.^2)));
YL=lagr1(x0,y0,X); YS=interp1(x0,y0,X,'spline'); YH=interp1(x0,y0,X, 'pchip'); yL=lagr1(x0,y0,x); yX=interp1(x0,y0,x); yS=interp1(x0,y0,x,'spline'); yH=interp1(x0,y0,x, 'pchip');
RL=abs((y-yL)./y);RS=abs((y-yS)./y); RH=abs((y-yH)./y);
RX=abs((y-yX)./y);RLj=abs(y-yL);mRLj=mean(RLj); RSj=abs(y-yS);
mRSj=mean(RSj);RHj=abs(y-yH);RXj=abs(y-yX);mRHj=mean(RHj);
mRXj=mean(RXj);mRL=mean(RL);mRX=mean(RX); mRS=mean(RS);mRH=mean(RH);
CZ=[x' y' yL' yX' yS' yH'],R=[x' RL' RX' RS' RH'], mR=[mRL' mRX' mRS' mRH']
Rj=[x' RLj' RXj' RSj' RHj'],mRj=[mRLj' mRXj' mRSj' mRHj'], plot(x0,y0,'bo',X,Y,'k-',X,YS,'rp', X,YH,'g>'),
legend('节点','被插值函数','三次样条函数','分段埃尔米特插值函数
98．
')
例6.7.9 设函数f(x)?1定义在区间[?5,5]上，节点(X(i),f(X (i)))的横坐1?x2标向量X的元素是首项a=-5，末项b=5，公差h=1.5的等差数列，构造分段三次样条函数Sn(x).把区间[?5,5]分成20等份，构成20个小区间，用限定端点约束条件和不限定端点约束条件的MATLAB程序计算各小区间中点xi处Sn(x)的值，并作出节点，插值点，
f(x)和Sn(x)的图形，并与分段埃尔米特插值函数的图形比较.
解 记节点的横坐标xi??5?ih,h?0.5,i?0,1,2,?,20,插值点xi?121?(xi?xi?1)，2i?0,1,2,?,19.
(1) 不限定端点约束条件，在MATLAB工作窗口输入程序
>>x0=-5:1.5:5; y0=1./(1+x0.^2); x1=-4.75:0.5:4.75;
x=-5:0.001:5;
y=1./(1+x.^2); S= splinetx1(x0,y0,x1,x,y)
title(' y=1/(1+x^2)及其三次样条插值函数,节点和插值点的图形')
运行后屏幕显示各小区间中点xi处Sn(x)的值，作出节点，插值点，f(x)和Sn(x)的图形(略).
'（2）限定端点约束条件，取dy1?dyn?Sn(?5)?0，在MATLAB工作窗口输入
程序
>>x0=-5:1.5:5;y0=1./(1+x0.^2);x1=-4.75:0.5:4.75;
x=-5:0.001:5;
y=1./(1+x.^2); S= splinetx (x0,y0,x1,x,y,0,0)
title(' y=1/(1+x^2)及其分段压紧三次样条函数,节点和插值点的图形')
运行后屏幕显示各小区间中点xi处Sn(x)的值（略），作出节点，插值点，f(x)和Sn(x)的图形(略).
如果调节端点约束条件或者增加节点的倍数，例如，在MATLAB工作窗口输入程序
>> x0=-5:0.5:5; y0=1./(1+x0.^2); x1=-4.75:0.5:4.75;
x=-5:0.001:5;
y=1./(1+x.^2); S= splinetx (x0,y0,x1,x,y,0,0)
title(' y=1/(1+x^2)及其增加节点后的不限定端点约束条件的三次样条函
数,节点和插值点的图形')
则运行后输出的图形中的三次样条函数与被插值函数的图像基本重合.
例6.7.10 设函数f(x)?0.5x?cosx定义在区间[?2?,2?]上，取n?7，按等距节点构造分段三次样条函数Sn(x)，用MATLAB程序计算各小区间中点xi处Sn(x)的值，分别作出局部和整体区间上的节点，插值点，f(x)，三次样条函数Sn(x)和分段三次埃尔米特插值函数Hn(x)的图形，并进行比较.
解 编写并保存名为sanci.m的M文件如下
subplot(2,2,1)
h=4*pi/7; x1=-2*pi:h:-2*pi+3*h;c1=0.5.*x1-cos(x1); xx1 =-2*pi+4*pi/14:h:-2*pi+pi/11+2*h; X1=-2*pi:0.001:-2*pi+3*h;
Y1=interp1(x1,c1,xx1, 'pchip'); YY1=interp1(x1,c1,X1,
'pchip');
pp1 = interp1 (x1,c1,xx1,'spline '); P1 = interp1 (x1,c1,X1,'spline '); cc1 =0.5.*X1-cos(X1);
plot(x1,c1,'bo',xx1,pp1,'k*',X1,P1,'k-',xx1,Y1,'rx',X1,
99．