1、单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,本资料仅供参考,不能作为科学依据。谢谢。本资料仅供参考,不能作为科学依据。本资料仅供参考,不能作为科学依据。谢谢。本资料仅供参考,不能作为科学依据。,第三章 插值法与最小二乘法,在实际问题中碰到函数 有些有解析表示式,但很复杂,有些只给出一些离散数据,这给我们求解函数值、导数值、零点、极值和积分值等带来了很多不便。对于这些情况,自然想法是,设法找到某个简单近似函数满足 。本章介绍两种方法,即插值法和最小二乘法。,1/57,3.1 插值法,3.1.1 多项式插值概念,在众多插值函数中,多项式是最简单最易计算
2、设,函数在区间,a,b,上连续,且在,n+1,个不一样点,上分别取值为 。,在多项式插值中,最基本、最简单问题是求一个次数不超出,代数多项式,(3-1),满足,其中 均为实数,(3-2),2/57,称为被插值函数;,称插值多项式;,条件(3-2)为插值条件;,为插值点。其,几何意义为,图3-1 插值几何意义,满足(3-2)式插值多项式是存在且唯一。原因是满足,3/57,(3-2)式多项式未知系数行列式为著名范德蒙德(Vardermonde)行列式,故解是存在且唯一。,3.1.2,拉格朗日插值多项式,在每个插值点 结构插值基函数 为一个,次多项式 ,且满足条件,(3-3),4/57,(3-4)
3、由(3-4)式和多项式 定义以及(3-2)式插值条件,我们有:,再由(3-4)式知,是,零点,故按因式定理必含有因式:,既是 次多项式,可知,5/57,再由条件 ,得,于是有,所以,拉格朗日插值多项式可写为,(3-5),6/57,例3-1,已知特殊角 ,用,近似值。,解:令,用 为节点,有,一次和二次拉格朗日插值公式求,7/57,用 为节点,有,8/57,拉格朗日插值多项式程序设计,计算公式,(3-5),程序为二重循环。由内循环(,j,循环)经过累乘求得基函数 ,然后经过外循环(,i,循环)得到插值结果,y,.,下列图描述了拉格朗日算法流程,对于该框图中相关参数说明以下:,N:,插值多项式次
4、数(插值点个数,-1,),Xi:,插值点,Yi:,插值点上函数值,LN:,拉格朗日插值结果,T:,存放累乘积,I:,外循环变量,J:,内循环变量,9/57,START,INPUT N,X,Xi,Yi(I=0,N),LN=0,I=0,N,T=1,J=0,N,J=I?,NO,T=(X-Xj)/(Xi-Xj)*T,END J,LN=LN+Yi*T,END I,OUTPUT,END,YES,10/57,拉格朗日算法程序,DIMENSION,X(I),Y(I),REAL,LN,READ(*,*)N,DO K=1,N,READ(*,*)X(I),Y(I),END DO,LN=0.0,DO I=0,N,T=
5、1,DO J=0,N,T=T*(X-X(J)/(X(I)-X(J),END DO,LN=LN+T*Y(I),END DO,WRITE(*,*)LN=,LN,11/57,1、,插值余项,对插值后,用 代替 必定会产生误差,其误差可表示为,称上式为插值多项式,插值余项,。,设 在区间 上有直到 阶导数,为,上 个互异节点,为满足条件,次插值多项式。,那么对于任何,有,(3-6),其中 且依赖于 。(普通取,最大值,3.1.3,插值多项式余项,12/57,3.1.4,牛顿插值多项式,在插值基函数 中,按这种方式来求 :由线性代数性质可知,一个不高于n次多项式可表示成,线性组合。在满足插值条件 情况下
6、可表示为(记为 ):,(3-6),上式称为,牛顿(,Newton,)插值多项式,。,13/57,3.1.4,牛顿插值多项式,1.1、,向前差分与牛顿向前插值多项式,取节点为等距离,则,其中,h,为步长。在两相临节点处函数值之差为,称为函数,f(x),在节点处以,h,为步长,一阶向前差分,(简称一阶差分)。又称一阶差分差分为二阶差分,记为,14/57,3.1.4,牛顿插值多项式,为了便于计算与一目了然,用表格描述为,表3.1,15/57,3.1.4,牛顿插值多项式,在,等距离节点情况下,能够用差分来确定,(3.6),式待定系数,并将牛顿插值,多项式加以简化。,在节点上,函数值是已知,所以,16
7、/57,3.1.4,牛顿插值多项式,将上式带入,(3.6),式,并令 ,则式,(3.6),可简化为,该式称为,牛顿向前插值公式,。其插值余项是,(3-7),(3-8),17/57,3.1.4,牛顿插值多项式,2、,向后差分与牛顿向后插值公式,a,、向后差分,b,、中心差分,18/57,3.1.4,牛顿插值多项式,3、,差商与牛顿插值多项式,用差商来确定,(3.6)式,中系数,一阶差商,二阶差商,依次类推,19/57,3.1.4,牛顿插值多项式,表3.差商表,20/57,3.1.4,牛顿插值多项式,x,4.0002,4.0104,4.0233,4.0294,f(x),0.6020817,0.60
8、31877,0.6045824,0.6052404,例3-2,给出列表函数f(x)=lg,x,值,以下表所表示,试用牛顿插值法求lg,4.01,值。,表3.3,解:,依据给定函数作差商表,如表3.3所表示。,x=4.01,x,1,=4.0002,x,2,=4.0233,x,3,=4.0294,21/57,3.1.4,牛顿插值多项式,lg,x,=f(x,0,)+f(x,0,x,1,)(x-x,0,)+f(x,0,x,1,x,2,)(x-x,0,)(x-x,1,),+f(x,0,x,1,x,2,x,3,)(x-x,0,)(x-x,1,)(x-x,2,),表3.4,x,4.0002,4.0104,4
9、0233,4.0294,0.6020817,0.6031877,0.6045824,0.6052404,0.108433,0.108116,0.107869,-0.013636,-0.013000,0.021781,f(x),一阶差商,二阶差商,三阶差商,将以上各插值点值、函数值和差商值代入插值公式得lg,4.01,=,0.6031443,22/57,3.1.5,Hermite插值多项式,前述插值多项式均只要求插值点函数值相等,为深入光滑函数曲线,有时还要求其导数值也相等,并等于已知值 次多项式 。该多项式 称为 Hermite插值多项式,它也是一个近似式。,求解Hermite插值多项式一个
10、最简单方法就是直接应用牛顿插值法。即已知节点 处函数值与导数值,时,把 视为 重节点,,并注意 个重节点 阶差商,23/57,3.1.5,Hermite插值多项式,例 3-2,已知函数,函数值、导数值以下,求其插值多项式及误差表示式。,-1,0,1,0,-4,-2,0,5,6,解:,24/57,3.1.5,Hermite插值多项式,作差商表以下。f(0,0)=y(0)/1!=0,f(0,0,0)=y”(0)/2!=3,f(1,1)=y(1)/1!=5.按差商表3.2有:,-1,0,0,0,1,1,0,-4,-4,-4,-2,-2,-4,0,0,2,5,4,3,2,3,-1,-1,1,0,2,1
11、25/57,3.1.6,三次样条插值多项式,1、,三次样条插值函数定义,对于一组已知数据 ,且 ,若函数 满足:,(1),、在每一个子区间 上 都能够用最高为三次多项式来表示。,、在 上二次连续可微。即 ,,连续。,、插值条件,则称为函数f(x)关于节点,三次样条插值函数,。,26/57,3.1.6,三次样条插值多项式,2、,边界条件问题,依据定义,在任意子区间上,三次样条函数可表示为,该式中有,4,个待定系数,故共有,4n,个待定系数。由定义中连续条件和插值条件可,列出以下方程:,27/57,3.1.6,三次样条插值多项式,上式共有,4n-2,个方程,所以还需要,2,个方程才能确定,4n,
12、个系数,这就要用到边界,条件。边界条件类型很多,常见有三种:,1、在边界上给定一阶导数值 ;,2、给定二阶导数值;(假如,,称其为,自然边界条件,3、,f(x)是以b-a为周期函数,则要求S(x)及其导数也是以b-a为周期函数,对应边界条件为,28/57,3.1.6,三次样条插值多项式,3、三弯矩方程,设S(x)在某一节点二阶导数为,(在力学上解释为:M,i,是细梁在截面处弯矩,因为该弯矩与相临两个弯矩相关,故称三弯矩方程),因为S(x)最高次数不超出三次,所以它二阶导数是线性函数。,令 ,于是,29/57,3.1.6,三次样条插值多项式,连续积分两次得,其中 为积分常数。由以下插值条件确定,
13、所以,(3-9),30/57,3.1.6,三次样条插值多项式,代入(3.9)式得,从上式可知,只要确定M值(M,0,M,1,M,n,共n+1个),就能够完全确定S(x)。在确定(3.9)式中待定系数时用了插值点条件,下面利用节点一阶导数连续条件,(3-10),(3-11),31/57,3.1.6,三次样条插值多项式,由(3.10)式得,于是,(3-12),(3-13),32/57,3.1.6,三次样条插值多项式,由(3.11)式、(3.12)式和(3.13)式得,若记,则方程可简写为,(3-14),33/57,3.1.6,三次样条插值多项式,即,上方程组有n-1个方程,但有n+1个未知数。所以
14、要确定M,i,还需用边界条件。,1、,第一个边界问题,边界条件,则,(3-15),34/57,3.1.6,三次样条插值多项式,整理得,其中,(3-16),35/57,3.1.6,三次样条插值多项式,2、第二种边界问题,边界条件,次边界条件表明,M,0,和,M,n,是已知,所以,(3.15),式直接改写为,(3-17),36/57,3.1.6,三次样条插值多项式,2、第二种边界问题,边界条件 a、,b、,并注意到 ,得,记,37/57,3.1.6,三次样条插值多项式,则方程可简写成,将上式和 与,(3.15),式合在一起,满足第1或第2或第3种边界条件三次样条插值函数,S(x),是存在且唯一。
15、3-18),38/57,3.1.6,三次样条插值多项式,例3-3,:,已知函数,y=f(x),函数值以下,x -1.5 0 1 2,y 0.125 -1 1 9,求三次样条函数,S(x),,使其满足边界条件 。,解:这是第一个边界条件问题,(3-14)式和(3-16)式个参数为,39/57,3.1.6,三次样条插值多项式,方程组为,求解后得,由,(3-10),式得,40/57,3.2,曲线拟合,3.2.1,问题描述,试验目标是寻找一些物理量之间关系,如一组试验数据(,x,i,y,i,),要寻找函数 一个近似表示式 ,这就是一个曲线拟合问题。,插值方法能够在一定程度上处理曲线拟合问题,但有显
16、著不足。首先,试验数据有误差,假如进行插值,则所求得曲线必须严格经过全部试验数据点,从而使得该曲线保留了试验误差;其次,试验数据往往很多,用插值法求得表示式所以缺乏实用性。,曲线拟合还有一个好坏标准问题,标准不一样决定着不一样拟合方法。惯用有最小二乘法,其基本思想是:使拟合曲线与试验点之差平方和最小。,41/57,3.2,曲线拟合,3.2.2,概述,试验数据表,试验编号,自变量 因变量,1 x,1,y,1,2 x,2,y,2,3 x,3,y,3,i x,i,y,i,m x,m,y,m,42/57,3.2,曲线拟合,设x和y之间为线性关系,则有:,经过试验数据 来确定(3-19)式中常数a和b,
17、即建立x与y关系,此为一个一元线性拟合问题(线性回归)。,问题:怎样用m组数据确定线性方程中a和b?见下列图。,(3-19),43/57,3.2,曲线拟合,怎样确定哪一条直线是最好?,a.,最小二乘法:使回归残差平方和最小。,b,.残差:试验数据y,i,与回归方程计算f(x,i,)之间差,用q,i,表示:,按最小二乘法定义:,(3-20),(3-21),44/57,3.2,曲线拟合,2.算法,由(,3-21,)式可知,,Q,是,a,和,b,函数,即,依据数学知识,要使,Q,最小,有:,将(,3-21,)式代入前两式:,45/57,3.2,曲线拟合,经过计算得到:,(3-22),(3-23),4
18、6/57,3.2,曲线拟合,由(3-22)和(3-23)式可验证:,3.方差分析,怎样衡量回归直线与试验数据之间吻合程度?或回归直线可信度是多少?,-这个问题由,方差分析来完成。,衡量标准:(几个统计量),a.残差平方和,b.回归平方和 其中,47/57,3.2,曲线拟合,c.,剩下标准差 其中n=1,d.相关系数,e.综合检验值,普通惯用是:相关系数R(,R:0=R=1;R越靠近1越好,)和综合检验值F(,F为显著性检验,它值越大越好,)。,48/57,3.2,曲线拟合,多元线性最小二乘法:自变量个数大于一,如:,当,Q/,b,i,=0时,得到一个线性方程组,用以求出 b,i,非线性最小二乘
19、法:将非线性方程回归转化为线性形式。,如 转化为 即相当于,若非线性方程不能回归转化为线性形式时,可采取后面求根法中迭代法。,49/57,3.2,曲线拟合,3.2.3 程序设计,程序功效是:由已知m组实测数据(x,i,y,i,),i=1,2,3,m,按最小二程序原理拟合多项式系数。其中用户能够依据详细问题预先确定拟合多项式次数,也能够按精度要求,使程序自动选择多项式次数。应该注意是,不论何种选择,其次数都不应该超出 m。,1、变量及数组使用说明,M:实测数据(x,i,y,i,)个数;,N:拟合多项式次数加1,例初值为2,即为一次多项式;,EPS:允许误差;,N0:用户直接选择多项式次数,=0时
20、表示按精度要求自动选择多项式次数;,50/57,3.2,曲线拟合,X(M):存放给定点值;,Y(M):存放给定点实测值;,F(M):存放给定点所拟合函数值;,S(2M):存放运算中形成正规方程组系数全部值;,A(M,M):存放正规方程组系数增广矩阵;,B(M):存放正规方程组常系数列向量;,Z(M):调用解方程组子程序得到解向量;,C(M):存放多项式系数列向量;,Q(M):存放残量之绝对值。,2、程序流程图,51/57,开始,输入 N0,M,X(I),Y(I),N0=0?,N=2,N=N0+1,计算S(K),K=1,2,2N-1,生成S(I+J-1)=A(I,J),I=1,N;J=1,N,计
21、算B(K)=,A(K,N+1),K=1,N,调用主元法解方程组子程序,解向量Z(I),X(I)=Z(I),I=1,N,计算拟合曲线在给定点X(I)上函数值 F(I),I=1,M,计算残量,Q(I)=ABS(F(I)-Y(I),寻找最大值Q(I,0,),N0=0?,Q(I,0,)/F(I,0,)EPS?,N=N+1,NQ0)THEN,510 Q0=Q(I),520 I0=I,54/57,3.2,曲线拟合,520 ENDIF 670 RETURN,530 END DO,540,IF,(N0=0),THEN,550,IF,(Q0/F(I0)=EPS),THEN,570 N=N+1,580,IF,(N
22、0,a,b为待定参数),使它能和以下观察值拟合。,x 1 2 3 4 5 6 7 8,y 15.3 20.5 27.4 36.6 49.1 65.6 87.8 117.6,2、已知函数有y=f(x)函数值以下,x 0 1 2 3 4,y -8 -7 0 19 56,求三次样条插值函数S(x),使满足边界条件:S(0)=0,S(4)=48,56/57,上机实习题,画出牛顿插值多项式程序结构流程图并用一个高级语言编制对应程序计算函数f(x)近似牛顿插值多项式N,5,(x),x=0.2+0.01k,k=020,x 0.2 0.24 0.28 0.32 0.36 0.40,y 0.1987 0.2377 0.2764 0.3146 0.3523 0.3894,57/57,






