数值分析第章——矩阵特征值问题计算ppt课件.ppt

1 第八章矩阵特征值问题计算 对n阶方阵A求数 和非零向量x 使其满足Ax x这样的 值称为矩阵A的特征值 非零向量x称为矩阵A的与特征值 相对应的一个特征向量 2 定义1设矩阵A B Rn n 若有可逆阵P 使则称A与B相似 定理1若矩阵A B Rn n且相似 则 1 A与B的特征值完全相同 2 若x是B的特征向量 则Px便为A的特征向量 8 1预备知识 3 定理2 设A Rn n具有完全的特征向量系 即存在n个线性无关 其中 i为A的特征值 P的各列为相应于 i的特征向量 的特征向量构成Rn的一组基底 则经相似变换可化A为 对角阵 即有可逆阵P 使 4 定理3 A Rn n 1 n为A的特征值 则 2 A的行列式值等于全体特征值之积 即 1 A的迹数等于特征值之和 即 5 定理4 6 定理5设A Rn n为对称矩阵 其特征值 1 2 n 则 1 对任意A Rn x 0 2 3 7 定理6 Gerschgorin圆盘定理 设A Rn n 则 表示以aii为中心 以半径为的复平面上的n个圆盘 2 如果矩阵A的m个圆盘组成的并集S 连通的 与其余 1 A的每一个特征值必属于下述某个圆盘之中 n m个圆盘不连接 则S内恰包含m个A的特征值 8 9 10 11 12 13 定理7 14 15 一幂法 1基本思想 任取非零向量v0 则可唯一表示为 设n阶矩阵A的特征值 满足 且其对应有n个线性无关的特征向量x1 x2 xn 即 8 2幂法和反幂法 16 则 17 其中 由假设条件 所以当k充分大时 有 从而 所以 18 即为矩阵A的对应特征值 1的近似特征向量 且 19 则对任意非零初始向量 下面的式子成立 定理 20 迭代公式 1 实质上是由矩阵A的乘幂Ak与非零向量v0相乘来构造向量序列 xk 从而计算主特征值及其对应的主特征向量 故称这种方法为幂法 21 则对任意非零初始向量 定理 按照下述方法构造的向量序列 则有 22 2 幂法实用计算公式 23 解取v0 0 0 1 T 则 24 直到k 8时的计算结果见下表 从而 25 二 幂法的加速 1 原点平移法 如果 是矩阵A的特征值 则对任意的实数p 矩阵A pE的特征值为 p 且A与A pE的特征向量相同 据此 如果要计算A的主特征值 1 只要选择合适的数p 使 1 p为矩阵A pE的主特征值 且 那么 对矩阵A pE应用幂法求其主特征值 1 p 收敛速度将会加快 这种通过求A pE的主特征值和特征向量 进而得到A的主特征值和特征向量的方法叫原点平移法 26 且使 显然 当 2 p n p 即P 2 n 2时 上式取最小值 如果希望计算 n 类似的讨论可知应选取p 1 n 1 2 则对任意实数p 矩阵A pE的主特征值为 1 p或 n p 要求 1 则选p使 27 例2用原点平移加速法求例1中矩阵A的主特征值与其对应的特征向量 解取p 2 5 做平移变换B A pE 则 对B应用幂法 仍取x0 0 0 1 T 则 28 迭代5步的计算结果见下表 0 5 1 0 7500 13 5000 6 7500 13 5000 10 1250 5 0 5 1 0 7500 13 5007 6 7503 13 5007 10 1256 4 0 5 1 0 7507 13 5179 6 76 13 5179 10 1406 3 0 5 1 0 7545 14 7 14 10 5625 2 0 5 1 0 875 4 2 4 3 5 1 k 可得到B的主特征值 1 13 5000特征向量v1 0 5 1 0 0 7500 T因此 A的主特征值为 1 1 p 11 0000 特征向量仍为v1 0 5 1 0 7500 T 29 30 31 设A为n阶实对称矩阵 称 为向量x的瑞利商 其中 x x xTx为内积 不难证明 对实对称矩阵A 如果其特征值满足 2 瑞利商加速 由幂法公式生成的xk的瑞利商满足 由此可见 R xk 比mk更快的收敛于 1 32 幂法的瑞利商加速迭代公式为 其中A为n阶实对称矩阵 对给定的误差限 当 mk mk 1 时 取 33 三 反幂法 反幂法是用于求非奇异矩阵A的按模最小的特征值和对应特征向量的方法 而结合原点平移法的反幂法则可以求矩阵A的任何一个具有先验了解的特征值和对应的特征向量 设矩阵A非奇异 其特征值 i i 1 2 n 满足 其相应的特征向量x1 x2 xn线性无关 则A 1的特征值为1 i 对应的特征向量仍为xi i 1 2 n 34 此时 A 1的特征值满足 因此 对A 1应用幂法 可求出其主特征值1 n 和特征向量xn uk 从而求得A的按模最小特征值 n 1 和对应的特征向量xn uk 这种方法称为反幂法 35 为了避免求A 1 可通过解线性方程组Avk uk 1得到yk 采用LU分解 即先对A进行LU分解A LU 此时反幂法的迭代公式为 36 对给定的误差 当 时 得 显然 反幂法的收敛速度取决于比值 比值越小 收敛越快 例3用反幂法求矩阵的对应于特征值的特征向量 37 解取解方程组 得 解方程组 得 38 与的对应分向量大体上成正比 所以对应于的特征向量为 39 QR算法也是一种迭代算法 是目前计算任意实的非奇异矩阵全部特征值问题的最有效的方法之一 该方法的基础是构造矩阵序列 并对它进行QR分解 由线性代数知识知道 若A为非奇异方阵 则A可以分解为正交矩阵Q与三角形矩阵R的乘积 即A QR 而且当R的对角线元素符号取定时 分解式是唯一的 8 4QR算法 40 若A为奇异方阵 则零为A的特征值 任取一数p不是A的特征值 则A pI为非奇异方阵 只要求出A pI的特征值 就很容易求出A的特征值 所以假设A为非奇异方阵 并不防碍讨论的一般性 设A为非奇异方阵 令 对进行QR分解 即把分解为正交矩阵与上三角形矩阵的乘积 作矩阵继续对进行QR分解 并定义 41 一般地 递推公式为 QR算法就是利用矩阵的QR分解 按上述递推公式构造矩阵序列 只要A为非奇异方阵 则由QR算法就完全确定这个矩阵序列具有如下性质 性质1所有都相似 它们具有相同的特征值 42 证明 因为 若令 则为正交阵 且有因此与A相似 他们具有相同的特征值 43 性质2的QR分解式为其中 证明用归纳法 显然当k 1时 有 假设有分解式 于是 44 因为 所以 定理1如果收敛于非奇异矩阵 为上三角形矩阵 则存在并且是上三角形矩阵 因为都是正交阵 所以也是正交阵 同样也是上三角阵 从而的QR分解式为 45 证明因为收敛 故下面极限存在 由于为上三角形矩阵 所以为上三角形矩阵 又因为所以存在 并且是上三角形矩阵 46 定理2 QR算法的收敛性 设A为n阶实矩阵 1 A的特征值满足 2 其中且设有三角分解式 L的单位下三角阵 U为上三角阵 则由QR算法得到的矩阵序列本质上收敛于上三角形矩阵 即满足 47 的极限不一定存在 证明因为矩阵决定的收敛性 又我们利用求 然后讨论的收敛性 由定理条件得 48 其中的 i j 元素为 于是 由假设 当i j时 故 设方阵X的QR的分解式为 由 49 由知 对充分大的k 非奇异 它应有唯一的QR分解式 并且 于是 但上三角阵的对角线元素不一定大于零 为此 引入对角阵 以便保证的对交线元素都是正数 50 从而得到的QR分解式 由的QR分解式的唯一性得到 从而 由于所以 51 从而 其中 于是 因为为上三角阵 为对角阵 且元素为1或 1 所以 52 的极限不一定存在 53 例1用QR算法求矩阵特征值 A的特征值为 1 4 1 2i 54 解令用施密特正交化过程将分解为 55 将与逆序相乘 求出 用代替A重复上面过程 计算11次得 56 由不难看出 矩阵A的一个特征值是4 另一个特征值是 1 其它两个特征值是方程 的根 求得为 57 上Hessenberg化 58 59 60 61 62 63 64 用正交变换化对称矩阵为对称三对角阵 65 带原点位移的QR算法 66 67 用单步QR方法计算上Hessenberg特征值 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 隐式QR算法 88 89 90 91 92 93 定理6 设A为n阶实对称矩阵 则必存在正交矩阵 P使得 其中 为A的n个特征值 证 设A的互不相同的特征值为 根据性质1和性质3知 对应于特征值 94 把它们标准正交化 特征向量组 特征向量共有n个 并有 其中 为A的n个特征值 这样的 特征值的特征向量是正交的 向量两两正交 以它们为列向量构成正交矩阵P 又由性质2知 A的属于不同 故这n个单位特征 95 由定理6可知 实对称矩阵的对角化问题 实质上是求正交矩阵P的问题 计算P的步骤如下 1 2 求出齐次线性方 标准正交的特征向量 求出实对称矩阵A的全部特征值 对于各个不同的特征值 对基础解系 这个向量组所含向量 96 3 4 则P为正交矩阵且使得 为对角阵 对角线上的元素为 相应特征向量的特征值 取 97 例13 解 得特征值 98 99 7 2对称QR方法 对称矩阵的三对角化 100 带原点位移的QR迭代 101 102 隐式QR迭代 103 104 定义Jacobi方法是用来计算实对称矩阵A的全部特征值及其相应特征向量的一种变换方法 Jacobi方法的基本思想是通过一系列的平面旋转矩阵构成的正交变换将对称矩阵逐步化为对角阵 从而得到A的全部特征值及其相应的特征向量 定理 1 如果n阶方阵A满足则称A为正交阵 7 3Jacobi方法 105 2 设A是n阶实对称矩阵 则A的特征值都是实数 并且有互相正交的n个特征向量 3 相似矩阵具有相同的特征值 4 设A是n阶实对称矩阵 P为n阶正交阵 则B PAP也是对称矩阵 5 n阶正交矩阵的乘积是正交矩阵 设A是n阶实对称矩阵 则必有正交矩阵P 使 106 由 6 可知 对于任意的阶实对称矩阵A 只要能求得一个正交阵P 使 则可得到A的全部特征值及其相应的特征向量 这就是雅克比方法的理论基础 其中的对角线元素是A的n个特征值 正交矩阵P的第i列是A的对应于特征值的特征向量 下面我们详细介绍雅可比方法 首先引进中的平面旋转变换 变换 107 记为 其中 108 X Y 称为平面内的平面旋转矩阵 容易得到如下性质 则称x 为中平面内的一个旋转变换 2 的主对角线元素中除第i个与第j个元素为外 其它元素均为1 非对角元素中除第i行第j列元素为 第j行第i列元素为外 其它元素均为零 1 为正交矩阵 109 3 只改变A的第i行第j行元素 AP只改变A的第i列第j列元素 所以只改变A的第i行 第j行 第i列 第j列元素 设A 为n阶实对称矩阵 为一非对角线元素 令则为实对称矩阵 且与A有同的特征值 通过直接计算知 110 111 当取满足关系式时 且 由于在正交相似变换下 矩阵元素的平方和不变 所以若用D A 表示矩阵A的对角线元素的平方和 用S A 表示A的非对角元素平方和 则由 11 式得 112 的非对角元素平方和A的非对角元素平方 且将事先选定的非对角元素消去了 这说明用对A作正交相似变换化为后 的对角线元素平方和比A的对角元素平方和增加了 和减少了 113 因此 只要我们逐次地用这种变换 就可以使得矩阵A的非对角线元素平方和趋于零 也即使得矩阵A逐步化为对角阵 这里需要说明一点 并不是对矩阵A的每一对非对角线非零元素进行一次这样的变换就能得到对角阵 因为在用变换消去的时候 只有第i行 第j行 第i列 第j列元素在变化 如果或为零 经变换后又往往不是零了 114 雅克比方法就是逐步对矩阵A进行正交相似变换 消去非对角线上的非零元素 直到将A的非对角线元素化为接近与零为止 从而求得A的去全特 征值 把逐次的正交相似变换矩阵乘起来 便是所要求的特征向量 Jacobi的计算步骤归纳如下 第一步在矩阵A的非对角元素中选取一个非零元素 一般来说 取绝对值最大的非对角线元素 115 第二步由公式求出 从而得平面旋转矩阵 第三步 的元素由公式 9 计算 第四步以代替A 重复第一 二 三步求出及继续重复这一过程 直到的非对交线元素全化为充分小 即小于允许误差 时为止 第五步的对角化元素为A的全部特征值的近似值 的第j列为对应于特征值 为的对角线上第j个元素 的特征向量 116 例1用雅克比方法求矩阵的特征值与特征向量 解首先i 1 j 2 由于故取所以 117 再取i 1 j 3 由 得 所以 118 继续做下去 直到非对角线元素趋于零 进行九次变换后 得 的对角线元素就是A的特征值 即 119 继续做下去 直到非对角线元素趋于零 进行九次变换后 得 的对角线元素就是A的特征值 即 120 用雅可比方法求得的结果都比较高 特别是求得的特征向量正交性很好 所以雅可比方法是求实对称矩阵的全部特征值及其对应特征向量的一个较好的方法 但由于上面介绍的雅可比