样条插值函数及应用

样条插值函数及应用

摘要

样条函数具有广泛的应用，是现代函数论的一个十分活跃的分支，是计算方法的主要基础和工具之一，由于生产和科学技术向前发展的推动以及电子计算机广泛应用的需要，人们便更多地应用这个工具，也更深刻的认识了它的本质。

在实际问题中所遇到许多函数往往很复杂，有些甚至是很难找到解析表达式的。根据函数已有的数据来计算函数在一些新的点处的函数值，就是插值法所需要解决的问题。

插值法是数值逼近的重要方法之一，它是根据给定的自变量值和函数值,求取未知函数的近似值。早在一千多年前,我国科学家就在研究历法时就用到了线性插值和二次插值。而在实际问题中,有许多插值函数的曲线要求具有较高的光滑性,在整个曲线中,曲线不但不能有拐点,而且曲率也不能有突变。因此，对于插值函数必须二次连续可微且不变号 ,这就需要用到三次样条插值。

关键词三次样条函数；插值法

目录

引言 (1)

第一章三次样条插值 (2)

1.1 样条插值函数简介 (2)

1.2 三次样条函数应用 (3)

第二章AMCM91A 估计水塔水流量 (5)

2.1 理论分析及计算 (6)

2.2运用MATLAB软件计算 (9)

参考文献 (14)

引言

样条函数具有广泛的应用，是现代函数论的一个十分活跃的分支，是计算方法的主要基础和工具之一，由于生产和科学技术向前发展的推动以及电子计算机广泛应用的需要，人们便更多地应用这个工具，也更深刻的认识了它的本质。上世纪四十年代，在研究数据处理的问题中引出了样条函数，例如，在1946年Schoenberg将样条引入数学,即所谓的样条函数，直到五十年代，还多应用于统计数据的处理方面，从六十年代起，在航空、造船、汽车等行业中，开始大量采用样条函数。

在我国，从六十年代末开始，从船体数学放样到飞机外形设计，逐渐出现了一个使用样，逐渐出现了一个使用样条函数的热潮，并推广到数据处理的许多问题中。

在实际生活中有许多计算问题对插值函数的光滑性有较高的要求，例如飞机机翼外形、发动机进、排气口都要求有连续的二阶导数，用三次样条绘制的曲线不仅有很好的光滑度，而且当节点逐渐加密时其函数值整体上能很好地逼近被插函数，相应的导数值也收敛于被插函数的导数值，不会发生“龙格现象”。

现在国内外学者对这方面的研究也越来越重视，根据我们的需要来解决不同的问题，而且函数的形式也在不断地改进，长期以来很多学者致力于样条插值的研究，对三次样条的研究已相当成熟。

第一章 三次样条插值

1.1 样条插值函数简介

在实际问题中所遇到许多函数)(x f 往往很复杂，有些甚至是很难找到解析表达式的。有时通过实验或者数值计算所得到的也只是一些离散的点()n i x i ...2,10，=上的函数值，即n i x f i i ...2,1,0),(y ==。

根据函数)(x f 已有的数据来计算函数)(x f 在一些新的点x 处的函数值，就是插值法所需要解决的问题。

插值法的基本思想就是，首先根据已有的函数值来构造一个简单的函数)(y x 作为)(x f 的近似表达式，然后用)(y x 来计算新的点上的函数值作为

)(x f 的近似值。通常可以选多项式函数作为近似函数)(y x ，因为多项式具有各阶导数，

求值也比见方便。常用的有Lagrange 插值、Newton 插值、Hermite 插值和样条插值。线性插值在分段点上仅连续而不可导，三次埃尔米特插值有连续的一阶导数，这样的光滑程度常不能满足物理问题的需要，样条函数可以同时解决这两个问题,使插值函数既是低阶分段函数,又是光滑的函数，并且只需在区间端点提供某些导数信息。

三次样条函数定义：

设在区间[]b ,a 上取1n +个节点b x x x x a n =<<<=...210，函数)(y x f =在各个节点处的函数值为,,...,1,0),(y i n i x f i ==若)(S x 满足：

(1) ;,...1,0,y )(S n i x i i ==

(2) 在区间[]b ,a 上，)(S x 具有连续的二阶导数；

(3) 在区间[]）（1,...,1,0x ,x 1-=-n i i i 上,)(S x 是x 三次的多项式； 则称)(S x 是函数)(y x f =在区间[]b ,a 上的三次样条插值函数。

由以上定义可以看出，虽然每个子区间上的多项式可以各不相同，但在相邻子区间的连接处却是光滑的。因此，样条插值也称为分段光滑插值。

从定义知要求出)(S x ，在每一个小区间[]）（1,...,1,0x ,x 1-=-n i i i 上确定4个待定系数，共有n 个小区间，故应有4n 个参数。

根据)(S x 在[]b ,a 上二阶导数连续，在节点()1-...2,1n i x i =出满足连续性条件

)0()0(S +=-i i x S x

)0()0(S ''+=-i i x S x )0()0(S ''''+=-i i x S x

共有3n-3个条件，再加上)(S x 满足插值条件;,...1,0,y )(S n i x i i ==共有4n-2个条件，因此还需要2个条件才能确定)(S x 。

通常可在区间端点上各加一个条件（称为边界条件），可根据实际问题的要求给定，通常有以下三种：

（1）已知端点的一阶导数值，即

'''

00')(S ,)(S n n f x f x ==

（2）俩端点的二阶导数已知，即

'''''

'00'')(S ,)(S n n f x f x ==

其特殊情况，0)(S )(S ''0''==n x x 称为自然边界条件。

（3）当)(x f 是以0x x n -为周期的函数时，则要求)(S x 也是周期函数。这时边界条件应满足

)0()0(S +=-i i x S x

)0()0(S ''+=-i i x S x )0()0(S ''''+=-i i x S x

而此时n y y =0。这样确定的样条函数)(S x ，称为周期函数。

1.2 三次样条函数应用

作函数)2sin()73(42x e x x y x ?+-=-在[]1,0取间隔为0.1的点图，用插值进行实验。 使用MATLAB 软件 程序代码如下： %产生原始数据 x=0:0.1:1;

y=(x.^2-3*x+7).*exp(-4*x).*sin(2*x);

%作图

subplot(1,2,1);

plot(x,y,x,y,'ro')

%待求插值点

xx=0:0.02:1;

yy=interp1(x,y,xx,'spline');

%作图

subplot(1,2,2)

plot(x,y,'ro',xx,yy,'b')

运行截图

图1.1 运行结果

第二章AMCM91A 估计水塔水流量

美国某洲的各用水管理机构要求各社区提供以每小时多少加仑计的用水率以及每天总的用水量,但许多社区并没有测量水流入或流出当地水塔的水量的设备,他们只能代之以每小时测量水塔中的水位,精度在0.5%以内,更为重要的是,无论什么时候,只要水塔中的水位下降到某一最低水位L时,水泵就启动向水塔重新充水至某一最高水位H,但也无法得到水泵的供水量的测量数据。因此，在水泵工作时，人们容易建立水塔中的水位与水泵工作时的用水量之间的关系。水泵每天向水塔充水一次或两次，每次约两小时。

表1 白某小镇某天的水塔水位

时间0 3316 6635 10619 13937 17921 21240 25223 28543

水位3175 3110 3054 2994 2947 2892 2850 2797 2752

时间32284 35935 39332 39435 43318 46636 49953 53936 57254

水位2697 水泵

工作

水泵

工作

3550 3445 3350 3260 3167 3087

时间60574 64554 68535 71854 75021 79154 82649 85968 89953

水位3012 2927 2842 2767 2697 水泵

工作

水泵

工作

3475 3397

时间93270

水位3340

单位：时间/秒；水位/0.01英尺

试估计在任何时刻，甚至包括水泵正在工作期间内，水从水塔流出的流量)

