计算方法插值方法

计算方法插值方法第一页，共八十页，编辑于2023年，星期一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.插值问题第二页，共八十页，编辑于2023年，星期一第三页，共八十页，编辑于2023年，星期一当选择代数多项式作为插值函数时，称为代数多项式插值问题:代数多项式插值问题:设函数y=f(x)在[a,b]有定义,且已知在n+1个点a≤x0<x1<……<xn≤b上的函数值y0,y1,……,yn.,要求一个次数不高于n的多项式

使满足插值原则称pn(x)为f(x)的n次插值多项式。

本章只讨论多项式插值与分段插值。因为多项式具有一些很好的特性，如它具有各阶导数，计算多项式的值比较方便，多项式四则运算后仍是多项式等等。第四页，共八十页，编辑于2023年，星期一

定理2

在n+1个互异基点处满足插值原则且次数不超过n的多项式pn(x)是存在并唯一的。证

其系数行列式

因此方程组存在唯一的解

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

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

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

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

（线性插值基函数）

第七页，共八十页，编辑于2023年，星期一

例1

已知

1011y100121x解

与精确值比较，这个结果有3位有效数字.的精确值为10.723805…,第八页，共八十页，编辑于2023年，星期一基函数的性质第九页，共八十页，编辑于2023年，星期一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)，是一个二次函数.第十页，共八十页，编辑于2023年，星期一先构造l0(x)因它有两个零点x1及x2故可表为其中c为待定系数，由条件l0(x0)=1,求得

于是，得同理可得第十一页，共八十页，编辑于2023年，星期一（拉格朗日二次插值多项式）于是求得显然，它满足第十二页，共八十页，编辑于2023年，星期一

例2

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

x100121144y101112解

这个结果同精确值比较，有4位有效数字.第十三页，共八十页，编辑于2023年，星期一

二次插值也称之为抛物插值。当三点(x0,y0),(x1,y1),(x2,y2)位于一条直线上时，显然插值函数的图形是直线。第十四页，共八十页，编辑于2023年，星期一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次插值基函数为第十五页，共八十页，编辑于2023年，星期一显然，它满足条件n次拉格朗日插值多项式

由lk(x)的定义，知

第十六页，共八十页，编辑于2023年，星期一若引入记号注意：

n次插值多项式pn(x)通常是次数为n的多项式，特殊情况次数可能小于n，例如，通过三点（x0,y0),（x1,y1),（x2,y2),的二次插值多项式p2(x)，如果三点共线，则y=p2(x)就是一直线，而不是抛物线，这时p2(x)是一次式。容易求得于是pn(x)又可改写成第十七页，共八十页，编辑于2023年，星期一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为待定系数，而第十八页，共八十页，编辑于2023年，星期一这样误差函数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).第十九页，共八十页，编辑于2023年，星期一

推论当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时，抛物线插值余项为

第二十页，共八十页，编辑于2023年，星期一计算框图：

第二十一页，共八十页，编辑于2023年，星期一第二十二页，共八十页，编辑于2023年，星期一1.5

牛顿插值公式2差商及其性质1具有承袭性的插值公式3差商形式的插值公式第二十三页，共八十页，编辑于2023年，星期一

先考察线性插值的点斜式表达式为：

由于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)第二十四页，共八十页，编辑于2023年，星期一以上论述表明，为了建立具有承袭性的插值公式，需要引进差商并研究性质。记c0=f(x0),从而有第二十五页，共八十页，编辑于2023年，星期一一阶差商定义为：

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

第二十六页，共八十页，编辑于2023年，星期一差商有如下的基本性质:

例如

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

k阶差商可表为函数值f(x0),f(x1),…,f(xk)的线性组合,即第二十七页，共八十页，编辑于2023年，星期一

定理4

若f(x)在[a,b]上存在n阶导数,且节点则n阶差商与导数关系如下证作函数pn-1(x)是过的插值多项式,由定理3的证明得证第二十八页，共八十页，编辑于2023年，星期一差商表

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个基点值。

注意：差商表中，对角线上的差商是构造牛顿型插值公式的重要数据。粗线框出的部分在计算机上可存入二维数组第二十九页，共八十页，编辑于2023年，星期一差商表的数据构成一个矩阵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]第三十页，共八十页，编辑于2023年，星期一

