线性注意力的线性代数前置知识

发表于 2026-08-28 更新于 2026-09-05 分类于技术

DeltaNet / KDA 这一系的 chunkwise 算法，数学上就建在三件事上：Householder 变换、WY 表示、以及把串行递推改写成闭式的 UT 变换。本文把这三层从定义推到可用于 kernel 的形式，并附两个可交互演示（广义 Householder 的 $\beta$ 系数、单位下三角求逆的前代法）。后续四篇 TileLang 实战可以直接回查本文的结论。

1. Householder 矩阵

1.1 定义

设 $\bm{u} \in \mathbb{R}^n$ 且 $\bm{u}^{\intercal}\bm{u} = 1$ （即 $\|\bm{u}\|_2 = 1$ ），定义

\mathbf{H} = \mathbf{I} - 2\bm{u}\bm{u}^{\intercal} \in \mathbb{R}^{n \times n}

1.1.1 几何意义与矩阵形式

把 $\bm{x}$ 沿 $\bm{u}$ 分解：分量为 $(\bm{u}^{\intercal}\bm{x})\bm{u}$ 。 $\mathbf{H}$ 的作用就是把这个分量减去两份：减一份落在法平面上（投影），减两份则穿过平面到对称位置，即关于法平面的镜像。平面内的分量不变，所以长度也不变。取 $\bm{u} = \bm{e}_1$ 、 $\bm{x} = [3, 2]^{\intercal}$ ：

写成式子就是 $\bm{x}' = \bm{x} - 2(\bm{u}^{\intercal}\bm{x})\bm{u}$ 。要把它整理成矩阵作用，关键是 $\bm{u}^{\intercal}\bm{x}$ 是标量，可与向量换序， $\bm{x}$ 于是变成右端公因子：

\bm{x}' = \bm{x} - 2(\bm{u}^{\intercal}\bm{x})\bm{u} = \bm{x} - 2\bm{u}(\bm{u}^{\intercal}\bm{x}) = \big(\mathbf{I} - 2\bm{u}\bm{u}^{\intercal}\big)\bm{x} = \mathbf{H}\bm{x}

$(\bm{u}^{\intercal}\bm{x})\bm{u}$ 里 $\bm{x}$ 夹在中间提不出来，换序后 $\bm{u}\bm{u}^{\intercal}$ 自然聚成一个矩阵。

1.2 通用形式：把归一化收进一个系数 $\tau$

$\mathbf{H} = \mathbf{I} - 2\bm{u}\bm{u}^{\intercal}$ 要求 $\|\bm{u}\|_2 = 1$ 。要支持任意非零 $\bm{u}$ ，只需除掉它的长度平方 $\bm{u}^{\intercal}\bm{u}$ ——外积 $\bm{u}\bm{u}^{\intercal}$ 对 $\bm{u}$ 是二次齐次的，除掉这个标量后结果只依赖 $\bm{u}$ 的方向，与长度无关：

\mathbf{H} = \mathbf{I} - \frac{2\bm{u}\bm{u}^{\intercal}}{\bm{u}^{\intercal}\bm{u}} = \mathbf{I} - \tau\,\bm{u}\bm{u}^{\intercal}, \qquad \tau := \frac{2}{\bm{u}^{\intercal}\bm{u}}

这是文献里最常见的形式。 $\|\bm{u}\|_2 = 1$ 时 $\tau = 2$ ，退回 §1.1 的定义。

这个系数的记法各家不一，下面几种都是同一件事：

写法	系数取值
$\mathbf{I} - \tau\bm{u}\bm{u}^{\intercal}$	$\tau = 2/(\bm{u}^{\intercal}\bm{u})$
$\mathbf{I} - \bm{u}\bm{u}^{\intercal}/\tau$	$\tau = (\bm{u}^{\intercal}\bm{u})/2$

记法不重要，知道 $\tau$ 只是外积前面的一个标量系数就够了——后面所有推导只用到「它是标量、可以随意搬动」这一条性质。

系数取什么值，几何上就完全不同。令 $k := \tau\,\bm{u}^{\intercal}\bm{u}$ （上表第一种记法下 $k = \beta\|\bm{u}\|^2$ ），则 $\mathbf{H}$ 沿 $\bm{u}$ 方向的特征值恰是 $1 - k$ ，而法平面内 $n-1$ 个方向的特征值恒为 1。于是三个整数点各对应一件事：

$k$	$\mathbf{H}$	作用
0	$\mathbf{I}$	恒等，什么都不做
1	$\mathbf{I} - \bm{u}\bm{u}^{\intercal}/(\bm{u}^{\intercal}\bm{u})$	正交投影到法平面（ $\mathbf{H}\bm{x}$ 落在垂足上）
2	$\mathbf{I} - 2\bm{u}\bm{u}^{\intercal}/(\bm{u}^{\intercal}\bm{u})$	镜像反射，等长、 $\mathbf{H}^2 = \mathbf{I}$

拖动下面的琥珀色手柄旋转 $\bm{u}$ （长度锁定，镜面随之转动），黑色手柄移动 $\bm{x}$ ，滑块调 $\beta$ 。 $k$ 从 0 连续走到 2 的过程，就是「不动 → 压到镜面上 → 穿过镜面」；继续加大到 $k > 2$ 会看到法向反向放大。这也是后文一直只要求「 $\tau$ 是标量」而不限定它取 2 的原因——DeltaNet 里的 $\beta_t \in (0,1)$ 就落在压缩那一段，从来不是真正的反射。

2. WY 表达式

结论先给： $p$ 个 Householder 反射的乘积，永远能写成单位矩阵减一个秩不超过 $p$ 的修正

\mathbf{H}_1\mathbf{H}_2\cdots\mathbf{H}_p = \mathbf{I} - \mathbf{W}\mathbf{Y}^{\intercal}, \qquad \mathbf{W},\,\mathbf{Y} \in \mathbb{R}^{n \times p}

证明用归纳构造，一列一列把反射「吸收」进 $\mathbf{W},\mathbf{Y}$ 。

2.1 基例 $p = 1$

$\mathbf{H}_1 = \mathbf{I} - 2\bm{u}_1\bm{u}_1^{\intercal} = \mathbf{I} - \mathbf{W}_1\mathbf{Y}_1^{\intercal}$ ，取

\mathbf{W}_1 = \big[2\bm{u}_1\big], \qquad \mathbf{Y}_1 = \big[\bm{u}_1\big]

均为 $n \times 1$ 。命题对 $p = 1$ 成立。

2.2 归纳步：一列一列「吸收」反射

归纳假设：已存在 $\mathbf{W}_k,\,\mathbf{Y}_k \in \mathbb{R}^{n \times k}$ 使 $\mathbf{H}_1\cdots\mathbf{H}_k = \mathbf{I} - \mathbf{W}_k\mathbf{Y}_k^{\intercal}$ 。要证 $\mathbf{H}_1\cdots\mathbf{H}_k\mathbf{H}_{k+1}$ 也能写成 $\mathbf{I} - \mathbf{W}_{k+1}\mathbf{Y}_{k+1}^{\intercal}$ ，其中 $\mathbf{W}_{k+1},\mathbf{Y}_{k+1} \in \mathbb{R}^{n \times (k+1)}$ 。

记新反射方向 $\bm{u} = \bm{u}_{k+1}$ ，并记

\bm{w} := 2\,(\mathbf{I} - \mathbf{W}_k\mathbf{Y}_k^{\intercal})\,\bm{u}

则把 $(\mathbf{I} - \mathbf{W}_k\mathbf{Y}_k^{\intercal})(\mathbf{I} - 2\bm{u}\bm{u}^{\intercal})$ 乘开后的四项可以合并。

\mathbf{I} - \mathbf{W}_k\mathbf{Y}_k^{\intercal} - 2\bm{u}\bm{u}^{\intercal} + 2\mathbf{W}_k\mathbf{Y}_k^{\intercal}\bm{u}\bm{u}^{\intercal} = \mathbf{I} - \mathbf{W}_k\mathbf{Y}_k^{\intercal} - \bm{w}\bm{u}^{\intercal}

