C ++犰狳没有正确解决条件差的矩阵

Question

关于犰狳内置的线性求解器，我有一个相对简单的问题。我是 C++ 的新手，但有其他语言的编码经验。我正在通过连续线性化解决流体流动问题，使用犰狳函数 Solve(A,b) 在每次迭代时获得解决方案。

我遇到的问题是我的矩阵非常病态。行列式大约是 10^-20，条件数是 75000。我知道这些条件很糟糕，但这就是我所拥有的。有谁知道是否可以将我的 A 矩阵和求解函数中的精度指定为双精度（可能是长双精度）？我知道犰狳中有双矩阵类，但我还没有找到任何更高精度的文档。

为了从另一个角度解决这个问题，我在 Mathematica 中编写了一些代码，LinearSolve 运行良好，程序收敛到正确答案。我的理由是 Mathematica 变量具有更高的精度，可以处理更高级别的舍入误差。

如果有人对此有任何见解，请告诉我。我知道还有其他方法可以处理条件不佳的矩阵（例如预处理和旋转），但我的工作更多是在物理方面，而不是在实际的数值解决方案中，所以我试图避免这种情况。

编辑：我只是将 Mathematica 版本的精度限制为小数点后 15 位，程序仍然收敛。这让我相信这不是一个精度可变的问题，而是方法的问题。

Answer 1

正如您所说的“您的工作更多是在物理学中”：我不会尝试提高准确性，而是使用 Moore-Penrose Pseudo-Inverse，它在 Armadillo 中可以通过函数pinv 获得。然后，您应该对参数tolerance 进行一些体验，以将其设置为合理的水平。

几何解释如下：不良条件数是由于行/列向量线性相关的事实。在物理学中，这种线性依赖通常有一个至少需要解释的起源。伪逆首先将矩阵投影到低维空间，在该空间中，通过删除所有奇异值小于参数tolerance 的奇异向量，向量“线性相关性较低”。结果矩阵有一个更好的条件数，这样标准逆矩阵的构造问题就更少了。