QR算法 · aiwiki.page

QR算法是数值线性代数中用于计算方矩阵的特征值的迭代方法，也可用于计算舒尔分解。其基本步骤是将矩阵分解为正交或酉因子 $Q$ 与上三角因子 $R$ ，然后将这两个因子按相反的顺序相乘。实际实现将这一原理与预先进行的矩阵约化、位移、隐式变换和降阶相结合。这些改进使QR算法成为求解一般稠密矩阵全部特征值问题的标准方法。(cs.cornell.edu)

基本迭代

设 $A_0=A$ 为一个 $n\times n$ 矩阵。无位移QR迭代为

A_k=Q_kR_k,\qquad A_{k+1}=R_kQ_k,

其中第一个等式是QR分解。对于实矩阵， $Q_k$ 是正交矩阵；对于复矩阵，它是酉矩阵。由于 $Q_k^*Q_k=I$ ，有

A_{k+1}=Q_k^*A_kQ_k,

其中 $^{*}$ 表示共轭转置。因此，每一步都是一次矩阵相似变换，在精确运算下保持特征值不变。(cs.cornell.edu)

若将各步变换累积为

Z_k=Q_0Q_1\cdots Q_{k-1},

则有 $A_k=Z_k^*AZ_k$ 。在适当的谱条件和初始子空间条件下，这一过程趋近于舒尔形。无位移迭代与正交迭代密切相关：它隐式地执行同时子空间迭代，而不是逐个求取特征向量。其收敛速度取决于特征值模长的比值；当模长相同时，收敛可能很慢，也可能无法得到三角形的极限矩阵。(cs.cornell.edu)

舒尔形与特征向量

对于复矩阵，所求的分解为

A=ZTZ^*,

其中 $Z$ 为酉矩阵， $T$ 为上三角矩阵。特征值就是 $T$ 的对角元素。对于实矩阵，使用实数运算可得到实舒尔形： $Z$ 为正交矩阵， $T$ 为准上三角矩阵，其中 $1\times1$ 块表示实特征值， $2\times2$ 块表示共轭复特征值对。(netlib.org)

$Z$ 的各列称为舒尔向量，一般并不分别对应特征向量，而是为适当的不变子空间提供标准正交基。要求得右特征向量，可利用三角或分块三角结构求解

Ty=\lambda y

然后计算 $x=Zy$ 。LAPACK提供了专门完成这一计算的例程。因此，更准确地说，QR算法是一种求舒尔形的方法，而非直接进行矩阵对角化的方法。(netlib.org)

海森伯格约化

若每次迭代都对一般稠密矩阵进行分解，每一步需要 $O(n^3)$ 次运算。因此，实际使用的QR算法首先将 $A$ 约化为上海森伯格矩阵，即上海森伯格矩阵：

H=U^*AU,\qquad h_{ij}=0\quad\text{当 }i>j+1.

这一预处理需要 $O(n^3)$ 次运算，通常采用豪斯霍尔德变换。恰当实现的QR步骤能够保持海森伯格结构，并将一次单移位或双移位扫描的运算量降至 $O(n^2)$ 。(cs.cornell.edu)

对于实对称矩阵，对称的海森伯格矩阵是三对角矩阵：只有主对角线及其两侧相邻对角线上的元素可能非零。这一更强的结构使迭代成本大幅降低，也支持专门的对称矩阵特征值求解算法。(cs.cornell.edu)

位移

带位移的QR步骤将基本迭代替换为

A_k-\mu_kI=Q_kR_k, \qquad A_{k+1}=R_kQ_k+\mu_kI,

其中 $I$ 是单位矩阵， $\mu_k$ 是选定的标量。相似变换恒等式仍然成立：

A_{k+1}=Q_k^*A_kQ_k.

接近某个特征值的位移能够加快该特征值从其余尚未收敛的矩阵部分中分离出来。因此，位移的选择是实际算法收敛的关键。(cs.cornell.edu)

