计算方法上机作业.docx
- 文档编号:10048056
- 上传时间:2023-02-08
- 格式:DOCX
- 页数:23
- 大小:192.56KB
计算方法上机作业.docx
《计算方法上机作业.docx》由会员分享,可在线阅读,更多相关《计算方法上机作业.docx(23页珍藏版)》请在冰豆网上搜索。
计算方法上机作业
计算方法与实习上机作业
班级:
软件s161
学号:
165840
姓名:
肖盈圣
P211
例1对
计算不定积分
MATLAB程序如下:
y0=log(6.0/5.0);
fprintf('y[%d]=%f\n',0,y0)
n=1;
while
(1)
y1=1.0/n-5*y0;
fprintf('y[%d]=%f\n',n,y1)
if(n>=40)break;
end
y0=y1;
n=n+1;
end
输出结果如下:
P216
例1求方程
在1.5附近的根。
C语言程序如下:
#include
#include
#defineeps5e-6
#definedelta1e-6
floatBisection(floata,floatb,float(*f)(float))
{
floatc,fc,fa=(*f)(a),fb=(*f)(b);
intn=1;
printf("二分次数\t\tc\t\tf(c)\n");
while
(1)
{
if(fa*fb>0){printf("不能用二分法求解");break;}
c=(a+b)/2,fc=(*f)(c);
if(fabs(fc) elseif(fa*fc<0){b=c;fb=fc;} else{a=c;fa=fc;} if(b-a printf("%d\t\t%f\t\t%f\n",n++,c,fc); } returnc; } floatf(floatx) { returnx*x*x+x*x-3*x-3; } voidmain() { floata=1,b=2; floatx; x=Bisection(a,b,f); printf("\n方程的根是: %f",x); } 输出结果如下: P219 例2求方程 在1.5附近的根。 C语言程序如下: #include #include #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) { returnx*x*x+x*x-3*x-3; } floatf1(floatx) { return3.0*x*x+2*x-3; } voidmain() { floatx0,y0; printf("请输入迭代初值x0\n"); scanf("%f",&x0); printf("x(0)=%f\n",x0); y0=Newton(f,f1,x0); printf("方程的根为: %f\n",y0); } 输出结果如下: P224 例1解方程组 C语言程序如下: #include #include voidmain() { voidColPivot(float*,int,float[]); inti; floatx[3]; floatc[3][4]={0.101,2.304,3.555,1.183,-1.347,3.712,4.623,2.137,-2.835,1.072,5.643,3.035}; ColPivot(c[0],3,x); for(i=0;i<=2;i++)printf("x[%d]=%f\n",i,x[i]); } voidColPivot(float*c,intn,floatx[]) { inti,j,t,k; 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 { 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)); } } 输出结果如下: P229 例3求解方程组 , C语言程序如下: #include voidmain() { floatx[4]; inti; floata[4][5]={1,2,-12,8,27,5,4,7,-2,4,-3,7,9,5,11,6,-12,-8,3,49}; voidDirectLU(float*,int,float[]); DirectLU(a[0],4,x); for(i=0;i<=3;i++)printf("x[%d]=%f\n",i,x[i]); } voidDirectLU(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)); } } 输出结果如下: P231 例4用雅可比迭代法解方程组 C语言程序如下: #include #include #defineeps1e-6 #definemax100 voidJacobi(float*a,intn,floatx[]) { inti,j,k=0; floatepsilon,s; floaty[3]; for(i=0;i while (1) { epsilon=0; k++; for(i=0;i { s=0; for(j=0;j { if(j==i)continue; s+=*(a+i*(n+1)+j)*x[j]; } y[i]=(*(a+i*(n+1)+n)-s)/(*(a+i*(n+1)+i)); epsilon+=fabs(y[i]-x[i]); } for(i=0;i if(epsilon {printf("迭代次数为: %d\n",k);return;} if(k>=max) {printf("迭代发散");return;} } free(y); } voidmain() { inti; floata[3][4]={5,2,1,8,2,8,-3,21,1,-3,-6,1}; floatx[3]; Jacobi(a[0],3,x); for(i=0;i<3;i++)printf("x[%d]=%f\n",i,x[i]); } 输出结果如下: P233 例5用高斯-赛德尔迭代解方程组 C语言程序如下: #include #include #defineN500 voidmain() { inti; floatx[3]; floatc[3][4]={8,-3,2,20,4,11,-1,33,6,3,12,36}; voidGaussSeidel(float*,int,float[]); GaussSeidel(c[0],3,x); for(i=0;i<=2;i++)printf("x[%d]=%f\n",i,x[i]); } voidGaussSeidel(float*a,intn,floatx[]) { inti,j,k=1; floatd,dx,eps; for(i=0;i<=n-1;i++)x[i]=0.0; while (1) { eps=0; for(i=0;i<=n-1;i++) { d=0; for(j=0;j<=n-1;j++) { if(j==i)continue; d+=*(a+i*(n+1)+j)*x[j]; } dx=(*(a+i*(n+1)+n)-d)/(*(a+i*(n+1)+i)); eps+=fabs(dx-x[i]); x[i]=dx; } if(eps<1e-6){printf("迭代次数: %d\n",k);return;} if(k>N) { printf("迭代发散\n"); return; } k++; } } 输出结果如下: P239 例1已知函数表 0.56160 0.56280 0.56401 0.56521 0.82741 0.82659 0.82577 0.82495 用三次拉格朗日插值多项式求x=0.5635时的函数近似值。 C语言程序如下: #include floatLagrange(floatx[],floaty[],floatxx,intn) { inti,j; floata[4],yy=0; for(i=0;i<=n-1;i++) { a[i]=y[i]; for(j=0;j<=n-1;j++) if(j! =i)a[i]*=(xx-x[j])/(x[i]-x[j]); yy+=a[i]; } returnyy; } voidmain() { floatx[4]={0.56160f,0.56280f,0.56401f,0.56521f}; floaty[4]={0.82741f,0.82659f,0.82577f,0.82495f}; floatxx=0.5635f,yy; yy=Lagrange(x,y,xx,4); printf("x=%f,y=%f\n",xx,yy); } 输出结果如下: P241 例2已知函数表 0.4 0.55 0.65 0.8 0.9 0.41075 0.57815 0.69675 0.88811 1.02652 用牛顿插值多项式求 和 。 C语言程序如下: #include #defineN4 voidDifference(floatx[],floaty[],intn) { floatf[5]; 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]; } return; } voidmain() { inti; floatb,varx=0.895f; floatx[N+1]={0.4f,0.55f,0.65f,0.8f,0.9f}; floaty[N+1]={0.41075f,0.57815f,0.69675f,0.88811f,1.02652f}; 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); } 输出结果如下: P249 例1用复化辛卜生法计算积分 。 C语言程序如下: #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); s2+=(*f)(a+k*h); } s=h/6*((*f)(a)+4*s1+2*s2+(*f)(b)); returns; } floatf(floatx) { return1/(1+x*x); } voidmain() { inti,n=2; floats; for(i=0;i<=2;i++) { s=Simpson(f,0,1,n); printf("s(%d)=%f\n",n,s); n*=2; } } 输出结果如下: P252 例3用变步长梯形法计算积分 。 C语言程序如下: #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)>1e-6); returnt2; } floatf(floatx) { return1/(1+x*x); } voidmain() { floats; s=AutoTrap(f,0,1); printf("T(%d)=%f\n",n,s); } 输出结果如下: P254 例5用龙贝格方法计算积分 。 C语言程序如下: #include #include floatf(floatx) { return1/(1+x*x); } floatRomberg(floata,floatb,float(*f)(float),floateps) { intn=1,k; floath=b-a,x,temp; floatT1,T2,S1,S2,C1,C2,R1,R2; T1=(b-a)/2*((*f)(a)+(*f)(b)); while (1) { temp=0; for(k=0;k<=n-1;k++) { x=a+k*h+h/2; temp+=(*f)(x); } T2=(T1+temp*h)/2; if(fabs(T2-T1) S2=T2+(T2-T1)/3; if(n==1){T1=T2;S1=S2;h/=2;n*=2;continue;} C2=S2+(S2-S1)/15; if(n==2){C1=C2;T1=T2;S1=S2;h/=2;n*=2;continue;} R2=C2+(C2-C1)/63; if(n==4){R1=R2;C1=C2;T1=T2;S1=S2;h/=2;n*=2;continue;} if(fabs(R2-R1) R1=R2;C1=C2;T1=T2;S1=S2;h/=2;n*=2; } } voidmain() { floateps=5e-6; printf("R=%f",Romberg(0,1,f,eps)); } 输出结果如下:
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 计算方法 上机 作业