矩阵特征值与特征向量的计算

1、第九章 矩阵特征值与特征向量的计算教学目的与要求:掌握用幂法和反幂法求矩阵特征值与特征向量的方法,了解 Jacobi 方法的适用范围和使用方法。 重点和难点:幂法和反幂法 教学内容:§1 幂法和反幂法一、 幂法幂法的基本思想是给定初始向量 (00x , 由迭代公式 产生向量序列(1( (0,1, 2, +=L k k x Ax k (k x :上述向量称为迭代向量。 (1(0(22(0( (0 =LLLLk k x Ax x A x x A x 于是由上式得(1 ( 1(01111( +k i u =nnk k k k i i i i i i x Ax A x A a u a 111

2、21112211+=+L k k k n n n a u a u a u设 ,由 10a 1(2,3, , i i n >=L 得 11lim 0+=k i i i k a u ,于是 121lim 0+=k ni i i k i a u故只要 k 充分大,就有 (1111111121+=+nk k k i i i i 1x a u a u a u 因此, 可以近似作为与 (1 +k x 1相应的特征向量。下面我们通过特征向量来计算特征值 1。用 ( k i x 表示 的第 i 个分量,由于( k x (1 1111(111( ( +k k i i k k i i x a x a u u

3、 ,所以 (11( (1,2, , +=L k i k ix i n x 上式这种由已知非零向量 及矩阵 (0x A 的乘幂 构造向量序列 kA ( k x 用来计算矩阵 A 按模最大的 特征值 1与对应的特征向量的方法称为 幂法 。例 1 用幂法的规范运算求矩阵 的按模最大的特征值及对应的特征向量。=1439A幂法的收敛速度取决于比值21,比值越小,收敛越快。反之,则很慢。此外,当矩阵 A 无 个线性 无关的特征量时,幂法收敛很慢,亦应考虑改用其它方法。n 二、反幂法分析反设 A 为 阶非奇异矩阵, ×n n , u 为 A 的特征值与相应的特征向量,即 =Au u ,由于,11=

4、A u u ,所以 1A 的特征值是 A 的特征值的倒数1,而相应的特征向量不变。如果 A 的特征值的次序为12n L ,则 1A 的特征值为11111L n n,因此,若对矩阵 1A 用幂法,即可计算出 1A 的按模最大的特征值n1,从而求得 A 的按模最小的特征值 n 。这就是反幂法的基本思想。因为反幂法的来源反幂法计算的主要步骤如下:1.对 A 进行三角分解 ; 2.求整数 LU A =r ,使得 ( ( 1max , k k rii nk rxx=x, 计算 (=k k x y;3. 解方程组(1+=k k Lz yUxz 例 2 用反幂法求矩阵 的最接近 =7232133312A 1

5、3的特征值,并求相应的特征向量。 §2 Jacobi方法Jacobi 方法的基本思想是:构造特殊的正交矩阵序列,通过正交变换,使 A 的非零的非对角元素逐次 化成零 ,并且使得非对角元素的平方和减小。从而使该矩阵近似为对角矩阵,得到全部特征值和特征向 量。一、 矩阵的旋转变换Jacobi 方法的关键是如何构造正交矩阵?先分析简单例子。 对二阶矩阵,只做依次正交变换,选择适当角 ,就可将 A 对角化。将这种思想推广到 维情况,设 n A 为 阶实对称矩阵,考虑 n nR 中平面旋转变换矩阵11cos sin 1( 1sin cos 11ij V =O L M O M L O它是将 n

6、阶单位阵中对角线上第 i 和第 j 个元换成 cos ,非对角元 和 分别换或 sin ij V ji V 和sin 而得到的,容易验证, ( ij V 是正交矩阵,若记 (1(1( T ij ij ijA V AV a = 则有 (122(122(1(1(1(1(1cos sin sin 2 sin cos sin 2 cos sin (, sin cos ii ii jj ij jj ii jj ij ik ki ik jk jk kj ik jk kl a a a a a a a a a a a a k i j a a a a a =+=+=+=+ (1(1(1(, , 1(sin 2c

7、os 22lk kl ij jk jj ii ij a a k l i j a a a a a =+如果 ,取 0ij a 使得 2tan 2(4ij ii jja a a =<则有 (1(10ij ji a a =,这样,就得到一个使 A 中非零的非对角元素 变成零的正交相似变换。对 重复上述过程,可得 ,如此继续下去,得到一个矩阵序列 =ij ji a a (1A (2A ( k A 。 可以证明由上述方法构造的旋转矩阵对 变换后, 就会使非对角元的平方和严格单调递减,而对角元的平方和单调递增。(k A 二、 J acobi 方法通过一系列旋转相似变换将 A 变成 (1+k A ,求

8、得 A 的全部特征值与特征向量的方法称为 Jacobi 方法。计算过程如下:(1令 ; (2求整数 (0, k k A =A j i , ,使得 (1, , maxk ij lml m n l ma =k a ; (3计算旋转矩阵 (4计算 ; (5计算 (1+k A(1(1 2( (k k lm l mE Aa +=, (令 =lk kl a A E 2(6若 (1( k E A+<,则 为特征值, 的各列为相应的特征 向量;否则, 返回 2,重复上述过程。(1 (1 (11122, , , k k k nna a a +L +k (0(1( k Q V V V =L 1k +数值计算

9、方法教案 -第九章 矩阵特征值与特征向量的计算定理 1 设 A 为 阶实对称阵,对 n A 用经典 Jacobi 方法得到序列 ( k A ,其中 (0=AA ,则( lim ( 0k k E A =例 3 用 Jacobi 方法求对称矩阵 的全部特征值。=612152224A 小结:1. 幂法的基本思想; 2. 幂法的应用; 3. 反幂法的基本思想及应用; 4. Jacobi方法作业:习题 9第 1, 2, 3题§3 方法QR 教学目的与要求:掌握用 QR 方法求矩阵特征值与特征向量的方法。 重点和难点:方法QR 教学内容:QR 方法是计算一般中小型矩阵全部特征值与特征向量的最有效

10、的方法之一。它是以矩阵正交三角分解为基础的一种变换方法。这里仅讨论实矩阵,并假定矩阵非奇异。一、矩阵的 QR 分解定理 9.2 任一非奇异实矩阵都可以分解成一个正交矩阵 Q 和一个上三角矩阵 R 乘积,而且当 R 的对 角元素为正时,分解是惟一的。对矩阵作 QR 分解的方法有多种,下面以 Schmit 正交化方法为例证明。 证明 设 A 为 阶非奇异实矩阵,将矩阵 n A 写成分块形式12=L n A a a a其中 ,因为 A 非奇异,所以 线性无关。=j a 12(, , , (1, 2, , j j nj a a a j n =L L 12L n a a a 取 111a a b =,

11、显然 11222b b a a b ><=, 12b b , 取 222/=b b b , 则 12121, , 0=b b b b 。 一般地,取11, =k kk k i i b a a b b i /(2,3, , =L k kk b b b k n 数值计算方法教案 -第九章 矩阵特征值与特征向量的计算则向量组 正交,且 12, , , L n b b b 1(1, 2, , =L k b k n 。式(9-13可改写成1111, , =+L k k k k k ka a b b a b b b b k 于是1211212121, , , , , , , , 0=L L L

12、 O Ln n n nna a b a b b A a a a b b b b b M QR =这就是用 Schmit 正交化方法对矩阵进行 分解的过程。QR 二、基本 QR 方法基本 方法的思想是利用矩阵的 QR 分解,通过迭代格式QR, 2, 1(1(L =+k Q R AR Q A k k k kk k 将 化成相似的上三角阵(或分块上三角阵 ,从而求出矩阵 (1=A A A 的全部特征值与特征向量。具体计算步骤为: 令 ,对 作 分解(1A=A (1A QR (111A Q R =令 (211A R Q =,作 分解QR (222A Q R = 重复上述过程,得迭代公式( (1(1,

13、2, k k kk k k A Q R k AR Q +=L 这样可以得到矩阵序列 ( k A ,并且矩阵是相似矩阵。事实上:由 (111A Q R =可得 ,于是11Q A R =1(211111A R Q Q AQ =数值计算方法教案-第九章 矩阵特征值与特征向量的计算 即A (2 与 A 相似。同理可得， A (k A(k = 2,3,L ，故它们有相同的特征值。有上述方法构造正交相似 矩阵的方法称为基本 QR 方法。 在一定条件下， 矩阵序列 A( k 收敛于一个上三角阵或分块上三角阵， 且对角块为 1×1 矩阵或 2×2 矩阵，其中 1×1 矩阵对角元素

14、为 A 的实特征值，2×2 矩阵含有 A 的一对复特征值。特别地，如果 A 是实 对称阵，则 A( k 收敛于对角矩阵。另外，当 K 充分大时， A 可以作为 A 的特征值的近似。 (k 的主对角元（或主对角子块的特征值）就 5 2 5 1 1 0 3 2 的特征值。 例 5 用基本 QR 方法求矩阵 A = 0 2 2 3 1 2 0 0 解 令A (1 = A ，对 A (1 作 QR 分解 A(1 = Q1 R1 0.9806 0.0377 0.1923 0.1038 5.0992 1.9612 5.4912 0.3922 0.1961 0.18804 0.8804 0.419

15、2 0.0000 2.0381 1.5852 2.5288 × = 0.0000 0.9813 2.5242 3.2736 0.1761 0.0740 0.0000 0.0000 0.0000 0.7822 0.3962 0.8989 0.0000 0.0000 0.0000 0.0000 然后，求得 A(2 = R1Q1 ，作 QR 分解 A(2 = Q2 R2 ，一直下去，得到 * * * 4.0000 1.8789 3.5910 * A12 = 1.3290 0.1211 * 1.0000 所以， A 的一个特征值为 4，一个特征值为 1 ，还有一对复特征值是方程 1.8789

16、 1.3290 2 3.5910 0.1211 =0 的根，即 2 1.8789 × 0.1211 1.3290 × ( 3.5910 = 0 ，得 1± 2i 。 事实上，矩阵 A 特征方程为 4 53 + 7 2 7 20 = 0 ，其特征值为 4,1 ， 1 ± 2i 。 小结：基本 QR 方法每次迭代都需作一次 QR 分解与矩阵乘法，计算量大，而且收敛速度慢。因此实 数值计算方法教案-第九章 矩阵特征值与特征向量的计算 际计算时，总是先用一系列相似变换将 A 化成拟上三角矩阵（称为上Hessenberg矩阵） ，然后对此矩阵用基 本 QR 方法。

17、这样可减少运算量。下面介绍一种化 A 为相似的拟上三角阵的方法：Householder变换。 作业：习题 9 第 8，12 题 数值实验： 1 程序设计基础知识： 求 A 的全部特征值、特征向量: EigensystemA 求 A 的数字特征值 EigenvaluseA 求 A 的特征向量组: EigenvectorsA 求 A 的特征多项式：DetA-x*IdentityMatrix3 ClearA,x A=1,2,1,-1,2,1,0,4,2; MatrixForm% EigenvaluesA EigenvectorsA; MatrixForm% DetA-x*IdentityMatrix

18、3 EigensystemA 2 1 0 2用幂法求矩阵 A = 0 2 1 按模最大特征值 1 和对应的特征向量 x1 0 1 2 取初始向量 V0=(0.5,0.5,1.1T A=2,-1,0,0,2,-1,0,-1,2; MatrixForm% vx=0.5,0.5,1.1; Dovy=A.vx;Printk," ",vy," ",vy1/vx1," ", vy2/vx2;vx=vy/MaxAbsvy,k,1,15 EigensystemA 数值计算方法教案-第九章 矩阵特征值与特征向量的计算 2 1 0 用反幂法求矩阵 A = 0 2 1 按模最小特征值 1 和对应的特征向量 x1 0 1 2 A=2,-1,0,0,2,-1,0,-1,2; MatrixForm% y=0,0,1.; Dox=LinearSolveA,y;Printk," y=x/MaxAbsx,k,1,20 EigensystemA ",x," ",y1/x1; 知识拓展及课题研究（适合研究生） ： 课题 1：利用极小值原理设计广义特征值算法，并评价其数值性质 课题 2：设