三次样条插值实验报告
数值分析实验报告-插值、三次样条

实验报告:牛顿差值多项式&三次样条问题:在区间[-1,1]上分别取n=10、20用两组等距节点对龙格函数21()25f x x作多项式插值及三次样条插值,对每个n 值,分别画出插值函数及()f x 的图形。
实验目的:通过编程实现牛顿插值方法和三次样条方法,加深对多项式插值的理解。
应用所编程序解决实际算例。
实验要求:1. 认真分析问题,深刻理解相关理论知识并能熟练应用;2. 编写相关程序并进行实验;3. 调试程序,得到最终结果;4. 分析解释实验结果;5. 按照要求完成实验报告。
实验原理:详见《数值分析 第5版》第二章相关内容。
实验内容:(1)牛顿插值多项式1.1 当n=10时:在Matlab 下编写代码完成计算和画图。
结果如下:代码:clear allclcx1=-1:0.2:1;y1=1./(1+25.*x1.^2);n=length(x1);f=y1(:);for j=2:nfor i=n:-1:jf(i)=(f(i)-f(i-1))/(x1(i)-x1(i-j+1));endendsyms F x p ;F(1)=1;p(1)=y1(1);for i=2:nF(i)=F(i-1)*(x-x1(i-1));p(i)=f(i)*F(i);endsyms PP=sum(p);P10=vpa(expand(P),5);x0=-1:0.001:1;y0=subs(P,x,x0);y2=subs(1/(1+25*x^2),x,x0);plot(x0,y0,x0,y2)grid onxlabel('x')ylabel('y')P10即我们所求的牛顿插值多项式,其结果为:P10(x)=-220.94*x^10+494.91*x^8-9.5065e-14*x^7-381.43*x^6-8.504e-14*x^5+123.36*x^4+2.0202e-1 4*x^3-16.855*x^2-6.6594e-16*x+1.0并且这里也能得到该牛顿插值多项式的在[-1,1]上的图形,并和原函数进行对比(见Fig.1)。
学会运用三次样条函数插值

实验目的:学会运用三次样条函数进行函数插值实验内容:给定插值条件如下:i 0 1 2 3 4 5 6 7X i8.125 8.4 9.0 9.485 9.6 9.959 10.166 10.2 Y i0.0774 0.099 0.280 0.60 0.708 1.200 1.800 2.177 作三次样条函数插值,取第一类边界条件Y0’=0.01087 Y7’=1001.画出插值曲线的图像。
程序如下:%将插值点的坐标值输入matlabx=[8.125 8.4 9.0 9.485 9.6 9.959 10.166 10.2];y=[0.0774 0.099 0.280 0.60 0.708 1.200 1.800 2.177];%调用三次样条函数,并给出第一类边界条件即 Y0’=0.01087 Y7’=100 cs = spline(x,[0.01087 y 100]);%给定x(0)到x(7)之间的作图点数,即1000xx=linspace(8.125,10.2,1000);%绘制图像plot(x,y,'o',xx,ppval(cs,xx),'-')grid on插值曲线:2.逆时针旋转座标轴45o保持(1.)中结点和边界条件的几何关系不变,再次作三次样条函数插值,画出插值曲线的图像。
程序如下:x=[8.125 8.4 9.0 9.485 9.6 9.959 10.166 10.2];y=[0.0774 0.099 0.280 0.60 0.708 1.200 1.800 2.177];%用坐标转换求出坐标轴旋转45度后插值点新的坐标值x1=x*cos(pi/4)+y*sin(pi/4);y1=-x*cos(pi/4)+y*sin(pi/4);%求出坐标轴转换后新的边界值Y0=(0.01087-1)/(0.01087+1);Y7=(100-1)/(100+1);cs = spline(x1,[Y0 y1 Y7]);xx=linspace(x1(1),x1(8),1000);plot(x1,y1,'o',xx,ppval(cs,xx),'-')grid on其中在求解边界条件用转换关系tan(45)θ-,θ为在原图像中边界点的切线与原x坐标轴的夹角差值曲线:。
Matlab实验报告六(三次样条与分段线性插值)范文

