Spa*_*r23 14 python octave scipy sparse-matrix eigenvector
我正在尝试计算几个 (5-500) 特征向量，这些特征向量对应于大型对称方形稀疏矩阵（高达 30000x30000）的最小特征值，其中非零值小于 0.1%。
我目前在 shift-invert 模式 (sigma=0.0) 下使用 scipy.sparse.linalg.eigsh，我通过有关该主题的各种帖子发现这是首选解决方案。但是，在大多数情况下，解决问题最多需要 1 小时。另一方面，如果我要求最大的特征值（在我的系统上为亚秒），则该函数非常快，这是文档中预期的。
由于我在工作中更熟悉 Matlab，所以我尝试在 Octave 中解决这个问题，使用 eigs (sigma=0) 在短短几秒钟内（低于 10 秒）给了我相同的结果。由于我想对包括特征向量计算在内的算法进行参数扫描，因此在 python 中也能获得这种时间增益。
我首先更改了参数（尤其是容差），但在时间尺度上并没有太大变化。
我在 Windows 上使用 Anaconda，但试图将 scipy 使用的 LAPACK/BLAS（这是一个巨大的痛苦）从 mkl（默认 Anaconda）切换到 OpenBlas（根据文档由 Octave 使用），但看不到变化表现。
我无法弄清楚，使用的 ARPACK 是否有改变（以及如何改变）？
我将以下代码的测试用例上传到以下 dropbox 文件夹：https ://www.dropbox.com/sh/l6aa6izufzyzqr3/AABqij95hZOvRpnnjRaETQmka?dl =0
在 Python 中
import numpy as np
from scipy.sparse import csr_matrix, csc_matrix, linalg, load_npz
M = load_npz('M.npz')
evals, evecs = linalg.eigsh(M,k=6,sigma=0.0)
在八度：
M=dlmread('M.txt');
M=spconvert(M);
[evecs,evals] = eigs(M,6,0);
任何帮助是appriciated ！
我根据评论和建议尝试了一些其他选项：
Octave：
eigs(M,6,0)并eigs(M,6,'sm')给我相同的结果：
[1.8725e-05 1.0189e-05 7.5622e-06 7.5420e-07 -1.2239e-18 -2.5674e-16]
而eigs(M,6,'sa',struct('tol',2))收敛到
[1.0423 2.7604 6.1548 11.1310 18.0207 25.3933]
快得多，但前提是容差值大于 2，否则它根本不会收敛并且值相差很大。
Python：
eigsh(M,k=6,which='SA')并且eigsh(M,k=6,which='SM')两者都不收敛（未达到收敛的 ARPACK 错误）。只eigsh(M,k=6,sigma=0.0)给出了一些特征值（大约一个小时后），对于最小的八度音阶是不同的（甚至发现了 1 个额外的小值）：
[3.82923317e-17 3.32269886e-16 2.78039665e-10 7.54202273e-07 7.56251500e-06 1.01893934e-05]
如果容差足够高，我也会从 中得到结果eigsh(M,k=6,which='SA',tol='1')，这些结果接近于其他获得的值
[4.28732218e-14 7.54194948e-07 7.56220703e-06 1.01889544e-05, 1.87247350e-05 2.02652719e-05]
再次使用不同数量的小特征值。计算时间仍然接近30分钟。虽然不同的非常小的值可能是可以理解的，因为它们可能代表 0 的倍数，不同的多重性让我感到困惑。
此外，SciPy 和 Octave 似乎存在一些根本差异，但我还无法弄清楚。
我首先想说的是，我不知道为什么你和@Bill 报告的结果是这样的。我只是想知道eigs(M,6,0)在 Octave 中是否对应于k=6 & sigma=0，或者也许是其他东西？
使用 scipy，如果我不提供 sigma，我可以通过这种方式在相当长的时间内得到结果。
import numpy as np
from scipy.sparse import csr_matrix
from scipy.sparse.linalg import eigsh
from time import perf_counter
M = np.load('M.npz')
a = csr_matrix((M['data'], M['indices'], M['indptr']), shape=M['shape'])
t = perf_counter()
b, c = eigsh(a, k=50, which='SA', tol=1e1)
print(perf_counter() - t)
print(b)
我完全不确定这是否有意义。
0.4332823531003669
[4.99011753e-03 3.32467891e-02 8.81752215e-02 1.70463893e-01
2.80811313e-01 4.14752072e-01 5.71103821e-01 7.53593653e-01
9.79938915e-01 1.14003837e+00 1.40442848e+00 1.66899183e+00
1.96461415e+00 2.29252666e+00 2.63050114e+00 2.98443218e+00
3.38439528e+00 3.81181747e+00 4.26309942e+00 4.69832271e+00
5.22864462e+00 5.74498014e+00 6.22743988e+00 6.83904055e+00
7.42379697e+00 7.97206446e+00 8.62281827e+00 9.26615266e+00
9.85483434e+00 1.05915030e+01 1.11986296e+01 1.18934953e+01
1.26811461e+01 1.33727614e+01 1.41794599e+01 1.47585155e+01
1.55702295e+01 1.63066947e+01 1.71564622e+01 1.78260727e+01
1.85693454e+01 1.95125277e+01 2.01847294e+01 2.09302671e+01
2.18860389e+01 2.25424795e+01 2.32907153e+01 2.37425085e+01
2.50784800e+01 2.55119112e+01]
我发现使用 sigma 并在适当的时间内获得结果的唯一方法是将 M 作为 LinearOperator 提供。我对这个东西不太熟悉，但据我了解，我的实现代表一个单位矩阵，如果在调用中未指定，则 M 应该是这个矩阵。原因是 scipy 将使用迭代求解器，而不是执行直接求解（LU 分解），这可能更适合。作为比较，如果您提供M = np.identity(a.shape[0])，它应该完全相同，那么 eigsh 需要很长时间才能产生结果。请注意，如果提供了此方法，则该方法不起作用sigma=0。但我不确定是否sigma=0真的那么有用？
import numpy as np
from scipy.sparse import csr_matrix
from scipy.sparse.linalg import eigs, eigsh, LinearOperator
from time import perf_counter
def mv(v):
return v
M = np.load('M.npz')
a = csr_matrix((M['data'], M['indices'], M['indptr']), shape=M['shape'])
t = perf_counter()
b, c = eigsh(a, M=LinearOperator(shape=a.shape, matvec=mv, dtype=np.float64),
sigma=5, k=50, which='SA', tol=1e1, mode='cayley')
print(perf_counter() - t)
print(np.sort(-5 * (1 + b) / (1 - b)))
同样，不知道它是否正确，但肯定与以前不同。如果有来自 scipy 的人的意见那就太好了。
1.4079377939924598
[3.34420263 3.47938816 3.53019328 3.57981026 3.60457277 3.63996294
3.66791416 3.68391585 3.69223712 3.7082205 3.7496456 3.76170023
3.76923989 3.80811939 3.81337342 3.82848729 3.84137264 3.85648208
3.88110869 3.91286153 3.9271108 3.94444577 3.97580798 3.98868207
4.01677424 4.04341426 4.05915855 4.08910692 4.12238969 4.15283192
4.16871081 4.1990492 4.21792125 4.24509036 4.26892806 4.29603036
4.32282475 4.35839271 4.37934257 4.40343219 4.42782208 4.4477206
4.47635849 4.51594603 4.54294049 4.56689989 4.58804775 4.59919363
4.63700551 4.66638214]