求解稀疏且高度病态的系统

计算科学 线性求解器 稀疏矩阵 拉帕克 条件数
2021-11-28 06:50:42

我打算解决 Ax = b 其中 A 是复数、稀疏、不对称和高度病态(条件数 ~ 1E+20)方阵或矩形矩阵。我已经能够在 LAPACK 中准确地使用 ZGELSS 解决系统问题。但是随着我系统中自由度的增加,在使用 ZGELSS 的 PC 上解决系统需要很长时间,因为没有利用稀疏性。最近我为同一系统尝试了 SuperLU(使用 Harwell-Boeing 存储),但结果对于条件数 > 1E+12 不准确(我不确定这是否是旋转的数值问题)。

我更倾向于使用已经开发的求解器。是否有一个强大的求解器可以快速(即利用稀疏性)和可靠地(鉴于条件数)解决我提到的系统?

3个回答

当你使用 ZGELSS 来解决这个问题时,你正在使用截断奇异值分解来规范这个极端病态的问题。的最小二乘解决方案,而是试图平衡找到一个最小化反对最小化. Ax=bxAxb

请注意,传递给 ZGELSS 的参数 RCOND 可用于指定应包含哪些奇异值,哪些奇异值应从解决方案的计算中排除。任何小于 RCOND*S(1) 的奇异值(S(1) 是最大的奇异值)都将被忽略。矩阵或右侧中的系数的噪声水平一无所知,所以很难说您是否使用过适量的正则化。 Ab

您似乎对使用 ZGELSS 获得的正则化解决方案感到满意,因此似乎正则化受截断 SVD 方法影响(在最小化 \| 的最小二乘解决方案中找到最小 |在与大于 RCOND*S(1)) 的奇异值相关联的奇异向量所跨越的解空间上对您来说是令人满意的。 xAxb

您的问题可以重新表述为“我怎样才能有效地获得这个大型、稀疏和非常病态的线性最小二乘问题的正则化最小二乘解?”

我的建议是使用迭代方法(例如 CGLS 或 LSQR)来最小化显式正则化最小二乘问题

minAxb2+α2x2

其中正则化参数被调整,以使阻尼最小二乘问题得到很好的调节,并且您对得到的正则化解决方案感到满意。 α

Jed Brown 已经在对问题的评论中指出了这一点,但是如果您的条件数很大,您在通常的双精度中实际上并不能做很多事情:在大多数情况下,您可能无法获得一位数的准确度您的解决方案,更糟糕的是,您甚至无法分辨,因为您无法准确评估与您的解决方案向量相对应的残差。换句话说:这不是您应该选择哪种线性求解器的问题——没有线性求解器可以为此类矩阵做一些有用的事情。

发生这种情况通常是因为您选择了不合适的基础。例如,如果您选择函数作为 Galerkin 方法的基础,则会得到此类病态矩阵。(这导致了 Hilbert 矩阵,众所周知,它的条件很差。)这种情况下的解决方案不是问哪个求解器可以求解线性系统,而是问是否有更好的基可以使用。我鼓励你也这样做:考虑重新表述你的问题,这样你就不会得到这些类型的矩阵。1,x,x2,x3,...

解决病态问题的最简单/最快的方法是提高计算的精度(通过蛮力)。另一种(但并非总是可能的)方法是重新制定您的问题。

您可能需要使用四倍精度(34 个十进制数字)。即使在课程中会丢失 20 位数字(由于条件编号),您仍然会得到 14 位正确的数字。

如果有任何兴趣,现在MATLAB 中也提供了四精度稀疏求解器。

(我是上述工具箱的作者)。