于是

\mathbf{H}_1\cdots\mathbf{H}_{k+1} = \mathbf{I} - \big[\,\mathbf{W}_k \;\; \bm{w}\,\big]\big[\,\mathbf{Y}_k \;\; \bm{u}\,\big]^{\intercal}

即 $\mathbf{W}_{k+1} = [\mathbf{W}_k,\ \bm{w}]$ ， $\mathbf{Y}_{k+1} = [\mathbf{Y}_k,\ \bm{u}_{k+1}]$ 。定理得证。 $\blacksquare$

3. compact WY：把 $\mathbf{W}$ 分解成 $\mathbf{Y}$ 与一个上三角 $\mathbf{T}$

由 §2.2 的递推回代可得 $\bm{w}_t = 2\mathbf{H}_1\cdots\mathbf{H}_{t-1}\bm{u}_t$ ，每根列都带着一个前缀乘积。

1989 年 Schreiber–Van Loan 的 compact WY 指出 $\mathbf{W}$ 一定能分解成 $\mathbf{Y}\mathbf{T}$ ，其中 $\mathbf{T}$ 是 $p \times p$ 上三角矩阵：

\mathbf{H}_1\cdots\mathbf{H}_p = \mathbf{I} - \mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal}, \qquad \mathbf{Y} = [\,\bm{u}_1, \dots, \bm{u}_p\,], \quad \mathbf{T} \in \mathbb{R}^{p \times p}\ \text{上三角}

两个 $n \times p$ 因子 $\mathbf{Y}$ 和 $\mathbf{Y}^{\intercal}$ 里于是只剩 $\bm{u}$ ， $p \ll n$ 时 $\mathbf{T}$ 的存储规模会很小。下面用 §1.2 的广义形式 $\mathbf{H} = \mathbf{I} - \tau\bm{u}\bm{u}^{\intercal}$ 写， $\bm{u}_t$ 不需归一， $\tau_t$ 也不受限制。

3.1 递推构造

\mathbf{T}_1 = [\tau_1], \qquad \mathbf{T}_{k+1} = \begin{bmatrix} \mathbf{T}_k & -\tau_{k+1}\,\mathbf{T}_k\,(\mathbf{Y}_k^{\intercal}\bm{u}_{k+1}) \\ \bm{0} & \tau_{k+1} \end{bmatrix}

证明与 §2.2 是同一个机制，只是新增的一列不往 $n$ 维因子里塞，而往 $\mathbf{T}$ 里塞。记 $\bm{u} = \bm{u}_{k+1}$ 、 $\tau = \tau_{k+1}$ ，乘开：

(\mathbf{I} - \mathbf{Y}_k\mathbf{T}_k\mathbf{Y}_k^{\intercal})(\mathbf{I} - \tau\bm{u}\bm{u}^{\intercal}) = \mathbf{I} - \mathbf{Y}_k\mathbf{T}_k\mathbf{Y}_k^{\intercal} - \tau\bm{u}\bm{u}^{\intercal} + \tau\,\mathbf{Y}_k\mathbf{T}_k(\mathbf{Y}_k^{\intercal}\bm{u})\,\bm{u}^{\intercal}

三个修正项正好是 $\big[\mathbf{Y}_k,\, \bm{u}\big]$ 夹一个 $(k+1) \times (k+1)$ 矩阵再夹 $\big[\mathbf{Y}_k,\, \bm{u}\big]^{\intercal}$ 的展开结果： $\mathbf{T}_k$ 占左上块， $\tau$ 占右下角，交叉项给出右上列（整体已提出一个负号，所以该元素带 $-$ ），左下为零：

= \mathbf{I} - \big[\mathbf{Y}_k,\, \bm{u}\big] \begin{bmatrix} \mathbf{T}_k & -\tau\mathbf{T}_k(\mathbf{Y}_k^{\intercal}\bm{u}) \\ \bm{0} & \tau \end{bmatrix} \big[\mathbf{Y}_k,\, \bm{u}\big]^{\intercal}

即 $\mathbf{Y}_{k+1} = [\mathbf{Y}_k,\ \bm{u}_{k+1}]$ 、 $\mathbf{T}_{k+1}$ 如上，定理得证。 $\blacksquare$

3.2 从递推读出的三个性质

$\mathbf{T}$ 上三角。左下块永远是 $\bm{0}$ ，上三角性逐层遗传。
对角线就是系数。 $\mathrm{diag}(\mathbf{T}) = (\tau_1, \dots, \tau_p)$ ，单位向量情形下全是 $2$ 。
与 §2 的关系是 $\mathbf{W} = \mathbf{Y}\mathbf{T}$ 。 $n \times p$ 的 $\mathbf{W}$ 换成了 $p \times p$ 的 $\mathbf{T}$ ，存储从 $O(np)$ 降到 $O(p^2)$ 。

3.3 $\mathbf{T}$ 的每个数就是一个外积的系数

§3.1 的递推说明了 $\mathbf{T}$ 怎么算，但没说 $\mathbf{T}$ 的每个数是什么。把 $\mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal}$ 硬展开一遍就清楚了，只需要三条基本式子。

式子一（行视角）：行向量乘矩阵 = 矩阵各行的线性组合。 $\mathbf{T}$ 是 $p \times p$ ，把 $\mathbf{Y}^{\intercal}$ 按行拆成 $[\bm{u}_1^{\intercal}; \cdots; \bm{u}_p^{\intercal}]$ ，用 $\mathbf{T}$ 的第 $l$ 行当系数：

(\mathbf{T}\mathbf{Y}^{\intercal})\ \text{的第 } l \text{ 行} = \sum_{j=1}^{p} \mathbf{T}_{lj}\,\bm{u}_j^{\intercal}

式子二（外积视角）： $\mathbf{A}\mathbf{B} = \sum_l (\mathbf{A}\ \text{第 } l \text{ 列})(\mathbf{B}\ \text{第 } l \text{ 行})$ ，即矩阵乘法是一堆秩一项之和。左乘 $\mathbf{Y} = [\bm{u}_1, \dots, \bm{u}_p]$ ：

\mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal} = \sum_{l=1}^{p} \bm{u}_l \Big(\sum_{j=1}^{p} \mathbf{T}_{lj}\bm{u}_j^{\intercal}\Big) = \sum_{i,j} \mathbf{T}_{ij}\,\bm{u}_i\bm{u}_j^{\intercal}

$\mathbf{T}$ 上三角意味着 $i > j$ 的项为零，展开只剩两类：

\mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal} = \underbrace{\tau_1\bm{u}_1\bm{u}_1^{\intercal} + \cdots + \tau_p\bm{u}_p\bm{u}_p^{\intercal}}_{\text{对角项}} + \underbrace{\mathbf{T}_{12}\bm{u}_1\bm{u}_2^{\intercal} + \mathbf{T}_{13}\bm{u}_1\bm{u}_3^{\intercal} + \cdots}_{\text{上三角项}}

结论： $\mathbf{T}_{ij}$ 就是外积 $\bm{u}_i\bm{u}_j^{\intercal}$ 的系数，整体前面还有一个负号： $\mathbf{H}_1\cdots\mathbf{H}_p = \mathbf{I} - \sum_{i \le j}\mathbf{T}_{ij}\,\bm{u}_i\bm{u}_j^{\intercal}$ 。

对角 $\mathbf{T}_{tt} = \tau_t$ ——第 $t$ 个反射自己的强度，与 §3.2 从递推读出的结论一致。
严格上三角 $\mathbf{T}_{ij}\ (i<j)$ ——反射 $i$ 与反射 $j$ 之间的累积耦合权重，即所有「先经过 $\bm{u}_i$ 再经过 $\bm{u}_j$ 」路径的加权和。

一句话总结：把 $\mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal}$ 按「行视角 + 外积分解」展开得到 $\sum_{i \le j}\mathbf{T}_{ij}\bm{u}_i\bm{u}_j^{\intercal}$ —— $\mathbf{T}$ 的第 $(i,j)$ 个数就是外积 $\bm{u}_i\bm{u}_j^{\intercal}$ 的系数。