对于对称三对角矩阵，广泛使用的一种选择是威尔金森位移：取右下角 $2\times2$ 主子矩阵中最接近右下角对角元素的特征值。适当处理两个特征值距离相等的情形后，这一策略通常能使末端次对角元素迅速趋于零。其局部收敛往往是三次的，但并非总能达到三次收敛；某些特殊的对称矩阵仅表现出严格的二次收敛。(arxiv.org)

弗朗西斯双移位

实非对称矩阵可能具有非实特征值，因此仅采用实数单移位策略并不足够。弗朗西斯双移位步骤将两个位移结合起来，这两个位移通常取自右下角的 $2\times2$ 块。如果它们是一对共轭数 $\mu,\bar{\mu}$ ，相应的二次多项式为

p(z)=(z-\mu)(z-\bar{\mu}) =z^2-2\operatorname{Re}(\mu)z+|\mu|^2.

其系数为实数，因此组合后的变换可以完全使用实数运算完成。从概念上说，它使用 $p(H)$ 的正交因子，但高效实现并不会显式构造完整的矩阵多项式。(cs.cornell.edu)

隐式QR与凸起追赶

隐式QR算法无需显式分解整个带位移矩阵，即可构造所需的相似变换。对于未约化的海森伯格矩阵，即所有次对角元素均非零的海森伯格矩阵，隐式Q定理说明：变换矩阵的第一列加上恢复海森伯格形的要求，就能确定该变换，所剩的自由度仅是符号或相位约定。(cs.cornell.edu)

在双移位步骤中，根据 $p(H)$ 的第一列构造一个小规模的豪斯霍尔德变换。从左右两侧施加该变换，会在海森伯格带的下方产生一小簇元素，称为凸起。随后的一系列局部变换将凸起向下移动，直至其移出矩阵。这一过程称为凸起追赶，它在保留海森伯格结构所带来的低运算成本的同时，实现了预期的带位移相似变换。(cs.cornell.edu)

降阶与终止

降阶将谱中已经收敛的部分从剩余问题中分离出来。如果某个海森伯格次对角元素变为精确的零，矩阵就会分裂成分块上三角形：

H= \begin{pmatrix} H_{11}&H_{12}\\ 0&H_{22} \end{pmatrix}.

此时，整个矩阵的特征值就是两个对角块的特征值。在浮点运算中，可以通过考虑数值尺度的判据，将足够小的次对角元素置零。生产级实现除了简单的绝对阈值外，还包含其他保障措施。(cs.cornell.edu)

对于实非对称问题，表示一对共轭复特征值的孤立 $2\times2$ 块就是有效的已收敛块，无需进一步变成对角形。算法成功终止时，会得到相应的三角或准三角舒尔形。当常规位移选择导致长时间停滞时，实际实现还会使用特殊位移。(netlib.org)

对称问题与计算成本

对于实对称矩阵，QR过程首先通过正交变换将其约化为对称三对角形。如果只需特征值，在三对角表示上进行一次隐式QR扫描的成本为 $O(n)$ 。假设每个被降阶分离出的特征值大致只需要有界次数的扫描，则求取三对角矩阵的全部特征值需要 $O(n^2)$ 次运算；但一般稠密对称矩阵的初始约化仍需要 $O(n^3)$ 次运算。累积计算全部特征向量会增加大量工作。(cs.cornell.edu)

对于一般稠密非对称矩阵，通常将初始海森伯格约化与后续QR迭代的总运算量估计为 $O(n^3)$ 。这一估计假定每个特征值或特征值对所需的平均扫描次数较少，并不意味着对所有输入都存在统一的次数上界。实际迭代次数取决于矩阵和位移策略。存储完整稠密矩阵需要 $O(n^2)$ 的空间。(netlib.org)

