计算方法1插值方法

第一章插值方法1.1问题的提法1.2拉格朗日插值公式1.5牛顿插值公式1.8样条插值1.6埃尔米特插值1.7分段插值法1.3插值余项1.9曲线拟合的最小二乘法1.1问题的提法函数y=f(x)给出一组函数值

x:x0x1x2……xny:y0y1y2……yn其中x0,x1,x2,…,xn是区间[a,b]上的互异点，要构造一个简单的函数p(x)作为f(x)的近似表达式，使满足

(插值原则、插值条件)

这类问题称为插值问题。p(x)-----f(x)的插值函数,

f(x)-----被插值函数x0,x1,x2,…,xn-----插值节点

求插值函数的方法称为插值法。

若x∈[a,b]，需要计算f(x)的近似值p(x),则称x为插值点。

1.插值问题当选择代数多项式作为插值函数时，称为代数多项式插值问题:代数多项式插值问题:设函数y=f(x)在[a,b]有定义,且已知在n+1个点a≤x0<x1<……<xn≤b上的函数值y0,y1,……,yn.,要求一个次数不高于n的多项式

使满足插值原则称pn(x)为f(x)的n次插值多项式。本章只讨论多项式插值与分段插值。因为多项式具有一些很好的特性，如它具有各阶导数，计算多项式的值比较方便，多项式四则运算后仍是多项式等等。定理2在n+1个互异基点处满足插值原则且次数不超过n的多项式pn(x)是存在并唯一的。证

其系数行列式

因此方程组存在唯一的解

，因此pn(x)存在并唯一。2插值多项式的存在唯一性1.2拉格朗日插值多项式1线性插值3一般情形2抛物插值1线性插值----n=1时的代数多项式插值已知f(x0)=y0，f(x1)=y1,x0≠x1y0y1yx0x1x要构造线性函数p1(x),使它满足插值条件

p1(x0)=y0,p1(x1)=y1.（线性插值多项式）

（拉格朗日线性插值多项式）

公式的结构：它是两个一次函数的线性组合

（线性插值基函数）

例1已知

1011y100121x解

与精确值比较，这个结果有3位有效数字.的精确值为10.723805…,基函数的性质2

抛物插值----n=2时的代数多项式插值已知f(x)在三个互异节点x0,x1,x2的函数值y0,y1,y2xx0x1x2yy0y1y2要构造次数不超过二次的多项式p2(x)

使满足插值条件

公式的构造：采用基函数方法构造p2(x)，先构造三个二次插值基函数

lj(x)(j=0,1,2)，使满足且lj(x)(j=0,1,2)，是一个二次函数.先构造l0(x)因它有两个零点x1及x2故可表为其中c为待定系数，由条件l0(x0)=1,求得

于是，得同理可得（拉格朗日二次插值多项式）于是求得显然，它满足

例2利用100,121和144的开方值求

x100121144y101112解

这个结果同精确值比较，有4位有效数字.二次插值也称之为抛物插值。当三点(x0,y0),(x1,y1),(x2,y2)位于一条直线上时，显然插值函数的图形是直线。3一般情形下面讨论通过n+1个节点x0<x1<……<xn的n次插值多项式pn(x),假设它满足条件

定义:若n次多项式lj(x)(j=0,1,…,n)在n+1个节点x0<x1<……<xn上满足条件

就称这n+1个n次多项式l0(x)，l1(x)，…,

ln(x)为节点x0，x1，…，xn

上的n次插值基函数。

类似n=1,n=2的情况，可得n次插值基函数为显然，它满足条件n次拉格朗日插值多项式

由lk(x)的定义，知

若引入记号注意：

n次插值多项式pn(x)通常是次数为n的多项式，特殊情况次数可能小于n，例如，通过三点（x0,y0),（x1,y1),（x2,y2),的二次插值多项式p2(x)，如果三点共线，则y=p2(x)就是一直线，而不是抛物线，这时p2(x)是一次式。容易求得于是pn(x)又可改写成1.3插值余项截断误差：

也称为插值多项式的余项。