本题是给出粗略等分点让你插入更多点用双线性插值法来作出更清晰的山区地貌图。
2.问题求解
x=0:400:2800;
y=0:400:2400;
z=[1430 1450 1470 1320 1280 1200 1080 940;
1450 1480 1500 1550 1510 1430 1300 1200;
2.分段线性插值与计算量与n无关;n越大,误差越小.
3.三次样条插值比分段线性插值更光滑。
4.‘linear’:分段线性插值;‘spline’:三次样条值。
【实验环境】
MatlabR2010b
二、实验内容
问题1对函数 ,x[-5,5],分别用分段线性插值和三次样条插值作插值(其中插值节点不少于20),并分别作出每种插值方法的误差曲线.
本次实验因为是我们课本没有的内容,心理上给了我很大的压力,幸好我们还能根据老师的课件以及例题去掌握这次实验所需要的各种插值法,但结果还好,两道题都做出来了。
plot(x,y,'*',x1,yl,'r',x1,y2,'b')
y0=1./(1+x1.^2);
y3=yl-y0;
y4=ys-y0;
holdon
plot(x1,y3,'y',x1,y4,'g')
3.结果
4误。
问题2山区地貌图在某山区(平面区域(0,2800)(0,2400)内,单位:米)测得一些地点的高程(单位:米)如表1,试作出该山区的地貌图.
1.分析问题
本题先取出少量的插值节点并作出图形,再用分段线性插值法和三次样条插值法做出更精确的图形,最后在作出误差曲线。
LAB03_三次样条插值实验

Lab03.三次样条插值实验【实验目的和要求】1.使学生深入理解三次样条插值法,深入进行程序设计能力训练;2.对第一与第二种边界条件,按三弯矩法,通过用Matlab 语言设计计算三次样条插值的程序,以提高学生程序设计的能力。
【实验内容】1.根据Matlab 语言特点,描述三次样条插值法。
2.对第一与第二种边界条件,按三弯矩法,用Matlab 语言设计计算三次样条插值的程序。
3对(1) 自然边界条件0)0.1()2.0(=''=''S S ;(2) 第一种边界条件55741.1)0.1( ,20271.0)2.0(='='S S . 输出用追赶法解出的弯矩向量),,(521M M M 和)1.02.0(i S + (i =0,1,…,8)的值,并画出)(x S y =的图形。
4.完成教材P45例8的计算,并将计算结果与Langrage 插值法计算的结果进行比较,由此说明三次样条插值的优越性。
【实验仪器与软件】1.CPU 主频在1GHz 以上,内存在128Mb 以上的PC ;2.Matlab 6.0及以上版本。
实验讲评:实验成绩:评阅教师:2011 年 月 日Lab03.三次样条插值实验一、算法描述1.根据Matlab语言特点,描述三次样条插值法.答:S(x) 在[x j, x j+1](j=1,2,⋯,n-1)上是三次多项式,于是S"(x)在[x j, x j+1] 上是一次多项式,如果S"(x) 在[x j,x j+1](j=1,2,⋯,n-1)两端点上的值已知,设S"(x j)=M j,S"(x j+1)=M j+1,则S"(x) 的表达式为:= ,其中h j =x j+1-x j,对S"(x) 进行两次积分,则得到1 个具有2二、程序设计2.对第一与第二种边界条件,按三弯矩法,用Matlab语言设计计算三次样条插值的程序。
Matlab实验报告六(三次样条与分段线性插值)

实验名称插值与拟合
所属课程数学软件与实验
实验类型综合型实验
专业信息与计算科学
班级
学号
姓名
指导教师
一、实验概述
【实验目的】
学会在matlab环境下使用几种不同的插值法和拟合两种方法构造函数依据已经知道的某些特殊点来推测实际问题中需要知道但又不便于测量出来的量。
【实验原理】
1.z=interp2(x0,y0,z0,x,y,’method’): 要求x0,y0单调;x, y可取为矩阵, 或x取行向量, y取为列向量, x,y的值分别不能超出x0,y0的范围。
2.分段线性插值与计算量与n无关;n越大, 误差越小.
3.三次样条插值比分段线性插值更光滑。
4.‘linear’ : 分段线性插值;‘spline’ : 三次样条
二、实验内容
问题1 对函数, x([-5,5], 分别用分段线性插值和三次样条插值作插值(其中插值节点不少于20), 并分别作出每种插值方法的误差曲线.
1180 1320 1450 1420 1400 1300 700 900];
mesh(x,y,z)
xi=0:20:2800;
yi=0:20:2400;
zi=interp2(x,y,z,xi',yi,'cubic');
mesh(xi,yi,zi)
3.结果
4.结论及分析
通过实验,结果正确,分析无误。
三、实验小结
1270 1500 1200 1100 1350 1450 1200 1150
1230 1390 1500 1500 1400 900 1100 1060
1180 1320 1450 1420 1400 1300 700 900
三次样条插值

三次样条插值的数值实验姓名: 王维滨 学号:0842011157 姓名: 李佳乐 学号:0842011034 姓名: 谢朝 学号:0842011062 姓名: 杨其荣 学号:08420110721.实验项目的性质和任务对三次样条插值进一步理解,并编写matlab 程序,实现这些功能2.算法设计和matlab 编程。
总前提:x i 第i 点的横坐标i y 第i 点的纵坐标i M ,,s 的记号,,,,,,23()()()()()()()2!3!i i i i i i i s x s x s x y s x x x x x x x =+-+-+- 有数值逼近书上的推导,我们令:111i i i i i x x u x x -+--=-,111i i i i i x x x x λ++--=-,11111111f (,,)i i i i i i i i i i i i i y y y y x x x x x x x x x +-+--++------=- 由于未知数的数目多于方程的个数,我们需要增加两个条件才能唯一确定一个分段三次函数1)D1的三次样条插值a .实验方案与原理:我们加上条件:,,,,11()(),()()n n s x f x s x f x ==我们建立三弯矩方程组:1211211111126(,,)26(,,),2,3.....126(,,)i i i i i i i i n n n n n M M f x x x u M M M f x x x i n M M f x x x λ-+-+--+=⎧⎪++==-⎨⎪+=⎩然后采用追赶法迭代求方程组,但是我们在程序中采用简单的方法(矩阵计算)直接求解降低编程难度,2)D2三次样条插值223123211i 11112111126(,,)26(,,),3,4.....226(,,)i i i i i i i n n n n n n n n M M f x x x u M u M M M f x x x i n u M M f x x x M λλλ-+-+----+-+=-⎧⎪++==-⎨⎪++=-⎩3)D3三次样条插值22321231i 111211126(,,)26(,,),3,4.....126(,,)n i i i i i i i n n n n n n n M M u M f x x x u M M M f x x x i n M u M M f x x x λλλ-+-+--+++=⎧⎪++==-⎨⎪++=⎩b .编写程序:见附录调用test,每次显示一个图像,关闭后按回车继续显示下一幅。
关于三次样条插值函数的学习报告