例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)=5第三十一页，共八十页，编辑于2023年，星期一3差商形式的插值公式

根据差商定义，把x看成[a,b]上一点，可得牛顿均差型线性插值多项式只要把后一式代入前一式，就得到第三十二页，共八十页，编辑于2023年，星期一我们称p

n(x)为牛顿均差插值多项式计算牛顿差商插值多项式的步骤：

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

余项公式

牛顿型插值多项式

p

n(x)显然满足插值条件，且就是次数不超过n次的插值多项式，系数第三十三页，共八十页，编辑于2023年，星期一

例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)第三十四页，共八十页，编辑于2023年，星期一用二次插值求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第三十五页，共八十页，编辑于2023年，星期一

不少实际问题不但要求在节点上函数值相等，而且还要求它的导数值也相等，满足这种要求的插值多项式就是Hermite插值多项式。这里不准备对Hermite插值作一般性论述,而仅仅研究两个具体问题.1.求作二次式p2(x)满足设用这一插值多项式p2(x)逼近某个取值解决这类问题与处理不带导数的插值一样,也有两种方法下面分别说明:1.6埃尔米特插值第三十六页，共八十页，编辑于2023年，星期一(1)基于承袭性,按牛顿插值特点,令

则不管系数c怎样取值,总有将c代入,得再用剩下的一个条件确定c,第三十七页，共八十页，编辑于2023年，星期一(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)第三十八页，共八十页，编辑于2023年，星期一若x0,x1为任意节点，那么可令x1-x0=h，不难验证，这时2.求作三次式p3(x)，使满足

仿照问题1的解法，记h=x1-x0,而令第三十九页，共八十页，编辑于2023年，星期一

仿照问题1的解法，导出其插值基函数为余项：对于问题1和问题2的插值余项分别是其中1，2均包含在由点x0,x1和x所界定的范围内第四十页，共八十页，编辑于2023年，星期一

例5

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

显然故其形式为(t)在（a,b）内有5个零点（重根算两个）。反复应用罗尔定理，的(4)(t)在（a,b）内至少有一个零点，故于是，余项公式可表为第四十二页，共八十页，编辑于2023年，星期一1高次插值的龙格现象2分段插值的概念3分段线性插值1.7分段插值法4分段三次插值第四十三页，共八十页，编辑于2023年，星期一1高次插值的龙格现象1901年龙格首先发现多项式插值有危险，他试图在区间[-5,5]内相等间隔的节点上用多项式对函数进行插值，却发现当插值多项式pn(x)的次数趋于无穷时，pn(x)在|x|<3.63内收敛，而在该区间之外发散。如下图所示。第四十四页，共八十页，编辑于2023年，星期一高次代数多项式插值的龙格（Runge)现象所以，七、八次以上的代数多项式插值很少使用。第四十五页，共八十页，编辑于2023年，星期一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次式.第四十六页，共八十页，编辑于2023年，星期一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第四十七页，共八十页，编辑于2023年，星期一在每个小区间[xi,xi+1]上有估计式而所以定理5当f(x)C2[a,b]且x[a,b]时，有估计式由此得知，S1(x)在[a,b]上一致收敛到f(x).第四十八页，共八十页，编辑于2023年，星期一

习题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进行分段线性插值。第四十九页，共八十页，编辑于2023年，星期一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]上可表示为第五十页，共八十页，编辑于2023年，星期一定理6若f

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

分段插值法的利弊

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

从定义知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三次样条插值第五十五页，共八十页，编辑于2023年，星期一常用的边界条件有三种:

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

2°两端的二阶导数值已知，即特别取时称为自然边界条件。3°当f(x)是以

xn-x0为周期的周期函数时，则要求S3(x)也是周期函数，端点要满足

这样确定的样条函数S3(x)，称为周期样条函数。第五十六页，共八十页，编辑于2023年，星期一若假定在节点xj处的值为由分段三次Hermite插值，可得及取边界条件

