系列笔记为上科学计算课上随意记下的笔记，可能有一些凌乱，勉强凑合着看。以记录和梳理为主，证明较少。
关于特征值的求解主要有两种方法，一类是基于矩阵相乘的迭代幂法，另一类是基于正交相似变换的方法。
幂法
幂法是求解矩阵主特征值（按模最大的特征值）及其对应特征向量的经典迭代方法，尤其适用于大型稀疏矩阵。其核心思想是通过构造向量序列的迭代过程逼近主特征方向，但收敛速度受次大特征值与主特征值的模长比值影响较大，且只能求取最大或最小、或最接近某个值的一个特征值。
幂法的基本原理
设矩阵
A
∈
R
n
×
n
A \in \mathbb{R}^{n \times n}
A∈Rn×n 有
n
n
n 个线性无关的特征向量
x
1
,
x
2
,
…
,
x
n
x_1, x_2, \dots, x_n
x1,x2,…,xn，对应特征值为
λ
1
,
λ
2
,
…
,
λ
n
\lambda_1, \lambda_2, \dots, \lambda_n
λ1,λ2,…,λn，且满足
∣
λ
1
∣
>
∣
λ
2
∣
≥
⋯
≥
∣
λ
n
∣
|\lambda_1| > |\lambda_2| \geq \cdots \geq |\lambda_n|
∣λ1∣>∣λ2∣≥⋯≥∣λn∣。任意初始向量
v
0
v_0
v0 可表示为：
v
0
=
∑
i
=
1
n
α
i
x
i
(
α
1
≠
0
)
v_0 = \sum_{i=1}^n \alpha_i x_i \quad (\alpha_1 \neq 0)
v0=∑i=1nαixi(α1=0). 通过迭代
v
k
=
A
v
k
−
1
v_k = A v_{k-1}
vk=Avk−1，得：
v
k
=
A
k
v
0
=
∑
i
=
1
n
α
i
λ
i
k
x
i
=
λ
1
k
(
α
1
x
1
+
∑
i
=
2
n
α
i
(
λ
i
λ
1
)
k
x
i
)
v_k = A^k v_0 = \sum_{i=1}^n \alpha_i \lambda_i^k x_i = \lambda_1^k \left( \alpha_1 x_1 + \sum_{i=2}^n \alpha_i \left( \frac{\lambda_i}{\lambda_1} \right)^k x_i \right)
vk=Akv0=i=1∑nαiλikxi=λ1k(α1x1+i=2∑nαi(λ1λi)kxi) 当
k
→
∞
k \to \infty
k→∞，
(
λ
i
λ
1
)
k
→
0
\left( \frac{\lambda_i}{\lambda_1} \right)^k \to 0
(λ1λi)k→0，故
v
k
v_k
vk 趋近于
λ
1
k
α
1
x
1
\lambda_1^k \alpha_1 x_1
λ1kα1x1。通过归一化（如取最大分量）可提取主特征向量
x
1
x_1
x1，而主特征值
λ
1
\lambda_1
λ1 通过相邻迭代向量分量的比值确定:
lim
n
→
∞
(
v
k
+
1
)
i
(
v
k
)
i
=
λ
1
\lim_{ n \to \infty } \frac{\left({v_{k+1}}\right)_{i}}{\left({v_{k}}\right)_{i}} = \lambda_{1}
n→∞lim(vk)i(vk+1)i=λ1
幂法在迭代过程中若不进行规范化处理，当主特征值
∣
λ
1
∣
≠
1
|\lambda_1| \neq 1
∣λ1∣=1时，向量序列
v
k
v_k
vk的模长会以指数速度增长（
∣
λ
1
∣
>
1
|\lambda_1|>1
∣λ1∣>1）或衰减（
∣
λ
1
∣
<
1
|\lambda_1|<1
∣λ1∣<1），导致数值溢出或精度丢失。为解决这一问题，需引入规范化操作，核心策略是通过对每次迭代后的向量进行尺度缩放，使其保持合理的数值范围。
u
k
=
v
k
∥
v
k
∥
∞
,
v
k
+
1
=
A
u
k
=
v
k
+
1
∥
v
k
∥
∞
u_k = \frac{v_k}{\|v_k\|_\infty},\: v_{k+1} = Au_{k} = \frac{v_{k+1}}{\|{v_{k}}\| _{\infty}}
uk=∥vk∥∞vk,vk+1=Auk=∥vk∥∞vk+1 从而：
u
k
→
x
1
max
(
x
1
)
max
(
v
k
)
→
λ
1
\begin{align*} u_{k}&\to \frac{x_{1}}{\max\: (x_{1})}\\ \max\: (v_{k})&\to \lambda_{1} \end{align*}
ukmax(vk)→max(x1)x1→λ1
幂法的加速策略
幂法的收敛速度由比值
r
=
∣
λ
2
/
λ
1
∣
r = |\lambda_2 / \lambda_1|
r=∣λ2/λ1∣ 决定。当
r
≈
1
r \approx 1
r≈1 时收敛缓慢，需引入加速技术。
原点平移法
构造矩阵
B
=
A
−
p
I
B = A - pI
B=A−pI，其特征值为
μ
i
=
λ
i
−
p
\mu_i = \lambda_i - p
μi=λi−p。选择适当的平移量
p
p
p，使得
∣
μ
1
∣
≫
∣
μ
i
∣
|\mu_1| \gg |\mu_i|
∣μ1∣≫∣μi∣（
i
≥
2
i \geq 2
i≥2），从而加速收敛。例如，当
A
A
A 的特征值为实数且满足
λ
1
>
λ
2
≥
⋯
≥
λ
n
\lambda_1 > \lambda_2 \geq \cdots \geq \lambda_n
λ1>λ2≥⋯≥λn，取
p
=
(
λ
2
+
λ
n
)
/
2
p = (\lambda_2 + \lambda_n)/2
p=(λ2+λn)/2 可最大化
∣
μ
1
∣
/
∣
μ
i
∣
|\mu_1| / |\mu_i|
∣μ1∣/∣μi∣，显著提升收敛速度。
Rayleigh商加速
设
A
A
A为实对称矩阵，对于任意非零向量，定义其Rayleigh商为：
R
A
(
x
)
:
=
x
T
A
x
x
T
x
R_{A}(x) := \frac{x^{T}Ax}{x^{T}x}
RA(x):=xTxxTAx 由于是实对称矩阵，其一定可以相似对角化，也就是有
n
n
n个相互独立的特征向量。于是
x
x
x可以写成：
x
=
∑
α
i
v
i
x =\sum \alpha_{i} v_{i}
x=∑αivi 从而很好证明：
1.
λ
n
⩽
(
A
x
,
x
)
(
x
,
x
)
⩽
λ
1
(
对于任何非零
x
∈
R
n
)
;
2.
λ
1
=
max
x
∈
R
n
(
A
x
,
x
)
(
x
,
x
)
;
3.
λ
n
=
min
x
∈
R
n
(
A
x
,
x
)
(
x
,
x
)
.
\begin{aligned} & 1.\quad \lambda_n\leqslant\frac{(Ax,x)}{(x,x)}\leqslant\lambda_1\quad(\text{对于任何非零 }x\in\mathbb{R}^n); \\ & 2.\quad \lambda_1=\underset{x\in\mathbb{R}^n}{\operatorname*{\operatorname*{max}}}\frac{(Ax,x)}{(x,x)}; \\ & 3.\quad \lambda_{n}=\min_{x\in\mathbf{R}^n}\frac{(Ax,x)}{(x,x)}. \end{aligned}
1.λn⩽(x,x)(Ax,x)⩽λ1(对于任何非零 x∈Rn);2.λ1=x∈Rnmax(x,x)(Ax,x);3.λn=x∈Rnmin(x,x)(Ax,x). 我们在采用幂法求特征值时，常常精度都是在
k
k
k阶。我们采用Rayleigh商，可以证明：
R
A
(
u
k
)
=
λ
1
+
O
(
(
λ
2
λ
1
)
2
k
)
R_{A}(u_{k}) = \lambda_{1} + O\left({(\frac{\lambda_{2}}{\lambda_{1}})^{2k}}\right)
RA(uk)=λ1+O((λ1λ2)2k)
Aitken加速
对线性收敛序列
{
a
k
}
\{a_k\}
{ak}，应用Aitken外推公式：
a
k
′
=
a
k
−
(
a
k
+
1
−
a
k
)
2
a
k
+
2
−
2
a
k
+
1
+
a
k
a_k' = a_k - \frac{(a_{k+1} - a_k)^2}{a_{k+2} - 2a_{k+1} + a_k}
ak′=ak−ak+2−2ak+1+ak(ak+1−ak)2 可消除线性误差项，加速逼近极限值。此方法无需额外矩阵运算，适用于幂法迭代序列的加速。
反幂法
反幂法用于求解矩阵的最小模特征值
λ
n
\lambda_n
λn 及对应特征向量，其核心是对
A
−
1
A^{-1}
A−1 应用幂法。因
A
−
1
A^{-1}
A−1 的特征值为
1
/
λ
n
≥
1
/
λ
n
−
1
≥
⋯
≥
1
/
λ
1
1/\lambda_n \geq 1/\lambda_{n-1} \geq \cdots \geq 1/\lambda_1
1/λn≥1/λn−1≥⋯≥1/λ1，反幂法迭代格式为：
A
v
k
=
μ
k
−
1
,
μ
k
=
v
k
∥
v
k
∥
A v_k = \mu_{k-1}, \quad \mu_k = \frac{v_k}{\|v_k\|}
Avk=μk−1,μk=∥vk∥vk 每次迭代需解线性方程组
A
v
k
=
μ
k
−
1
A v_k = \mu_{k-1}
Avk=μk−1，通常通过LU分解提高效率。
结合原点平移法，可计算接近某给定值
p
p
p 的特征值：对
B
=
A
−
p
I
B = A - pI
B=A−pI 应用反幂法，若
p
≈
λ
j
p \approx \lambda_j
p≈λj，则
B
−
1
B^{-1}
B−1 的主特征值为
1
/
(
λ
j
−
p
)
1/(\lambda_j - p)
1/(λj−p)，从而快速逼近
λ
j
\lambda_j
λj 及其特征向量。
QR算法
Householder方法以矩阵的正交相似变换为基础，可以用来求解矩阵全部的特征值和特征向量。
QR算法的核心是对矩阵进行QR分解：
A
=
:
A
1
=
Q
1
R
1
A=:A_{1} = Q_{1}R_{1}
A=:A1=Q1R1 其中
Q
1
,
R
1
Q_{1},R_{1}
Q1,R1分别是正交矩阵和上三角矩阵。这个过程中有三种典型方法，分别是Householder方法，Schimidt正交化，Givens矩阵法。三种分解方法的时间复杂度从小到大依次为： 1.Householder方法 (2/3*n^3) 2.Gram-Schmidt正交化 (3/3*n^3) 3.Givens矩阵 (4/3*n^3) 复杂度分析来自知乎回答QR分解的三种实现方法（Gram schmidt等）各自有什么优势和劣势？ - vonZooming的回答 - 知乎，本人未验证。
但是，如果在进行QR分解前使用Householer方法将
A
A
A正交相似化为Hessenberg矩阵，其QR分解的复杂度一下就降到了
O
(
n
2
)
O(n^{2})
O(n2)，大大加快了分解速度。
进行QR分解之后进行迭代：
A
k
=
Q
k
R
k
⟹
A
k
+
1
:
=
R
k
Q
k
=
Q
k
+
1
R
k
+
1
A_{k} = Q_{k}R_{k}\implies A_{k+1}:=R_{k}Q_{k} = Q_{k+1}R_{k+1}
Ak=QkRk⟹Ak+1:=RkQk=Qk+1Rk+1 其中有着一些性质：
A
k
+
1
=
Q
k
T
A
k
Q
k
=
Q
~
k
T
A
1
Q
~
k
A_{k+1} = Q_{k}^{T}A_{k}Q_{k} = \tilde{Q}_{k}^{T}A_{1}\tilde{Q}_{k}
Ak+1=QkTAkQk=Q~kTA1Q~k 其中
Q
~
k
:
=
∏
i
=
1
k
Q
k
\tilde{Q}_{k}:= \prod_{i=1}^{k}Q_{k}
Q~k:=∏i=1kQk。此外还有：
A
k
=
Q
~
k
R
~
k
A^{k} = \tilde{Q}_{k}\tilde{R}_{k}
Ak=Q~kR~k 在这样不断迭代下，可以证明最终会收敛到
s
c
h
u
r
schur
schur标准型，于是特征值就可以从中直接读出。具体实现细节见下：
1. 正交变换：Householder变换与Givens变换
作用：通过正交变换（保持向量长度和夹角）对矩阵进行数值稳定的相似变换。
Householder变换：
反射变换，通过一个超平面反射将向量中的多个元素变为零。用途：批量零化矩阵的某一行或一列（例如将矩阵化为上Hessenberg形式）。公式：
H
=
I
−
2
v
v
T
v
T
v
H = I - 2\frac{vv^T}{v^Tv}
H=I−2vTvvvT，其中
v
v
v 为构造的反射向量。 Givens变换：
平面旋转，通过旋转消去特定位置的元素（如次对角线元素）。用途：逐次零化单个元素（例如QR分解中的上三角化）。公式：
G
(
i
,
j
,
θ
)
G(i,j,\theta)
G(i,j,θ) 在
(
i
,
j
)
(i,j)
(i,j) 平面的旋转矩阵。
联系：两者均为正交矩阵，用于保持矩阵特征值不变的前提下简化矩阵结构。
2. 上Hessenberg矩阵
定义：次对角线以下元素全为零（或接近零），形如：
A
=
[
∗
∗
⋯
∗
∗
∗
⋯
∗
0
∗
⋱
⋮
0
0
⋯
∗
]
A = \begin{bmatrix}* & * & \cdots & * \\* & * & \cdots & * \\ 0 & * & \ddots & \vdots \\ 0 & 0 & \cdots & * \end{bmatrix}
A=
∗∗00∗∗∗0⋯⋯⋱⋯∗∗⋮∗
构造方法：
对任意方阵
A
A
A，通过 Householder变换 进行相似变换（左乘
H
H
H 和右乘
H
T
H^T
HT），逐步将每列的次对角线以下元素零化。 意义：上Hessenberg矩阵保留了原矩阵的特征值，但结构更简单，使得后续QR算法的计算量从
O
(
n
3
)
O(n^3)
O(n3) 降为
O
(
n
2
)
O(n^2)
O(n2)。
3. 实Schur分解
目标：将实矩阵
A
A
A 分解为准上三角矩阵（实Schur形式）：
A
=
Q
T
Q
T
A = Q T Q^T
A=QTQT 其中
Q
Q
Q 为正交矩阵，
T
T
T 为准上三角矩阵（对角块为1×1或2×2实块，对应实特征值或共轭复特征对）。 关键点：
若
A
A
A 已为上Hessenberg矩阵，QR算法可高效逼近实Schur形式。2×2对角块处理复数特征值，避免引入复数运算。
4. QR算法
核心思想：通过迭代将矩阵收敛到实Schur形式，从而提取特征值。 步骤：
预处理：用Householder变换将
A
A
A 化为上Hessenberg矩阵
H
H
H。迭代：
对
H
H
H 进行QR分解：
H
=
Q
R
H = QR
H=QR（用Givens变换处理次对角线元素）。计算新矩阵：
H
′
=
R
Q
H' = RQ
H′=RQ（RQ乘积保持上Hessenberg结构）。重复直到收敛为准上三角矩阵（实Schur形式）。 加速技巧：
位移策略（Wilkinson位移）：加速对角块收敛。隐式QR算法：避免显式构造
Q
Q
Q，直接通过Givens旋转更新矩阵。
嘉兴这条老街，即将大变样！
再见，快乐大本营：那些年我们一起追过的经典综艺回顾