.176.第六讲线性方程组定常选代方法不断重复以上步骤，于是，我们就得到一个序列(1),(2满足以下递推关系(k+1) = 2(k) + M-1(b - Ar(k)), k = 1, 2, ...这就构成了一个选代方法，由于每次选代的格式是一样的，因此称为定常选代通常，构造一个好的定常选代，需要考虑以下两点：(1）以M为系数矩阵的线性方程组必须要比原线性方程组更容易求解；(2）M应该是A的一个很好的近似，且送代序列【]要收敛常用的定常选代方法是基于矩阵分裂的迭代方法，主要包括：·Jacobi方法·Gauss-Seidel (G-S)方法·超松弛（SOR，SuccessiveOver-Relaxation）方法·对称超松弛（SSOR，SymmetricSOR）方法·加速超松弛（AOR，AcceleratedOver-Relaxation）方法6.2.1矩阵分裂送代方法首先给出矩阵分裂的定义定义6.1（矩阵分裂MatrixSplitting）设AERnxn非奇异，称(6.6)A-M-N为A的一个矩阵分裂,其中 M 非奇异考虑线性方程组Ar =b,(6.7)其中AERnxn非奇异.选代方法的基本思想：给定一个选代初始值r(O)，通过一定的选代格式生成一个选代序列(())=0使得lim 2(K) = ,会 A-1b.给定一个矩阵分裂（6.6)，则原方程组（6.7)就等价于Mr=N+b.于是我们就可以构造出以下的选代格式Mr(+1)= Nr(k) +b, k=0,1,...或(k+1) = M-1N(k) + M-1b ≤ Ga(k) + g, k = 0, 1,..,(6.8)其中G=M-1N称为送代矩阵.这就是基于矩阵分裂（6.6）的选代方法.易知，选取不同的M，就可以构造出不同的选代方法
· 176 · 第六讲 线性方程组定常迭代方法 不断重复以上步骤, 于是, 我们就得到一个序列 x (1), x(2), . . . , , x(k) , . . . . 满足以下递推关系 x (k+1) = x (k) + M−1 (b − Ax(k) ), k = 1, 2, . . . . 这就构成了一个迭代方法. 由于每次迭代的格式是一样的, 因此称为 定常迭代. 通常, 构造一个好的定常迭代, 需要考虑以下两点: (1) 以 M 为系数矩阵的线性方程组必须要比原线性方程组更容易求解; (2) M 应该是 A 的一个很好的近似, 且迭代序列 {xk} 要收敛. 常用的定常迭代方法是基于矩阵分裂的迭代方法, 主要包括: • Jacobi 方法 • Gauss-Seidel (G-S) 方法 • 超松弛 (SOR, Successive Over-Relaxation) 方法 • 对称超松弛 (SSOR, Symmetric SOR) 方法 • 加速超松弛 (AOR, Accelerated Over-Relaxation) 方法 6.2.1 矩阵分裂迭代方法 首先给出 矩阵分裂 的定义. 定义 6.1 (矩阵分裂 Matrix Splitting) 设 A ∈ R n×n 非奇异, 称 A = M − N (6.6) 为 A 的一个矩阵分裂, 其中 M 非奇异. 考虑线性方程组 Ax = b, (6.7) 其中 A ∈ R n×n 非奇异. 迭代方法的基本思想: 给定一个迭代初始值 x (0) , 通过一定的迭代格式生成一 个迭代序列 {x (k)}∞ k=0, 使得 lim k→∞ x (k) = x∗ ≜ A −1 b. 给定一个矩阵分裂 (6.6), 则原方程组 (6.7) 就等价于 Mx = Nx + b. 于是我们就可以构造出以下的迭代 格式 Mx(k+1) = Nx(k) + b, k = 0, 1, . . . , 或 x (k+1) = M−1Nx(k) + M−1 b ≜ Gx(k) + g, k = 0, 1, . . . , (6.8) 其中 G = M−1N 称为 迭代矩阵. 这就是基于矩阵分裂 (6.6) 的迭代方法. 易知, 选取不同的 M, 就可以 构造出不同的迭代方法
6.2定常选代方法.177.6.2.2Jacobi送代将矩阵A分裂为A=D-L-U.其中D为A的对角线部分，-L和-U分别为A的严格下三角和严格上三角部分在矩阵分裂A=M-N中取M=D,N=L+U,则可得Jacobi送代方法(k+1) = D-1(L + U)(k) + D-1b, k = 0, 1,2, ....(6.9)对应的迭代矩阵为G) = D-1(L + U).写成分量形式即为T1(k)(k+1)i=1,2,...,n.aijLiaiij=l.j+i+由于Jacobi 选代中α(k+1)的更新顺序与i无关,即可以按顺序i=1,2,,n计算,也可以按顺序i=n,n.一1,·.,2.1计算，或者乱序计算.因此Jacobi选代非常适合并行计算算法6.1.Jacobi选代1; Given an initial guess r(0)2:while not converge dofor i=l to n do3:_(k+1)4:Ciend for5:6: end while我们有时也将Jacobi迭代格式写为(k+1) = 2(k) + D-1(b - Ar(k) = 2(k) + D-1rk, k = 0, 1, 2, .*-其中rk会b一Az(）是次选代后的残量.这表明，(k+1)是通过对(k）做一个修正得到的下面给出求解二维离散Poisson方程（6.5）的Jacobi选代方法算法6.2.求解二维离散 Poisson方程的 Jacobi 选代方法I: Given an initial guess w(0)2:while not convergedofor i=1 to N do3:for j=l to N do4:+)=(f +1 +1++++1) /45:ui.s6:end for
6.2 定常迭代方法 · 177 · 6.2.2 Jacobi 迭代 将矩阵 A 分裂为 A = D − L − U, 其中 D 为 A 的对角线部分, −L 和 −U 分别为 A 的严格下三角和严格上三角部分. 在矩阵分裂 A = M − N 中取 M = D, N = L + U, 则可得 Jacobi 迭代 方法: x (k+1) = D−1 (L + U)x (k) + D−1 b, k = 0, 1, 2, . . . . (6.9) 对应的迭代矩阵为 GJ = D−1 (L + U). 写成分量形式即为 x (k+1) i = 1 aii bi − ∑n j=1,j̸=i aijx (k) j , i = 1, 2, . . . , n. † 由于 Jacobi 迭代中 x (k+1) i 的更新顺序与 i 无关, 即可以按顺序 i = 1, 2, . . . , n 计算, 也可以按顺序 i = n, n − 1, . . . , 2, 1 计算, 或者乱序计算. 因此 Jacobi 迭代非常适合并行计算. 算法 6.1. Jacobi 迭代 1: Given an initial guess x (0) 2: while not converge do 3: for i = 1 to n do 4: x (k+1) i = ( bi − ∑n j=1,j̸=i aijx (k) j ) / aii 5: end for 6: end while 我们有时也将 Jacobi 迭代格式写为 x (k+1) = x (k) + D−1 (b − Ax(k) ) = x (k) + D−1 rk, k = 0, 1, 2, . . . , 其中 rk ≜ b − Ax(k) 是 k 次迭代后的残量. 这表明, x (k+1) 是通过对 x (k) 做一个修正得到的. 下面给出求解二维离散 Poisson 方程 (6.5) 的 Jacobi 迭代方法. 算法 6.2. 求解二维离散 Poisson 方程的 Jacobi 迭代方法 1: Given an initial guess v (0) 2: while not converge do 3: for i = 1 to N do 4: for j = 1 to N do 5: u (k+1) i,j = ( h 2fi,j + u (k) i+1,j + u (k) i−1,j + u (k) i,j+1 + u (k) i,j−1 ) /4 6: end for
第六讲线性方程组定常选代方法.178.end for7:8:end while6.2.3Gauss-Seidel选代在分裂A=M-N中取M=D-L,N=U,即可得Gauss-Seidel(G-S)送代方法(k+1) = (D - L)-1U(K) + (D - L)-1b)(6.10)对应的选代矩阵为GGs = (D - L)-1U将G-S选代改写为Dr(h+1) = Lr(k+1) + Ur() + b,即可得分量形式i-1n1(k+1)(k+1)(kV7=1.2..b:aijrQijs.naiij=1j=i+1算法6.3.Gauss-Seidel选代1: Given an initial guess r(0)2:while not converge dofori=ltondo3:A(k+1)(k+1)4:Liai=+endfor5:6:end while与Jacobi选代方法类似，我们也可以给出求解问题（6.5)的G-S选代方法G-S方法的主要优点是充分利用了已经获得的最新数据但在G-S方法中，未知量的更新必须按自然顺序进行，因此不适合并行计算下面我们介绍一种适合并行计算的更新顺序：红黑排序，即将二维网格点依次做红黑记号，如右图所示在计算过程中，对未知量的值进行更新时，我们可以先更新红色节点，此时所使用的只是黑色节点的数据，然后再更新黑色节点，这时使用的是红色节点的数据.于是我们得到下面红黑排序G-S选代方法由于在更新红点时，各个点之间是相互独立的，因此可以并行计算，同样，在更新黑点时，各个点之间也是相互独立的，因此也可以并行计算
· 178 · 第六讲 线性方程组定常迭代方法 7: end for 8: end while 6.2.3 Gauss-Seidel 迭代 在分裂 A = M − N 中取 M = D − L, N = U, 即可得 Gauss-Seidel (G-S) 迭代 方法: x (k+1) = (D − L) −1Ux(k) + (D − L) −1 b. (6.10) 对应的迭代矩阵为 GGS = (D − L) −1U. 将 G-S 迭代改写为 Dx(k+1) = Lx(k+1) + Ux(k) + b, 即可得分量形式 x (k+1) i = 1 aii bi − ∑ i−1 j=1 aijx (k+1) j − ∑n j=i+1 aijx (k) j , i = 1, 2, . . . , n. 算法 6.3. Gauss-Seidel 迭代 1: Given an initial guess x (0) 2: while not converge do 3: for i = 1 to n do 4: x (k+1) i = 1 aii ( bi − i ∑−1 j=1 aijx (k+1) j − ∑n j=i+1 aijx (k) j ) 5: end for 6: end while 与 Jacobi 迭代方法类似, 我们也可以给出求解问题 (6.5) 的 G-S 迭代方法. G-S 方法的主要优点是充分利用了已经获得的最新数据. 但在 G-S 方法中, 未知量的更新必须按自然顺序进行, 因此不适合并行计算. 下面我们介绍一种适合并行计算的更新顺序: 红黑排序, 即将 二维网格点依次做红黑记号, 如右图所示. 在计算过程中, 对未知量的值进行更新时, 我们可以先更新红色 节点, 此时所使用的只是黑色节点的数据, 然后再更新黑色节点, 这 时使用的是红色节点的数据. 于是我们得到下面红黑排序 G-S 迭代 方法. 由于在更新红点时, 各个点之间是相互独立的, 因此可以并行 计算. 同样, 在更新黑点时, 各个点之间也是相互独立的, 因此也可 以并行计算
6.2定常选代方法.179 .算法6.4.求解二维离散Poisson方程的红黑排序G-S选代方法1: Given an initial guess w(0)2:while not converge dofor（i，j）为红色节点do3:(k+1)(k)()(k)(h2fi + u+1, + u-1, + u3+1 + u3-1)4:ui.j1end for5:for(i,i）为黑色节点do6:ug+1)(k+1)+ u(*+1)(k+1) + u(k+1)(h2fij +uif)+u+7:-114end for8:9:end while6.2.4SOR送代在G-S方法的基础上，我们可以通过引入一个松弛参数w来加快收敛速度.这就是SOR(SuccessiveOverrelaxation)方法[78].该方法的基本思想是将G-S方法中的第k+1步近似解与第k步近似解做一个加权平均，从而给出一个新的近似解，即(k+1) = (1 - w)() + w (D-1(L(+1) + U(k) + D-1b)(6.11)整理后即为(k+1) = (D - wL)-1 (1 - w)D + wU) r(k) + w(D - wL)-lb,(6.12)其中w称为松弛参数（relaxationparameter)，当w=1时，SOR即为G-S方法，当w<1时，称为低松弛(underrelaxation）方法，当w>1时，称为超松弛（overrelaxation方法.在大多数情况下，当w>1时会取得比较好的收敛效果十SOR方法曾经在很长一段时间内是科学计算中求解线性方程组的首选方法SOR的选代矩阵为GsOR = (D - wL)-1 (1 - w)D +wU),对应的矩阵分裂为1-wp+U.M=N-D-Lww由（6.12)可知SOR迭代的分量形式为i-1nUr(++1) = (1 - w)2(k)(k+1)KbC;aiaiij=1i=i+i-12wr(k)(k+1)kVb,-aiiaijraiij=1j=i
6.2 定常迭代方法 · 179 · 算法 6.4. 求解二维离散 Poisson 方程的红黑排序 G-S 迭代方法 1: Given an initial guess v (0) 2: while not converge do 3: for (i, j) 为红色节点 do 4: u (k+1) i,j = 1 4 ( h 2fi,j + u (k) i+1,j + u (k) i−1,j + u (k) i,j+1 + u (k) i,j−1 ) 5: end for 6: for (i, j) 为黑色节点 do 7: u (k+1) i,j = 1 4 ( h 2fi,j + u (k+1) i+1,j + u (k+1) i−1,j + u (k+1) i,j+1 + u (k+1) i,j−1 ) 8: end for 9: end while 6.2.4 SOR 迭代 在 G-S 方法的基础上, 我们可以通过引入一个松弛参数 ω 来加快收敛速度. 这就是 SOR (Successive Overrelaxation) 方法 [78]. 该方法的基本思想是将 G-S 方法中的第 k + 1 步近似解与第 k 步近似解做一 个加权平均, 从而给出一个新的近似解, 即 x (k+1) = (1 − ω)x (k) + ω ( D−1 (Lx(k+1) + Ux(k) ) + D−1 b ) . (6.11) 整理后即为 x (k+1) = (D − ωL) −1 ((1 − ω)D + ωU) x (k) + ω(D − ωL) −1 b, (6.12) 其中 ω 称为松弛参数 (relaxation parameter) . 当 ω = 1 时, SOR 即为 G-S 方法, 当 ω < 1 时, 称为低松弛 (under relaxation) 方法, 当 ω > 1 时, 称为超松弛 (over relaxation) 方法. 在大多数情况下, 当 ω > 1 时会取 得比较好的收敛效果. † SOR 方法曾经在很长一段时间内是科学计算中求解线性方程组的首选方法. SOR 的迭代矩阵为 GSOR = (D − ωL) −1 ((1 − ω)D + ωU), 对应的矩阵分裂为 M = 1 ω D − L, N = 1 − ω ω D + U. 由 (6.12) 可知 SOR 迭代的分量形式为 x (k+1) i = (1 − ω)x (k) i + ω aii bi − ∑ i−1 j=1 aijx (k+1) j − ∑n j=i+1 aijx (k) j = x (k) i + ω aii bi − ∑ i−1 j=1 aijx (k+1) j − ∑n j=i aijx (k) j
第六讲线性方程组定常选代方法.180.算法6.5.求解线性方程组的SOR选代方法1: Given an initial guess r(0) and parameter w2:while not converge dofor i=l to n do3:2(*+1) = (1 - w) (k) +4:aiend for5:6:end while+SOR方法最大的优点是引人了松弛参数w:通过选取适当的w就可以大大提高方法的收敛速度但是SOR方法最大的难点就是如何选取最优的参数6.2.5SSOR送代方法将SOR方法中的L和U相互交换位置，则可得选代格式(k+1) = (D - wU)-1 (1 - w)D + wL) r(k) + w(D - wU)-1b将这个选代格式与SOR相结合，就可以得到下面的两步选代方法(++) = (D -wL)-[(1 -w)D + wUjr(k) + w(D -wL)-1b,(+1) = (D - wU)-[(1 - w)D + wL) r(k+) + w(D - wU)-1b.这就是对称超松弛（SSOR）迭代方法，相当于将L与U同等看待，交替做两次SOR代消去中间选代向量(+),可得(k+1) = GssORz(k) + g,其中送代矩阵GssOR = (D - wU)-1[(1 - w)D +wL](D -wL)-1[(1 -w)D +wU]对应的矩阵分裂为1[D-w(L+U)+w?LD-1]M :w(2-w(D-wL)D-1(D-wU),7w(2-w)u(2-[(1-w)D+wl)D-[(1-w)D+U],N=十对于某些特殊问题，SOR方法不收敛，但仍然可能构造出收敛的SSOR方法
· 180 · 第六讲 线性方程组定常迭代方法 算法 6.5. 求解线性方程组的 SOR 迭代方法 1: Given an initial guess x (0) and parameter ω 2: while not converge do 3: for i = 1 to n do 4: x (k+1) i = (1 − ω)x (k) i + ω aii ( bi − i ∑−1 j=1 aijx (k+1) j − ∑n j=i+1 aijx (k) j ) 5: end for 6: end while † SOR 方法最大的优点是引入了松弛参数 ω: 通过选取适当的 ω 就可以大大提高方法的收敛速度. 但是 SOR 方法最大的难点就是如何选取最优的参数. 6.2.5 SSOR 迭代方法 将 SOR 方法中的 L 和 U 相互交换位置, 则可得迭代格式 x (k+1) = (D − ωU) −1 ((1 − ω)D + ωL) x (k) + ω(D − ωU) −1 b. 将这个迭代格式与 SOR 相结合, 就可以得到下面的两步迭代方法 x (k+ 1 2 ) = (D − ωL) −1 [ (1 − ω)D + ωU] x (k) + ω(D − ωL) −1 b, x (k+1) = (D − ωU) −1 [ (1 − ω)D + ωL] x (k+ 1 2 ) + ω(D − ωU) −1 b. 这就是 对称超松弛 (SSOR ) 迭代方法, 相当于将 L 与 U 同等看待, 交替做两次 SOR 迭代. 消去中间迭代向量 x (k+ 1 2 ) , 可得 x (k+1) = GSSORx (k) + g, 其中迭代矩阵 GSSOR = (D − ωU) −1 [ (1 − ω)D + ωL] (D − ωL) −1 [ (1 − ω)D + ωU] . 对应的矩阵分裂为 M = 1 ω(2 − ω) [ D − ω(L + U) + ω 2LD−1U ] = 1 ω(2 − ω) (D − ωL)D−1 (D − ωU), N = 1 ω(2 − ω) [ (1 − ω)D + ωL] D−1 [ (1 − ω)D + ωU] . † 对于某些特殊问题, SOR 方法不收敛, 但仍然可能构造出收敛的 SSOR 方法