计算方法 实验一 插值Word格式文档下载.docx
- 文档编号:22550466
- 上传时间:2023-02-04
- 格式:DOCX
- 页数:16
- 大小:156.22KB
计算方法 实验一 插值Word格式文档下载.docx
《计算方法 实验一 插值Word格式文档下载.docx》由会员分享,可在线阅读,更多相关《计算方法 实验一 插值Word格式文档下载.docx(16页珍藏版)》请在冰豆网上搜索。
Ø
拉格朗日插值
插值基函数:
:
牛顿插值
牛顿插值公式:
Pn(x)=a0+a1(x-x0)+…+an(x-x0)(x-x1)…(x-xn-1)
自然样条函数
4)实验内容
联系所学的知识,按照实验题目要求进行编程以及对数值的分析,能够更加深刻的理解拉格朗日插值,多项式插值,和自然样条插值方法及程序设计算法。
5)实验结果
6)实验结果分析与小结
1.使用拉格朗日插值进行计算时,首先计算插值基函数,再构造拉格朗日插值多项式对数据进行计算。
进行准确性分析,采用事后误差估计。
1910年的实际人口约为91772000,判断插值计算得到的1965年和2002年的人口数据准确性。
~~
共有7组数据则取前6组得L(1965)≈178340.5L(2002)≈-3964442.345
取后6组得L(1965)=193082.5
则R(x)≈
R(1965)≈[55/(-60)]*(193081.5-178340.5)≈-13512.5
x=2002时,所得结果明显偏离事实甚远。
可见拉格朗日插值的求值效果不理想。
2.N(1965)≈193096.277344与拉格朗日插直球的结果相近
N(2002)≈15823.279616仍然偏离事实
3.S(1910)≈88219.0
S(1965)≈191873.430622与前两个插值求得的结果相近
S(2002)≈157868.777837比前两个插值求得的结果靠近与事实
在自然样条插值计算的编程中,我用到了Courant法求解线性方程组的方法,主要是用来解决Mn的求值问题。
4.比较以上实验结果,并加以分析可得,自然样条插值相对其他更好。
但是,脱离了给出的x的区间,在两端延伸的区间内进行数值计算所得到的数值依然不够精确,不够理想。
只有在已知[Xi,Xi+1]的区间内计算,才能够得到比较满意的数值。
程序详见附录。
附录:
编程语言采用c语言,运行环境为win-tc。
1.拉格朗日插值
#include<
stdio.h>
#defineMAX_N20
typedefstructtagPOINT
{
doublex;
doubley;
}POINT;
intmain()
intn;
inti,j;
POINTpoints[MAX_N+1];
doublediff[MAX_N+1];
doublex,tmp,newton=0;
printf("
\nInputnvalue:
"
);
scanf("
%d"
&
n);
if(n>
MAX_N)
{
TheinputnislargerthenMAX_N,pleaseredefinetheMAX_N.\n"
return1;
}
if(n<
=0)
Pleaseinputanumberbetween1and%d\n"
MAX_N);
Nowinputthe(x_i,y_i),i=0,...,%d:
\n"
n);
for(i=0;
i<
=n;
i++)
%lf%lf"
points[i].x,&
points[i].y);
Nowinputthexvalue:
%lf"
x);
i++)diff[i]=points[i].y;
n;
for(j=n;
j>
i;
j--)
diff[j]=(diff[j]-diff[j-1])/(points[j].x-points[j-1-i].x);
tmp=1;
newton=diff[0];
tmp=tmp*(x-points[i].x);
newton=newton+tmp*diff[i+1];
newton(%f)=%f\n"
x,newton);
getch();
return0;
}
2.牛顿插值
doublex,tmp=0,lagrange=0,tx,ty;
getch();
diff[i]=0;
tx=1;
ty=1;
for(j=0;
j<
j++)
if(i!
=j)
tx=tx*(x-points[j].x);
ty=ty*(points[i].x-points[j].x);
diff[i]=tx/ty;
{
tmp=points[i].y*diff[i];
printf("
%f"
tmp);
lagrange+=tmp;
lagrange(%f)=%f\n"
x,lagrange);
3.自然样条函数
math.h>
doublex,tmp,sa=0,sb=0,s,t[3];
doubleh[MAX_N],v[MAX_N],u[MAX_N],M[MAX_N+1],d[MAX_N];
/*courant中所用变量*/
intn2;
inti2,j2,k2;
intmi2;
doublemx2,tmp2;
staticdoublea2[MAX_N][MAX_N],b2[MAX_N],x2[MAX_N],y2[MAX_N];
staticdoublel2[MAX_N][MAX_N],u2[MAX_N][MAX_N];
h[i]=points[i+1].x-points[i].x;
for(i=1;
v[i]=h[i]/(h[i]+h[i-1]);
u[i]=1-v[i];
tmp=0;
tmp+=points[i-1].y/((points[i-1].x-points[i].x)*(points[i-1].x-points[i+1].x));
tmp+=points[i].y/((points[i].x-points[i-1].x)*(points[i].x-points[i+1].x));
tmp+=points[i+1].y/((points[i+1].x-points[i].x)*(points[i+1].x-points[i-1].x));
d[i]=tmp*6;
/*线性方程组求解解放在M[i]中,M[0]=M[n]=0*/
n2=n-1;
M[0]=0;
M[n]=0;
for(i2=0;
i2<
n2;
i2++)
for(j2=0;
j2<
j2++)
a2[i2][j2]=0;
if(i2-1>
=0)a2[i2][i2-1]=u[i2+1];
a2[i2][i2]=2;
if(i2+1<
n2)a2[i2][i2+1]=v[i2+1];
b2[i2]=d[i2+1];
/*courant*/
for(i2=0;
i2++)u2[i2][i2]=1;
for(k2=0;
k2<
k2++)
for(i2=k2;
l2[i2][k2]=a2[i2][k2];
=k2-1;
l2[i2][k2]-=(l2[i2][j2]*u2[j2][k2]);
for(j2=k2+1;
u2[k2][j2]=a2[k2][j2];
u2[k2][j2]-=(l2[k2][i2]*u2[i2][j2]);
u2[k2][j2]/=l2[k2][k2];
y2[i2]=b2[i2];
=i2-1;
y2[i2]-=(l2[i2][j2]*y2[j2]);
y2[i2]/=l2[i2][i2];
for(i2=n2-1;
i2>
=0;
i2--)
x2[i2]=y2[i2];
for(j2=i2+1;
x2[i2]-=(u2[i2][j2]*x2[j2]);
i2++){M[n2+1]=x2[i2];
/*printf("
%f\t"
M[n2+1]);
*/}
i=0;
while(x>
=points[i].x&
&
n)i++;
if(i>
0)i--;
sa=points[i+1].x-x;
sb=x-points[i].x;
s=M[i]*sa*sa*sa/(h[i]*6)+M[i+1]*sb*sb*sb/(h[i]*6)+
(points[i].y-M[i]*h[i]*h[i]/6)*sa/h[i]+
(points[i+1].y-M[i+1]*h[i]*h[i]/6)*sb/h[i];
S(%f)=%f\n"
x,s);
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 计算方法 实验一 插值 实验