如何在 scipy 中实现 ILU 预调节器？

Question

对于 scipy.sparse.linalg 中的迭代求解器，例如 bicg、gmres 等，可以选择为矩阵 A 添加预条件子。但是，documentation 并不清楚我应该给什么作为预处理器。如果我使用 ilu = sp.sparse.linalg.spilu(A)，ilu 不是任何矩阵，而是包含许多东西的对象。

有人为 Python 2.7 问了一个类似的问题 here，但我不适合我（Python 3.7，scipy 版本 1.1.0）

所以我的问题是如何为那些迭代算法合并不完整的 LU 预处理器？

Answer 1

作为预处理器，bicg or gmres接受

稀疏矩阵
密集矩阵
线性运算符

在您的情况下，预条件子来自因式分解，因此它必须作为线性运算符传递。

因此，您可能希望根据通过 spilu 获得的 ILU 分解显式定义线性运算符。沿线的东西：

sA_iLU = sparse.linalg.spilu(sA)
M = sparse.linalg.LinearOperator((nrows,ncols), sA_iLU.solve)

此处，sA 是 CSC 格式的稀疏矩阵，M 现在将成为您将提供给迭代求解器的预条件子线性运算符。

一个完整的例子based on the question you mentioned:

import numpy as np
from scipy import sparse
from scipy.sparse import linalg

A = np.array([[ 0.4445,  0.4444, -0.2222],
              [ 0.4444,  0.4445, -0.2222],
              [-0.2222, -0.2222,  0.1112]])
sA = sparse.csc_matrix(A)

b = np.array([[ 0.6667], 
              [ 0.6667], 
              [-0.3332]])

sA_iLU = sparse.linalg.spilu(sA)
M = sparse.linalg.LinearOperator((3,3), sA_iLU.solve)

x = sparse.linalg.gmres(A,b,M=M)

print(x)

备注：

我实际上使用的是密集矩阵作为示例，而在您的情况下从具有代表性的稀疏矩阵开始更有意义。
线性运算符 M 的大小是硬编码的。
ILU 无论如何都没有配置，但使用默认值。
这几乎是所建议的 in the comments to that aforementioned answer，但是，我不得不做一些小的调整以使其 Python3 兼容。

如何在 scipy 中实现 ILU 预调节器？

How to implement ILU precondioner in scipy?

python

numpy

linear-algebra

scipy

sparse-matrix