盖尔圆定理与特征值估计：从范数到相似变换的实用指南

做矩阵的人十个里有九个在特征值上栽过跟头。中小型矩阵求全特征值numpy一行就能出结果可矩阵一大或者参数带点误差精确求值不仅贵而且容易给你一种虚假的安全感。真正常用的做法反而是先估算特征值的大致范围——上界、下界、落在复平面的哪个区域——用很小的代价换一颗定心丸。这也是矩阵论里“特征值估计”这一支存在的意义。本文会沿着最实用的路线走先看几种朴素的上下界估计再重点讲盖尔圆Gershgorin圆盘定理包括它怎么用、怎么证明、怎么通过相似变换调整圆盘大小最后用一个三阶矩阵完整跑一遍手算和代码验证。适合正在学矩阵论、数值线性代数的学生也适合工程上想快速判断矩阵性质、又不想动不动上特征值求解器的从业者。1. 为什么“只求范围”常常比“精确计算”更实用1.1 现实场景里的一堆“范围判断”特征值的问题在工程里几乎都是“范围问题”。判断一个矩阵可不可逆本质是问0是不是特征值判断一个迭代法收不收敛本质是问谱半径最大特征值模是不是小于1判断一个动力系统稳不稳定本质是问所有特征值的实部是不是都小于0判断一个协方差矩阵能不能直接求逆本质是问最小特征值是不是远离0。这些问题的共同特征是只要知道特征值落在哪个区间答案就出来了根本不需要每个特征值的具体位置。举几个实际场景。第一是PageRank和网络分析特征值特别是谱半径决定了迭代算法的收敛速度但网络矩阵动辄上亿维你不可能先把全部特征值求出来再决定用哪个算法通常是先用估计判断一下谱半径的量级。第二是预处理器的调试很多迭代法收敛慢不是矩阵本身的问题而是谱分布不好先估一下谱范围能帮你决定要不要换预处理器省下的时间不是一点点。第三是控制系统的稳定性分析连续系统要求状态矩阵的特征值都在左半平面如果盖尔圆并集完全落在左半平面那么不需要求解器也能给出稳定结论。这些场景的共同点是要的不是“第几个特征值等于几”而是“它们大概在哪个区域”。区域性判断通常比精确计算便宜一个数量级以上。1.2 精确求特征值有多贵、多脆弱为什么不能每次都老老实实把特征值算出来答案有两个贵而且脆弱。先说贵。稠密n阶矩阵求全部特征值经典QR算法是O(n^3)量级。n1000时大约要10^9次浮点运算勉强还能忍n10000时就是10^12次单机跑起来得按小时算。而且很多工程矩阵不是拿来“求一次特征值”就结束的参数稍微一变就得重算。相比之下算一遍行和、列和、或者做一次简单的相似变换成本只有O(n^2)甚至更少。这个差距在大型问题上完全是两个世界。再说脆弱。特征值对矩阵扰动非常敏感这不是危言耸听。数值线性代数里著名的Wilkinson矩阵就是典型反例一个看着人畜无害的矩阵只要某个元素抖动一点点特征值就剧烈变化条件数可以达到10^13以上。工程中的矩阵都是从测量、离散化、估计中来的本身带误差算出来的特征值再精确也只是“某个近似矩阵”的特征值。所以在动手精算之前先用估计把握一下大体范围不仅省钱还能给结果上一道保险——万一求解器返回一个超出估计范围的答案第一反应不应该是怀疑估计而应该是怀疑求解过程出了问题。2. 先用范数和迹拿一个粗糙的上下界2.1 谱半径与诱导范数最朴素的估计是所有特征值的模都不超过任意一个诱导范数。也就是ρ(A) max_i |λ_i| ≤ ||A||。证明一句话设Ax λxx≠0取范数得到|λ|·||x|| ||Ax|| ≤ ||A||·||x||约去||x||就得到|λ|≤||A||。因为对任意特征值成立所以谱半径不超过范数。实操中最常用的两个诱导范数是行和范数和列和范数∞-范数是每行绝对值之和的最大值1-范数是每列绝对值之和的最大值。两个都算一下取小者当上界用。比如后面会反复用到的矩阵A [[10, 9, 0], [1, 2, 0], [0, 1, 3]]行和分别是19、3、4所以||A||_∞ 19列和分别是11、12、3所以||A||_1 12。谱半径ρ(A) 11所以1-范数给的上界12已经比较靠谱∞-范数19就比较松。这个差异提醒我们不同范数估计的松紧程度差别很大实际操作时多算几个再挑。2.2 特征值平方和约束Frobenius范数第二个常用约束是Schur不等式Σ_{i1}^n |λ_i|² ≤ ||A||F² Σ{i,j} |a_ij|²。理解它不需要硬背公式。对任意方阵A做Schur分解A UTU^H其中U是酉矩阵T是上三角矩阵。酉变换不改变Frobenius范数所以||A||_F ||T||_F。T的对角元就是A的特征值但T还有非对角元于是||T||F² Σ|λ_i|² Σ{ij}|T_ij|² ≥ Σ|λ_i|²。这个不等式的用途有两个。一个是验证特征值数量如果n个特征值模的平方和已经被一个较大的特征值占了很多其余特征值的模就必须很小。另一个是算上界|λ_i| ≤ ||A||_F对每个特征值都成立。对于上面那个A||A||_F² 100811419 196所以每个特征值的模不超过14。比1-范数的12松但比∞-范数的19紧。可见这几个界各有优劣实践中叠加使用。另外还有一个经典约束是特征值的乘积等于行列式和等于迹。这两个通常不单独用来做界但当你已经确定了n−1个特征值范围时可以用它们推断剩下的那个的范围。比如两个特征值实部都知道在某个区间第三个的实部就可以用迹来反推。2.3 实部虚部的另一个界Hermitian部分范数给的是模的界但很多问题关心的是实部。比如稳定性判断要实部小于0正定性判断要实部大于0。这时候可以用矩阵的Hermitian部分。记H (A A^H)/2。对A的任意特征值λ和特征向量x有x^H A x λ x^H x取实部得到 (x^H H x)/(x^H x) Re(λ)。因为H是Hermitian矩阵Rayleigh商的范围就在H的最小特征值和最大特征值之间所以λ_min(H) ≤ Re(λ) ≤ λ_max(H)。类似地令K (A − A^H)/(2i)可以得到虚部的界λ_min(K) ≤ Im(λ) ≤ λ_max(K)。对那个例子AH [[10,5,0],[5,2,0.5],[0,0.5,3]]它的特征值大致在[-0.4, 12.4]这个范围所以A所有特征值实部都在[-0.4, 12.4]内。真实特征值1、3、11确实都落在里面。这个方法在小矩阵上可以直接用numpy算出H的特征值在中大规模矩阵上也可以用Lanczos等方法只求两端特征值成本可控。这里有一个实操心得Hermitian部分的界通常不比盖尔圆紧但它给出的区间是“硬”的而且对复特征值同样有效。盖尔圆擅长处理模与位置Hermitian部分擅长处理实部两者配合起来往往能给出比单独使用任何一个都更完整的图像。3. 盖尔圆把特征值“圈”在圆盘里3.1 定理内容与一行证明Gershgorin圆盘定理是我在工程里用得最多的特征值估计工具。设A是一个n阶矩阵记第i行除对角元外的绝对值之和为R_i Σ_{j≠i} |a_ij|定义第i个盖尔圆为D_i { z ∈ C : |z − a_ii| ≤ R_i }。定理说A的所有特征值都落在这些圆盘的并集里。证明很短也很能体现这类问题的套路。取A的任意特征值λ以及对应的特征向量x。设|x_k|是x所有分量中绝对值最大的那个显然|x_k|0。把特征方程Ax λx的第k行写出来λ x_k Σ_{j1}^n a_kj x_j移项得到(λ − a_kk) x_k Σ_{j≠k} a_kj x_j。两边取模再用三角不等式放缩|λ − a_kk| |x_k| ≤ Σ_{j≠k} |a_kj| |x_j| ≤ Σ_{j≠k} |a_kj| |x_k| R_k |x_k|。约去|x_k|得到|λ − a_kk| ≤ R_k也就是λ落在第k个盖尔圆里。这里的关键是取特征向量的最大分量那一行来展开因为只有最大分量才能保证|x_j| ≤ |x_k|从而完成放缩。3.2 圆盘连成一片时特征值个数怎么算盖尔圆定理还有更强的一个版本如果k个盖尔圆组成的集合与其余n−k个圆互不相交那么这个连通分量里恰好包含A的k个特征值按代数重数计。这个结论的直观理解很重要。把A写成A D B其中D是对角部分B是非对角部分。考虑连续路径A(t) D tBt从0变到1。t0时A(0)D特征值就是对角元a_ii每个特征值正好落在自己的圆心上t1时A(1)A就回到原矩阵。随着t增大特征值是连续变化的而任意时刻每个特征值都落在“当前时刻的盖尔圆并集”里。当前时刻第i个盖尔圆的半径是t·R_i圆心始终是a_ii。如果某个特征值想从一个连通分量“跑”到另一个它必须穿过两个连通分量之间的空隙但空隙里没有任何盖尔圆覆盖特征值不可能出现在那里所以每个连通分量里的特征值数量从头到尾保持不变。t0时每个圆正好对应一个特征值因此最终每个连通分量中有几个圆就有几个特征值。这个性质比原定理实用得多。它给了你一个“计数”能力不仅能说特征值在这个并集里还能说这个区域里有几个特征值。当你看到两个圆交叠在一起、跟其他圆分开时就可以肯定这两个圆里一共只有两个特征值不多不少。3.3 列和半径、相似变换与常用推论除了行和半径还可以用列和半径。因为A^T与A有相同的特征值复矩阵要考虑A^H特征值是共轭关系而A^T的第i行非对角元绝对值之和正是A第i列的非对角元之和C_i Σ_{j≠i} |a_ji|。于是可以用C_i定义“列盖尔圆”。行圆并集和列圆并集都包含所有特征值所以它们的交集也包含所有特征值。实际计算中行圆和列圆往往能把范围压缩不少。比行列圆更进一步的是相似变换。对任意正对角矩阵D diag(d_1, ..., d_n)矩阵D^{-1}AD与A有相同的特征值。这个变换不改变特征值但是会改变非对角元的比例变换后第i行第j列元素变为a_ij·d_j/d_i所以新的盖尔圆半径是R_i Σ_{j≠i} |a_ij| · d_j / d_i。圆心不变还是a_ii。这就给了我们一个调节手段想让第i个圆变小就调整d_i和其他d_j的比例。这是盖尔圆工程应用里最核心的技巧后面实操部分专门演示。由盖尔圆定理可以直接得到两个常用推论。第一如果A严格对角占优也就是对每个i都有|a_ii| R_i那么原点不在任何盖尔圆里0不是特征值所以A可逆。这是数值线性代数里判断可逆性最廉价的方法。第二如果A是实矩阵且它的盖尔圆关于实轴对称那么复特征值必然成共轭对出现对应地出现在对称的两个圆里如果某个圆盘关于实轴不对称那它里面不可能含有非实数的共轭特征值。4. 实操手算、调圆、代码验证一条龙4.1 手算一个三阶矩阵的盖尔圆下面用矩阵B [[10, 9, 0], [1, 2, 0], [0, 1, 3]]来完整跑一遍。先算行半径第一行非对角元素绝对值之和是9第二行是1第三行是1所以行盖尔圆是D1: |z − 10| ≤ 9即实轴上区间[1, 19]D2: |z − 2| ≤ 1区间[1, 3]D3: |z − 3| ≤ 1区间[2, 4]。三个圆在实轴上的投影分别是[1,19]、[1,3]、[2,4]可以看到它们完全连成一片覆盖了[1,19]这一整段。也就是说仅靠行圆只能判断特征值都在[1,19]里但它到底是靠近1还是靠近19完全说不清。再用列半径算一次第一列非对角元绝对值之和是1第二列是10第三列是0所以列盖尔圆是E1: |z − 10| ≤ 1区间[9, 11]E2: |z − 2| ≤ 10区间[-8, 12]E3: |z − 3| ≤ 0就是单点3。列圆并集是[-8, 12]和行圆并集[1, 19]取交集得到所有特征值都在[1, 12]内。这个范围已经比单独用行圆好不少但还是不够精细。真实特征值其实是1、3、11其中11恰好被列圆E1锁定在[9,11]里3则被E3锁定为单点。这个例子说明行圆、列圆各算一遍把两个并集叠起来就已经能排除很多不可能的位置了。4.2 用相似变换把圆盘“收紧”到孤立接着做相似变换。目标是让第一行那个大圆缩小以便把特征值11单独隔离出来。第一行的大半径来源于a_12 9要让它对第一行的贡献变小就要让d_2/d_1尽量小。我取D diag(1, 0.5, 1)也就是d_1 1, d_2 0.5, d_3 1。计算变换后D^{-1}BD的元素第(1,2)元素变成9×0.5/1 4.5第(2,1)元素变成1×1/0.5 2第(3,2)元素变成1×0.5/1 0.5。于是新的行半径是R_1 4.5所以D1 [5.5, 14.5]R_2 |2| 0 2所以D2 [0, 4]R_3 0.5所以D3 [2.5, 3.5]。变换后在实轴上三个区间是D1 [5.5,14.5]D2 [0,4]D3 [2.5,3.5]。D2和D3明显重叠连成一个连通分量[0,4]D1与它们之间有空隙4到5.5之间没有圆覆盖。按照连通分量计数规则D1单独一个圆所以[5.5,14.5]内恰好有一个特征值D2和D3共两个圆所以[0,4]内恰好有两个特征值。真实特征值1、3、11完美对应1和3在[0,4]11在[5.5,14.5]。把几种估计方式放在一起看方法D1区间D2区间D3区间并集/连通情况行圆[1, 19][1, 3][2, 4]连成 [1, 19]列圆[9, 11][-8, 12][3, 3]连成 [-8, 12]与行圆取交集———[1, 12]相似变换[5.5, 14.5][0, 4][2.5, 3.5]两组{D2,D3} 与 {D1}这个例子把相似变换的价值体现得很明显原本三个圆连成一整段几乎给不出分区信息调一下d的尺度两个分量就分开了而且带着“个数”的结论。4.3 用Python快速完成全套验证手算容易出错实际工程里我都是写一个几行的脚本直接跑。下面这段代码可以在任何装了numpy的Python环境里运行。import numpy as np def gersh_circles(A): A np.asarray(A, dtypefloat) centers np.real(np.diag(A)) row_r np.sum(np.abs(A), axis1) - np.abs(np.diag(A)) col_r np.sum(np.abs(A), axis0) - np.abs(np.diag(A)) return centers, row_r, col_r A np.array([ [10., 9., 0.], [1., 2., 0.], [0., 1., 3.] ]) centers, row_r, col_r gersh_circles(A) print(圆心:, centers) print(行半径:, row_r) print(行区间:, list(zip(centers - row_r, centers row_r))) print(列半径:, col_r) print(列区间:, list(zip(centers - col_r, centers col_r))) D np.diag([1.0, 0.5, 1.0]) B np.linalg.inv(D) A D centers2, row_r2, col_r2 gersh_circles(B) print(相似变换后行半径:, row_r2) print(相似变换后行区间:, list(zip(centers2 - row_r2, centers2 row_r2))) print(真实特征值:, np.linalg.eigvals(A))输出大致是圆心: [10. 2. 3.] 行半径: [9. 1. 1.] 行区间: [(1.0, 19.0), (1.0, 3.0), (2.0, 4.0)] 列半径: [ 1. 10. 0.] 列区间: [(9.0, 11.0), (-8.0, 12.0), (3.0, 3.0)] 相似变换后行半径: [4.5 2. 0.5] 相似变换后行区间: [(5.5, 14.5), (0.0, 4.0), (2.5, 3.5)] 真实特征值: [11. 3. 1.]建议每次调试时把真实特征值打印出来对照确认估计没有算错。如果想可视化可以再加几行matplotlib画圆每行圆心画一个圆再把真实特征值画在复平面里一眼就能看出估计区域是否覆盖了所有特征值。5. 常见翻车点与实操心得5.1 盖尔圆不是“一圆一特征值”初学盖尔圆最容易犯的错误是把定理反向使用看到第i个圆就以为里面一定有特征值或者看到某个圆比较大就以为里面特征值很多。定理只保证每个特征值至少落在某一个圆里并不保证每个圆里都有特征值。一个简单的反例是2阶矩阵[[0, 10], [0.1, 0]]特征值是±1。行圆有两个第一个圆心0半径10第二个圆心0半径0.1。第二个圆 |z|≤0.1 里没有任何特征值因为±1都在外面。但它完全被第一个圆包含所以也不违反定理。这个例子提醒我们在做判断时要用“并集覆盖”和“连通分量计数”而不是“逐个圆对应”。5.2 相似变换的尺度怎么选相似变换不是随便选的。我的经验是分三步走。第一步先看哪些行半径太大找出是哪些非对角元造成的。第二步针对性地调整d_i和d_j的比例。想缩小第i个圆就减小所有与它相连的d_j/d_i也就是要么增大d_i要么减小其他d_j。第三步接受代价被你调小的圆可能变大。比如实操那个例子里为了缩小第一行的大圆我把d_2调小结果第二行的半径从1涨到了2但这个代价换来的是连通分量被拆开非常值。如果矩阵规模不大完全可以枚举几组d值看哪一组能把关心的区域孤立出来。如果矩阵规模很大有一种理论选法对非负矩阵B取|a_ij|、对角置0取d为B的Perron特征向量可以让所有圆盘同时收缩到最小。这个结论非常漂亮但工程上多数时候不必追求全局最优只要把关心的特征值区域“抠”出来就够了。5.3 结合实部虚部界做更稳的判断光靠盖尔圆有时候会误判。特别是当盖尔圆横跨虚轴或横跨0点的时候你无法确定特征值到底在左边还是右边。这时候建议叠加Hermitian部分的实部估计和反Hermitian部分的虚部估计。具体做法是先算A的H和K用numpy求H的最小、最大特征值得到实部硬范围再用盖尔圆看整体位置和连通分量如果盖尔圆范围很大、实部范围很小就说明特征值集中在某条竖直带内虚部可能较大。我在工程里的习惯是所有矩阵一上来先把这四个量全部算一遍——行圆、列圆、H的端特征值、K的端特征值。这几步成本都很低但能把特征值的分布大致框在一个很小的区域内之后再决定要不要做精确特征值分解。5.4 大型稀疏矩阵场景下的使用建议最后聊一下大规模问题。对n很大但稀疏的矩阵盖尔圆估计的计算量只有O(nnz)也就是非零元个数几乎可以忽略。这使它成为大型迭代算法里非常理想的预筛工具。实际用法之一是判断稀疏矩阵是否可逆。如果行圆并集不包含原点可以直接说可逆如果原点落在某个大圆里则无法下结论需要再用列圆或相似变换进一步判断。用法之二是给迭代法的收敛性提供参考。很多迭代算法要求谱半径严格小于1如果盖尔圆并集完全落在以原点为中心的单位圆内那谱半径一定小于1收敛性就有保障。用法之三是在做预处理的时候估计条件数的量级。条件数是最大奇异值与最小奇异值之比虽然不是特征值但估计思路类似。大规模场景下做相似变换要小心因为D虽然是对角阵变换后矩阵的稀疏结构不变但数值范围可能发生很大变化。我通常只在中小型矩阵上用相似变换精调大规模矩阵只做行圆、列圆和实部虚部估计图的就是省事和稳。上面这套流程概括起来就是先范数框模再Hermitian部分框实部然后用盖尔圆定位和计数最后用相似变换把模糊的地方调清楚。我自己的体会是特征值估计这活儿工具不难难的是结合具体问题判断“哪个界够用就行”。很多场合你根本不需要知道特征值精确到小数点后几位只要知道它不在0附近、不在虚轴附近、谱半径没有超过1就已经能做决定了。如果哪次估计结果恰好把你关心的区域框得很准那种感觉比跑通一个复杂算法还舒服。