插值与拟合方法建模0727
- 格式:ppt
- 大小:991.50 KB
- 文档页数:59
数值计算方法插值与拟合数值计算方法在科学计算和工程应用中起着重要的作用,其中插值和拟合是其中两个常用的技术。
插值是指通过已知的离散数据点来构造出连续函数或曲线的过程,拟合则是找到逼近已知数据的函数或曲线。
本文将介绍插值和拟合的基本概念和常见的方法。
一、插值和拟合的基本概念插值和拟合都是通过已知数据点来近似表达未知数据的方法,主要区别在于插值要求通过已知数据点的函数必须经过这些数据点,而拟合则只要求逼近这些数据点。
插值更加精确,但是可能会导致过度拟合;拟合则更加灵活,能够通过调整参数来平衡拟合精度和模型复杂度。
二、插值方法1. 线性插值线性插值是一种简单的插值方法,通过已知数据点构造出线段,然后根据插值点在线段上进行线性插值得到插值结果。
2. 拉格朗日插值拉格朗日插值是一种基于多项式插值的方法,通过已知数据点构造出一个多项式,并根据插值点求解插值多项式来得到插值结果。
3. 分段线性插值分段线性插值是一种更加灵活的插值方法,通过将插值区间分成若干小段,然后在每个小段上进行线性插值。
三、拟合方法1. 最小二乘法拟合最小二乘法是一种常用的拟合方法,通过最小化实际观测点和拟合函数之间的残差平方和来确定拟合函数的参数。
2. 多项式拟合多项式拟合是一种基于多项式函数的拟合方法,通过选择合适的多项式次数来逼近已知数据点。
3. 曲线拟合曲线拟合是一种更加灵活的方法,通过选择合适的曲线函数来逼近已知数据点,常见的曲线包括指数曲线、对数曲线和正弦曲线等。
四、插值与拟合的应用场景插值和拟合在实际应用中具有广泛的应用场景,比如图像处理中的图像重建、信号处理中的滤波器设计、金融中的风险评估等。
五、插值与拟合的性能评价插值和拟合的性能可以通过多种指标进行评价,常见的评价指标包括均方根误差、相关系数和拟合优度等。
六、总结插值和拟合是数值计算方法中常用的技术,通过已知数据点来近似表达未知数据。
插值通过已知数据点构造出连续函数或曲线,拟合则找到逼近已知数据的函数或曲线。
数学建模插值与拟合数据插值与拟合插值与插值函数:已知由(可能未知或⾮常复杂)产⽣的⼀批离散数据,且个互异插值节点,在插值区间内寻找⼀个相对简单的函数,使其满⾜下列插值条件:再利⽤已求得的计算任⼀⾮插值节点的近似值,这就是插值。
其中称为插值函数,称为被插函数。
最⼩⼆乘拟合:已知⼀批离散的数据,互不相同,寻求⼀个拟合函数,使与的误差平⽅和在最⼩⼆乘意义下最⼩。
在最⼩⼆乘意义下确定的称为最⼩⼆乘拟合函数。
1)Lagrange插值法a.待定系数法:假设插值多项式,利⽤待定系数法即可求得满⾜插值条件的插值函数。
关键在于确定待定系数。
b.利⽤基函数的构造⽅法⾸先构造个满⾜条件:的次插值基函数,再将其线性组合即可得如下的Lagrange插值多项式:其中c.Lagrange插值余项注:上述两种构造⽅法所得的Lagrange插值多项式是⼀样的,即满⾜插值条件的Lagrange插值多项式是唯⼀的。
2)分段线性插值作分段线性插值的⽬的在于克服Lagrange插值⽅法可能发⽣的不收敛性缺点。
所谓分段线性插值就是利⽤每两个相邻插值节点作线性插值,即可得如下分段线性插值函数:其中特点:插值函数序列具有⼀致收敛性,克服了⾼次Lagrange插值⽅法的缺点,故可通过增加插值节点的⽅法提⾼其插值精度。
但存在于节点处不光滑、插值精度低的缺点。
3)三次样条插值三次样条插值的⽬的在于克服Lagrange插值的不收敛性和提⾼分段线性插值函数在节点处的光滑性。
所谓三次样条插值⽅法就是在满⾜下列条件:a.b.在每个⼦区间上是三次多项式的三次样条函数中寻找满⾜如下插值条件:以及形如等边界条件的插值函数的⽅法。
特点:三次样条插值函数序列⼀致收敛于被插函数,因此可通过增加节点的⽅法提⾼插值的精度。
4)插值⽅法的Matlab实现⼀维数据插值MATLAB中⽤函数interp1来拟合⼀维数据,语法是YI = INTERP1(X,Y,XI,⽅法)其中(X,Y)是已给的数据点,XI 是插值点,其中⽅法主要有'linear' -线性插值,默认'pchip' -逐段三次Hermite插值'spline' -逐段三次样条函数插值其中最后⼀种插值的曲线⽐较平滑例:x=0:.12:1; x1=0:.02:1;y=(x.^2-3*x+5).*exp(-5*x).*sin(x);plot(x,y,'o'); hold on;y1=interp1(x,y,x1,'spline');plot(x1,y1,':')如果要根据样本点求函数的定积分,⽽函数⼜是⽐较光滑的,则可以⽤样条函数进⾏插值后再积分,在MATLAB 中可以编写如下程序:function y=quadspln(x0,y0,a,b)f=inline(‘interp1(x0,y0,x,’’spline’’)’,’x’,’x0’,’y0’);y=quadl(f,a,b,1e-8,[],x0,y0);现求six(x)在区间[0,pi]上的定积分,只取5点x0=[0,0.4,1,2,pi];y0=sin(x0);I=quadspln(x0,y0,0,pi)结果得到的值为 2.01905,精确值为2⼆元函数插值:MATLAB中⽤函数interp2来拟合⼆维⽹格(X,Y)上的数据Z,语法是YI = INTERP2(X,Y, Z,XI, YI,⽅法)其中(X,Y,Z)是已给的数据点,(XI,YI)是插值点坐标,其中⽅法主要有'linear' -线性插值,默认'pchip' -逐段三次Hermite插值'spline' -逐段三次样条函数插值其中最后⼀种插值的曲⾯⽐较平滑例:[x,y]=meshgrid(-3:.6:3,-2:.4:2);z=(x.^2-2*x).*exp(-x.^2-y.^2-x..*y);[x1,y1]=meshgrid(-3:.2:3,-2:.2:2);%⽣成⽹格,x1和y1均为同样size的矩阵z1=interp2(x,y,z,x1,y1,’spline’); %z1是矩阵,size 和x1,y1相同surf(x1,y1,z1);axis([-3,3,-2,2,-0.7,1.5]);-33如果数据不是在⽹格上取的,则可⽤函数griddata 来解决语法是YI = griddata(X,Y, Z ,XI, YI ,‘v4’)其中(X , Y ,Z )是已给的数据点,(XI ,YI )是插值点坐标,其中除了⽅法‘v4’外还有 'linear' -线性插值,默认 'cublc' -逐段三次Hermite 插值 'nearest' 其中‘v4’⽅法⽐较好例x=-3+6*rand(200,1); %⽣成随机点的x坐标向量xy=-2+4*rand(200,1); %⽣成随机点的y坐标向量yz=(x.^2-2*x).*exp(-x.^2-y.^2-x.*y); % 上述点的样本值向量z[x1,y1]=meshgrid(-3:.2:3,-2:.2:2); %⽣成⽹格,x1和y1均为同样size的矩阵z1=griddata(x,y,z,x1,y1,’v4’);surf(x1,y1,z1);axis([-3,3,-2,2,-0.7,1.5]);⽣成的图类似上图。
插值和拟合【1 】试验目标:懂得数值剖析建模的办法,控制用Matlab进行曲线拟合的办法,懂得用插值法建模的思惟,应用Matlab一些敕令及编程实现插值建模.试验请求:懂得曲线拟合和插值办法的思惟,熟习Matlab相干的敕令,完成响应的演习,并将操纵进程.程序及成果记载下来.试验内容:一.插值1.插值的根本思惟·已知有n +1个节点(xj,yj),j = 0,1,…, n,个中xj互不雷同,节点(xj, yj)可算作由某个函数 y= f(x)产生;·结构一个相对简略的函数y=P(x);·使P经由过程全体节点,即 P (xk) = yk,k=0,1,…, n ;·用P (x)作为函数f ( x )的近似.2.用MA TLAB作一维插值盘算yi=interp1(x,y,xi,'method')注:yi—xi处的插值成果;x,y—插值节点;xi—被插值点;method—插值办法(‘nearest’:最临近插值;‘linear’:线性插值;‘spline’:三次样条插值;‘cubic’:立方插值;缺省时:线性插值).留意:所有的插值办法都请求x是单调的,并且xi不克不及够超出x的规模.演习1:机床加工问题机翼断面下的轮廓线上的数据如下表:x 0 3 5 7 9 11 12 13 14 15y 0用程控铣床加工机翼断面的下轮廓线时每一刀只能沿x偏向和y偏向走异常小的一步.表3-1给出了下轮廓线上的部分数据但工艺请求铣床沿x偏向每次只能移动单位.这时需求出当x 坐标每转变单位时的y 坐标. 试完成加工所需的数据,画出曲线. 步调1:用x0,y0两向量暗示插值节点;步调2:被插值点x=0:0.1:15; y=y=interp1(x0,y0,x,'spline'); 步调3:plot(x0,y0,'k+',x,y,'r')grid on答:x0=[0 3 5 7 9 11 12 13 14 15 ]; y0=[0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6 ]; x=0:0.1:15;y=interp1(x0,y0,x,'spline'); plot(x0,y0,'k+',x,y,'r') grid on0510150.511.522.53.用MA TLAB 作网格节点数据的插值(二维)z=interp2(x0,y0,z0,x,y,’method’)注:z—被插点值的函数值;x0,y0,z0—插值节点;x,y—被插值点;method—插值办法(‘nearest’:最临近插值;‘linear’:双线性插值; ‘cubic’:双三次插值;缺省时:双线性插值).留意:请求x0,y0单调;x,y可取为矩阵,或x取行向量,y取为列向量,x,y的值分离不克不及超出x0,y0的规模.4.用MA TLAB作散点数据的插值盘算cz =griddata(x,y,z,cx,cy,‘method’)注:cz—被插点值的函数值;x,y,z—插值节点;cx,cy—被插值点;method—插值办法(‘nearest’:最临近插值;‘linear’:双线性插值; ‘cubic’:双三次插值;'v4‘:Matlab供给的插值办法;缺省时:双线性插值).演习2:航行区域的警示线某海域上频仍地有各类吨位的船只经由.为包管船只的航行安然,有关机构在低潮时对水深进行了测量,下表是他们供给的测量数据:水道水深的测量数据x 129.0140.0 103.5 88.0 185.5 195.0 105.5y 7.5 141.5 23.0 147.0 22.5 137.5 85.5z 4 8 6 8 6 8 8x157.5 107.5 77.0 81.0 162.0 162.0 117.5y -6.5 -81.0 3.0 56.5 -66.5 84.0 -33.5z 9 9 8 8 9 4 9个中(x, y)为测量点,z为(x, y)处的水深(英尺),水深z是区域坐标(x, y)的函数z= z (x, y),船的吨位可以用其吃水深度来反应,分为4英尺.英尺.5英尺和英尺 4 档.航运部分要在矩形海域(75,200)×(-50,150)上为不合吨位的航船设置警示标识表记标帜.请依据测量的数据描写该海域的地貌,并绘制不合吨位的警示线,供航运部分应用. x=[129 140 103.5 88 185.5 195 105.5 157.5 107.5 77 81 162 162 117.5];y=[7.5 141.5 23 147 22.5 137.5 85.5 -6.5 -81 3 56.5 -66.5 84 -33.5];z=[-4 -8 -6 -8 -6 -8 -8 -9 -9 -8 -8 -9 -4 -9];cx=75:0.5:200;cy=-70:0.5:150;cz=griddata(x,y,z,cx,cy','cubic');meshz(cx,cy,cz),rotate3dxlabel('X'),ylabel('Y'),zlabel('Z')%pausefigure(2),contour(cx,cy,cz,[-5 -5]);grid on,hold onplot(x,y,'+')xlabel('X'),ylabel('Y')200XYZXY80100120140160180200-60-40-20020406080100120140演习3:估量水塔的水流量—93,请绘出三次样条插值曲线,并盘算一天的总的用水量. 解:t0=[0.46,1.38,2.4,3.41,4.43,5.44,6.45,7.47,8.45,11.49,12.49,13.42,14.43,15.44,16.37,17.38,18.49,19.50,20.40,24.43,25.32];v0=[11.2,9.7,8.6,8.1,9.3,7.2,7.9,7.4,8.4,15.6,16.4,15.5,13.4,13.8,12.9,12.2,12.2,12.9,12.6,11.2,3.5]; t=0:0.1:26; y=interp1(t0,v0,t,'spline'); plot(t0,v0,'k+',t,y,'r') grid on0510********-10-55101520二.曲线拟合已知一组(二维)数据,即平面上 n 个点(xi,yi) i=1,…n, 追求一个函数(曲线)y=f(x), 使 f(x) 在某种准则下与所稀有据点最为接近,即曲线拟合得最好.最经常应用的办法是线性最小二乘拟合 1.多项式拟合⏹对给定的数据(xj,yj),j = 0,1,…, n;⏹拔取恰当阶数的多项式,如二次多项式g(x)=ax^2+bx+c;⏹使g(x)尽可能逼近(拟合)这些数据,但是不请求经由给定的数据(xj,yj); 2.多项式拟合指令1)多项式f(x)=a1xm+ …+amx+am+1拟合指令:a=polyfit(x,y,m)a:输出多项式拟合系数a[a1,a2,…,am];x,y:输出长度雷同的数组;m:多项式的次数. 2)多项式在x处的值y的盘算敕令:y=polyval(a,x)演习4:对下面一组数据作二次多项式拟合写出拟合敕令:plot(x,y,'k+',x,z,'r')作出数据点和拟合曲线:0.10.20.30.40.50.60.70.80.91写出拟合的二次多项式:0317.01293.208108.9)(2-+-=x x x f3.可化为多项式的非线性拟和曲线改直是工程中又一经常应用的断定曲线情势的办法,很多罕有的函数都可以经由过程恰当的变换转化为线性函数.(1)幂函数 by ax c =+ln ln ln y c a b x -=+(2)指数函数 xy ab c =+ln ln ln y c a x b -==(3)抛物函数 2,(0)y ax bx c x =++≠b ax xcy +=- 演习5:完成教材P93页的习题5的第一小题. x0=[0,300,600,1000,1500,2000];x=0:100:2000;y0=[0.9689,0.9322,0.8969,0.8519,0.7989,0.7491];y=interp1(x0,y0,x,'spline');plot(x0,y0,'k+',x,y,'r')grid on0200400600800100012001400160018002000。
插值与拟合方法在实际中,常常要处理由实验或测量所得到的一批离散数据.插值与拟合方法就是要通过这些数据去确定某一类已知函数的参数或寻找某个近似函数,使所得到的近似函数与已知数据有较高的拟合精度.插值问题:要求这个近似函数(曲线或曲面)经过所已知的所有数据点.通常插值方法一般用于数据较少的情况.数据拟合:不要求近似函数通过所有数据点,而是要求它能较好地反映数据的整体变化趋势。
共同点:插值与拟合都是根据实际中一组已知数据来构造一个能够反映数据变化规律的近似函数的方法,由于对近似要求的准则不同,因此二者在数学方法上有很大的差异.插值问题的一般提法:已知某函数)(x f y =(未知)的一组观测(或试验)数据),,2,1)(,(n i y x ii⋅⋅⋅=,要寻求一个函数)(x φ,使iiy x =)(φ),,2,1(n i ⋅⋅⋅=,则)()(x f x ≈φ.实际中,常常在不知道函数)(x f y =的具体表达式的情况下,对于i x x =有实验测量值iy y =),,2,1,0(n i ⋅⋅⋅=,寻求另一函数)(x φ使满足:)()(i i i x f y x ==φ),,2,1,0(n i ⋅⋅⋅=称此问题为插值问题,并称函数)(x φ为)(x f 的插值函数,nx x x x ,,,,21⋅⋅⋅称为插值节点,),,2,1,0()(n i y x ii⋅⋅⋅==φ称为插值条件,即)()(iiix f y x ==φ),,2,1,0(n i ⋅⋅⋅=,则)()(x f x ≈φ.(1) 拉格朗日(Lagrange )插值设函数)(x f y =在1+n 个相异点nx x x x ,,,,21⋅⋅⋅上的函数值为ny y y y ,,,,21⋅⋅⋅,要求一个次数不超过n 的代数多项式nnnx a x a x a a x P +⋅⋅⋅+++=221)(使在节点i x 上有),,2,1,0()(n i y x P ii n ⋅⋅⋅==成立,称之为n 次代数插值问题,)(x P n称为插值多项式.可以证明n 次代数插值是唯一的.事实上: 可以得到j n j n i i j in y x x xx x P j i ∑∏==⎥⎥⎦⎤⎢⎢⎣⎡⎪⎪⎭⎫⎝⎛--=≠00)()( 当1=n 时,有二点一次(线性)插值多项式:101001011)(y x x x x y x x x x x P --+--=当n =2时,有三点二次(抛物线)插值多项式:2120210121012002010212))(())(())(())(())(())(()(y x x x x x x x x y x x x x x x x x y x x x x x x x x x P ----+----+----=(2)牛顿(Newton ) 插值牛顿插值的基本思想:由于)(x f y =关于二节点10,x x 的线性插值为)()()()()()()()()(00101000010101x x x x x f x f x p x x x x x f x f x f x p ---+=---+= 假设满足插值条件)2,1,0()()(2===i x p y x f iii的二次插值多项式一般形式为))(()()(1212x x x x c x x c c x p --+-+= 由插值条件可得⎪⎩⎪⎨⎧=--+-+=-+=)())(()()()()(21202202101011000x f x x x x c x x c c x f x x c c x f c 可以解出⎪⎪⎪⎩⎪⎪⎪⎨⎧------=--==020101121220101100)()()()()()(),(x x x x x f x f x x x f x f c x x x f x f c x f c所以))(()())(()()(10211020102x x x x c x p x x x x c x x c c x p --+=--+-+=类似的方法,可以得到三次插值多项式等,按这种思想可以得到一般的牛顿插值公式.函数的差商及其性质对于给定的函数)(x f ,用),,,(10n x x x f ⋅⋅⋅表示关于节点nx x x ,,,1⋅⋅⋅的n 阶差商,则有一阶差商:01011)()(),(x x x f x f x x f --=,121221)()(),(x x x f x f x x f --= 二阶差商:021021210),(),(),,(x x x x f x x f xx x f --=n 阶差商:0110211),,,(),,,(),,,(x x x x x f x x x f x x x f n n n n -⋅⋅⋅-⋅⋅⋅=⋅⋅⋅-差商有下列性质:(1)差商的分加性:∑∏=≠=-=⋅⋅⋅nk nk j j j kk n x xx f xx x f 0)(01)()(),,,(.(2)差商的对称性:在),,,(1nx x x f ⋅⋅⋅中任意调换jix x ,的次序其值不变.牛顿插值公式: 一次插值公式为))(,()()(01001x x x x f x f x p -+=二次插值公式为))()(,,()())()(,,())(,()()(1021011021001002x x x x x x x f x p x x x x x x x f x x x x f x f x p --+=--+-+=于是有一般的牛顿插值公式为)())()(,,,()()())()(,,,())()(,,())(,()()(11010111010102100100----⋅⋅⋅--⋅⋅⋅+=-⋅⋅⋅--⋅⋅⋅+⋅⋅⋅+--+-+=n n n n n n x x x x x x x x x f x p 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 p可以证明:其余项为))(())()(,,,,()(11010n n n x x x x x x x x x x x x f x R --⋅⋅⋅--⋅⋅⋅=-实际上,牛顿插值公式是拉格朗日插值公式的一种变形,二者是等价的.另外还有著名的埃尔米特(Hermite )插值等.(3)样条函数插值方法样条,实质上就是由分段多项式光滑连接而成的函数,一般称为多项式样条.由于样条函数的特殊性质,决定了样条函数在实际中有着重要的应用.样条函数的一般概念定义 设给定区间],[b a 的一个分划b x x x a n=<⋅⋅⋅<<=∆1:,如果函数)(x s 满足条件:(1) 在每个子区间),,2,1](,[1n i x x ii ⋅⋅⋅=-上是k 次多项式; (2) )(x s 及直到k -1阶的导数在],[b a 上连续.则称)(x s 是关于分划△的一个k 次多项式样条函数,nx x x ,,,1⋅⋅⋅称为样条节点,121,,,-⋅⋅⋅n x x x 称为内节点,nx x ,0称为边界节点,这类样条函数的全体记作),(k S P∆,称为k 次样条函数空间.若),()(k S x s P∆∈,则)(x s 是关于分划△的k 次多项式样条函数.k 次多项式样条函数的一般形式为∑∑=-=+-+=ki n j k j jii k x x k i x x s 011)(!!)(βα其中),,1,0(k i i=α和)1,,2,1(-=n j jβ均为任意常数,而)1,,2,1(,0,)()(-=⎪⎩⎪⎨⎧<≥-=-+n j x x x x x x x x jj kj kj在实际中最常用的是2=k 和3的情况,即为二次样条函数和三次样条函数. 二次样条函数:对于],[b a 上的分划b x x x a n=<⋅⋅⋅<<=∆1:,则)2,()(!2!2)(11222102∆βαααP n j j jS x x x x x s ∈-+++=∑-=+其中)1,2,1(,0,)()(22-=⎪⎩⎪⎨⎧<≥-=-+n j x x x x x x x x j j j j . 三次样条函数:对于],[b a 上的分划b x x xa n =<⋅⋅⋅<<=∆10:,则)3,()(!3!3!2)(1133322103∆βααααP n j j jS x x x x x x s ∈-++++=∑-=+其中)1,2,1(,0,)()(33-=⎪⎩⎪⎨⎧<≥-=-+n j x x x x x x x x jjj j .1 二次样条函数插值)2,()(2∆∈P S x s 中含有2+n 个待定常数,故应需要2+n 个插值条件,因此,二次样条插值问题可分为两类:问题(1):已知插值节点ix 和相应的函数值),,2,1,0(n i y i⋅⋅⋅=,以及端点0x (或n x )处的导数值0'y (或ny '),求)2,()(2∆∈PS x s 使得⎩⎨⎧'=''='⋅⋅⋅==))(()(),,2,1,0()(20022n n i i y x s y x s n i y x s 或(5.1)问题(2):已知插值节点ix 和相应的导数值),,2,1,0(n i y i⋅⋅⋅=',以及端点0x (或n x )处的函数值0y (或ny ),求)2,()(2∆∈P S x s 使得⎩⎨⎧==⋅⋅⋅='='))(()(),,2,1,0()(20022n n i i y x s y x s n i y x s 或(5.2)事实上,可以证明这两类插值问题都是唯一可解的.对于问题(1),由条件(5.1)⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧'=+='==-+++==++==++=∑-=00210211222102121211112020201002)(,,3,2,)(2121)(21)(21)(y x x s n j y x x x x x s yx x x s y x x x s j j i i j i jj j ααβααααααααα 引入记号T n ),,,,,(11210-=ββααα X 为未知向量,T nn y y y y ),,,,(10'= C 为已知向量, ⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡---=-0010)(21)(21211)(212110211211021212212222211200x x x x x x x x x x x x x x x n n n n n A 于是,问题转化为求方程组C AX =的解Tn ),,,,,(1121-=ββααα X 的问题,即可得到二次样条函数的)(2x s 的表达式.对于问题(2)的情况类似.2.三次样条函数插值由于)3,()(3∆∈P S x s 中含有3+n 个待定系数,故应需要3+n 个插值条件,因此可将三次样条插值问题分为三类: 问题(1):已知插值节点jx 和相应的函数值),,2,1,0(n j y j⋅⋅⋅=,以及两个端点0x ,n x 处的导数值0'y ,ny ',求)3,()(3∆∈PS x s 使满足条件⎪⎩⎪⎨⎧='='⋅⋅⋅==),0()(),,1,0()(33n j y x s n j y x s j j j j(5.3)问题(2):已知插值节点jx 和相应的函数值),,2,1,0(n j y j⋅⋅⋅=,以及两个端点0x ,nx 处的二阶导数值0y '',n y '',求)3,()(3∆∈PS x s 使满足条件⎪⎩⎪⎨⎧=''=''⋅⋅⋅==),0()(),,1,0()(33n j y x s n j y x s j j j j(5.4)问题(3):类似地,求)3,()(3∆∈PSx s 使满足条件⎪⎩⎪⎨⎧=+=-==)2,1,0)(0()0(),,1,0()(0)(3)(33k x s x s n j y x s k n k j j(5.5)这三类插值问题的条件都是3+n 个,可以证明其解都是唯一的〔8〕.一般的求解方法可以仿照二次样条的情况处理方法,在这里给出一种更简单的方法.仅依问题(1)为例,问题(2)和问题(3)的情况类似处理.由于在)3,()(3∆PS x s ∈区间],[b a 上是一个分段光滑,且具有二阶连续导数的三次多项式,则在子区间],[1+j jx x 上)(3x s ''是线性函数,记),,,1,0)((3n j x s d jj =''=为待定常数.由拉格朗日插值公式可得nj x x h h x x d h x x d x s j j j jj j jj j ,,1,0,,)(1113=-=-+-=''+++显然jjj h d dx s -='''+13)(在],[1+j jx x上为常数.于是在],[1+j j x x 上有31233)(6)(2))(()(j jjj j j j j j x x h d d x x d x x x s y x s --+-+-'+=+(5.6)则当1+=j x x 时,由(5.6)式和问题(1)的条件得121231362)()(+++=-++'+=j j jj j j j j j j y h d d h d h x s y x s故可解得)2(6)(113+++--='j j j jjj j d d h h y y x s(5.7)将(5.7)式代入(5.6)式得)1,,1,0](,[,)(6)(2)()2(6)(1312113-=∈--+-+-⎥⎥⎦⎤⎢⎢⎣⎡+--+=++++n j x x x x x h d d x x d x x d d h h y y y x s j j j jj j j jj j j j j j j j(5.8) 在],[1j j x x-上同样的有),,2,1](,[,)(6)(2)()2(6)(131112111111113n j x x x x x h d d x x d x x d d h h y y y x s j j j j j j j j j j j j j j j j =∈--+-+-⎥⎥⎦⎤⎢⎢⎣⎡+--+=------------(5.9) 根据)(3x s的一阶导数连续性,由(5.9)式得)()2(6)0(311113j j j j j j j j x s d d h h y y x s '=++-=-'---- 结合(5.7)式整理得⎪⎪⎭⎫ ⎝⎛---+=++++--+-+----11111111162j j j j j j j j j j j j j j j j j h y y h y y h h d h h h d d h h h 引入记号⎪⎪⎭⎫ ⎝⎛---+=+=--+--111116,j j j j j j j j j j j j j h y y h y y h h c h h h a ,111--+=-j j j j h h h a .则)1,,2,1(,2)1(11-==++-+-n j c d a d d a j j j j j j(5.10)再由边界条件:nny x s y x s '=''=')(,)(33得⎪⎪⎩⎪⎪⎨⎧⎪⎪⎭⎫ ⎝⎛--'=+⎪⎪⎭⎫ ⎝⎛'--=+----111100010106262n n n n n n n h y y y h d d y h y y h d d(5.11)联立(5.10),(5.11)式得方程组C D A =⋅(5.12)其中⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡---=----2121212112112200n n n n a a a a a aA ,⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=-n n d d d d 110 D ,⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎪⎪⎭⎫ ⎝⎛--'⎪⎪⎭⎫ ⎝⎛'--=----111110001066n n n n n n hy y y h c c y h y y h C 由方程组(6.12)可以唯一解出),,1,0(n j d j=,代入(5.8)式就可以得三次样条函数)(3x s 的表达式.B样条函数插值方法磨光函数实际中的许多问题,往往是既要求近似函数(曲线或曲面)有足够的光滑性,又要求与实际函数有相同的凹凸性,一般插值函数和样条函数都不具有这种性质.如果对于一个特殊函数进行磨光处理生成磨光函数(多项式),则用磨光函数构造出样条函数作为插值函数,既有足够的光滑性,而且也具有较好的保凹凸性,因此磨光函数在一维插值(曲线)和二维插值(曲面)问题中有着广泛的应用.由积分理论可知,对于可积函数通过积分会提高函数的光滑度,因此,我们可以利用积分方法对函数进行磨光处理.定义 若)(x f 为可积函数,对于0>h ,则称积分⎰+-=22,1)(1)(hx h x h dt t f h x f为)(x f 的一次磨光函数,h 称为磨光宽度.同样的,可以定义)(x f 的k 次磨光函数为)1()(1)(22,1,>=⎰+--k dt t f h x f hx h x h k h k事实上,磨光函数)(,x f h k 比)(x f 的光滑程度要高,且当磨光宽度h 很小时)(,x f h k 很接近于)(x f .等距B样条函数对于任意的函数)(x f ,定义其步长为1的中心差分算子δ如下:⎪⎭⎫ ⎝⎛--⎪⎭⎫ ⎝⎛+=2121)(x f x f x f δ在此取0)(+=x x f ,则002121+++⎪⎭⎫ ⎝⎛--⎪⎭⎫ ⎝⎛+=x x x δ是一个单位方波函数(如图5-1),记0)(+=Ωx x δ.并取1=h ,对)(0x Ω进行一次磨光得++++-+++-+++--+-+=-=⎥⎥⎦⎤⎢⎢⎣⎡⎪⎭⎫ ⎝⎛--⎪⎭⎫ ⎝⎛+==⎰⎰⎰⎰)1(2)1(2121)()(11212100212101x x x dt t dt t dt t t dt t x x xx x x x x x ΩΩ显然)(1x Ω是连续的(如图5-2).)(1x Ωo1-1/2 0 1/2 x -1 0 1 x 图5-1图5-2类似地可得到k 次磨光函数为kk j jk j k j k x k C x ++=+⎪⎭⎫ ⎝⎛-++-=Ω∑21!)1()(11 实际上,可以证明:)(x kΩ是分段k 次多项式,且具有1-k 阶连续导数,其k 阶导数有2+k个间断点,记为)1,,2,1,0(21+⋅⋅⋅=+-=k j k j x j.从而可知)(x kΩ是对应于分划+∞<<⋅⋅⋅<<<-∞∆+110:k x x x 的k 次多项式样条函数,称之为基本样条函数,简称为k 次B样条.由于样条节点为)1,,2,1,0(21+⋅⋅⋅=+-=k j k j xj是等距的,故)(x k Ω又称为k 次等距B样条函数.对于任意函数)(x f 的k 次磨光函数,由归纳法可以得到 [4,8] :⎪⎭⎫⎝⎛+≤≤--Ω=⎰∞+∞--22)()(1)(1,h x t h x dt t f htx h x f k h k 特别地,当1)(=x f 时,有1)(11⎰+∞∞--=-dt htx hk Ω,从而1)(⎰+∞∞-=dx x k Ω,且当k ≥1时有递推关系⎥⎦⎤⎢⎣⎡⎪⎭⎫ ⎝⎛-Ω⎪⎭⎫ ⎝⎛---⎪⎭⎫ ⎝⎛+Ω⎪⎭⎫ ⎝⎛++=Ω--212121211)(11x x k x k x k x k k k一维等距B样条函数插值等距B样条函数与通常的样条如下的关系: 定理设有区间],[b a 的均匀分划nab h n j jh x x j -=⋅⋅⋅=+=),,,1,0(:0∆,则对任意 k 次样条函数),()(k S x S p k ∆∈都可以表示为B样条函数族1021-=-=⎭⎬⎫⎩⎨⎧⎪⎭⎫⎝⎛+---n j k j k k j h x x Ω的线性组合[14].根据定理 5.1,如果已知曲线上一组点()jjy x ,,其中),,1,0,0(0n j h jh x x j⋅⋅⋅=>+=,则可以构造出一条样条磨光曲线(即为B样条函数族的线性组合)⎪⎭⎫⎝⎛--=∑--=j h x x c x S n kj k j k 01)(Ω 其中)1,,1,(-⋅⋅⋅+--=n k k j c j为待定常数.用它来逼近曲线,既有较好的精度,又有良好的保凸性.实际中,最常用的是3=k 的情况,即一般形式为⎪⎭⎫ ⎝⎛--=∑+-=j h x x c x S n j j 01133)(Ω 其中3+n 个待定系数)1,,0,1(+⋅⋅⋅-=n j c j可以由三类插值条件确定.由插值条件(5.3)得()()()⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧'=-'='==-='=-'='∑∑∑+-=+-=+-=n n j j n i n j j i n j j y j n c h x S ni y j i c x S y j c h x S 113311330113031)(,,1,0,)(1)(ΩΩΩ(5.13)注意到)(3x Ω的局部非零性及其函数值:61)1(,32)0(33=±=ΩΩ,当2≥x 时0)(3=x Ω;且由)21()21()(223--+='x x x ΩΩΩ知,21)1(,0)0(33=±'='ΩΩ,当2≥x 时0)(3='x Ω.则(5.13)中的每一个方程中只有三个非零系数,具体的为⎪⎩⎪⎨⎧'=+-==++'=+-+-+--n n n i i i i y h c c n i y c c c y h c c 2,,1,0,6421111011(5.14)由方程组(5.14)容易求解出)1,,0,1(+⋅⋅⋅-=n j c j,即可得到三次样条函数)(3x S 表达式.类似地,由插值条件(5.4)得待定系数的)1,,0,1(+⋅⋅⋅-=n j c j所满足的方程组为⎪⎩⎪⎨⎧''=+-==++''=+-+-+--nn n n i i i i y h c c c n i y c c c y h c c c 21111021012,,1,0,642(5.15)由插值条件(5.5)得待定系数的)1,,0,1(+⋅⋅⋅-=n j cj所满足的方程组为⎪⎪⎩⎪⎪⎨⎧==++=-+---=-++-=-+-+-+-+--+--+--ni y c c c c c c c c c c c c c c c c c c c i i i i n n n n n n n n ,,1,0,640)()(2)(0)(0)(0)()(4)(1111011111111011(5.16)方程组(5.15),(5.16)也都是容易求解的.注:上述等距B样条插值公式也适用于近似等距的情形,但在端点0x 和n x 处误差可能较大,实际应用时,为了提高在端点0x 和nx 处的精度,可以适当向左右延拓几个节点.二维等距B样条函数插值设有空间曲面),(y x f z =(未知),如果已知二维等距节点()()τj y ih x y x ji++=0,,)0,(>τh 上的值为),,2,1,0;,,2,1,0(m j n i z ij⋅⋅⋅=⋅⋅⋅=,则相应的B样条磨光曲面的一般形式为⎪⎭⎫ ⎝⎛--⎪⎭⎫⎝⎛--=∑∑--=--=j y y i h x x c y x s l m lj k ij n ki τΩΩ0011),( 其中),,2,1,0;,,2,1,0(m j n i c ij⋅⋅⋅=⋅⋅⋅=为待定常数,l k ,可以取不同值,常用的也是2,=l k 和3的情形.这是一种具有良好保凸性的光滑曲面(函数),在工程设计中是常用的,但只能使用于均匀分划或近似均匀分划的情况.(4) 最小二乘拟合方法最小二乘拟合方法的思想:由于一般插值问题并不总是可解的(即当插值条件多于待定系数的个数时,其问题无解),同时,问题的插值条件本身一般是近似的,为此,只要求在节点上近似地满足插值条件,并使它们的整体误差最小,这就是最小二乘拟合法.最小二乘拟合方法可以分为线性最小二乘拟合方法和非线性最小二乘拟合方法.线性最小二乘拟合方法设{}m k kx 0)(=φ是一个线性无关的函数系,则称线性组合∑==mk k k x a x 0)()(φφ为广义多项式.如三角多项式:∑∑==+=mk k mk kkx b kx ax 0sin cos )(φ.设由给定的一组测量数据),(iiy x 和一组正数),,2,1(n i w i⋅⋅⋅=,求一个广义多项式∑==mk k k x a x 0)()(φφ使得目标函数[]21)(∑=-=ni i i i y x w S φ(5.17)达到最小,则称函数)(x φ为数据),,2,1)(,(n i y x ii⋅⋅⋅=关于权系数),,2,1(n i w i⋅⋅⋅=的最小二乘拟合函数,由于)(x φ关于待定系数ia 是线性的,故此问题又称为线性最小二乘问题. 注意:这里{}m k kx 0)(=φ可根据实际来选择,权系数iw 的选取更是灵活多变的,有时可选取1=i w ,或nw i 1=,对于nw i1=,则相应问题称为均方差的极小化问题.最小二乘拟合函数的求解要使最小二乘问题的目标函数(5.17)达到最小,则由多元函数取得极值的必要条件得),,2,1,0(0m k a Sk==∂∂ 即),,2,1,0(0)()(10m k x y x a w i k ni i m k i k k i ⋅⋅⋅⋅==⎥⎦⎤⎢⎣⎡-∑∑==φφ 亦即),,2,1,0()()()(001m k x y w a x x w n i i k i i j mj n i i k i j i ⋅⋅⋅⋅==⎥⎦⎤⎢⎣⎡∑∑∑===φφφ(5.18)是未知量为ma a a a ,,,,21⋅⋅⋅的线性方程组,称(5.18)式为正规方程组.实际中可适当选择函数系{}m k kx 0)(=φ,由正规方程组解出ma a a a ,,,,210⋅⋅⋅,于是可得最小二乘拟合函数∑==mk kk x a x 0)()(φφ.一般线性最小二乘拟合方法将上面一元函数的最小二乘拟合问题推广到多元函数,即为多维线性最小二乘拟合问题.假设已知多元函数),,,(21nx x x f y ⋅⋅⋅=的一组测量数据);,,,(21iniiiy x x x ⋅⋅⋅),,2,1(m i ⋅⋅⋅=和一组线性无关的函数系{}N k nk x x x 021),,,(=⋅⋅⋅φ,求函数∑=⋅⋅⋅=⋅⋅⋅Nk n k k n x x x a x xx 02121),,,(),,,(φφ对于一组正数mw w w ,,,21⋅⋅⋅,使得目标函数[]2121),,,(∑=⋅⋅⋅-=mi ni i i i i x x x y w S φ达到最小.其中待定系数N a a a a,,,,210⋅⋅⋅由正规方程组),,2,1,0(),(),(0N k y a Nj k j k j⋅⋅⋅==∑=φφφ确定,此处ini i i k mi i k ni i i k mi ni i i j i k j y x x x w y x x x x x x w ),,,(),(),,,(),,,(),(21121121⋅⋅⋅=⋅⋅⋅⋅⋅⋅=∑∑==φφφφφφ注:上面的函数φ关于ia 都是线性的,这就是线性最小二乘拟合问题,对于这类问题的正规组总是容易求解的.如果φ关于ia 是非线性的,则相应的问题称为非线性最小二乘拟合问题.非线性最小二乘拟合方法假设已知多元函数),,,(21nx x x f y ⋅⋅⋅=的一组测量数据);,,,(21iniiiy x x x ⋅⋅⋅),,2,1(m i ⋅⋅⋅=,要求一个关于参数),,2,1,0(N j a j⋅⋅⋅=是非线性的函数),,,;,,,(1021Nn a a a x x x ⋅⋅⋅⋅⋅⋅=φφ对一组正数mw w w ,,,21⋅⋅⋅使得目标函数[]21102110),,,;,,,(),,,(∑=⋅⋅⋅⋅⋅⋅-=⋅⋅⋅mi N ni i i i i N a a a x x x y w a a a S φ达到最小,则称之为非线性最小二乘问题.这类问题属于无约束的最优化问题,一般问题的求解是很复杂的,通常情况下,可以采用共轭梯度法、最速下降法、拟牛顿法和变尺度法等方法求解.实例:黄河小浪底调水调沙问题问题的提出2004年6月至7月黄河进行了第三次调水调沙试验,特别是首次由小浪底、三门峡和万家寨三大水库联合调度,采用接力式防洪预泄放水,形成人造洪峰进行调沙试验获得成功.整个试验期为20多天,小浪底从6月19日开始预泄放水,直到7月13日恢复正常供水结束.小浪底水利工程按设计拦沙量为75.5亿立方米,在这之前,小浪底共积泥沙达14.15亿吨.这次调水调试验一个重要目的就是由小浪底上游的三门峡和万家寨水库泄洪,在小浪底形成人造洪峰,冲刷小浪底库区沉积的泥沙.在小浪底水库开闸泄洪以后,从6月27日开始三门峡水库和万家寨水库陆续开闸放水,人造洪峰于29日先后到达小浪底,7月3日达到最大流量2700立方米/每秒,使小浪底水库的排沙量也不断地增加.下面是由小浪底观测站从6月29日到7月10日检测到的试验数据:表5-1: 试验观测数据单位:水流为立方米/每秒,含沙量为公斤/立方米·84··85·注:以上数据主要是根据媒体公开报道的结果整理而成的,不一定与真实数据完全相符.现在,根据试验数据建立数学模型研究下面的问题:(1) 给出估算任意时刻的排沙量及总排沙量的方法;(2) 确定排沙量与水流量的变化关系.模型的建立与求解对于问题(1),根据所给问题的试验数据,要计算任意时刻的排沙量,就要确定出排沙量随时间变化的规律,可以通过插值来实现.考虑到实际中排沙量应该是随时间连续变化的,为了提高精度,我们采用三次B样条函数进行插值.下面构造三次B样条函数)(x S y =.由试验数据,时间是每天的早8点和晚8点,间隔都是12个小时,共24个点)24,,2,1(⋅⋅⋅=i t i.为了计算方便,令)23,,,1,0(122128⋅⋅⋅=+⎥⎦⎤⎢⎣⎡⋅+-=i i t x i i(5.19)则it 对应于)23,,1,0(1⋅⋅⋅=+=i i x i.于是以)23,,1,0(⋅⋅⋅=i x i为插值节点(等距),步长1=h .其相应的排沙量为)23,,1,0(⋅⋅⋅=i y i 对应关系如下表:·86·表5-2: 插值数据对应关系单位:排沙量为公斤函数)(x S y =所满足的条件为 (1)23,,1,0,)(⋅⋅⋅==i y x S ii;(2) 3500)(,56400)(2223222323231212-=--≈'='=--≈'='x x y y x S y x xy yx S y .取)(x S 的三次B样条函数一般形式为∑-=⎪⎭⎫⎝⎛--=24103)(j j j h x x c x S Ω·87·其中)24,,1,0,1(⋅⋅⋅-=j cj为待定常数,1=h .在这里⎪⎪⎪⎩⎪⎪⎪⎨⎧≥<<+-+-≤+-=Ω2,021,342611,3221)(23233x x x x x x x x x且易知⎪⎪⎪⎩⎪⎪⎪⎨⎧≥±===Ω2,01,610,32)(3x x x x和⎪⎪⎩⎪⎪⎨⎧≥±===Ω'2,01,210,0)(3x x x x 根据B样条函数的性质,)(x S ''在[]23,x x 上连续,则有()∑-=--'='='2413)(j jj xx c x S y Ω由插值条件(1),(2)可得到下列方程组()()()⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧'=-'=''=-'='⋅⋅⋅==-=∑∑∑-=-=-=23241323024130241323)()(23,,1,0,)(y j c x S y j c x S i y j i c x S j j j j i j j i ΩΩΩ 即⎪⎩⎪⎨⎧'=+-'=+-⋅⋅⋅==++-+-23242311112223,,1,0,64y c c y c c i y c c c i i i i 将232324112,2y c c y c c '+='-=-代入前24个方程中的第一个和最后一个,便可得到方程组F AC =,其中·88·⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⋅⋅⋅⋅⋅⋅=⨯232102424,421410141014124c c c c C A ,⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡'-'+=3400048000684000458400266626232322100 y y y y y y F显然A 为满秩阵,方程组F AC =一定有解,用消元法求解可得问题的解为56044.39830=c , 4117111.2031=c , 2159510.7882=c , 9189845.6433=c ,1203106.6364=c , 8239727.8115=c ,8249182.1166=c , 1263543.7217=c ,9287842.9988=c , 2302284.2839=c ,4317419.86810=c , 1304836.24311=c ,3307635.15912=c ,6305423.11913=c ,2270672.36214=c ,4240287.43115=c ,0154177.91216=c ,4103000.92017=c ,99818.406218=c , 43725.454719=c ,49279.775020=c ,32155.445221=c , 2098.444222=c ,7450.777923=c ,-450.777924311.2034,2232324011='+=='-=-y c c y c c . 将)24,,1,0,1(⋅⋅⋅-=j c j代入()∑-=--==24131)(j jj x c x S y Ω(5.20)即得排沙量的变化规律.由(5.19)和(5.20)式可得到第i 时间段(12小时为一段)内,任意时刻]12,0[∈t 的排沙量.则总的排沙量为()dt j t c dx x S Y j j⎰∑⎰-=--Ω==284824132411)(经计算可得1110844.1⨯=Y 吨,即从6月29日至7月10日小浪底水库排沙总量大约为1.844亿吨,此与媒体报道的排沙量基本相符.对于问题(2),研究排沙量与水量的关系,从试验数据可以看出,开始排沙量是随着水流量的增加而增长,而后是随着水流量的减少而减少.显然,变化规律并非是线性的关系,为此,我们问题分为两部分,从开始水流量增加到最大值2720立方米/每秒(即增长的过程)为一段,从水流量的最大值到结束为第二段,分别来研究水流量与排沙量的关系.具体数据如表5-3和5-4.表5-3: 第一阶段试验观测数据 单位:水流为立方米/每秒,含沙量为公斤/立方米表5-4: 第二阶段试验观测数据单位:水流为立方米/每秒,含沙量为公斤/立方米对于第一阶段,由表5-3用Matlab作图(如图5-3)可以看出其变化趋势,我们用多项式作最小二乘拟合.·90··91·图5-3设拟合函数为∑==mk kk x a x 1)(φ确定待定常数),,1,0(m k ak=使得211111102])([∑∑∑===⎥⎦⎤⎢⎣⎡-=-=i i i m k k i k i i y x a y x S φ有最小值.于是可以得到正规方程组为m k x y a x mj i k i i j i j k i ,,1,0,0111111 ==⎪⎭⎫⎝⎛∑∑∑===+ 当3=m 时,即取三次多项式拟合,则3,2,1,0,1113111321112111110111==⎪⎭⎫⎝⎛+⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛∑∑∑∑∑==+=+=+=k x y a x a x a x a x i k i i i k i i k i i k i i k i求解可得73321108423.1,103172.1,3.1784,-2492.9318--⨯=⨯-===a a a a .于是可得拟合多项式为332213)(x a x a x a a x +++=φ,最小误差为847.72=S ,拟合效果如图所示.·92·图:三次拟合效果,带*号的为拟合曲线.类似地,当4=m 时,即取四次多项式拟合,则正规方程组为4,3,2,1,0111411143111321112111110111==⎪⎭⎫⎝⎛+⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛∑∑∑∑∑∑==+=+=+=+=k x y a x a x a x a x a x i ki i i k i i k i i k i i k i i k i求解可得104633210109312.1,1094.1,102626.7,12.0624,-7434.6557---⨯-=⨯=⨯-===a a a a a 于是可得拟合多项式为443322104)(x a x a x a x a a x ++++=φ,最小误差为102.66=S ,拟合效果如图5-5所示.图5-5:四次拟合效果,带*号的为拟合曲线.从上面的三次多项式拟合和四次多项拟合效果来看,差别不大.基本可以看出排沙量与水流量的关系.图5-6:第二段三·93··94· 次多项式拟合效果对于第二阶段,由表5-4可以类似地处理.我们用线性最小二乘法作三次和四多项式拟合.拟合效果如图5-6和5-7所示,最小误差分别为5.459=S 和1.236=S . 从拟合效果来看,显然四次多项式拟合要比三次多项式拟合好的多.图5-7:第二段四次多项式拟合效果。
插值和数据拟合一、 插值方法问题:已知n+1个节点(x j ,y j )(j=0,1,…,n),a=x 0<x 1<…< x n =b ,求任一插值点x*处的插值y*方法:构造一个相对简单的函数y=f(x),使得f 通过所有节点,即f(x j )= y j ,再用y=f(x)计算x*的值。
1. 拉格朗日多项式插值设f(x)是n 次多项式,记作1110()n n n n n L x a x a x a x a --=++++要求对于节点(,)j j x y 有(),0,1,,n j j L x y j n ==将n+1个条件带入多项式,就可以解出多项式的n+1个系数。
实际上,我们有n 次多项式011011()()()()()()()()()i i n i i i i i i i n x x x x x x x x l x x x x x x x x x -+-+----=----满足1,()0,,,0,1,,i j i jl x i j i j n =⎧=⎨≠=⎩则0()()nn i i i L x y l x ==∑就是所要的n 次多项式,称为拉格朗日多项式。
由拉格朗日多项式计算的插值称为拉格朗日插值。
一般来讲,并不是多项式的阶数越高就越精确,一般采用三阶、二阶或一阶(线性)多项式,对相邻点进行分段插值。
2. 样条插值在分段插值时,会造成分段点处不光滑,如果要求在分段点处光滑,即不仅函数值相同,还要一阶导数和二阶导数相同,则构成三阶样条插值。
一般用于曲线绘制,数据估计等。
例 对21,[5,5](1)y x x =∈-+,用n=11个等分节点做插值运算,用m=21个等分插值点作图比较结果。
见inter.m 程序二、 曲线拟合 三、 给药方案 1. 问题一种新药用于临床必须设计给药方案,在快速静脉注射的给药方式下,就是要确定每次注射剂量多大,间隔时间多长.我们考虑最简单的一室模型,即整个机体看作一个房室,称为中心室,室内血液浓度是均匀的.注射后浓度上升,然后逐渐下降,要求有一个最小浓度1c 和一个最大浓度2c .设计给药浓度时,要使血药浓度保持在1c ~2c 之间.2. 假设(1)药物排向体外的速度与中心室的血药浓度成正比,比例系数是k(>0),称为排出速度.(2)中心室血液容积为常数V ,t=0的瞬间注入药物的剂量为d ,血药浓度立即为dV. 3. 建模设中心室血药浓度为c(t),满足微分方程(0)dckc dtd c V=-=用分离变量法解微分方程,有()ktd c te V-=(*) 4. 方案设计每隔一段时间τ,重复注入固定剂量D ,使血药浓度c(t)呈周期变化,并保持在1c ~2c 之间.如图:设初次剂量加大到D 0,易知0221,D Vc D Vc Vc ==-,2121()11ln[],()()ln c Vc t t t c t c k d k c τ=-=-= 那么,当12,c c 确定后,要确定给药方案0{,,}D D τ,就要知道参数V 和k .5. 由实验数据做曲线拟合确定参数值已知1210,25(/)c c g ml μ==,一次注入300mg 药物后,间隔一定ln lndc kt V=- 记12ln ,,lndy c a k a V==-=,则有 12y a t a =+求解过程见medicine_1.m得120.2347, 2.9943a a =-=,由d=300(mg)代入算出k=0.2347,V=15.02(L) 从而有0375.5(),225.3(), 3.9()D mg D mg τ===小时四、 口服给药方案 1. 问题口服给药相当于先有一个将药物从肠胃吸收入血液的过程,可简化为一个吸收室,一个中心室,记t 时刻,中心室和吸收室的血液浓度分别是1()()c t c t 和,容积分别是V ,V1,中心室的排除速度为k ,吸收速度为k1,且k,k1分别是中心室和吸收室血液浓度变化率与浓度的比例系数,t=0口服药物的剂量为d ,则有11111,(0)dc dk c c dt V =-= (1) 111,(0)0V dckc k c c dt V=-+= (2) 解方程(1)有111()k td c te V -=代入方程(2)有111()()k t kt k d c t e e V k k--=--其中三个参数1,,dk k b V=,可由下列数据拟合得到:(非线性拟合)。
数学建模精品教材-第九章插值与拟合第九章插值与拟合插值:求过已知有限个数据点的近似函数。
拟合:已知有限个数据点,求近似函数,不要求过已知数据点,只要求在某种意义下它在这些点上的总偏差最小。
插值和拟合都是要根据一组数据构造一个函数作为近似,由于近似的要求不同,二者的数学方法上是完全不同的。
而面对一个实际问题,究竟应该用插值还是拟合,有时容易确定,有时则并不明显。
§1 插值方法下面介绍几种基本的、常用的插值:拉格朗日多项式插值、牛顿插值、分段线性插值、Hermite 插值和三次样条插值。
1.1 拉格朗日多项式插值1.1.1 插值多项式用多项式作为研究插值的工具,称为代数插值。
其基本问题是:已知函数 f x 在区间[a,b]上n +1个不同点x ,x , L,x 处的函数值 y f x i 0,1, L,n,求一个0 1 n i i至多n次多项式nx a +a x + L +a x (1)n 0 1 n使其在给定点处与 f x同值,即满足插值条件 x f x y i 0,1, L,n(2) n i i ix称为插值多项式,x i 0,1, L,n称为插值节点,简称节点,[a,b]称为插值区n i间。
从几何上看,n次多项式插值就是过n +1个点 x , f x i 0,1, L,n,作一条i i多项式曲线 y x近似曲线 y f x。
nn次多项式(1)有n +1个待定系数,由插值条件(2)恰好给出n +1个方程2 na +a x +a x + L +a x y0 1 0 2 0 n 0 02 na +a x +a x + L +a x y0 1 1 2 1 n 1 1(3)L L L L L L L L L L L L2 na +a x +a x + L +a x y0 1 n 2 n n n n 记此方程组的系数矩阵为A,则2 n1 x x L x0 0 02 n1 x x L x1 1 1 detAL L L L L L L2 n1 x x L xn n n是范德蒙特Vandermonde行列式。
Matlab数学建模学习笔记——插值与拟合⽬录插值与拟合插值和拟合的区别图⽚取⾃知乎⽤户yang元祐的回答插值:函数⼀定经过原始数据点。
假设f(x)在某区间[a,b]上⼀系列点上的值y_i=f(x_i),i=0,1,\dots,n。
插值就是⽤较简单、满⾜⼀定条件的函数\varphi(x)去代替f(x)。
插值函数满⾜条件\varphi(x_i)=y_i,i=0,1,\dots,n拟合:⽤⼀个函数去近似原函数,不要求过已知数据点,只要求在某种意义下它在这些点上的总偏差最⼩。
插值⽅法分段线段插值分线段插值就是将每两个相邻的节点⽤直线连起来,如此形成的⼀条折线就是就是分段线性插值函数,记作I_n(x),它满⾜I_n(x_i)=y_i,且I_n(x)在每个⼩区间[x_i,x_{i+1}]上是线性函数(i=0,1\dots,n-1)。
I_n(x)可以表⽰为I_n(x)=\sum_{i=0}^n y_il_i(x),其中l_i(x)= \begin{cases} \frac{x-x_{i-1}}{x_i-x_{i-1}},&x\in [x_{i-1},x_i],i \neq 0,\\ \frac{x-x_{i+1}}{x_i-x_{i+1}},&x\in [x_i,x_{i+1}],i \neq n,\\ 0,&其他 \end{cases}I_n(x)有良好的收敛性,即对x\in [a,b],有\lim _{n \rightarrow \infin}I_n(x)=f(x)⽤I_n(x)计算x点的插值的时候,只⽤到x左右的两个点,计算量与节点个数n⽆关。
但是n越⼤,分段越多,插值误差越⼩。
拉格朗⽇插值多项式朗格朗⽇(Lagrange)插值的基函数为\begin{aligned} l_i(x)&=\frac{(x-x_0)\cdots(x-x_{i-1})(x-x_{i+1})\cdots(x-x_n)}{(x_i-x_0)\cdots(x_i-x_{i-1})(x_i-x_{i+1})\cdots(x_i-x_n)}\\ &= \prod_{j=0\\j\neq i}^{n} \frac{x-x_j}{x_i -x_j},i=0,1,\cdots,n。
在数值计算中,插值和拟合是两种常用的方法,用于通过已知数据点推测未知数据点的数值。
插值是一种通过已知数据点构建一个函数,以便在这些数据点之间进行预测。
而拟合是一种将一个函数与已知数据点进行匹配,以便预测未知数据点的数值。
插值的目标是通过经过已知数据点的连续函数来准确地估计未知数据点的数值。
最简单的插值方法是线性插值,它假设两个相邻数据点之间的函数值是线性变化的。
线性插值可以用于计算两个已知数据点之间的任何位置的函数值。
如果我们有更多的数据点,可以使用更高阶的插值方法,如二次插值或三次插值。
这些方法使用多项式来表示数据点之间的函数,以便更准确地预测未知数据点。
然而,插值方法并不总是最理想的选择。
在某些情况下,通过已知数据点精确地构建一个连续函数是不可能的。
这可能是因为数据点之间的差异太大,或者数据点的数量太少。
在这种情况下,拟合方法可以提供更好的预测结果。
拟合的目标是找到一个函数,使其与已知数据点的误差最小。
最常用的拟合方法是最小二乘拟合,它通过最小化数据点的残差的平方和来找到最佳拟合函数。
最小二乘拟合可以用于各种不同的函数类型,如线性拟合、多项式拟合、指数拟合等。
根据数据点的分布和特性,我们可以选择适当的拟合函数来获得最准确的预测结果。
在实际应用中,插值和拟合方法经常同时使用。
例如,在地理信息系统中,我们可能需要通过已知地点的气温数据来估计未知地点的气温。
我们可以使用插值方法来构建一个连续函数,以便在已知地点之间预测未知地点的气温。
然后,我们可以使用拟合方法来匹配这个连续函数与其他已知数据点,以提高预测的准确性。
插值和拟合方法在科学、工程、金融等各个领域都有广泛的应用。
在科学研究中,它们可以用于数据分析和预测,以帮助我们理解和解释实验结果。
在工程中,它们可以用于控制系统设计、信号处理和机器学习等领域。
在金融领域,它们可以用于市场预测和风险管理等重要任务。
总而言之,插值和拟合是数值计算中常用的方法,用于通过已知数据点推测未知数据点的数值。
数值分析中的插值与拟合插值和拟合是数值分析中常用的技术,用于估计或预测数据集中缺失或未知部分的数值。
在本文中,我们将讨论插值和拟合的概念、方法和应用。
一、插值插值是通过已知数据点之间的连续函数来估计中间数据点的数值。
插值方法可以根据不同的数据和需求选择合适的插值函数,常用的插值方法包括拉格朗日插值、牛顿插值和埃尔米特插值。
1.1 拉格朗日插值拉格朗日插值是一种基于多项式的插值方法。
通过已知的n个数据点,可以构建一个n-1次的插值多项式。
这个多项式通过已知数据点上的函数值来准确地经过每一个点。
1.2 牛顿插值牛顿插值方法也是一种多项式插值方法,通过差商的概念来构建插值多项式。
差商是一个递归定义的系数,通过已知数据点的函数值计算得出。
牛顿插值可以通过递推的方式计算出插值多项式。
1.3 埃尔米特插值埃尔米特插值是一种插值方法,适用于已知数据点和导数值的情况。
它基于拉格朗日插值的思想,通过引入导数信息来逼近数据的真实分布。
埃尔米特插值可以更准确地估计数据点之间的值,并且可以保持导数的连续性。
二、拟合拟合是通过一个模型函数来逼近已知数据点的数值。
拟合方法旨在找到最适合数据集的函数形式,并通过最小化误差来确定函数的参数。
常见的拟合方法包括最小二乘法、多项式拟合和曲线拟合。
2.1 最小二乘法最小二乘法是一种常用的拟合方法,通过最小化数据点到拟合函数的误差平方和来确定最佳拟合曲线或曲面。
最小二乘法适用于线性和非线性拟合问题,可以用于拟合各种类型的非线性函数。
2.2 多项式拟合多项式拟合是一种基于多项式函数的拟合方法。
通过多项式的线性组合来近似已知数据集的数值。
多项式拟合可以通过最小二乘法或其他优化算法来确定拟合函数的系数。
2.3 曲线拟合曲线拟合是一种用曲线函数来逼近已知数据点的拟合方法。
曲线函数可以是非线性的,并且可以根据数据的特点进行选择。
曲线拟合可以通过优化算法来确定拟合函数的参数。
三、应用插值和拟合在数值分析中有广泛的应用。
插值与拟合方法插值和拟合是数学中常用的方法,用于根据已知数据点的信息,推断出未知数据点的数值或函数的形式。
插值和拟合方法是经典的数学问题,应用广泛,特别是在数据分析、函数逼近和图像处理等领域。
1.插值方法:插值方法是通过已知数据点的信息,推断出两个已知数据点之间的未知数据点的数值。
插值方法的目的是保证插值函数在已知数据点处与实际数据值一致,并且两个已知数据点之间的连续性良好。
最常用的插值方法是拉格朗日插值法和牛顿插值法。
拉格朗日插值法根据已知数据点的横纵坐标,构造一个多项式函数,满足通过这些数据点。
拉格朗日插值法可以用于任意次数的插值。
牛顿插值法是使用差商的概念进行插值。
差商是指一个多项式在两个数据点之间的斜率。
牛顿插值法通过迭代计算得到与已知数据点一致的多项式。
插值方法的优点是可以精确地经过已知数据点,但是在两个已知数据点之间的插值部分可能会出现震荡现象,从而导致插值结果不准确。
2.拟合方法:拟合方法是通过已知数据点的信息,找出一个函数或曲线,使其能够最好地拟合已知数据点。
拟合方法的目标是寻找一个函数或曲线,尽可能地逼近已知数据点,并且能够在未知数据点处进行预测。
最常用的拟合方法是最小二乘法。
最小二乘法是通过求解最小化残差平方和的问题来进行拟合。
残差是指已知数据点与拟合函数的差异。
最小二乘法的目标是找到一个函数,使得所有数据点的残差平方和最小。
拟合方法的优点是可以得到一个光滑的函数或曲线,从而可以预测未知数据点的数值。
但是拟合方法可能会导致过拟合问题,即过度拟合数据点,导致在未知数据点处的预测结果不准确。
除了最小二乘法,还有其他的拟合方法,如局部加权回归和样条插值等。
局部加权回归是一种基于最小二乘法的拟合方法,它通过赋予不同的数据点不同的权重,来实现对未知数据点的预测。
样条插值是一种基于多项式插值的拟合方法,它将整个数据集分段拟合,并且在分段部分保持连续性和光滑性。
总结:插值和拟合方法是数学中的经典方法,用于根据已知数据点的信息,推断出未知数据点的数值或函数的形式。