03-基本迭代方法

1 .第三讲 基本迭代方法 • 定常迭代法 • 收敛性分析 • 正则分裂 • 交替方向迭代法 • 加速技巧 

2 .迭代方法 求解线性方程组的一般方法: 直接法 和 迭代法. 迭代法: 定常迭代法 和 子空间迭代法. 快速算法通常是与问题本身 结构密切相关, 一般只适用 快速算法: 特殊问题的特殊结构, FFT, 多重网格, 快速多极子, 等等. 于某些特定问题的求解. 注记 在实际应用中, 这些方法往往结合使用, 如混合方法, 预处理方法等. 本讲主要内容: 在本讲中, 如果没有特别指 出, 一般假定 A 是非奇异的 定常迭代法 (基本迭代法, 矩阵分裂迭代法), 收敛性方向, 加速技巧. 实矩阵. 2/98 

3 .1 定常迭代法 1.1 矩阵分裂与定常迭代 1.2 Jacobi 迭代 1.3 Gauss-Seidel 迭代 1.4 SOR 迭代 1.5 SSOR 迭代 1.6 AOR 与 SAOR 迭代 1.7 Richardson 迭代 1.8 块迭代方法 定常迭代法有时也称为经典迭代法, 基本迭代法 或 不动点迭代法. 3/98 

4 .迭代法基本想法 当直接求解 Ax = b 比较困难时, 我们可以求解一个近似等价方程 组 M x = b , 其中 M 是对 A 的某种意义下的近似. 设 M x = b 的解为 x(1) . 则它与原方程的解 x∗ = A−1 b 之间的差满足 ( ) A x∗ − x(1) = b − Ax(1) 如果 x(1) 已经满足精度要求, 则可以停止计算, 否则需要进行修正. 4/98 

5 .设修正量为 ∆x, 则 ∆x 满足方程 A ∆x = b − Ax(1) 直接求解该方程比较困难, 因此我们还是求解近似方程 M ∆x = b − Ax(1) . 得到一个近似的修正量 ∆˜ x. 于是修正后的近似解为 x = x(1) + M −1 (b − Ax(1) ) x(2) = x(1) + ∆˜ 如果 x(2) 满足精度要求, 则停止计算, 否则继续按以上方式进行修正. 5/98 

6 .不断重复以上过程, 于是我们就得到一个序列 x(1) , x(2) , . . . , , x(k) , . . . . 满足以下递推关系 x(k+1) = x(k) + M −1 (b − Ax(k) ) 这就构成了一个迭代方法. 由于每次迭代的格式是一样的, 因此称为 定常迭代. 如何构造 好 的定常迭代方 法 (1) 以 M 为系数矩阵的线性方程组比原线性方程组更容易求解 (2) M 应该是 A 的一个很好的近似, 且迭代序列 {xk } 收敛 6/98 

7 .常见的定常迭代方法 • Jacobi 迭代 • Gauss-Seidel (G-S) 迭代 • 超松弛 (SOR, Successive Over-Relaxation) 迭代 • 对称超松弛 (SSOR, Symmetric SOR) 迭代 • 加速超松弛 (AOR, Accelerated Over-Relaxation) 迭代 • 交替方向 (ADI) 迭代和对称与斜对称 (HSS) 迭代 关键技术 矩阵分裂 7/98 

8 .1.1 矩阵分裂与定常迭代 定义 1 (矩阵分裂, Matrix Splitting) 设 A ∈ Rn×n 非奇异, 称 A=M −N (3.1) 为 A 的一个 矩阵分裂 , 其中 M 非奇异. 8/98 

9 .给定一个矩阵分裂 (3.1), 则原方程组 Ax = b 就等价于 M x = N x + b . 于是我们就可以构造出以下的迭代格式 M x(k+1) = N x(k) + b , k = 0, 1, . . . , 或 x(k+1) = M −1 N x(k) + M −1 b ≜ Gx(k) + g , k = 0, 1, . . . , (3.2) 其中 G ≜ M −1 称为 迭代矩阵. N 这就是基于矩阵分裂 A = M − N 的迭代方法. 选取不同的 M , 就可以构造出不同的迭代方法. 9/98 

10 .1.2 Jacobi 迭代 记 A = D − L − U , 其中 D 为 A 的对角部分, −L 和 −U 为 A 的严 格下三角和严格上三角部分. 取 M = D, N = L + U , 则可得 Jacobi 迭代 方法: x(k+1) = D−1 (L + U )x(k) + D−1 b , k = 0, 1, 2, . . . . (3.3) 对应的迭代矩阵为 GJ = D−1 (L + U ) ( ) (k+1) 1 ∑ n (k) 分量形式: xi = bi − aij xj , i = 1, 2, . . . , n. aii j=1,j̸=i 10/98 