4. UT 变换：把递推变成一次三角求解

WY、compact WY、UT 是同一个恒等式的三个阶段，改动的只有「右端各因子怎么算」：

\mathbf{H}_1\cdots\mathbf{H}_p = \mathbf{I} - \mathbf{W}\mathbf{Y}^{\intercal} = \mathbf{I} - \mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal}

每一步消掉上一步留下的一个具体瓶颈：

阶段	得到	消掉的瓶颈	新引入的代价
WY（§2，1987）	$p$ 个秩一修正聚成一个秩 $p$ 修正	反射只能逐个应用、无法块化	$\mathbf{W}$ 额外占 $O(np)$ ； $\bm{w}_t$ 依赖 $\bm{w}_{i<t}$ ，串行
compact WY（§3，1989）	$\mathbf{W} = \mathbf{Y}\mathbf{T}$	存储 $O(np) \to O(p^2)$ （ $\mathbf{Y}$ 本来就要存）	$\mathbf{T}$ 逐列递推，第 $k{+}1$ 列等第 $k$ 列， $p$ 步串行
UT（本节，Joffrain 等，2006）	$\mathbf{T}$ 的闭式	串行递推 $\to$ 1 次 Gram 矩阵 $+$ 1 次单位三角求解	多一次 $O(p^2)$ 回代

UT 变换要证的只有一件事：递推的串行依赖是实现引入的，不是数学必需—— $\mathbf{T}$ 的全部信息就是 $\mathbf{G} = \mathbf{Y}^{\intercal}\mathbf{Y}$ 加一组 $\tau$ ，而 $\mathbf{G}$ 的 $p^2$ 个内积互不依赖。仍用 §1.2 的广义形式 $\mathbf{H}_t = \mathbf{I} - \tau_t\bm{u}_t\bm{u}_t^{\intercal}$ ， $\tau_t$ 任意，推导对一般秩一修正同样成立，不要求正交性。

4.1 推导：把 WY 递推堆成一个方程

复习 §2.2 的归纳步：吸收第 $t$ 个反射时， $\mathbf{W}$ 的新列是定义出来的（§2.2 取 $\tau=2$ ，这里换成广义系数）：

\bm{w}_t := \tau_t\big(\mathbf{I} - \mathbf{W}_{t-1}\mathbf{Y}_{t-1}^{\intercal}\big)\,\bm{u}_t

由归纳假设 $\mathbf{I} - \mathbf{W}_{t-1}\mathbf{Y}_{t-1}^{\intercal} = \mathbf{H}_1\cdots\mathbf{H}_{t-1}$ ，它有等价写法 $\bm{w}_t = \tau_t\,\mathbf{H}_1\cdots\mathbf{H}_{t-1}\,\bm{u}_t$ （新列 = 前缀乘积作用在新的 $\bm{u}_t$ 上再乘系数）。把这个定义式乘开，得一条只含内积的递推：

\bm{w}_t = \tau_t\Big(\bm{u}_t - \sum_{i<t}\bm{w}_i\,\big(\bm{u}_i^{\intercal}\bm{u}_t\big)\Big) \;\Longleftrightarrow\; \bm{w}_t + \sum_{i<t}\tau_t\,\mathbf{G}_{it}\,\bm{w}_i = \tau_t\bm{u}_t, \qquad \mathbf{G} := \mathbf{Y}^{\intercal}\mathbf{Y}

（ $\mathbf{G}_{it} = \bm{u}_i^{\intercal}\bm{u}_t$ ，即第 $(i,t)$ 元是 $\bm{u}_i$ 与 $\bm{u}_t$ 的内积。 $\mathbf{D} := \mathrm{diag}(\tau_1,\dots,\tau_p)$ ， $\mathrm{diag}$ 表示以所列标量为对角元、其余位置全为零的对角矩阵，形状与 $\mathbf{G}$ 同为 $p\times p$ ，保证 $\mathbf{G}\mathbf{D}$ 可乘。）把后文用到的矩阵展开写清楚：

\mathbf{G} = \mathbf{Y}^{\intercal}\mathbf{Y} = \begin{bmatrix} \bm{u}_1^{\intercal}\bm{u}_1 & \cdots & \bm{u}_1^{\intercal}\bm{u}_p\\ \vdots & \ddots & \vdots\\ \bm{u}_p^{\intercal}\bm{u}_1 & \cdots & \bm{u}_p^{\intercal}\bm{u}_p \end{bmatrix} \in \mathbb{R}^{p \times p}, \qquad \mathbf{D} = \mathrm{diag}(\tau_1,\dots,\tau_p) = \begin{bmatrix} \tau_1 & & \\ & \ddots & \\ & & \tau_p \end{bmatrix} \in \mathbb{R}^{p \times p}

\mathbf{G}\,\mathbf{D} = \begin{bmatrix} \tau_1\,\bm{u}_1^{\intercal}\bm{u}_1 & \cdots & \tau_p\,\bm{u}_1^{\intercal}\bm{u}_p\\ \vdots & \ddots & \vdots\\ \tau_1\,\bm{u}_p^{\intercal}\bm{u}_1 & \cdots & \tau_p\,\bm{u}_p^{\intercal}\bm{u}_p \end{bmatrix}

\mathrm{striu}(\mathbf{G}\mathbf{D}) = \begin{bmatrix} 0 & \tau_2\,\bm{u}_1^{\intercal}\bm{u}_2 & \cdots & \tau_p\,\bm{u}_1^{\intercal}\bm{u}_p\\ & 0 & \ddots & \vdots\\ & & \ddots & \tau_p\,\bm{u}_{p-1}^{\intercal}\bm{u}_p\\ & & & 0 \end{bmatrix}, \qquad \mathrm{striu}(\mathbf{G})\ \text{同形，去掉各 } \tau_t

要点一条：右乘对角阵 $\mathbf{D}$ 是列缩放， $(\mathbf{G}\mathbf{D})_{it} = \mathbf{G}_{it}\,\tau_t$ ，第 $t$ 列统一乘 $\tau_t$ ——其 $(i,t)$ 元正是上面递推里 $\bm{w}_i$ 前的系数。由递推出发推导时缩放天然是列向的，因此全文只用 $\mathrm{striu}(\mathbf{G}\mathbf{D})$ 一种写法， $\mathrm{striu}$ 总是作用在缩放后的整个矩阵上。

（左乘 $\mathbf{D}$ 的行缩放 $\mathbf{D}\mathbf{G}$ 是等价的另一种缩放，经恒等式 $\mathbf{D}\big(\mathbf{I}+\mathrm{striu}(\mathbf{G}\mathbf{D})\big) = \big(\mathbf{I}+\mathrm{striu}(\mathbf{D}\mathbf{G})\big)\mathbf{D}$ 可互换；本文不引入第二种记号。）

每条方程只往回看 $i < t$ ，是严格三角依赖。逐列核对三项：

位置	第 $t$ 列	对应递推里的
$\mathbf{W}$	$\bm{w}_t$	待求项
$\mathbf{W}\,\mathrm{striu}(\mathbf{G}\mathbf{D})$	$\sum_{i<t}\mathbf{G}_{it}\tau_t\,\bm{w}_i$	左端第二项
$\mathbf{Y}\mathbf{D}$	$\tau_t\bm{u}_t$	右端

中间那行需要列视角： $(\mathbf{W}\mathbf{S})$ 的第 $t$ 列 $= \mathbf{W}\,\bm{s}_t$ ，是 $\mathbf{W}$ 的所有列以 $\mathbf{S}$ 第 $t$ 列元素为系数的线性组合。而 $\mathrm{striu}(\mathbf{G}\mathbf{D})$ 的第 $t$ 列在 $i \ge t$ 处为零，组合实际只含 $\bm{w}_{i<t}$ ——「新列只由旧列线性组合而成」就是递推串行结构的矩阵语言。

矩阵乘的另两种视角顺带补齐（ $\mathbf{A}\in\mathbb{R}^{m\times k}$ ， $\mathbf{S}\in\mathbb{R}^{k\times n}$ ）：

