摘要

1. 为了评估矩阵求逆方法的相对有效性，对计算的逆矩阵中的误差有一个先验界限是有用的。本文确定了一些最有效的直接方法的误差范围。为了充分说明我们所使用的技术，我们对浮点计算和定点计算分别进行了一些分析。在所有情况下，都假定计算是使用$t$ 个二进制位的精度执行，但应该了解，在同时具有固定和浮点设备的计算机上，定点数中允许的位数要大于浮点数尾数的位数。用于分析浮点计算的技术本质上是[8]中的技术，并且假设读者已经熟悉该论文。

2. 用向量范数和矩阵范数来表示误差界是最方便的，除了明确指出相反的情况外，我们自始至终都使用欧几里得向量范数和谱矩阵范数。 为方便起见，在第9节中给出了这些规范的主要性质。

在最近的一篇论文中，我们分析了在求解方程时产生的舍入误差的影响

通过高斯消元法与“尺寸旋转”。我们证明了计算出的解是

我们可以得到$E$ 和$\delta b$ 的范围，前提是对连续简化矩阵的最大元素与$A$ 的最大元素之比做了假设。本文讨论的所有方法都依赖于将原始矩阵$A^{\left(1\right)}$ 转化为矩阵$A^{\left(2\right)},A^{\left(3\right)}，\cdot \cdot \cdot ，$$A^{\left(k\right)}$ 的连续变换，使得每个$A^{\left(3\right)}$ 等价于$A^{\left(1\right)}$ ，最终$A^{\left(k\right)}$ 是三角形。在任何方法的先验分析中，一个基本要求是对数量$r$ 有一个界限，定义为

对于基于高斯消元的方法，只有在每个阶段使用尺寸旋转时，我们才能够获得一般非对称矩阵的合理上限。这个界限远不是明确的，但至少导致了一个有用的结果，并且证明了限制$R$ 大小的因素。对于一些在实践中很重要的矩阵类，给出了$R$ 的一个更小的界限，并且这个界限即使在采用有限形式的大小旋转时仍然成立。最后，我们证明了对于正定矩阵$R\leqq 1$ ，即使交换“尺寸旋转”也可以。

为了与这种高斯类型的方法进行比较，我们对基于通过正交矩阵的等效变换将$A$ 简化为三角形的方法进行了分析。对于这些，很容易证明对于任何矩阵$R<\sqrt{n}$ ，这样大小的控制就得到了保证。由于对$R$ 大小的限制，一般矩阵的正交约简得到的误差范围小于消除方法得到的误差范围。 然而，正交约简比消除需要更多的计算，并且由于在实践中$R$ 的值远远低于我们获得的界限，因此我们认为使用正交变换很少，如果有的话，是合理的。

3. 我们并不声称我们给出的先验界限在任何方面都是“最好的”，并且在许多情况下，我们知道它们可以通过更详细的论证而有所减少。 然而，我们认为错误分析的本质简单性可能会被过多的细节所掩盖。 在任何情况下，界限都不考虑统计效应，在实践中，这总是通过一个比详细阐述论点所能达到的更为重要的因素来减少真实误差。我们还提出了一些论点，这些论点部分解释了为什么在实践中计算的逆运算比预期的要准确得多，即使考虑到统计效应也是如此。

误差的所有先验界限都明确包含因子$\Vert A^{-1}\Vert$ 。这意味着在得到计算的逆之前，我们几乎不可能利用这些结果。 分析的性质使得如何使用计算的逆来获得其误差的后验上限变得很明显，但结果的主要价值是表明给定方法在计算中给定精度的适用范围， 因此应根据这一点进行判断。

高斯消元法中$R$ 的上限

4.当一般矩阵通过高斯消元法简化为三角形时，我们首先推导出$R$ 的上限，在每个阶段选择整个剩余方阵中最大模数的元素作为关键元素。 我们将此称为“完全”旋转尺寸，而不是在每个阶段选择前导列中的最大元素，我们称之为“部分”旋转尺寸。

用$A^{\left(n\right)}$ 表示原始矩阵，用$A^{\left(n-1\right)},A^{\left(n-2\right)},\dots ,A^{\left(1\right)}$ 表示简化矩阵(使$A^{\left(x\right)}$ 是$r$ 阶矩阵)，用$p_r$ 表示$A^{\left(r\right)}$ 的关键元素的模，使得

然后我们得到

另一方面，Hadamard 定理给出了

因为$A^{\left(r\right)}$ 的每一列的长度都以$r^{½}{p}_r$ 为界。 我们写作

但保留

取关系式 (4.4) 和 (4.5) 的对数并写作

我们得到

当$r=2,3,\cdot \cdot \cdot n-1$ 时，将( 4.7 )除以$r\left(r-1\right)$ ，( 4.8 )再除以$\left(n-1\right)$ ，观察到

给出

其中

将$|\det A^{\left(n\right)}|$ 代入(4.3)(其中$r=n$ )得到

因此

这给出了最后一个基准值与第一个基准值之比的上限。特别是，如果$p_n\leqq 1$ ，则$p_1\leqq g(n)$ 。

另一方面，如果我们写

因此，$A^{\left(n\right)}$ 的所有列的长度都以$l$ 为界限，则

因此根据(4.11)

特别是当$l\leqq 1$ 时

函数$f(n)$ 和$g(n)$ 的增长相对缓慢，这很容易证明

在表 1 中，我们给出了代表$n$ 值的$f(n)$ 和$g(n)$ 的值。

标准化的两种形式将在本文中使用，分别称为规范化(I)和规范化(II)。(4.15)和(4.18)对于任何以适当方式归一化的矩阵都成立，实际上它们的真实性决不依赖于(4.20)和(4.21)中的左侧不等式。然而，在我们看来，只有当原始矩阵的所有行和列都具有可比较的范数时，旋转尺寸才是合理的策略。

对于规范化 (I)，这可以通过缩放所有行和列来实现，以便它们至少包含一个满足$\frac{1}{2}\leqq a_{ij}|<1$ 的元素。对于规范化 (II)，这之后应该对列进行额外的缩放，因此满足（4.21）。F. L. Bauer 建议将此类矩阵称为“平衡”矩阵。 在实践中，没有必要确保精确的平衡，而只是没有完全由“小”元素组成的行或列。 严格来说，在选择枢轴之前，应在还原的每个阶段进行平衡，但这些后来的平衡并不重要。

5. 不等式 (4.15) 肯定不有效。 因为只有当 (4.4) 中的等式对所有$r$ 都成立时，等式才能成立。 如果这是真的，我们可以证明，正如我们对$r=n$ 所做的那样，那么

因此

除了小$r$ 之外，这个比值仅略大于1。这意味着关键序列$p_n,p_{n-1},\cdot \cdot \cdot$ 起初增长相当缓慢。现在，等式(4.4)只有在$A^{\left(r\right)}$ 的所有元素都等于$\pm p_r$ 且所有列正交时才成立。对于某些$r$ 值，例如$r=3$ ，这显然是不可能的。进一步，如果$a\genfrac{}{}{0ex}{}{\left(r\right)}{ij}=\pm p_r$ ，则所有$a\genfrac{}{}{0ex}{}{\left(r-1\right)}{ij}$ 为$\pm 2p_r$ 或为零。因此，如果在任何阶段等号在(4.4)中成立，那么$p_{r-1}/p_r=2$ 或$A^{\left(r-1\right)}$ 都为空，在任何一种情况下(5.2)都无法满足。

这些考虑表明，支点的最小上限可能比我们给出的要小得多。根据我们的经验，最后一个主元远远低于给定的界限，并且在实践中还没有遇到$p_l/p_n$ 等于8的矩阵。请注意，我们肯定有

因为$p_r$ 是$\left(n+1-r\right)$ 阶矩阵的归约的最后一个支点，但如果$p_l$ 接近最大值，则初始增长率必须远低于对应于 (5.3) 的增长率 。

正定矩阵

6. 有一些特殊类型的矩阵，可以给出最大主元的更好界限，令人惊讶的是，这些矩阵都不需要完全的主元来确定大小。