11 .算法 1.1 Jacobi 迭代方法 1: Given an initial guess x(0) 2: while not converge do 3: for i = 1 to n do ( ) (k+1) 1 ∑ i−1 (k) ∑n (k) 4: xi = bi − aij xj − aij xj aii j=1 j=i+1 (k+1) Jacobi 迭代中 xi 的更新 5: end for 顺序与 i 无关, 因此非常适合 6: end while 并行计算 Jacobi 迭代格式也可以写为 x(k+1) = x(k) + D−1 (b − Ax(k) ), k = 0, 1, 2, . . . , 即 x(k+1) 是 x(k) 的一个修正. 11/98 

12 .1.3 Gauss-Seidel 迭代 取 M = D − L, N = U , 即可得 Gauss-Seidel (G-S) 迭代 方法: x(k+1) = (D − L)−1 U x(k) + (D − L)−1 b (3.4) 对应的迭代矩阵为 GGS = (D − L)−1 U 将 G-S 迭代改写为 Dx(k+1) = Lx(k+1) + U x(k) + b, 即可得分量形式 ( ) (k+1) 1 ∑ i−1 (k+1) ∑ n (k) xi = bi − aij xj − aij xj aii j=1 j=i+1 12/98 

13 .算法 1.2 Gauss-Seidel 迭代方法 1: Given an initial guess x(0) 2: while not converge do 3: for i = 1 to n do ( ) (k+1) 1 ∑ i−1 (k+1) ∑n (k) 4: xi = bi − aij xj − aij xj aii j=1 j=i+1 5: end for 6: end while 注记 G-S 迭代的主要优点是 充分利用了已经获得的最新数据. 在 G-S 迭代中, 未知量的更新必须按自然顺序进行, 不适合并行. 13/98 

14 .1.4 SOR 迭代 在 G-S 迭代中, 如果对修正量进行适当调整, 可能会加快收敛速度, 即 ( ) (k+1) (k) ω ∑ i−1 (k+1) ∑n (k) xi = xi + bi − aij xj − aij xj aii j=1 j=i 这就是 SOR (Successive Overrelaxation) 迭代方法. 其中 ω 为 松弛参数: • 当 ω = 1 时, SOR 即为 G-S 迭代, • 当 ω < 1 时, 称为低松弛 (under relaxation) 迭代, • 当 ω > 1 时, 称为超松弛 (over relaxation) 迭代. 14/98 

15 .SOR 迭代方法 (cont.) 事实上, 也可以将 SOR 看成是将 G-S 迭代中的第 k + 1 步近似解与 第 k 步近似解做加权平均所得到, 即 ( ) x(k+1) = (1 − ω)x(k) + ω D−1 (Lx(k+1) + U x(k) ) + D−1 b . (3.5) 整理后即为 注意向量加权与分量加权的 x(k+1) = (D − ωL)−1 ((1 − ω)D + ωU ) x(k) + ω(D − ωL)−1 b, (3.6) 区别. SOR 迭代曾经在很长一段时间内是求解线性方程组的首选方法 15/98 

16 .SOR 迭代方法 (cont.) SOR 的迭代矩阵为 ( ) GSOR = (D − ωL)−1 (1 − ω)D + ωU 对应的矩阵分裂为 1 1−ω M= D − L, N= D+U ω ω 16/98 

17 .算法 1.3 求解线性方程组的 SOR 迭代方法 1: Given an initial guess x(0) and parameter ω 2: while not converge do 3: for i = 1 to n do ( ) (k+1) (k) ω ∑ i−1 (k+1) ∑n (k) 4: xi = xi + bi − aij xj − aij xj aii j=1 j=i 5: end for 6: end while 注记 SOR 迭代最大的优点是引入了松弛参数 ω: 通过选取适当的 ω 就可 以大大提高方法的收敛速度. 但是 SOR 迭代最大的难点就是如何 选取最优的参数. 17/98 

18 .1.5 SSOR 迭代 将 SOR 迭代中的 L 和 U 相互交换位置, 则可得迭代格式 x(k+1) = (D − ωU )−1 ((1 − ω)D + ωL) x(k) + ω(D − ωU )−1 b. 将这个迭代格式与 SOR 相结合, 就可以得到下面的两步迭代方法  x(k+ 12 ) = (D − ωL)−1 [(1 − ω)D + ωU ]x(k) + ω(D − ωL)−1 b, [ ] x(k+1) = (D − ωU )−1 (1 − ω)D + ωL x(k+ 12 ) + ω(D − ωU )−1 b. 这就是 对称超松弛 (SSOR ) 迭代方法. 1 消去中间迭代向量 x(k+ 2 ) , 可得 x(k+1) = GSSOR x(k) + g, 其中迭代矩阵 [ ] [ ] GSSOR = (D − ωU )−1 (1 − ω)D + ωL (D − ωL)−1 (1 − ω)D + ωU . 18/98 

19 .对应的矩阵分裂为 1 [ ] M= D − ω(L + U ) + ω 2 LD−1 U ω(2 − ω) 1 = (D − ωL)D−1 (D − ωU ), ω(2 − ω) 1 [ ] [ ] N= (1 − ω)D + ωL D−1 (1 − ω)D + ωU . ω(2 − ω) 注记 对于某些特殊问题, SOR 迭代不收敛, 但仍然可能构造出收敛的 SSOR 迭代. 注记 一般来说, SOR 迭代的渐进收敛速度对参数 ω 比较敏感, 但 SSOR 的 收敛速度对参数 ω 不太敏感. 19/98 

