数学建模与数学实验：ch5 数值分析法建模

1、1,拟合与插值,第5章 数值分析法建模,作业： 习题5：1，2 上机实验3：习题5：3,2,拟合函数,已知函数f在离散样本点集M上的函数值y1,y2,yn，拟合是指选取含参数的函数形式f(x;a1,a2 ,am) ，通过调整待定参数a1,a2 ,am, 使得在样本点集上的函数值与实际值误差达到最小。 如果选取线性函数f，则称线性拟合或者线性回归，否则称为非线性拟合或者非线性回归。 拟合特点：函数不要求经过样本数据点。因此当数据与真实值有误差时，一般用数据拟合,3,插值函数,插值是指已知某函数的在离散点集M上的函数值或者导数信息，通过求解函数形式，使得函数在给定离散点上满足一定约束条件。插值函数

2、特点：必须经过所有样本数据点。当数据是准确无误差时，一般用插值,4,曲线拟合步骤,1.确定经验公式形式,2.计算经验公式中的参数,3.检验经验公式有效性,5,利用已知的结论,曲线改直技术,作散点图法,多项式近似,如何确定经验公式形式,6,最小二乘法原理就是找一组参数 使最小二乘误差 最小,数学模型,曲线拟合问题最常用的解法线性最小二乘法的基本思路,选定一组基函数 r1(x), r2(x), rn(x), nm, 选取拟合函数形式为 f(x)=a1r1(x)+a2r2(x)+ +anrn(x) ，其中 a1,a2, an 为待定系数,正规方程,特例：多项式拟合,模型求解,Matlab求解,7,差

3、分与差商概念 一阶前向差分 二阶前向差分 m阶前向差分 一阶差商 二阶差商 m阶差商,8,差商与导数的联系,若y=f(x)在a,b上m阶可微，且 则,若节点为等距分布时，有 ，h为节点距，则,定拟合多项式的最佳阶数： 差分表/插商表波动最小原则 等距节点用差分表，不等距节点用差商表,9,差分表,10,xk f(xk) 一阶差商 二阶差商 三阶差商 n 阶差商,差商表,11,与拟合有关的MATLAB 函数,polyfit:多项式拟合 poly2sym:由多项式系数向量得多项式符号表达式 polyval:计算多项式函数在指定处的函数值 poly:计算过固定点的多项式 lsqcurvefit, ls

