基础算法

Jacobi-Davidson算法的基本形式在算法4.13中给出。稍后，我们将描述带有重启和其他策略的更复杂的变体。

在该算法的每次迭代中，计算厄米矩阵A的特征对(\theta, {u})，对应于A的最大特征值。迭代过程在残差的范数 A{u}-\theta {u}低于给定阈值\epsilon时终止。

要应用此算法，我们需要指定一个起始向量{v}_0和一个容差\epsilon。完成后，将得到最大特征值 \lambda=\lambda_{\max}(A)及其对应的特征向量 {x}={x}_{\max}的近似值。计算的特征对 (\widetilde\lambda, \widetilde{x}) 满足 \Vert A\widetilde{x} - \widetilde\lambda\widetilde{x} \Vert\leq \epsilon。

现在，我们将描述一些实现细节，参考算法4.13中的相应阶段。

(1)

这是过程的初始化阶段。搜索子空间在每次迭代中通过向量{t}扩展，我们从给定的向量t=v_0开始这一过程。理想情况下，该向量应在所需特征向量的方向上有显著分量。除非对所需特征向量有所了解，否则从随机向量开始可能是个好主意。这使得我们有信心认为所需特征向量在起始向量中具有非零分量，这对检测特征向量是必要的。

(3)-(5)

这代表了新向量t相对于集合{v}_1,\ldots, {v}_{{m}-1}的修正Gram-Schmidt过程。如果{m}=1，这是一个空循环。设{t}_{in}表示正交化开始前的向量，{t}_{out}表示阶段(3)-(5)完成后的向量。建议（参见[96]）如果 \Vert{t}_{out}\Vert _2/\Vert{t}_{in}\Vert _2\leq \kappa，其中\kappa是一个适度的常数，例如\kappa=.25，则重复Gram-Schmidt过程一次。这保证了在相对意义上，正交性的损失被限制在1/\kappa倍机器精度内。修正Gram-Schmidt正交化与迭代改进的模板在算法4.14中给出。

(7)-(9)

在这一阶段，计算矩阵{M}\equiv {V}^\ast A{V}的上三角部分的第{m}列。矩阵{V}表示列向量为{v}_j的n乘{m}矩阵，{V^A}同样。

(10)

计算{m}\times {m}厄米矩阵{M}的最大特征对，其上三角部分包含元素{M}_{i,j}，可以使用LAPACK中的适当例程完成（参见第4.2节）。

(12)

向量{u^A}可以按此处描述进行更新，也可以重新计算为{u^A}=A{u}，取决于哪种成本更低。选择在于{m}次更新和另一次与A的乘法；如果A平均每行少于{m}个非零元素，则通过A{u}计算更可取。如果{u^A}计算为A{u}，则无需存储向量v^A_j。

(14)

如果 \Vert A{u}-\theta {u}\Vert _2\leq \epsilon，则算法终止。在这种情况下，A有一个特征值\lambda，满足 \vert\lambda -\theta\vert\leq \epsilon。对于相应的归一化特征向量，如果\lambda是简单的且与其他特征值分离良好，则角度的界限相似。这也导致\lambda的更严格的界限；参见（4.5）或第4.8节。

收敛到 \lambda\neq\lambda_{\max}(A) 可能发生，但通常不太可能。例如，如果 v_0\perp {x}_{\max}，或者 选定的\theta非常接近另一个特征值 \lambda\neq\lambda_{\max}(A)，这可能发生。对于任何迭代求解器，特别是如果 \epsilon取值不够小（例如，大于机器精度的平方根），这种情况可能发生。

(15)

扩展向量{t}的近似解可以通过Krylov求解器计算，例如MINRES或SYMMLQ。使用左或右预条件化时，必须选择适用于非对称系统的Krylov求解器（如GMRES、CGS或Bi-CGSTAB），因为预条件化后的算子通常不对称。带有左预条件化的Krylov子空间方法的近似解模板在算法4.15中给出。右预条件化的情况，稍微更昂贵，在算法4.16中给出模板。关于迭代Krylov子空间求解器，参见[41]。近似解必须与{u}正交，但如果从初始猜测{t}_0=0开始，使用Krylov求解器时自动满足这一条件。在大多数情况下，不需要将校正方程求解到高精度；在第{m}次迭代中，相对精度为2^{-m}似乎足够。建议对迭代求解器的迭代步数设置限制。

Davidson[99]建议取 {t}=(\mathrm{\rm diag}(A)-\theta I)^{-1}r，但在此情况下{t}不与{u}正交。此外，对于对角矩阵，这种选择会导致停滞，这说明了这种方法的问题。

为了限制存储，算法可以在某个适当的值 {m}={m}_{\max} 处终止，并以{u}的最新值作为v_0重新启动。我们将在第4.7.3节中描述一个具有更复杂重启策略的Jacobi-Davidson算法变体。

请注意，大多数计算密集型操作，即成本与n成正比的操作，可以轻松并行化。此外，多向量更新可以通过适当的Level 2 BLAS例程执行（参见第10.2节）。

在接下来的小节中，我们将描述Jacobi-Davidson算法的更复杂变体。在第4.7.3节中，我们将介绍一个允许重启的变体，如果希望保持所涉及子空间的维度有限，这将很方便。该变体也适用于在发现特征对后进行重启，以定位下一个特征对。该技术基于收缩。生成的算法设计用于计算给定矩阵的几个最大或最小特征值。

在第4.7.4节中，我们将描述一个适用于计算A内部特征值的Jacobi-Davidson方法变体。