(t

f，并估计一天的总用水量，表1中给出了某个真实小镇某一天的真实数据。

表1中给出了从第一次测量开始的以秒为单位的时刻，以及该时刻的高度单位为百分之一英尺的水塔中水位的测量值，例如，3316秒后，水塔中的水位达到31.10英尺。水塔是一个垂直圆形柱体，高为40英尺，直径57英尺，通常当水塔的水位降至27.00英尺时水泵开始向水塔充水，而当水塔的水位升至35.50英尺时水泵停止工作。

2.1 理论分析及计算

1. 水塔充水时间的确定 (1) 第一次充水时间的确定

当时间t=32284秒时，水位26.97英尺,约低于最低水位27英尺，因此可作为第一次开始充水时刻。

当t=39435秒时，水塔水位35.5英尺，恰为最高水位，因此可作为第一次充水的结束时刻。充水时间为dt =(39435-32284)/3600=1.9864小时，也接近充水时间2小时。

(2) 第二次充水时间的确定

当时间t=75021秒时，水位26.97英尺,约低于最低水位27英尺，因此可作为第二次开始充水时刻。

当t=82649秒时，水泵在工作，但充水时间达到dt =(8264-7502)/3600=2.1189小时；但下一时刻t=85968时，水塔水位34.75英尺，低于最高水位35.50 英尺。

因此可将t=82649秒作为第二次充水的结束时刻,且该时刻水位为最大充水高度35.50 英尺。

2.计算各时刻塔内水的体积

单位转换为1英尺=0.3048米, 1升=1/3.785411加仑 体积计算公式为4/.2h d v π=

表2 不同时刻水体积表

其中(1)表示第一段开始，(2) 表示第二段开始，(3) 表示第三段开始； 单位：时间(小时)，水流量：加仑/小时；

时间 水体积 时间 水体积 时间 水体积

0 (1) 606125 8.9678 514872 19.0375 542554 0.9211 593716 10.9542 (2) 677715 19.9594 528236 1.8431 583026 12.0328 657670 20.8392 514872 2.9497 571571 12.9544 639534 22.9581 (3) 677715 3.8714 562599 13.8758 622352 23.8800 663397 4.9781 552099 14.9822 604598 24.9869 648506 5.9000 544081 15.9039 589325 25.9083 637625 7.0064 533963 16.8261 575008 7.9286 525372 17.9317 558781

3.计算各时刻点的水流量(加仑/小时) 水流量公式为：

dt

t dv t f )

()(=

以上25个时刻处的水流量采用差分的方法得到，共分三段分别处理。 差分公式为：

(1) 对每段前两点采用向前差分公式

)

(243)(12

1i i i i i i t t V V V t f --+-=

+++

(2) 对每段最后两点采用向后差分公式

12

134()2()

i i i i i i V V V f t t t ----+=

-

(3) 对每段中间点采用中心差分公式

)

(1288)(12

112i i i i i i i t t V V V V t f -+-+-=

+--++

得到各点水流量表

表3 不同时刻水流量表

其中(1)表示第一段开始，(2) 表示第二段开始，(3) 表示第三段开始： 单位：时间(小时)，水流量：加仑/小时； 4．用三次样条拟合流量数据

对表3中25个时刻点的流量数据采用三次样条插值得到一条光滑曲线，作为任意时刻的流量曲线,见图2.1。

时间 水流量 时间 水流量 时间 水流量 0 (1) 14404 8.9678 11023 19.0375 16652 0.9211 11182 10.954(2) 19469 19.9594 14495 1.8431 10063 12.0328 20195 20.8392 14648 2.9497 11012 12.9544 18941 22.958(3) 15220 3.8714 8798 13.8758 15903 23.8800 15263 4.9781 9991 14.9822 18055 24.9869 13711 5.9000 8124 15.9039 15646 25.9083

9634 7.0064 10161 16.8261 13742 7.9286

8487

17.9317

14962

图2.1 水塔流量图

其中*表示数据点，实线为样条曲线 5.一天总用水量计算 一天流水总量计算: 方法1:直接积分法:

?==24

1)332986()(加仑dt t f S

方法2:分段计算法

第一次充水前用水912535148726061251=-=V (加仑)

第一次充水后第二次充水前用水1628435148726777152=-=V (加仑) [22.9581，23.88]期间用水143186633976777153=-=V (加仑)

第一次充水期间用水: ?

==9542.109678.84)30326()(加仑dt t f V 第二次充水期间用水: ?

==9581

.228392.205)(31605)(加仑dt t f V [23.88，24]期间用水: ?

==24

88

.236)(1524)(加仑dt t f V

总共用水)(3318696

1

2加仑==∑=i i V S

两种方法结果相差%34.01

2

1=-=S S S err 6.水泵水流量计算

第一次充水期间水塔体积增加

)(1628435148726777151加仑=-=?V

充水时间:

)(9864.19678.89542.101小时=-=?t

第一次充水期间水泵平均流量

)/(97246)(1

9542

.109678

.811小时加仑=?+?=

?

t dt

t f V p

第二次充水期间水塔体积增加

)(1628435148726777152加仑=-=?V

充水时间:

)(1189.28392.209581.222小时=-=?t

第二次充水期间水泵平均流量

)/(91769)(2

9581

.228392

.2022小时加仑=?+?=

?

t dt

t f V p

则整个充水期间水泵平均流量

)/(945072/)(21小时加仑=+=p p p

2.2运用MATLAB 软件计算

1.依据理论分析，编写程序代码。 MATLAB 程序代码

c=0.3048; %1英尺等于0.3048米 p=1.0/3.785; %1升=1/3.785411加仑 d=57*c; h=31.75*c;

v=pi*d*d*h/4*1000*p;

data=[0,3175;3316,3110;6635,3054;10619,2994;13937,2947;17921,2892;

21240,2850;25223,2797;28543,2752;32284,2697;39435,3550;

43318,3445;46636,3350;49953,3260;53936,3167;57254,3087;

60574,3012;64554,2927;68535,2842;71854,2767;75021,2697;

82649,3550;85968,3475; 89953,3397;93270,3340]; %原始数据

t=data(:,1)/3600; %计算时间(小时为单位)

v=pi*d*d*data(:,2)/100*c/4*1000*p; %计算体积

%计算差分

n=length(v);

f=zeros(n,1); %存储差分值

n1=10; %计算第一段

for i=1:n1

if i<=2 %前两点采用向前差分

f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2))/(2*(t(i+1)-t(i)));