4、qnonlin非线性最小二乘拟合 fit,fittype一般非线性拟合 cftool 拟合工具箱,作多项式f(x)=anxn+ +a1x+a0拟合,可利用Matlab命令,a=polyfit(x,y,n,yy=polyval(a,xx).计算多项式a在xx处的值,12,当精确函数 y = f(x) 非常复杂或未知时，在区间a,b上一系列节点 x0 xm 处测得函数值 y0 = f(x0), , ym = f(xm)，由此构造一个简单易算的 近似函数 g(x) ，满足条件 g(xj) = f(xj) (j = 0, m) (*) 这类问题称为“插值问题,插值问题,g(x) 称为f(x) 的插值函

5、数,一般取多项式函数,x0 xm称为插值节点, 插值节点互不相同,条件(*)称为插值条件，区间a,b称为插值区间,13,x0,x0,x2,xm-1,xm,f(x,g(x,精确函数 y = f(x) 非常复杂或未知,g(x) 称为f(x) 的插值函数,节点 x0 xm称为插值节点,14,基本思想：在n次多项式空间Pn中找一组合适的基函数 0(x),1(x), n(x),使,pn(x)=a0 0(x) +a1 1(x) +an n(x,不同的基函数的选取导致不同的插值方法,无论是从理论和计算的角度,还是从应用的 角度看,多项式都是最简单的函数,因此,多项式 插值是最基本的插值方法,15,根据插值条

6、件，应有,拉格朗日(Lagrange)插值,16,称为拉格朗日插值基函数,拉格朗日插值多项式公式也可以通过构造插值基函数方法直接得到,其中li(x) 为n次多项式,拉格朗日(Lagrange)插值,17,拉格朗日(Lagrange)插值,特别地,两点一次(线性)插值多项式,三点二次(抛物线)插值多项式,18,拉格朗日多项式 插值的这种振荡现象 叫 Runge现象,采用拉格朗日多项式插值：选取不同插值节点个数n+1，其中n为插值多项式的次数，当n分别取2,4,6,8,10时，绘出插值结果图形,例,1ch.m,19,Newtons Interpolation基本原理,Lagrange 插值虽然易算

7、，但若要增加一个节点时，全部基函数 li(x) 都需要重新计算,能否重新在Pn中寻找新的基函数？希望每加一个节点时，只附加一项上去即可,20,选取1,x - x0 , (x - x0)(x - x1),(x-x0)(x-x1) (x-xn-1)构成Pn的一组基函数,Newtons Interpolation基本原理,利用插值条件Nn(xj)=f(xj), j=0,1，n代入上式，得关于Ak (k=0,1，n)的线性代数方程组,当xj 互异时，系数矩阵 非奇异，方程有唯一解,21,Lagrange插值与Newton插值的异同点,两者都是通过给定n+1个互异的插值节点，求一条n次多项式曲线近似地表

8、示待插值的函数曲线,Lagrange插值和Newton法插值的结果和余项都是一致的，因为都是利用n次多项式插值,区别： Lagrange插值法在求每个基函数的时候要用到所有结点，因此如果需要再多加一个结点的话，需要重新求出函数才可，而这需要大工作量，于是数学家们就发明了Newton法。 Lagrange插值法是通过构造n+1个n次基函数，作线性组合而得到 Newton法插值是通过求各阶差商，递推得到公式 f(x)=f(x0)+(x-x0)fx0,x1+(x-x0)(x-x1)fx0,x1,x2.+(x-x0).(x-xn-1)fx0,x1.xn,22,在数学上，光滑程度的定量描述是：函数(曲线

9、)的k阶导数存在且连续，则称该曲线具有k阶光滑性。 光滑性的阶次越高，则越光滑。是否存在较低次的分段多项式达到较高阶光滑性的方法？三次样条插值就是一个很好的例子,三次样条插值,23,Cubic Spline Lagrange Interpolation Interpolation,24,三次样条插值,f(x)为被插值函数,插值条件,连接条件,25,yi=interp1(x，y，xi，method,nearest ：最邻近插值linear ： 线性插值； spline ： 三次样条插值； cubic ： 立方插值。 缺省时： 线性插值,插值方法要求x是单调的，且xi不能够超过x的范围,用MATL

10、AB作一维插值计算,扩展材料,26,1. 拉格朗日插值:自编程序,2. 分段线性插值:已有程序 y=interp1(x0,y0,x,3. 三次样条插值:已有程序 y=interp1(x0,y0,x,spline) 或 y=spline(x0,y0,x,注意：几乎所有的插值方法都要求x0单调(例外：cscvn)， 并且x不能够超过x0的范围,interp2,interp3,interpn多元函数插值,27,pp=spline (x,y) ：样条函数的表示（PP结构） pp=cscvn(x,y): 自然样条函数的表示（不要求x单调！） yi=ppval(pp,xx)或 yi=spline(x,y,

11、xx)：样条函数求值 fprime = fnder(pp), fprime = fnder(pp,dorder) 样条函数求导 inters = fnint(pp) , intgrf = fnint(pp,value)样条函数积分 fnplt(pp)样条函数绘图,样条函数的相关MATLAB命令,28,机床加工问题,用程控铣床加工机翼断面的下轮廓线时 每一刀只能沿x方向和y方向走非常小的一步。 表3-1给出了下轮廓线上的部分数据 但工艺要求铣床沿x方向每次只能移动0.1单位. 这时需求出当x坐标每改变0.1单位时的y坐标。 试完成加工所需的数据，画出曲线,29,求解机床加工问题,x0=0 3 5

12、 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 on,30,例1：每隔1小时测量一次温度，测得的12个温度依次为：5，8，9，15，25，29，31，30，22，25，27，24。试估计每隔6分钟的温度值,Temp.m,hours=1:12; temps=5 8 9 15 25 29 31 30 22 25 27 24; h=1:0.1:12; t=interp1(hours,temp

13、s,h,spline); %(直接输出数据将是很多的) plot(hours,temps,+,h,t,hours,temps,r:) %作图 xlabel(Hour),ylabel(Degrees Celsius,31,要求x0,y0单调；x，y可取为矩阵，或x取行向量，y取为列向量，x,y的值分别不能超出x0,y0的范围,z=interp2(x0,y0,z0,x,y,method,nearest 最邻近插值 linear 双线性插值 cubic 双三次插值 缺省时, 双线性插值,用MATLAB作网格节点数据的插值(二维,32,cz =griddata（x，y，z，cx，cy，method,要

14、求cx取行向量，cy取为列向量,nearest 最邻近插值 linear 双线性插值 cubic 双三次插值 v4 Matlab提供的插值方法 缺省时, 双线性插值,用MATLAB作散点数据的插值计算,33,航行区域的警示线,某海域上频繁地有各种吨位的船只经过。 为保证船只的航行安全，有关机构在低潮时对水深进行了测量,下表是测量数据： 表3 水道水深的测量数据 x129.0 140.0 103.5 88.0 185.5 195.0 105.5 y 7.5 141.5 23.0 147.0 22.5 137.5 85.5 z 4 8 6 8 6 8 8 x157.5 107.5 77.0 81.

15、0 162.0 162.0 117.5 y-6.5 -81.0 3.0 56.5 -66.5 84.0 -33.5 z 9 9 8 8 9 4 9,34,航行区域的警示线,其中（x, y）为测量点，z为（x, y）处的水深(英尺)，水深z是区域坐标（x, y）的函数z= z (x, y）， 船的吨位可以用其吃水深度来反映，分为 4英尺、4.5英尺、5英尺和 5.5英尺 4 档。 航运部门要在矩形海域（75，200）（50，150）上为不同吨位的航船设置警示标记。 请根据测量的数据描述该海域的地貌，并绘制不同吨位的警示线，供航运部门使用,35,Matlab求解,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