视角	式子	读法	本文实例
列	$(\mathbf{A}\mathbf{S})$ 第 $j$ 列 $= \mathbf{A}\bm{s}_j$	$\mathbf{A}$ 各列的线性组合	上面的逐列核对
行	$(\mathbf{A}\mathbf{S})$ 第 $i$ 行 $= \bm{a}_i^{\intercal}\mathbf{S}$	$\mathbf{S}$ 各行的线性组合	$\mathbf{G}$ 第 $i$ 行 $= \bm{u}_i^{\intercal}\mathbf{Y}$
外积	$\mathbf{A}\mathbf{S} = \sum_{j=1}^{k}\bm{a}_j\bm{s}_j^{\intercal}$	$k$ 个秩一项之和	$\mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal} = \sum_t \bm{w}_t\bm{u}_t^{\intercal}$

行视角那一行就是「Gram 矩阵一行 = 一个向量对全体的内积表」。外积视角完全不走线性组合的路子，且直接给出一个关键事实：

\mathbf{I}-\mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal} = \mathbf{I}-\sum_{t=1}^{p}\bm{w}_t\bm{u}_t^{\intercal}

右端是 $p$ 个秩一修正之和，与单个 $\mathbf{H}_t = \mathbf{I}-\tau_t\bm{u}_t\bm{u}_t^{\intercal}$ 的秩一结构同构。连乘的「积」变成了秩一项的「和」——这就是 $n \times p$ 的 $\mathbf{W},\mathbf{Y}$ 能承载全部信息的直观内容。

于是 $p$ 条递推合成一个矩阵方程：

\mathbf{W}\big(\mathbf{I} + \mathrm{striu}(\mathbf{G}\mathbf{D})\big) = \mathbf{Y}\mathbf{D}

4.2 解出 $\mathbf{T}$ ：闭式

系数矩阵 $\mathbf{M} := \mathbf{I} + \mathrm{striu}(\mathbf{G}\mathbf{D})$ 是单位上三角：对角元为 1（ $\bm{w}_t$ 完整出现在自己的方程里）， $i > t$ 处为零。对角恒非零，故必可逆、无需选主元——「UT」的 U 就是 unit triangular。于是

\mathbf{W} = \mathbf{Y}\,\underbrace{\mathbf{D}\,\mathbf{M}^{-1}}_{=:\ \mathbf{T}} \qquad\Longrightarrow\qquad \mathbf{H}_1\cdots\mathbf{H}_p = \mathbf{I} - \mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal}, \qquad \mathbf{T} = \mathbf{D}\,\big(\mathbf{I} + \mathrm{striu}(\mathbf{G}\mathbf{D})\big)^{-1}

对右边这个 $\mathbf{T}$ 求一次逆可以顺手验证它与递推自洽： $\mathbf{T}^{-1} = \mathbf{M}\mathbf{D}^{-1} = \mathbf{D}^{-1} + \mathrm{striu}(\mathbf{G}\mathbf{D})\,\mathbf{D}^{-1} = \mathbf{D}^{-1} + \mathrm{striu}(\mathbf{G})$ （列缩放再逐列除回， $\tau_t$ 恰好消掉），即

\mathbf{T}^{-1} = \mathbf{D}^{-1} + \mathrm{striu}(\mathbf{G})

这一步同时回收前两节： $\mathbf{W} = \mathbf{Y}\mathbf{T}$ 不再是 §3 的构造性分解，而是方程的解； $\mathbf{T}$ 上三角（单位上三角之逆仍上三角）、 $\mathrm{diag}(\mathbf{T}) = \mathrm{diag}(\mathbf{D}) = (\tau_1,\dots,\tau_p)$ ，与 §3.1 的递推逐位一致。逐列递推只是解这个三角方程的一种（串行的）实现方式。

求逆本身便宜： $\mathbf{N} := \mathrm{striu}(\mathbf{G}\mathbf{D})$ 严格三角因而幂零（ $\mathbf{N}^p = \bm{0}$ ），Neumann 级数 $(\mathbf{I}+\mathbf{N})^{-1} = \mathbf{I}-\mathbf{N}+\mathbf{N}^2-\cdots$ 在 $p$ 项内精确截断；数值上就是一次回代， $O(p^2)$ 。这两句话其实是同一件事，§4.5 展开说明。

4.3 为什么这个形式适合并行

$\mathbf{G} = \mathbf{Y}^{\intercal}\mathbf{Y}$ 的 $p^2$ 个内积互不依赖，一次 GEMM 算完，是 Tensor Core 最擅长的形状；
剩下的三角求解只在 $p \times p$ 规模上进行， $p$ 取 64 或 128 时整个 $\mathbf{T}$ 放得进一个线程块的 shared memory；
对照 §3 的递推：第 $k+1$ 列必须等第 $k$ 列，只能用到单个线程。

分块也能直接拼： $\mathbf{Y} = [\mathbf{Y}_1, \mathbf{Y}_2]$ 时

\mathbf{T} = \begin{bmatrix} \mathbf{T}_1 & -\mathbf{T}_1(\mathbf{Y}_1^{\intercal}\mathbf{Y}_2)\mathbf{T}_2 \\ \bm{0} & \mathbf{T}_2 \end{bmatrix}

块内 $\mathbf{T}_1,\mathbf{T}_2$ 各自并行算，跨块耦合项一次 GEMM 拼上——序列建模里的 chunkwise 并行就是这个结构：块内各算 $\mathbf{T}$ ，块间靠耦合项拼接。

方向约定：上三角对应 $\mathbf{H}_1\cdots\mathbf{H}_p$ 这个次序，反序乘积得到的是下三角那一半；等价地，把 $\bm{u}_t$ 按行堆成 $\mathbf{K}$ （而非按列堆成 $\mathbf{Y}$ ）也会翻成下三角。两件事都只是记号，但弄反会得到错的矩阵，是实现里最容易踩的坑。

4.4 前代法求逆：为什么「一次回代」就够

§4.2 说求 $\mathbf{M}^{-1}$ 「数值上就是一次回代， $O(p^2)$ 」，这一节把这句话摊开。下面统一用单位下三角的 $\mathbf{L}$ 来写（对角恒为 1，严格下三角是数据），上三角情形把顺序倒过来即可。

第一步：把求逆变成解 $n$ 个方程组。 求 $\mathbf{X} = \mathbf{L}^{-1}$ 等价于解 $\mathbf{L}\mathbf{X} = \mathbf{I}$ 。把两边按列分块， $\mathbf{X} = [\bm{x}_1,\dots,\bm{x}_n]$ 、 $\mathbf{I} = [\bm{e}_1,\dots,\bm{e}_n]$ ，矩阵方程裂成 $n$ 个互不相干的方程组：

\mathbf{L}\bm{x}_j = \bm{e}_j, \qquad j = 1,\dots,n

第二步：下三角让每一行只剩一个新未知量。 逐行看 $\mathbf{L}\bm{x}_j = \bm{e}_j$ ：第 1 行只含 $x_{1j}$ （系数 $l_{11}=1$ ），直接解出；第 $i$ 行含 $x_{1j},\dots,x_{ij}$ 共 $i$ 个未知量，而前 $i-1$ 个在处理前面几行时已经求出，于是第 $i$ 行同样只剩 $x_{ij}$ 一个新未知量。从上往下顺推、全程不需要消元，这就是前代法（forward substitution）。

第三步： $\bm{e}_j$ 的稀疏性给出显式递推。 $\bm{e}_j$ 只有第 $j$ 个分量为 1，代进去三段都能直接读出：

位置	方程	结果
$i < j$	右端为 0，逐行顺推	$x_{ij} = 0$ （逆仍是下三角）
$i = j$	$0 + 1\cdot x_{jj} = 1$	$x_{jj} = 1$ （逆仍是单位对角）
$i > j$	$\sum_{k<i} l_{ik}x_{kj} + x_{ij} = 0$	$x_{ij} = -\big(l_{ij} + \sum_{k=j+1}^{i-1} l_{ik}x_{kj}\big)$

