Q^* T Q的稳定性

算法4.9所示的收缩正交变换构造显然是稳定的（即，它在各分量上相对精确地表示了在精确算术中得到的变换）。毫无疑问，相似变换Q^* T Q能够完全保留矩阵T的特征值的数值精度。然而，在锁定和/或清除过程中，三对角形式的保持程度存在一个严重的问题。为了确保如果T y = y \theta，则T_+ \equiv Q^* T Q是对称的且在数值上是三对角的（即，非三对角带外的所有元素相对于\Vert T \Vert都是微小的），需要对基本算法进行修改。

T_+ = Q^* T Q为三对角形式依赖于表达式中的项g y^* T R消失：

Q^* T Q = L^*TR + g y^*T R + \theta e_1 e_1^*.

然而，进一步检查发现：

e_1^* Q = e_1^* L + (e_1^* y) g^* = \eta_1 g^*,

其中\eta_1是y的第一个分量。因此，\Vert g \Vert = \frac{1}{\vert \eta_1 \vert}。当y的第一个分量很小时，可能会出现数值困难。具体来说，y^* T = \theta y^*，因此在精确算术中 y^* T R = 0。然而，在有限精度下，计算出的\mathrm{fl}(y^* T) = \theta y^* + z^*。误差z相对于\Vert T \Vert将在\epsilon_M的量级上，但：

\Vert g y^*T R \Vert = \Vert g\Vert \cdot \Vert z^* R \Vert =\frac{1}{\vert \eta_1 \vert} \Vert z^* R \Vert,

因此，这一项可能相当大，如果\eta_1 = O(\epsilon_M)，它可能达到O(1)的量级。这是一个严重的问题，并且在实践中如果不进行我们提出的修改，这种情况就会发生。

解决办法是引入逐步可接受的扰动和对向量y的重缩放，以同时强制满足以下条件：

Q^* y = e_1 \ \ \mathrm{and}\ \ y^* T Q = \theta e_1^*

在有限精度下具有足够的准确性。为了实现这一点，我们将设计一个方案，以实现：

y^* T q_j = 0 \ \ \mathrm{for} \ \ j > 1

在数值上。如[420]所示，这一修改足以在数值上建立q_i^* T q_j = 0，相对于\Vert T \Vert，对于i > j + 1。

基本思想如下：如果在第j步，计算出的量\mathrm{fl}(y^* T q_j)不够小，则通过用一个数\phi缩放向量y_j和用一个数\psi缩放分量\eta_{j+1}来调整它，使其足够小。在进行q_{j+1}的计算之前，进行这种重缩放，我们有y_j \leftarrow y_j \phi和\hat{y}_j \leftarrow \hat{y}_j \psi，其中y^* = [ y_j^* , \hat{y}_j^*]。当然，\Vert y \Vert不应因这种缩放而改变，因此这是必需的。这给出了以下方程组来确定\phi和\psi：如果\vert y_j^* T_j \hat{q}_j + \rho_j \beta_j \eta_{j+1} \vert > \epsilon_M \tau_{j+1}，

\begin{aligned} (\tau_j \phi)^2 + (1 - \tau_j^2)\phi^2 &= 1, \\ y_j^* T_j \hat{q}_j\phi + \rho_j \beta_j \eta_{j+1}\psi &= \pm \epsilon_M \tau_{j+1}. \end{aligned}

如果\rho_j在\sqrt{\epsilon_M}的量级上，这种缩放可以被吸收到\rho_j中，而不改变y，也不影响Q的数值正交性。当y被修改时，之前计算的q_i, 2 \le i \lt j,不需要改变。在第j步之后，向量y_j在后续步骤中仅被重新缩放，定义q_i, 2 \le i \le j,的公式对y_j的缩放是不变的。有关详细信息，请参阅[420]。

算法4.10实现了这一方案来计算可接受的Q ¹。在实践中，这种变换可以在不存储Q的情况下就地计算和应用，以获得T \leftarrow Q^* T Q和V \leftarrow VQ。然而，这种实现相当微妙，构造Q的细节要求避免额外的存储。该代码与上述描述略有不同，因为向量y仅在所有列2到m确定后才被重新缩放到单位范数。

有几个实现问题：

(3)

这里显示的扰动y(i)避免了特征向量y中精确零初始项的问题。理论上，当T未简化时，这种情况不应该发生，但如果\theta在初始向量中弱表示，这种情况可能在数值上发生。有一种更简洁的实现方法，不修改零项。这是最简单（也是最粗糙）的修正。

(10)

一旦\tau_j = \Vert y_j\Vert足够大，就不需要进一步的修正。删除这个if-clause将恢复到算法4.9所示的未修改的Q计算。

(16)

这显示了几种可能的修改y的方法之一，以实现Q^* T Q的非三对角带外元素在数值上微小的目标。更复杂的策略会将尽可能多的缩放吸收到对角元素Q(j,j) = \rho中。这里，不缩放\rho而是缩放y的分支设计为使{\rm fl}(1 + (\rho \psi)^2) = {\rm fl}(1 + \rho^2) = 1。

1. 目前，有一种更为优秀的算法用于构建和应用这些变换；参见 D. Sorensen, Deflation for Implicitly Restarted Arnoldi Methods, CAAM Technical Report, TR 98-12, Rice University, 1998 (revised August 2000).