§3计算实对称矩阵特征值的Jacobi方法

前面介绍计算矩阵的部分特征值和特征向量的一些方法 现介绍的Jacobi方法是计算一个实对称矩阵的全部特征值和特征向量的方法 设是任一n阶实对称矩阵 则必存在一个直交相似变换矩阵U将它化为一个对角阵 即其中就是A的n个特征值 而U的第j列向量是与相应的特征值向量 但是这种直交相似变换矩阵U必须用一个无限迭代过程求得 也就是说 必须通过一系列的直交相似变换把矩阵A化为一个对角阵 3计算实对称矩阵特征值的Jacobi方法 Jacobi方法就是通过一系列特殊的直交相似变换矩阵 Givens平面旋转 把矩阵A化为一个对角阵 用R p q 表示如下形式的阶矩阵 3 1 3 1Givens平面旋转矩阵 它的元素除 p q 外 其余元素与n阶单位阵相应位置元素相同 通常称这种矩阵R p q 为Givens平面旋转矩阵 为旋转角 容易验证 R p q 是一个直交阵 有若令 3 2 则矩阵B和A中 只有p q二行和p q二列元素有区别 它们之间有如下关系 3 3 为使bpq 0 必须即旋转角应满足关系式 3 4 常将限制在区间上 若 则可取 3 5 称A的 p q 位置元素为旋转主元 用表示原来给定的n阶实对称矩阵 即令其中 古典的Jacobi方法首先选取矩阵主对角线上方的绝对值最大的元素作为旋转主元 根据 3 4 或 3 5 式选取使矩阵的 p q 元素等于零 即 设此法已进行k 1步 得到相似矩阵 其主对角线上方绝对值最大的元素为 则第k步将选取作为旋转主元 取Givens平面旋转矩阵R p q 中的旋转角满足 3 2Jacobi方法及其收敛性 或并令 3 6 一般地 经过有限次上述变换不可能把A化为一个对角阵 这是因为在中的元素 但在中 可能变成非零元素 因此 Jacobi方法是一种迭代法 然而 我们将证明 当时 3 7 其中是矩阵A的n个特征值 记 3 8 其中的主对角元全为零 据 3 6 和 3 3 式 有 从而有 3 9 由于是中绝对值最大的元素 因此即 3 10 据 3 9 和 3 10 式便有当时 这就证明 当时 趋于一个对角阵 一 旋转主元的选取在古典的Jacobi方法中 每一步变换之前都要寻查的主对角线上方绝对值最大的元素 以它作为旋转主元 在计算机上这样寻查一遍比较费时间 实际应用Jacobi方法时 常作一些修改 以节省寻查旋转主元的时间 首先 逐次选取 1 2 1 3 1 n 2 3 2 n n l n 元素为旋转主元 确定旋转角和Givens平面旋转矩阵依次消去上三角部分的 p q 元素 从 1 2 到 n 1 n 称为一轮 做完第一轮后 再按 1 2 1 3 n 1 n 的次序做第二轮 第三轮 在每一轮消元时 都给定一个控制量 称为消元容限 如果 p q 元素的绝对值小于 就跳过这一步 容限的值 逐轮减小 最后可取 3 3实用的Jacobi方法及其计算步骤 接近计算机所能表示的最小正数作为容限 但是 没有统一的方法来确定消元容限 常用的一种方法是在第一轮取容限第二轮取容限其中是做完第一轮消元后得到的的 i j 元素 仿此继续选取第三轮消元容限 或者 更简单些 在第一轮取以后逐轮取 平面旋转矩阵的旋转角应满足关系 为了避免分母出现零的危险 我们把这个关系式改写成记 3 11 则 3 12 在很小时 就跳过这一步 因此实际上在 3 12 中不会发生分母为零的情形 为了计算和 可令 3 13 二 平面旋转矩阵参数和的计算 由于 所以t满足二次方程 3 14 解得 3 15 为取 3 15 中绝对值的较小者 则应取 3 16 如果 那么据 3 12 式计算b将会产生较大的误差 因此 取 3 17 此时 t近似地满足方程 3 14 最后 按下列公式计算和 据 3 3 3 18 和 3 19 式 并注意到便可得到计算的元素的公式如下 3 20 在时 3 21 3 22 3 23 三 元素的计算 设逐次所用的平面旋转矩阵为 则据 3 6 式有令 3 24 则 3 25 设可以看成是一个对角阵 非主对角元素都接近于零 则的主对角元素是矩阵A的特征值的近似 据 3 25 式可得从而的第j列向量就是矩阵A的特征值所对应的特征向量 并且得到的特征向量系是一个标准直交系 记 四 特