[理学]数值方法第2版答案.doc

C语言编程习题第二章习题2-25. 用二分法编程求 6x4 -40x2+9=0 的所有实根。#include <stdio.h>#include <math.h>#define N 10000double A,B,C;double f(double x) return (A*x*x*x*x+B*x*x+C);void BM(double a,double b,double eps1,double eps2) int k; double x,xe;double valuea = f(a);double valueb = f(b);if (valuea > 0 && valueb > 0 | valuea <0 && valueb < 0) return;printf("Finding root in the range: %.3lf, %.3lfn", a, b); for(k=1;k<=N;k+) x=(a+b)/2; xe=(b-a)/2; if(fabs(xe)<eps2 | fabs(f(x)<eps1) printf("The x value is:%gn",x); printf("f(x)=%gnn",f(x); return; if(f(a)*f(x)<0) b=x; else a=x; printf("No convergence!n");int main() double a,b,eps1,eps2,step,start; printf("Please input A,B,C:n"); scanf("%lf %lf %lf",&A,&B,&C); printf("Please input a,b, step, eps1,eps2:n"); scanf("%lf %lf %lf %lf %lf",&a,&b,&step,&eps1,&eps2); for (start=a; (start+step) <= b; start += step) double left = start;double right = start + step;BM(left, right, eps1, eps2); return 0;运行：Please input A,B,C:6 -40 9Please input a,b, step, eps1,eps2:-10 10 1 1e-5 1e-5Finding root in the range: -3.000, -2.000The x value is:-2.53643f(x)=-0.00124902Finding root in the range: -1.000, 0.000The x value is:-0.482857f(x)=0.00012967Finding root in the range: 0.000, 1.000The x value is:0.482857f(x)=0.00012967Finding root in the range: 2.000, 3.000The x value is:2.53643f(x)=-0.00124902有时若把判别语句if(fabs(xe)<eps2 | fabs(f(x)<eps1)改为if(fabs(xe)<eps2 && fabs(f(x)<eps1)会提高精度，对同一题运行结果：Finding root in the range: -3.000, -2.000The x value is:-2.53644f(x)=-4.26496e-007Finding root in the range: -1.000, 0.000The x value is:-0.482861f(x)=-7.3797e-006Finding root in the range: 0.000, 1.000The x value is:0.482861f(x)=-7.3797e-006Finding root in the range: 2.000, 3.000The x value is:2.53644f(x)=-4.26496e-007习题2-35. 请用埃特金方法编程求出x=tgx在4.5（弧度）附近的根。#include <stdio.h>#include <math.h>#define N 100#define PI 3.1415926void SM(double x0,double eps) int k; double x; double x1,x2; for(k=1;k<=N;k+) x1=sin(x0)/cos(x0); x2=sin(x1)/cos(x1); x=x0-(x1-x0)*(x1-x0)/(x2-2*x1+x0); if(fabs(x-x0)<eps) printf("The x value is:%gn",x); return; x0=x; printf("No convegence!n");int main() double eps,x0; printf("Please input eps,x0:n"); scanf("%lf %lf",&eps,&x0); SM(x0,eps); return 0;运行:Please input eps,x0:1e-5 4.5The x value is:4.49341习题2-411请编出用牛顿法求复根的程序，并求出 P(z)=z4-3z3+20z2+44z+54=0 接近于z0=2.5+4.5i 的零点。#include <cstdio>#include <cmath>#define MAX_TIMES 1000typedef struct double real, image; COMPLEX;COMPLEX Aa5=54,0,44,0,20,0,-3,0,1,0;COMPLEX Bb4=44,0,40,0,-9,0,4,0;COMPLEX zero = 0, 0;double eps1=1e-6;double eps2=1e-6;COMPLEX multi(COMPLEX a,COMPLEX b) COMPLEX result; result.real = a.real * b.real - a.image * b.image; result.image = a.image * b.real + a.real * b.image; return result;COMPLEX Div(COMPLEX a,COMPLEX b) COMPLEX z3; double s; s=(b.real*b.real)+(b.image*b.image); z3.real=b.real; z3.image=-b.image; z3=multi(a,z3); z3.real=z3.real/s; z3.image=z3.image/s; return z3;COMPLEX add(COMPLEX a,COMPLEX b) COMPLEX result; result.real = a.real + b.real; result.image = a.image + b.image; return result;COMPLEX subtract(COMPLEX a, COMPLEX b) COMPLEX result; result.real = a.real - b.real; result.image = a.image - b.image; return result;COMPLEX times(COMPLEX z,int n) int i; COMPLEX result=1, 0; for (i=0;i<n;i+) result=multi(result,z); return result;double distance(COMPLEX a, COMPLEX b) return sqrt(a.real - b.real) * (a.real - b.real) + (a.image - b.image) * (a.image - b.image);double complex_abs(COMPLEX a) return sqrt(a.real * a.real + a.image * a.image);COMPLEX f(COMPLEX x) int i; COMPLEX result=zero; for (i=0;i<5;i+) result=add(result,multi(Aai,times(x,i); return result;COMPLEX ff(COMPLEX x) int i; COMPLEX result=zero; for (i=0;i<4;i+) result=add(result,multi(Bbi,times(x,i); return result;int main() COMPLEX z0,z1,result; double x,y; int k; printf("please input x,yn"); scanf("%lf %lf",&x,&y); z0.real=x; z0.image=y; for(k=0; k<MAX_TIMES; k+) z1 = subtract(z0, Div(f(z0), ff(z0); result = f(z1); if (distance(z0,z1)<eps1 | complex_abs(result)<eps2) printf("The root is: z=(%.3f + %.3fi), f(z)=(%e + %ei)n", z1.real,z1.image, result.real, result.image); return 0; z0 = z1; printf("No convergence!n"); return 0;运行:please input x,y2.5 4.5The root is: z=(2.471 + 4.641i), f(z)=(-1.122705e-007 + -1.910245e-007i)习题2-62. 请编程用劈因子法求高次方程 x4+ x 3+5 x 2+4 x +4=0的所有复根。#include<cstdio>#include<cmath>int main()float b20,c20;float delta;float eps;int print=0;float fru0,fru1,s0,s1,frv0,frv1;float deltau,deltav;float u,v;float a20;int n,j;float r0,r1;float r;int i, k;printf("Please input max exponent:n");scanf("%d",&n);printf("Please input epslion:n");scanf("%f",&eps);printf("Please input the coefficient:n");for(i=0;i<=n;i+) scanf("%f,",&ai);for(j=0;j<=n;j+) printf("+%.3fX%d",aj,n-j);printf("n"); for(k=1;k<=2;k+) u=an-1/an-2;v=an/an-2;if(k=2) u=0;v=4;while(1) b0=a0;b1=a1-u*b0;for(j=2;j<=n;j+) bj=aj-u*bj-1-v*bj-2;r0=bn-1;r1=bn+u*bn-1;c0=b0;c1=b1-u*b0;for(j=2;j<=n-2;j+)cj=bj-u*cj-1-v*cj-2;s0=cn-3;s1=cn-2+u*cn-3;fru0=u*s0-s1;fru1=v*s0;frv0=-s0;frv1=-s1;deltau=-(frv1*r0-frv0*r1)/(fru0*frv1-fru1*frv0);deltav=-(fru1*r0-fru0*r1)/(frv0*fru1-frv1*fru0);u=u+deltau;v+=deltav;delta=u*u-4*v;if(fabs(deltau)<eps)&&(fabs(deltav)<eps)print=1;printf("found roots:");r=-u/2;i=(int) sqrt(fabs(delta)/2;if(delta<0) printf("%.3f+%.3fit",r,fabs(double)i);printf("%.3f-%.3fin",r,fabs(double)i);else printf("%.3ftt",r+i);printf("%.3fn",r-i);if(n=1) printf("Found root :%.3fn",a1/a0);if(k=2) return 0;if(k=1&&print) break; /* end while */ /* end for */return 0;运行:Please input max exponent:4Please input epslion:1e-5Please input the coefficient:1 1 5 4 4+1.000X4+1.000X3+5.000X2+4.000X1+4.000X0found roots:-0.500+0.000i-0.500-0.000ifound roots:0.000+2.000i0.000-2.000i习题3-15. 请编程用列全主元高斯约当消去法求矩阵A, B的逆矩阵。 #include <cstdio>#include <cmath>#define MAX 10int ROW;static float AMAX+12*MAX+1;void putout()int i,j;printf("nThe reverse matrix is:n");for(i=1;i<=ROW;i+) for(j=ROW+1;j<=2*ROW;j+) printf("%.3f ",Aij);printf("n");int get()float maxMAX+1;int count_rMAX+1;float tmp2*MAX+1;float pass1,pass2;int i,j,k;for(i=1;i<=ROW;i+) maxi=Aii;count_ri=i;for(k=i;k<=ROW;k+) if(maxi<Aki) maxi=Aki;count_ri=k; /* find the max one*/if(maxi=0) return -1;if(count_ri!=i) for(k=1;k<=2*ROW;k+)tmpk=Acount_rik;Acount_rik=Aik;Aik=tmpk; /*change the rows*/pass1=Aii;for(k=1;k<=2*ROW;k+) Aik/=pass1;for(j=1;j<=ROW;j+) if(j=i) continue;else pass2=Aji;for(k=1;k<=2*ROW;k+) Ajk=Ajk-pass2 *Aik; /* get simple*/putout();return 0;int main()int i,j;float tmp=-1;float p;printf("please input the number of ROW:");scanf("%d",&ROW);printf("n");if(ROW>=MAX) printf("the number of row is too big!n");return 0;printf("Please input the matrix A:n");for(i=1;i<=ROW;i+) for(j=1;j<=ROW;j+) scanf("%f,",&p);Aij=p;if(tmp<fabs(Aij) tmp=fabs(Aij);if(tmp=0) printf("No inverse!");return 0;for(i=1;i<=ROW;i+) Aii+ROW=1;return get();运行：please input the number of ROW:3Please input the matrix A:-3 8 52 -7 41 9 -6The reverse matrix is:0.026 0.396 0.285 0.068 0.055 0.094 0.106 0.149 0.021please input the number of ROW:4Please input the matrix A:2 1 -3 -13 1 0 7-1 2 4 -21 0 -1 5The reverse matrix is:-0.047 0.588 -0.271 -0.9410.388 -0.353 0.482 0.765-0.224 0.294 -0.035 -0.471-0.035 -0.059 0.047 0.294习题3-24. 编程用追赶法解下列三对角线方程组。#include <stdio.h>#include <math.h>#define ROW 5static float AROW+1=0,0,0,0,0,0;static float BROW+1=0,0,0,0,0,0;static float CROW+1=0,0,0,0,0,0;static float XROW+1=0,0,0,0,0,0;int main()int i;printf("Please input the A:n");for(i=2;i<=ROW;i+) scanf("%f,",&Ai);printf("Please input the B:n");for(i=1;i<=ROW;i+) scanf("%f,",&Bi);printf("Please input the C:n");for(i=1;i<=ROW-1;i+) scanf("%f,",&Ci);printf("Please input the F:n");for(i=1;i<=ROW;i+) scanf("%f,",&Xi);C1=C1/B1;for(i=2;i<=ROW-1;i+)Ci=Ci/(Bi-Ai*Ci-1);X1=X1/B1;for (i=2; i<=ROW; i+)Xi=(Xi-Ai*Xi-1)/(Bi-Ai*Ci-1);for(i=ROW-1;i>=1;i-)Xi=Xi-Ci*Xi+1;printf ("The value is : n");for(i=1;i<=ROW;i+) printf ("x%d=%.3fn",i,Xi);return 0;运行：Please input the A:-1 -1 -1 -1Please input the B:4 4 4 4 4 Please input the C:-1 -1 -1 -1Please input the F:100 200 200 200 100The value is : x1=46.154x2=84.615x3=92.308x4=84.615x5=46.154习题4-31. 用高斯赛德尔方法编程解下列线性方程组，要求当时迭代终止。#include<cstdio>#include<cmath>float x6;float y6;float eps;float a66=4,-1,0,-1,0,0,-1,4,-1,0,-1,0,0,-1,4,0,0,-1,-1,0,0,4,-1,0,0,-1,0,-1,4,-1,0,0,-1,0,-1,4;float b6=0,5,0,6,-2,6;int gs()int i,k,j;float s;for(i=0;i<6;i+) yi=xi=0;for(k=0;k<20;k+) for(i=0;i<6;i+) s=0;for(j=0;j<6;j+) if(j!=i) s=s+aij*yj;yi=(bi-s)/aii;for(i=0;(i<6)&&(abs(yi-xi)<eps);i+) ;if(i>=5) break;for(j=0;j<6;j+) xj=yj;if(k>20) printf("No convergence./n");return -1;return 1;int main()int tag,i,m;printf("nthe max circles=");scanf("%d",&m);printf("eps=");scanf("%lf",&eps);tag=gs();if(tag>0) for(i=0;i<=5;i+) xi=yi;printf("x(%d)=%.3en",i+1,xi);return 0;运行：nthe max circles=10eps=1e-5x(1)=1.000e+000x(2)=2.000e+000x(3)=1.000e+000x(4)=2.000e+000x(5)=1.000e+000x(6)=2.000e+000习题4-4#include <cmath>#include <cstdio>int agsdl(double a10, double b,int n, double x,double eps,double w,int m) int i,j,pcount; double p,t,s,q; for(i=0;i<10;i+) xi=0; p=eps+1.0; for (pcount=0;p>=eps&&pcount<m;pcount+) p=0.0;for (i=0; i<n; i+) t=xi; s=0.0;for (j=0; j<n; j+)if (j!=i) s=s+aij*xj;xi=(bi-s)/aii;q=w*fabs(xi-t);if (q>p) p=q; return pcount;int main() int i,j,n,m,count=0;double eps,w3;static double a1010;static double x10,b10;printf("nthe max circles=");scanf("%d",&m);printf("n=");scanf("%d",&n);printf("nmatrix A:n"); for(i=0;i<n;i+)for(j=0;j<n;j+) scanf("%lf",&aij); printf("nvector b:n");for(i=0;i<n;i+) scanf("%lf",&bi);printf("eps=");scanf("%lf",&eps);for(i=0;i<3;i+) printf("w(%d)=",i); scanf("%lf",&wi);for(j=0;j<3;j+) printf("w(%d)=%.3fn",j,wj);if(agsdl(a,b,n,x,eps,wj,m)<m) for(i=0;i<n;i+) printf("x(%d)=%.3e ",i,xi);printf("n");elseprintf("failn");return 0;运行：the max circles=10n=3matrix A:4 -1 0-1 4 -10 -1 4vector b:1 4 -3eps=1e-5w(0)=1.03w(1)=1w(2)=1.1w(0)=1.030x(0)=5.000e-001 x(1)=1.000e+000 x(2)=-5.000e-001 w(1)=1.000x(0)=5.000e-001 x(1)=1.000e+000 x(2)=-5.000e-001 w(2)=1.100x(0)=5.000e-001 x(1)=1.000e+000 x(2)=-5.000e-001 第5章习题5-12. 编程求下列矛盾方程组的解:#include<fstream.h>#include<stdio.h>#include<stdlib.h>double eps=1.0e-5;/error rangedouble *a; /save matrix Adouble *ab;/save augmented matrix ABdouble *ta;/save matrix ATAdouble *b; /save vectors Bdouble *tb; /save ATBdouble *ca;/save matrix Adouble *cb; /save vectors Bdouble *x; /save result vectors Xint *fr;/row flagint *fc;/column flagint m,cm;/row num of matrix Aint n;/column num of matrix Aint r;void init();int G_J();int d_0(int);void change();void g_a_m();void show();void show1();void main()int t;init();t=1;while(1)r=G_J();if (m=r) break;if (d_0(r)=0) break;change();t=0;if (n=r)if (t)cout<<"该线性方程组是非矛盾线性方程组，有唯一解："<<endl;show();elsecout<<"该线性方程组是矛盾线性方程组，有唯一最小二乘解："<<endl;show();elseif (t)cout<<"该线性方程组是非矛盾线性方程组，有普遍解："<<endl;show1();elsecout<<"该线性方程组是矛盾线性方程组，有普遍最小二乘解："<<endl;show1();void init() /open the data file-ab.txtifstream fin("ab.txt");if(!fin)cout<<"Can't open file ab.txt to read!"<<endl;exit(-1);/read the row and column information of matrix Afin>>m>>n;cm=m;/require room of matrix Aa=(double *)malloc(m+1)*sizeof(double *);int i,j;for(i=1;i<=m;i+)ai=(double *)malloc(n+1)*sizeof(double);/read matrix A's data from ab.txtfor(i=1;i<=m;i+)for(j=1;j<=n;j+)fin>>aij;/init matrix Bb=(double*)malloc(m+1)*sizeof(double);for(i=1;i<=m;i+)fin>>bi;/close the data filefin.close();/backup matrix A and Bca=a; cb=b;return;/generate augmented matrixvoid g_