计算方法与实习报告.docx
- 文档编号:9622108
- 上传时间:2023-02-05
- 格式:DOCX
- 页数:21
- 大小:63.78KB
计算方法与实习报告.docx
《计算方法与实习报告.docx》由会员分享,可在线阅读,更多相关《计算方法与实习报告.docx(21页珍藏版)》请在冰豆网上搜索。
计算方法与实习报告
实习题1
1.用2种不同的顺序计算
,分析其误差的变化。
思路:
采用迭代法用C++来编程实现,一种用正序即n从1到10000取值,另一种用逆序
(一)正序:
采用递增的算法,将n从1到10000进行
的累加求和。
1)程序代码为:
#include
#include
usingnamespacestd;
voidmain()
{
floatsum=0;
intn=1;
while(n<=10000){
sum=sum+(float)1/(n*n);
n++;
}
cout<<"Thesumthatfrommintomaxis:
"< } (2)运行输出的结果为: Thesumthatfrommaxtominis: =1.64473 (二)逆序: 采用递减的算法,将n从10000到1进行 的累加求和。 (1)程序代码为: #include #include usingnamespacestd; voidmain() { floatsum=0; intn=10000; while(n>0){ sum=sum+(double)1/(n*n); n--; } cout<<"Thesumthatfrommaxtominis: "< } (2)运行输出的结果为: Thesumthatfrommintomaxis: 1.64483 结果分析: 由运行的结果可以看出正序的结果为1.64473,逆序的结果为1.64483,而题目所给出的精确值为1.644834。 可以发现从小到大求和比从大到小求和更接近实际值,因为在从大到小的求解过程中,会出现大数“吃掉”小数的现象,从而造成较大误差。 实习题2 1: 用牛顿法求下列方程的根: 思路: 给定初始值x0, 为根的容许误差, 为|f(x)|的容许误差,N为迭代次数的容许值。 在满足根的容许误差,函数值的容许误差及迭代次数条件下,计算x1=x0-f(x0)/f'(x0)。 以此进行下次循环。 (1)程序代码为: #include #include usingnamespacestd; #defineN100 #defineeps1e-6 #defineeta1e-8 floatNewton(float(*f)(float),float(*f1)(float),floatx0) { floatx1,d; intk=0; do { x1=x0-(*f)(x0)/(*f1)(x0); if(k++>N||fabs((*f1)(x1)) { printf("\nNewton迭代发散"); break; } d=fabs(x1)<1? x1-x0: (x1-x0)/x1; x0=x1; printf("x(%d)=%f\t",k,x0); } while(fabs(d)>eps&&fabs((*f)(x1))>eta); returnx1; } floatf(floatx) { returnfloat(x*x-exp(x)); } floatf1(floatx) { returnfloat(2*x-exp(x)); } voidmain(){ inti=0; floatx=0; for(i=0;i<100;i++){ x=x-A(x)/B(x); cout<<"x="< i++; } } (2)运行过程: x=-1 x=-0.733044 x=-0.703808 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 x=-0.703467 (3)运行结果为: x=-0.703467 结果分析: 得出的结果为x=-0.703467,由运行过程可知在此值时结果已经稳定,且经检验结果精确到了小数点后第六位。 实习题3 1: 用列主元高斯消去法解方程组: (1)x1+x2+3x4=4 2x1+x2–x3+x4=1 3x1–x2–x3+3x4=–3 –x1+2x2+3x3–x4=4 思路: 将方程用增广矩阵[A|b]=(aij)n (n+1)表示,分别进行消元和回代过程 消元过程: 对k=1,2,.....,n-1 (1)选主元 (2)如果 =0,则矩阵A奇异,程序结束;否则执行(3) (3)如果 ,则交换第k行与第 行对应元素位置 (4)消元,对i=k+1,........,n计算 ,对j=k+1,..........,n+1计算 回代过程 (1)若 ,则矩阵A奇异,程序结束;否则执行 (2) (2) ;对i=n-1,......2,1计算 (1)程序代码为: #include #include #include usingnamespacestd; voidmain() { inti; floatx[4]; floatc[4][5]={1,1,0,3,4, 2,1,-1,1,1, 3,-1,-1,3,-3, -1,2,3,-1,4}; voidDirect(float*,int,float[]); Direct(c[0],4,x); for(i=0;i<=3;i++)printf("x[%d]=%f\n",i,x[i]); } voidDirect(float*c,intn,floatx[]) { inti,j,k,t; floatp; for(i=0;i<=n-2;i++) { k=i; for(j=i+1;j<=n-1;j++) if(fabs(*(c+j*(n+1)+i))>(fabs(*(c+k*(n+1)+i))))k=j; if(k! =i) for(j=i;j<=n;j++) { p=*(c+i*(n+1)+j); *(c+i*(n+1)+j)=*(c+k*(n+1)+j); *(c+k*(n+1)+j)=p; } for(j=i+1;j<=n-1;j++) { p=(*(c+j*(n+1)+i))/(*(c+i*(n+1)+i)); for(t=i;t<=n;t++)*(c+j*(n+1)+t)-=p*(*(c+i*(n+1)+t)); } } for(i=n-1;i>=0;i--) { for(j=n-1;j>=i+1;j--) (*(c+i*(n+1)+n)-=x[j]*(*(c+i*(n+1)+j))); x[i]=*(c+i*(n+1)+n)/(*(c+i*(n+1)+i)); } } (2)运行结果为: x[0]=-1.333333 x[1]=2.333333 x[2]=-0.333333 x[3]=1.000000 4: 编写用追赶法解三对角线性方程组的程序,并解下列方程组。 (2)Ax=b其中 思路: 应用高斯消去法求解时,消元过程可以简化,得到同解的三角方程组,进行回代 消去过程为: 回代过程为 (1)程序代码为: #include voidmain() { floatx[10]; inti; floata[10][11]={-4,1,0,0,0,0,0,0,0,0,-27, 1,-4,1,0,0,0,0,0,0,0,-15, 0,1,-4,1,0,0,0,0,0,0,-15, 0,0,1,-4,1,0,0,0,0,0,-15, 0,0,0,1,-4,1,0,0,0,0,-15, 0,0,0,0,1,-4,1,0,0,0,-15, 0,0,0,0,0,1,-4,1,0,0,-15, 0,0,0,0,0,0,1,-4,1,0,-15, 0,0,0,0,0,0,0,1,-4,1,-15, 0,0,0,0,0,0,0,0,1,-4,-15}; voidDirect(float*,int,float[]); Direct(a[0],10,x); for(i=0;i<=9;i++)printf("x[%d]=%f\n",i,x[i]); } voidDirect(float*u,intn,floatx[]) { inti,r,k; for(r=0;r<=n-1;r++) { for(i=r;i<=n;i++) for(k=0;k<=r-1;k++) *(u+r*(n+1)+i)-=*(u+r*(n+1)+k)*(*(u+k*(n+1)+i)); for(i=r+1;i<=n-1;i++) { for(k=0;k<=r-1;k++) *(u+i*(n+1)+r)-=*(u+i*(n+1)+k)*(*(u+k*(n+1)+r)); *(u+i*(n+1)+r)/=*(u+r*(n+1)+r); } } for(i=n-1;i>=0;i--) { for(r=n-1;r>=i+1;r--) *(u+i*(n+1)+n)-=*(u+i*(n+1)+r)*x[r]; x[i]=*(u+i*(n+1)+n)/(*(u+i*(n+1)+i)); } } (2)运行结果为: x[0]=8.705758 x[1]=7.823033 x[2]=7.586371 x[3]=7.522452 x[4]=7.503439 x[5]=7.491305 x[6]=7.461785 x[7]=7.355835 x[8]=6.961556 x[9]=5.490389 实习题4 2: 按下列数据: xi0.300.420.500.580.660.72 yi1.044031.084621.118031.156031.198171.23223 作五次插值,并求x1=0.46,x2=0.55,x3=0.60时的函数近似值。 思路: (1)输入 令 (2)对i=0,1,2,....,n计算 (1)程序代码为: #include #include usingnamespacestd; #defineN5 voidDifference(floatx[],floaty[],intn) { float*f=newfloat[n+1]; intk,i; for(k=1;k<=n;k++) { f[0]=y[k]; for(i=0;i f[i+1]=(f[i]-y[i])/(x[k]-x[i]); y[k]=f[k]; } deletef; return; } voidmain() { inti; floatb,varx; floatx[N+1]={0.30,0.42,0.50,0.58,0.66,0.72}; floaty[N+1]={1.04403,1.08462,1.11803,1.15603,1.19817,1.23223}; Difference(x,y,N); b=y[N]; for(i=N-1;i>=0;i--) b=b*(varx-x[i])+y[i]; printf("Nn(%f)=%f",varx,b); } 在main函数中: 分别另varx为0.46,0.55,0.60可得: (2)运行结果为: =0.46时, =1.100724; =0.55时, =1.141271; =0.60时, =1.166194; 实习题5 1: 试分别用抛物线 和指数曲线 拟合下列数据: xi11.522.533.544.5 yi33.479.50122.65159.05189.15214.15238.65252.50 xi55.566.577.58 yi267.55280.50296.65301.40310.40318.15325.15 比较2个拟合函数的优劣。 (1)用抛物线 拟合: a.程序代码为: matlab代码为: >>A=[11.522.533.544.555.566.577.58]; B=[33.479.5122.65159.05189.15214.15238.65252.5267.55280.5296.65301.4310.4318.15325.15]; S=zeros(1,5);T=zeros(3,1); >>fork=1: 5 S(k)=sum(A.^(k-1)); end >>fork=1: 3 T(k)=sum(A.^(k-1).*B); end >>C=zeros(3,3); >>a=zeros(3,1); >>fori=1: 3 forj=1: 3 C(i,j)=S(i+j-1); end end >>a=C\T; >>fork=1: 3 disp(a(k)); end b.运行结果为: a[0]=-45.3333 a[1]=94.2302 a[2]=-6.1316 c.拟合曲线: (2)用指数曲线 拟合 将等式两边取自然对数得: lny=lna+bx 记y=lnya=lnay=a+bx 得到新数据: x11.522.533.544.555.566.577.58 y3.514.384.815.075.245.375.475.535.595.645.695.715.745.765.78 a.程序代码为: #include #include #include voidmain() { voidColPivot(float*,int,float[]); inti; floatx[2]; floatc[2][3]={15,67.5,79.29, 67.5,373.5,373.5}; ColPivot(c[0],2,x); for(i=0;i<=1;i++)printf("x[%d]=%f\n",i,x[i]); } voidColPivot(float*c,intn,floatx[]) { inti,j,k,t; floatp; for(i=0;i<=n-2;i++) { k=i; for(j=i+1;j<=n-1;j++) if(fabs(*(c+j*(n+1)+i))>(fabs(*(c+k*(n+1)+i))))k=j; if(k! =i) for(j=i;j<=n;j++) { p=*(c+i*(n+1)+j); *(c+i*(n+1)+j)=*(c+k*(n+1)+j); *(c+k*(n+1)+j)=p; } for(j=i+1;j<=n-1;j++) { p=(*(c+j*(n+1)+i))/(*(c+i*(n+1)+i)); for(t=i;t<=n;t++)*(c+j*(n+1)+t)-=p*(*(c+i*(n+1)+t)); } } for(i=n-1;i>=0;i--) { for(j=n-1;j>=i+1;j--) (*(c+i*(n+1)+n)-=x[j]*(*(c+i*(n+1)+j))); x[i]=*(c+i*(n+1)+n)/(*(c+i*(n+1)+i)); } } b.运行结果为: x[0]=4.208903 x[1]=0.239355 c.拟合曲线为: 分析: 带入数据进行检验,我们可以发现抛物线 与指数曲线 拟合相比,抛物线更加精确,而指数曲线具有较大的误差,因为在指数曲线拟合时,将指数函数按照一次函数来进行处理的,而拟合时,拟合函数的复杂度越高越准确,所以抛物线的拟合方法更准确。 实习题6 用复化梯形公式和复化Simpson公式计算积分 和 观察n为多少时所得近似值具有6位有效数字。 (1)用复化梯形公式计算积分 a.程序代码为 #include #include intn; floatAutoTrap(float(*f)(float),floata,floatb) { inti; floatx,s,h=b-a; floatt1,t2=h/2*((*f)(a)+(*f)(b)); n=1; do { s=0; t1=t2; for(i=0;i<=n-1;i++) { x=a+i*h+h/2; s+=(*f)(x); } t2=(t1+s*h)/2; n*=2; h/=2; } while(fabs(t2-t1)>0.5e-6); returnt2; } floatf(floatx) { returnsqrt(1+cos(x)*cos(x)); } voidmain() { floats; doublepi=3.141592653; s=AutoTrap(f,0,pi/2); printf("T(%d)=%f\n",n,s); } b.运行结果为: T[8]=1.910099 结果分析: 运行结果表示n为8时所得近似值具有6位有效数字。 (2)复化Simpson公式计算积分 a.程序代码为 #include #include floatSimpson(float(*f)(float),floata,floatb,intn) { intk; floats,s1,s2=0; floath=(b-a)/n; s1=(*f)(a+h/2); for(k=1;k<=n-1;k++) { s1+=(*f)(a+k*h+h/2); for(k=1;k<=n-1;k++) { s1+=(*f)(a+k*h+h/2); s2+=(*f)(a+k*h); } s=h/6*((*f)(a)+4*s1+2*s2+(*f)(b)); returns; } } floatf(floatx) { if(x==0)return1; elsereturnsin(x)/x; } voidmain() { doublepi=3.141592653; inti,n=2; floats; for(i=0;i<=2;i++) { s=Simpson(f,0,pi/4,n); printf("s(%d)=%f\n",n,s); n*=2; } } b.运行结果为: s (2)=1.005897 s(4)=0.887991 s(8)=0.824189 计算方法上机实验心得 通过用计算方法这门课的上机实验学习我对于c++和matlab程序的拥有一定掌握,并能完成一些在课本上的题目,如牛顿迭代法,列主元高斯消去法,雅可比迭代法,高斯赛德尔迭代法,插值法,复化梯形公式,复化Simpson公式等等。 在求解题目的过程中发现自己的编程能力逐步得到了提高,并且能够发现并解决一些简单的程序问题,同时加深了我们对课本上的知识的理解和认识,可以说是使得我们能够学以致用,加强了理论与实际的联系,我们从中获益良多。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 计算方法 实习 报告