第4章 福大科学工程与计算-插值法

科学和工程计算第4章插值法插值法插值法是一种古老的数学方法，早在一千多年前的隋唐时期定制历法时就广泛应用了二次插值。刘焯将等距节点的二次插值应用于天文计算。插值理论却是在17世纪微积分产生后才逐步发展起来的，Newton插值公式理论是当时的重要成果。由于计算机的使用以及航空、造船、精密仪器的加工，插值法在理论和实践上都得到进一步发展，获得了广泛的应用。引言拉格朗日插值均差与牛顿插值公式埃尔米特插值分段低次插值三次样条插值插值法求拟合这组数据的多项式能否存在一个性能优良、便于计算的函数一、插值问题------(1)这就是插值问题,(1)式为插值条件,其插值函数的图象如图整体误差的大小反映了插值函数的好坏为了使插值函数更方便在计算机上运算,一般插值函数都使用代数多项式和有理函数。x0x1x2x3x4xP(x)

f(x)Lagrange插值多项式为了求得便于使用的简单的插值多项式P(x),我们先讨论n=1的情形要求线性插值多项式L1(x),使它满足：L1(x)的几何意义就是通过这两点的直线；n=2的情况，假定插值节点为考虑通过n+1个节点……？n

1希望找到li(x)，i=0,…,n

使得

li(xj)=ij