定理3设区间[a,b]含有节点x0，x1，…，xn,而函数f(x)在[a,b]内有连续的直到n+1阶导数，且f(xi)=yi(i=0,1,…,n)已给，则当x∈[a,b]时，对于满足插值条件的n次插值多项式pn(x)，成立其中ξ∈(a,b)且依赖于x，

证

仅需要考察插值点x不是插值节点xi的情形，否则插值余项公式显然成立。令式中c为待定系数，而这样误差函数R(t)=f(t)-g(t)至少有n+2个互异的零点x,x0,x1,…,xn。根据罗尔定理,R(n+1)(ξ)=0依此类推，可知R(n+1)(t)在(a,b)内至少有一个零点ξ，使R”(t)在(a,b)内至少有n个互异零点R’(t)在(a,b)内至少有n+1个互异的零点另一方面直接对R(x)求导知上式令t=ξ,并将c的表达式代入，得

据此稍加整理，就得到余项表达式，证毕由于xi(i=0,1,…,n)都是

（t)的零点，根据插值条件有此外，若取则又有g(x)=f(x).

推论当f(x)是次数不超过n的多项式时，其n次插值多项式就是f(x)本身。证明因为f(n+1)(x)=0，从而Rn(x)≡0，pn(x)≡f(x)。误差公式的用法：，则截断误差估计为

如果设min{x0,x1,…,xn}=a，max{x0,x1,…,xn}=b当插值点x∈(a,b)时称为内插，否则称为外插。特别，当n=1时，线性插值余项为

当n=2时，抛物线插值余项为

计算框图：

1.5牛顿插值公式2差商及其性质1具有承袭性的插值公式3差商形式的插值公式先考察线性插值的点斜式表达式为：

由于1具有承袭性的插值公式

p1(x

)=p0(x)+c1(x-x0)其中，修正项的系数显然，不管系数c2如何取值p2(x)均能满足

p2(x0)=f(x0),p2(x1)=f(x1)再利用剩下的一个条件p2(x2)=f(x2)来确定c2,结果有

可看作是零次插值多项式，上式表明再修正p1(x)以进一步得到抛物线插值多项式p2(x

),为此令p2(x

)=p1(x)+c2(x-x0)(x-x1)以上论述表明，为了建立具有承袭性的插值公式，需要引进差商并研究性质。记c0=f(x0),从而有一阶差商定义为：

一般地,n阶差商定义为为统一起见，定义零阶差商为函数值本身，即f(xi),2差商及其性质二阶差商定义为一阶差商的差商：

差商有如下的基本性质:

例如这个性质可用归纳法证明.这个性质也表明,差商具有对称性,即任意改变节点的次序后其值不变。即

k阶差商可表为函数值f(x0),f(x1),…,f(xk)的线性组合,即

定理4若f(x)在[a,b]上存在n阶导数,且节点则n阶差商与导数关系如下证作函数pn-1(x)是过的插值多项式,由定理3的证明得证差商表

xif(xk)1阶2阶3阶4阶x0f(x0)

x1f(x1)f(x0,x1)

x2f(x2)f(x1,x2)f(x0,x1,x2)

x3f(x3)f(x2,x3)f(x1,x2,x3)f(x0,x1,x2,x3)

x4f(x4)f(x3,x4)f(x2,x3,x4)f(x1,x2,x3,x4)f(x0,x1,x2,x3,x4)

┊┊┊┊┊┊……计算规律：任一个k(≥1)阶差商的数值等于一个分式的值，其分子为所求差商左侧的数减去左上侧的数，分母为所求差商同一行最左边的基点值减去由它往上数第k个基点值。

注意：差商表中，对角线上的差商是构造牛顿型插值公式的重要数据。粗线框出的部分在计算机上可存入二维数组差商表的数据构成一个矩阵F：F00=f(x0)F10=f(x1),F11=f[x0,x1]F20=f(x2),F21=f[x1,x2],F22=f[x0,x1,x2]F30=f(x3),F31=f[x2,x3],F32=f[x1,x2,x3],F33=f[x0,x1,x2,x3]Fi,j-1=f[xi-j+1,…,xi]Fi-1,j-1=f[xi-j,,…,xi-1]计算机上计算均差表的公式

