슬라이드 1

Size: px

Start display at page:

Download "슬라이드 1"

侨桥桑
5 years ago
Views:

1 08-09 年度第一学期计算方法 (B) 童伟华管理科研楼 05 室中国科学技术大学数学科学学院

2 第七章计算矩阵的特征值与特征向量

3 特征值与特征向量在实际工程计算中, 经常会遇到特征值和特征向量的计算, 如 : 机械结构或电磁振动中的固有值问题 ; 物理学中的各种临界值等定义 : 设有阶矩阵 A, 若有数 λ及非零向量 v满足方程 Av = λv, 则称 λ 为 A的特征值, v为属于特征值 λ的特征向量线性代数中特征值和特征向量的计算步骤 : 计算矩阵的特征多项式 : det( λi A) = λ + + ( ) det( A) ( A λ i Iv ) = 0, i =,,, 求解齐次线性方程组的解空间 : 困难 : 求解高次多项式的根 3

4 特征值与特征向量求解特征值与特征向量的方法按求解方法的特性可以分为 : 直接法 : 可用来计算矩阵的全部特征值和特征向量, 一般用于稠密矩阵 3, 计算代价为 O ( ) 且对矩阵的元素是相对不敏感的 ( 注意 : 此处的直接法必须作迭代, 因为求特征值问题等价于求多项式的根, 而后者不可能存在非迭代的方法如果经验表明, 用一个固定的迭代次数收敛 ( 几乎 ) 从不失败, 则称该方法是直接的 ) 迭代法 : 一般只提供特征值和特征向量的一个子集的近似, 一般用于稀疏矩阵或者只适合执行矩阵 - 向量乘法的矩阵为了得到几个足够精确的特征值可能需要运行足够长的时间, 收敛性强烈的依赖于矩阵的元素求解特征值与特征向量的方法按矩阵的对称性可以分为 : 非对称特征值问题 : 非对称矩阵的特征值问题对称特征值问题 : 对称矩阵的特征值问题, 不仅内容非常丰富和优美, 而且在实际问题中经常遇到 4

5 幂法在实际问题中, 矩阵按模最大的特征值往往起重要的作用, 譬如矩阵的谱半径决定了迭代矩阵是否收敛幂法 : 计算按模最大特征值及相应特征向量的数值方法要求 : 矩阵 A具有完备的特征向量系, 即 A 有个线性无关的特征向量在实际问题中, 常遇到的实对称矩阵或具有个互不相同特征值的矩阵就具有这种性质 5

6 幂法设矩阵的特征值和特征向量如下 : 特征值 : 特征向量 : 幂法可以求基本思想 : 不断利用矩阵向量乘法, 分析得到的向量序列, 计算出矩阵按模最大特征值及相应特征向量取初值 x = A x (0) Ax ( + ) ( ) λ λ λ v, v,, v λ, v, 做迭代 x Ax x α v α v α v ( + ) (0) (0) =, = ( + ) x = A ( αv+ αv + + αv) = αav+ αav + + αav = αλ v+ αλ v+ + αλ v A是非亏损的, 即特征值的几何重数等于代数重数 6

7 幂法情形 : 若 λ > λ 收敛速度取决于 : λ, 则有 ( ) λ λ x = λ αv+ αv + + αv λ λ ( α v ) ( α v ) ( ) x λ α 0, ( + ) + x λ λ x v x ( + ) ( ) ( + ) / x λ λ 7

8 幂法情形 : 若 λ = λ > λ3 λ, λ = λ ( ) λ x = λ αv+ ( ) αv + + αv λ ( ( ) ) α 0 x λ α v + α v, ( ) (+ ) ( ) + x / x λ v x + λx, + + x / x λ v x λx 其它情形 ( ) ( ) ( ) ( ) ( ) ( ) 8

9 幂法在幂法中, 我们构造的序列 ( ) λ λ x = λ αv+ αv + + αv λ λ + 可以看出, 当 0, λ < ( ) x, λ > ( ) 为避免 x 分量过大 ( 上溢 ) 或过小 ( 下溢 ), 在实际运算中采用规范运算 : ( + ) ( ) x = Ay y = x / x ( + ) ( + ) ( + ) 9

10 幂法不难发现 : 从而有 ( ) y = Ay / x = = A y / x x ( ) ( ) ( ) (0) ( ) () ( ) y = y = Ay / Ay y ( ) (0) (0) ( ) = λ λ λ αv+ α α λ v + + λ v λ λ λ αv+ α + + α λ v λ v 0