第三行就是全部计算量所在。「上三角恒为 0、对角恒为 1、下三角按这个递推」三条合起来就是完整算法，它也顺手证明了单位下三角矩阵在乘法与求逆下封闭。

第四步：顺序必须列优先、列内自上而下。 看清递推的依赖：算 $x_{ij}$ 只用到 $\mathbf{L}$ 的第 $i$ 行（已知）和本列已求出的 $x_{kj}\ (k<i)$ ，完全不碰其他列。所以 $n$ 列天然可并行，而每列内部是一条长度 $n-j$ 的串行链。反过来若按行回代，会用到其他列尚未求出的元素，推不下去——这就是「按列回代」四个字的全部含义。

用 $\mathbf{L} = \begin{bmatrix}1&0&0\\2&1&0\\3&4&1\end{bmatrix}$ 走一遍，每步只填一格：

步	方程	求得
①	第 1 列， $i=j=1$ ： $1\cdot x_{11} = 1$	$x_{11} = 1$
②	第 1 列， $i=2$ ： $2\cdot 1 + x_{21} = 0$	$x_{21} = -2$
③	第 1 列， $i=3$ ： $3\cdot 1 + 4\cdot(-2) + x_{31} = 0$	$x_{31} = 5$
④	第 2 列： $x_{12} = 0$ （ $i<j$ ）， $x_{22} = 1$ （对角）	$x_{12}=0,\ x_{22}=1$
⑤	第 2 列， $i=3$ ： $3\cdot 0 + 4\cdot 1 + x_{32} = 0$	$x_{32} = -4$
⑥	第 3 列： $x_{13}=x_{23}=0$ ， $x_{33}=1$	该列无需回代

\mathbf{L}^{-1} = \begin{bmatrix}1&0&0\\-2&1&0\\5&-4&1\end{bmatrix}

验一下第 3 行： $3\cdot 1 + 4\cdot(-2) + 1\cdot 5 = 0$ 、 $3\cdot 0 + 4\cdot 1 + 1\cdot(-4) = 0$ 、 $1\cdot 1 = 1$ ，非对角全 0、对角全 1。

下面这个分步演示把上表的六步画成了矩阵填格，可以点「下一步」逐格看依赖关系：

计算量：第 $j$ 列约 $(n-j)(n-j+1)/2$ 次乘加， $n$ 列合计 $\approx n^3/6$ ，存储 $O(n^2)$ 。对比通用求逆的 $n^3$ 量级，单位三角结构省掉一半以上，而且不需要选主元。

最后一点实践上的便宜：要用的往往不是 $\mathbf{L}^{-1}$ 本身，而是 $\mathbf{L}^{-1}\mathbf{R}$ 。这时把上面的递推直接作用在 $\mathbf{R}$ 的行上，

\Delta[i,:] = \mathbf{R}[i,:] - \sum_{k<i} l_{ik}\,\Delta[k,:]

就能跳过 $\mathbf{L}^{-1}$ 这个 $n\times n$ 中间量。GPU kernel 里普遍这样写：省一块 shared memory，代价是那条 $n$ 步串行链暴露在关键路径上。

4.5 有限 Neumann 级数：前代法的矩阵语言

上一节结尾那条串行链是可以拆掉的，工具就是 §4.2 一笔带过的 Neumann 级数。有限 Neumann 级数与前代消去不是两个算法，是同一个算法的两种记法：级数是它的全局矩阵形式，代入是它的局部分量形式。写成矩阵形式的好处很实在——串行代入链能换成 $\log$ 级深度的矩阵乘。

符号沿用 §4.4 的单位下三角 $\mathbf{L}$ ，记 $\mathbf{N} := \mathbf{I}-\mathbf{L}$ （严格下三角），方程 $\mathbf{L}\bm{x}=\bm{b}$ 改写成 $(\mathbf{I}-\mathbf{N})\bm{x}=\bm{b}$ 。（§4.2 的 $\mathbf{N}$ 是严格上三角，转置即可互换，不影响下面任何一步。）

幂零性给出精确截断。 $\mathbf{N}^k$ 的非零元只出现在第 $k$ 条次对角线以下——每乘一次 $\mathbf{N}$ 就把非零区往左下推一格，推到 $k=n$ 时整块出界：

\mathbf{N}^{n} = \bm{0} \qquad\Longrightarrow\qquad \mathbf{L}^{-1} = (\mathbf{I}-\mathbf{N})^{-1} = \mathbf{I}+\mathbf{N}+\mathbf{N}^{2}+\cdots+\mathbf{N}^{n-1}

验证只要一行裂项相消： $(\mathbf{I}-\mathbf{N})\sum_{k=0}^{n-1}\mathbf{N}^{k} = \mathbf{I}-\mathbf{N}^{n} = \mathbf{I}$ 。

这里和解析情形有个本质差别：标量级数 $(1-x)^{-1}=\sum x^{k}$ 要求 $|x|<1$ ，矩阵情形要求谱半径 $\rho(\mathbf{N})<1$ ；而严格三角矩阵特征值全为 0， $\rho(\mathbf{N})=0$ 。所以收敛性无条件成立，且级数是有限的——不存在「迭代到什么精度」的问题，第 $n$ 项之后是恒等式，不是近似。

逐项累加就是逐行代入。 把部分和写成 Horner 形式：

\bm{x}^{(0)} = \bm{b}, \qquad \bm{x}^{(k)} = \bm{b} + \mathbf{N}\bm{x}^{(k-1)} = \sum_{j=0}^{k}\mathbf{N}^{j}\bm{b}

关键观察是 $\mathbf{N}$ 的前 $k$ 行在算到 $\mathbf{N}^{k}$ 时已经贡献完毕，所以每次迭代恰好「定死」一个分量：

迭代	表达式	本步定死的分量	前代法对应步骤
$k=0$	$\bm{x}^{(0)}=\bm{b}$	$x_1 = b_1$ （ $\mathbf{N}$ 第 1 行全 0，之后不再改动）	解第 1 行
$k=1$	$\bm{x}^{(1)}=\bm{b}+\mathbf{N}\bm{b}$	$x_2$ （只依赖已定死的 $x_1$ ）	解第 2 行
$k$	$\bm{x}^{(k)}=\bm{b}+\mathbf{N}\bm{x}^{(k-1)}$	$x_{k+1}$	解第 $k+1$ 行
$k=n-1$	级数展开完毕	$x_n$	解第 $n$ 行

两者连中间量都一致： $x_i$ 一旦在第 $i-1$ 次迭代定死就不再变动，其值等于前代法处理第 $i$ 行时写下的那个数。信息传播方向相同（自上而下）、步数相同（ $n$ ）、结果逐位相同。

用 §4.4 那个 $\mathbf{L} = \begin{bmatrix}1&0&0\\2&1&0\\3&4&1\end{bmatrix}$ 核对， $\mathbf{N} = \mathbf{I}-\mathbf{L}$ ：

\mathbf{N} = \begin{bmatrix}0&0&0\\-2&0&0\\-3&-4&0\end{bmatrix},\quad \mathbf{N}^{2} = \begin{bmatrix}0&0&0\\0&0&0\\8&0&0\end{bmatrix},\quad \mathbf{N}^{3} = \bm{0}

\mathbf{I}+\mathbf{N}+\mathbf{N}^{2} = \begin{bmatrix}1&0&0\\-2&1&0\\5&-4&1\end{bmatrix} = \mathbf{L}^{-1}

与 §4.4 前代法逐格填出的结果完全相同。留意 $x_{31}=5$ 的来路：前代法写成 $-(3 + 4\cdot(-2)) = 5$ ，级数里则拆成 $\mathbf{N}$ 贡献的 $-3$ 与 $\mathbf{N}^{2}$ 贡献的 $+8$ 。前代法中「代入已求出的 $x_{21}$ 」这个动作，在级数里就是多乘一次 $\mathbf{N}$ 。

代价对照：换的不是 FLOP，是深度。

	前代法	Neumann（逐项累加）	Neumann（平方倍增）
