我目前正在使用 P1 FEM 离散化求解谐波方程。生成的矩阵是 SPD 且相当稀疏，因此我使用预条件共轭梯度 (CG) 求解器来找到解决方案。我试图通过 Jacobi 预条件子和不完整的 Cholesky 预条件子（使用 : IC(0) 的稀疏模式）来加速收敛。两者都非常适合大幅减少 CG 迭代的次数，但是即使与 Jacobi 预处理迭代 (90) 相比，Choleksy 预处理迭代 (30) 少 3 倍，但它仍然慢了一点。我已将其缩小到求解三角系统，这似乎解释了每次迭代的运行时间增加了 3 倍（我假设这是因为它们本质上是编译器无法并行的）。我的问题是是否有办法加快速度。我在 SSOR 的上下文中阅读了 Eisenstat 技巧，在 Saad 的书中提到了以下内容：
此外，对于某些类型的矩阵，IC(0) 预条件子可以用这种方式表示，其中 D 可以通过递推公式获得。
虽然我还没有找到更多的细节。根据我的分析，不完全 Cholesky 分解的计算在运行时可以忽略不计。我可以做任何相关的优化吗？用填充来研究不完整的 Cholesky 分解是一个好主意吗？或者也许专注于代数多重网格预处理器会更好？
编辑：我在类似矩阵的对角线上添加了一张放大图像，类似于我所拥有的：
原始系统是 1600 万乘 1600 万元素，所以我无法存储它的完整图像，但我已经测试了几个较小的变体和所有非零条目都聚集在对角线周围，这是可以预料的，因为 FEM 网格中的顶点是从左到右和从上到下（空间上）枚举的，所以至少有一些排序，尽管我怀疑它是最优的（Dirichlet 节点也被取出，因为它们只影响右侧）。
编辑 2：这似乎是一个已知的瓶颈，但它的解决方案绝不简单：https : //arxiv.org/pdf/1908.00741.pdf 在尝试 RCM 如何重新排序之后，我可能会改为研究 multigrid会影响运行时间。
编辑 3：我实施了 Neil Lindquist 的建议，但是这些建议并没有在运行时产生显着的改进，但它们确实大大减少了内存需求和内存访问。我假设对于我的问题，我受到计算而不是内存访问的限制，主要是三角求解缺乏并行化。为此，我想我应该参考链接的论文，或者改用域分解或多重网格预处理器（因为两者都可以并行化）。