elseif i<=n1-2 %采用三点中心差分公式

f(i)=abs(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2))/(12*(t(i+1)-t(i)));

elseif i>=n1-1

f(i)=abs(3*v(i)-4*v(i-1)+v(i-2))/(2*(t(i)-t(i-1)));

end

end

n2=21; %计算第二段

for i=n1+1:n2

if i<=n1+2 %前两点采用向前差分

f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2))/(2*(t(i+1)-t(i)));

elseif i<=n2-2

f(i)=abs(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2))/(12*(t(i+1)-t(i)));

elseif i>=n2-1

f(i)=abs(3*v(i)-4*v(i-1)+v(i-2))/(2*(t(i)-t(i-1)));

end

end

n3=25; %计算第三段

for i=n2+1:n3

if i<=n2+2 %前两点采用向前差分

f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2))/(2*(t(i+1)-t(i)));

elseif i<=n3-2

f(i)=abs(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2))/(12*(t(i+1)-t(i))); elseif i>=n3-1

f(i)=abs(3*v(i)-4*v(i-1)+v(i-2))/(2*(t(i)-t(i-1)));

end

end

plot(t,f,'r*'); %画原始点图

tmin=min(t);

tmax=max(t);

tt=tmin:0.1:tmax; %获得离散的时间点,用于作样条曲线

ff=spline(t,f,tt); %计算三次样条插值

hold on

plot(tt,ff,'b'); %画样条曲线

xlabel('时间(小时)');

ylabel('流量(加仑/小时)');

title('水塔流量图');

hold off

dt=0.05;

t2=0.5:dt:24.5; %获得离散的时间点,用于积分

nn=length(t2);

f2=spline(t,f,t2);

s=(f2(1)+f2(nn)+2*sum(f2(2:nn-1)))*dt/2 ;%计算24小时用水量,采用复化梯形公式fprintf('(全部积分法)1天总水流量s= %8.2f\n',s);

v10=v(11)-v(10); %第一次水塔增加的水

dt1=t(11)-t(10);

tp=t(10):dt:t(11); %第一次充水其间流出的水

nn=length(tp);

yp=spline(t,f,tp); %计算三次样条插值

v11=(yp(1)+yp(nn)+2*sum(yp(2:nn-1)))*dt/2 ;

v1=v10+v11; %第一次充水总量

p1=v1/dt1; %第一次充水的平均水流量

v20=v(22)-v(21); %第二次水塔增加的水

dt2=t(22)-t(21);

tp1=t(21):dt:t(22); %第二次充水其间流出的水

nn=length(tp1);

yp1=spline(t,f,tp1); %计算三次样条插值

v21=(yp1(1)+yp1(nn)+2*sum(yp1(2:nn-1)))*dt/2 ;

v2=v20+v21; %第二次充水总量

p2=v2/dt2; %第二次充水的平均水流量

p=(p1+p2)/2; %两次充水平均水流量

fprintf('两次充水平均水流量p= %8.2f\n',p);

vv1=v(1)-v(10); %第一次充水前总流量

vv2=v(11)-v(21); %两次充水间总流量

vv3=v(22)-v(23); %t=83649到85968其间流量

%第一次充水期间流量

ta=t(10):dt:t(11); %获得离散的时间点,用于积分

nn=length(ta);

fa=spline(t,f,ta);

s1=(fa(1)+fa(nn)+2*sum(fa(2:nn-1)))*dt/2 ;

%第二次充水期间流量

tb=t(21):dt:t(22); %获得离散的时间点,用于积分

nn=length(tb);

fb=spline(t,f,tb);

s2=(fb(1)+fb(nn)+2*sum(fb(2:nn-1)))*dt/2 ;

%t=85968到86400流量

tc=t(23):dt:24; %获得离散的时间点,用于积分

nn=length(tc);

fc=spline(t,f,tc);

s3=(fc(1)+fc(nn)+2*sum(fc(2:nn-1)))*dt/2 ;

ss=vv1+vv2+vv3+s1+s2+s3;

fprintf('(部分积分法)1天总水流量ss= %8.2f\n',ss);

%误差分析

err=abs((s-ss)/s);

fprintf('两种计算法总量相对误差%6.2f%%\n',err*100);

2.运行代码，得出结果。

运行结果如下

图2.2 运行结果

参考文献

[1] 施吉林,刘淑珍.计算机数值方法[M].北京:高等教育出版社，1999.

[2] 常庚哲.关于三次样条函数的两点注记[J].数学的实践与认识，1979.

[3] 林成森编著.数值计算方法[M].科学出版社,2005.

[4] 李庆扬著.数值分析[M].清华大学出版社,2010.

样条插值函数与应用

样条插值函数及应用

摘要 样条函数具有广泛的应用，是现代函数论的一个十分活跃的分支，是计算方法的主要基础和工具之一，由于生产和科学技术向前发展的推动以及电子计算机广泛应用的需要，人们便更多地应用这个工具，也更深刻的认识了它的本质。 在实际问题中所遇到许多函数往往很复杂，有些甚至是很难找到解析表达式的。根据函数已有的数据来计算函数在一些新的点处的函数值，就是插值法所需要解决的问题。 插值法是数值逼近的重要方法之一，它是根据给定的自变量值和函数值,求取未知函数的近似值。早在一千多年前,我国科学家就在研究历法时就用到了线性插值和二次插值。而在实际问题中,有许多插值函数的曲线要求具有较高的光滑性,在整个曲线中,曲线不但不能有拐点,而且曲率也不能有突变。因此，对于插值函数必须二次连续可微且不变号 ,这就需要用到三次样条插值。 关键词三次样条函数；插值法

目录 引言 0 第一章三次样条插值 (1) 1.1 样条插值函数简介 (1) 1.2 三次样条函数应用 (2) 第二章AMCM91A 估计水塔水流量 (4) 2.1 理论分析及计算 (5) 2.2运用MATLAB软件计算 (8) 参考文献 (13)