乘加次数	$\approx n^{3}/6$	$\approx n^{4}$	$\approx n^{3}\log_2 n$
串行深度	$n$ 步	$n$ 步	$\lceil\log_2 n\rceil$ 步
主体算子	AXPY / GEMV（level-2）	GEMM	GEMM
选主元	不需要	不需要	不需要

逐项累加没有意义（深度没降、算术还涨），真正可用的是平方倍增：

\sum_{k=0}^{2^{m}-1}\mathbf{N}^{k} = (\mathbf{I}+\mathbf{N})(\mathbf{I}+\mathbf{N}^{2})(\mathbf{I}+\mathbf{N}^{4})\cdots(\mathbf{I}+\mathbf{N}^{2^{m-1}})

每乘一个因子就把已累加的项数翻倍， $\lceil\log_2 n\rceil$ 步覆盖全部 $n$ 项。 $n = 64$ 时：前代法 64 步串行、约 4.4 万次乘加；倍增法 6 步、约 157 万次乘加——多做 36 倍算术，换深度从 64 降到 6。小规模三角求解本来就填不满 SM，算术单元是闲着的，这笔账在 Tensor Core 上通常划算。

块化才是实际写法。 纯倍增把 $\mathbf{N}$ 的三角稀疏性浪费掉了（ $\mathbf{N}^{2}$ 起非零区迅速收缩，却仍按稠密 GEMM 算）。工程上取两项截断的分块形式：把 $\mathbf{L}$ 二分为 $\begin{bmatrix}\mathbf{L}_{11}&\bm{0}\\ \mathbf{L}_{21}&\mathbf{L}_{22}\end{bmatrix}$ ，则

\mathbf{L}^{-1} = \begin{bmatrix}\mathbf{L}_{11}^{-1}&\bm{0}\\ -\mathbf{L}_{22}^{-1}\mathbf{L}_{21}\mathbf{L}_{11}^{-1}&\mathbf{L}_{22}^{-1}\end{bmatrix}

这是分块层面的「一次前代」：两个对角块各自递归求逆（互不依赖，可并行），耦合项两次 GEMM 拼上。结构与 §4.3 的 $\mathbf{T}$ 分块式同构——那里拼的是三角矩阵乘法，这里拼的是求逆，都是「块内独立、块间一次耦合」。递归到块足够小（典型 16~32）时切回 §4.4 的逐格前代，这是 cuBLAS / CUTLASS 里 TRSM 与三角求逆的常见路线。

回到 §4.2 的 $\mathbf{T} = \mathbf{D}(\mathbf{I}+\mathbf{N})^{-1}$ ： $p = 64$ 或 128 的块内三角求逆正好落在这个尺度上，两三层分块递归就能把它压平成几次小 GEMM。这也是 §4 全节的落点——UT 变换把递推消成三角方程，Neumann 级数再把三角方程消成矩阵乘，两步都在做同一件事：拿算术换深度。

一句话总结：严格三角矩阵幂零，Neumann 级数在 $n$ 项精确截断；逐项累加与逐行代入是同一过程的两种记法——前者把它写成矩阵乘之和，于是能用平方倍增把 $n$ 步串行压成 $\lceil\log_2 n\rceil$ 步 GEMM。

5. 闭式推导

前文 §4 的 UT 变换是「闭式」的一个典型例子：把一段串行递推（dlarft 逐列构造 T）改写成只依赖原始输入的显式表达（T = D(I + striu(GD))⁻¹）。这个思路不限于 Householder 块化，下面用两个更小的例子再演示同样的手法，方便对比理解：先是一阶线性递推（§5.1，工具是「除以累积乘积」），再是二阶的 Fibonacci（§5.2，工具是「打包成向量 + 对角化」）。

5.1 一阶线性递推式的闭式推导

问题：给定

X_n = a_n X_{n-1} + b_n, \qquad X_0 \text{ 已知}

其中 $a_n, b_n$ 是已知系数序列。目标是跳过逐次迭代，直接写出 $X_n$ 的表达式。

思路：除过去，把递推变成求和

定义累积乘积

A_n := \prod_{i=1}^{n} a_i = a_n a_{n-1} \cdots a_1

逐项观察模式：

$n = 1$ ：

X_1 = a_1 X_0 + b_1

两边除以 $A_1 = a_1$ ，两项各除各的：

\frac{X_1}{A_1} = X_0 + \frac{b_1}{A_1}

$n = 2$ ： 把 $X_1$ 整个代入（先摊开，不化简）：

\begin{aligned} X_2 &= a_2(a_1 X_0 + b_1) + b_2 \\ &= a_2 a_1 X_0 + a_2 b_1 + b_2 \end{aligned}

除以 $A_2 = a_2 a_1$ ： $X_0$ 项的系数被整除消掉， $b_2$ 是新来的、身上没背 $a$ ， $b_1$ 背了一个 $a_2$ 正好被消掉：

\frac{X_2}{A_2} = X_0 + \frac{b_1}{a_1} + \frac{b_2}{a_2 a_1} = X_0 + \frac{b_1}{A_1} + \frac{b_2}{A_2}

$n = 3$ ： 代入 $X_2$ ：

X_3 = a_3 a_2 a_1 X_0 + a_3 a_2 b_1 + a_3 b_2 + b_3

除以 $A_3 = a_3 a_2 a_1$ ： $b_1$ 背的两个 $a$ 、 $b_2$ 背的一个 $a$ 全部对消：

\frac{X_3}{A_3} = X_0 + \frac{b_1}{A_1} + \frac{b_2}{A_2} + \frac{b_3}{A_3}

一般 $n$ ： 归纳即得

\frac{X_n}{A_n} = X_0 + \sum_{i=1}^{n} \frac{b_i}{A_i}

记 $S_n = \sum_{i=1}^{n} b_i/A_i$ ，得闭式：

\boxed{\,X_n = A_n\,(X_0 + S_n)\,}

5.2 二阶线性递推：Fibonacci 闭式

§5.1 的递推每步只回看一步。Fibonacci $F_n = F_{n-1} + F_{n-2}$ 要回看两步，「除以累积乘积」那招直接用不上了。但整套气质和 §4 的 UT 推导是一样的：每一步都被上一步逼出来，全程只问三个递进的问题——

依赖「前两步」怎么变成依赖「前一步」？——打包成向量；
一整条 $n$ 步的链能不能缩成一次运算？——常系数；
缩出来的那个对象里，耦合怎么拆开？——找它「退化成数乘」的方向。

注意这里不先猜 $F_n = \lambda^n$ 再回代验证（那是答案在前、理由在后的事后验算），特征方程会被第 3 问自己逼出来。

5.2.1 打包：二阶标量递推 → 一阶向量递推

（a）打包动机. 问题在于 $F_n$ 同时看 $F_{n-1}$ 和 $F_{n-2}$ 。但如果把相邻两个值打包成一个状态，新状态就只依赖旧状态：

\bm{v}_n := \begin{pmatrix} F_n \\ F_{n-1} \end{pmatrix}

（b）矩阵表示. 代入递推 $F_n = F_{n-1}+F_{n-2}$ ，逐行看：

\begin{pmatrix} F_n \\ F_{n-1} \end{pmatrix} = \begin{pmatrix} F_{n-1}+F_{n-2} \\ F_{n-1} \end{pmatrix} = \underbrace{\begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}}_{\mathbf{A}} \begin{pmatrix} F_{n-1} \\ F_{n-2} \end{pmatrix} \qquad\Longrightarrow\qquad \bm{v}_n = \mathbf{A}\,\bm{v}_{n-1}

（c）初值与效果. 这一步零猜测：「依赖两步」被换成了「状态里装两个数、每步只依赖上一步」——依赖长度 $2 \to 1$ ，代价是状态从标量变二维向量。初值 $F_0=0,\ F_1=1$ 给

\bm{v}_1 = \begin{pmatrix} F_1 \\ F_0 \end{pmatrix} = \begin{pmatrix} 1 \\ 0 \end{pmatrix}

（d）数值验证. 每一步直接乘出来：

