scipy.sparse.linalg.lobpcg¶

scipy.sparse.linalg.lobpcg(A, X, B=None, M=None, Y=None, tol=None, maxiter=None, largest=True, verbosityLevel=0, retLambdaHistory=False, retResidualNormsHistory=False)[源代码]¶

局部最优挡路预条件共轭梯度法

LOBPCG是求解大型对称正定(SPD)广义特征问题的预处理特征解算器。

参数

A{稀疏矩阵、稠密矩阵、线性运算符}: 问题的对称线性算子，通常是稀疏矩阵。通常称为“刚度矩阵”。
Xndarray、Float32或Float64: 的初始近似 k 特征向量(非稀疏)。如果 A 有 shape=(n,n) 然后 X 应该有形状 shape=(n,k) 。
B{密集矩阵，稀疏矩阵，线性运算符}，可选: 广义本征问题中的右侧算子。默认情况下， B = Identity 。通常被称为“质量矩阵”。
M{密集矩阵，稀疏矩阵，线性运算符}，可选: 预调剂到 A ；默认情况下 M = Identity 。 M 应该近似于 A 。
Yndarray、Float32或Float64，可选: n×大小约束矩阵(非稀疏)，大小<n迭代将在Y的列空间的B正交互补中执行。Y必须是满秩的。
tol标量，可选: 求解器公差(停止标准)。默认值为 tol=n*sqrt(eps) 。
maxiter整型，可选: 最大迭代次数。默认值为 maxiter = 20 。
largest布尔值，可选: 如果为True，则求解最大的特征值，否则求解最小的特征值。
verbosityLevel整型，可选: 控制解算器输出。默认值为 verbosityLevel=0 。
retLambdaHistory布尔值，可选: 是否返回特征值历史记录。默认值为False。
retResidualNormsHistory布尔值，可选: 是否返回历史残差规范。默认值为False。

退货

wndarray: 数组 k 特征值
vndarray: 一组 k 特征向量。 v 形状与 X 。
lambdasndarray列表，可选: 特征值历史，如果 retLambdaHistory 是真的。
rnormsndarray列表，可选: 剩余规范的历史，如果 retResidualNormsHistory 是真的。

注意事项

如果两者都有 retLambdaHistory 和 retResidualNormsHistory 为True，则返回元组的格式如下 (lambda, V, lambda history, residual norms history) 。

在以下内容中 n 表示矩阵大小，并且 m 所需特征值的数量(最小或最大)。

LOBPCG代码在内部解决大小的本征问题 3m 在每次迭代中调用“标准”密集特征解析器，因此如果 m 与之相比还不够小 n ，调用LOBPCG代码没有意义，而应该使用“标准”的特征解析器，例如，在本例中使用Numpy或Scipy函数。如果调用LOBPCG算法 5m > n ，它很可能会在内部中断，因此代码尝试调用标准函数。

不是那样的 n 对于LOBPCG来说应该很大，而不是比率 n / m 应该很大。它是您称之为LOBPCG的 m=1 和 n=10 ，但它还是有效的 n 很小。该方法适用于非常大的 n / m [4].

衔接速度基本上取决于两个因素：

求特征值与特征值的睡觉的相对分离程度。你可以试着改变 m 让这一切变得更好。
这个问题的条件有多好。这可以通过使用适当的预处理来改变。例如，杆振动测试问题(在测试目录下)对于大型 n 因此，除非使用有效的预处理，否则收敛将很慢。对于这个特定的问题，一个好的简单的预条件函数将是一个线性解 A ，这很容易编码，因为A是三对角线的。

参考文献

1: A.V.Knyazev(2001年)提出的最优预处理特征解析器：局部最优挡路预处理共轭梯度法。暹罗科学计算杂志23，第2期，第517-541页。 DOI:10.1137/S1064827500366124
2: A.V.Knyazev，I.Lashuk，M.E.阿根廷人和E.Ovchinnikov(2007年)，挡路在HYPRE和PETSC中的局部最优预条件特征值Xolver(BLOPEX)。 arXiv:0705.2626
3: A.V.Knyazev的C语言和MATLAB实现：https://bitbucket.org/joseroman/blopex
4: S·Yamada，T·Imamura，T·Kano和M.Machida(2006)，地球模拟器上量子多体问题精确数值方法的高性能计算。刊登在2006年ACM/IEEE超级计算会议论文集上。 DOI:10.1145/1188455.1188504

示例

求解 A x = lambda x 有约束和预处理。

>>> import numpy as np
>>> from scipy.sparse import spdiags, issparse
>>> from scipy.sparse.linalg import lobpcg, LinearOperator
>>> n = 100
>>> vals = np.arange(1, n + 1)
>>> A = spdiags(vals, 0, n, n)
>>> A.toarray()
array([[  1.,   0.,   0., ...,   0.,   0.,   0.],
       [  0.,   2.,   0., ...,   0.,   0.,   0.],
       [  0.,   0.,   3., ...,   0.,   0.,   0.],
       ...,
       [  0.,   0.,   0., ...,  98.,   0.,   0.],
       [  0.,   0.,   0., ...,   0.,  99.,   0.],
       [  0.,   0.,   0., ...,   0.,   0., 100.]])

约束条件：

>>> Y = np.eye(n, 3)

对于特征向量的初始猜测，应该有线性独立的列。列维度=请求的特征值的数量。

>>> rng = np.random.default_rng()
>>> X = rng.random((n, 3))

本例中A的倒数预处理器：

>>> invA = spdiags([1./vals], 0, n, n)

前置条件函数必须由函数定义：

>>> def precond( x ):
...     return invA @ x

预条件函数的自变量x是内部的矩阵 lobpcg ，从而使用矩阵-矩阵乘积 @ 。

预处理器函数作为 LinearOperator ：

>>> M = LinearOperator(matvec=precond, matmat=precond,
...                    shape=(n, n), dtype=float)

现在让我们解决矩阵A的特征值问题：

>>> eigenvalues, _ = lobpcg(A, X, Y=Y, M=M, largest=False)
>>> eigenvalues
array([4., 5., 6.])

请注意，Y中传递的向量是3个最小特征值的特征向量。返回的结果与那些结果正交。

scipy.sparse.linalg.eigsh

scipy.sparse.linalg.svds