用于查找大型但结构化矩阵的最小特征对的预处理器

计算科学 宠物 线性求解器 稀疏矩阵 本征系统 预处理
2021-12-08 18:42:23

我试图找到与大矩阵的第二小的特征值相对应的特征向量。是图拉普拉斯算子,具有以下结构:,其中D是对角矩阵,SB是大型稀疏矩阵。S的维度为4,000,000 \times 10,000,000,但只有大约40,000,000个非零条目。B具有相似的规模和稀疏性。所以我可以快速执行矩阵向量乘法:Lv = Dv - S(B(S^\intercal v))(4,000,000×4,000,000)LLL=DSBSDSBS4,000,000×10,000,00040,000,000BLv=DvS(B(Sv))

目前我正在使用 Scipy(它调用在 ARPACK 中实现的 Arnoldi 算法)来查找L的最小特征值和相应的特征向量。我不是直接找到L的最小特征对,而是找到M^{-1}L的最大特征对,其中M = L + cI(添加cI和反转不会改变特征向量。)为了计算M^{-1}v,其中v是任意向量,我使用共轭梯度算法。M1M=L+cIcIM1vv

但是这种方法很慢:Arnoldi 算法通常需要数百次迭代才能收敛,而每次 Arnoldi 迭代需要数百次共轭梯度 (CG) 迭代。由于我正在使用相同的设计矩阵M解决数百个 CG 问题,因此我认为预处理可以提高效率。但是什么预处理器适合这个问题呢?我已经读到存在可证明的良好的拉普拉斯算子稀疏近似(“超解析器”),但是可以在不破坏其稀疏性的情况下对它们进行移位和反转,以形成M1的良好近似吗?或者是否有利用M结构的预处理技术?

非常感谢。

3个回答

如果您只对最小特征值感兴趣,则应用于矩阵的共轭梯度法在相当少的步数后为您提供了一个很好的近似值,并且您不必求解任何线性系统。详细信息在 Y. Saad 关于迭代方法的书中,但这里有一个简短的总结:L

  1. 根据在共轭梯度 (cg) 迭代中计算的系数,您可以构造cg 步数维数的这个矩阵是矩阵到当前 Krylov 空间的 Galerkin 投影。TL

  2. 当 cg 步数增加时,的两个极值特征值接近的特征值。TL

  3. 由于很小且是三对角线,因此计算其特征值相当便宜,例如使用 LAPACK。T

  4. 一旦你有了特征值,计算相应的特征向量就相当便宜了。

作为一个侧节点:deal.II中实现的 cg 求解器有一个选项,允许您在每一步中计算的谱,从而观察特征值的收敛。T

假设您的整个矩阵是正定的(它绝对是对称的),那么我建议将代数多重网格(AMG)方法作为预处理器。他们自己计算稀疏矩阵的层次结构。如果您已经在使用 PETSc,请查看 hyprepreconditioner。使用它可能会迫使您实际乘以矩阵的元素,但这可能是值得的。(hypre 也有可能接受未组装的精细矩阵。)

您绝对应该尝试 Jacobi-Davidson 方法。这种方法的一个优点是不需要只需要一个近似的Jacobi-Davidson 方法的收敛性可以通过使用预处理器来改进。如果可以显式计算,那么您可以尝试对该矩阵进行“0”级不完全 Cholesky 分解。可能你需要一个小的转变来避免不完整的分解在零或负枢轴上崩溃。您也可以使用没有前置条件的 Jaobi-Davidson,但在这种情况下收敛可能不太好。存在针对对称矩阵优化的 Jacobi-Davidson 的 matlab 代码:JDCGM1vM1vL=DSBST