数值分析第一次大作业(宋宝瑞).doc

1、猜想ln(cond(Hn)与n之间在n较小时满足线性关系，在n较大时则不稳定，因此，设n分别为10,20,50,100，如下图所示，可以看出在n14时，函数曲线开始趋于平稳并不断波动，波动幅度较小，并有上升趋势。程序：N=100;y=zeros(N,1);for n=1:N Hn=hilb(n); condn=cond(Hn); y(n)=log(condn);endx=1:N;subplot(224);plot(x,y);xlabel(n);ylabel(ln(cond(Hn);title(ln(cond(Hn)-n(n=100);grid ;2、可以看出，log(cond(Hnpre)/cond(Hn)在y=0上下波动。对比整个曲线图，可以发现，经过预处理，条件数有所下降，但下降值不大。程序：x=1:100;y=zeros(100,1);y1=zeros(100,1);y2=zeros(100,1);for n=1:100 Hn=hilb(n); inv_D=zeros(n); for i=1:n inv_D(i,i)=1/sqrt(Hn(i,i); end Hnpre=inv_D*Hn*inv_D; y(n)=log(cond(Hnpre)/cond(Hn); y1(n)=log(cond(Hnpre); y2(n)=log(cond(Hn);end%plot(x,y1);plot(x,y,r,x,y1,g,x,y2,k);legend(y,y1,y2);xlabel(n);ylabel(log(cond(Hnpre)/cond(Hn);title(log(cond(Hnpre)/cond(Hn)-n(n=100);grid on;3、 设b所有值为1，nx1x2x34 -4.0000E+00-4.0000E+00-4.0000E+006.0000E+016.0000E+016.0000E+01-1.8000E+02-1.8000E+02-1.8000E+021.4000E+021.4000E+021.4000E+025 5.0000E+005.0000E+005.0000E+00-1.2000E+02-1.2000E+02-1.2000E+026.3000E+026.3000E+026.3000E+02-1.1200E+03-1.1200E+03-1.1200E+036.3000E+026.3000E+026.3000E+026 -6.0000E+00-6.0000E+00-6.0000E+002.1000E+022.1000E+022.1000E+02-1.6800E+03-1.6800E+03-1.6800E+035.0400E+035.0400E+035.0400E+03-6.3000E+03-6.3000E+03-6.3000E+032.7720E+032.7720E+032.7720E+037 7.0000E+007.0000E+007.0000E+00-3.3600E+02-3.3600E+02-3.3600E+023.7800E+033.7800E+033.7800E+03-1.6800E+04-1.6800E+04-1.6800E+043.4650E+043.4650E+043.4650E+04-3.3264E+04-3.3264E+04-3.3264E+041.2012E+041.2012E+041.2012E+048 -8.0000E+00-8.0000E+00-8.0000E+005.0400E+025.0400E+025.0400E+02-7.5600E+03-7.5600E+03-7.5600E+034.6200E+044.6200E+044.6200E+04-1.3860E+05-1.3860E+05-1.3860E+052.1622E+052.1622E+052.1622E+05-1.6817E+05-1.6817E+05-1.6817E+055.1480E+045.1480E+045.1480E+049 9.0000E+009.0000E+009.0001E+00-7.2000E+02-7.2000E+02-7.2000E+021.3860E+041.3860E+041.3860E+04-1.1088E+05-1.1088E+05-1.1088E+054.5045E+054.5045E+054.5045E+05-1.0090E+06-1.0090E+06-1.0090E+061.2613E+061.2613E+061.2613E+06-8.2368E+05-8.2368E+05-8.2368E+052.1879E+052.1879E+052.1879E+0510 -9.9980E+00-9.9987E+00-1.0001E+019.8983E+029.8989E+029.9005E+02-2.3756E+04-2.3758E+04-2.3761E+042.4021E+052.4022E+052.4025E+05-1.2611E+06-1.2612E+06-1.2613E+063.7833E+063.7835E+063.7839E+06-6.7260E+06-6.7263E+06-6.7269E+067.0006E+067.0008E+067.0015E+06-3.9378E+06-3.9380E+06-3.9383E+069.2370E+059.2373E+059.2380E+0511 1.0927E+011.0988E+011.0992E+01-1.3124E+03-1.3188E+03-1.3192E+033.8410E+043.8580E+043.8587E+04-4.7823E+05-4.8015E+05-4.8022E+053.1396E+063.1513E+063.1515E+06-1.2060E+07-1.2102E+07-1.2102E+072.8484E+072.8575E+072.8576E+07-4.1864E+07-4.1989E+07-4.1990E+073.7293E+073.7398E+073.7398E+07-1.8420E+07-1.8469E+07-1.8469E+073.8689E+063.8786E+063.8785E+0612 -9.6349E+00-1.2768E+01-1.1739E+011.4311E+031.8176E+031.6833E+03-5.1063E+04-6.3363E+04-5.9039E+047.7783E+059.4708E+058.8711E+05-6.3022E+06-7.5528E+06-7.1069E+063.0320E+073.5851E+073.3868E+07-9.1788E+07-1.0728E+08-1.0171E+081.7935E+082.0753E+081.9736E+08-2.2572E+08-2.5889E+08-2.4688E+081.7662E+082.0099E+081.9214E+08-7.8125E+07-8.8290E+07-8.4591E+071.4921E+071.6757E+071.6087E+07可以看出，当n9时，三种方法计算出的结果相差越来越大。程序：for n=12; b=ones(n,1); Hn=hilb(n); D=zeros(n); for i=1:n D(i,i)=sqrt(Hn(i,i); end Hnpre=inv(D)*Hn*inv(D); x1=inv(Hn)*b; x=inv(Hnpre)*inv(D)*b; x2=inv(D)*x; x3=Hnb; y=x1 x2 x3;End4、 取零向量为初始解时，对于所有的n值，只需迭代两次便得到所要求的结果，而取所有数为1作初始解时，则有如下结果456789101112k100110011001100110011001100110011001x0.9981 0.9992 1.0007 1.0019 1.0028 1.0033 1.0032 1.0027 1.0018 0.0208 0.0043 -0.0131 -0.0252 -0.0329 -0.0344 -0.0297 -0.0202 -0.0078 -0.0489 0.0061 0.0463 0.0679 0.0728 0.0618 0.0388 0.0097 -0.0205 0.0312 -0.0328 -0.0374 -0.0238 -0.0024 0.0187 0.0347 0.0438 0.0464 0.0240 -0.0296 -0.0485 -0.0458 -0.0299 -0.0085 0.0129 0.0312 0.0338 -0.0165 -0.0411 -0.0458 -0.0374 -0.0221 -0.0045 0.0455 -0.0035 -0.0308 -0.0409 -0.0393 -0.0308 0.0521 0.0055 -0.0226 -0.0362 -0.0397 0.0542 0.0110 -0.0166 -0.0323 0.0541 0.0147 -0.0119 0.0535 0.0177 0.0533 可以看出计算结果不收敛，因此我们可以得出结论：用高斯-赛德尔迭代法解此方程组并不是对任意的初始向量都收敛。Hilbert 矩阵是病态的，对于对原Hilbert 矩阵进行预处理后的新Hilbert 矩阵的条件数在一定范围内呈波动的变化规律。从整体上观察，对于大多数的n，进行预处理后的Hilbert 矩阵地条件数有明显的降低，这就说明预处理改善了大多数Hilbert 矩阵的性质。用迭代法求解方程时，如果系数矩阵为Hilbert 矩阵，则求解结果的误差较大，根据迭代法基本定理可知用高斯-赛德尔迭代法解系数矩阵为Hilbert 矩阵时是不收敛的。当n越大时，Hn的病态越严重。程序：for n=4:12; Hn=hilb(n); b=Hn(:,1); x k=Gauss_Seidel_iterative(Hn,b)EndGauss_Seidel_iterative.mfunction x,k=Gauss_Seidel_iterative(A,b)% 用Gauss_Seidel迭代求解线性方程组,矩阵A是方阵 x0=ones(1,length(b); % 赋初值 tol=10(-9); % 给定误差界 N=1000; % 给定最大迭代次数 n,n=size(A); % 确定矩阵A的阶 k=1; % 迭代过程 while k=N x(1)=(b(1)-A(1,2:n)*x0(2:n)/A(1,1); for i=2:n x(i)=(b(i)-A(i,1:i-1)*x(1:i-1)-A(i,i+1:n)*x0(i+1:n)/A(i,i); end if max(abs(x-x0)=tol fid = fopen(G_S_iter_result.txt, wt); fprintf(fid,n*用Gauss_Seidel迭代求解线性方程组的输出结果*nn); fprintf(fid,迭代次数: %d次nn