11 幂法情形 : 若 y λ > λ λ, 则有 α v ± λ < ( ) λ ( αv) ( ) αv y λ ( αv) α v λ > αv x = Ay = λ y ( + ) ( ) ( ) x = λy λ ( + ) ( ) ( ) A) 若 λ > 0, 则 { y } 收敛, λ, 0, 0 B) 若 λ < 0 ( ) ( ), 则 { y }{, y + } 分别收敛于互为反号的向 ( + ) ( ) 量, λ x, v y x, v y ( + ) ( )

12 幂法情形 : 若 y ( ) 令 = λ 其它情形 = λ > λ3 λ, λ = λ λ λ α + ( ) α + + α λ λ λ α + ( ) α + + α λ { ( ) } { ( ) y, y } x = Ax ( + ) ( ), 则有 v v v v v v 分别收敛于两个向量, 且不是互为反号, 则 v x + λ x v x x ( + ) ( ) ( + ) ( ) λ λ x / y, λ = λ ( + ) ( )

若 x是 A的特征向量, 那么 R( x) 的梯度为零向量 ; 反之, 若 R( x) = 0 且, 则是特征向量, R x 是相对应的特征值 x 0 x ( )

13 幂法 (Rayleigh 商 ) 设 A R, 对 x R 定义 R ( x) T = x Ax T xx 称为矩阵 A在 x 的 Rayleigh 商性质 : 对于给定的 x, 最小二乘问题 mi Ax αx α 有显示解, 即 α = R( x) T 不难发现 R( x) = ( Ax R( x) x)/ x x, 因此 : 若 x是 A的特征向量, 那么 R( x) 的梯度为零向量 ; 反之, 若 R( x) = 0 且, 则是特征向量, R x 是相对应的特征值 x 0 x ( ) 当趋于但未必等于特征向量时, R( x) 就是一个自然的特征值估计利用函数 R( x) 的光滑性, 可以得出一个重要结论 : R( x) R( q) = O( x q ), x q x 3

14 幂法幂法迭代算法 4

15 反幂法基本思想 : 计算按模最小特征值及相应特征向量的数值方法矩阵与特征值之间的关系 : Av = λv A v = v λ A : λ λ λ > 0 A A : A µ µ µ 反幂法迭代 : ( + ) ( ) x = A y ( + ) ( + ) ( + ) y = x / x A + + x = A y Ax = y 为避免求, 可用分解法求解线性方程组 : ( ) ( ) ( ) ( ) 5

16 反幂法设 λ 是矩阵特征值的近似, 带原点位移的反幂法 i A λ i λi 迭代进行如下 : 若 ( + ) ( ) ( ) λ x = A ii y ( + ) ( + ) ( + ) y = x / x λi λi < λj λi, j i, 则迭代收敛用途 : 用于求矩阵 A的最接近给定初值 λ 的特征值 λ i i及相应的特征向量 6

17 反幂法带原点位移的反幂法迭代算法 7

18 实对称矩阵的 Jacobi 方法特征值与特征向量的计算 : 要充分考虑矩阵的特性, 譬如对称性正定性稀疏性等, 来选用合适的算法定理 : 若矩阵使得 A 是对称的, 则存在正交阵 Q λ λ Q AQ =, λi, i =,,, λ 基本思想 : 对于一般的矩阵, 不可能直接找到 Q, 而是构造一系列特殊形式的正交阵对 A进行正交变换, 使得对角元比重逐渐增加, 非对角元变小 8

19 实对称矩阵的 Jacobi 方法 Gives 旋转变换 : cosθ siθ Q( pq,, θ ) = siθ cosθ p 列 q 列 9

20 实对称矩阵的 Jacobi 方法设矩阵, 则 T A = ( a ij ) B = Q (,, θ) AQ(,, θ ) = ( ) bip = bpi = api cosθ aqi si θ, i pq, biq = bqi = api siθ + aqi cos θ, i pq, bpp = app cos θ + aqq si θ apq si θ bqq = app si θ + aqq cos θ + apq si θ app aqq bpq = bqp = apq cosθ + si θ pq pq b ij 目标 : a pp aqq bpq = bqp = a pq cos θ + si θ = 0 0

21 实对称矩阵的 Jacobi 方法若记 aqq app s =, t = taθ a 则有 t pq t + ts = 按模的较小根, s 0 0 =, s = 0

22 实对称矩阵的 Jacobi 方法若记则有 cosθ = siθ = + t t + t c d bip = bpi = ca pi daqi, i p, q biq = bqi = da pi + caqi, i p, q bpp = a pp ta pq bqq = aqq + ta pq bpq = bqp = 0 b = a a = + ij ij pq i j i j b ii aii apq i= i=