关于三次样条插值函数的学习报告三次样条插值函数是一种广泛应用于数值分析领域的插值方法,用于逼近一组已知数据点构成的函数。
在这篇学习报告中,我将介绍三次样条插值函数的定义、原理、应用及其优缺点,并通过实际例子说明其如何在实际问题中使用。
一、三次样条插值函数的定义三次样条插值函数是指用分段三次多项式对一组已知数据点进行插值的方法。
具体来说,对于已知数据点$(x_0,y_0),(x_1,y_1),...,(x_n,y_n)$,三次样条插值函数会在每相邻两个数据点之间构造一个三次多项式,使得这些多项式在相应的数据点上满足插值条件,并且在相邻两个多项式之间满足一定的连续性条件。
二、三次样条插值函数的原理三次样条插值函数的原理是利用三次多项式在每个数据点上的取值和导数值来确定三次多项式的系数,从而构造出满足插值条件和连续性条件的插值函数。
具体来说,对于每个相邻的数据点$(x_i,y_i),(x_{i+1},y_{i+1})$,我们可以构造一个三次多项式$S_i(x)$,满足以下条件:1.$S_i(x_i)=y_i$,$S_i(x_{i+1})=y_{i+1}$,即在数据点上满足插值条件;2.$S_i'(x_{i+1})=S_{i+1}'(x_{i+1})$,$S_i''(x_{i+1})=S_{i+1}''(x_{i+1})$,即在数据点上满足连续性条件。
通过求解上述条件,可以得到每个相邻数据点之间的三次多项式$S_i(x)$,从而得到整个插值函数。
三、三次样条插值函数的应用三次样条插值函数在数值分析领域有广泛的应用,尤其在曲线拟合、数据逼近等问题中起到重要作用。
例如,当我们需要根据已知的离散数据点绘制平滑的曲线图形时,可以使用三次样条插值函数来进行插值,从而得到更加连续和光滑的曲线。
另外,在信号处理、图像处理等领域也常常会用到三次样条插值函数。
例如,在数字图像处理中,我们需要对像素点进行插值以得到更高分辨率的图像,三次样条插值函数可以很好地满足这个需求,使图像更加清晰和真实。
三次样条插值报告

