高等教育东北大学数值分析实验报告.docx

1、高等教育东北大学数值分析实验报告数值分析实验报告课题一 迭代格式的比较一、 问题提出设方程f(x)=x- 3x 1=0 有三个实根 x=1.8793 , x=-0.34727 ,x=-1.53209现采用下面三种不同计算格式，求 f(x)=0的根 x 或x 1、 x = 2、 x = 3、 x = 二、要求1、编制一个程序进行运算，最后打印出每种迭代格式的敛散情况；2、用事后误差估计来控制迭代次数，并且打印出迭代的次数；3、初始值的选取对迭代收敛有何影响；4、分析迭代收敛和发散的原因。三、目的和意义1、通过实验进一步了解方程求根的算法；2、认识选择计算格式的重要性；3、掌握迭代算法和精度控制；

2、4、明确迭代收敛性与初值选取的关系。程序代码：#include#include int k=1,k1=1,k2=1,k3=1; float x1,x2,x3; float one(float x0) x1=(3*x0+1)/(x0*x0); return(x1); float two(float x0) x2=(pow(x0,3)-1)/3; return(x2); float three(float x0) x3=pow(3*x0+1,0.33333); return(x3);main() float x,x0; printf(输入x0=); scanf(%f,&x); x0=x; x1=o

3、ne(x0); printf(第一个公式迭代结果: n); while(fabs(x0-x1)1e-5) printf(x1=%6.5fn,x1);x0=x1;x1=one(x0);k1+; printf(x1=%6.5f n,x1); printf(k1=%in,k1);x0=x;x2=two(x0);printf(第二个公式迭代结果: n);while(fabs(x0-x2)1e-5)printf(x2=%6.5fn,x2);x0=x2;x2=two(x0);k2+; printf(k2=%in,k2); x0=x; x3=three(x0); printf(第三个公式迭代结果: n);

4、while(fabs(x0-x3)1e-5) printf(x3=%6.5fn,x3); x0=x3; x3=three(x0); k3+; printf(x3=%6.5fn,x3); printf(k3=%in,k3);scanf(%);实验结果：四、程序运行结果讨论和分析：对于第一种迭代格式，收敛区间-8.2 -0.4，在该收敛区间内迭代收敛于-1.53209，只能求得方程的一个根；对于第二种迭代格式，收敛区间-1.5 1.8，在该收敛区间内迭代收敛于-0.34730，同样只能求得方程的一个根；对于第三种迭代格式，收敛区间-0.3 +)，在该收敛区间内迭代收敛于1.87937，只能求得方程

5、的一个根；由以上结果很容易发现，初值的选取对迭代敛散性有很大影响。以第一种迭代格式为例，当初值大于等于-0.3时，迭代格式发散；当初值小于等于-8.3时，迭代格式也发散；只有初值在-0.3和-8.3之间时，迭代格式才收敛于1.53209。其他迭代格式也有这样的性质，即收敛于某个数值区间，超出这个区间迭代格式就是发散的，这就是所谓迭代格式的收敛性。对于不同迭代格式在不同区间具有不同的敛散性的原因，我认为可以从一下两方面理解：1、迭代法是一种逐次逼近法，其基本思想是将隐式方程归结为一组显式的计算公式，就是说，迭代过程实质上是个逐步显式化的过程。2、我们可以用几何图像来更好地理解迭代过程。由图可知，

6、在某些区间选取的初始值随着迭代次数的增加会越来越逼近精确值，即收敛于精确值，而在另外一些区间选取的初始值随着迭代次数的增加却离精确值越来越远，即不会收敛于一个确定值。课题二 线性方程组的直接算法一、问题提出给出下列几个不同类型的线性方程组，请用适当算法计算其解。 1、设线性方程组=x= ( 1, -1, 0, 1, 2, 0, 3, 1, -1, 2 ) 2、设对称正定阵系数阵线方程组= x = ( 1, -1, 0, 2, 1, -1, 0, 2 ) 3、三对角形线性方程组 = x= ( 2, 1, -3, 0, 1, -2, 3, 0, 1, -1 ) 二、要求1、 对上述三个方程组分别利

