MATLAB 牛顿插值法例题与程序.docx

题目一：多项式插值
某气象观测站在8：00（AM）开始每隔10分钟对天气作如下观测，用三次多项式插值函数（Newton）逼近如下曲线，插值节点数据如上表，并求出9点30分该地区的温度（x=10）。
x12345678y22.523.324.421.7025.228.524.825.4二、数学原理
假设有n+1个不同的节点及函数在节点上的值（x，y），……（x，y），插值多项式有如下形式：
（1）
其中系数（i=0,1,2……n）为特定系数，可由插值样条（i=0,1,2……n）确定。
根据均差的定义，把x看成[a,b]上的一点,可得
f(x)=f（）+f[]（）
f[x,]=f[]+f[x,]（）
……
f[x,，…x]=f[x,，…x]+f[x,，…x](x-x)
综合以上式子，把后一式代入前一式，可得到：
f(x)=f[]+f[]（）+f[]（）（）+
…+f[x,，…x]（）…(x-x)+f[x,，…x,]=N（x）+
其中
N（x）=f[]+f[]（）+f[]（）（）+
…+f[x,，…x]（）…(x-x)（2）
=f(x)-N（x）=f[x,，…x,]（3）
=（）…(x-x)
Newton插值的系数（i=0,1,2……n）可以用差商表示。一般有
[]（k=0,1,2，……，n）（4）
把（4）代入（1）得到满足插值条件N（i=0,1,2，……n）的n次Newton插值多项式
N（x）=f（）+f[]（）+f[]（）（）+……+f[]（）（）…（）.
其中插值余项为：

介于之间。
三、程序设计
function[y,A,C,L]=newdscg(X,Y,x,M)
%y为对应x的值，A为差商表，C为多项式系数，L为多项式
%X为给定节点，Y为节点值，x为待求节点
n=length(X);m=length(x);%n为X的长度
fort=1:m
z=x(t);A=zeros(n,n);A(:,1)=Y';
s=0.0;p=1.0;q1=1.0;c1=1.0;
forj=2:n
fori=j:n
A(i,j)=(A(i,j-1)-A(i-1,j-1))/(X(i)-X(i-j+1));
end
q1=abs(q1*(z-X(j-1)));c1=c1*j;
end
C=A(n,n);q1=abs(q1*(z-X(n)));
fork=(n-1):-1:1
C=conv(C,poly(X(k)));
d=length(C);C(d)=C(d)+A(k,k);
end
y(k)=polyval(C,z);%输出y值
end
L(k,:)=poly2sym(C);%输出多项式

>>symsM,X=[1,3,5,7];Y=[22.5,24.4,25.2,24.8];x=10;
>>[y,A,C,L]=newdscg(X,Y,x,M)

y=
21.7313
A=

22.5000000
24.40000.950000
25.20000.4000-0.13750
24.8000-0.2000-0.1500-0.0021

C=

-0.0021-0.11871.452121.1688

L=

-x^3/480-(19*x^2)/160+(697*x)/480+3387/160
结果分析和讨论
对于不超过三次的插值多项式，x如果选取1,3,5,7这三个点能够得到较好的三次插值多项式L=-0.0021x^3-0.1187x^2+1.4521x+21.1688。当x=10时，也即9点30分时的温度为21.7317度，结果分析知此值应是偏小的。对于选取不同的插值节点，能够得到不同的插值多项式，误差也不尽相同。
完成题目的体会与收获
对于牛顿插值法有了更深的了解，合理选择插值节点很重要。加深了对其原理的认识，学会了牛顿插值法的matlab编程，对matlab计算方法更加熟悉。通过完成这道题使我受益匪浅。