三次样条插值多项式实验的目的及意义:为了取得理想结果:在不增加更多的插值条件下,能够求得一个插值多项式,既有良好的逼近效果,又有好的光滑性,引进三次样条插值 多项式。
如果已知函数y=f(x)在节点a=x0<x1<…<xn=b 处的函数值和导数值:()i i x f y =,i=0,1,2,…,n如果S(x)满足条件:1. S(x)是一个分段的三次多项式且()i i y x S =;2. S(x)在[a,b]具有二阶连续导数。
则称S(x)是三次样条插值函数。
S(x)的具体形式为:()()()()⎪⎪⎩⎪⎪⎨⎧∈∈∈=-]12,121,01,[,...............][,][,n n n x x x x s x x x x s x x x x s x s其中()x S i 在[]i i x x ,1-上是三次多项式()iiiiid x c x b x a x S +++=23由插值条件()ii y x S =,i=0,1,2,…,n ,得n+1个条件。
边界条件一:()()nn y x S y x S '',''00== 边界条件二:()()nn y x S y x S '''',''''00==数学公式:()()2211133[2]()[2()]()i i i i i i i i i i ih x x x x h x x x x H x y y h h ---+-----=++2211122()()()()i i i i i i i ix x x x x x x x m m h h -------+ 算法描述:Step1:输入未知数X 及(xi,yi),i=0,1,…,n ; Step2:计算步长H[i]; Step3:计算[][]()⎪⎪⎪⎩⎪⎪⎪⎨⎧+=-=+=-+++i i i i i i i ii i i i i x x f x x f u g u hh h ,,311111λλλStep4:根据边界条件,求解相应的方程得到m0,m1,…, mn Step5:判断X 属于[]i i x x ,1-,i=1,2,…,n 中的哪一个Step6:计算()x s y i i ≈Step7:输出y. 程序原代码如下: #include "stdio.h" #define N 5 void main() { int i,k; float X,s,y0,yn;float a[N][N+1],h[N],u[N],v[N],g[N],m[N],p[N],q[N],w[N];printf("please input X:"); //X 为未知数的大小scanf("%f",&X);printf("please input x:"); //输入x的大小for(i=0;i<N;i++)scanf("%f",&a[i][0]);printf("please input y:"); //输入y的大小for(i=0;i<N;i++)scanf("%f",&a[i][1]);for(i=1;i<N;i++)h[i]=a[i][0]-a[i-1][0]; //计算步长for(i=1;i<N;i++){v[i]=h[i+1]/(h[i]+h[i+1]);u[i]=1-v[i];g[i]=3*u[i]*(a[i+1][1]-a[i][1])/h[i+1]+3*v[i]*(a[i][1]-a[i-1][1])/h[i]; }printf("\t(1)已知边界条件1\n");printf("\t(2)已知边界条件2\n");printf("请选择边界条件序号:");scanf("%d",&k);if(k==1){printf("请输入y0和yn的一阶导:"); //输入边界条件一scanf("%f%f",&m[0],&m[N-1]);p[0]=0; //用追赶法求解m[N]q[0]=0;g[1]=g[1]-v[1]*m[0];g[N-2]=g[N-2]-u[N-2]*m[N-1];for(i=1;i<N;i++){w[i]=2-u[i]*p[i-1];p[i]=v[i]/w[i];q[i]=(g[i]-u[i]*q[i-1])/w[i];}m[N-2]=q[N-2];for(i=N-3;i>0;i--)m[i]=q[i]-p[i]*m[i+1];printf("输出m[i]的值:\n");for(i=0;i<N;i++)printf("%f\n",m[i]);for(i=1;i<N;i++) //计算最终结果if(X>a[i-1][0]&&X<a[i][0])s=(h[i]+2*(X-a[i-1][0]))*(X-a[i][0])*(X-a[i][0])*a[i-1][1]/(h[i]*h[i]*h[i]) +(h[i]-2*(X-a[i][0]))*(X-a[i-1][0])*(X-a[i-1][0])*a[i][1]/(h[i]*h[i]*h[i])+ (X-a[i-1][0])*(X-a[i][0])*(X-a[i][0])*m[i-1]/(h[i]*h[i])+(X-a[i-1][0])*(X-a[i-1][0])*(X-a[i][0])*m[i]/(h[i]*h[i]);printf("s(%f)=%f\n",X,s);}if(k==2){printf("请输入y0和yn的二阶导:"); //输入边界条件二scanf("%f%f",&y0,&yn);g[0]=3*(a[1][1]-a[0][1])/h[1]-h[1]*y0/2;g[N-1]=3*(a[N-1][1]-a[N-2][1])/h[N-1]+h[N-1]*yn/2;q[0]=g[0];u[0]=1;v[N-1]=1;w[0]=2;for(i=1;i<N;i++){w[i]=2-v[i]*u[i-1]/w[i-1];q[i]=g[i]-v[i]*q[i-1]/w[i-1];}m[N-1]=q[N-1]/w[N-1];for(i=N-2;i>=0;i--)m[i]=(q[i]-u[i]*m[i+1])/w[i];printf("输出m[i]的值:\n");for(i=0;i<N;i++)printf("%f\n",m[i]);for(i=1;i<N;i++)if(X>=a[i-1][0]&&X<=a[i][0])s=(h[i]+2*(X-a[i-1][0]))*(X-a[i][0])*(X-a[i][0])*a[i-1][1]/(h[i]*h[i]*h[i]) +(h[i]-2*(X-a[i][0]))*(X-a[i-1][0])*(X-a[i-1][0])*a[i][1]/(h[i]*h[i]*h[i])+ (X-a[i-1][0])*(X-a[i][0])*(X-a[i][0])*m[i-1]/(h[i]*h[i])+(X-a[i-1][0])*(X-a[i-1][0])*(X-a[i][0])*m[i]/(h[i]*h[i]);printf("s(%f)=%f\n",X,s);}}数值计算:已知y=f(x)的如下数值求三次样条插值函数S(x),满足条件1.s’(0)=0,s’(4)=482.s’’(0)=0,s’’(4)=24Please input X:2.5Please input x:0 1 2 3 4 Please input y:-8 -7 0 19 56(1)已知边界条件1(2)已知边界条件2请选择边界条件的序号:1请输入y0和yn的一阶导:0 48 0.0000003.00000012.00000027.00000048.000000s(2.500000)=7.625000press any key tocontinue请选择边界条件的序号:2请输入y0和yn的二阶导:0 24 -0.0000003.00000011.99999927.00000248.0000007.625000press any key tocontinue s(2.500000)=7.625000 对计算结果进行评价分析:()()443845h M x S x f ≤-三次样条插值函数与三次Hermite 插值函数相比,不仅光滑度有提高,而且要求求解时还不需要增加内节点处的导数值,因此比较实用。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
三次样条插值实验报告专业班级学号姓名一、 实验内容和要求1、阅读上面的文字和程序,试运行,检验程序和上面叙述的正确性。
2、阅读上面的 MATLAB 程序;查资料,了解各 MATLAB 语句及命令。
3、画程序流程图,理解并描述算法。
4、修改上面的程序,能根据给定数据点,求(1)自然样条插值,边界 S ′′(a )=0, S ′′(b )=0 (2)第二类边界条件,S ′′(a )和S ′′(b )是确定的。
5、使用上面的程序,根据数据点(0,1),(1,0),(2,0),(3,1),(4,2), (5,2),和(6,1),求三种不同的三次样条插值,其中S ′(0)=−0.6,S ′(6)=-1.8,S ′′(0)=1,S ′′(6)=−1;S ′′(0),S ′′(6)=0.在同一坐标系中,画出这 3 个三次样条插值和这些数据点。
6、写实验报告(实验内容+算法描述+程序+写成分段函数的结果描述+截图)。
二、 算法说明定义:设有N+1个点 ,其中 。
如果存在N 个三次多项式 ,系数为 ,满足如下性质: (1)(2)01n a x x x b =<<<=()k S x 1[,],0,1,,1k k x x x k N +∈=-0[,]Nk k k x y =,0,1,2,3(),(),(),()k k k k S x S x S x S x 23,0,1,2,3()S ()()()()k k k k k k k k S x x s s x x s x x s x x ==+-+-+-()k k S x y =0,1,,k N =(3)则称函数 为三次样条函数。
因为是分段三次多项式,所以 在区间 上是分段线性的。
(1)用 代入上式,可得(2)将(2)式积分两次,会引入两个积分常数,并得到(3) 将 代入上(3),并使用 ,可分别得到包含的方程: (4)求解这两个方程,求出 ,而且将这些值代入方程(3)中,可得到三次多项式方程:(5)表达式5可以简化成只包含未知系数的形式。
为求这些值,必须使用(5)式的导数,即(6) 1lim ()lim ()kkk k x x x x S x S x +--→→=1,2,,1k n =-1lim ()lim ()kkkk x x x x S x S x +--→→''=1,2,,1k n =-1lim ()lim ()kkkk x x x x S x S x +--→→''''=1,2,,1k n =-()S x 1111()S (x )S (x )k kkk k k k k k x x x x S x x x x x ++++--''''''=+--()S x ()S x ''[,]a b 111(),()kkk k kk kM S x M S x h x x +++''''===-和11()M k kkk k k k x x x x S x M h h ++--''=+1,0,1,, 1.k k x x x k N +≤≤=-33111()()()()()66k k k k k k k k k kkMM S x x x x x p x x q x x h h +++=-+-+-+-1,k k x x +和+1+1()()k k k k k k S x y S x y ==和k k p q 和2211,66k k k k k k k k k k M M y h p h y h q h ++=+=+k k p q 和3311111()()()()()6666k k k k k k k kk k k k k k k kkM M y M hy M h S x x x x x x x x x h h h h +++++⎛⎫⎛⎫=--+-+--+-- ⎪⎪⎝⎭⎝⎭{}k M 221111()()()226k k k k k k k kkk k k k k k k M M y M h y m hS x x x x x h h h h h ++++⎛⎫'=--+---+- ⎪⎝⎭在处计算(6),并简化结果可得到 ,其中 (7) 同理,在式(6)中用 并计算在 处的解,可得 (8)利用节点处一阶导函数连续及方程(6)、(7),可得到包含的重要关系式(9)其中方程组(9)中的未知数是要求的值 ,其他的项可通过数据点集 进行简单的计算得到的常量。
因此方程(9)是包含N+1个未知数,具有N-1个线性方程的不定方程组。
所以需要另外两个方程组才能求解,即边界条件。
如果已知,则 100001113(())23(())2N N N N N MM d S x h M M S x d h ---'=--'=--(10)(11) 根据(9)(10)求出 后,可利用下面的公式计算 的样条系数 。
k x 1()36k k k k k k k M m S x h h d +'=--+1k k k k y y d h +-=11,(),k k k S x -'-代替可得到k x 11111()36k k k k k k k Mm S x h h d -----'=++11,,k k k M M M -+11112()k k k k k k k k h M h h M h M μ---++++=16(),1,2,, 1.k k k d d k N μ-=-=-{}0Nk M {}0(,)Nk k k x y =0(),()n S x S x ''{}0Nk M ()k S x 111112222322222111N N N N N N N N N b c M a b c M a b C M a b M μμμμ---------⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦,0k k s y =,1,1(2)6k k k k k h M M s d ++=-,,22k k Ms =,1,36k k k k M M s h +-=,1k k k k y y d h +-=,1k k k h x x +=-三、 源程序csfit1:第一类边界条件function S=csfit1(X,Y,dx0,dxn)N=length(X)-1; H=diff(X); D=diff(Y)./H; A=H(2:N-1);B=2*(H(1:N-1)+H(2:N)); C=H(2:N); C=H(2:N); U=6*diff(D);B(1)=B(1)-H(1)/2;U(1)=U(1)-3*(D(1)-dx0); B(N-1)=B(N-1)-H(N)/2;U(N-1)=U(N-1)-3*(dxn-D(N)); for k=2:N-1temp=A(k-1)/B(k-1); B(k)=B(k)-temp*C(k-1); U(k)=U(k)-temp*U(k-1); EndM(N)=U(N-1)/B(N-1); for k=N-2:-1:1M(k+1)=(U(k)-C(k)*M(k+2))/B(k); EndM(1)=3*(D(1)-dx0)/H(1)-M(2)/2; M(N+1)=3*(dxn-D(N))/H(N)-M(N)/2;{},k jsfor k=0:N-1S(k+1,1)=(M(k+2)-M(k+1))/(6*H(k+1));S(k+1,2)=M(k+1)/2;S(k+1,3)=D(k+1)-H(k+1)*(2*M(k+1)+M(k+2))/6;S(k+1,4)=Y(k+1);endcsfit2:第二类边界条件function S=csfit2(X,Y,dx0,dxn)N=length(X)-1;H=diff(X);D=diff(Y)./H;A=H(2:N-1);B=2*(H(1:N-1)+H(2:N));C=H(2:N);C=H(2:N);U=6*diff(D);U(1)=U(1)-dx0;U(N-1)=U(N-1)-dxn;for k=2:N-1temp=A(k-1)/B(k-1);B(k)=B(k)-temp*C(k-1);U(k)=U(k)-temp*U(k-1);EndM(N)=U(N-1)/B(N-1);for k=N-2:-1:1M(k+1)=(U(k)-C(k)*M(k+2))/B(k);EndM(1)=dx0;M(N+1)=dxn;for k=0:N-1S(k+1,1)=(M(k+2)-M(k+1))/(6*H(k+1));S(k+1,2)=M(k+1)/2;S(k+1,3)=D(k+1)-H(k+1)*(2*M(k+1)+M(k+2))/6; S(k+1,4)=Y(k+1);end画图:x1=0:.01:1;y1=polyval(S1(1,:),x1-X(1));x2=1:.01:2;y2=polyval(S1(2,:),x2-X(2));x3=2:.01:3;y3=polyval(S1(3,:),x3-X(3));x4=3:.01:4;y4=polyval(S1(4,:),x4-X(4));x5=4:.01:5;y5=polyval(S1(5,:),x5-X(5));x6=5:.01:6;y6=polyval(S1(6,:),x6-X(6));>> plot(x1,y1,x2,y2,x3,y3,x4,y4,x5,y5,x6,y6,X,Y,'.') >> hold on>> x1=0:.01:1;y1=polyval(S2(1,:),x1-X(1));x2=1:.01:2;y2=polyval(S2(2,:),x2-X(2));x3=2:.01:3;y3=polyval(S2(3,:),x3-X(3));x4=3:.01:4;y4=polyval(S2(4,:),x4-X(4));x5=4:.01:5;y5=polyval(S2(5,:),x5-X(5));x6=5:.01:6;y6=polyval(S2(6,:),x6-X(6));>> plot(x1,y1,x2,y2,x3,y3,x4,y4,x5,y5,x6,y6,X,Y,'.') >> hold on>> x1=0:.01:1;y1=polyval(S3(1,:),x1-X(1));x2=1:.01:2;y2=polyval(S3(2,:),x2-X(2));x3=2:.01:3;y3=polyval(S3(3,:),x3-X(3));x4=3:.01:4;y4=polyval(S3(4,:),x4-X(4));x5=4:.01:5;y5=polyval(S3(5,:),x5-X(5));x6=5:.01:6;y6=polyval(S3(6,:),x6-X(6));>> plot(x1,y1,x2,y2,x3,y3,x4,y4,x5,y5,x6,y6,X,Y,'.')四、实验结果第一类边界条件:第二类边界条件:自然边界条件:图像:五、说明与分析实验结果显示,对同一组数据,不同边界条件,会导致形成曲线有一定差距,但总体趋势不变。