\begin{aligned} \bm{v}_2 &= \mathbf{A}\binom{1}{0} = \binom{1}{1} = \binom{F_2}{F_1} \quad \checkmark \\ \bm{v}_3 &= \mathbf{A}\binom{1}{1} = \binom{2}{1} = \binom{F_3}{F_2} \quad \checkmark \\ \bm{v}_4 &= \mathbf{A}\binom{2}{1} = \binom{3}{2} = \binom{F_4}{F_3} \quad \checkmark \end{aligned}

5.2.2 常系数 ⇒ 整条链缩成一次幂

$\mathbf{A}$ 与 $n$ 无关（这就是「常系数」），所以往回一路代：

\bm{v}_n = \mathbf{A}\,\bm{v}_{n-1} = \mathbf{A}^2\bm{v}_{n-2} = \cdots = \mathbf{A}^{n-1}\,\bm{v}_1

对比 §5.1 的变系数情形：那里每步的乘子 $a_t$ 各不相同，只能保留前缀积 $A_t = a_1\cdots a_t$ ；这里每步乘的是同一个 $\mathbf{A}$ ，前缀积坍缩成一个幂。于是整个问题变成一件事：

F_n = \big[\mathbf{A}^{n-1}\bm{v}_1\big]_1 = \Big[\mathbf{A}^{n-1}\binom{1}{0}\Big]_1

5.2.3 困难定位： $\mathbf{A}^{n-1}$ 难算，因为 $\mathbf{A}$ 不是对角的

如果 $\mathbf{A}$ 是对角阵 $\mathrm{diag}(\lambda_1,\lambda_2)$ ，那 $\mathbf{A}^{n-1}$ 就是 $\mathrm{diag}(\lambda_1^{\,n-1},\lambda_2^{\,n-1})$ ，每个分量各自求幂， $O(1)$ 写完。而 $\mathbf{A} = \binom{1\ \ 1}{1\ \ 0}$ 不是对角的——两个分量每步互相搅。

顺序的问题来了：能不能换个坐标系，让 $\mathbf{A}$ 变成对角的？ 一个矩阵在某组基下对角，等价于这组基的每个成员被 $\mathbf{A}$ 作用后只是伸长：

\mathbf{A}\bm{v} = \lambda\,\bm{v}, \qquad \bm{v} \neq \bm{0}

这不是假设解的形状，而是「要达成对角化，方向必须满足什么」的必要条件。移项：

(\mathbf{A}-\lambda\mathbf{I})\bm{v} = \bm{0}

非零解存在 $\iff$ 系数矩阵奇异 $\iff$ 行列式为零。代进 $\mathbf{A}$ ：

\det\begin{pmatrix} 1-\lambda & 1 \\ 1 & -\lambda \end{pmatrix} = (1-\lambda)(-\lambda) - 1 = \lambda^2 - \lambda - 1 = 0

那个熟悉的二次方程出现了——但这次是被「让 $\mathbf{A}$ 数乘化的方向必须存在」逼出来的，不是猜 $F_n=\lambda^n$ 猜出来的。解它：

\lambda = \frac{1\pm\sqrt5}{2} \quad\Longrightarrow\quad \varphi := \frac{1+\sqrt5}{2} \approx 1.618, \qquad \psi := \frac{1-\sqrt5}{2} \approx -0.618

5.2.4 把两个方向具体求出来

$\lambda=\varphi$ 时解 $(\mathbf{A}-\varphi\mathbf{I})\bm{v}=\bm{0}$ ：取第二行 $x - \varphi y = 0$ ，即 $x=\varphi y$ ，令 $y=1$ 得方向 $\bm{v}_\varphi=(\varphi,\ 1)^{\intercal}$ 。直接验证它确实只被伸长——注意用 $\varphi^2=\varphi+1$ （这正是 $\varphi$ 满足的方程）：

\mathbf{A}\binom{\varphi}{1} = \binom{\varphi+1}{\varphi} = \binom{\varphi^2}{\varphi} = \varphi\binom{\varphi}{1}

同理 $\psi^2=\psi+1$ 给出 $\bm{v}_\psi=(\psi,\ 1)^{\intercal}$ ，且 $\mathbf{A}\bm{v}_\psi=\psi\,\bm{v}_\psi$ 。 $\varphi\neq\psi$ 所以两个方向不共线，张成整个平面——这就是我们要的对角化基。

5.2.5 换坐标：把初值拆进这两个方向

$\binom{1}{0}$ 在这组基下写成组合 $\binom{1}{0}=c_\varphi\binom{\varphi}{1}+c_\psi\binom{\psi}{1}$ 。按分量写出就是普通二元一次方程组：

\begin{cases} \varphi\,c_\varphi + \psi\,c_\psi = 1 & \text{(第一行)} \\ c_\varphi + c_\psi = 0 & \text{(第二行)} \end{cases}

第二行给 $c_\psi=-c_\varphi$ ，代进第一行： $c_\varphi(\varphi-\psi)=1$ 。而 $\varphi-\psi=\sqrt5$ ，所以

c_\varphi = \frac{1}{\sqrt5}, \qquad c_\psi = -\frac{1}{\sqrt5}

5.2.6 传播：换基之后，每步只乘自己的数

$\bm{v}_1$ 拆成两个方向后， $\mathbf{A}$ 对每个方向单独作用：

\begin{aligned} \bm{v}_n = \mathbf{A}^{n-1}\bm{v}_1 &= \frac{1}{\sqrt5}\,\mathbf{A}^{n-1}\binom{\varphi}{1} -\frac{1}{\sqrt5}\,\mathbf{A}^{n-1}\binom{\psi}{1} \\ &= \frac{1}{\sqrt5}\,\varphi^{\,n-1}\binom{\varphi}{1} -\frac{1}{\sqrt5}\,\psi^{\,n-1}\binom{\psi}{1} \end{aligned}

两个方向的耦合在换基那一刻就没了：原坐标里每步是 $\mathbf{A}$ （两分量纠缠），新坐标里每步只是乘 $\varphi$ 或乘 $\psi$ （两分量各自独立）。

5.2.7 取第一分量，闭式一步到位

F_n = \big[\bm{v}_n\big]_1 = \frac{1}{\sqrt5}\,\varphi^{\,n-1}\cdot\varphi - \frac{1}{\sqrt5}\,\psi^{\,n-1}\cdot\psi = \boxed{\,\frac{\varphi^{\,n} - \psi^{\,n}}{\sqrt5}\,}

数值核对（到这一步才验算，且只是检查、不是方法）：

\begin{aligned} n=1 &: \frac{\varphi-\psi}{\sqrt5} = \frac{\sqrt5}{\sqrt5} = 1 \quad \checkmark \\ n=2 &: \frac{\varphi^2-\psi^2}{\sqrt5} = \frac{(\varphi-\psi)(\varphi+\psi)}{\sqrt5} = \frac{\sqrt5\cdot 1}{\sqrt5} = 1 \quad \checkmark \quad (\varphi+\psi=1) \\ n=3 &: \frac{\varphi^3-\psi^3}{\sqrt5} \approx \frac{4.236 - (-0.236)}{2.236} = 2 \quad \checkmark \end{aligned}

和 §4 对上一句：UT 变换是把「 $p$ 步串行递推」换成「一次 Gram 矩阵 + 一次三角求解」，这里是把「 $n$ 步串行迭代」换成「一次对角化 + 两次求幂」。手法都是换一个让耦合消失的表示，只不过 UT 里那个表示是上三角 $\mathbf{T}$ ，这里是特征基。

5.3 积分因子视角

回到 §5.1 的一阶递推。令 $Y_n = X_n / A_n$ ，代入原递推：

\begin{aligned} A_n Y_n &= a_n A_{n-1} Y_{n-1} + b_n \\ &= A_n Y_{n-1} + b_n \quad \text{（因为 $a_n A_{n-1} = A_n$）} \\ \Rightarrow\quad Y_n - Y_{n-1} &= \frac{b_n}{A_n} \end{aligned}