引言 样条函数具有广泛的应用，是现代函数论的一个十分活跃的分支，是计算方法的主要基础和工具之一，由于生产和科学技术向前发展的推动以及电子计算机广泛应用的需要，人们便更多地应用这个工具，也更深刻的认识了它的本质。上世纪四十年代，在研究数据处理的问题中引出了样条函数，例如，在1946年Schoenberg将样条引入数学,即所谓的样条函数，直到五十年代，还多应用于统计数据的处理方面，从六十年代起，在航空、造船、汽车等行业中，开始大量采用样条函数。 在我国，从六十年代末开始，从船体数学放样到飞机外形设计，逐渐出现了一个使用样，逐渐出现了一个使用样条函数的热潮，并推广到数据处理的许多问题中。 在实际生活中有许多计算问题对插值函数的光滑性有较高的要求，例如飞机机翼外形、发动机进、排气口都要求有连续的二阶导数，用三次样条绘制的曲线不仅有很好的光滑度，而且当节点逐渐加密时其函数值整体上能很好地逼近被插函数，相应的导数值也收敛于被插函数的导数值，不会发生“龙格现象”。 现在国内外学者对这方面的研究也越来越重视，根据我们的需要来解决不同的问题，而且函数的形式也在不断地改进，长期以来很多学者致力于样条插值的研究，对三次样条的研究已相当成熟。

Matlab中插值函数汇总和使用说明

MATLAB中的插值函数 命令1：interp1 功能：一维数据插值（表格查找）。该命令对数据点之间计算内插值。它找出一元函数f(x)在中间点的数值。其中函数f(x)由所给数据决定。 x：原始数据点 Y：原始数据点 xi：插值点 Yi：插值点 格式 (1) yi = interp1(x,Y,xi) 返回插值向量yi，每一元素对应于参量xi，同时由向量x 与Y 的内插值决定。参量x 指定数据Y 的点。若Y 为一矩阵，则按Y 的每列计算。yi 是阶数为length(xi)*size(Y,2)的输出矩阵。 (2) yi = interp1(Y,xi) 假定x=1:N，其中N 为向量Y 的长度，或者为矩阵Y 的行数。 (3) yi = interp1(x,Y,xi,method) 用指定的算法计算插值： ’nearest’：最近邻点插值，直接完成计算； ’linear’：线性插值（缺省方式），直接完成计算； ’spline’：三次样条函数插值。对于该方法，命令interp1 调用函数spline、ppval、mkpp、umkpp。这些命令生成一系列用于分段多项式操作的函数。命令spline 用它们执行三次样条函数插值； ’pchip’：分段三次Hermite 插值。对于该方法，命令interp1 调用函数pchip，用于对向量x 与y 执行分段三次内插值。该方法保留单调性与数据的外形； ’cubic’：与’pchip’操作相同； ’v5cubic’：在MATLAB 5.0 中的三次插值。 对于超出x 范围的xi 的分量，使用方法’nearest’、’linear’、’v5cubic’的插值算法，相应地将返回NaN。对其他的方法，interp1 将对超出的分量执行外插值算法。 (4)yi = interp1(x,Y,xi,method,'extrap') 对于超出x 范围的xi 中的分量将执行特殊的外插值法extrap。 (5)yi = interp1(x,Y,xi,method,extrapval) 确定超出x 范围的xi 中的分量的外插值extrapval，其值通常取NaN 或0。

对样条函数及其插值问题的一点认识

对样条函数及其插值问题的一点认识 样条函数是计算数学以及计算机辅助设计几何设计的重要工具。1946年，I. J. Schoenberg 著名的关于一元样条函数的奠定性论文“Contribution to the problem of application of equidistant data by analytic functions ”发表，建立了一元样条函数的理论基础。自此以后，关于样条函数的研究工作逐渐深入。随着电子计算机技术的不断进步，样条函数的理论以及应用研究得到迅速的发展和广泛的应用。经过数学工作者的努力，已经形成了较为系统的理论体系。 所谓（多项式）样条函数，乃指具有一定光滑性的分段（分片）多项式。一元n 次且n －1阶连续可微的样条函数具有如下的表示式： 1()()()()N n n j j j s x p x c x x x +==+--∞<<+∞∑[] 011,00,01,,...,,(1),...,(),,...,,n n n n N n N N u un u u u u x x x x x S x x x x ++++ +≥??=??

Matlab中插值函数汇总和使用说明.

