数值分析实验报告1——Hilbert矩阵的求解
矩阵分析与数值分析实验报告
《矩阵分析与数值分析》实验报告院系：姓名：学号：所在班号：任课老师：一．设错误！未找到引用源。
，分别编制从小到大和从大到小的顺序程序计算错误！未找到引用源。
并指出有效位数。
程序如下：function sum3j=input('请输入求和个数 "j":');A=0;B=0;double B;double A;for n=2:jm=n^2-1;t=1./m;A=A+t;enddisp('从小到大：')s=Afor n=j:-1:2m=n^2-1;t=1./m;B=B+t;enddisp('从大到小:')s=B运行结果：>> sum3请输入求和个数 "j":100从小到大：s =0.740049504950495从大到小:s =0.740049504950495>> sum3请输入求和个数 "j":10000从小到大：s =0.749900004999506从大到小:s =0.749900004999500>> sum3请输入求和个数 "j":1000000从小到大：s =0.749999000000522从大到小:s =0.749999000000500二、解线性方程组1．分别Jacobi 迭代法和Gauss-Seidel 迭代法求解线性方程组。
⎪⎪⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎭⎫ ⎝⎛----000121001210012100124321x x x x 迭代法计算停止的条件为：6)()1(3110max -+≤≤<-k j k j j x x 。
解：（1）Jacobi 迭代法程序代码： function jacobi(A, b, N) clc;clear;A=[-2 1 0 0;1 -2 1 0;0 1 -2 1;0 0 1 -2]; b=[-1 0 0 0]'; N=100;n = size(A,1); D = diag(diag(A)); L = tril(-A,-1); U = triu(-A,1); Tj = inv(D)*(L+U); cj = inv(D)*b; tol = 1e-06; k = 1;format longx = zeros(n,1); while k <= Nx(:,k+1) = Tj*x(:,k) + cj;disp(k); disp('x = ');disp(x(:,k+1)); if norm(x(:,k+1)-x(:,k)) < toldisp('The procedure was successful')disp('Condition ||x^(k+1) - x^(k)|| < tol was met after k iterations') disp(k); disp('x = ');disp(x(:,k+1)); break endk = k+1; end结果输出The procedure was successfulCondition ||x^(k+1) - x^(k)|| < tol was met after k iterations 60 x =0.799998799067310.599998427958700.399998056850090.19999902842505（2）Gauss-Seidel迭代法程序代码：function gauss_seidel(A, b, N)clc;clear;A=[-2 1 0 0;1 -2 1 0;0 1 -2 1;0 0 1 -2];b=[-1 0 0 0]';N=100;n = size(A,1);D = diag(diag(A));L = tril(-A,-1);U = triu(-A,1);Tg = inv(D-L)*U;cg = inv(D-L)*b;tol = 1e-06;k = 1;x = zeros(n,1);while k <= Nx(:,k+1) = Tg*x(:,k) + cg;disp(k); disp('x = ');disp(x(:,k+1));if norm(x(:,k+1)-x(:,k)) < toldisp('The procedure was successful')disp('Condition ||x^(k+1) - x^(k)|| < tol was met after k iterations') disp(k); disp('x = ');disp(x(:,k+1));breakendk = k+1;end结果输出The procedure was successfulCondition ||x^(k+1) - x^(k)|| < tol was met after k iterations31x =0.799999213979350.599998971085610.399999167590770.199999583795392. 用Gauss列主元消去法、QR方法求解如下方程组：⎪⎪⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎭⎫ ⎝⎛---017435222331325212214321x x x x (1)Gauss 列主元消去法 程序代码：function x=Gaussmain(A,b) clc;clear; format longA=[1 2 1 2;2 5 3 -2;-2 -2 3 5;1 3 2 3]; b=[4 7 -1 0]'; N=length(A); x=zeros(N,1); y=zeros(N,1); c=0; d=0;A(:,N+1)=b; for k=1:N-1 for i=k:4if c<abs(A(i,k))d=i;c=abs(A(i,k)); end endy=A(k,:);A(k,:)=A(d,:); A(d,:)=y; for i=k+1:N c=A(i,k);for j=1:N+1A(i,j)=A(i,j)-A(k,j)*c/A(k,k); end end endb=A(:,N+1);x(N)=b(N)/A(N,N); for k=N-1:-1:1x(k)=(b(k)-A(k,k+1:N)*x(k+1:N))/A(k,k); end结果输出 ans =18.00000000000000 -9.571428571428576.00000000000000-0.42857142857143（2）QR方法：程序代码function QR(A,b)clc;clear;format longA=[1 2 1 2;2 5 3 -2;-2 -2 3 5;1 3 2 3];b=[4 7 -1 0]';[Q,R]=qr(A)y=Q'*b;x=R\y结果输出Q =-0.31622776601684 -0.04116934847963 -0.75164602800283 0.57735026918962 -0.63245553203368 -0.49403218175557 -0.15032920560056 -0.57735026918963 0.63245553203368 -0.74104827263336 -0.22549380840085 -0.00000000000000 -0.31622776601684 -0.45286283327594 0.60131682240226 0.57735026918963 R =-3.16227766016838 -6.00832755431992 -0.94868329805051 2.84604989415154 0 -2.42899156029822 -4.65213637819829 -4.15810419644272 0 0 -0.67648142520255 -0.52615221960200 0 0 0 4.04145188432738 x =17.99999999999989-9.571428571428515.99999999999997-0.42857142857143三、非线性方程的迭代解法1．用Newton迭代法求方程()06cos22x=-++=-xexf x的根，计算停止的条件为：6110-+<-kkxx；编程如下：function newton(f,df,x,a,a0)syms xf=input('please enter your equation:') a0=input('please enter you x(0):');df=diff(f)e=1e-6;a1=a0+1;N=0;while abs(a1-a0)>ea=a0-subs(f,a0)/subs(df,a0); a1=a0; a0=a; N=N+1; endfprintf('a=%0.6f',a) N运行结果： >> newtonplease enter your equation:exp(x)+2^(-x)+2*cos(x)-6 f =exp(x)+2^(-x)+2*cos(x)-6 please enter you x(0):2df =exp(x)-2^(-x)*log(2)-2*sin(x) a=1.829384 N =42．利用Newton 迭代法求多项式07951.2954.856.104.5x 234=+-+-x x x的所有实零点，注意重根的问题。
Hilbert代数方程组的数值解法
Hilbert 代数方程组的数值解法1 Hilbert 矩阵的条件数和矩阵的阶数的关系编制Matlab 程序：clear,clc format long %format short for n=1:14; A=hilb(n); %A=rand(n);condA(n)=cond(A,inf); endplot(log10(condA)) 得出图形图1：Hilbert 矩阵和阶次的关系由图可见a.Hilbert 矩阵的2-条件数的对数与阶次近似正比;b.阶次超过12后，计算困难，舍入误差会导致计算不准 结论：Hilbert 矩阵的2-条件数随阶次指数增长，关系大概是： Log 10(Cond(A))= 1.33791720780254(n-1)0246810121424681012141618nl o g 10(c o n d A )条件数的对数-阶次2.各种求解方法的对比2.1 直接法：Gauss 消去方法（程序见附录，下同）在双精度型变量下，即eps=2.220446049250313e-016，对于直接法Gauss 消去方法，求解Hilbert 矩阵方程组，在阶次n=13时，解得: X=[1 0.99997 1.0013 0.97881 1.1906 -0.035441 4.6171 -7.3967 14.089 -12.542 9.9162 -2.3817 1.5624]出现明显错误，可见Gauss 消去方法对病态问题比较敏感。
求解阶次不高。
2.2 Jacobi 迭代法很遗憾，jakobi 迭代法的迭代矩阵的谱半径随阶次线性上升（程序见附录2），普遍大于1，计算结果发送，无法迭代出满意的结果。
需放弃图2：Jacobi 迭代矩阵的谱半径2.3 Gauss-Seidel 迭代与SOR 迭代求解先分析SOR 迭代矩阵的谱半径和松弛因子w 的关系0510152025510152025阶次谱半径Jacobi 迭代收敛性分析图3：SOR 迭代矩阵的谱半径和松弛因子w 的关系可以发现SOR 迭代和G-S 迭代的谱半径都普遍接近于1，因此松弛因子的选择对计算结果的影响不大。
数值分析计算方法实验报告
end;
end;
X=x;
disp('迭代结果:');
X
format short;
输出结果：
因为不收敛，故出现上述情况。
4.超松弛迭代法：
%SOR法求解实验1
%w=1.45
%方程组系数矩阵
clc;
A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13]
b=[10,5,-2,7]'
b=[10,5,-2,7]'
[m,n]=size(A);
if m~=n
error('矩阵A的行数和列数必须相同');
return;
end
if m~=size(b)
error('b的大小必须和A的行数或A的列数相同');
return;
end
if rank(A)~=rank([A,b])
error('A矩阵的秩和增广矩阵的秩不相同,方程不存在唯一解');
3．实验环境及实验文件存档名
写出实验环境及实验文件存档名
4.实验结果及分析
输出计算结果，结果分析和小结等。
解：1.高斯列主元消去法：
%用高斯列主元消去法解实验1
%高斯列主元消元法求解线性方程组Ax=b
%A为输入矩阵系数,b为方程组右端系数
%方程组的解保存在x变量中
format long;
A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13]
return;
end
c=n+1;
A(:,c)=b;
for k=1:n-1
数值分析Hilbert矩阵病态线性方程组的求解Matlab程序
（Hilbert矩阵）病态线性方程组的求解理论分析表明，数值求解病态线性方程组很困难。
考虑求解如下的线性方程组的求解1Hx=b,期中H是Hilbert矩阵,H=（h j）四，h j=,।,j=1,2，…，nij-11,估计矩阵的2条件数和阶数的关系2 .对不同的n,取x=（1,1，…，1）w n,分别用Gauss消去，Jacobi迭代，Gauss-seidel迭代，SOR迭代和共轲梯度法求解，比较结果。
3 .结合计算结果，试讨论病态线性方程组的求解。
第1小题：condition.m。
蟆1小题程序t1=20;%>数n=20x1=1:t1;y1=1:t1;fori=1:t1H=hilb(i);y1(i)=log(cond(H));endplot(x1,y1);xlabel('阶数n');ylabel('2-条件数的对数(log(cond(H))');title('2-条件数的对数(10g(cond(H))与阶数n的关系图’);t2=200;%J>数n=200x2=1:t2;y2=1:t2;fori=1:t2H=hilb(i);y2(i)=log(cond(H));endplot(x2,y2);xlabel('阶数n');ylabel('2-条件数的对数(log(cond(H))');title('2-条件数的对数(10g(cond(H))与阶数n的关系图’);画出Hilbert矩阵2-条件数的对数和阶数的关系n=200时n=20时马KiMn的美方四*--wficEJ-帛=«』从图中可以看出，1)在n小于等于13之前，图像近似直线log(cond(H))-1.519n-1.8332)在n大于13之后，图像趋于平缓，并在一定范围内上下波动，同时随着n的增加稍有上升的趋势第2小题：soke.m%蹴2小题主程序N=4000;xGauss=zeros(N,1);xJacobi=zeros(N,1);xnJ=zeros(N,1);xGS=zeros(N,1);xnGS=zeros(N,1);xSOR=zeros(N,1);xnSOR=zeros(N,1);xCG=zeros(N,1);xnCG=zeros(N,1);forn=1:N;x=ones(n,1);t=1.1;%初始值偏差x0=t*x;%迭代初始值e=1.0e-8;%^定的误差A=hilb(n);b=A*x;max=100000000000;%可能最大的迭代次数w=0.5;%SORt代的松弛因子G=Gauss(A,b);[J,nJ]=Jacobi(A,b,x0,e,max);[GS,nGS]=G_S(A,b,x0,e,max);[S_R,nS_R]=SOR(A,b,x0,e,max,w);[C_G,nC_G]=CG(A,b,x0,e,max);normG=norm(G'-x);xGauss(n)=normG;normJ=norm(J-x);nJ;xJacobi(n)=normJ;xnJ(n)=nJ;normGS=norm(GS-x);nGS;xGS(n)=normGS;xnGS(n)=nGS;normS_R=norm(S_R-x);nS_R;xSOR(n)=normS_R;xnSOR(n)=nS_R;normC_G=norm(C_G-x);nC_G;xCG(n)=normC_G;xnCG(n)=nC_G;endGauss.m%Gauss消去法functionx=Gauss(A,b)n=length(b);l=zeros(n,n);x=zeros(1,n);脸肖去过程fori=1:n-1forj=i+1:nl(j,i)=A(j,i)/A(i,i);fork=i:nA(j,k)=A(j,k)-l(j,i)*A(i,k);endb(j)=b(j)-l(j,i)*b(i);endend%回代过程x(n)=b(n)/A(n,n);fori=n-1:-1:1c=A(i,:).*x;x(i)=(b(i)-sum(c(i+1:n)))/A(i,i);endJacobi.m%Jacobi迭代,x0表示迭代初值,e表示允许误差(迭代停止条件),n表示迭代次数,m可能最大的迭代次数function[x,n]=Jacobi(A,b,x0,e,m)n=length(A);D=diag(diag(A));U=-triu(A,1);L=-tril(A,-1);B=D\(L+U);f=D\b;x=B*x0+f;n=1;whilenorm(x-x0)>ex0=x;x=B*x0+f;n=n+1;ifn>mdisp('Jacobi迭代次数过多，迭代可能不收敛‘);break;endendG_S.m%Gauss-Seidel迭代,x0表示迭代初值，e表示允许误差(迭代停止条件),n表示迭代次数,m可能最大的迭代次数function[x,n]=G_S(A,b,x0,e,m)n=length(A);D=diag(diag(A));U=-triu(A,1);L=-tril(A,-1);B=(D-L)\U;f=(D-L)\b;x=B*x0+f;n=1;whilenorm(x-x0)>ex0=x;x=B*x0+f;n=n+1;ifn>mdisp('Gauss-Seidel迭代次数过多，迭代可能不收敛’)；break;endendSOR.m%SOR超松弛迭代，x0表示迭代初值,e表示允许误差(迭代停止条件),n表示迭代次数,m可能最大的迭代次数,w松弛因子function[x,n]=SOR(A,b,x0,e,m,w)n=length(A);D=diag(diag(A));U=-triu(A,1);L=-tril(A,-1);B=(D-w*L)\((1-w)*D+w*U);f=(D-w*L)\b*w;x=B*x0+f;n=1;whilenorm(x-x0)>ex0=x;x=B*x0+f;n=n+1;ifn>mdisp('SOR超松弛迭代次数过多，迭代可能不收敛’)；break;endendCG.m%C献朝梯度法，x0表示迭代初值，e表示允许误差(迭代停止条件),n表示迭代次数,m可能最大的迭代次数function[x,n]=CG(A,b,x0,e,m)r=b-A*x0;p=r;alpha=(r'*r)/(p'*(A*p));x=x0+alpha*p;r1=b-A*x;n=1;whilenorm(r1)>ebelta=(r1'*r1)/(r'*r);p=r1+belta*p;r=r1;x0=x;alpha=(r'*r)/(p'*(A*p));x=x0+alpha*p;r1=b-A*x;n=n+1;ifn>mdisp('CG共朝梯度法迭代次数过多，迭代可能不收敛’);break;endend第2小题选择不同的阶数n,取x=(l,l,…RI分别使用Gauss消去，Jacobi迭代，Gauss-seidel迭代,SOR迭代和共挽梯度法(CG)求解，比较结果。
数值分析_线性方程组迭代解法Hilbert矩阵
0.999427 1.000368 1.000068 0.999538 1.003545 0.999834 1.000724 1.002775 0.99582 x 9 0.99888 0.9983 0.99933 0.995748 0.999679
0.996358 0.999739
1.002141 1.001106 1.001204 1.001784 1.005742 1.00315 1.00155 1.00181 1.002206 1.00062
1.000066 0.999977 1.000001 1.000041 0.998534 1.000054 0.999712 0.999029 1.005874 1.000672 1.00175 1.004447
0.996128 0.998886 0.998141 0.994858 10 x 0.993805 0.998912 0.998328 0.998514 0.999227 1.000283 1.000399 1.001125 1.005241 1.001468 1.001894 1.002573 1.006835 1.001551 1.001869 1.002117
数值分析课程上机报告
数值分析第二次上机实习报告
——线性方程组迭代解法
一、问题描述
设 Hn = [hij ] ∈ Rn×n 是 Hilbert 矩阵, 即 hij = 对 n = 2,3,4，…15, 1 i + j −1
1 x ∈ R n×n ，及 bn = H n x ，用 SOR 迭代法和共轭梯度法来求解，并与直 取= 1
四、计算结果及其分析
1. 超松弛迭代法分析 令初始向量 x0=u(1,1,……)T，给定不同的初始向量与松弛因子，计算不同情 况下解的情况。计算结果如下。
数值求解Hilbert病态线性方程组
病态线性代数方程组的求解理论的分析表明，求解病态的线性代数方程组是困难的。
考虑方程组Hx = b 的求解，其中H 为Hilbert 矩阵，n n ij h H ⨯=)(，11-+=j i h ij ，n j i ,...,2,1,=1. 估计Hilbert 矩阵2-条件数与阶数的关系；2. 选择问题的不同维数，分别用Gauss 消去法，Jacobi 迭代，GS 迭代和SOR 迭代求解，比较结果；3. 讨论病态问题求解的算法。
解： 1、取Hilbert 矩阵阶数最高分别为n=20和n=100。
采用Hilbert 矩阵的2-条件数作为实验的比较对象，画出的曲线如下图所示：lg(())n cond H n从图中可以看出，在n ≤13之前，图像接近直线，在n ＞13之后，图像趋于平缓，在一定的范围内上下波动。
为了比较图像的线性部分，作出一条直线与已知曲线进行比较。
比较直线的关系式为：833.1519.1))(lg(-=n H cond n ，结果下图所示。
nl g (c o n d (H n ))lg(cond(Hn))~n 关系图从图2中可以看出，当n 较小时，n H cond n ~))(lg(之间近似满足线性关系。
当n 继续增大到100时，n H cond n ~))(lg(关系图下图所示：从图中可以看出，图像的走势符合在n=20时的猜想，在n 大于一定的值之后，图像趋于平缓，且在一定范围内震荡，同时又有一定上升趋势，但上升速度很慢。
2、选择不同的阶数n ，设方程组的精确解为xz=(1,1,…,1)T进行计算，用四种方法解x_Guass1、x_Jacobi1、x_GS1、x_SOR1对比表如下nl g (c o n d (H n ))lg(cond(Hn))~n关系图nl g (c o n d (H n ))lg(cond(Hn))~n 关系图Gauss消去法求解：选择问题的阶数为3~8时，用Gauss消去法求得的解与精确解一致，当阶数为9~14时，解开始出现偏差，而且n越大，偏差越大。
数值分析实验报告1
显然该多项式的全部根为1,2,…，20共计20个，且每个根都是单重的。现考虑该多项式的一个扰动
其中是一个非常小的数。这相当于是对（1。1）中的系数作一个小的扰动.我们希望比较(1。1）和（1。2）根的差别,从而分析方程(1。1）的解对扰动的敏感性.
实验内容：为了实现方便，我们先介绍两个Matlab函数:“roots”和“poly”。
（3)取矩阵阶数n=20或者更大，重复上述实验过程，观察记录并分析不同的问题及消去过程中选择不同的主元时计算结果的差异，说明主元素的选取在消去过程中的作用。
（4)选取其他你感兴趣的问题或者随机生成矩阵，计算其条件数。重复上述实验，观察记录并分析实验结果。
实验过程：
程序：
建立M文件:
function x=gauss（n,r）
实验总结：
利用MATLAB来进行病态问题的实验,虽然其得出的结果是有误差的，但是可以很容易的得出对一个多次的代数多项式的其中某一项进行很小的扰动，对其多项式的根会有一定的扰动的，所以对于这类病态问题可以借助于MATLAB来进行问题的分析。
学号：06450210
姓名：万轩
实验二插值法
实验2.1(多项式插值的振荡现象)
学号:06450210
姓名：万轩
实验五解线性方程组的直接方法
实验5。1（主元的选取与算法的稳定性）
问题提出：Gauss消去法是我们在线性代数中已经熟悉的.但由于计算机的数值运算是在一个有限的浮点数集合上进行的，如何才能确保Gauss消去法作为数值算法的稳定性呢？Gauss消去法从理论算法到数值算法，其关键是主元的选择。主元的选择从数学理论上看起来平凡，它却是数值分析中十分典型的问题.
其中若变量a存储n+1维的向量,则该函数的输出u为一个n维的向量。设a的元素依次为，则输出u的各分量是多项式方程
数值分析实验报告(包括高斯消去、二分法、牛顿迭代法)
for k=1:N
x=(a+b)/2;
fx=feval(f,x);fa=feval(f,a);
if abs((b-a)/2)<e || abs(fx)<e
disp('the number of iterations is');k
f=input('please enter a function:f(x)=');
x0=input('please enter the initial value:x0=');
e=input('please enter error:e=');
N=input('please enter the largest number of iterations:N=');
disp('the approximate solution is');x
disp('f(x) is');fx
disp('the number of iterations is');k
return
else
x0=x;
end
end
end
disp('The maximum number of iterations is reached, stop calculation');
开课学院、实验室：实验时间：2014年1月1日
课程
名称
数值分析基础性实验
实验项目
名称
数值计算算法及实现
数值分析——线性方程组直接解法Hilbert矩阵
数值分析第一次上机实习报告——线性方程组直接解法一、问题描述设 H n = [h ij ] ∈ R n ×n 是 Hilbert 矩阵, 即11ij h i j =+- 对n = 2,3,4，…13,(a) 取11n n x R ⨯⎛⎫ ⎪=∈ ⎪ ⎪⎝⎭，及n n b H x =，用Gauss 消去法和Cholesky 分解方法来求解n n H y b =，看看误差有多大.(b) 计算条件数：2()n cond H(c) 使用某种正则化方法改善(a)中的结果.二、方法描述1. Gauss 消去法Gauss 消去法一般用于系数矩阵稠密且没有任何特殊结构的线性方程组。
设H =[h ij ]，y = (y 1,y 2,…，y n )T . 首先对系数矩阵H n 进行LU 分解，对于k=1,2，…n,交替进行计算：1111),,1,,1(),1,2,,k kj kj kr rj r k ik ik ir rk r kk u h l u j k k n l a l u i k k n u -=-=⎧=-=+⎪⎪⎨⎪=-=++⎪⎩∑∑…… 由此可以得到下三角矩阵L=[l ij ]和上三角矩阵U=[u ij ]. 依次求解方程组Ly=b 和Ux=y ，111,1,2,,1(),,1,,1i i i ir r r n i i ir r r i ii y b l y i n x y u x i n n u -==+⎧=-=⎪⎪⎨⎪=-=-⎪⎩∑∑…… 即可确定最终解。
2. Cholesky 分解法对于系数矩阵对称正定的方程组n n H y b =，存在唯一的对角元素为正数的下三角矩阵L ，使得H=LL T 。
因此，首先对矩阵H n 进行Cholesky 分解，即1122111()1()j jj jj jk k j ij ij ik jk k jj l h l l h l l l -=-=⎧=-⎪⎪⎨⎪=-⎪⎩∑∑ 1,i j n =+… L 的元素求出之后，依次求解方程组Ly=b 和L T x=y ，即1111111(),2,3,i i i ik k k ii b y l y b l y i n l -=⎧=⎪⎪⎨⎪=-=⎪⎩∑… 11(),1,2,n n nn n i i ki k k i nn y x l x y l x i n n l =+⎧=⎪⎪⎨⎪=-=--⎪⎩∑…1 由此求得方程组n n H y b =的解。
Hilbert矩阵病态线性代数方程组的求解
实验一病态线性代数方程组的求解1.估计Hilbert矩阵2-条件数与阶数的关系运行tiaojianshu.m 输入m=10 可以得到如下表的结果2.选择不同维数，分别用Guass消去（LU分解），Jacobi迭代，GS 迭代，SOR迭代求解，比较结果。
说明：Hx=b，H矩阵可以由matlab直接给出，为了设定参考解，我们先设x为分量全1的向量，求出b,然后将H和b作为已知量，求x，与设定的参考解对比。
对于Jacobi迭代，GS迭代，SOR迭代，取迭代初值x0为0向量，迭代精度eps=1.0e-6，迭代次数<100000, SOR迭代中w=1.2和0.8分别计算。
a. n=5b. n=8c. n=10d. n=15取不同的n值，得到如下结果:对于Guass法，可以看出来，随着n的增大，求解结果误差变大，这是因为随着n增大，系数矩阵的条件数变大，微小的扰动就容易造成很大的误差。
最后得不到精确解。
对于Jacobi迭代，计算结果为Inf，说明是发散的。
对于GS迭代和SOR迭代，结果是收敛的，但是可以看出迭代次数比较多，并且对于不同维数GS和SOR收敛速度不一样，有时候GS快，有时SOR快。
对SOR取不同的w迭代速度也不一样，存在一个最优的松弛因子w。
并且可以知道，迭代次数多少跟初值x0也有关系。
3.讨论病态问题求解的算法。
通过上面的实验分析，可以看出，求解病态矩阵的时候要小心，否则可能得不到所要求的精确度。
可以采用高精度运算，用双倍多倍字长，使得由于误差放大而损失若干有效数字位之后，还能保留一些有效位。
另外可以通过对原方程作某些预处理，降低系数矩阵的条件数，因为cond(aA)=cond(A),所以不能通过将每一个方程乘上相同的常数来达到这个目标，可考虑将矩阵的每一行和每一列分别乘上不同的常数，亦即找到可逆的对角阵D1和D2将方程组化为D1AD2y=D1b，x=D2y这称为矩阵的平衡问题，但是这样计算量比原问题本身要多。
数值分析实验报告--实验6--解线性方程组的迭代法
1 / 8数值分析实验六：解线性方程组的迭代法2016113 张威震1 病态线性方程组的求解1.1 问题描述理论的分析表明，求解病态的线性方程组是困难的。
实际情况是否如此，会出现怎样的现象呢？实验内容：考虑方程组Hx=b 的求解，其中系数矩阵H 为Hilbert 矩阵，,,1(),,,1,2,,1i j n n i j H h h i j n i j ⨯===+-这是一个著名的病态问题。
通过首先给定解（例如取为各个分量均为1）再计算出右端b 的办法给出确定的问题。
实验要求：（1）选择问题的维数为6，分别用Gauss 消去法、列主元Gauss 消去法、J 迭代法、GS 迭代法和SOR 迭代法求解方程组，其各自的结果如何？将计算结果与问题的解比较，结论如何？（2）逐步增大问题的维数（至少到100），仍然用上述的方法来解它们，计算的结果如何？计算的结果说明了什么？（3）讨论病态问题求解的算法1.2 算法设计首先编写各种求解方法的函数，Gauss 消去法和列主元高斯消去法使用实验5中编写的函数myGauss.m 即可，Jacobi 迭代法函数文件为myJacobi.m ，GS 迭代法函数文件为myGS.m ，SOR 方法的函数文件为mySOR.m 。
1.3 实验结果1.3.1 不同迭代法球求解方程组的结果比较选择H 为6*6方阵，方程组的精确解为x* = (1, 1, 1, 1, 1, 1)T ，然后用矩阵乘法计算得到b ，再使用Gauss 顺序消去法、Gauss 列主元消去法、Jacobi 迭代法、G-S 迭代法和SOR 方法分别计算得到数值解x1、x2、x3、x4，并计算出各数值解与精确解之间的无穷范数。
Matlab 脚本文件为Experiment6_1.m 。
迭代法的初始解x 0 = (0, 0, 0, 0, 0, 0)T ，收敛准则为||x(k+1)-x(k)||∞<eps=1e-6，SOR方法的松弛因子选择为w=1.3，计算结果如表1。
数值分析实习报告
数值分析实习报告姓名：***学号：***班级：***1.题目取h=0.1,利用Euler公式求解dy/dx=y-2*x/y (0<=x<=1)y(0)=12.思路利用左矩形公式得到的公式yn+1=yn+hf(xn,yn)进行迭代，在进行迭代的过程使用for循环，让n从1取到10，每迭代一次输出一个x与y，同时x加0.1再进行下一次迭代，一直到循环结束。
3.程序clear; y=1, x=0, %初始化for n=1:10y=1.1*y-0.2*x/y, x=x+0.1,end4.运行结果y =1 x =0y =1.1000 x =0.1000y =1.1918 x =0.2000y =1.2774 x =0.3000y =1.3582 x =0.4000y =1.4351 x = 0.5000y =1.5090 x = 0.6000y =1.5803 x = 0.7000y =1.6498 x = 0.8000y =1.7178 x =0.9000y =1.7848 x =1.00001.题目取h=0.1,利用Euler公式求解dy/dx=y-2*x/y (0<=x<=1)y(0)=12.思路利用右矩形公式得到的公式yn+1=yn+hf(xn+1,yn+1)进行迭代，在进行迭代的过程使用for 循环，让n从1取到10，每迭代一次输出一个x与y，同时x加0.1再进行下一次迭代，一直到循环结束。
3.程序clear; y0=1, x0=0, %初始化for n=1:10yp=y0+0.1*(y0-2*x0/y0);x=x0+0.1;x0=xyp=y0+0.1*(yp-2*x/yp)y0=ypend4.运行结果y =1 x =0y =1.0918 x =0.1000y =1.1763 x =0.2000y =1.2546 x =0.3000y =1.3278 x =0.4000y =1.3964 x =0.5000y =1.4609 x =0.6000y =1.5216 x =0.7000y =1.5786 x =0.8000y =1.6321 x =0.9000y =1.6819 x =1.00001.题目取h=0.1,利用梯形法求解dy/dx=y-2*x/y (0<=x<=1)y(0)=12.思路利用右矩形公式得到的公式yn+1=yn+hf(xn+1,yn+1)进行迭代，在进行迭代的过程使用for 循环，让n从1取到10，每迭代一次输出一个x与y，同时x加0.1再进行下一次迭代，一直到循环结束。
Hilbert矩阵病态性分析
Hilbert矩阵病态性分析一、问题叙述Hilbert矩阵是著名的病态矩阵，用迭代法解矩阵方程时，如果系数矩阵是Hilbert矩阵，则求解结果误差较大。
本文就研究Hilbert矩阵的病态性，和Hilbert 矩阵的阶数与迭代求解误差大小的关系。
二、问题分析MATLAB中有专门的Hilbert矩阵及其准确逆矩阵的生成函数，从而可以通过MATLAB求解出系数矩阵是Hilbert矩阵的矩阵方程比较准确的解，再根据条件数估计求解误差的大小。
三、实验程序及注释m=input('input m:='); %输入矩阵的阶数N=[m];for k=1:length(N)n=N(k); %矩阵的阶H=hilb(n); %产生n阶Hilbert矩阵disp(H) %输出n阶Hilbert矩阵Hi=invhilb(n); %产生完全准确的n阶逆Hilbert矩阵b=ones(n,1); %生成n阶全1向量x_approx=H\b; %利用左除H求近似解x_exact=Hi*b; %利用准确逆Hilbert矩阵求准确解ndb=norm(H*x_approx-b);nb=norm(b);ndx=norm(x_approx - x_exact);nx=norm(x_approx);er_actual(k)=ndx/nx; %实际相对误差K=cond(H); %计算Hilbert矩阵的条件数er_approx(k)=K*eps; %最大可能的近似相对误差er_max(k)=K*ndb/nb; %最大可能的相对误差enddisp('Hilbert矩阵阶数'),disp(N)format short edisp('实际误差 er_actual'),disp(er_actual),disp('')disp('近似的最大可能差 er_approx'),disp(er_approx),disp('') disp('最大可能误差 er_max'),disp(er_max),disp('')四、实验数据结果及分析程序运行后，输入矩阵阶数4时的输出为：input m:=41.0000 0.5000 0.3333 0.25000.5000 0.3333 0.2500 0.20000.3333 0.2500 0.2000 0.16670.2500 0.2000 0.1667 0.1429Hilbert矩阵阶数4实际误差er_actual1.0284e-013近似的最大可能误差er_approx3.4447e-012最大可能误差er_max4.7732e-011当输入不同的矩阵阶数时，误差大小如表1所示。
数值求解Hilbert病态线性方程组
病态线性代数方程组的求解理论的分析表明，求解病态的线性代数方程组是困难的。
考虑方程组Hx = b 的求解，其中H 为Hilbert 矩阵，n n ij h H ⨯=)(，11-+=j i h ij ，n j i ,...,2,1,=1. 估计Hilbert 矩阵2-条件数与阶数的关系；2. 选择问题的不同维数，分别用Gauss 消去法，Jacobi 迭代，GS 迭代和SOR 迭代求解，比较结果；3. 讨论病态问题求解的算法。
解： 1、取Hilbert 矩阵阶数最高分别为n=20和n=100。
采用Hilbert 矩阵的2-条件数作为实验的比较对象，画出的曲线如下图所示：lg(())n cond H n从图中可以看出，在n ≤13之前，图像接近直线，在n ＞13之后，图像趋于平缓，在一定的围上下波动。
为了比较图像的线性部分，作出一条直线与已知曲线进行比较。
比较直线的关系式为：833.1519.1))(lg(-=n H cond n ，结果下图所示。
nl g (c o n d (H n ))lg(cond(Hn))~n 关系图从图2中可以看出，当n 较小时，n H cond n ~))(lg(之间近似满足线性关系。
当n 继续增大到100时，n H cond n ~))(lg(关系图下图所示：从图中可以看出，图像的走势符合在n=20时的猜想，在n 大于一定的值之后，图像趋于平缓，且在一定围震荡，同时又有一定上升趋势，但上升速度很慢。
2、选择不同的阶数n ，设方程组的精确解为xz=(1,1,…,1)T进行计算，用四种方法解x_Guass1、x_Jacobi1、x_GS1、x_SOR1对比表如下nl g (c o n d (H n ))lg(cond(Hn))~n关系图nl g (c o n d (H n ))lg(cond(Hn))~n 关系图Gauss消去法求解：选择问题的阶数为3~8时，用Gauss消去法求得的解与精确解一致，当阶数为9~14时，解开始出现偏差，而且n越大，偏差越大。
《数值分析》实验报告2
《数值分析》实验报告一、问题的提出求解线性方程组的迭代法，即是用某种极限过程去逐步逼近线性方程组的精确解的过程，迭代法是解大型稀疏矩阵方程组的重要方法。
二、实验名称运用MATLAB编程实现雅可比（Jacobi）迭代和高斯-赛德尔（Gauss-Seidel）迭代。
三、实验目的1、熟悉了解雅可比（Jacobi）迭代和高斯-赛德尔（Gauss-Seidel）迭代的算法。
2、学习MATLAB软件的功能。
四、基本原理五、实验环境操作环境：Windows10实验平台：Matlab7.1软件六、试验设计1、jacobi迭代法（1）算例：课本p54页例1（2）程序清单Jacobi迭代法的MATLAB函数文件Jacobi.m如下：function [y,n]=jacobi(A,b,x0,eps)if nargin==3eps=1.0e-6;elseif nargin<3errorreturnendD=diag(diag(A)); %求A的对角矩阵L=-tril(A,-1); %求A的下三角阵U=-triu(A,1); %求A的上三角阵B=D\(L+U);f=D\b;y=B*x0+f;n=1; %迭代次数while norm(y-x0)>=epsx0=y;y=B*x0+f;n=n+1;end（3）实验结果及分析:>> A=[8,-3,2;4,11,-1;6,3,12];>> b=[20,33,36]';>> [x,n]=jacobi(A,b,[0,0,0]',1.0e-6)x =3.00002.00001.0000n =162、Gauss-seidel迭代法（1）算例：课本p54页例1（2）程序清单:Gauss-Serdel迭代法的MATLAB函数文件gauseidel.m如下：function [y,n]=gauseidel(A,b,x0,eps)if nargin==3eps=1.0e-6;elseif nargin<3errorreturnendD=diag(diag(A)); %求A的对角矩阵L=-tril(A,-1); %求A的下三角阵U=-triu(A,1); %求A的上三角阵G=(D-L)\U;f=(D-L)\b;y=G*x0+f;n=1; %迭代次数while norm(y-x0)>=epsx0=y;y=G*x0+f;n=n+1;end（3）实验结果及分析:>> A=[8,-3,2;4,11,-1;6,3,12];>> b=[20,33,36]';>> [x,n]=gauseidel(A,b,[0,0,0]',1.0e-6)x =3.00002.00001.0000n =9七、结果说明：高斯-赛德尔迭代法比雅可比迭代法收敛得快一些（达到相同精度所需迭代次数较少）。
Hilbert矩阵
希尔伯特矩阵的注记（2009-11-21）一、希尔伯特矩阵的来源及矩阵的正定性设f (x )∈C [0，1]，求n 次多项式P (x ) = a 0 + a 1x + a 2 x 2 + …… + a n x n使得⎰-=12)]()([dx x f x P L取最小值。
其中，P (x )称为函数f (x )的最佳平方逼近多项式。
由极值必要条件，对⎰∑=-=1210])([),,,(dx x f x a a a a L nj j j n中的变量求导数，得⎰∑⎰⎰∑=+=-=-=∂∂1110)(22])([2nj iji j n j jj iidx x f x dx xa dx x f xa x a L令其为零，得方程组∑⎰===++n j ij n i dxx f x a j i 01),,1,0()(11令⎰=1)(dx x f x b ii ，将方程组写为矩阵形式⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡+++++n n b b b a a a n n n n n1010)12/(1)2/(1)1/(1)2/(13/12/1)1/(12/11 方程组的系数矩阵就是著名的Hilbert 矩阵。
多项式可表示为内积形式⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡=n n x x a a a x P 1][)(10 所以⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡=n nn n a a a x x x x a a a x P 10102]1[1][)]([积分，得⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡+++++=⎰n n a a a n n n n n a a a dx x P101012)12/(1)2/(1)1/(1)2/(13/12/1)1/(12/11][)]([由此可知，希尔伯特矩阵是对称正定矩阵。
二、希尔伯特矩阵的病态性希尔伯特矩阵是著名的病态矩阵。
工程数值分析实验报告(3篇)
第1篇一、实验目的本次实验旨在通过数值分析的方法，对工程实际问题进行建模、求解和分析。
通过学习数值方法的基本原理和算法，提高解决实际工程问题的能力。
二、实验内容1. 线性方程组的求解2. 矩阵特征值与特征向量的计算3. 函数插值与曲线拟合4. 数值微分与积分三、实验步骤1. 线性方程组的求解（1）编写程序实现高斯消元法、克劳斯消元法和列主元素法（2）设计输入界面，用户输入增广矩阵的行和列，填写系数及常数项（3）分别运用三种方法求解线性方程组，比较求解结果的正确性、数值稳定性和计算效率2. 矩阵特征值与特征向量的计算（1）编写程序实现幂法、QR算法和逆幂法（2）设计输入界面，用户输入矩阵的行和列，填写矩阵元素（3）分别运用三种方法计算矩阵的特征值与特征向量，比较求解结果的准确性和计算效率3. 函数插值与曲线拟合（1）编写程序实现拉格朗日插值、牛顿插值和样条插值（2）设计输入界面，用户输入函数的自变量和函数值，选择插值方法（3）分别运用三种方法进行函数插值，比较插值结果的准确性和光滑性4. 数值微分与积分（1）编写程序实现有限差分法、龙格-库塔法和辛普森法（2）设计输入界面，用户输入函数的导数或积分的上下限，选择数值方法（3）分别运用三种方法进行数值微分和积分，比较求解结果的准确性和计算效率四、实验结果与分析1. 线性方程组的求解通过实验，我们发现列主元素法在求解线性方程组时具有较好的数值稳定性，计算效率也较高。
而高斯消元法和克劳斯消元法在处理大型稀疏矩阵时存在一定的困难。
2. 矩阵特征值与特征向量的计算实验结果表明，QR算法和逆幂法在计算矩阵特征值与特征向量时具有较高的准确性和计算效率。
幂法在处理大型稀疏矩阵时表现出较好的性能。
3. 函数插值与曲线拟合在函数插值和曲线拟合实验中，样条插值方法具有较好的准确性和光滑性。
拉格朗日插值和牛顿插值方法在处理简单函数时表现良好，但在处理复杂函数时可能存在精度问题。
《数值分析》课程实验报告范文
《数值分析》课程实验报告范文《数值分析》课程实验报告姓名：学号：学院：机电学院日期：2022年某月某日目录实验一函数插值方法1实验二函数逼近与曲线拟合5实验三数值积分与数值微分7实验四线方程组的直接解法9实验五解线性方程组的迭代法15实验六非线性方程求根19实验七矩阵特征值问题计算21实验八常微分方程初值问题数值解法24实验一函数插值方法一、问题提出对于给定的一元函数的n+1个节点值。
试用Lagrange公式求其插值多项式或分段二次Lagrange插值多项式。
实验二函数逼近与曲线拟合一、问题提出从随机的数据中找出其规律性，给出其近似表达式的问题，在生产实践和科学实验中大量存在，通常利用数据的最小二乘法求得拟合曲线。
在某冶炼过程中，根据统计数据的含碳量与时间关系，试求含碳量与时间t的拟合曲线。
t(分)051015202530354045505501.272.162.863.443.874.154.374.51 4.584.024.64二、要求1、用最小二乘法进行曲线拟合；2、近似解析表达式为；3、打印出拟合函数，并打印出与的误差，；4、另外选取一个近似表达式，尝试拟合效果的比较；5、某绘制出曲线拟合图。
三、目的和意义1、掌握曲线拟合的最小二乘法；2、最小二乘法亦可用于解超定线代数方程组；3、探索拟合函数的选择与拟合精度间的关系四、实验步骤：第一步先写出线性最小二乘法的M文件functionc=lpoly(某,y,m)n=length(某);b=zero(1:m+1);f=zero(n,m+1); fork=1:m+1f(:,k)=某.^(k-1);enda=f＇某f;b=f＇某y＇;c=a＼b;c=flipud(c);第二步在命令窗口输入：>>lpoly([0,5,10,15,20,25,30,35,40,45,50,55],[0,1.27,2.16,2.86,3.44,3.87,4.15,4.37,4.51,4.58,4.02,4.64],2)回车得到：an=-0.00240.20370.2305即所求的拟合曲线为y=-0.0024某2+0.2037某+0.2305在编辑窗口输入如下命令：>>某=[0,5,10,15,20,25,30,35,40,45,50,55];>>y=-0.0024某某.^2+0.2037某某+0.2305;>>plot(某,y)命令执行得到如下图五、实验结论分析复杂实验数据时,常采用分段曲线拟合方法。
数值分析希尔伯特病态线性方程组
病态线性代数方程组的求解理论的分析表明，求解病态的线性代数方程组是困难的。
考虑方程组Hx = b 的求解，其中H 为Hilbert 矩阵，n n ij h H ⨯=)(，11-+=j i h ij ，n j i ,...,2,1,=1. 估计Hilbert 矩阵2-条件数与阶数的关系；2. 选择问题的不同维数，分别用Gauss 消去法，Jacobi 迭代，GS 迭代和SOR 迭代求解，比较结果；3. 讨论病态问题求解的算法。
解： 1、取Hilbert 矩阵阶数最高分别为n=20和n=100。
采用Hilbert 矩阵的2-条件数作为实验的比较对象，画出的曲线如下图所示：lg(())n cond H n从图中可以看出，在n ≤13之前，图像接近直线，在n ＞13之后，图像趋于平缓，在一定的范围内上下波动。
为了比较图像的线性部分，作出一条直线与已知曲线进行比较。
比较直线的关系式为：833.1519.1))(lg(-=n H cond n ，结果下图所示。
nl g (c o n d (H n ))lg(cond(Hn))~n 关系图从图2中可以看出，当n 较小时，n H cond n ~))(lg(之间近似满足线性关系。
当n 继续增大到100时，n H cond n ~))(lg(关系图下图所示：从图中可以看出，图像的走势符合在n=20时的猜想，在n 大于一定的值之后，图像趋于平缓，且在一定范围内震荡，同时又有一定上升趋势，但上升速度很慢。
2、选择不同的阶数n ，设方程组的精确解为xz=(1,1,…,1)T进行计算，用四种方法解x_Guass1、x_Jacobi1、x_GS1、x_SOR1对比表如下nl g (c o n d (H n ))lg(cond(Hn))~n关系图nl g (c o n d (H n ))lg(cond(Hn))~n 关系图表所示。
Gauss 消去法求解：选择问题的阶数为3~8时，用Gauss 消去法求得的解与精确解一致，当阶数为9~14时，解开始出现偏差，而且n 越大，偏差越大。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
数值分析课程实验报告
题目：病态线性方程组的求解
理论分析表明，数值求解病态线性方程组很困难。
考虑求解如下的线性方程组的求解
Hx = b ，期中H 是Hilbert 矩阵，()ij n n H h ⨯=，1
1
ij h i j =+-，i ，j = 1,2，…，n
1. 估计矩阵的2条件数和阶数的关系
2. 对不同的n ，取(1,1,
,1)n
x =∈
，分别用Gauss 消去，Jacobi 迭代，Gauss-seidel 迭
代，SOR 迭代和共轭梯度法求解，比较结果。
3. 结合计算结果，试讨论病态线性方程组的求解。
解答过程
1.估计矩阵的2-条件数和阶数的关系
矩阵的2-条件数定义为：1
222
()Cond A A A
-=⨯，将Hilbert 矩阵带入有：
1222
()Cond H H H -=⨯
调用自编的Hilbert_Cond 函数对其进行计算，取阶数n = 50，可得从1阶到50阶的2-条件数，以五位有效数字输出，其中前10项见表1。
表1.前十阶Hilbert 矩阵的2-条件数
从表1可以看出，随着阶数每递增1，Hilbert 矩阵的2-条件数都至少增加一个数量级，但难以观察出明显的相依规律。
故考虑将这些数据点绘制在以n 为横轴、Cond （H ）2为纵轴的对数坐标系中（编程用Hilbert_Cond 函数同时完成了这个功能），生成结果如图1。
图1.不同阶数下Hilbert矩阵的2-条件数分布
由图可见，当维数较小时，在y-对数坐标系中Cond（H）2与n有良好的线性关系；但n超过10后，线性趋势开始波动，n超过14后更是几乎一直趋于平稳。
事实上，从n = 12开始，系统便已经开始提出警告：“Warning: Matrix is close to singular or badly scaled.Results may be inaccurate.”。
也就是说，当n较大时，H矩阵已经接近奇异，计算结果可能是不准确的。
通过查阅相关资料，我找到了造成这种现象的原因：在matlab中，用inv函数求条件数过大的矩阵的逆矩阵将是不可靠的。
而调用系统自带的专门对Hilbert矩阵求逆的invhilb（n）函数则不存在这个问题，生成结果如图2。
图2. 修正后的不同阶数下Hilbert矩阵的2-条件数分布
简便起见，取n 不大于10的前十项进行线性拟合，结果如图3。
C o n d (H )
n
图3.1~10阶Hilbert 矩阵2-条件数的线性拟合
由拟合结果知，Cond （H ）2与n 的关系为：
2lg () 1.65183 1.47863Cond H n
=-+
其线性相关系数r=0.99985，可见二者具有较好的线性关系。
对上式稍作变形得：
1.651831.478632()10n Cond H -+=
这就是对Hilbert 矩阵的-2条件数与其阶数n 的关系估计。
可见Hilbert 矩阵的2-条件数会随其阶数n 的增加呈指数增大趋势，因此当n 较大时Hilbert 矩阵将是严重病态的，甚至导致matlab 中inv 求逆运算失真。
2.对不同的n ，采用各种方法求解方程
调用自编的Calc 主函数（其中包括的Hilbert 函数以及b_函数可创建出对应阶数的H 矩阵以及向量b ，Gauss_Cal 函数、Jacobi_Cal 函数、Gauss_Seidel_Cal 函数、SOR_Cal 函数（该函数自动寻找最优松弛因子，然后以最优因子进行求解）以及CG_Cal 函数则可完成各
自方法的求解），分别取n = 2,5,10,20,50，对于迭代法设定终止计算精度为10
10ε-=，所得
计算结果以16位有效数字输出，分别见表2~表6。
表2.n=2的计算结果
从表2可看出，n=2时，四种迭代法都能够收敛，迭代次数最大为e+2量级（J法），最小仅要2次（CG法），并且五种解法都能给出非常精确的结果，最大误差为e-10量级（GS 法）。
表3.n=5的计算结果
从表3可以看出，n=5时，J法已经不收敛，其余解法依然收敛（但值得注意的是GS 法以及SOR法的谱半径以及非常接近1，达到了收敛的边缘），最大迭代次数达到e+6量级（GS法），最小需要7次（CG法），计算精度依然较高，最大误差为e-6量级。
SOR法比
GS法节省了较多的迭代次数。
表4.n=10的计算结果
从表4可以看出，n=10时，各种解法的收敛性与n=5时相同（GS法与SOR法的谱半径进一步趋近1），最大迭代次数达e+8量级（SOR法），最小需要8次（CG法），计算精度已经较低，最大误差达到e-3量级。
此时有一个异常现象：SOR法的迭代次数不但不比GS法少，反而超出好几倍。
但根据计算出的两种算法下的迭代矩阵谱半径，可以看出SOR 法为0.999999999871931，而GS法为0.999999999997045，在最优松弛因子之下的SOR法确实具有更小的迭代矩阵谱半径，因此应当更快收敛。
经检查，计算程序的编制应该没有错误，但计算结果确实与理论不符，希望老师指点迷津！
表5.n=20的计算结果
从表5可以看出，n=20时，各种解法的收敛性依然没变（但GS法以及SOR法的谱半径在matlab的十六位显示精度下已经等于1），最大迭代次数e+8量级（SOR法），最小需要10次（CG法），Gauss消元法的计算结果已失去意义（误差e+2量级），迭代法的计算误差均为e-3量级。
且此时SOR的迭代次数依然高于GS。
表6.n=50的计算结果
从表6可以看出，n=50时，计算结果的特点与n=20相似，不同之处在于Gauss消元法的误差进一步增大（e+3量级）。
3.综合讨论
根据上一部分的计算结果，可以总结出求解病态矩阵的两个特点：
1.收敛性差，收敛速度慢。
对于本例，当阶数n大于2以后J法就不收敛了，这是收敛性差的一个体现。
GS法和SOR法虽然一直保持收敛（事实上这是由于Hilbert矩阵是正定对称的，所以决定了GS法以及SOR法一定收敛），但随着n的增大它们的谱半径迅速趋近于1（n=20时matlab的16位显示精度已经无法显示出它们的谱半径与1的差别），根据理论我们知道，趋近于1的谱半径决定了计算的收敛速度将会非常慢，之前计算结果中的迭代次数N也很好地验证了这一点。
2.计算精度低，阶数较高时误差惊人。
根据理论，线性方程组系数矩阵的条件数的直观意义就是初始数据扰动的相对误差的放大倍数。
由于数值计算时的舍入误差等因素，初始数据的微小扰动难以避免，因此解病态方程组时，这种微小的初始相对误差将被病态矩阵巨大的条件数所放大，最后造成计算结果的巨大偏差。
对于本例，n=20时Gauss消元法的计算误差竟然达到e+2量级，使得计算结果完全失去意义，而其他几种迭代法的计算误差看似能够接受（e-3量级），但仔细分析上一节
的计算结果不难发现：从n=2到n=5再到n=10，几种迭代法的迭代次数迅速增大，且计算精度显著降低；而从n=10到n=20再到n=50，迭代次数和计算精度却趋于稳定，由于矩阵的条件数一直呈指数增大，因此这种计算结果的趋势明显是不合理的。
究其原因，我认为和一开始图1所示的2-条件数的计算失真是类似的。
当矩阵的条件数非常大时，inv指令的求逆结果不可靠，从而造成图1中2-条件数趋于平稳的那个阶段，而在GS迭代以及SOR迭代中也需要对矩阵进行求逆运算，因此很有可能又是因为inv对高条件数矩阵求逆的偏差造成最终计算结果的失真。
因此，我预计实际上当n=20时，GS法和SOR法的迭代次数以及计算误差要远远高于n=10的值，同理n=50时的结果又要远远大于n=20的值。
不难看出，造成上述两个特点的本质因素还是因为病态矩阵的高条件数，无论收敛性、收敛速度还是计算精度等，都与矩阵的条件数直接相关，因此改善病态矩阵的求解的根本方法还是在于降低矩阵的条件数。
因此我尝试着对Hilbert矩阵进行了一些预处理，并发现在某些预处理方法下条件数的阶数被有效减小，因此预处理是求解病态方程组的一个十分有效且必要的环节。
此外，若给定条件数并横向比较上一节中的不同求解方法，可以看出CG法的表现是最令人满意的，当矩阵阶数n不太大时（不超过10），CG法都能以非常少的迭代次数（只需要几次）以及较高的计算精度完成计算任务，相比之下，GS法以及SOR法的上千万次迭代次数显得非常吃力。
此外，事实上n较小时Gauss消元法的计算结果也比较理想，具有很快的计算速度和较高的计算精度。
我们知道，在传统意义上Gauss消元法的最大弊端在于n非常大时计算量大得令人难以接受（n三次方量级），而迭代法的单步计算量则为n 平方量级，但在求解病态问题时，n只能比较小，且迭代法的迭代次数却非常大，所以造成迭代法的计算速度反而远慢于Gauss消元法。
综上，求解病态方程组时首先应该对系数矩阵进行预处理，尽可能降低它的条件数；其次可以选择共轭梯度（CG）法对其进行求解。
如果阶数n不大，Gauss消元法也是比较有效的。