最重要的一类是正定矩阵。 我们将证明，当一个正定矩阵通过高斯消元减少而没有任何大小的旋转时，任何减少的矩阵中的元素都不会超过原始矩阵中的最大元素。 现在可以方便地调用原始矩阵$A^{\left(1\right)}$ 和简化矩阵$A^{\left(2\right)},A^{\left(3\right)},\cdot \cdot \cdot ,A^{\left(20\right)}$ 以使$A^{\left(r\right)}$ 是有序的

我们首先证明$A^{\left(2\right)}$ 是正定的。

第一个约简中使用的乘数$m_{il}$ 定义为

并且$a\genfrac{}{}{0ex}{}{\left(1\right)}{i1}$ 不为零，因为$A^{\left(1\right)}$ 是正定的。 $A^{\left(2\right)}$ 的元素由下式给出

根据$A^{\left(1\right)}$ 的对称性，我们得到

如果$A^{\left(2\right)}$ 不是正定的，则存在$x_2,x_3,\dots ,x_n$ 的非空集，使得${\sum }_{a_{ij}}^{\left(2\right)}{x}_i{x}_j\leqq 0$ 。如果我们写

然后有了这个$x_1,x_2,\cdot \cdot \cdot x_n$ 我们得到

这是不可能的，因为$A^{\left(1\right)}$ 是正定的。 我们还有进一步

$a\genfrac{}{}{0ex}{}{\left(2\right)}{ii}$ 是正的，因为$A^{\left(2\right)}$ 是正定的。因此，$A^{\left(2\right)}$ 的对角线项不大于$A^{\left(1\right)}$ 的对角线项，因为正定矩阵的最大元素位于对角线上

同样地，$A\left(3\right)\cdot \cdot \cdot A\left(n\right)$ 是正定的。 特别是如果$|a\genfrac{}{}{0ex}{}{\left(1\right)}{ij}\leqq 1$ （所有的$i,j$ ）, 那么这对所有$a\genfrac{}{}{0ex}{}{\left(r\right)}{ij}$ 都是正确的。

乘数可能很大，例如对于矩阵

因此，高斯消元要么在浮点运算中执行，要么需要以类似于通常在反向替换中使用的方式引入比例因子。

三角分解

7. 对于正定对称矩阵的对称 Cholesky 分解，我们可以证明一个更令人满意的结果。 我们证明如果$A$ 是正定的并且$|a_{ij}|\leqq 1$ ，则存在一个下三角矩阵$L$ 使得

证明使用归纳法。假设对于高达$\left(n-1\right)$ 阶的矩阵是正确的。如果$A$ 是$n$ 阶正定矩阵，我们可以将其划分为以下形式

其中$B$ 是阶$\left(n-1\right)$ 的正定矩阵，$a$ 是阶$\left(n-1\right)$ 的向量。 通过假设有一个$L$ 满足

如果我们选择$l$ 使得

既然$L$ 是三角形且非奇异的，那么显然可以这样做

其中

两边都取行列式

现在$\det A$ 是正的，因为$A$ 是正定的； 因此$l\genfrac{}{}{0ex}{}{2}{nn}$ 是正的并且$l_{nn}$ 是实数并且可以被认为是正的。 如果$l$ 的分量用 $l_{n1},l_{n2},...,l_{n,n-1}$ 表示，则由方程 (7.7)

并且每个$l_{ni}$ ，满足$|l_{ni}|\leqq 1$ 。 等式 (7.6) 和 (7.7) 现在给出了$A$ 的适当分解，并且三角矩阵的所有元素都在$\pm 1$ 之间。 这样做的结果是可以在定点算术中执行$L{L}^T$ 分解。 顺便说一下，我们已经证明了$L$ 中任意一行元素的平方和都小于$1$ ，这个结果在后面会用到。

$R<g\left(n\right)$ 的其他矩阵

8. Hessenberg 形式的矩阵与特征值问题有关。我们现在考虑使用归一化$l$ 减少上 Hessenberg 矩阵，使用“部分”旋转来确定大小。在每个阶段，简化矩阵的前列中只有两个元素，其中一个元素位于尚未修改的行中。由于即使对大小进行部分旋转也可以确保乘数在$\pm 1$ 之间，因此很容易看出第$r$ 个缩减矩阵中的最大元素位于$\pm r$ 之间，因此最大的关键值为$n$ ，这只能通过最后一个主元来实现。进一步，如果我们称最终的三角矩阵为$U$ , $M$ 是$A$ 乘以$U$ 得到的矩阵，因此

很明显，$M$ 是下三角，其所有元素的模量都小于或等于1。如果我们用通常的方式写，

那么$L=M^{-1}$ ，并且我们已经证明了$L^{-1}$ 的所有元素都在$\pm 1$ 之间。

Hessenberg 矩阵的一种特殊类型是三对角矩阵。 根据大致相同的论点，很明显$U$ 的最大元素现在位于$\pm 2$ 之间，并且$L^{-1}$ 的所有元素再次位于$\pm 1$ 之间。 此外${\Sigma }_j|u_{ij}|\leqq 3$ 。

最后，我们考虑具有主导对角元素的矩阵，使得

如果我们在不旋转大小的情况下执行消除，那么我们通常在第一阶段

因此

此外

因此，对角线项在所有$A^{\left(r\right)}$ 中以完全相同的方式占主导地位，并且从(8.5)开始，$A^{\left(r\right)}$ 中任意列元素的模的和随着$r$ 的增大而减小。因此

因此这个比值$R$ 肯定小于$2$ 。通过类似的论证，我们看到，如果我们把约简写成以下形式

那么$M$ (或$L^{-1}$ )中的所有元素都小于$1$ 。

即使对角元素不是严格占主导地位的，最后一种类型的行为也很常见。 例如，$C_n$ 形式的矩阵，其中

在高斯消元过程中保持这个特征，即使在任何阶段都不满足 (8.3) 时，对角线也会逐渐占主导地位。

在实际中遇到的许多矩阵并不严格属于我们所提到的任何特殊类别，但在性质上与它们相似。根据我们的经验，当使用大小的“完全”主元时，关键序列即使不是单调递减的，通常呈总体下降趋势。对于病态矩阵，这通常是非常明显的。这并不奇怪，因为如(4.11)所示，最后一个主元$p_1$ 满足

当$A^{\left(n\right)}$ 是病态时，$|\det A^{\left(n\right)}|$ 通常比$p_n^n$ 小得多。

特殊矩阵的行为可能给人的印象是“部分”旋转大小与“完全”旋转一样有效。 以下示例表明它可能还不够。 考虑以下类型的矩阵$B_n$ ：

如果我们使用部分主元，那么主元总是在左上角，序列是$+1,+1,\cdot \cdot \cdot +1,{\left(-2\right)}^{n-1}$ 。如果使用完全主元，则序列为$+1,-2,+2,+2,\cdot \cdot \cdot ,+2$ 。我们稍后将展示如何修改$B_n$ ，从而得到一个条件良好的矩阵，对其大小进行部分主元变换后得到完全不准确的逆矩阵。尽管完全旋转并不总是最好的策略，但我们还没能构造出一个例子来证明它是一个非常糟糕的策略，无论如何，第16节的分析给出了误差的上限。

这个例子也说明了平衡的重要性。如果将$B_n$ 的最后一列乘以$2^{-n}$ ，那么部分主元法给出的关键选择与完全主元法相同，结果将是灾难性的。(见第29节)

谱范数的一些基本性质

9. 许多结果将需要重复使用，在这里组装它们很方便。我们首先给出我们使用的谱范数的属性。 我们假设$A$ 自始至终都是实数。

（A）定义：$\Vert A\Vert =\max [\lambda \left(A{A}^T\right){]}^{\frac{1}{2}}$

（B）$\Vert A\Vert =\max \Vert Ax\Vert /\Vert x\Vert ,\Vert x\Vert \ne 0$ ，因此$\Vert Ax\Vert \leqq \Vert A\Vert \Vert x\Vert$

（C）$\Vert A\Vert$ 满足矩阵范数的一般关系。

（D）如果$A^{\left(r\right)}$ 是$A$ 的任意主子矩阵，则$\Vert A^{\left(r\right)}\Vert \leqq \Vert A\Vert$ (根据(B))

