《数模插值与拟合建模知识分享》由会员分享，可在线阅读，更多相关《数模插值与拟合建模知识分享（133页珍藏版）》请在金锄头文库上搜索。
1、第八章 插值与拟合建模,重庆邮电大学 数理学院 沈世云,一.插值,主要内容,2、掌握用数学软件包求解插值问题。,1、了解插值的基本内容。,1一维插值,2二维插值,3实验作业,拉格朗日插值,分段线性插值,三次样条插值,一 维 插 值,一、插值的定义,二、插值的方法,三、用Matlab解插值问题,返回,返回,二维插值,一、二维插值定义,二、网格节点插值法,三、用Matlab解插值问题,最邻近插值,分片线性插值,双线性插值,网格节点数据的插值,散点数据的插值,返回,2.3 拉格朗日插值,1、线性插值,假设已知,在区间,上的两个结点,和它们的函数值,求一个一次多项式,，使得多项式,在结点上满足条件,这
2、种插值方法称为线性插值方法（也称两点插值）。,可以求出：,分段线性插值,计算量与n无关; n越大，误差越小.,2、抛物插值,已知,在区间,上的三个结点,和它们的函数值,求一个次数不超过2的多项式,，使得它在结点上满足条件,这种插值方法称为抛物线插值法，,可求出：,3、n次拉格朗日插值,假设取区间,上的,个结点,，并且已知函数,在这此点的函数值,现在求一个次数不超过,的多项式,，使得满足条件,这种插值方法称为,次多项式插值（或称代数插值），,利用拉格朗日插值插值方法可得,上述多项式称为 n次拉格朗日（Lagrange）插值多项式，,函数,称为拉格朗日插值基函数。,当n=1,2时，n次拉格朗日（L
3、agrange）插值多项式即为线性插值多项式和抛物插值多项式。,比分段线性插值更光滑。,在数学上，光滑程度的定量描述是：函数(曲线)的k阶导数存在且连续，则称该曲线具有k阶光滑性。 光滑性的阶次越高，则越光滑。是否存在较低次的分段多项式达到较高阶光滑性的方法？三次样条插值就是一个很好的例子。,三次样条插值,三次样条插值,g(x)为被插值函数。,例1 已知函数发f(x)的函数表如下：,求其拉格朗日插值多项式，并求,的近似值。,解 由于给出了4个插值结点，所以可做出次数不超过3的拉格朗日插值多项式。,将上列4式代入n=3 的拉格朗日插值公式，可得所要求的插值多项式为,将x=2.5代入可得f(2.5
4、)的近似值为1.8496。,拉格朗日插值法适合节点较少的情况，当节点较多的大范围高次插值的逼近效果往往并不理想且当插值结点增加时，计算越来越繁。为了提高精度和减少计算，还有牛顿插值法下、三次样条插值等，具体可参阅有关书籍。,用MATLAB作插值计算,一维插值函数：,yi=interp1(x，y，xi，method),nearest ：最邻近插值linear ： 线性插值； spline ： 三次样条插值； cubic ： 立方插值。 缺省时： 分段线性插值。,注意：所有的插值方法都要求x是单调的，并且xi不能够超过x的范围。,例：在1-12的11小时内，每隔1小时测量一次温度，测得的温度依次为
5、：5，8，9，15，25，29，31，30，22，25，27，24。试估计每隔1/10小时的温度值。,To MATLAB (temp),hours=1:12; temps=5 8 9 15 25 29 31 30 22 25 27 24; h=1:0.1:12; t=interp1(hours,temps,h,spline); (直接输出数据将是很多的) plot(hours,temps,+,h,t,hours,temps,r:) %作图 xlabel(Hour),ylabel(Degrees Celsius),例2 在用外接电源给电容器充电时，电容器两端的电压V将会随着充电时间t发生变化，已
6、知在某一次实验时，通过测量得到下列观测值，分别用拉格朗日插值法、分段线性插值法、三次样条插值法画出V随着时间t变化的曲线图，分别计算当时间t=7s时，三种插值法各自算得电容器两端电压的近似值。,解 由于MATLAB没有提供现成的拉格朗日插值命令，我们可以编写一个函数lglrcz.m来完成，其他两种插值法可用现成的命令。,用MATLAB软件进行三种插值计算的程序为szczqx.m。,程序lglrcz.m： function y=lglrcz(x0,y0,x) n=length(x0)； m=length(x) ； for i=1:m z=x(i) ； s=0.0； for k=1:n p=1.0
7、； for j=1:n if j=k p=p*(z-x0(j)/(x0(k)-x0(j) ； end end s=p*y0(k)+s； end y(i)=s； end,程序为azczf.m t=1,2,3,4,6.5,9,12 ； v=6.2,7.3,8.2,9.0,9.6,10.1,10.4 ； t0=0.2:0.1:12.5； lglr=lglrcz(t,v,t0) ； laglr=lglrcz(t,v,7) ； fdxx=interp1(t,v,t0) ； fendxx=interp1(t,v,7) ； scyt=interp1(t,v,t0,spline) ； sancyt=inter
8、p1(t,v,7,spline) plot(t,v,*,t0,lglr,r,t0,fdxx,g,t0,scyt,b) gtext(lglr) gtext(fdxx) gtext(scyt),执行结果是 laglr=9.52988980716254 fendxx=9.70000000000000 sancyt=9.67118039327444,图形如图1所示。,图中曲线lglr表示拉格朗日插值曲线，scyt表示三次样条插值曲线，fdxx表示分段线性插值曲线。,案例: 水流量的估计,美国某州的用水管理机构要求各社区提供以每小时多少加仑计的用水量以及每天所用的总水量。许多社区没有测量流入或流出水塔水
9、量的装置，只能代之以每小时测量水塔中的水位，其误差不超过5%。需要注意的是，当水塔中的水位下降到最低水位L时，水泵就自动向水塔输水直到最高水位H，此期间不能测量水泵的供水量，因此，当水泵正在输水时不容易建立水塔中水位和用水量之间的关系。水泵每天输水一次或两次，每次约2小时。,2.1 问题,试估计任何时刻（包括水泵正在输水时间）从水塔流出的水流量f(t)，并估计一天的总用水量。已知该水塔是一个高为40ft（英尺），直径为57ft（英尺）的正圆柱，表5-1给出了某个小镇一天水塔水位的真实数据，水位降至约27.00ft水泵开始工作，水位升到35.50ft停止工作。（注：1ft（英尺）=0.3048m
10、（米）,表1 某小镇某天水塔水位,2.2 问题分析,本实验所指流量可视为单位时间内流出水的体积。由于水塔是正圆柱形，横截面积是常数，所以在水泵不工作时段，流量很容易根据水位相对时间的变化算出。问题的难点在于如何估计水泵供水时段的流量。,水泵供水时段的流量只能靠供水时段前后的流量经插值或拟合得到。作为用于插值或拟合的原始数据，我们希望水泵不工作时段的流量越准确越好。这此流量大体上可由两种方法计算，一是直接对表12-1中的水量用数值微分算出各时段的流量，用它们拟合其它时刻或连续时间的流量；二是先用表中数据拟合水位一时间函数，求导数即可得到连续时间的流量。,有了任何时刻的流量，就不难计算一天的总用水
11、量。其实，水泵不工作时段的用水量可以由测量记录直接得到，由表1中下降水位乘以水搭的截面积就是这一时段的用水量。这个数值可以用来检验数据插值或拟合的结果。,在具体给出本问题的解答之前，先介绍一个简单的数据插值方法。,2.2 问题求解,为了表示方便，我们将2.1节问题中所给表1中的数据全部化为国际标准单位（表2），时间用小时（h），高度用米（m）：,表12-2 一天内水塔水位记录,1模型假设,（1）流量只取决于水位差，与水位本身无关，故由物理学中Torriceli定律：小孔流出的液体的流速正比于水面高度的平方根。题目给出水塔的最低和最高水位分别是8.1648m,和10.7352m,（设出口的水位为
12、零），因为sqrt,，约为1，所以可忽略水位对流速的影响。,（2）将流量看作是时间的连续函数，为计算简单，不妨将流量定义成单位时间流出水的高度，即水位对时间变化率的绝对值（水位是下降的），水塔截面积为,（m2），得到结果,后乘以s即可。,2流量估计方法,首先依照表2所给数据，用MATLAB作出时间水位散点图（图2）。,下面来计算水箱流量与时间的关系。,根据图1.，一种简单的处理方法为，将表2中的数据分为三段，然后对每一段的数据做如下处理：,设某段数据,，相邻数据中点的平均流速用下面的公式（流速=（右端点的水位右端点的水位）/区间长度）：,每段数据首尾点的流速用下面的公式计算：,用以上公式求得时
13、间与流速之间的数据如表3。,由表3作出时间流速散点图如图3。,1）插值法,由表12-3，对水泵不工作时段1，2采取插值方法，可以得到任意时刻的流速，从而可以知道任意时刻的流量。,我们分别采取拉格朗日插值法，分段线性插值法及三次样条插值法；对于水泵工作时段1应用前后时期的流速进行插值，由于最后一段水泵不工作时段数据太少，我们将它与水泵工作时段2合并一同进行插值处理（该段简称混合时段）。,我们总共需要对四段数据（第1，2未供水时段，第1供水时段，混合时段）进行插值处理，下面以第1未供水时段数据为例分别用三种方法算出流量函数和用水量（用水高度）。,下面是用MATLAB实现该过程的程序。,t=0,0.
14、46,1.38,2.395,3.41,4.425,12.44,6.45,7.465,8.45,8.97； v=29.89,21.74,18.48,16.22,16.30,15.32,13.04,15.45,13.98,16.35,19.29； t0=0:0.1:8.97；,lglr=lglrcz(t,v,t0)； /*注：lglrcz为一函数，程序同lglrcz.m*/ lglrjf=0.1*trapz(lglr) fdxx=interpl(t,v,t0)； fdxxjf=0.1*trapz(fdxx) scyt=interpl(t,v,t0,spline)； sancytjf=0.1*tra
15、pz(scyt) plot(t,v,*,t0,lglr,r,t0,fdxx,g,t0,scyt,b) gtext(lglr) gtext(fdxx) gtext(scyt),运行结果 lglrjf=145.6231 ； fdxxjf=147.1430； sancytjf=1412.6870,图4是对第1未供水段数据用三种不同方法得到的插值函数图，,图中曲线lglr、fdxx和scyt分别表示用拉格朗日插值法，分段线性插值法及三次样条插值法得到的曲线。,由表2知，第1未供水时段的总用水高度为146(=968822)，可见上述三种插值方法计算的结果与实际值（146）相比都比较接近。考虑到三次样条插
16、值方法具有更加良好的性质，建议采取该方法。,其他三段的处理方法与第1未供水时段的处理方法类似，这里不再详细叙述，只给出数值结果和函数图像（图5图7），图中曲线标记同图4。,图5 第一供水段时间流速示意图,图6 第2未供水段时间流速示意图,图7 混合时段时间流速示意图,图8是用分段线性及三次样条插值方法得到的整个过程的时间流速函数示意图。,表4 各时段及一天的总用水量（用水高度）,表5是对一天中任取的4个时刻分别用3种方法得到的水塔水流量近似值。,注：拉格朗日插值法分段线性插值法三次样条插值法,2）拟合法,（1）拟合水位时间函数,从表12-2中的测量记录看，一天有两次供水时段和三次未供水时段，分别对第1，2未供水时段的测量数据直接作多项式拟合，可得到水位函数（注意，根据多项式拟合的特点，此处拟合多项式的次数不宜过高，一般以36次为宜）。对第3未供水时段来说，数据过少不能得到很好的拟合。,设t，h分别为已输入的时刻和水位测量记录（由表12.2提供，水泵启动的4个时刻