告: Matlab中插值函数汇总和使用说明收藏 命令1 interp1 功能一维数据插值(表格查找。该命令对数据点之间计算内插值。它找出一元函数f(x在中间点的数值。其中函数f(x由所给数据决定。x:原始数据点 Y:原始数据点 xi:插值点 Yi:插值点 格式 (1yi = interp1(x,Y,xi 返回插值向量yi,每一元素对应于参量xi,同时由向量x 与Y 的内插值决定。参量x 指定数据Y 的点。 若Y 为一矩阵,则按Y 的每列计算。yi 是阶数为length(xi*size(Y,2的输出矩阵。 (2yi = interp1(Y,xi 假定x=1:N,其中N 为向量Y 的长度,或者为矩阵Y 的行数。 (3yi = interp1(x,Y,xi,method 用指定的算法计算插值: ’nearest’:最近邻点插值,直接完成计算; ’linear’:线性插值(缺省方式,直接完成计算;

’spline’:三次样条函数插值。对于该方法,命令interp1 调用函数spline、ppval、mkpp、umkpp。这些命令生成一系列用于分段多项式操作的函 数。命令spline 用它们执行三次样条函数插值; ’pchip’:分段三次Hermite 插值。对于该方法,命令interp1 调用函数p chip,用于对向量x 与y 执行分段三次内插值。该方法保留单调性与数据的外形; ’cubic’:与’pchip’操作相同; ’v5cubic’:在MATLAB 5.0 中的三次插值。 对于超出x 范围的xi 的分量,使用方法’nearest’、’linear’、’v5cubic’的插值算法,相应地将返回NaN。对其他的方法,interp1 将对超出的分量执行外插值算法。 (4yi = interp1(x,Y,xi,method,'extrap' 对于超出x 范围的xi 中的分量将执行特殊的外插值法extrap。 (5yi = interp1(x,Y,xi,method,extrapval 确定超出x 范围的xi 中的分量的外插值extrapval,其值通常取NaN 或0。 例1 1.>>x = 0:10; y = x.*sin(x; 2.>>xx = 0:.25:10; yy = interp1(x,y,xx; 3.>>plot(x,y,'kd',xx,yy 复制代码 例2 1.>> year = 1900:10:2010;

三次样条插值课后题集

例1 设)(x f 为定义在[0,3]上的函数，有下列函数值表： 且2.0)('0=x f ，1)('3-=x f ，试求区间[0,3]上满足上述条件的三次样条插值函数)(x s 本算法求解出的三次样条插值函数将写成三弯矩方程的形式： ) ()6()() 6()(6)(6)(211123 13 1j j j j j j j j j j j j j j j j x x h h M y x x h h M y x x h M x x h M x s -- + -- + -+ -= +++++其中，方程中的系数 j j h M 6， j j h M 61+， j j j j h h M y )6(2- ， j j j j h h M y ) 6(211++- 将由Matlab 代码中的变量Coefs_1、Coefs_2、Coefs_3以及Coefs_4的值求出。 以下为Matlab 代码： %============================= % 本段代码解决作业题的例1 %============================= clear all clc % 自变量x 与因变量y ，两个边界条件的取值 IndVar = [0, 1, 2, 3]; DepVar = [0, 0.5, 2, 1.5];

LeftBoun = 0.2; RightBoun = -1; % 区间长度向量，其各元素为自变量各段的长度h = zeros(1, length(IndVar) - 1); for i = 1 : length(IndVar) - 1 h(i) = IndVar(i + 1) - IndVar(i); end % 为向量μ赋值 mu = zeros(1, length(h)); for i = 1 : length(mu) - 1 mu(i) = h(i) / (h(i) + h(i + 1)); end mu(i + 1) = 1; % 为向量λ赋值 lambda = zeros(1, length(h)); lambda(1) = 1; for i = 2 : length(lambda) lambda(i) = h(i) / (h(i - 1) + h(i));

三次样条插值方法的应用

CENTRAL SOUTH UNIVERSITY 数值分析实验报告

三次样条插值方法的应用 一、问题背景 分段低次插值函数往往具有很好的收敛性，计算过程简单，稳定性好，并且易于在在电子计算机上实现，但其光滑性较差，对于像高速飞机的机翼形线船体放样等型值线往往要求具有二阶光滑度，即有二阶连续导数，早期工程师制图时，把富有弹性的细长木条（即所谓的样条）用压铁固定在样点上，在其他地方让他自由弯曲，然后沿木条画下曲线，称为样条曲线。样条曲线实际上是由分段三次曲线并接而成，在连接点即样点上要求二阶导数连续，从数学上加以概括就得到数学样条这一概念。下面我们讨论最常用的三次样条函数及其应用。 二、数学模型 样条函数可以给出光滑的插值曲线（面），因此在数值逼近、常微分方程和偏微分方程的数值解及科学和工程的计算中起着重要的作用。 设区间[]b ,a 上给定有关划分b x x n =<<<= 10x a ，S 为[]b ,a 上满足下面条件的函数。 ● )(b a C S ,2∈； ● S 在每个子区间[]1,+i i x x 上是三次多项式。 则称S 为关于划分的三次样条函数。常用的三次样条函数的边界条件有三种类型： ● Ⅰ型 ()()n n n f x S f x S ''0'',==。 ● Ⅱ型 ()()n n n f x S f x S ''''0'''',==，其特殊情况为()()0''''==n n x S x S 。 ● Ⅲ型 ()() 3,2,1,0,0==j x S x S n j j ，此条件称为周期样条函数。 鉴于Ⅱ型三次样条插值函数在实际应用中的重要地位，在此主要对它进行详细介绍。 三、算法及流程 按照传统的编程方法，可将公式直接转换为MATLAB 可是别的语言即可；另一种是运用矩阵运算，发挥MATLAB 在矩阵运算上的优势。两种方法都可以方便地得到结果。方法二更直观，但计算系数时要特别注意。这里计算的是方法一的程序，采用的是Ⅱ型边界条件，取名为spline2.m 。 Matlab 代码如下： function s=spline2(x0,y0,y21,y2n,x) %s=spline2(x0,y0,y21,y2n,x) %x0,y0 are existed points,x are insert points,y21,y2n are the second

matlab插值(详细 全面)

Matlab中插值函数汇总和使用说明 MATLAB中的插值函数为interp1，其调用格式 为： yi= interp1(x,y,xi,'method') 其中x，y为插值点，yi为在被插值点xi处的插值结果；x,y为向量， 'method'表示采用的插值方法，MATLAB提供的插值方法有几种： 'method'是最邻近插值， 'linear'线性插值； 'spline'三次样条插值； 'cubic'立方插值．缺省时表示线性插值 注意：所有的插值方法都要求x是单调的，并且xi不能够超过x的范围。 例如：在一天24小时内，从零点开始每间隔2小时测得的环境温度数据分别为 12，9，9，10，18 ，24，28，27，25，20，18，15，13， 推测中午12点（即13点）时的温度． x=0:2:24; y=[12 9 9 10 18 24 28 27 25 20 18 15 13]; a=13; y1=interp1(x,y,a,'spline') 结果为： 27.8725 若要得到一天24小时的温度曲线，则： xi=0:1/3600:24; yi=interp1(x,y,xi, 'spline'); plot(x,y,'o' ,xi,yi)

命令1 interp1 功能一维数据插值（表格查找）。该命令对数据点之间计算内插值。它找出一元函数f(x)在中间点的数值。其中函数f(x)由所给数据决定。 x：原始数据点 Y：原始数据点 xi：插值点 Yi：插值点 格式 (1)yi = interp1(x,Y,xi) 返回插值向量yi，每一元素对应于参量xi，同时由向量x 与Y 的内插值决定。参量x 指定数据Y 的点。 若Y 为一矩阵，则按Y 的每列计算。yi是阶数为length(xi)*size(Y,2)的输出矩阵。(2)yi = interp1(Y,xi) 假定x=1:N，其中N 为向量Y 的长度，或者为矩阵Y 的行数。 (3)yi = interp1(x,Y,xi,method) 用指定的算法计算插值： ’nearest’：最近邻点插值，直接完成计算； ’linear’：线性插值（缺省方式），直接完成计算； ’spline’：三次样条函数插值。对于该方法，命令interp1 调用函数spline、ppval、mkpp、umkpp。这些命令生成一系列用于分段多项式操作的函数。命令spline 用它们执行三次样条函数插值； ’pchip’：分段三次Hermite 插值。对于该方法，命令interp1 调用函数pchip，用于对

关于三次样条插值函数的学习报告(研究生)资料

学习报告—— 三次样条函数插值问题的讨论 班级：数学二班 学号：152111033 姓名：刘楠楠

样条函数： 由一些按照某种光滑条件分段拼接起来的多项式组成的函数；最常用的样条函数为三次样条函数，即由三次多项式组成，满足处处有二阶连续导数。 一、三次样条函数的定义： 对插值区间[,]a b 进行划分，设节点011n n a x x x x b -=<< <<=，若 函数2()[,]s x c a b ∈在每个小区间1[,]i i x x +上是三次多项式，则称其为三次样条函数。如果同时满足()()i i s x f x = (0,1,2)i n =，则称()s x 为()f x 在 [,]a b 上的三次样条函数。 二、三次样条函数的确定： 由定义可设：101212 1(),[,] (),[,]()(),[,] n n n s x x x x s x x x x s x s x x x x -∈??∈?=???∈?其中()k s x 为1[,]k k x x -上的三次 多项式，且满足11(),()k k k k k k s x y s x y --== (1,2,,k n = 由2()[,]s x C a b ∈可得：''''''()(),()(),k k k k s x s x s x s x -+-+== 有''1()(),k k k k s x s x -++= ''''1()(),(1 ,2,,1)k k k k s x s x k n -+ +==-， 已知每个()k s x 均为三次多项式，有四个待定系数，所以共有4n 个待定系数，需要4n 个方程才能求解。前面已经得到22(1)42n n n +-=-个方程，因此要唯一确定三次插值函数，还要附加2个条件，一般上，实际问题通常对样条函数在端点处的状态有要求，即所谓的边界条件。 1、第一类边界条件：给定函数在端点处的一阶导数，即 ''''00(),()n n s x f s x f == 2、第二类边界条件：给定函数在端点处的二阶导数，即

数值分析作业-三次样条插值

数值计算方法作业 实验4.3 三次样条差值函数 实验目的： 掌握三次样条插值函数的三弯矩方法。 实验函数: dt e x f x t ? ∞ -- = 2 221)(π 实验内容： （1） 编程实现求三次样条插值函数的算法，分别考虑不同的边界条件； （2） 计算各插值节点的弯矩值； （3） 在同一坐标系中绘制函数f(x)，插值多项式，三次样条插值多项式的曲线 比较插值结果。 实验4.5 三次样条差值函数的收敛性 实验目的： 多项式插值不一定是收敛的，即插值的节点多，效果不一定好。对三次样条插值函数如何呢？理论上证明三次样条插值函数的收敛性是比较困难的，通过本实验可以证明这一理论结果。 实验内容： 按照一定的规则分别选择等距或非等距的插值节点，并不断增加插值节点的个数。 实验要求： （1） 随着节点个数的增加，比较被逼近函数和三样条插值函数的误差变化情 况，分析所得结果并与拉格朗日插值多项式比较； （2） 三次样条插值函数的思想最早产生于工业部门。作为工业应用的例子，考

虑如下例子：某汽车制造商根据三次样条插值函数设计车门曲线，其中一 算法描述： 拉格朗日插值： 错误!未找到引用源。 其中错误!未找到引用源。是拉格朗日基函数，其表达式为：() ∏ ≠=--=n i j j j i j i x x x x x l 0) ()( 牛顿插值： ) )...()(](,...,,[.... ))(0](,,[)0](,[)()(1102101210100----++--+-+=n n n x x x x x x x x x x f x x x x x x x f x x x x f x f x N 其中????? ?? ?? ?????? --=--= --= -)/(]),...,[],...,[(]...,[..],[],[],,[)()(],[01102110x x x x x f x x x f x x x f x x x x f x x f x x x f x x x f x f x x f n n n n i k j i k j k j i j i j i j i 三样条插值： 所谓三次样条插值多项式Sn(x)是一种分段函数，它在节点Xi(a

三次样条插值自然边界条件

例：已知一组数据点，编写一程序求解三次样条插值函数满足 并针对下面一组具体实验数据 0.25 0.3 0.39 0.45 0.53 0.5000 0.5477 0.6245 0.6708 0.7280 求解，其中边界条件为. 1）三次样条插值自然边界条件源程序： function s=spline3(x,y,dy1,dyn) %x为节点，y为节点函数值，dy1，dyn分别为x=0.25,0.53处的二阶导 m=length(x);n=length(y); if m~=n error('x or y输入有误') return end h=zeros(1,n-1); h(n-1)=x(n)-x(n-1); for k=1:n-2 h(k)=x(k+1)-x(k); v(k)=h(k+1)/(h(k+1)+h(k)); u(k)=1-v(k); end g(1)=3*(y(2)-y(1))/h(1)-h(1)/2*dy1; g(n)=3*(y(n)-y(n-1))/h(n-1)+h(n-1)/2*dyn; for i=2:n-1 g(i)=3*(u(i-1)*(y(i+1)-y(i))/h(i)+v(i-1)*(y(i)-y(i-1))/h(i-1)); end for i=2:n-1; A(i,i-1)=v(i-1); A(i,i+1)=u(i-1); end A(n,n-1)=1; A(1,2)=1; A=A+2*eye(n); M=zhuigf(A,g); %调用函数，追赶法求M fprintf('三次样条（三对角）插值的函数表达式\n'); syms X;

for k=1:n-1 fprintf('S%d--%d:\n',k,k+1); s(k)=(h(k)+2*(X-x(k)))./h(k).^3.*(X-x(k+1)).^2.*y(k)... +(h(k)-2*(X-x(k+1)))./h(k).^3.*(X-x(k)).^2.*y(k+1)... +(X-x(k)).*(X-x(k+1)).^2./h(k).^2*M(k)+(X-x(k+1)).*... (X-x(k)).^2./h(k).^2*M(k+1); end s=s.'; s=vpa(s,4); %画三次样条插值函数图像 for i=1:n-1 X=x(i):0.01:x(i+1); st=(h(i)+2*(X-x(i)))./(h(i)^3).*(X-x(i+1)).^2.*y(i)... +(h(i)-2.*(X-x(i+1)))./(h(i)^3).*(X-x(i)).^2.*y(i+1)... +(X-x(i)).*(X-x(i+1)).^2./h(i)^2*M(i)+(X-x(i+1)).*... (X-x(i)).^2./h(i)^2*M(i+1); plot(x,y,'o',X,st); hold on End plot(x,y); grid on %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %调用的函数： %追赶法 function M=zhuigf(A,g) n=length(A); L=eye(n); U=zeros(n); for i=1:n-1 U(i,i+1)=A(i,i+1); end U(1,1)=A(1,1); for i=2:n L(i,i-1)=A(i,i-1)/U(i-1,i-1); U(i,i)=A(i,i)-L(i,i-1)*A(i-1,i); end Y(1)=g(1); for i=2:n Y(i)=g(i)-L(i,i-1)*Y(i-1); end M(n)=Y(n)/U(n,n); for i=n-1:-1:1 M(i)=(Y(i)-A(i,i+1)*M(i+1))/U(i,i);

Matlab中插值函数汇总和使用说明

告： Matlab中插值函数汇总和使用说明收藏 命令1 interp1 功能一维数据插值（表格查找）。该命令对数据点之间计算内插值。它找出一元函数f(x)在中间点的数值。其中函数f(x)由所给数据决定。x：原始数据点 Y：原始数据点 xi：插值点 Yi：插值点 格式 (1)yi = interp1(x,Y,xi) 返回插值向量yi，每一元素对应于参量xi，同时由向量x 与Y 的内插值决定。参量x 指定数据Y 的点。 若Y 为一矩阵，则按Y 的每列计算。yi 是阶数为length(xi)*size(Y,2)的输出矩阵。 (2)yi = interp1(Y,xi) 假定x=1:N，其中N 为向量Y 的长度，或者为矩阵Y 的行数。 (3)yi = interp1(x,Y,xi,method) 用指定的算法计算插值： ’nearest’：最近邻点插值，直接完成计算； ’linear’：线性插值（缺省方式），直接完成计算； ’spline’：三次样条函数插值。对于该方法，命令interp1 调用函数spline、ppval、mkpp、umkpp。这些命令生成一系列用于分段多项式操作的函

数。命令spline 用它们执行三次样条函数插值； ’pchip’：分段三次Hermite 插值。对于该方法，命令interp1 调用函数p chip，用于对向量x 与y 执行分段三次内插值。该方法保留单调性与数据的外形； ’cubic’：与’pchip’操作相同； ’v5cubic’：在MATLAB 5.0 中的三次插值。 对于超出x 范围的xi 的分量，使用方法’nearest’、’linear’、’v5cubic’的插值算法，相应地将返回NaN。对其他的方法，interp1 将对超出的分量执行外插值算法。 (4)yi = interp1(x,Y,xi,method,'extrap') 对于超出x 范围的xi 中的分量将执行特殊的外插值法extrap。 (5)yi = interp1(x,Y,xi,method,extrapval) 确定超出x 范围的xi 中的分量的外插值extrapval，其值通常取NaN 或0。 例1 1.>>x = 0:10; y = x.*sin(x); 2.>>xx = 0:.25:10; yy = interp1(x,y,xx); 3.>>plot(x,y,'kd',xx,yy) 复制代码 例2 1.>> year = 1900:10:2010; 2.>> product = [75.995 91.972 105.711 12 3.203 131.669 150.697 179.323 203.212 226.505

三次样条插值函数matlab程序绝不坑爹

x0=[0 0.9211 1.8431 2.9497 3.8714 4.9781 5.9 7.0064 7.9286 8.9678 10.9542 12.0328 12.9544 13.8758 14.9822 15.9039 16.8261 17.9317 19.0375 19.9594 20.8392 22.9581 23.88 24.9869 25.9083]; >> >> y0=[14405 11180 10063 11012 8797 9992 8124 10160 8488 11018 19469 20196 18941 15903 18055 15646 13741 14962 16653 14496 14648 15225 15264 13708 9633]; >> x=0:0.1:25.9; >> y1=interp1(x0,y0,x,'spline'); >> pp1=csape(x0,y0); %样条插值工具箱函数 y2=ppval(pp1,x); %计算x对应的y值 pp2=csape(x0,y0,'second'); y3=ppval(pp2,x); xydata=[x',y1',y2',y3'] subplot(1,2,1) plot(x0,y0,'+',x,y1) title('Spline1') subplot(1,2,2) plot(x0,y0,'+',x,y2) title('Spline2') dx=diff(x); dy=diff(y2); dy_dx=dy./dx; dy_dx0=dy_dx(1) ytemp=y2(13<=x&x<=15); ymin=min(ytemp); xmin=x(y2==ymin); xymin_1315=[xmin,ymin]

数值分析作业-三次样条插值

数值计算方法作业 实验 三次样条差值函数 实验目的： 掌握三次样条插值函数的三弯矩方法。 实验函数: dt e x f x t ? ∞ -- = 2 221)(π 实验内容： （1） 编程实现求三次样条插值函数的算法，分别考虑不同的边界条件； （2） 计算各插值节点的弯矩值； （3） 在同一坐标系中绘制函数f(x)，插值多项式，三次样条插值多项式的曲线 比较插值结果。 实验 三次样条差值函数的收敛性 实验目的： 多项式插值不一定是收敛的，即插值的节点多，效果不一定好。对三次样条插值函数如何呢理论上证明三次样条插值函数的收敛性是比较困难的，通过本实验可以证明这一理论结果。 实验内容： 按照一定的规则分别选择等距或非等距的插值节点，并不断增加插值节点的个数。 实验要求： （1） 随着节点个数的增加，比较被逼近函数和三样条插值函数的误差变化情 况，分析所得结果并与拉格朗日插值多项式比较； （2） 三次样条插值函数的思想最早产生于工业部门。作为工业应用的例子，考 虑如下例子：某汽车制造商根据三次样条插值函数设计车门曲线，其中一

k x 012345678910 k y k y' 算法描述： 拉格朗日插值： 其中是拉格朗日基函数，其表达式为： () ∏ ≠ = - - = n i j j j i j i x x x x x l ) ( ) ( 牛顿插值： ) )...( )( ]( ,... , , [ .... ) )( ]( , , [ )0 ]( , [ ) ( ) ( 1 1 2 1 1 2 1 1 - - - - + + - - + - + = n n n x x x x x x x x x x f x x x x x x x f x x x x f x f x N 其中 ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? - - = - - = - - = - ) /( ]) ,... , [ ] ,... , [ ( ] ... , [ . . ] , [ ] , [ ] , , [ ) ( ) ( ] , [ 1 1 2 1 1 x x x x x f x x x f x x x f x x x x f x x f x x x f x x x f x f x x f n n n n i k j i k j k j i j i j i j i 三样条插值： 所谓三次样条插值多项式Sn(x)是一种分段函数，它在节点Xi(a

matlab-牛顿插值法-三次样条插值法

(){}2 1 ()(11),5,10,20: 1252 1()1,(0,1,2,,)()2,(0,1,2,,)()()2 35,20:1100 (i i i i n n k k k Newton f x x n x f x x i i n f x n x y i n Newton N x S x n x k y f x =-≤≤=+=-+====-+ =L L 题目：插值多项式和三次样条插值多项式。 已知对作、计算函数在点处的值； 、求插值数据点的插值多项式和三次样条插值多项式；、对计算和相应的函数值),()() (1,2,,99)4:()max ()()max ()n k n k n k n k n k n k k k N x S x k E N y N x E S y S x ==-=-L 和； 、计算，； 解释你所得到的结果。 算法组织: 本题在算法上需要解决的问题主要是：求出第二问中的Newton 插值多项式 )(x N n 和三次样条插值多项式()n S x 。如此，则第三、四问则迎刃而解。计算两 种插值多项式的算法如下： 一、求Newton 插值多项式)(x N n ，算法组织如下： Newton 插值多项式的表达式如下： )())(()()(110010--???--+???+-+=n n n x x x x x x c x x c c x N 其中每一项的系数c i 的表达式如下： 1102110) ,,,(),,,(),,,(x x x x x f x x x f x x x f c i i i i i -???-???= ???=- 根据i c 以上公式，计算的步骤如下： ?? ??? ?? ?????+??????? ???????????----) ,,,,(1) ,,,(),,,,(),(,),,(2)(,),(),(11101111011010n n n n n n n n x x x x f n x x x f x x x f n x x f x x f x f x f x f 、计算、计算、计算、计算 二、求三次样条插值多项式)(x S n ，算法组织如下：

三次样条插值函数

一．介绍

二．程序框图 开始 输入未知数X及 (xi,yi),i=0,1,…,n 计算步长H[i] 计算λ、μ、d 根据边界条件，求 解相应的方程得到 M1,…, Mn 将M代入原方程， 得到分段函数 结束

三．源码 syms h n=9;%插入节点数，可以根据题目更改 h=2/(n+1); u=0.5; v=0.5; f=inline('1/(1+25*x.^2)');%输入函数，这个也可以根据题目更改g=inline('3/h*((c-b)/h-(b-a)/h)','a','b','c','h'); for i=1:n+2 x(1)=-1; x(i+1)=x(i)+2/(n+1); y(i)=f(x(i)); end for i=1:n d(i)=g(y(i),y(i+1),y(i+2),h); end A=zeros(n,n); for i=1:n A(i,i)=2; end for i=1:n-1 A(i,i+1)=u; A(i+1,i)=v; end B=zeros(n,1); for i=1:n B(i,1)=d(i) end C=inv(A)*B for i=1:n M(i)=C(i,1); end x=(-1:h/50:1); k=1./(1+25*x.^2); cs=spline(x,k); plot(x,k,'r.'); hold on; ezplot('1/(1+25*x^2)',[-1 1]); title('三次样条插值曲线和f(x)曲线') 四．结果

系数矩阵 弯矩M

分段函数不同次幂x对应的系数 三次样条插值函数与原函数

三次样条插值多项式matlab

三次样条插值多项式 ——计算物理实验作业四 陈万物理学2013级 主程序： clear,clc; format rat x = [1,4,9,16,25,36,49,64]; y = [1,2,3,4,5,6,7,8]; f1 = ; fn = 1/16; [a,b,c,d,M,S] = spline(x,y,f1,fn); 子程序1： function [a,b,c,d,M,S]=spline(x,y,f1,fn) % 三次样条插值函数 % x是插值节点的横坐标 % y是插值节点的纵坐标 % u是插值点的横坐标 % f1是左端点的一阶导数 % fn是右端点的一阶导数 % a是三对角矩阵对角线下边一行 % b是三对角矩阵对角线 % c是三对角矩阵对角线上边一行 % S是插值点的纵坐标

n = length(x); h = zeros(1,n-1); deltay = zeros(1,n); miu = zeros(1,n-1); lamda = zeros(1,n-1); d = zeros(1,n-1); for j = 1:n-1 h(j) = x(j+1)-x(j); deltay(j) = y(j+1)-y(j); end % 得到h矩阵 for j = 2:n-1 sumh = h(j-1) + h(j); miu(j) = h(j-1) / sumh; lamda(j) = h(j) / sumh; d(j) = 6*( deltay(j)/h(j)-(deltay(j-1)/h(j-1)))/sumh; end % 根据第一类边界条件，作如下规定 lamda(1) = 1; d(1) = 6*(deltay(1)/h(1)-f1)/h(1); miu(1) = 1; d(n) = 6*(fn-deltay(n-1)/h(n-1))/h(n-1);

三次样条插值函数的Matlab代码

并针对下面一组具体实验数据 求解，其中边界条件为. 解：Matlab计算程序为： clear clc x=[0.25 0.3 0.39 0.45 0.53] y=[0.5000 0.5477 0.6245 0.6708 0.7280] n=length(x); for i=1:n-1 h(i)=x(i+1)-x(i); end for i=1:n-2 k(i)=h(i+1)/(h(i)+h(i+1)); u(i)=h(i)/(h(i)+h(i+1)); end for i=1:n-2 gl(i)=3*(u(i)*(y(i+2)-y(i+1))/h(i+1)+k(i)*(y(i+1)-y(i))/h(i)); end g0=3*(y(2)-y(1))/h(1); g00=3*(y(n)-y(n-1))/h(n-1); g=[g0 gl g00]; g=transpose(g) k1=[k 1]; u1=[1 u]; Q=2*eye(5)+diag(u1,1)+diag(k1,-1) m=transpose(Q\g) syms X; for i=1:n-1 p1(i)=(1+2*(X-x(i))/h(i))*((X-x(i+1))/h(i))^2*y(i); p2(i)=(1-2*(X-x(i+1))/h(i))*((X-x(i))/h(i))^2*y(i+1); p3(i)=(X-x(i))*((X-x(i+1))/h(i))^2*m(i); p4(i)=(X-x(i+1))*((X-x(i))/h(i))^2*m(i+1); p(i)=p1(i)+p2(i)+p3(i)+p4(i); p(i)=simple(p(i)); end s1=p(1) s2=p(2) s3=p(3) s4=p(4) for k=1:4

数值分析作业-三次样条插值教学提纲

数值分析作业-三次样 条插值

数值计算方法作业 实验4.3 三次样条差值函数 实验目的： 掌握三次样条插值函数的三弯矩方法。 实验函数: dt e x f x t ? ∞ -- =2 221)(π 求f(0.13)和f(0.36)的近似值 实验内容： （1） 编程实现求三次样条插值函数的算法，分别考虑不同的边界条件； （2） 计算各插值节点的弯矩值； （3） 在同一坐标系中绘制函数f(x)，插值多项式，三次样条插值多项式的曲 线比较插值结果。 实验4.5 三次样条差值函数的收敛性 实验目的： 多项式插值不一定是收敛的，即插值的节点多，效果不一定好。对三次样条插值函数如何呢？理论上证明三次样条插值函数的收敛性是比较困难的，通过本实验可以证明这一理论结果。

实验内容： 按照一定的规则分别选择等距或非等距的插值节点，并不断增加插值节点的个数。 实验要求： （1） 随着节点个数的增加，比较被逼近函数和三样条插值函数的误差变化情 况，分析所得结果并与拉格朗日插值多项式比较； （2） 三次样条插值函数的思想最早产生于工业部门。作为工业应用的例子， 考虑如下例子：某汽车制造商根据三次样条插值函数设计车门曲线，其中一段数据如下： k x 0 1 2 3 4 5 6 7 8 9 10 k y 0.0 0.79 1.53 2.19 2.71 3.03 3.27 2.89 3.06 3.19 3.29 k y ' 0.8 0.2 算法描述： 拉格朗日插值： 其中 是拉格朗日基函数，其表达式为：() ∏ ≠=--=n i j j j i j i x x x x x l 0 ) ()( 牛顿插值： ) )...()(](,...,,[.... ))(0](,,[)0](,[)()(1102101210100----++--+-+=n n n x x x x x x x x x x f x x x x x x x f x x x x f x f x N