（E）$\Vert a{b}^T\Vert =\Vert a\Vert \Vert b\Vert$ (根据(B))。

（F）如果$A$ 是对称的，那么$\Vert A\Vert =\max |{\lambda }_i\left(A\right)|$ 。

（G）如果A是正定的，$\Vert A\Vert =\max {\lambda }_i\left(A\right)={\lambda }_i$ ，$\Vert A^{-1}\Vert =\max 1/{\lambda }_i\left(A\right)=1/{\lambda }_n$ 。

（H）$\Vert A\Vert =\Vert A^T\Vert$ 。

（I）$\Vert A{A}^T\Vert =\Vert A{\Vert }^2$ 。

（J）$n{\max }_{i,j}|a_{ij}|\geqq \Vert A\Vert \geqq {\max }_{i,j}|a_{ij}|$ 。

下面的结果是误差分析的基础，将反复使用。

这些结果是由恒等式得出的

${\left(I+A\right)}^{-1}$ 存在于$\left(Q\right)$ 之后。 写作${\left(I+A\right)}^{-1}=R$ ,

因为$(1-\Vert A\Vert )$ 是正数，所以得到(i)。(9.4)给出

从而得到(ii)；从(ii)我们推断

从而得到(iii)。

(S) 如果$A$ 是非奇异的并且$\Vert A^{-1}E\Vert <1$ ，那么$\left(A+E\right)$ 是非奇异的，并且我们有

因此

其中$F=A^{-1}E$ ，使得$\Vert F\Vert <1$ 。

证明如下。我们可以写作

因此

${\left(I+F\right)}^{-1}$ 存在于$\left(Q\right)$ 之后。 这给出

这给出了 (i)，并且 (ii) 是一个直接的结果。

仅当${\left(A+E\right)}^{-1}$ 可被视为$A^{-1}$ 的合理近似时才有价值。 通常我们不会知道乘积$A^{-1}E$ 有什么特别之处，因此必须使用不等式来简化$F$

这样做的必要性很可能导致我们估计

远高于其真实价值。 当进行替换 (9.12) 时，(i) 变为

我们将左边的表达式称为${\left(A+E\right)}^{-1}$ 中的相对误差。 请注意，我们得到的这个相对误差的界限包含$\Vert A^{-1}\Vert$ 本身作为一个因素，而绝对误差的界限包含$\Vert A^{-1}{\Vert }^2$ 。这是否代表真实的事态取决于关键取决于 (9.12) 是否是一个相当尖锐的不等式。

表2中的示例说明了最后一点。我们给出了一个病态矩阵$A$ ，并展示了在$A$ 中两种不同扰动$E_1$ 和$E_2$ 对逆矩阵的影响。这两种扰动具有相同的数量级，但$E_{{}_1}$ 导致了$\Vert A^{-1}{\Vert }^2\Vert E_1\Vert$ 阶扰动和$\Vert A^{-1}\Vert \Vert E_2\Vert$ 阶扰动。当$E$ 可以用这种形式表示时，扰动E的行为将与$E_2$ 的行为相同

其中$G$ 与$E$ 的数量级相同。稍后我们将看到，在误差分析中会遇到这种扰动。 在 (S) 的应用中，我们通常会从头开始限制自己，例如

因为否则我们将无法证明计算出的逆函数是真实逆函数的合理近似。

在这个条件下，(i)和(ii)变为

(T) 最终结果将在几乎所有分析中以某种形式使用。 如果 A 是非奇异的并且

其中$a、b和c$ 是标量，如果$a\Vert A^{-1}\Vert <1$

如果右边用$k\Vert A^{-1}\Vert$ 表示

然后在表格中写下(9.17)

并对双方进行规范。假设所有$a\Vert A^{-1}\Vert$ 和 ($b\Vert A^{-1}\Vert +c$ ) 明显小于单位，(i) 和 (iii) 表明$Y$ 是$A^{-1}$ 的良好近似。 在应用中，$P$ 总是出现在（9.17）中，但$Q$ 或$R$ 或两者都可能不存在。

我们再次强调，如果$\Vert A^{-1}\Vert \Vert P\Vert$ 是$\Vert A^{-1}P\Vert$ 的一个很差的近似，那么结果将是悲观的，如果$P$ 可以用$AS$ 的形式表示，其中$S$ 为与我们$P$ 相同的数量级（当然，这句话中隐含了一个归一化的$A$ ）。

三角方程组的解

10. 我们分析的所有反演方法都需要求解具有系数三角矩阵的方程组。 这些集合在定点算术中的解通常需要某种形式的缩放，所以我们首先考虑浮点解。 取上三角方程组

通常，$x_r$ 是使用关系从计算的$x_{r+1},\dots ,x_n$ 计算出来的

在[8]的符号中。 从中给出的内积结果，我们得到

其中每个$0$ 表示一个值，通常每个因子都不同，满足$|\theta |\leqq 1$ 。

因此计算的解满足

其中，对于$n=5$ ，

我们可以写作

对于$n$ 组方程的解

对于计算出的$X$ ，我们得到

因此

如果我们仅假设 U 具有归一化 (I)，那么我们可以通过替换所有 u1 来获得 V 的上层纠缠。 通过统一。 我们得到

其中我们将$V$ 的第$r$ 列中的所有倍数替换为$\left(n+3-r\right)$ 。

我们推断$\Vert V{\Vert }_E\sim \frac{n^2}{\sqrt{12}}$ ，事实上

（这里和其他地方的结果都成立，当$n<10$ 时保持的结果并不比$n\geqq 10$ 时的结果更糟，但我们不太关心低阶矩阵。）因此

当$B=I$ 时，我们发现计算出的$X$ 恰好是三角形的，并且没有对应于$\Vert B\Vert$ 的项现在出现在 (10.13) 中。 这是因为逆的第 $r$ 列是$U_x=e_r$ 的解，因此 (10.5) 中$\delta b$ 对应的误差是零向量。 除了对应于$e_n$ 的解之外，我们不需要整个优化矩阵$V$ ，但我们不会试图考虑这一点。因此，对于$UX=I$ 的计算解，我们有，

以及第 9 节的 (T)，

只有当$(0.4)2^{-t}{n}^{5/2}\Vert U^{-1}\Vert$ 明显小于1时，我们才能保证$X$ 是一个好的逆矩阵。

另一方面，如果我们知道（正如我们经常会知道的那样）

那么我们得到

求逆的结果是

11. 定点运算的相应结果不太确定，因为它们取决于运算工具的细节以及我们是否有关于解大小的先验信息。 为了明确起见，我们假设标量积可以在双精度累加器中精确累加，并且可以划分到这个双精度累加器中。 所有桌面机器和许多高速计算机都提供了这些设施。

如果我们事先知道计算解在$2^{k-1}$ 和$2^k$ 之间，那么我们可以写作

求解

变量$y_r$ , 由计算出的$y_{r+1},\cdot \cdot \cdot ,y_n$ 通过以下关系确定

其中$\hslash$ 表示我们取括号中表达式的定点结果。分子是精确累加的，唯一的误差是除法的舍入误差。(我们默认标量积的中间值不需要引入比最终值更大的比例因子。在ACE中，这个困难可以通过精确累加$2^{-p}$ 倍的标量积来避免，其中$2^p$ 大于可以处理的最大阶矩阵。因此，在标量积中不会出现溢出，引入的舍入误差($2^{p-2t}$ )远低于其他的水平。我们不再提及这一点。)因此

给出

应用于$UX=B$ 的解，我们得到一个$X$ 满足

其中

当$B=I$ ，我们得到

其中 (11.9) 得到满足，给出

在这一点上，我们默认$k\geqq 0$ 。 如果$k<0$ ，则不需要比例因子，并且 (11.9) 变为$\Vert C\Vert <n{2}^{-t}$ 。 对于反演问题$k$ 必然满足$k\geqq 0$ ，因为$X$ 的一个元素为($1/u_{nn}$ )且不小于1。 (11.5) 的最严格应用给出

