Jacobi G-S SOR迭代法在matlab中的例子

4 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
x0=zeros(20,1),运行程序，得到i=38。 Jacobi迭代法绘图程序： function y=jacobian1(A,b,x0) residue=zeros(1000,1); if size(b,2)~=1 b=b'; end if size(x0,2)~=1 x0=x0'; end D=diag(diag(A)); B=D\(D-A); g=D\b; x1=B*x0+g; i=1; while(norm(x1-x0,inf)>1e-10) residue(i)=norm(x1-x0,inf); x0=x1; x1=B*x1+g; i=i+1; end semilogy(residue);
（k）
已经算出，用它代替 x1（k-1） ，其他分量仍用 xi（k-1） 。类似
的，计算 xl（k）时，因 x1（k） ， ， ，xl-1（k）都已算出，用它们代替 x1（k-1） ， ， ， ，xl-1（k-1） ，其他分量仍用 xk-1 的分量，于是有 Xk=D-1Lxk+D-1Uxk-1+g， k=1，2， ， ， （6）
10 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
end D=diag(diag(A)); L=-tril(A-D); U=-triu(A-D); M=(D-w*L)\((1-w)*D+w*U); g=w*(D-w*L)\b; x1=M*x0+g; i=1; while(norm((x1-x0),2)>1.0e-6) residue(i)=norm(x1-x0,inf); x0=x1; x1=M*x0+g; i=i+1; end semilogy(residue); y=x1; i 由公式 wopt=
3 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
if size(b,2)~=1 b=b'; end if size(x0,2)~=1 x0=x0'; end D=diag(diag(A)); B=D\(D-A); g=D\b; x1=B*x0+g; i=1; while(norm((x1-x0),2)>0.0001) x0=x1; x1=B*x1+g; i=i+1; end y=x0; i-1 其中，A是一个20阶矩阵。由于Jacobi迭代法需满足条件A及 2D-A均正定，故可取A=rand（20），令A=A*A’，然后令 A=A+100*eye （20） ， 使A满足严格对角占优。 算得cond （A ， inf）=3.1223，说明该矩阵条件数较好。令b=ones(20,1)，
B <1 或者 B
1
2、 设 B 是 jacobi 迭代的迭代矩阵。 若 则 G-S 迭代收敛。
<1，
3、 若系数矩阵 A 对称， 对角线元素 aii>0 （i=1， 2， ， ， n） ， 则 jacobi 迭代收敛的充分必要条件是 A 和 2D-A 都正定；若 A 正定，则 G-S 迭代收敛。 三、MATLAB 中的矩阵实验 Jacobi 矩阵算法： function y=jacobian(A,b,x0)
[数学基地班 赵晨晓 2011301000007]
可以看出图像弯曲不平，性质没有上一条好。
可以看出， 只需要 6 次迭代， 就能满足要求。 因为 w 较小时， 收敛稳定性高，w 较大时，收敛速度较快。 另一方面，当矩阵条件数没有那么好的时候，例如可令 A=rand（20） ，A=A*A’， cond(A,inf)= 2.8346e+04，此时， 由于 jacobi 迭代条件较为苛刻，jacobi 迭代图像如下所示：
5 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
y=x0; i 可 得 到 此 矩 阵 的 收 敛 图 像 ：
G-S迭代法算法： function y=GS(A,b,x0) if size(b,2)~=1 b=b'; end if size(x0,2)~=1 x0=x0'; end
9 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
L=-tril(A-D); U=-triu(A-D); M=(D-w*L)\((1-w)*D+w*U); g=w*(D-w*L)\b; x1=M*x0+g; i=1; while(norm((x1-x0),2)>1.0e-6) residue(i)=norm(x1-x0,inf); x0=x1; x1=M*x0+g; i=i+1; end y=x1; i SOR 迭代法绘图算法： function y=SOR(A,b,x0,w) residue=zeros(1000,1); if size(b,2)~=1 b=b'; end if size(x0,2)~=1 x0=x0';
13 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
这说明此时 jacobi 迭代法不收敛； 而 G-S 迭代法图像如下所示：
14 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
这说明它需要迭代 3000 多次才能达到预期目标。 而对于 SOR 迭代法，由于ρ（B）过大，它不收敛。 综上，可以得到结论： （1） 当矩阵收敛时， G-S 迭代法收敛速度总比 Jacobi 迭代法收敛速度要快。而选取恰当 w 后，SOR 迭代法速度比 G-S 迭代法收敛速度 还要快。 （2） 当矩阵条件数较小时，即此时矩阵性质比较 好，三种算法的收敛速度都比较快。 （3） 以上启示我们在进行矩阵迭代之前，最好对 其进行预处理，以使其条件数较小，性质优 良。
6 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
D=diag(diag(A)); L=-tril(A-D); U=-triu(A-D); M=(D-L)\U; g=(D-L)\b; x1=M*x0+g; i=1; while(norm((x1-x0),2)>1.0e-6) x0=x1; x1=M*x0+g; i=i+1; end y=x1; i-1 同样的 A，得到 i=8。 C-G 迭代法绘图算法： function y=GS(A,b,x0) residue=zeros(1000,1); if size(b,2)~=1 b=b'; end if size(x0,2)~=1
1 1 1 ( B)* ( B)
利用函数 vrho，代入计算，得到
wopt=0.9，此时收敛图像为：
11 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
可以看出这是一条较为光滑且下降较快的曲线。作为对比， 当 w=1.5 时，其图像为：
12 / 15
[数值分析报告]
其中，D 为严格下三角阵，L 为严格上三角阵，那个（1）可 写为 x=Bx+g （4）
其中，B=D-1(L+U),g=D-1b.若给定初始向量 X0=（x1（0） ， x2 （ 0 ） ， ， ， ， x n（ 0 ） ）T，并代入（4）的右端，就 可以计算出一个新的向量 x1，即 x1=Bx0+g；
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
关于 Jacobi G-S SOR 方法的收敛速度比较 摘要： 本论文主要通过判断比较 Jacobi G-S SOR 三种
计算方法的收敛速度，来分析三种算法的优劣。辅助软件为 MATLAB。 随着计算技术的发展，计算机的存储量日益增大，计算速 度也迅速提高，直接法在计算机上可以求解的线性方程组的 规模也越来越大， 但直接法大多需要对系数矩阵 A 进行分解， 因此一般不能保证 A 的稀疏性。而实际应用中，特别是偏微 分方程的数值求解时，常常遇到的就是大型稀疏线性方程的 求解问题。因此寻求能够保持稀疏性的有效算法就成为数值 线性代数中的一个非常重要的研究课题。 一、三种迭代方法内容简介 Jacobi 迭代 考虑非奇异线性代数方程组 Ax=bHale Waihona Puke 令 A=D-L-U （1） （2 ）
1 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
再把 x1 代入（4）右端，又可得到一个向量 x2；以此类推有 Xk=Bxk-1+g，k=1，2，3， ， ， （5 ）
这被称为 jacobi 迭代格式。B 称为 jacobi 迭代的迭代矩阵，g 称为常数项。 Gauss-Seidel 迭代 假设不按 jacobi 格式，而是在计算 xk 的第一个分量用 xk-1 的各个分量计算。但当计算 xk 的第二个分量 x2（k）时，因 x1
k+1
-xk，则有
xk+1= x+xk。
Xk+1 可以看作是在向量 xk 上加上修正项 x 而得到的。 若
2 / 15
[数值分析报告]
[数学基地班 赵晨晓 2011301000007]
修正项的前面加上一个参数ω，便得到松弛法的计算公式 Xk+1=xk+ω x=（1-ω）xk+ω（D-1Lxk+1+D-1Uxk+D-1b） 其中ω叫做松弛因子。当ω>1 时叫做超松弛；ω<1 时叫 做低松弛； ω=1 即为 G-S 迭代。 我们把超松弛迭代简称为 SOR。
可得到此矩阵的收敛图像：
8 / 15
[数值分析报告]