一般有Fi,j=f[xi-j,xi-j+1,…,xi-1,xi]

例3已知函数y=f(x)的观测数据如表，试构造差商表，并求f(2,4,5)及f(2,4,5,6)的值。

x02456f(x)159-413解

n=4,构造差商表

xif(xi)1阶2阶3阶4阶0245621159-4132-13170-515-15f(2,4,5)=-5f(2,4,5,6)=53差商形式的插值公式根据差商定义，把x看成[a,b]上一点，可得牛顿均差型线性插值多项式只要把后一式代入前一式，就得到我们称p

n(x)为牛顿均差插值多项式计算牛顿差商插值多项式的步骤：（1）作差商表（2）根据公式计算牛顿型插值多项式（表中对角线上各差商值就是pn(x)的各项系数）。

余项公式

牛顿型插值多项式

p

n(x)显然满足插值条件，且就是次数不超过n次的插值多项式，系数例4已知函数y=f(x)的观测数据如上例，试用全部节点构造牛顿插值多项式，并用二次插值求f(3)的近似值。

解用全部基点时，n=4，先作差商表,见上例。

p4(x)=f(0)+f(0,2)(x-0)+f(0,2,4)(x-0)(x-2)+f(0,2,4,5)(x-0)(x-2)(x-4)+f(0,2,4,5,6)(x-0)(x-2)(x-4)(x-5)xif(xi)1阶2阶3阶4阶0245622-1317159-4130-515-151=1+2x-x(x-2)(x-4)+x(x-2)(x-4)(x-5)用二次插值求f(3)时n=2，x=3，作内插取

x0=2,x1=4,x2=5f(3)≈p2(3)=f(2)+f(2,4)(3–2)+f(2,4,5)(3-2)(3-4)=7-5(3-2)(3-4)=12

不少实际问题不但要求在节点上函数值相等，而且还要求它的导数值也相等，满足这种要求的插值多项式就是Hermite插值多项式。这里不准备对Hermite插值作一般性论述,而仅仅研究两个具体问题.1.求作二次式p2(x)满足设用这一插值多项式p2(x)逼近某个取值解决这类问题与处理不带导数的插值一样,也有两种方法下面分别说明:1.6埃尔米特插值(1)基于承袭性,按牛顿插值特点,令则不管系数c怎样取值,总有将c代入,得再用剩下的一个条件确定c,(2)用基函数方法，为简化计算，先设x0=0,x1=1,而令均为二次式，它们分别满足其中基函数下面我们将利用以上条件来分别确定由条件0(1)=0知x=1是0(x)的一个根，故可令将0(x)的其它两个条件代入，得方程组解得a=b=-1,于是有0(x)=1-x2同理可得1(x)=x2，0(x)=x(1-x)若x0,x1为任意节点，那么可令x1-x0=h，不难验证，这时2.求作三次式p3(x)，使满足仿照问题1的解法，记h=x1-x0,而令仿照问题1的解法，导出其插值基函数为余项：对于问题1和问题2的插值余项分别是其中1，2均包含在由点x0,x1和x所界定的范围内

例5求满足的插值多项式及其余项表达式。解：由给定条件，可确定次数不超过3的插值多项式，由于此多项式通过点故其形式为其中A为待定常数，可由条件确定，通过计算得:为了求余项R(x)=f(x)-P(x)的表达式,设其中K(x)为待定函数，构造

显然故其形式为(t)在（a,b）内有5个零点（重根算两个）。反复应用罗尔定理，的(4)(t)在（a,b）内至少有一个零点，故于是，余项公式可表为1高次插值的龙格现象2分段插值的概念3分段线性插值1.7分段插值法4分段三次插值1高次插值的龙格现象1901年龙格首先发现多项式插值有危险，他试图在区间[-5,5]内相等间隔的节点上用多项式对函数进行插值，却发现当插值多项式pn(x)的次数趋于无穷时，pn(x)在|x|<3.63内收敛，而在该区间之外发散。如下图所示。高次代数多项式插值的龙格（Runge)现象所以，七、八次以上的代数多项式插值很少使用。2分段插值的概念