12.如果事先不知道比例因子，随着计算的进行，$X$ 的元素会在必要时被密封。如果$X$ 的第一个计算元素需要整个计算所需的最大比例因子，那么分析就是我们刚才给出的。如果使用了大小的完全旋转，这通常会很容易。否则情况将是这样。假设在即将计算$x_r$ 的阶段，当前比例因子为$2^{k_{r+1}}$ ，我们可以将变量$x_i$ 的当前密封值表示为$x\genfrac{}{}{0ex}{}{\left(r+1\right)}{i}$（$i=r+1,\cdot \cdot \cdot ,n$ ） ，因此

下一步是计算$x\genfrac{}{}{0ex}{}{\left(r\right)}{r}$ 并在必要时引入增加的比例因子$2^{k_r}$ ，我们得到

现在可能有必要在任何或所有后续阶段引入更高的比例因子，以便通过四舍五入$x_r,x_{r+1},x_{r+2},\cdot \cdot \cdot ,x_n$ 获得这些变量的最终值$x\genfrac{}{}{0ex}{}{\left(1\right)}{r}\cdot \cdot \cdot x\genfrac{}{}{0ex}{}{\left(1\right)}{n}$ 到逐渐变少的数字。 可能发生的最坏情况是发生了 ($k_1-k_{r+1}$ ) 连续舍入。 如果那我们写作

那么

因此，最终值满足

右边的误差项在$(n+1-r)2^{k_1-t}$ 范围内。这比在整个过程中使用最终比例因子$2^{k_t}$ 时获得的因子$\left(n+1-r\right)$ 更差。尽管如此，我们认为逐步引入比例因子通常会得到更精确的结果，特别是对于条件非常病态的矩阵。我们可以通过考虑矩阵$A$ 的求逆运算来说明这一点，如下所示，运算到五位数。

当使用递进的比例因子计算逆矩阵的最后一列时，我们对第一个实例的最后一个元素得到884.96。它随后被四舍五入为885，然后是$9\times 10^2$ ，最后一列是$10^2\left[74552,975,9\right]$ 。如果我们从一开始就引入正确的比例因子，我们将计算出最后一列$10^2\left[75818,992,9\right]$ 。 很明显，前者的答案更好，尽管它的残差大约是后者的 100 倍。

矩阵乘法

13.对于浮点矩阵乘法，我们写作

根据浮点标量积的结果，我们得到

与${ϵ}_i$ 的通常界限。 因此，对于 4 阶矩阵，我们通常有：

因此自然

且

在本文中唯一的应用中，我们有$B=A^T$ 且

在表格中写作

我们看到，计算出的乘积的相对误差的界与乘积无关。请注意，即使是计算出来的$C$ 也是完全对称的，因此我们只需要构造$C$ 的上半部分。

对于$A{A}^T$ 的定点计算，即使比例因子是逐步引入的，我们也有

其中$2^{k_1}$ 是最终的比例因子。如果先计算对角线上的元素，则会在初始阶段确定最终的比例因子，因为最大的元素在对角线上。由式(13.7)给出

并且因此

因为$\Vert A{A}^T\Vert \geqq \max |{\left(A{A}^T\right)}_{ij}|\geqq \frac{1}{2}\cdot 2^{k_1}$ 。 这是最令人满意的误差范围。

高斯消元法的一般误差分析

14.  现在考虑一个一般的矩阵$A$ 的求逆，通过使用一个完整的关键策略的高斯消元法将其还原为三角形式。这种策略只是决定了$A$ 的行和列的重新排序，如果我们认为$A$ 以这种方式被重新排序，那么在每个阶段主元素都在左上角，它将简化符号，而不会失去任何通用性。应该强调的是，这种重新排序是事先不知道的。

我们首先给出一个与所使用的算法类型无关的一般分析。我们用$a\genfrac{}{}{0ex}{}{\left(r\right)}{ij}$ 表示第$r$ 个矩阵$A^{\left(r\right)}$ 的计算元素，用$m_{ij}$ 表示计算的乘数。$\left(i,j\right)$ 位置的元素历史记录如下。

其中$ϵ\genfrac{}{}{0ex}{}{（k）}{ij}$ 是根据计算出的$a\genfrac{}{}{0ex}{}{\left(k-1\right)}{ij},m_{i,k-1}$ 和$a\genfrac{}{}{0ex}{}{\left(k-1\right)}{k-1,j}$ 计算$a\genfrac{}{}{0ex}{}{\left(k-1\right)}{ij}-m_{i,k-1}a\genfrac{}{}{0ex}{}{\left(k-1\right)}{k-1,j}$ 时产生的误差。元素$a\genfrac{}{}{0ex}{}{\left(i\right)}{ij}$ 是第$i$ 个关键行的元素，不会发生进一步变化。对方程 (14.1) 求和并消去公共元素，

(ii)$i>j$ 。 我们现在有与 (14.1) 相同的方程，

下一个值$a\genfrac{}{}{0ex}{}{\left(j+1\right)}{ij}$ 被指定为零，并且在后续的归约过程中保持为零。元素$a\genfrac{}{}{0ex}{}{\left(j\right)}{ij}$ 用于计算$m_{ij}$

其中${\eta }_{ij}$ 是除法中的舍入误差。因此

其中

将(14.1)-(14.3)式求和并加上(14.5)，

注意，方程(14.2)和(14.7)给出了最终三角矩阵的元素$a\genfrac{}{}{0ex}{}{\left(r\right)}{rs}$ 、元素$m_{ra}$ 和添加了扰动的元素$a\genfrac{}{}{0ex}{}{\left(1\right)}{rs}$ 之间的关系。 将$L$ 表示由 $m_{ij}$ 形成的下三角矩阵，增加单位对角线，将$U$ 表示由关键行形成的上三角矩阵，(14.2) 和 (14.7) 给出

其中$E^{\left(k\right)}$ 是由$ϵ\genfrac{}{}{0ex}{}{\left(k\right)}{ij}$ 构成的矩阵。注意，它的第1行到$\left(k-1\right)$ 为空，第1列到$\left(k-2\right)$ 为空。

我们可以这样描述这个结果。“计算得到的矩阵$L$ 和$U$ 是通过正确计算$A^{\left(1\right)}+E$ 得到的矩阵。此外，$A^{\left(k\right)}$ 是通过精确计算得到的矩阵

浮点数中的高斯消元法

15. 我们现在假设矩阵是按比例缩放的，因此所有元素都略小于$1/g(n)$ [见第4节]。严格来说，对于浮点计算来说，这种扩展是无关的，但它有助于与其他地方的定点计算进行比较和误差分析。我们还假设矩阵是平衡的，尽管没有利用这一事实。我们认为这是必要的，以证明完整的关键战略。我们知道，对于精确计算，$A^{\left(k\right)}$ 的元素将被限制在$g(k)/g(n)$ 之间，因此每个$A^{\left(k\right)}$ 的每个元素都小于1。我们假设计算出的$A^{\left(k\right)}$ 也是如此。

我们可能会觉得最后一段的结果保证了这一点，因为最终的$L$ 和$U$ 确实对应于$A+F^{\left(n\right)}$ 的精确计算。不幸的是，每个$A^{\left(k\right)}$ 对应于对应的$A^{\left(1\right)}+F^{\left(k\right)}$ 的精确计算，因此，尽管我们可以声称对$A^{\left(1\right)}+F^{\left(n\right)}$ 进行了精确的高斯消元，但不能说我们已经选择了对应于这个固定矩阵的每个阶段的最大主元。我们希望稍后再讨论这一点，但目前这仍然是一个假设。我们明确地假设

那么

请注意，$|m_{ij}|\leqq 1$ 是确定的，因为我们在每个阶段都将计算出的$A^{\left(k\right)}$ 的最大元素作为基准值。

现在我们得到

因此

说明，在我们的假设 (15.1) 下，$\left(A+E\right)$ 的所有元素都小于$1/g(n)$ 。

因为我们假设所有$a\genfrac{}{}{0ex}{}{\left(k\right)}{ij}$ 都小于1，所以很明显

因此，矩阵$E^{\left(k\right)}$ 的边界通常是$n=5，k=3$ ，

并且$E$ 为

使用欧几里得范数作为近似值，我们得到

