第五讲 对称特征值问题 Jacobi迭代方法 2 Rayleigh商迭代方法 3 对称QR迭代方法 4 分而治之法 5 对分法和反选代 6 奇异值分解 7 扰动分析 现代数位分析(数值线性代数),潘建瑜 http://math.ecnu.edu.cn/~jypan
第五讲 对称特征值问题 1 Jacobi 迭代方法 2 Rayleigh 商迭代方法 3 对称 QR 迭代方法 4 分而治之法 5 对分法和反迭代 6 奇异值分解 7 扰动分析 现代数值分析(数值线性代数), 潘建瑜 http://math.ecnu.edu.cn/~jypan
秦 对称特征值问题的常用算法有(直接法) ·Jacobi迭代:最古老,收敛速度较慢,但精度较高,适合并行计算。 ·Rayleigh商迭代:一般具有三次收敛性,但需要解方程组。 ·对称Q迭代:对于对称三对角矩阵,若只计算特征值,则速度最快(运算 量为O(m2),如果还需要计算特征向量,则运算量约为6m3。 ·分而治之法(Divide-and-Conquer):同时计算特征值和特征向量的一种快速 算法.基本思想是将大矩阵分解成小矩阵,然后利用递归思想求特征值.在 最坏的情形下,运算量为O(n3).在实际应用中,平均为O(2.3).如果使用 快速多极子算法(MM)后,理论上的运算量可降低到O(n log?n),其中p 是一个较小的整数,这使得分而治之算法成为目前计算对称三对角矩阵的 特征值和特征向量的首选方法。 http://math.ecnu.edu.cn/~jypan 2/77
对称特征值问题的常用算法有 (直接法) • Jacobi 迭代: 最古老, 收敛速度较慢, 但精度较高, 适合并行计算。 • Rayleigh 商迭代: 一般具有三次收敛性, 但需要解方程组。 • 对称 QR 迭代: 对于对称三对角矩阵, 若只计算特征值, 则速度最快 (运算 量为 O(n 2 )), 如果还需要计算特征向量, 则运算量约为 6n 3。 • 分而治之法 (DivideandConquer) : 同时计算特征值和特征向量的一种快速 算法. 基本思想是将大矩阵分解成小矩阵, 然后利用递归思想求特征值. 在 最坏的情形下, 运算量为 O(n 3 ). 在实际应用中, 平均为 O(n 2.3 ). 如果使用 快速多极子算法 (FMM) 后, 理论上的运算量可降低到 O(n logp n), 其中 p 是一个较小的整数, 这使得分而治之算法成为目前计算对称三对角矩阵的 特征值和特征向量的首选方法。 http://math.ecnu.edu.cn/~jypan 2/77
秦 ·对分法和反迭代:对分法主要用于求解对称三对角矩阵在某个区间中的特 征值,运算量约为O(k),其中k为所需计算的特征值的个数.反迭代用于 计算特征向量,在最佳情况下,即特征值“适当分离”时,运算量约为O(k), 但在最差情况下,即特征值成串地紧靠在一起时,运算量约为O(k2),而且 不能保证特征向量的精度(虽然实际上它几乎是精确的)。 除了Jacobi迭代和Rayleigh商迭代外,其余算法都需要先将对称矩阵三 对角化,这个过程大约需花费3”3的工作量,如果需要计算特征向量的 话,则运算量约为 http://math.ecnu.edu.cn/~jypan 3/77
• 对分法和反迭代: 对分法主要用于求解对称三对角矩阵在某个区间中的特 征值, 运算量约为 O(kn), 其中 k 为所需计算的特征值的个数. 反迭代用于 计算特征向量, 在最佳情况下, 即特征值 “适当分离” 时, 运算量约为 O(kn), 但在最差情况下, 即特征值成串地紧靠在一起时, 运算量约为 O(k 2n), 而且 不能保证特征向量的精度 (虽然实际上它几乎是精确的)。 ✍ 除了 Jacobi 迭代和 Rayleigh 商迭代外, 其余算法都需要先将对称矩阵三 对角化, 这个过程大约需花费 4 3 n 3 的工作量, 如果需要计算特征向量的 话, 则运算量约为 8 3 n 3 . http://math.ecnu.edu.cn/~jypan 3/77
秦 1|Jacobi达代 基本思想 通过一系列的Jacobi旋转将A正交相似于一个对角矩阵,即: A(0)=A,A(k1)=JA(k)J:=0,1;..., 且A收敛到一个对角矩阵,其中为Jacobi旋转,即Givens变换: jk cosk·-sin Ok Jk =G(ik,jk,0k)= ik sin0k·cos0g jk http://math.ecnu.edu.cn/~jypan 4/77
1 Jacobi 迭代 基本思想 通过一系列的 Jacobi 旋转 将 A 正交相似于一个对角矩阵,即: A (0) = A, A(k+1) = J ⊺ kA (k)Jk, k = 0, 1, . . . , 且 A(k) 收敛到一个对角矩阵, 其中 Jk 为 Jacobi 旋转, 即 Givens 变换: Jk = G(ik, jk, θk) = ik jk I cos θk · · · − sin θk . . . . . . . . . sin θk · · · cos θk I ik jk http://math.ecnu.edu.cn/~jypan 4/77
秦 引理设A∈R2x2是对称矩阵,则存在Givens变换G∈R2x2使得GTAG为 对角阵 (板书) http://math.ecnu.edu.cn/~jypan 5/77
引理 设 A ∈ R 2×2 是对称矩阵, 则存在 Givens 变换 G ∈ R 2×2 使得 G ⊺AG 为 对角阵. (板书) http://math.ecnu.edu.cn/~jypan 5/77
为了使得A收敛到一个对角矩阵,其非对角元素必须趋向于0. 秦 记o(A)为所有非对角元素的平方和,即 oiA=∑a喝=A呢-∑ 我们的目标就是使得of(A)尽快趋向于0. 引理设A=[aijnxn∈Rnxm是对称矩阵,A=[liilnxn=JTAJ,J= G(i,j,),其中0的选取使得a=ai=0,则 of(A)=ofA)-2a喝 (板书) http://math.ecnu.edu.cn/~jypan 6/77
为了使得 A(k) 收敛到一个对角矩阵, 其非对角元素必须趋向于 0. 记 off(A) 为所有非对角元素的平方和, 即 off(A) = X i̸=j a 2 ij = ∥A∥ 2 F − Xn i=1 a 2 ii , 我们的目标就是使得 off(A) 尽快趋向于 0. 引理 设 A = [aij ]n×n ∈ R n×n 是对称矩阵, Aˆ = [ˆaij ]n×n = J ⊺AJ, J = G(i, j, θ), 其中 θ 的选取使得 aˆij = ˆaji = 0, 则 off(Aˆ) = off(A) − 2a 2 ij (板书) http://math.ecnu.edu.cn/~jypan 6/77
算法1.1 Jacobi迭代算法 1:Given a symmetric matrix AE Rnxn 2:if eigenvectors are desired then set J=I and shift=1 4:end if 5:while not converge do 6: choose an index pair (i,j)such that aij 7: T=(aii-ajj)/(2aij) 8: t=sign(T)/(+v1+T2) %计算tan0 9: c=l/W1+tz,s=c~t%计算Givens变换 10: A=G(i,j,0)AG(i,j,0) %实际计算时不需要做矩阵乘积 1: if shift =1 then 12: J=J·G(i,j,) 13: end if 14:end while
算法 1.1 Jacobi 迭代算法 1: Given a symmetric matrix A ∈ R n×n 2: if eigenvectors are desired then 3: set J = I and shif t = 1 4: end if 5: while not converge do 6: choose an index pair (i, j) such that aij ̸= 0 7: τ = (aii − ajj )/(2aij ) 8: t = sign(τ )/(|τ | + √ 1 + τ 2) % 计算 tan θ 9: c = 1/√ 1 + t 2, s = c · t % 计算 Givens 变换 10: A = G(i, j, θ) ⊺AG(i, j, θ) % 实际计算时不需要做矩阵乘积 11: if shif t = 1 then 12: J = J · G(i, j, θ) 13: end if 14: end while
秦 a的选取问题 直观选取方法:所有非对角元素中绝对值最大的一个→经典Jacobi算法 算法l.2经典Jacobi迭代算法 1:Given a symmetric matrix AE Rnxn 2:while off(A)>tol do 3: choose (i,j)such that laijl maxk lakil 4: T=(ai-ai)/(2a) 5: t=sign(T)/(+v1+72) 6: c=1/W1+t2,s=c.t 7: A=G(i,j,0)TAG(i,j,0) 8: J=J.G(i,j,0)%if eigenvectors are desired 9:end while http://math.ecnu.edu.cn/~jypan 8/77
aij 的选取问题 直观选取方法: 所有非对角元素中绝对值最大的一个 =⇒ 经典 Jacobi 算法 算法 1.2 经典 Jacobi 迭代算法 1: Given a symmetric matrix A ∈ R n×n 2: while off(A)> tol do 3: choose (i, j) such that |aij | = maxk̸=l |akl| 4: τ = (aii − ajj )/(2aij ) 5: t = sign(τ )/(|τ | + √ 1 + τ 2) 6: c = 1/√ 1 + t 2, s = c · t 7: A = G(i, j, θ) ⊺AG(i, j, θ) 8: J = J · G(i, j, θ) % if eigenvectors are desired 9: end while http://math.ecnu.edu.cn/~jypan 8/77
可以证明,经典Jacobi算法至少是线性收敛的. 秦 定理对于经典Jacobi算法1.2,有 omu)s)of),D) 故k步迭代后,有 4o≤(-)4m)=(-)4 http://math.ecnu.edu.cn/~jypan 9/77
可以证明, 经典 Jacobi 算法至少是线性收敛的. 定理 对于经典 Jacobi 算法 1.2, 有 off(A (k+1)) ≤ 1 − 1 N off(A (k) ), N = n(n − 1) 2 . 故 k 步迭代后, 有 off(A (k) ) ≤ 1 − 1 N k off(A (0)) = 1 − 1 N k off(A). http://math.ecnu.edu.cn/~jypan 9/77
事实上,经典Jacobi算法最终是二次局部收敛的. 秦 定理经典Jacobi算法1.2是N步局部二次收敛的,即对足够大的k,有 oA+)=O(o(4)) 么经典Jacobi算法的每一步都要寻找绝对值最大的非对角元,比较费时. http://math.ecnu.edu.cn/-jypan 10/77
事实上, 经典 Jacobi 算法最终是二次局部收敛的. 定理 经典 Jacobi 算法 1.2 是 N 步局部二次收敛的, 即对足够大的 k, 有 off A (k+N) = O off2 A (k) . ✍ 经典 Jacobi 算法的每一步都要寻找绝对值最大的非对角元, 比较费时. http://math.ecnu.edu.cn/~jypan 10/77