Scipy稀疏反转或溶解导致UMFPACK_ERROR_OUT_OF_MEMORY

Question

我正在尝试按如下方式转换大的(150000,150000)稀疏矩阵：

import scipy as sp
import scipy.sparse.linalg as splu

#Bs is a large sparse matrix with shape=(150000,150000)

#calculating the sparse inverse
iBs=splu.inv(Bs)

导致以下错误消息：

Traceback (most recent call last):
    iBs=splu.inv(Bs)
  File "/usr/lib/python2.7/dist-packages/scipy/sparse/linalg/dsolve/linsolve.py", line 134, in spsolve
autoTranspose=True)
  File "/usr/lib/python2.7/dist-packages/scipy/sparse/linalg/dsolve/umfpack/umfpack.py", line 603, in linsolve
self.numeric(mtx)
  File "/usr/lib/python2.7/dist-packages/scipy/sparse/linalg/dsolve/umfpack/umfpack.py", line 450, in numeric
umfStatus[status]))
RuntimeError: <function umfpack_di_numeric at 0x7f2c76b1d320> failed with UMFPACK_ERROR_out_of_memory

我重新调整程序以简单地求解线性微分方程组：

import numpy as np

N=Bs.shape[0]

I=np.ones(N)

M=splu.spsolve(Bs,I)

并且我再次遇到相同的错误

我在具有16 GB RAM的计算机上使用此代码，然后将其移至具有32 GB RAM的服务器上，仍然无济于事。

以前有没有人遇到过？

Answer 1

[首先让我说，应该在http://scicomp.stackexchange.com上最好提出这个问题，因为那里有大量的计算科学和数值线性代数专家。

让我们从基础开始：从不反转稀疏矩阵，这完全没有意义。请参见MATLAB Central上的discussion，尤其是Tim Davis的comment。

简而言之：没有用于对矩阵进行数值反转的算法。每当您尝试数值计算NxN矩阵的逆时，实际上就可以解决N个线性系统，其中N个rhs向量对应于单位矩阵的列。]

换句话说，当您计算时

from scipy.sparse import eye
from scipy.sparse.linalg import (inv, spsolve)

N = Bs.shape[0]
iBs = inv(Bs)
iBs = spsolve(Bs, eye(N))

最后两个语句（inv(eye)和spsolve(Bs, eye(N))）是等效的。请注意，当您提出错误的假设时，单位矩阵（eye(N)）不是一个向量（np.ones(N)）。

这里的意思是，矩阵逆在数值线性代数中很少有用：Ax = b的解不是通过inv（A）* b计算的，而是通过专门的算法计算的。

转到您的特定问题，对于大型稀疏方程组，没有black-box求解器。只有对矩阵问题的结构和性质有充分的了解，才可以选择正确的求解器类别。矩阵的属性又是您要解决的问题的结果。例如。当通过FEM离散椭圆PDE系统时，最终会得到对称的正稀疏代数方程组。一旦知道了问题的性质，就可以选择正确的解决方案。

在您的情况下，您尝试使用通用直接求解器，而不对方程重新排序。众所周知，这会生成填充，从而破坏iBs函数第一阶段中spsolve矩阵的稀疏性（应该将其分解）。请注意，全双精度150000 x 150000矩阵需要大约167 GB的内存。为了减少因式分解过程中的填充，有很多技术可以对方程进行重新排序，但是您提供的信息不足以提供明智的提示。

[很抱歉，您应该考虑在http://scicomp.stackexchange.com上重新提出问题，清楚说明您要解决的问题，以便对矩阵的结构和性质有所了解。

Answer 2

稀疏数组仅将矩阵的非零条目放入内存。现在假设您进行了反转。这意味着矩阵的几乎所有条目都变为非零。稀疏矩阵经过内存优化。

您可以在稀疏矩阵上应用一些操作，而不会丢失“ spare”属性：

• 加法，仅添加一个常数就可以使稀疏矩阵保持稀疏。