矩阵与数值分析

《矩阵与数值分析》由会员上传分享，免费在线阅读，更多相关内容在学术论文-天天文库。

1、2013级工科硕士研究生《矩阵与数值分析》课程数值实验题目N1Q6一、设分别编制从小到大和从大到小的顺序程序分别计算7=2J-1o10000,J1000000并指山两种方法计算结果的有效位数。Matlab程序如下：function[$i,sd]=S(N)formatlong;si=0;sd=0;forj=N:-l:2si=1.0e6/(jA2-l)+si;endforj=2:Nsd=1.0e6/(jA2-l)+sd;endend在matlab命令窗口中输入：[si,sd]=S(10000)运行结果：si=7.499000049995000e+005sd=7.49900004

2、9994994e+005在matlab命令窗口中输入：[si,sd]=S(1000000)运行结果：si=7.499990000005000e+005sd=7.499990000005200e+0051结果分析:si为从大到小的顺序求和的值，sd为从小到大的顺序求和的值。当N分别为10000和1000000时，si分别为7.499000049995000e+005和7.499990000005000e+005,可以看出这两个数的有效值均为13位；而sd分别为7.499000049994994e+005和7.499990000005200e+005,这两个数的有效值均为16位

3、。这就出现了我们在矩阵理论课上所学的“大数吃小数”的问题。为了使结果更为精确我们必须避免在四则运算屮出现“大数吃小数”的情况，应该按从小到大的顺序进行求和。二、解线性方程组1.分別利用Jacobi迭代法和Gauss-Seidel迭代法求解线性方程组/U=/>,其屮常向fi为维随机生成的列向y:,系数矩阵A具有如下形式其中1-、z2-1Tn^+2/^,w-lH-l为《-1阶矩阵，L为n-1阶单位矩阵，迭代法计算停止的条件-12为：k+r'l

4、=jz(n)fori=l:n-l%此处n赋值即可，如n=100forj=l:n-lif(i==j)T(ij)=2;endif(i==(j+D)T(i,j)=-1；endif(i==(j_l))T(ij)=-1；endendendI=eye(n-1);k=l;formm=l:(n-l)A(k:(k+n-2),k:(k+n-2))=T+2*I;k=k+n-l;endk=l;forxx=l:(n-2)A(k:(k+n-2)，(k+n-l):(k+2*n-3))=-I;A((k+n-l):(k+2*n-3),k:(k+n-2))=-I;k=k+n-l;endb=randn((n-1

5、)A2,1);xO=zeros((n-l)A2,1);Jacobi迭代法Matlab程序如下：functionn=jacobi(A,b,xO)eps=1.0e-10;M=10000;D=diag(diag(A));L=-tril(A，-l);U=-triu(A,l);B=D(L+U);f=Db;x=B*xO+f;n=l;whilenorm(x-xO)>=epsx0=x;x=B*xO+f;n=n+l;if(n〉=M)dispCWarning:迭代次数太多，可能不收敛return;endendGauss-Seidel迭代法Matlab程序如下：functionn=gause

6、idel(A，b，xO)eps=1.0e-10;%求A的对角矩阵%求A的下三角阵%求A的上三角阵%迭代次数M=10000;D=diag(diag(A));L=-tril(A,-l);U=-triu(A,l);G=(D-L)U;f=(D-L)b;x=G*xO+f;n=l;whilenorm(x-xO)〉=epsx0=x;x=G*xO+f;n=n+l;if(n>=M)dispCWarning:return;end%求A的对角矩阵%求A的下三角阵%求A的上三角阵%迭代次数迭代次数太多，可能不收敛end在matlab命令窗口中输入:[A,b,x0J=jz(10);J10=jac

7、obi(A，b，x0)G1O=gauseidel(A,b,xO)[A,b,x0]=jz(20);J20=jacobi(A,b,x0)G20=gauseidel(A,b,x0)[A,b,xO]=jz(3O);J30=jacobi(A,b,x0)G3O=gauseidel(A,b,xO)运行结果：J10=436；GIO=214J20=1810；G20=913J30=3990；G30=2001结果分析：N=10且M=1000吋，Jacobi迭代法和Gauss—seidel迭代法的迭代次数分别为436和214；N=20且M