计算方法实验报告.docx
- 文档编号:25421454
- 上传时间:2023-06-08
- 格式:DOCX
- 页数:19
- 大小:357.70KB
计算方法实验报告.docx
《计算方法实验报告.docx》由会员分享,可在线阅读,更多相关《计算方法实验报告.docx(19页珍藏版)》请在冰豆网上搜索。
计算方法实验报告
计算方法实验报告
姓名:
耿乐
学号:
913000710013
指导老师:
李建良
实验日期:
2014.11
实验一数值积分的数值实验
一、摘要:
本文对函数
,x
[-1,1],取n+1个基点,
=0,1,2,…,n,
步长
求n次插值多项式
计算最大相对误差,并将
和
画在坐标系中。
二、算法:
根据节点的特征,本文采用了Newton插值法来生成插值函数。
程序算法如下:
最后,将原函数和插值多项式图像画在同一坐标系中进行比较,并计算相应n值时的最大相对误差。
三、源程序
#include
#include
#definen4
intmain()
{
inti,j,t,k;
doubleco[n+1]={0},x[n+1],f[n+1],c[n+1][n+1]={0};
doubletable[n+1][n+1]={0},temp[n+1]={0};
doubleh=2.0/n,max,a,te;
doublefun(doublex);
doublecha(double*co,doublex);
for(i=0;i<=n;i++)
{
x[i]=-1.0+i*h;
table[i][0]=fun(x[i]);
}
for(j=1;j<=n;j++)
{
for(i=0;i<=n-j;i++)
{
table[i][j]=(table[i+1][j-1]-table[i][j-1])/(x[i+j]-x[i]);
}
f[j]=table[0][j];
}//生成差商表
printf("差商表为:
\n");
for(i=0;i<=n;i++)
{
printf("%8.4f\t",x[i]);
for(j=0;j<=n;j++)
printf("%8.4f",table[i][j]);
printf("\n");
}
printf("\n");
c[0][0]=1.0;
for(k=1;k<=n;k++)
{
c[k][0]=0;
for(i=0;i<=k;i++)
{
c[k][i+1]=c[k-1][i];
}
for(j=0;j { c[k][j]=c[k][j]-c[k-1][j]*x[k-1]; } }//求(x-x0)(x-x1)···(x-xi)展开系数 /*for(i=0;i<=n;i++) { for(j=0;j<=n;j++) { printf("%8.4f",c[i][j]); } printf("\n"); }*/ for(i=0;i<=n;i++) { f[i]=table[0][i]; }//取出差商表第一行 printf("\n"); for(i=0;i<=n;i++) { for(j=0;j<=n;j++) { co[i]+=f[j]*c[j][i]; } }//求各次项系数 printf("插值函数为: \nPn(x)="); for(i=n;i>=0;i--) { if(co[i]>=-0.00000001&&co[i]<=0.00000001) continue; printf("+(%8.5fx^%d)",co[i],i); } printf("\n"); printf("最大相对误差为: \n"); max=0; for(a=-1.0;a<=1.0;a+=0.01) { te=1.0-cha(co,a)/fun(a); If(te<0) te=-te; if(te>max) max=te; }//求最大相对误差 printf("%f\n",max); return0; } doublefun(doublex)//原函数 { return(1.0/(1.0+25.0*x*x)); } doublecha(double*co,doublex)//插值函数 { doublet=0; inti; for(i=n;i>=0;i--) { t+=co[i]*pow(x,i); } returnt; } 4、实验结果及分析 1、n=2 最大相对误差为6.0072 2、n=4 最大相对误差为7.4483 3、n=6 最大相对误差为12.9112 4、n=8 最大相对误差为21.4013 5、n=10 最大相对误差为32.5478 5、实验结论 可以看出,插值次数较低时,插值函数与原函数差异明显,相对误差较大;而插值次数增大时,仅在零点附近较小的区域,插值函数与原函数较为接近,而其他部分,尤其是靠近区域两侧边缘的区间上,误差极大。 这一现象称为龙格现象。 因此,可采用分段Hermite插值法进行改进。 实验二线性方程组求解的数值实验 一、摘要 本文对从2到9的每一个n值,解n阶方程组 在这里A和b如下定义: 其中 ,最后解释其中现象。 二、算法 由于本案例的特殊性——系数矩阵为对称矩阵,可以采用平方根法来求解。 同时由于系数矩阵中元素的特殊分布,在采用高斯消去法时即是在采用列主元消去法。 两种算法如下: 1)高斯消去法 2)平方根法 三、源程序 1)高斯消去法 #include #defineN9 intmain() { doublea[N+1][N+1],b[N+1],t,x[N+1]; inti,j,k; doublep(intx); for(i=1;i<=N;i++) { for(j=1;j<=N;j++) { a[i][j]=1.0/(i+j-1); } b[i]=p(N+i-1)-p(i-1); }//生成增广矩阵 printf("增广矩阵\n"); for(i=1;i<=N;i++) { for(j=1;j<=N;j++) { printf("%-7.4f",a[i][j]); } printf("\t%-7.4f\n",b[i]); } printf("\n"); for(k=1;k<=N;k++) { t=a[k][k]; for(j=k;j<=N;j++) { a[k][j]=a[k][j]/t; } b[k]=b[k]/t;//每次化三角前把首个元素化为1 for(i=k+1;i<=N;i++) { t=a[i][k]/a[k][k]; for(j=k;j<=N;j++) { a[i][j]=a[i][j]-t*a[k][j]; } b[i]=b[i]-t*b[k]; }//化下三角 } printf("化简后的增广矩阵\n"); for(i=1;i<=N;i++) { for(j=1;j<=N;j++) { printf("%-6.2f",a[i][j]); } printf("\t%-8.4f\n",b[i]); } printf("\n"); for(i=N;i>=1;i--) { t=0; for(j=i+1;j<=N;j++) { t+=a[i][j]*x[j]; } x[i]=(b[i]-t)/a[i][i]; }//倒序求解 printf("方程组解\n"); for(i=1;i<=N;i++) { printf("x%d=%25f\n",i,x[i]); } return0; } doublep(intx) { doublet; t=14.0+N*(12+3*N); t=-7+N*N*t; t=x*x*t+2.0; t=(x*x*t)/24.0; returnt; } 2)平方根法 #include #defineN9 intmain() { inti,j,s; doubletemp; doublep(intx); doublea[N+1][N+1],b[N+1],x[N+1],y[N+1]; doublelow[N+1][N+1]={0},d[N+1]; printf("增广矩阵\n"); for(i=1;i<=N;i++) { for(j=1;j<=N;j++) { a[i][j]=1.0/(i+j-1); } b[i]=p(N+i-1)-p(i-1); }//生成增广矩阵 for(i=1;i<=N;i++) { for(j=1;j<=N;j++) { printf("%-6.3f",a[i][j]); } printf("\t"); printf("%-6.3f\n",b[i]); } for(j=1;j<=N;j++) { low[j][j]=1; temp=0; for(s=1;s { temp+=low[j][s]*low[j][s]*d[s]; } d[j]=a[j][j]-temp; for(i=j+1;i<=N;i++) { temp=0; for(s=1;s { temp+=low[i][s]*low[j][s]*d[s]; } low[i][j]=(a[i][j]-temp)/d[j]; } } y[1]=b[1]; for(i=2;i<=N;i++) { temp=0; for(s=0;s {temp+=low[i][s]*y[s];} y[i]=b[i]-temp; } x[N]=y[N]/d[N]; for(i=N-1;i>=1;i--) { temp=0; for(s=i+1;s<=N;s++) {temp+=low[s][i]*x[s];} x[i]=y[i]/d[i]-temp; } printf("\n方程组解为\n"); for(i=1;i<=N;i++) {printf("x%d=%25f\n",i,x[i]);} return0; } doublep(intx) { doublet; t=14.0+N*(12+3*N); t=-7+N*N*t; t=x*x*t+2.0; t=(x*x*t)/24.0; returnt; } 四、计算结果 以下给出两种方法下n为2,5,9的计算结果: 1)n=2 2)n=5 3)n=9 6、实验结论 由计算可以看出,随着n的增大,两种计算结果差别变得非常大,用其他软件计算解时会发现以上解的误差会比较大,这是因为随着n增大,系数矩阵的条件数变大,微小的扰动就容易造成很大的误差,可以认为该方程组是病态的最后得不到精确解。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 计算方法 实验 报告