这些结果都代表了极端的上限，即使在减少过程中发生了最大可能的增长，也无法获得。通常情况下，随着过程的进行，元素的大小有所减小，则舍入误差肯定小于$\left(|a\genfrac{}{}{0ex}{}{\left(k-1\right)}{ij}|+|a\genfrac{}{}{0ex}{}{\left(k-1\right)}{k-1,j}|\right)2^{-t}$ 在每个阶段，逐渐变小。例如，如果元素在每个阶段都减小2倍，则 E 边界的最后一个元素变为$2\times 2^{-t}/g(n)$ 而不是 $\left(n-1\right)\times 2^{-t}$ 。

仍然存在由于$A^{\left(k\right)}$ 为零而导致进程过早中断的可能性。只有当$A+F^{\left(k\right)}$ 正好是奇异值时才会发生这种情况。现在$\Vert F^{\left(k\right)}\Vert \leqq \left(0.41\right)2^{-t}{n}^2$ ，因为 F(k) 的每个元素肯定受我们对 E 的限制的限制。因此，如果$\left(0.41\right)2^{-t}{n}^2\Vert A^{-1}\Vert <1$ ，则不会发生分解。如果我们要保证一个好的逆，我们在任何情况下都需要一个更强的条件。

16. 求逆运算通过求解完成

在这两个步骤中执行哪一个

将上述结果应用于三角方程组的浮点解，我们最终得到：

其中根据（15.8）

根据（10.14）

根据（10.13）

因此

或

我们需要对右边每一项的范数有满意的界限。 很明显，除非$2^{-t}{n}^{7/2}\Vert A^{-1}\Vert$ 明显小于 1，否则我们无法保证合理的逆。 因此我们假设

而这个条件当然可以保证在还原过程中不会出现故障。 我们进一步假设$n\geqq 10$ 。 联系 (16.4) 意味着$A$ 是非奇异的。 等式 (a) 和 (b) 给出

现在根据（a）

因为$U$ 的所有元素都小于1。

根据（c）和（d），

应用这些结果并记住$\Vert L\Vert \leqq 0.8n$ ，因为它的所有元素都被1约束，我们得到

我们现在已经建立了第9节的类型(T)结果

因此

我们可以很容易地推导出元素以1为界限的矩阵的结果。如果$B$ 是这样的矩阵，那么我们使用$A=B/g\left(n\right)$ 。通常情况下，我们会除以大于$g(n)$ 的最小的 2 次方，这样就不会涉及四舍五入。如果$X$ 是$A$ 的逆矩阵，那么将$Z=X/g\left(n\right)$ 作为$B$ 的逆矩阵，然后我们得到

因此

乍一看，术语$n^{7/2}$ 的出现令人失望。 但是，我们必须记住，Goldstine 和 von Neumann 对于正定对称矩阵的结果包含因子$n^2\Vert A\Vert \Vert A^{-1}\Vert$ ，并且$\Vert A\Vert$ 可能与$n$ 一样大。 在我们分析的相应位置，我们用$0.8n$ 替换了$\Vert L\Vert$ 和$\Vert U\Vert$ ，这在我们的结果中贡献了一个因子$n$ 。 因此，我们只有附加因子$n^{\frac{1}{2}}$ ，这个结果是针对一般矩阵的。

因子$g(n)$ 的存在是不可避免的，但请注意它应该被$\max \left[1,R\right]$ 代替，其中$R$ 是增长率，在每个特定情况下，在实践中$R$ 通常小于 1。 对于我们的每一种特殊类型的矩阵，我们都有一个$R$ 的先验上限，它比$g(n)$ 小得多。

定点计算

17. 我们不希望对所有类型的不动点计算进行重复分析，我们只介绍一种特定的类型。在描述这一点之前，我们来证明它的选择是合理的。

求解已执行误差分析的非对称方程的替代方法是

(i)将问题转化为正定矩阵[3]的$A{\Lambda }^T$ 求逆问题，

(ii) 通过不同类型的正交三角化 (§21-23)。

这些方法都不需要少于$2n^3$ 次乘法，而我们刚刚描述的方法需要$n^3$ 次。 如果我们准备使用$2n^3$ 乘法，则以下策略可提供较高的准确性。

(I)对上述矩阵求逆，从而确定交换和所需的比例因子。

(II)利用(I)中获得的信息，在定点重复计算，尽可能累积标量积。

在 (II) 中，我们直接从$LU=A$ 确定$L$ 和$U$ ，其中$\Lambda$ 的行和列按 (I) 确定的顺序排列。 我们应该计算 $½L$ 以防某些乘数由于不同的舍入误差而刚好大于 1。 由于这项方法，我们得到

因为$½L$ 和$U$ 的每个元素都是通过除以精确累加的标量积来直接确定的。 类似地，在求解$LY=I$ 和$UX=Y$ 时，比例因子是预先已知的。 我们给出的一些边界将根据$ma{x}_{ij}|x_{ij}|$ 大于或小于$ma{x}_{ij}|y_{ij}|$ 而不同。如果是这样，前一种情况下的边界在左侧给出，而后者在右侧给出。 我们得到

给出

前提是$2^{-t}{n}^2\Vert A^{-1}\Vert \leqq 0.1$ 。 （常数只保留一位小数。）在这两种情况下，我们都得到第 9 节类型 (T) 的结果，其中

在实践中，第二种情况并不常见，因为$X$ 的最大元素通常比Y的最大元素大，特别是当$A$ 处于病态状态时。对于一些特殊的矩阵，我们已经知道${\max }_{i,j}|y_{ij}|={\max }_{ij}|l\genfrac{}{}{0ex}{}{-1}{ij}|$ 等于1。在任何情况下，${\max }_{i,j}|y_{{}_{i,j}}|\leqq 2^{n-1}$ ，因此对于小矩阵$\Vert Y\Vert$ 不能很大。

因此我们通常会有

在另一种情况下，我们甚至不能保证$X$ 是一个近似的逆，除非$n^3{2}^{-t}\Vert {\Lambda }^{-1}\Vert \leqq 0.1$ ，然后我们得到

对称ChoIesky分解的误差分析

18. 我们省略了正定矩阵的高斯消元分析，现在可以在没有关键策略的浮点运算中完成。 分析类似于上一节的分析，除了我们知道关键元素的大小没有增长，并且$L$ 和$U$ 的范数与$A$ 的范数更特别相关。我们转而采用在其他情况下很重要的$L{L}^T$ 分解。

我们考虑对正定矩阵$A$ 进行$L$ 的定点计算。我们证明如果

那么$L$ 的所有计算元素都满足$|l_{ij}|<1$ 和

证明使用的是双重归纳法。假设当我们计算$L$ 的前$s$ 行和$\left(s+1\right)$ 第$r$ 行元素时，计算得到的元素都是1，并且是由矩阵精确分解($A+E_{sr}$ )得到的。这里的$E_{sr}$ 是一个对称矩阵，除了它的$s$ 阶主矩阵和它的$\left(s+1\right)$ 第一行和列的前$r$ 个元素中有边界为(18.2)的元素外，其他元素为零。对于$n=5,s=3,r=1$ ，我们通常有：

由式(18.1)和我们的假设，$A+E_{rs}$ 是一个元素以1为界的正定矩阵。

现在考虑确定$l_{s+1,r+1}$ 。 计算它的方程式是

如果这个商被精确地计算出来，它将是对应于$\left(A+E_{rs}\right)$ 的精确分解的$\left(s+1,r+1\right)$ 元素，并且由于这是一个正定矩阵，其元素以1为界，所以精确商以1为界。 因此，计算出的$l_{s+1,r+1}$ 也受单位限制，并且

因此

这样就完成了我们对元素$l_{s+1,s+1}$ 的归纳。 这是根据关系计算的

同样，如果精确计算平方根，它将小于1，因为它是$A+E_{s+1,s}$ 的精确分解中的一个元素，它是正定的，元素以单位为界。

我们不想深入研究平方根子程序。在 ACE 上，平方根是通过双精度累加器除以牛顿过程计算的，在四舍五入时要注意，计算值是正确四舍五入的值

我们将在我们的分析中假设这一点，尽管我们有