针对对称问题的专用替代算法包括分治法特征值求解器、雅可比方法，以及结合逆迭代的二分法。根据所需的特征值或特征向量，这些方法可能比QR算法更合适。(cs.cornell.edu)

数值稳定性与局限

QR算法的主要数值优势在于使用正交或酉变换。正确实现的隐式QR迭代，在适当控制降阶的情况下，具有后向稳定性：计算得到的舒尔分解与输入矩阵受到微小扰动后的精确分解非常接近。这是关于数值稳定性的结论，并不保证每个计算出的特征值都具有很小的相对误差。(cs.cornell.edu)

特征值和特征向量的精度还取决于原问题的敏感性，可用适当的条件数来量化。因此，LAPACK提供了专门用于估计特征值—特征向量对及不变子空间条件性的工具。稳定的算法无法消除输入矩阵本身固有的敏感性。(netlib.org)

在海森伯格约化之前，可以选择对非对称问题进行矩阵平衡，以重新调整其数值尺度。这可能提高精度，但不加区分地进行平衡也可能增大特征向量误差；缩放及其逆操作都需要谨慎处理。(arxiv.org)

实现与相关方法

LAPACK将非对称矩阵的计算分为若干阶段：xGEHRD 将矩阵约化为海森伯格形，xHSEQR 计算其特征值，并可选择计算舒尔形，xTREVC 等例程则求取特征向量。前缀 x 表示一组适用于不同数值类型和精度的例程。(netlib.org)

高级实现使用小凸起多移位QR，在处理多个位移的同时，将更新组织为高效的分块矩阵运算。积极提前降阶通过检查矩阵右下角的一个窗口，在常规次对角元素判据能够检测到收敛之前，就识别出已收敛的特征值。这些技术提高了实际性能，但并未改变底层的相似变换原理。(netlib.org)

相关的QR迭代也用于奇异值分解，在矩阵约化为双对角形之后执行。它们通过从左右两侧施加变换来计算奇异值，而不是通过相似变换计算特征值。相关的QZ算法则通过广义舒尔形处理矩阵对的广义特征值问题。(cs.cornell.edu)

历史

约翰·G. F. 弗朗西斯与薇拉·库布拉诺夫斯卡娅在1959—1961年前后独立提出了QR方法。弗朗西斯于1959年10月29日提交了他的第一篇理论论文；库布拉诺夫斯卡娅于1960年7月5日提交了最初的概要。他们的早期论文确立了QR作为求解矩阵全部特征值问题的方法，其中弗朗西斯的工作尤其与隐式双移位构造密切相关。(atm.org.uk)

后续发展主要集中在可靠的位移策略、更严格的降阶判据、多移位扫描以及高效的分块实现上。LAPACK 3.1纳入了结合积极提前降阶的小凸起多移位QR，将原始方法发展为结构复杂得多的生产级特征值求解器。(netlib.org)

参考来源

CS 4220: Numerical Analysis — March 6, 2026cs.cornell.edu
A Parallel QR Algorithm for the Nonsymmetric Eigenvalue Problemnetlib.org
LAPACK: CHSEQRnetlib.org
LAPACK Users' Guidenetlib.org
CS 6210: Matrix Computations — October 29, 2025cs.cornell.edu
CS 6210: Matrix Computations — November 5, 2025cs.cornell.edu
The Asymptotics of Wilkinson's Iteration: Loss of Cubic Convergencearxiv.org
LAPACK 3.1 xHSEQR: Tuning and Implementation Notes on the Small Bulge Multi-shift QR Algorithm with Aggressive Early Deflationnetlib.org
A Parallel Implementation of the Nonsymmetric QR Algorithmnetlib.org
LAPACK Working Note 13netlib.org
On Matrix Balancing and Eigenvector Computationarxiv.org
The QR Algorithm: 50 Years Later—Its Genesis by John Francis and Vera Kublanovskaya and Subsequent Developmentsatm.org.uk