算法

在直接迭代变体中， 我们用A进行乘法运算并解B的线性系统； 这对应于C=B^{-1}A。它计算一个 基V_j，其中矩阵束\{A,B\}由一个 实对称三对角矩阵表示，

T_{j}=\begin{bmatrix} \alpha_1 & \beta_1 & & \\ \beta_1 & \alpha_2 & \ddots & \\ & \ddots & \ddots & \beta_{j-1}\\ & & \beta_{j-1} & \alpha_j \end{bmatrix}, \tag{5.8}

满足基本递归关系，

AV_j=BV_{j}T_{j}+re_j^{\ast}. \tag{5.9}

与标准厄米特情况（4.10）相比较。 现在基V_j是B-正交的，

V_j^{\ast}BV_j=I_{j}, \tag{5.10}

且矩阵T与A的一个截面是相合的，

V_j^{\ast}AV_j=T_{j}. \tag{5.11}

我们通过引入一个辅助基来简化B-正交化的描述，

W_j=BV_j, \tag{5.12}

它是B^{-1}-正交的，

W_j^{\ast}B^{-1}W_j=I_{j}, \tag{5.13}

并且对于它，W^{\ast}V=I。

正如标准情况一样， 计算T_{j}的特征解，

T_{j}s_i^{(j)}=s_i^{(j)}\theta_i^{(j)},

并得到一个 Ritz 值 \theta_i^{(j)} 和一个 Ritz 向量，

x_i^{(j)}=V_js_i^{(j)}, \tag{5.14}

近似于束的特征对（5.1）。 其残差，

\begin{aligned} r_i^{(j)} &= Ax_i^{(j)}-Bx_i^{(j)}\theta_i^{(j)}=AV_js_i^{(j)}- BV_js_i^{(j)}\theta_i^{(j)} \\ &= (AV_j-BV_jT_j)s_i^{(j)} = r e^*_j s_i^{(j)} = Bv_{j+1}\beta_js_{j,i}^{(j)}, \end{aligned}

是B-正交于由V_j张成的Krylov空间。

我们可以像在标准厄米特情况中那样估计残差的范数（4.13），但现在使用B^{-1}-范数会更好，得到

\Vert r_i^{(j)}\Vert _{B^{-1}}=\vert\beta_js_{j,i}^{(j)}\vert\;, \tag{79}

利用了事实 (Bv_{j+1})^{\ast}B^{-1}(Bv_{j+1})=v_{j+1}^{\ast}Bv_{j+1}。 在衡量收敛性时自然使用B^{-1}-范数； 参见[353, Chap. 15]。

与标准情况一样，我们需要 监控T的次对角元素\beta_j， 以及其特征向量的最后一个元素s_{j,i}^{(j)}。一旦这个 乘积变小，我们就可以标记一个特征值为收敛， 而不必实际执行矩阵-向量 乘法（5.14）。我们保存这个 操作直到步骤j，那时估计（5.15） 指示收敛。

我们得到以下算法。

算法 5.4 GHEP的Lanczos方法
(1)  取初始向量 q，计算 r = Bq,\; \beta_0 = \vert q^*r\vert^{1/2}
(2)  for j = 1,2,\dots,
(3)      w_j = r/\beta_{j-1}, v_j = q/\beta_{j-1}
(4)      r = Av
(5)      r = r - w_{j-1}\beta_{j-1}
(6)      \alpha_j = v_j^*r
(7)      r = r - w_j \alpha_j
(8)      如果需要，重正交化
(9)      解 Bq = r 得到 q
(10)     \beta_j = \vert q^* r\vert^{1/2}
(11)     计算 T_j 的特征分解 T_j = S \Theta^{(j)} S^*
(12)     如果收敛，停止
(13) end for
(14) 计算特征向量 x_i^{(j)} = V_js_i^{(j)}，对于每个收敛的 i

让我们一步一步地评论这个算法：

(1)

如果有一个好的想要的特征向量的猜测，就用它作为 起始q。在其他情况下选择一个随机方向，例如， 由正态分布的随机数组成。注意 \alpha_1=q^{\ast}Aq/(q^{\ast}Bq)是起始向量的Rayleigh商，而\beta_1测量其 残差的B^{-1}-范数（5.15）。

(4)

这是大矩阵A出现的地方。 任何执行矩阵-向量乘法的例程都可以使用。

(8)

只需重新正交化基V或W中的一个， 而不是两者。 选择与标准情况相同：

1. 完全： 现在我们希望使v基向量B-正交， 计算
   r=r-B(V_j(V_j^{\ast}r))
   并重复直到向量r与基V_j正交。 我们必须为每次重新正交化应用一次B的矩阵-向量乘法，并且必须使用经典 的Gram-Schmidt过程。

   如果我们保存两个基V_j和W_j并减去W_j的列的倍数，

   r=r-W_j(V_j^{\ast}r),
   直到正交性获得，几乎总是只需一次。 现在我们可以使用改进的Gram-Schmidt正交化。

2. 选择性：仅在必要时重新正交化， 必须监控正交性，如 §4.4.4中所述 对于标准情况。注意现在对称矩阵 V_j^{\ast}W_j=V_j^{\ast}BV_j涉及其中，一些向量v_k 必须由相应的向量w_k替换。
3. 局部：用于巨大的矩阵，当存储 整个基V_j很困难时。仅在寻找一个或两个极端特征值时才建议使用。我们确保w_{j+1}与v_{j-1}和 v_j正交，通过减去 r=r-w_{j-1}(v_{j-1}^{\ast}r); r=r-w_j(v_j^{\ast}r)一次 在这一步。

(9)

在这里我们需要解一个正定矩阵B的系统。这在标准情况中 （4.1）是不需要的。

(11)

对于每个步骤j，或在适当的间隔， 计算对称三对角矩阵 T_{j}（5.8）的特征值 \theta_i^{(j)}和特征向量s_i^{(j)}。与标准情况相同的过程。

(12)

当找到一个足够大的基V_j时，算法停止， 使得三对角矩阵T_{j}（5.8）的特征值 \theta_i^{(j)}给出束（5.1）的所有 特征值的良好近似。

如果基 V_j不是完全B-正交的，估计（5.15）对于残差 可能过于乐观。 那么Ritz向量x_i^{(j)}（5.14）可能具有范数 小于1，我们必须用以下估计替换，

\Vert r_i^{(j)}\Vert=\vert\beta_js_{j,i}^{(j)}\vert/\Vert V_js_i^{(j)}\Vert _B\;.

(14)

只有在步骤（5.4）中的测试指示它们 已经收敛时，才计算原始矩阵束（5.1）的特征向量。然后基V_j用于矩阵-向量 乘法以得到 特征向量（5.14），

x_i^{(j)}=V_js_i^{(j)},

对于每个标记为收敛的i。