正文牛顿插值法

牛顿插值法

摘要：值法利用函数f (x)在某区间中若干点的函数值，作出适当的特定函数，在这些点上取已知值，在区间的其他点上用这特定函数的值作为函数f (x)的近似值。如果这特定函数是多项式，就称它为插值多项式。利用插值基函数很容易得到拉格朗日插值多项式，公式结构紧凑，在理论分析中甚为方便，但当插值节点增减时全部插值基函数均要随之变化，整个公式也将发生变化，这在实际计算中是很不方便的，为了克服这一缺点，提出了牛顿插值。

牛顿插值通过求各阶差商，递推得到的一个公式：f(x)=f[x0]+f[x0,x1](x-x0)+f[x0,x1,x2](x-x0)(x-x1)+...f[x0,...xn](x-x 0)...(x-xn-1)+Rn(x)

关键词：牛顿插值法流程图程序实现

一、插值法的由来

在许多实际问题及科学研究中，因素之间往往存在着函数关系，然而，这种关系经常很难有明显的解析表达，通常只是由观察与测试得到一些离散数值。有时，即使给出了解析表达式，却由于表达式过于复杂，不仅使用不便，而且不易于进行计算与理论分析。解决这类问题的方法有两种：一种是插值法,另一种是拟合法。插值法是一种古老的数学方法，它来自生产实践，早在一千多年前，我国科学家在研究历法上就应用了线性插值与二次插值，但它的基本理论却是在微积分产生之后才逐渐完善的，其应用也日益增多，特别是在计算机软件中，许多库函数，如等的计算实际上归结于它的逼近函数的计算。逼近函数一般为只含有算术运算的简单函数，如多项式、有理分式（即多项式的商）。在工程实际问题当中，我们也经常会碰到诸如此类的函数值计算问题。被计算的函数有时不容易直接计算，如表达式过于复杂或者只能通过某种手段获取该函数在某些点处的函数值信息或者导数值信息等。因此，我们希望能用一个“简单函数”逼近被计算函数，然后用该简单函数的函数值近似替代被计算函数的函数值。这种方法就叫插值逼近或者插值法。

逐次线性插值法优点是能够最有效地计算任何给定点的函数值，而不需要写出各步用到的插值多项式的表达式。但如果解决某个问题时需要插值多项式的表达式，那么，它的这个优点就成了它的缺点了。能不能根据插值条件构造一个插值多项式，它既有具体的表达式，又很容易用它计算任何点的函数值呢？牛顿插值法能作到这一点。

二、牛顿插值法的概念

牛顿插值多项式的表达式

设

)