===niiinyxlxL0)()(li(x)每个li

有n

个根x0…

xi…xnLagrangePolynomial与有关，而与无关节点f然后令例1:解:且在例1中,如果只给出两个节点169和225,也可以作插值多项式,即1次Lagrange插值多项式,有两个插值基函数,这种插值方法称为Lagrange线性插值,也可以在n+1个节点中取相邻的两个节点作线性插值例2.解:Lagrange插值基函数为Lagrange线性插值多项式为所以Lagrange插值多项式的缺点：插值基函数计算复杂高次插值的精度不一定高高次插值通常优于低次插值但绝对不是次数越高就越好，嘿嘿……三、插值余项Remainder满足不会完全成立因此,插值多项式存在着截断误差,那么我们怎样估计这个截断误差呢？令设其中近似函数误差根据Rolle定理,再由Rolle定理,依此类推由于因此所以定理2.Lagrange型余项余项表达式只有在f(x)的高阶导数存在时才能应用。设则例3:解:练习1：练习2：第4章插值和拟合均差与牛顿插值公式Lagrange插值多项式的缺点我们知道,Lagrange插值多项式的插值基函数为理论分析中很方便，但是当插值节点增减时全部插值基函数就要随之变化，整个公式也将发生变化，这在实际计算中是很不方便的；Lagrange插值虽然易算，但若要增加一个节点时，全部基函数li(x)都需重新算过。解决由线性代数的知识可知,任何一个n次多项式都可以表示成共n+1个多项式的线性组合那么，是否可以将这n+1个多项式作为插值基函数呢？显然，多项式组线性无关，因此，可以作为插值基函数基函数有再继续下去待定系数的形式将更复杂。。。。。。为此引入差商的概念差商(亦称均差)/*

divideddifference*/1阶差商

/*the1stdivideddifferenceoffw.r.t.xi

andxj

*/2阶差商定义2.11101010111010],,...,[],,...,[],,...,[],...,,[],...,[++--+++--=--=kkkkkkkkkkkxxxxxfxxxfxxxxxfxxxfxxf(k+1)阶差商差商的计算方法(表格法):规定函数值为零阶差商差商表差商具有如下性质:Warning:myheadisexploding…Whatisthepointofthisformula?差商的值与xi

的顺序无关！Newton插值公式Newton插值公式及其余项例：已知x=1,4,9的平方根为1,2,3，利用牛顿基本差商公式求的近似值。解：从而得二阶牛顿基本差商公式为因此计算得的近似值为kxkf(xk)一阶差商二阶差商三阶差商四阶差商0123412345147863301-1-1/3-2-3/2-1/61/24二、代数插值多项式的存在唯一性且满足--------(2)--------(3)--------(4)上述方程组的系数行列式为n+1阶范德蒙行列式由Cramer法则,线性方程组(4)有唯一解定理1.

则满足插值条件的插值多项式存在且唯一.注：若不将多项式次数限制为n

，则插值多项式不唯一。例如也是一个插值多项式，其中可以是任意多项式。练习2.已知函数在点的函数值，求其三次插值多项式。解：其三次插值多项式即为函数本身：（3）已知100,121,144的平方根，计算115平方根近似值x100121144y101112§4.3埃尔米特插值/*HermiteInterpolation*/§3HermiteInterpolation

求Hermite多项式的基本步骤：写出相应于条件的hi(x)、hi(x)的组合形式；对每一个hi(x)、hi(x)找出尽可能多的条件给出的根；根据多项式的总阶数和根的个数写出表达式；根据尚未利用的条件解出表达式中的待定系数；最后完整写出H(x)。习题1:求不超过3次的多项式H(x),使它满足插值条件习题2:求不超过2次的多项式H(x),使它满足插值条件§4分段低次插值/*piecewisepolynomialapproximation*/RememberwhatIhavesaid?IncreasingthedegreeofinterpolatingpolynomialwillNOTguaranteeagoodresult,sincehigh-degreepolynomialsareoscillating.例：在[5,5]上考察的Ln(x)。取

-5

-4

-3

-2

-1

0

1

2

3

4

5

-0.5

0

0.5

1

1.5

2

2.5

n

越大，端点附近抖动越大，称为Runge现象Ln(x)f(x)分段低次插值§4PiecewisePolynomialApproximation

分段线性插值

/*piecewiselinearinterpolation*/在每个区间上，用1阶多项式

(直线)逼近f(x):记，易证：当时，一致失去了原函数的光滑性。

分段Hermite插值

/*Hermitepiecewisepolynomials*/给定在上利用两点的y及y’构造3次Hermite函数导数一般不易得到。Howcanwemakeasmoothinterpolationwithoutaskingtoomuchfromf?Headache…三次样条插值早期工程师制图时，把富有弹性的细长木条（所谓样条）用压铁固定在样点上，在其它地方让它自由弯曲，然后画下长条的曲线，称为样条曲线。它实际上是由分段三次曲线并接而成，在连接点即样点上要求二阶导数连续，从数学上加以概括就得到数学样条这一概念。

三次样条插值样条本质上是一段一段的三次多项式拼合而成的曲线在拼接处,不仅函数是连续的,且一阶和二阶导数也是连续的一、三次样条插值函数定义1.------(1)注：三次样条与分段Hermite插值的根本区别在于S(x)自身光滑，不需要知道f的导数值（除了在2个端点可能需要）；而Hermite插值依赖于f在所有插值点的导数值。f(x)H(x)S(x)要求出S(x)，则在每个小区间上要确定4个待定系数，共有n个小区间，所以应确定4n个参数。共(n+1)+(3n-3)=4n-2个条件，因此还需要两个条件才能确定S(x)可在区间端点a,b上各加一个条件（边界条件），具体要根据实际问题要求给定；已知两端的一阶导数值2.两端的二阶导数已知其特殊情况为3.当f(x)是为周期的周期函数时，则要求S(x)也是周期函数，这时边界条件应满足：这样确定的样条函数S(x)称为周期样条函数；加上任何一类边界条件(至少两个)后一般使用第一、二类边界条件,常用第二类边界条件样条插值函数的建立即或可直接利用分段三次Hermit插值，只要假定可得加以整理后可得------(10)由条件------(11)由于以上两式相等,得基本方程组如果问题要求满足第一类(一阶)边界条件:即------(12)基本方程组化为n-1阶方程组将上式化为矩阵形式------(13)------(14)这是一个三对角方程组如果问题要求满足第二类(二阶自然)边界条件:由(11)式,可知------(15)----(16)------(17)------(18)与基本方程组(12)联合,并化为矩阵形式,得-----(19)(19)式与(14)一样,都是三对角方程组,解是唯一的；例1.对于给定的节点及函数值解:由(12)式可得由(19)式得基本方程组将上述结果代入(10)式定理.

最后，介绍一个有用的结果小结曲线拟合当函数只在有限点集上给定函数值，要在包含该点击的区间上用公式给出函数的简单表达式，这些都涉及到在区间[a,b]上用简单函数逼近已知复杂函数的问题，这就是函数逼近问题。插值法就是函数逼近问题的一种拟解决的问题：计算复杂的函数值已知有限点集上的函数值，给出在包含该点集的区间上函数的简单表达式函数逼近——对函数类A中给定的函数f(x)，记作要求在另一类简单的便于计算的函数类B中求函数使p(x)与f(x)的误差在某种度量意义下最小。逼近问题函数逼近曲线拟合实例：考察某种纤维的强度与其拉伸倍数的关系,下表是实际测定的24个纤维样品的强度与相应的拉伸倍数是记录:纤维强度随拉伸倍数增加而增加并且24个点大致分布在一条直线附近必须找到一种度量标准来衡量什么曲线最接近所有数据点(1)仍然是已知x1…xm

;y1…ym,求一个简单易算的近似函数P(x)

f(x)。但是①

m

很大；②

yi

本身是测量值，不准确，即yi

f(xi)这时没必要取P(xi)=yi,

而要使P(xi)yi

总体上尽可能小。使误差在某种度量意义下最小常见做法：

使最小/*minimaxproblem*/

太复杂使最小不可导，求解困难使最小/*Least-Squaresmethod*/最小二乘法的基本概念一般使用称为平方误差从而确定(1)中的待定系数注意(1)式是一条直线因此将问题一般化仍然定义平方误差我们选取的度量标准是---------(2)---------(3)法方程组由可知因此可假设因此求最小二乘解转化为二次函数由多元函数取极值的必要条件得即---------(4)引入记号则由内积的概念可知---------(5)---------(6)显然内积满足交换律方程组(4)便可化为---------(7)将其表示成矩阵形式-----(8)并且其系数矩阵为对称阵所以法方程组的系数矩阵非奇异,即根据Cramer法则,法方程组有唯一解即是的最小值所以因此作为一种简单的情况,基函数之间的内积为平方误差例1.回到本节开始的实例，从散点图可以看出纤维强度和拉伸倍数之间近似与线性关系故可选取线性函数为拟合函数，其基函数为建立法方程组根据内积公式，可得法方程组为解得平方误差为拟合曲线与散点的关系如右图:例2.求拟合下列数据的最小二乘解x=.24.65.951.241.732.012.232.522.772.99y=.23-.26-1.10-.45.27.10-.29.24.561解:从数据的散点图可以看出因此假设拟合函数与基函数分别为6.7941-5.347563.2589-5.34755.1084-49.008663.2589-49.00861002.5

1.6163-2.382726.7728通过计算,得法方程组的系数矩阵及常数项矩阵为Go!用Gauss列主元消去法,得

-1.0410-1.26130.030735拟合的平方误差为图象如图例3.在某化学反应里,测得生成物浓度y%与时间t的数据如下，试建立y关于t的经验公式x=1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16y=4.00,6.40,8.00,8.80,9.22,9.50,9.70,9.86,10.00,10.20,10.32,10.42,10.50,10.55,10.58,10.60解:具有图示的图形的曲线很多，本题特提供两种形式例：xy(xi,yi),i=1,2,…,m方案一：设baxxxPy+=)(求a

和b

使得最小。=-+=miiiiybaxxba12)(),(jButhey,thesystemofequationsforaandbisnonlinear!Takeiteasy!Wejusthavetolinearizeit…线性化

/*linearization*/：令，则bXaY+就是个线性问题将化为后易解a

和b。),(iiYX),(iiyx方案二：设xbeaxPy/)(-=(a>0,b>0)线性化：由可做变换xbay-lnlnbBaAxXyY-====,ln,1,lnBXAY+就是个线性问题将化为后易解A

和B),(iiYX),(iiyx两边取对数,得得即为拟合函数基函数为解法方程组得平方误差为用最小二乘法得即无论从图形还是从平方误差考虑在本例中指数函数拟合比双曲线拟合要好平方误差为定义权函数：①

离散型/*discretetype*/根据一系列离散点拟合时，在每一误差前乘一正数wi

，即误差函数

，这个wi

就称作权