显然这样造出来的三次多项式S3(x)不论叁数mj如何选取，均可直接验证所有节点不但连续，而且还有连续一阶导数，于是问题就归结到如何选取叁数mj使S3(x)在每个节点上也有连续的二阶导数。为此我们利用条件来确定mj。第五十七页，共八十页，编辑于2023年，星期一对

S3(x)求二阶导数得于是在每个小区间[xi,xi+1]的左右两端分别有第五十八页，共八十页，编辑于2023年，星期一为了保证S(x)在节点

xi处有连续的二阶导数，应有即令整理化简为第五十九页，共八十页，编辑于2023年，星期一如采用第一种边界条件

m1,m1,…,mn-1的n-1个方程，称为基本方程组方程为只含第六十页，共八十页，编辑于2023年，星期一其系数矩阵

它的系数矩阵的非零元素集中在三条对角线上而被称作是三对角型的，求解这类方程组的有效方法是追赶法。第六十一页，共八十页，编辑于2023年，星期一求三次样条插值函数的计算步骤归纳为:

步1输入初始数据3计算步骤与例题步2i从0到n-1计算步3i从1到n-1计算

步4用追赶法解三对角方程求出mi（i=1,…,n-1)步5计算S3(x)在若干点上的值，并打印结果。第六十二页，共八十页，编辑于2023年，星期一第六十三页，共八十页，编辑于2023年，星期一x0123y0000端点条件为：m0=1，m3=0.求三次样条插值函数的分段表达式，并求f(1.5)。

求m1,m2

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

给定函数表

解

第六十四页，共八十页，编辑于2023年，星期一利用公式求得三次样条函数如下：因为1.5[1,2],所以

第六十五页，共八十页，编辑于2023年，星期一

例给定区间[0,3]上3个点的函数值

f(0)=0,f(1)=2,f(3)=4,试求数

a,b,c,d,使函数

S(x)为给定点上的三次样条插值函数。

解设

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

得由得由得由得求得

第六十六页，共八十页，编辑于2023年，星期一即第六十七页，共八十页，编辑于2023年，星期一1.9曲线拟合的最小二乘法1直线拟合2多项式拟合3观测数据的修匀第六十八页，共八十页，编辑于2023年，星期一曲线拟合的问题：设函数y=f(x)在n个互异点的观测数据为

求一个简单的近似函数φ(x)，使之“最好”地逼近f(x)，而不必满足插值原则。称函数y=φ(x)为经验公式或拟合曲线。

xix1x2…..xnyiy1y2…..yn

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

第六十九页，共八十页，编辑于2023年，星期一例子：（注意它与插值法的不同）第七十页，共八十页，编辑于2023年，星期一1.直线拟合.

假设所给数据点(xi,yi)(i=1,2,…,N)的分布大致成一直线，我们求作拟合直线尽可能地从所给数据点(xi,yi)附近通过，就是说，近似地成立这里，数据点的数目N》2，，因此，拟合直线的构造，本质上是个解超定（矛盾）方程组的代数问题设表示按拟合直线

y=a+bx求得的近似值，两者之差称为残差。显然，残差的大小是衡量拟合好坏的重要标志。第七十一页，共八十页，编辑于2023年，星期一

构造拟合曲线可以采用下列三种准则之一：分析以上三种准则,(1)和(2)两种由于含有绝对值运算,不便于实际应用,最常用的是准则(3)称作曲线拟合的最小二乘法.按最小二乘法,作直线拟合应使总误差(1)使残差的最大绝对值为最小:(2)使残差的绝对值之和为最小:(3)使残差的平方和为最小:为最小.即解方程组第七十二页，共八十页，编辑于2023年，星期一正则方程组为

例已知实验数据表，试用最小二乘法求经验公式拟合这组数据。解作散点图，容易看出数据点接近一条直线，因此设经验公式为2112840y2468x正则方程组为解得得经验公式为y=-12.5+6.55xN=4,第七十三页，共八十页，编辑于2023年，星期一2

多项式拟合

使总误差

假设所给数据点(xi,yi)(i=1,2,…,N)的分布大致成一曲线，我们求作m(m<<N)次多项式为最小.即解方程组得即有这个关于系数aj的线性方程组通常称为正则方程组。（*）第七十四页，共八十页，