7、用Gauss顺序消去法与Gauss列主元消去法；平方根法与改进平方根法；追赶法求解（选择其一）；2、 应用结构程序设计编出通用程序；3、 比较计算结果，分析数值解误差的原因；4、 尽可能利用相应模块输出系数矩阵的三角分解式。三、目的和意义1、通过该课题的实验，体会模块化结构程序设计方法的优点；2、运用所学的计算方法，解决各类线性方程组的直接算法；3、提高分析和解决问题的能力，做到学以致用；4、通过三对角形线性方程组的解法，体会稀疏线性方程组解法的特点。程序代码：1.Gsuss列主元消去法#include#include#includeusing namespace std;int main()

8、 int n, i,j,k,m; coutn; float *A=new double*(n+1); for(i=1;i=n;i+) Ai=new doublen+1; float *b=new doublen+1; float *x=new doublen+1; float l; float temp1,temp2,temp3; cout输入系数矩阵A:endl; for(i=1;i=n;i+) for(j=1;jAij; cout输入向量b:; for(i=1;ibi; coutendl; for(k=1;kn;k+) temp1=abs(Akk); m=k; for(i=k;i=n;i+

9、)/找最大值的列主元 if(temp1abs(Aik) temp1=abs(Aik);m=i;/m是确定的最列主元的行标 if(temp1=0) coutno unique solution!endl; exit(0); if(m!=k)/换行 for(j=1;j=n;j+) temp2=Akj; Akj=Amj; Amj=temp2; temp3=bk; bk=bm; bm=temp3; for(i=k+1;i=n;i+)/消元 l=Aik/Akk; for(j=k+1;j=n;j+) Aij=Aij-l*Akj; bi=bi-l*bk; if(Ann=0) coutno unique so

10、lution!=1;i-) float sum=0; for(j=i+1;j=10;j+) sum=sum+Aij*xj; xi=(bi-sum)/Aii; cout输出结果向量x:endl; for(i=1;i=10;i+) coutxiendl; system(pause); return 0;2.平方根法#include#include#includeusing namespace std;int main() int n,i,j,k,m; coutn; double *A=new double*(n+1); for(i=1;i=n;i+) Ai=new doublen+1; doubl

11、e *b=new doublen+1; double *x=new doublen+1; double *y=new doublen+1; cout输入系数对称正定矩阵A:endl; for(i=1;i=n;i+) for(j=1;jAij; cout输入向量b:; for(i=1;ibi; coutendl; for(k=1;k=n;k+) double sum=0; for(m=1;m=k-1;m+) sum=sum+pow(Akm,2.0); sum=Akk-sum; Akk=sqrt(sum); for(i=k+1;i=n;i+) double temp1=0; for(m=1;m=k

12、-1;m+) temp1=temp1+Aim*Akm; temp1=Aik-temp1; Aik=temp1/Akk; double temp2=0; for(m=1;m=1;k-) double temp3=0; for(m=k+1;m=n;m+) temp3=temp3+Amk*xm; xk=(yk-temp3)/Akk; cout输出结果向量x:endl; for(i=1;i=n;i+) coutxiendl; system(pause); return 0;3.追赶法#include#include#includeusing namespace std;int main() int n

13、,i; coutn; double *a=new doublen+1; double *c=new doublen+1; double *d=new doublen+1; double *b=new doublen+1; double *x=new doublen+1; double *y=new doublen+1; cout输入系数矩阵A数据：endl; for(i=1;iai; for(i=1;ici; for(i=1;idi; cout输入b :endl; for(i=1;ibi; for(i=1;i=n-1;i+) ci=ci/ai; ai+1=ai+1-di+1*ci; cout输

14、出解向量a:endl; for(i=1;i=n;i+)coutaiendl; cout输出解向量c:endl; for(i=1;i=n;i+)coutciendl; y1=b1/a1; for(i=2;i=n;i+) yi=(bi-di*yi-1)/ai; cout输出解向量y:endl; for(i=1;i=n;i+)coutyi=1;i-) xi=yi-ci*xi+1; cout输出解向量x:endl; for(i=1;i=n;i+)coutxiendl; system(pause); return 0;四、 程序运行结果分析 在方法的选择上存在一定的误差，可以选一些更准确的方法求解；程序

15、中对变量的类型设定，若设成double型，结果可以更精确；计算机在做运算时，会根据需要对中间结果进行舍入，这也会对最终结果有影响； 课题三 线性方程组的迭代法一、问题提出1、设线性方程组=x= ( 1, -1, 0, 1, 2, 0, 3, 1, -1, 2 ) J迭代算法：程序设计流程图：源程序代码：#include#include #include void main() float a5051,x150,x250,temp=0,fnum=0; int i,j,m,n,e,bk=0; printf(使用Jacobi迭代法求解方程组:n); printf(输入方程组的元：nn=); scan

16、f(%d,&n); for(i=1;in+1;i+) x1i=0; printf(输入方程组的系数矩阵:n); for(i=1;in+1;i+) j=1; while(jn+1) scanf(%f,&aij); j+; printf(输入方程组的常数项:n); for(i=1;in+1;i+) scanf(%f,&ain+1); printf(n); printf(请输入迭代次数:n); scanf(%d,&m); printf(请输入迭代精度:n); scanf(%d,&e); while(m!=0) for(i=1;in+1;i+) for(j=1;jn+1;j+) if (j!=i) t

17、emp=aij*x1j+temp; x2i=(ain+1-temp)/aii; temp=0; for(i=1;itemp) temp=fnum; if(temp=pow(10,-4) bk=1; for(i=1;in+1;i+) x1i=x2i; m-; printf(原方程组的解为:n); for(i=1;in+1;i+) if(x1i-x2i)=e|(x2i-x1i)=e) printf(x%d=%7.4f ,i,x1i); 运行结果：2、设对称正定阵系数阵线方程组= x = ( 1, -1, 0, 2, 1, -1, 0, 2 ) GS迭代算法：#include#include#inc

18、ludeconst int m=11;void main() int choice=1; while(choice=1) double amm,bm,e,xm,ym,w,se,max; int n,i,j,N,k; coutGauss-Seidol迭代法endl; coutn; for(i=1;i=n;i+) cout请输入第i个方程的各项系数：; for(j=1;jaij; cout请输入各个方程等号右边的常数项:n; for(i=1;ibi; coutN; coute; for(i=1;i=n;i+) xi=0; yi=xi; k=0; while(k!=N) k+; for(i=1;i=

19、n;i+) w=0; for(j=1;j=n;j+) if(j!=i) w=w+aij*yj; yi=(bi-w)/double(aii); max=fabs(x1-y1); for(i=1;imax) max=se; if(maxe) coutendl; for(i=1;i=n;i+) coutxi=yiendl; break; for(i=1;i=n;i+) xi=yi; if(k=N) cout迭代失败！endl; choice=0;3、三对角形线性方程组 = x= ( 2, 1, -3, 0, 1, -2, 3, 0, 1, -1 ) SOR方法：# include # include

20、 #include/*定义全局变量*/float *a; /*存放A矩阵*/float *b; /*存放b矩阵*/float *x; /*存放x矩阵*/float p; /*精确度*/float w; /*松弛因子*/int n; /*未知数个数*/int c; /*最大迭代次数*/int k=1; /*实际迭代次数*/*SOR迭代法*/void SOR(float xk) int i,j; float t=0.0; float tt=0.0; float *xl; xl=(float *)malloc(sizeof(float)*(n+1); for(i=1;in+1;i+) t=0.0;

21、tt=0.0; for(j=1;ji;j+) t=t+aij*xlj; for(j=i;jn+1;j+) tt=tt+aij*xkj; xli=xki+w*(bi-t-tt)/aii; t=0.0; for(i=1;in+1;i+) tt=fabs(xli-xki); tt=tt*tt; t+=tt; t=sqrt(t); for(i=1;ic) if(tp) k+; SOR(xk); else printf(nReach the given precision!n); printf(nCount number is %dn,k); /*程序*开始*/void main() int i,j;p

22、rintf(SOR方法n);printf(请输入方程个数:n);scanf(%d,&n);a=(float *)malloc(sizeof(float)*(n+1); for(i=0;in+1;i+) ai=(float*)malloc(sizeof(float)*(n+1);printf(请输入三对角矩阵:n);for(i=1;in+1;i+) for(j=1;jn+1;j+) scanf(%f,&aij);for(i=1;in+1;i+) for(j=1;jn;j+)b=(float *)malloc(sizeof(float)*(n+1);printf(请输入等号右边的值:n);for(i=1;in+1;i+) scanf(%f,&bi);x=(float *)malloc(sizeof(float)*(n+1);printf(请输入初始的x:);for(i=1;in+1;i+) scanf(%f,&xi);printf(请输入精确度:);scanf(%f,&p);printf(请输入迭代次数:);sca