所谓分段插值就是将插值函数逐段多项式化。设已知节点a≤x0<x1<……<xn≤b上的函数值y0,y1,……,yn.，求函数SK(x)满足：1ºSk(x)C[a,b],2ºSk(xk)=yk,(k=0,1,…,n)3ºSk(x)在每个小区间[xk,xk+1]上是k次式，

则称Sk(x)为分段k次式.3分段线性插值设给定f(x)函数值(xi,yi),i=0,1,…,n,求在分划：a≤x0<x1<……<xn≤b下的一次式S1(x)，满足条件S1(xi)=yi,i=0,1,…,n

由定义可知S1(x)在每个小区间[xi,xi+1]上可表示为或表示为式中hi=xi+1-xi,而0(x)=1-x，1(x)=x在每个小区间[xi,xi+1]上有估计式而所以定理5当f(x)C2[a,b]且x[a,b]时，有估计式由此得知，S1(x)在[a,b]上一致收敛到f(x).

习题29对函数

解

xi+1/2f(xi+1/2)S

(

xi+1/2

)相对误差0.51.52.53.54.50.750000.033650.800000.307690.137930.075470.047060.350000.150000.079410.048640.137500.055200.062500.05520进行分段线性插值。4分段三次Hermite插值若在节点xk(k=0,1,…,n)上除已知函数值yk外还已知导数值,就可构造一个导数连续的分段插值函数S3(x)，它满足：1ºS3(x)C1[a,b],2ºS3(xi)=yi,3ºS3(x)在每个小区间[xi,xi+1]上是三次多项式，根据两点三次Hermite插值多项式可知，S3(x)在每个小区间[xi,xi+1]上可表示为定理6若f

(x)C4[a,b],且x[a,b],有估计式

分段插值法的利弊

分段插值法是一种显式算法，其算法简单，而且收敛性能得到保证，只要节点间距充分小，分段插值法总能获得所要求的精度，而不会象高次插值那样发生龙格现象。分段插值法的另一个重要特点是它的局部性质，如果修改某个数据，那么插值曲线仅仅在某个局部范围内受到影响，而代数插值却会影响整个插值区间。我们看到，和分段线性插值相比较，分段三次Hermite插值虽然改善了精度，但这种插值要求给出节点上的导数值，所要提供的信息“太多”，同时它的光滑性也不高（只有连续的一阶导数）。改进这种插值以克服其缺点，这就导致了所谓三次样条插值。1样条函数的概念2三次样条插值3计算步骤与例题1.8样条插值1样条函数的概念

样条插值的思想：逐段选取适当的低次多项式，按一定的光滑性要求连接起来构成插值函数。

定义若函数Sk(x)Ck-1[a,b],且在每个小区间[xj,xj+1]上是k次多项式，其中a=x0<x1<……<xn=b是给定节点，则称Sk(x)是节点x0,x1,x2,…,xn上的k次样条函数.若在节点xj上给定函数值yj=f(xj)(j=0,1,…,n),并成立S(xj)=yj,j=0,1,…,n.(*).

则称Sk(x)为三次样条插值函数.

称xoy平面上的点(xi,yi)(i=0,1,…,n)为样点。

从定义知S3(x)在每个小区间[xj,xj+1]上是三次多项式，要确定4个待定系数，给定n+1个样点(xi,yi)(i=0,1,…,n)，共有n个小区间，需要确定4n叁数。在定义中，已指定了3n-3个条件，即所以，一般需在区间[a,b]端点a=x0,b=xn补充2个边界条件。

再加上共有4n-2个条件2三次样条插值常用的边界条件有三种:

1°已知两端的一阶导数值，即

2°两端的二阶导数值已知，即特别取时称为自然边界条件。3°当f(x)是以xn-x0为周期的周期函数时，则要求S3(x)也是周期函数，端点要满足

这样确定的样条函数S3(x)，称为周期样条函数。若假定在节点xj处的值为由分段三次Hermite插值，可得及取边界条件