对于任何合理的$k$ ，$A$ 的条件不需要从 (18.1) 做太多改变。 由于计算值是正确舍入的值，

对于$1\times 1$ 矩阵，结果显然是正确的，因此归纳完成。注意，如果有任意一个$i$ 很小(并且如果$A$ 是病态的，几乎可以肯定这是正确的)，则误差矩阵($L{L}^T-A$ )的对应元素远小于$2^{-t-1}$ 。换句话说，当$A$ 的元素不能精确地用$t$ 个二进制数字表示时，分解产生的误差可能远小于对$A$ 的元素进行初始舍入所产生的误差。$L{L}^T$ 分解是应用于对称矩阵的 Rutishauser的$L-R$ 技术的基础，也是当$A$ 为正定时计算$\det (B-\lambda A)$ 和 $\det (AB-\lambda I)$ 的零点的基础。 因此，$L{L}^T$ 与$A$ 的接近程度具有重要的实际意义。

逆运算的完成

19.现在可以通过求解来完成反演

然后计算$X{X}^T$ 。在这部分工作中，我们不能再避免缩放，我们给出了浮点计算的分析，因为这一阶段所犯的错误并不重要。利用三角矩阵求逆的结果，并利用$L$ 中任意一行元素的平方和不大于一的事实(第7节)，我们得到(10.18)：

计算出的逆$Y$ 满足

根据(13.5)。我们只需要计算$Y$ 的上半部分，因为即使是计算出来的$Y$ 也是完全对称的。从$A$ 计算$Y$ 所涉及的乘法总数为$\frac{1}{2}{n}^3$ ，因此该方法充分利用了对称性。

连同等式 (19.2) 和 (19.3)

使我们能够评估 Y 的准确性。由于$L$ 和$L^{-1}$ 的范数与$A$ 和$A^{-1}$ 的范数之间的密切关系，因此直接计算相对误差比评估$\Vert AY-I\Vert$ 更简单；我们在第 20 节对此进行评论。如果我们假设它会简化评估

其中第一个或第二个更具限制性取决于$\Vert A^{-1}\Vert$ 的大小。在任何情况下，第一种方法都足以确保在计算$L$ 时不发生故障，并且其他任何地方都不会出现故障。

我们将A的对称性和(19.2)到(19.5)的关系的推论归纳在公式(19.6)中：

因此

条件 (19.5) 允许我们以多种方式表达$\Vert F{\Vert }^2$ 中的项，但使用与$H$ 的主要成分相同的形式更方便。最后我们有

结合(19.6b)、(19.8)和(19.9)，$Y$ 的相对误差为：

其中右边的三项分别来自$L、X和Y$ 的计算。 通过在$X和Y$ 的计算中使用定点运算，我们可以将最后两项中的$n^2$ 减少到$n$ ，但是如果$A$ 非常病态，以至于$\Vert A^{-1}\Vert$ 很大，那么无论如何，三角化的误差是最重要的，而计算$X{X}^T$ 的误差是最小的。

我们可以先解$LX=I$ ，然后再计算$X^{T}X$ ，而不是$L^{T}X=I$ 。在这两种情况下，相对误差的界是相同的，当然，计算出的逆函数是不同的。

残差矩阵

20.在大多数误差分析中，我们处理的是残差矩阵$\left(AY-I\right)$ 。值得注意的是，如果我们求解$L^{T}X=I$ ，我们可以获得比求解$LX=I$ 更好的矩阵边界。此外，计算的各个部分对残差矩阵的贡献与对相对误差的贡献并不相同。

使用$L^{T}X=\left(I+F\right)$ ，我们有

由$E、F、G$ 产生的分量范数的界限的数量级分别为

包含$F$ 和$F^T$ 的项做出相同的贡献，而$F{F}^T$ 项的贡献可以忽略不计。 请注意，$G$ 对残差的贡献可能很容易是最大的，即使它对逆本身的影响最小也是如此。

另一方面，如果我们使用

那么F和G是有界的，但是

包含$F^T$ 和$F$ 的项的数量级现在是

因此，如果A是病态的，那么前一项可能是严重的。然而，除了特殊构造的矩阵，这两项之间的差异在实际中并不会发生。这是由三角矩阵逆的性质决定的，将在第27-28节讨论。

正交化简为三角式

21. 我们现在考虑通过与正交矩阵的预乘将矩阵化简为三角式。 在这种类型的方法中，原始矩阵被连续变换为 ，其中除了前一个矩阵中的那些之外，序列的每个成员在对角线下方至少有一个额外的零。 这些零可以一次产生一个 (Givens)，或者一次产生整个子柱 (Householder)。 无论该方法的细节如何，误差分析都具有某些共同特征，正如我们现在将要展示的那样。

我们可以这样定义数学过程。从

生成

其中

矩阵$R_r$ 是连续正交矩阵。

我们现在考虑实际的过程，用$A_r$ 表示第$r$ 个计算矩阵。与此对应的$A_r$ (注意，这不是前面阶段中通过精确计算得到的矩阵)，该技术定义了一个精确的正交矩阵$R_r$ 。计算得到的矩阵${\bar{R}}_r$ 与此不同，我们可以写成

然后，我们将$A_r$ 乘以$\overline{R_r}$ ，并进一步计算误差，使计算出的$A_{r+1}$ 和$B_{r+1}$ 满足：

请注意，我们将那些本应为零的元素取为恰好等于零，并且$F_r$ 必须包括该假设中隐含的误差。 如果$F_r$ 要小，那么${\bar{R}}_r{A}_r$ 必须至少有小元素来代替那些应该为零的元素。

我们现在可以用两种不同的方法进行分析。

(A)我们可以写作

并假设我们可以找到$\Vert H_r\Vert$ 和$\Vert K_r\Vert$ 的界限

其中$h$ 和$k$ 是$n$ 和$r$ 的简单函数。

最后一组方程是

一个三角形集，计算出的$X$ 满足

预乘$R_1^T{R}_2^T\cdot \cdot \cdot R\genfrac{}{}{0ex}{}{T}{s-1}$ ，我们有

其中

因为所有的$R_i$ 和$S_i$ 都是正交的，我们得到

(B)我们可以使用${\bar{R}}_r$ ，并假设$F_r$ 和$G_r$ 具有这种形式的边界

将(21.9)先乘$\bar{R}\genfrac{}{}{0ex}{}{-1}{1}\bar{R}\genfrac{}{}{0ex}{}{-1}{2}\cdot \cdot \cdot \bar{R}\genfrac{}{}{0ex}{}{-1}{s-1}$ 我们得到

其中

现在$S_i$ 不再是完全正交的，但假设我们可以证明$\Vert E_r\Vert \leqq a<1$ 。那么我们有

因为$R\genfrac{}{}{0ex}{}{-1}{r}$ 是正交的。因此，我们就有了与(21.12)相同形式的关系：

现在观察，我们可以允许$a$ 的值大大大于$2^{-t}$ 而不会损失太多。 例如假设$a=10^{-4}$ 。 现在在 Householder 方法中$s=n$ ，因此即使$n=10^4$ ，因子$1/{\left(1-a\right)}^{s-1}$ 也仅等于$e$ 。 因此我们看到，尽管将$A^r$ 与${\bar{R}}_r$ 预乘时产生的误差应保持尽可能低（$2^{-t}$ 阶），并且确切的${\bar{R}}_r{A}_r$ 应在应产生零的位置具有$2^{-t}$ 阶的元素 ，我们可以容忍矩阵${\bar{R}}_r$绝不是正交的工作精度。 我们对${\bar{R}}_r$ 的真正要求是它的范数和$\bar{R}\genfrac{}{}{0ex}{}{-1}{r}$ 的范数都应该合理地接近统一。 (B)分析显然更深入。

这与在对称特征问题中使用正交变换的情况形成了对比。在这里，我们与${\bar{R}}_r$ 做前乘，用$\bar{R}\genfrac{}{}{0ex}{}{T}{T}$ 做后乘，我们假设得到的矩阵仍然是对称的。除非${\bar{R}}_r$ 是准确正交的，否则这不是一个准确的相似变换。