())(())(()(1

1

1

2

1

x

x x x x x c x x x x c x x c c N n n

n

--???--???+--+-+=

问题是如何根据插值条件

()y x N i

i

n

=,i=0,1,2???n

来计算待定系数c c c c n ???210,,？ 由

)

()(00

x y

x N f n

==

知， )

(00

0x y

c f ==。 由 ()(1

1

1

x y

x N f n

==

知

y x x c c 1

1

1

)(=-+

因而

[]

x x x x x x x x y y c f

f f 1

1010

1

11

,)

()(?--=

-

-=，

其中 []x x f 10, 称为函数f(x)在x x 10,点的一阶商。

由

)

()(22

2

x y

x N f n

==

知

y x x x x c x x c c 2

1

2

2

2

1

1

))(()(=--+-+

因而

()()()]

,,[)

(],[],[)

(]

,[)

)(()](,[)](,[)

)(()](,[))(()](,[2

1

2

1

2

1

2

1

12

1

21

2

2

2

1

1

1

1

2

1

2

2

2

1

1

1

2

1

2

2

2

1

2

2

x x x x x x x x x x x x x x x y y x x x x x x x x x x x x y y x x x x x x x x y y y y x x x x x x x x y y c f f f f f f f f ?--=----=

-----+-=

-----+-=

-----=

其中],,[210x x x f 称为函数f (x)在x x x 210,,点的二阶差商。实际上，它是一阶差商的差商。一般地，如果已知一阶差商],[],,[11x x x x i i i i f f +-， 那么就可以计算二

阶差商 x x

x x x x x x x i i i i i i i i i f f f 1

1

1111_]

,[],[],,[-+-++--=

类似于上述过程不断地推导下去，可得

],

,,[)

(]

,,[],[],

,,,,[)

(]

,,,[],,,[],

,,,[)

(]

,,[],,[210012102144321004321043214

3210032103213x x x x x x x x x x x x x c x x x x x x x x x x x x x x x c

x x x x x x x x x x x x c n n n n f f f f f f f f f ????-???-???=

?

????--=

?--=-

其中，

],,,,[3210x x x x f ,,,,[43210x x x x x f ],,,,,,[543210x x x x x x f 分别称为函数

f (x)在相应点处的三阶差商，四阶差商和n 阶差商。实际上，

c

c c c n

???2

1

,, 的计算可通过以下简易地构造函数的差商来完成。 x

c x f 0

0)(=

x

1

)

(1x f

c x x f 1

10],[=

x

2

)

(2x f

]

,[21x x f

c x x x f 2210],,[=

x

3

)

(3x f

]

,[32x x f

]

,1[32

,x x x f

c x x x x f 33210],,,[=

x

4

)

(4x f

]

,[43x x f

]

,,[432x x x f

]

,,1[432

,x x x x f

c x x x x x f 4

43210],,,,[=

. . .

. . .

. . .

. . .

. . .

. . .

按上述方式构造插值多项式的方法叫做牛顿插值法。根据插值多项式的惟一性知，其截断误差与拉格朗日插值法相同， 即：

)

()()!

1(1

π1)

1(x n n n n f

R +++=

ξ

但也可以表示成差商形式。这是因为以x x x x n ???210,,为节点的多项式

)

(],[)()(11101

x x x x N

N

n n n

n f x x ∏???+=

+++

从而 )(],[)()()(111101111x x x x x N x N x n n n n n n n n f f +++++++∏???+==

于是)(x N n 的截断误差可表为

)

(],,,[)(1110x f x n n n

x x x x R ∏++???=

顺便指出，因为牛顿插值多项式具有性质：

)

())(](,[)()(121101

x x x x x x x x x N

N

n n n n

f x x ---???--???+=

所以，类似于逐次线性插值法，也可以把上述和式中的第二项

)

())(](,[12110x x x x x x x x x n n f --???--???

看成是估计

)

(1

x N

n - 的一种实用误差估计式。

与差商概念密切联系的另一个概念是差分，它是指在等距节点上函数值的差。所谓等距节点，是指对给定的常数h （称为步长），节点

)2,1,0(,0

n i ih x x i

???=+=称f

x x

k

i i f f ?

?-+)()(1

为

x i

处的一阶向前差

分；

称 f

x x i

i i f f ?

?--)()(1

为

x i

处的一阶向后差分；

称 f

x x i

h i h i f f δ

?-+)1()(22

为

x i

处的中心差分。一阶差分的差分称为二阶差分，

即 f

f f

i

i

i ???

-?

+2

1

称为

x i

处的二阶向前差分。

一般地，m 阶向前和向后差分可定义如下：

????????=?

????

==?-???-??-+--+-3,2,,111111m m m f f f f f f i m i m i i m i m i

三、牛顿插值法的实现

1、【算法】

步骤1：输入节点（xj，yj），精度ξ，计值点xx，f0→p，1→T，1→i；

步骤2：对k=1，2，……，i依次计算k阶均差

f[xi-k,xi-k+1,…,xi] = (f[xi-k+1,…,xi]- f[xi-k,…,xi])/( xi -xi-k )

步骤3：（1）、若| f[x1,…,xi]- f[x0,…,xi-1]|< ξ，则p为最终结果Ni-1(x)，余项Ri-1= f[x0,…,xi](xx-xi-1)T。

（2）、否则(xx-xi-1)*T→T，p+ f[x0,…,xi]*T→p，转步骤4。

步骤4：若i

2、【流程图】

i+1→i

YES

NO

开始 输出ξ，xx ，n 及（x j ，y j ） f 0→g 0，f i →q 0 k= 1，2，……，i

（q k-1－g k-1）（x i －x i-k ）→q k

k

|g i-1-q i-1|<ξ

q i (x i -x i-1)T →R 输出p ，r ，i STOP 1 (x i -x i-1)T →T p+q i *T →p

k = 1，2，……，i

YES i+1→ i NO

q k →g k

k

k

q i (xx-x n-1)*T →R 输出P ，R ，n

STOP 2

3、【程序清单】#include"stdio.h"

#define n 4//牛顿插值的次数

void main()

{

float a[n+1][n+2]={0},s=0,t=1,x;

int i,j;

printf("请输入xi及yi的值//要求先输入xi再输入yi然后输入下一组\n");

for(i=0;i

for(j=0;j<2;j++)

scanf("%f",&a[i][j]);

for(j=1;j

for(i=j;i

a[i][j+1]=(a[i][j]-a[i-1][j])/(a[i][0]-a[i-j][0]);

printf("输出xi，yi及各阶均差\n");

for(i=0;i

{

for(j=0;j

printf("%6.5f ",a[i][j]);

printf("\n");

}

printf("输出牛顿插值表达式\n");

printf("N%d(x)=",n);

for(i=0;i

{

printf("%6.5f",a[i][i+1]);

for(j=0;j

printf("(x-%3.2f)",a[j][0]);

if(i==n)

break;

printf("+");

}

printf("\n");

printf("输入插值点x=");

scanf("%f",&x);

for(i=0;i

for(j=0;j

t*=(x-a[j][0]);s+=a[i][i+1]*t;

}printf("N%d(%4.3f)=%6.5f\n",n,x,s);} 4.【程序实现】

参考文献:

Richard L. Burden, J. Douglas Faires, Numerical Analysis

(Seventh Edition), Brooks Pub. Co.,2001.

2. 蔡大用，白峰杉. 高等数值分析. 清华大学出版社，北京，1998.

3. 邓建中，刘之行. 计算方法（第二版）.西安交通大学出版社，2001.

4. 韩旭里. 数值分析. 中南大学出版社，2003.

致谢

本文得以顺利完成，非常感谢我的指导教师。从论文的选题直到论文的最终完成，他都给予我尽心尽力的指导。老师严谨的治学态度深深地影响着我，对我今后的学习，工作，生活必将产生影响。借此机会，特向老师表示最诚挚的感谢。

感谢南昌工程学院理学系的所有领导和老师。他们严谨的学风，渊博的知识，诲人不倦的品格一直感染和激励着我不断上进，使我大学四年的时光充实而有意义，“自强不息，格物致知”，在这里我所学到的一切，必将使我受益终生。

在本论文的写作中，我也参照了大量的著作和文章，许多学者的科研成果及写作思路给了我很大的启发，在此向这些学者们表示由衷的感谢，感谢我的家人，同学，朋友对我的大力支持，他们的无私奉献，关爱和支持使我能够继续去追求自己的人生理想和目标。感谢所有关心，帮助和支持我的人。