显然这样造出来的三次多项式S3(x)不论叁数mj如何选取，均可直接验证所有节点不但连续，而且还有连续一阶导数，于是问题就归结到如何选取叁数mj使S3(x)在每个节点上也有连续的二阶导数。为此我们利用条件来确定mj。对S3(x)求二阶导数得于是在每个小区间[xi,xi+1]的左右两端分别有为了保证S(x)在节点xi处有连续的二阶导数，应有即令整理化简为如采用第一种边界条件

m1,m1,…,mn-1的n-1个方程，称为基本方程组方程为只含其系数矩阵

它的系数矩阵的非零元素集中在三条对角线上而被称作是三对角型的，求解这类方程组的有效方法是追赶法。求三次样条插值函数的计算步骤归纳为:

步1输入初始数据3计算步骤与例题步2i从0到n-1计算步3i从1到n-1计算步4用追赶法解三对角方程求出mi（i=1,…,n-1)步5计算S3(x)在若干点上的值，并打印结果。x0123y0000端点条件为：m0=1，m3=0.求三次样条插值函数的分段表达式，并求f(1.5)。

求m1,m2

的方程组形为方程组化为解得例

给定函数表

解

利用公式求得三次样条函数如下：因为1.5[1,2],所以

例给定区间[0,3]上3个点的函数值f(0)=0,f(1)=2,f(3)=4,试求数a,b,c,d,使函数S(x)为给定点上的三次样条插值函数。

解设

根据定义，由得d=0，故则由

得由得由得由得求得

即1.9曲线拟合的最小二乘法1直线拟合2多项式拟合3观测数据的修匀曲线拟合的问题：设函数y=f(x)在n个互异点的观测数据为求一个简单的近似函数φ(x)，使之“最好”地逼近f(x)，而不必满足插值原则。称函数y=φ(x)为经验公式或拟合曲线。

xix1x2…..xnyiy1y2…..yn通常选择函数类型的做法：描出散点图，再根据专业知识和经验来选择φ(x)的类型。

例子：（注意它与插值法的不同）1.直线拟合.

假设所给数据点(xi,yi)(i=1,2,…,N)的分布大致成一直线，我们求作拟合直线尽可能地从所给数据点(xi,yi)附近通过，就是说，近似地成立这里，数据点的数目N》2，，因此，拟合直线的构造，本质上是个解超定（矛盾）方程组的代数问题设表示按拟合直线y=a+bx求得的近似值，两者之差称为残差。显然，残差的大小是衡量拟合好坏的重要标志。构造拟合曲线可以采用下列三种准则之一：分析以上三种准则,(1)和(2)两种由于含有绝对值运算,不便于实际应用,最常用的是准则(3)称作曲线拟合的最小二乘法.按最小二乘法,作直线拟合应使总误差(1)使残差的最大绝对值为最小:(2)使残差的绝对值之和为最小:(3)使残差的平方和为最小:为最小.即解方程组正则方程组为

例已知实验数据表，试用最小二乘法求经验公式拟合这组数据。解作散点图，容易看出数据点接近一条直线，因此设经验公式为2112840y2468x正则方程组为解得得经验公式为y=-12.5+6.55xN=4,2多项式拟合

使总误差

假设所给数据点(xi,yi)(i=1,2,…,N)的分布大致成一曲线，我们求作m(m<<N)次多项式为最小.即解方程组得即有这个关于系数aj的线性方程组通常称为正则方程组。（*）(1)由已知数据画出函数粗略的图形——散点图，确定拟合多项式的次数m;(2)列表计算(3)写出正则方程组，求出a0,a1,…,am;(4)写出拟合多项式多项式拟合的一般方法可归纳为以下几步：定理7设节点x0,x1,…,xn互异,则正则方程组(*)有唯一解.定理8设aj(j=0,1,…,m)为正则方程组(*)的解,则必为满足最小二乘的拟合多项式.例

已知一组观测数据表，试用最小二乘法求一个多项式拟合这组数据。解作散点图，可以看出这些点接近一条抛物线，因此设所求的多项式为其正则方程组为得c1=4.7143,c2=-2.7857,c3=0.5000521