23 实对称矩阵的 Jacobi 方法 Jacobi 方法 : 取 pq, 使得 a pq = max aij, 连续实施 Gives i j 变换, 即 ( + ) T ( ) A = Q ( pq,, θ) A Q( pq,, θ) 定理 : 若对称矩阵 A, 则 ( + ) A diag{ λ,, λ } 优点 : 特征值计算结果精度一般比较高, 特征向量的正交性也较好, 适合于并行计算缺点 : 当计算稀疏矩阵的特征值时,Gives 变换后不能保持原矩阵稀疏的性质用途 : 计算对称稠密矩阵的全部特征值及相应的特征向量 3

24 实对称矩阵的 Jacobi 方法实对称矩阵的 Jacobi 迭代算法 4

25 QR 方法 Householder 矩阵 : 若 v, v =, 则矩阵 T H = I vv 称为 Householder 矩阵, 且满足以下的性质 : () det( H) = ; () H 为对称正交矩阵 ; T (3) 若 xy,, x = y, 记 v =, H = I vv, y x 则 Hx = y Householder 变换 : 镜像变换的推广 y x 5

26 QR 方法 (QR 分解定理 ): 设 A 是列满秩的, 则矩阵 A能分解为 A = QR的形式, 其中 Q 是酉矩阵, R 是非奇异上三角阵推论 : 设是列满秩的, 则矩阵能唯一分解为 A A A = QR的形式, 其中 Q 是正交阵, R 是非奇异上三角阵对于非满秩的矩阵或一般形式的 m 阶矩阵, 有类似的结论 QR 方法 : 一种求解完全特征值问题的算法, 属于变换方法, 与 Gram-Schmidt 正交化过程密切相关 6

27 QR 方法 Householder 变换 : 设, 如何构造及使得 Hx = αe? x x e x x e 取 v = T,, 则有 0 x v H H = I vv Hx = x e = αe 计算 QR 分解的算法 : 设 Aααα = [,,, ] T 第一步 : 取 vαe= ( α )/ ρ, Q = I vv, 则 QAααααe = [ α e,,, ], α = sig( a ), ρ = α () () ( ) ( ) ( ) ( ) 第步 : 记 αα = [ a,,, T ] T, 取 α, I vαe= ( αq )/ ρ, = ( ) ( ) T I vv () ( ) ( + ) ( + ) = ααα a + α ai i + α + i= Q QA [ e, e e,, e e,,, ] 7

28 QR 方法不断利用 Householder 变换 : () (3) ( ) α a a3 a (3) ( ) α a3 a ( ) Q QA = α3 a = R 3 α A = ( Q Q ) R = ( Q Q ) R = ( Q Q ) R = QR 3 算法时间复杂度 : O(4 / 3) 8

29 QR 方法设矩阵的特征值满足 : 交序列 : A A = A > > > QR = RQ + 则当时, 有 ( ) aij 0, j< i ( ) aii λi, i =,,, 有以下关系 : λ λ λ A = Q A Q = Q Q AQ Q Q = Q Q, R = R R A Q R + + = + +, 作以下相似正收敛于上三角矩阵, 对角线元素即为矩阵的特征值 9

30 QR 方法 QR 方法的收敛性证明较困难, 参见 Matrix Computatios,4 th ed., p368 p37. QR 方法的改进带位移的 QR 算法实矩阵的双重步 QR 算法广义特征值问题 QZ 算法 (QR 算法的推广形式 ) 30

31 特征值与特征向量三对角阵二分法 Sturm 序列方法分而治之方法高阶稀疏矩阵 Laczos 方法 Aroldi with restartig Implicit restartig Jacobi-Davidso 方法 SVD 分解 Jacobi 方法 QR 方法 Golub-Kaha 方法 3

Microsoft PowerPoint - Eng-math-lecture14.ppt [Compatibility Mode]

Microsoft PowerPoint - Eng-math-lecture14.ppt [Compatibility Mode] -- 第讲一特征值与特征向量的概念定义设是阶矩阵如果数和维非零列向量 x 使关系式 x x 成立那末这样的数称为方阵的特征值非零向量 x称为的对应于特征值的特征向量说明特征向量 x 特征值问题是对方阵而言的阶方阵的特征值就是使齐次线性方程组 ( E x 有非零解的值即满足方程 E 的都是矩阵的特征值 // // E a a a a a a a a a