大矩阵的近似谱

计算科学 线性代数 本征系统 图论 近似
2021-11-26 00:23:02

我想计算一个大型稀疏矩阵(数十万行)的所有特征值)。这很难。

我愿意接受一个近似值。是否有近似方法可以做到这一点?

虽然我希望对这个问题有一个一般性的答案,但我也会对以下具体案例的答案感到满意。我的矩阵是一个大图的归一化拉普拉斯算子。特征值将介于 0 和 2 之间,其中大量聚集在 1 附近。

4个回答

如果您的图形是无向的(正如我所怀疑的那样),则矩阵是对称的,并且您无法比 Lanczsos 算法做得更好(如果需要稳定性,可以进行选择性重新正交化)。由于全谱由 100000 个数字组成,我猜你主要对谱密度感兴趣。

要获得近似的谱密度,请取 100 维左右的领先 Krylov 子空间的谱,并将其离散密度替换为平滑版本。

领先的 Krylov 谱将具有几乎分辨好的孤立特征值(如果存在的话),近似于非孤立谱末端的特征值,并且在两者之间有些随机,其分布的累积分布函数类似于真实谱的分布. 如果维度增长,它将以精确的算术收敛到它。(如果你的算子是无限维的,情况仍然如此,你会得到真实光谱密度函数在连续光谱上的积分。)

Arnold Neumaier 的答案在Lin Lin、Yousef Saad 和 Chao Yang (2016) 的论文“大矩阵的近似光谱密度”的第 3.2 节中进行了更详细的讨论

还讨论了其他一些方法,但本文末尾的数值分析表明,Lanczos 方法优于这些替代方法。

如果您可以考虑不是特征值但在某种意义上仍然可以告诉您有关频谱的函数的事物,那么我认为您应该查看莱斯大学的 Mark Embree 的一些工作。

这是表征光谱的另一种方法。

给定一个特征值问题 $\mathbf{A} v_k = \lambda_k v_k$(假设实对称 $\mathbf{A}$ 和分离的特征值;尽管后者可能不是必需的),让我们尝试估计拖尾谱密度 $$ S(\omega) = \sum_k \frac{\pi^{-1}\sigma}{\sigma^2 + (\lambda_k - \omega)^2} = \frac{\sigma}{\pi} \mathop {\mathrm{Tr}} [\sigma^2 + (\omega - \mathbf{A})^2]^{-1} $$ 击中后例如http://dx.doi.org/10.1016/0377- 0427(96)00018-0Avk=λkvk (assume real symmetric A and separated eigenvalues; although the latter is probably not necessary), let's attempt to estimate the smeared spectral density

S(ω)=kπ1σσ2+(λkω)2=σπTr[σ2+(ωA)2]1
After hitting e.g. 在文献检索中,我们知道迹的无偏蒙特卡罗估计是 $$ S(\omega) = \frac{\sigma}{\pi}\langle z^T [\sigma^2 + (\omega - \mathbf{A})^2]^{-1} z \rangle $$ 其中随机向量 $z$ 的每个条目包含 $+1$ 或 $-1$,每个条目的概率为 0.5。对于给定的 $\sigma$ 和 $\omega$,逆积 $[\sigma^2 + (\omega - \mathbf{A})^2]^{-1} z$ 可以用共轭来计算梯度方法,或稀疏 LU 在 $[\omega + i \sigma - \mathbf{A}]^{-1} [\omega - i \sigma - \mathbf{A}]^{-1}$ 上以最小化填充-在。这也允许估计大型矩阵的 $S(\omega)$。在实践中,似乎 CG 解决方案不需要非常准确,计算平均值也不需要很多向量。这可能取决于问题。
S(ω)=σπzT[σ2+(ωA)2]1z
where each entry of the random vector z contains either +1 or 1 with probability 0.5 for each. For given σ and ω, the inverse product [σ2+(ωA)2]1z can be computed for example with the conjugate gradient method, or sparse LU on [ω+iσA]1[ωiσA]1 to minimize fill-in. This allows estimation of S(ω) also for large matrices. In practice, it seems the CG solution doesn't need to be very accurate, and neither are many vectors necessary in computing the average. This may depend on the problem.

以上似乎比类似涂抹的 Krylov 谱密度更均匀地加权了部分光谱——尝试 diag(linspace(0, 1, 150000)) ——尽管也许有办法纠正这个问题?这有点类似于伪光谱方法,但结果表明点附近的(涂抹)特征值数,而不是到最近的特征值的反比距离。ω, rather than the inverse distance to the nearest eigenvalue.

编辑:计算上述数量的一种更好的替代方法是计算切比雪夫矩(通过与上述类似的随机评估),然后从中重建谱密度。这既不需要矩阵求逆,也不需要对每个 $\omega$ 进行单独计算。请参阅 http://theorie2.physik.uni-greifswald.de/downloads/publications/LNP_chapter19.pdf和其中的参考资料。ω. See