这就变成了一个纯差分方程，直接求和即得 $Y_n = X_0 + \sum_{i=1}^{n} b_i / A_i$ 。 $A_n = \prod a_i$ 就是离散版本的「积分因子」，和常微分方程里 $e^{\int p(x)\,dx}$ 的作用一模一样。

5.4 与 Householder 块化的类比

	一阶线性递推（§5.1）	Fibonacci（§5.2）	Householder 块化（§4）
递推式	$X_n = a_n X_{n-1} + b_n$	$F_n = F_{n-1} + F_{n-2}$	$\bm{w}_t = \tau_t\bm{u}_t - \mathbf{W}\,\mathrm{striu}(\cdot)_t$
消依赖的手法	两边除以累积乘积 $A_n = \prod a_i$	打包成向量 + 换到特征基	把 $p$ 条递推堆成一个矩阵方程
得到的方程	$Y_n - Y_{n-1} = b_n/A_n$ （纯差分）	$\bm{v}_n = \mathbf{A}^{n-1}\bm{v}_1$ ， $\mathbf{A}\bm{v}=\lambda\bm{v}$	$\mathbf{W}(\mathbf{I} + \mathrm{striu}(\mathbf{G}\mathbf{D})) = \mathbf{Y}\mathbf{D}$
闭式	$X_n = A_n(X_0 + S_n)$	$F_n = (\varphi^n - \psi^n)/\sqrt5$	$\mathbf{T} = \mathbf{D}(\mathbf{I} + \mathrm{striu}(\mathbf{G}\mathbf{D}))^{-1}$
原始数据	$\{a_i\}, \{b_i\}, X_0$	$\mathbf{A}$ 与 $\bm{v}_1$	$\mathbf{G} = \mathbf{Y}^{\intercal}\mathbf{Y}$ 与 $\{\tau_t\}$
收尾操作	一次前缀求和	两次标量求幂	一次 GEMM + 一次三角求解

三者的共同结构：把一个串行依赖链（每一步等前一步结果）变成一个只依赖原始数据的闭式表达式。 差别只在「消依赖」用的工具——除以累积乘积、换到特征基、把递推堆成三角方程；本质都是找一个让递推项对消掉的表示。

5.5 启示：闭式推导的一般步骤

从以上三个例子可以抽象出闭式推导的通用框架：

先把依赖长度压到 1 — 若一步要回看多步（Fibonacci 回看两步），先把相邻状态打包成向量，用「状态变宽」换「依赖变短」。这一步纯粹是记账，不需要任何猜测。
识别依赖链 — 哪一步必须等前一步？（一阶递推里 $X_{n-1}$ 是 $X_n$ 的前置；dlarft 里第 $k+1$ 列必须等第 $k$ 列）
把递推摊开，找出那个「整除因子」／让耦合消失的表示 — 展开前几项，看被累积起来的到底是什么。一阶递推里是 $\prod a_i$ ，Fibonacci 里是 $\mathbf{A}$ 的特征基，UT 里是那个单位上三角 $\mathbf{M}$ 。
两边同除／换基／堆成方程 — 依赖项被这个因子对消掉，剩下的是不含中间量的关系式。
解出闭式并验证自洽 — 一阶递推求和即得；Fibonacci 取回第一分量后代 $n=1,2,3$ 核对；UT 里对 $\mathbf{T}$ 求逆回到 $\mathbf{T}^{-1} = \mathbf{D}^{-1} + \mathrm{striu}(\mathbf{G})$ ，与 §3.1 的递推逐位对上。

第 2 步是全部难点所在：能不能找到那个因子／表示，取决于递推的系数结构。线性（每步只做一次乘加）是关键——非线性递推一般没有这种因子，也就没有闭式。另外第 0 步和第 2 步之间有个耐人寻味的呼应：打包成向量把「时间上的长依赖」换成了「空间上的分量耦合」，而对角化再把分量耦合彻底解开。 序列建模里把 $\mathbf{A}$ 换成数据相关的转移矩阵，这两步就分别对应 state expansion 和 chunkwise 并行化。

6. 总结

一个恒等式贯穿全文： $\prod_{t=1}^{p}(\mathbf{I} - \tau_t\bm{u}_t\bm{u}_t^{\intercal}) = \mathbf{I} - \mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal}$ ， $\mathbf{T}$ 是 $p \times p$ 上三角。三个表达式只是同一个恒等式的不同写法与不同算法。
$\mathbf{T}$ 的每个数是一个外积的系数（§3.3）。展开成 $\sum_{i \le j}\mathbf{T}_{ij}\bm{u}_i\bm{u}_j^{\intercal}$ 后，对角 $\mathbf{T}_{tt} = \tau_t$ 是单个反射的强度，上三角 $\mathbf{T}_{ij}$ 是反射之间的累积耦合。上三角性来自外积乘法 $(\bm{u}\bm{v}^{\intercal})(\bm{x}\bm{y}^{\intercal}) = (\bm{v}^{\intercal}\bm{x})\bm{u}\bm{y}^{\intercal}$ 只保留「前 $\times$ 后」方向，是乘积次序的必然结果。
两次不同方向的优化：

阶段	形式	优化对象	代价
WY（1987）	$\mathbf{I} - \mathbf{W}\mathbf{Y}^{\intercal}$	能块化了	$\mathbf{W}$ 额外占 $O(np)$
compact WY（1989）	$\mathbf{I} - \mathbf{Y}\mathbf{T}\mathbf{Y}^{\intercal}$	存储： $O(np) \to O(p^2)$	$\mathbf{T}$ 逐列递推， $p$ 步串行
UT（2006）	$\mathbf{T} = \mathbf{D}\,(\mathbf{I} + \mathrm{striu}(\mathbf{G}\mathbf{D}))^{-1}$	依赖图： $p$ 步串行 $\to$ 1 次 GEMM + 1 次三角求解	多一次 $O(p^2)$ 回代

第二步不是「想办法并行化递推」，而是证明根本不需要递推： $\mathbf{T}$ 的信息量不超过 Gram 矩阵 $\mathbf{G} = \mathbf{Y}^{\intercal}\mathbf{Y}$ 加一组 $\tau$ ，逐列递推只是一种实现方式带进来的人为依赖。

单位三角求逆就是一次前代（§4.4）。 $\mathbf{L}\mathbf{X} = \mathbf{I}$ 按列裂成 $n$ 个独立方程组， $\bm{e}_j$ 的稀疏性直接给出三条规则： $i<j$ 恒为 0、 $i=j$ 恒为 1、 $i>j$ 走递推 $x_{ij} = -(l_{ij} + \sum_{k=j+1}^{i-1}l_{ik}x_{kj})$ 。合计 $\approx n^3/6$ 次乘加、不需选主元；列间可并行、列内 $n$ 步串行。实践中往往连 $\mathbf{L}^{-1}$ 都不物化，把递推直接作用在右端项的行上。
前代法与有限 Neumann 级数是同一个算法的两种记法（§4.5）。 $\mathbf{N} = \mathbf{I}-\mathbf{L}$ 严格三角故幂零， $\mathbf{L}^{-1} = \sum_{k=0}^{n-1}\mathbf{N}^{k}$ 精确成立（不是近似， $\rho(\mathbf{N})=0$ 无收敛条件）；Horner 迭代 $\bm{x}^{(k)} = \bm{b}+\mathbf{N}\bm{x}^{(k-1)}$ 每步恰好定死一个分量，与前代法逐行代入逐位一致。写成级数的价值在于可以用平方倍增 $\prod_{j}(\mathbf{I}+\mathbf{N}^{2^{j}})$ 把 $n$ 步串行压成 $\lceil\log_2 n\rceil$ 步 GEMM： $n=64$ 时深度 $64 \to 6$ ，代价是算术量涨约 36 倍。

一句可迁移的启示：碰到串行递推时，先问这个依赖是数学必需还是实现引入的。 把递推展开看它到底依赖哪些量（这里是 $p^2$ 个内积），如果那些量本身互不依赖，递推就一定能改写成方程一次解掉。序列建模里叫 chunking，数值线代里叫 blocking，同一个 trick 隔了三十年在两个领域各自被发明了一次。

参考：