20 .算法 1.4 SSOR 方法 1: Given an initial guess v (0) and parameter ω 2: while not converge do 3: for i = 1 to n do ( ) (k+ 12 ) (k) ω ∑ i−1 (k+ 21 ) ∑n (k) 4: xi = xi + bi − aij xj − aij xj aii j=1 j=i 5: end for 6: for i = n to 1 do ( ) (k+1) (k+ 12 ) ω ∑ i−1 (k+1) ∑ n (k+ 12 ) 7: xi = xi + bi − aij xj − aij xj aii j=1 j=i 8: end for 9: end while 20/98 

21 .1.6 AOR 与 SAOR 迭代 加速超松弛 (AOR , Accelerated Over-Relaxation ) 迭代方法: x(k+1) = (D−γL)−1 ((1 − ω)D + (ω − γ)L + ωU ) x(k) +ω(D−γL)−1 b, 其中 γ 和 ω 为松弛参数. 对应的矩阵分解为 1 1 M= (D − γL), N= [(1 − ω)D + (ω − γ)L + ωU ] ω ω 迭代矩阵为 [ ] GAOR = (D − γL)−1 (1 − ω)D + (ω − γ)L + ωU 21/98 

22 .AOR 迭代方法 (cont.) (1) 当 γ = ω 时, AOR 迭代即为 SOR 迭代; (2) 当 γ = ω = 1 时, AOR 迭代即为 G-S 迭代; (3) 当 γ = 0, ω = 1 时, AOR 迭代即为 Jacobi 迭代. 注记 AOR 迭代中含有两个参数. 因此在理论上, 通过选取合适的参数, AOR 迭代会收敛得更快. 但也是因为含有两个参数, 使得参数的选 取变得更加困难, 因此较少使用. 注记 与 SSOR 类似, 我们也可以定义 SAOR 迭代. 22/98 

23 .1.7 Richardson 迭代 Richardson 迭代是一类形式非常简单的方法, 其迭代格式为 x(k+1) = x(k) + ω(b − Ax(k) ), k = 0, 1, 2, . . . . 它可以看作是基于以下矩阵分裂的迭代方法: 1 1 M= I, N= I − A. ω ω 对应的迭代矩阵为 GR = I − ωA. 23/98 

24 .定理 1 设 A ∈ Rn×n 是对称正定矩阵, λ1 和 λn 分别是 A 的最大和 最小特征值, 则 Richardson 迭代方法收敛当且仅当 1 0<ω< . λ1 另外, Richardson 迭代的最优参数为 2 ω∗ = arg min ρ(GR ) = , ω λ1 + λn 即当 ω = ω∗ 时, 迭代矩阵的谱半径达到最小, 且有    1 − ωλn if ω ≤ ω∗  λ − λ 1 n κ(A) − 1 ρ(GR ) = = if ω = ω∗   λ1 + λn κ(A) + 1   ωλ1 − 1 if ω ≥ ω∗ . 24/98 

25 .注记 如果在每次迭代时取不同的参数, 即 x(k+1) = x(k) + ωk (b − Ax(k) ), k = 0, 1, 2, . . . , 则每次迭代的格式就不一样了, 因此不再是定常迭代, 而是 非定常 (Nonstationary ) 迭代. 此时称为 非定常 Richardson 迭代. 25/98 

26 .1.8 块迭代方法 将 A 写成分块形式:   A11 A12 · · · A1p   A21 A22 · · · A2p  A=  .. .. .. ..  ,  . . . .  Ap1 Ap2 · · · App 26/98 

27 .• 块 Jacobi 迭代: (k+1) ∑ p (k) Aii xi = bi − Aij xj , i = 1, 2, . . . , p. j=1,j̸=i • 块 Gauss-seidel 迭代: (k+1) ∑ i−1 (k+1) ∑ p (k) Aii xi = bi − Aij xj − Aij xj , i = 1, 2, . . . , p. j=1 j=i+1 • 块 SOR 迭代: ( ) ∑ i−1 ∑ p + ωA−1 (k+1) (k) (k+1) (k) xi = (1 − ω)xi ii bi − Aij xj − Aij xj , j=1 j=i+1 i = 1, 2, . . . , p. 27/98 

28 .2 收敛性分析 2.1 迭代法的收敛性 2.2 不可约对角占优矩阵 2.3 对称正定矩阵 2.4 相容次序 28/98 

29 .2.1 迭代法的收敛性 定义 2 (迭代法的收敛性) 对于迭代方法 x(k+1) = Gx(k) + g, 如果存在 x∗ , 使得对任意的初始向量 x(0) , 都有 lim x(k) = x∗ , k→∞ 则称迭代格式 (3.2) 是收敛的, 否则就称为发散的. 注记 基于矩阵分裂的迭代方法, 其收敛性取决于迭代矩阵的谱半径. 29/98