22. 如果分析要给出令人满意的结果，那么根据我们使用 (A) 或 (B)，我们应该能够为$H$ 和$K$ 或$F$ 和$G$ 的范数获得令人满意的界限是至关重要的。 前者意味着为$E_r$ 找到非常好的边界，但后者在这方面要求不是很高。 现在，在高斯消元法中，找到边界的主要障碍是减少矩阵的元素大小逐渐增加的可能性，因为通常当前舍入误差与元素的当前大小成正比。 对于正交（或接近正交）矩阵，这个困难不再存在。$A_1$ 的每一列的欧几里德范数通过精确的正交变换得以保留。 因此，如果最初每列中元素的平方和小于 1，则这在所有后来的矩阵中仍然成立。

在实践中，在定点工作时，我们必须缩放$A_1$ ，使舍入误差的累积不会导致任何$A_r$ 的任何列的范数超过1。出于同样的原因，我们必须取$B_1=\frac{1}{2}I$ 。对于浮点数，我们不需要采取这些额外的预防措施，但在估计$H_r$ 和$K_r$ 或$F_r$ 和$G_r$ 的边界时，需要处理可能出现的规范的轻微增加。

对Givens Reduction的详细分析

23. 我们现在对Givens Reduction做一个简要的分析。 我们在浮点中工作，因为在定点工作中积累标量积所带来的收益并不大。 我们取$A_1$ 原本进行归一化 (II)。

考虑归约中的一个典型阶段，在该阶段，$r$ 列已变为三角形且元素 $\left(r+2,r+1\right);\cdot \cdot \cdot$ ； 第 $\left(r+1\right)$ 列中的$\left(k,r+1\right)$ 已设为零，因此通常对于$n=5、r=1、k=3$ ，当前矩阵的形式为

下一步是通过在$\left(r+1,k+1\right)$ 平面中旋转使$\left(r+1,k+1\right)$ 元素为零。 这仅修改行$\left(r+1\right)$ 和$\left(k+1\right)$ 中的元素，并保持前面的零不变。在这一阶段，我们假设$A$ 和$B$ 所有列的欧氏范数有界

其中$p$ 是到目前为止已经完成的变换次数，并在下一阶段的分析中验证这个假设。

旋转矩阵的分量可以根据这些关系计算出来

其中所有$|{ϵ}_i|\leqq 2^{-t}$ 并且我们假设了一个平方根子程序，它给出的误差小于$2^t$ 中的 2 个部分。 因此计算的$\cos \theta$ 和$\sin \theta$ 满足

在上一节的表示法中，我们有一些很粗略的不等式，

与$b$ 和$c$ 的大小无关。 因此，在第 21 节的符号中，我们可以取$a=6.2^{-t}$ ，很明显，偏离正交性并不重要。

现在我们转向F和G矩阵的确定。除了行$\left(r+1\right)$ 和$\left(k+1\right)$ 之外，这些都是空值，并且对于$F$ ，这两行从$1$ 到$r$ 的位置都是空值。

$F$ 的元素$\left(r+1,k+1\right)$ 必须特殊处理，但$F$ 和$G$ 的每个其他非零元素的边界由下式给出

或由

其中$x$ 和$y$ 是当前$A或B$ 的某列第$\left(r+1\right)$ 和$\left(k+1\right)$ 个位置的已有元素。对于${\delta }_1$ ，忽略$2^{-2t}$ ，

同样的结果也适用于${\delta }_2$ 。这使我们能够证明关于转换后的$A和B$ 的列范数的原始假设。因为此转换前后任何列的范数之比肯定介于$\left(1\pm 9\times 2^{-t}\right)$ 之间。这已经足够覆盖$\bar{R}$ 中的误差以及乘以它所产生的${\delta }_1$ 和${\delta }_2$ 的误差了。所有转换后的列的所有范数都位于$(1\pm 9.2^{-t}{)}^{\frac{n\left(n-1\right)}{2}}$ 之间。 这些范数将是统一的，因为我们在任何情况下都需要$n^2\times 2^{-t}$ 明显小于统一才能得到有用的逆。$\bar{R}A$ 的$\left(k+1,r+1\right)$ 元素的模量的精确值为

这肯定小于$2^{-t}$ 。 将其替换为零时，我们犯了一个错误，其边界小于我们为行$\left(r+1\right)$ 和$\left(k+1\right)$ 的其他元素获得的边界。 结合当前$F$ 和$G$ 的这些结果，

汇总这些结果并取${\sum }_1^r{S}^{\frac{1}{2}}=\frac{2}{3}{r}^{3/2}$ ，${\sum }_1^r{S}^{3/2}=\frac{2}{5}{r}^{5/2}$ ，根据（21.12）和（21.17）我们得到

我们用$1/{\left(1-a\right)}^{\frac{n\left(n-1\right)}{2}}$ 替换所有地方的$1/{\left(1-a\right)}^2$ 。我们暂时写

现在$M$ 是求解最终三角方程组的误差

因此，从(10.10)和(10.17)可知，计算得到的$X$ 满足

当然在那里

由于$A_s$ 和$B_s$ 每列的范数小于$\left(1+d\right)$ ，因此有

这是标准形式的结果。对于一个相当好的逆函数，我们需要

并假设$n\geqq 10$ 这肯定给出

我们将进一步要求

如果我们假设这一点，当$n>10$ 时我们肯定有

从第9节(T)我们得到

结果是一个归一化(II)的矩阵。对于归一化(I)，按照通常的方式，我们立即得到：

更危险的项是 (2.1)$2^{-t}{n}^3\Vert A^{-1}\Vert$ ，这主要来自在每个阶段计算$\bar{R}A$ 时产生的舍入误差。

我们不能通过定点工作获得太多，因为没有阶段可以积累任何大小的标量积。 在计算$\cos \theta$ 和$\sin \theta$ 时必须小心。最令人满意的技术是精确地累加$\left(b^2+c^2\right)$ ，然后按最大$2^{2k}$ 缩放它，其中$2^{2k}\left(b^2+c^2\right)<1$ 。 然后我们从$2^{k}b$ 和$2^{k}c$ 计算$\cos \theta$ 和$\sin \theta$ 。 计算出的矩阵比浮点计算更接近于正交矩阵，而且在计算$\left(x\cos \theta +y\sin \theta \right)$ 和$\left(-x\sin \theta +y\cos \theta \right)$ 时仅产生一个舍入误差。 精度的提高最多相当于一个小的常数因子，我们不能获得$\sqrt{n}$ 的因子。

Householder的减少

24. 只有$\left(n-1\right)$ 转换涉及到Householder减少。在第一步，我们乘以一个形式为$(I-2{\omega }_r{\omega }_r^T)$ 的矩阵，其中$w_r^T$ 是这种形式

这会在第$r$ 列的对角线下方产生零。 它保留$A和B$ 的前$\left(r-1\right)$ 行和$A$ 的前$\left(r-1\right)$ 列不变，因此$A$ 的${\left(n-r+1\right)}^2$ 个元素被修改， B 的$n\left(n-r+1\right)$ 个元素被修改。在进入下一阶段之前，我们被迫对这些元素进行舍入，这样我们就没有希望获得比$F_r$ 和$G_r$ 更清晰的边界

这给出了$N和P$ ，

即使我们让定义在第311页的$a$ 等于0。因此，如果我们充分小心地计算$(I-2w_r{w}_r^T)A_r$ 和$(I-2w_r{w}_r^T)B_r$ ，我们希望得到的最多是$\sqrt{n}$ 的因子，也许是一个小的常数因子。

如果我们使用传统的浮点运算，我们甚至无法得到$\sqrt{n}$ 因子。如前所述，很容易看出偏离正交性是轻微的，因此不重要，我们可以把注意力集中在计算$A_{r+1}$ 本身所产生的误差上。

我们假设我们有一个精确的$w_{{}_r}$ 并估计在计算 $A_r-2w_r{w}_r^{{}^T}{A}_r$ 时产生的错误。第一阶段是典型的，所以我们处理$A_1-2w_1{w}_1^{{}^T}{A}_1$ 。 我们有

$q_{1k}$ 可以写成这种形式

其中

现在$\Sigma x_i^2=1$ 和${\Sigma }_{i}a\genfrac{}{}{0ex}{}{2}{ik}=1$ （后者大约在后期阶段持续存在），但我们不能得出比这更明确的结论

在计算$\hslash \left(w_1{w}_1^T{A}_1\right)$ 时，仅来自$q_1$ 的贡献就是$\Vert w_1\Vert \Vert q_1\Vert$ ，我们对此唯一的限制是$n^{3/2}{2}^{-t}$ 。 类似地，在后面的阶段，我们有$r^{3/2}{2}^{-t}$ 量级的贡献，它们一起给出$n^{5/2}{2}^{-t}$ 项。 在实践中，在浮点计算中，将方程式表示为方程式 (24.7) 中所示的形式很方便（第一步是典型的）。

这种形式避免了不必要的平方根，但不会导致误差下界。

在定点运算中，我们确实可以实现因子$n^{\frac{1}{2}}$ 的增益，但只有在具有良好运算能力的计算机上才值得这样做。 为确保准确确定$w$ ，应准确累加$S_2$ ，然后乘以$2^{2k}$ ，其中$k$ 是$2^{2k}{S}^2$ 不大于零的最大整数。 然后可以使用$2^k{a}_{i1}$ 而不是$a_{i1}$ 来计算向量 $w$ 。 这给出了$w$ 的非常准确的确定，我们只需要关心计算中的错误 ($A_1-2w_1{w}_1^T{A}_1$ )。 在形成$p_1=\hslash \left(w_1^T{A}_1\right)$ 时，我们可以累加标量积，并且由于$A_1$ 的列具有以单位为界的范数，因此$w_1^T{A}_1$ 的元素小于单位。 如果我们写

那么$q_1$ 的所有元素都小于$2^{-t-1}$ ，因此$\Vert q_1\Vert \leqq \sqrt{n}{2}^{-t-1}$ 。 虽然$2w_1{w}_1^T{A}_1$ 的元素可能大于 1，但我们知道$A_1-2w_1{w}_1^T{A}_1$ 的元素没有统一，因此前者的溢出不会造成任何危害。 为了始终安全，我们必须最初对$A$ 进行缩放，以便其列充分小于 1，从而使$A_r$ 的列始终小于 1。 如果足够小心，我们可以获得计算出的逆$X$ 满足

其中$a和b$ 是统一的。 对于规范化(I) 的矩阵，这给出

非对称矩阵求逆方法的比较

25. Givens 方法需要四倍于 Gaussian 反演的乘法运算，因此在单精度 Givens 反演和双精度 Gaussian 之间进行公平比较。 在此基础上，高斯解要好得多。 另一方面，Householder 方法只需要两倍于高斯反演的乘法。 然而，我们已经看到 [(17.7) 和第 9 节，(T)]，如果我们准备执行两次高斯消元，我们可以期望逆满足

这比 Householder 更有利，即使在计算时非常小心，前提是“R”（枢轴比）小于 $\sqrt{n}$ 。 根据我们的经验，这几乎总是正确的，并且有许多重要的矩阵类型，我们事先知道它是正确的。

比较高斯消元法和Householder，当两者都以浮点数执行时，前者的相对误差为$2^{-t}{n}^{7/2}r\Vert A^{-1}\Vert$ ，后者的相对误差为$2^{-t}{n}^3\Vert A^{-1}\Vert$ 。在我们看来，因素$n^{7/2}$ 的存在源于分析的紧迫性，而不是源于任何真正的缺点，Householder 的案例很少。

我们进行误差分析的另一种方法是Goldstine和von Neumann的基于正定矩阵$A{A}^T$ 求逆的方法。然而，这种方法需要$2n^3$ 次乘法，并且在误差方面，它用$\Vert A^{-1}{\Vert }^2$ 代替$\Vert A^{-1}\Vert$ 。如果$A$ 完全是病态的，则此方法没有任何可取之处。

对结果的一般评论

26. 我们获得的界限在所有情况下都是严格的上限。 一般而言，舍入误差的统计分布将大大减少出现在相对误差中的$n$ 的函数。 我们可能期望在每种情况下，该函数都应该被替换为不大于其平方根并且通常明显更小的东西。 例如，在形成浮点标量积时，我们使用了结果

然后用它们的上限$\left(n-i\right)2^{-t}$ 替换了${ϵ}_i$ 。 现在要达到这些界限，不仅所有单独的舍入误差都必须具有最大值，而且$a_i{b}_i$ 必须具有特殊分布。

另一方面，总是出现在相对误差中的项$\Vert A^{-1}\Vert$ 通常不受任何此类考虑。 有一些简单的标准来判断相对误差是有用的，这是由以下考虑提供的。 考虑对$A$ 的元素中$2^{-t}$ 阶扰动的逆的影响。我们有

如果这些扰动不是特别相关，那么我们可以写

并期望右侧是左侧的合理近似值。

因此我们有

我们期望能够选择某个$E$ ，其右侧不会被严重高估。如果求逆的方法给出了$X$ 的逆

我们可以说，逆矩阵中的误差是由$A$ 的单个元素中的$\left(1/n\right)f(n)2^{-t}$ 扰动引起的误差的数量级。可以认为，计算出的$X$ 将是某个矩阵$\left(A+F\right)$ 的精确逆，其中$F$ 的元素为$\left(1/n\right)f(n)2^{-t}$ ，但这很少是正确的。实际上，$X$ 的逆矩阵元素通常与$A$ 的逆矩阵元素有$\left(1/n\right)f(n)2^{-t}\Vert A^{-1}\Vert$ 的数量级不同。因为如果$A$ 是病态的，$\Vert A^{-1}\Vert$ 将非常大，所以最后的结果有点令人惊讶。然而，一般来说，$X$ 的第$r$ 列确实是$\left(A+E_r\right)$ 的精确逆的第$r$ 列，其中$E$ 具有阶数为$\left(1/n\right)f(n)2^{-t}$ 的元素，但它要求每一列都有不同的$E_r$ 。对于一个病态的$2\times 2$ 矩阵，读者可以很容易地验证这些注释。它们对于对称矩阵和非对称矩阵同样成立。

计算逆的高精度

27. 通常的经验是，即使通过统计效应对误差进行了最慷慨的估计，计算出来的逆函数也比预期的要准确得多。对这种现象不太可能有一个简单的解释，但下面这个影响非常明显的例子确实对一些主要原因提供了相当多的见解。我们考虑希尔伯特矩阵五阶对称部分的求逆运算，涉及8个小数。矩阵是正定的，$\Vert A^{-1}\Vert$ 大约为$10^6$ 。它的条件非常糟糕，将原始矩阵的元素截断为 8 位小数已经影响了逆矩阵的第二个数字。 这是意料之中的，因为$\Vert A^{-1}\Vert \Vert E\Vert \Vert A^{-1}\Vert$ 是有序的

自然，我们将计算得到的逆矩阵与截断矩阵的逆矩阵进行比较，而不是与精确的希尔伯特段的逆矩阵进行比较!我们用$H$ 表示确切的希尔伯特段，用$\bar{H}$ 表示截断元素的矩阵，用$X$ 表示计算出的$\bar{H}$ 的逆矩阵。现在我们从表3中看到

另一方面，对于通过第19节中的技术计算的$X$ ，我们可以看到

在求解过程中所产生的舍入误差的总影响远小于初始截断误差的影响。

有两个主要的现象可以显著说明解决方案的显著准确性。

(i)我们在第18节中看到$L{L}^T=\bar{H}+E$ , 并且$E$ 的典型元素以$2^{-t-1}|l_{ii}$ 为界。如果某些元素$|l_{ii}|$ 远小于1，则$E$ 相应元素的边界远小于$2^{-t}$ 。现在，因为$\bar{H}$ 是病态的，所以$l_{ii}$ 很小，并且随着$i$ 的增加而逐渐变小。图中显示了矩阵$E$ ，我们可以看到，当我们移动到右下角时，它的元素也逐渐变小。这里它们远小于$\frac{1}{2}\cdot 10^{-8}$ (对于8进制计算，我们可能期望的最佳值)。请注意，$L{L}^T$ 如此接近$\bar{H}$ ，正是因为$\bar{H}$ 是病态的!