大矩阵对角化工具

Question

我想计算扩散核，其中涉及对一个大矩阵A取exp(b*A)。为了玩弄b的值，我想将A对角化（这样exp(A)就可以快速运行）。

我的矩阵大约是25k x 25k，但非常稀疏 - 只有约60k个值是非零的。Matlab的“eigs”函数以及Octave的“eig”和R的“eigen”都会因内存不足而运行失败。是否有一种工具可以找到大型稀疏矩阵的分解？

不知道这是否相关，但A是邻接矩阵，因此它是对称的，并且是满秩的。

- Xodarap

Answer 1

你尝试过在Matlab中使用SVD和针对稀疏矩阵的“svds”吗？

编辑：还有一件事，由于维度很大，请不要进行完整秩的SVD，使用一个小秩（例如500）即可使你的解决方案适合内存。这会削减小的特征值及其向量，因此对精度影响不大。

- Yin Zhu

---

@ A = USV'。S是一个方阵，其特征值在其对角线上。 - Yin Zhu

@duffymo：请参考维基百科关于奇异值分解和特征值分解之间关系的内容。无论如何，OP正在询问A的对角化，而svd可以实现这一点。 - vad

我在我的电脑上尝试了这个，但它显示“超出最大变量大小”。虽然我的电脑不是顶级配置，所以问题可能不在Matlab上。 - Xodarap

@ Xodarap。Matlab的svds并不是最好的选择。你也可以尝试http://tedlab.mit.edu/~dr/SVDLIBC/。当你开始时，首先通过尝试小矩阵来确保SVD能够胜任你的工作。 - Yin Zhu

@Yin - 这是一个很好的想法，但你知道会丢失多少信息吗？我觉得由于这是一个邻接矩阵，而且行之间相互独立，只找到10%的特征值（或其他）将会丢失90%的信息。 - Xodarap

@Xodarap。非常少。http://en.wikipedia.org/wiki/Singular_value_decomposition#Matrix_approximation您还可以在小矩阵上进行一些实验，以查看近似值和精确值之间的差异。 - Yin Zhu

Answer 2

如果您有访问 64 位计算机并且已编译具有 64 位支持的 Octave，则可能可以解决此问题。

此外，我不知道您在什么平台上运行所有内容，但在基于 UNIX 的系统中，您可以使用 ulimit 来增加用户进程允许的最大堆栈大小。

例如，您可以运行

ulimit -u unlimited

这样做可以确保您的进程没有内存限制等问题。总的来说，这不是一个好主意，因为您可能会遇到失控的进程，从而完全拖垮您的机器。建议尝试以下方法：

ulimit -s [stacksize]

为增加栈大小限制。

- Marcus P S

Answer 3

你考虑过以下的性质吗： exp(A*t) = L^(-1) {(sI-A)^(-1)} 其中L^(-1)是拉普拉斯变换的反演？- 假设你可以反演(sI-A)

- gd047

---

我认为反拉普拉斯变换的时间复杂度是O(n^4)，非常慢。而且容易出错。请参考http://www.cs.cornell.edu/cv/ResearchPDF/19ways+.pdf。 - Xodarap

Answer 4

Octave有splu函数，可以对稀疏矩阵进行LU分解。我不确定它是否能处理25k x 25k的矩阵，但值得一试。

另外，如果您的矩阵结构如下：A = [B zeros;zeros C]，那么您可以分别对B和C进行对角化，然后将它们组合成一个矩阵。我猜您可以对eig做类似的操作。

- vad

---

愚蠢的问题，但是一旦我有了LU分解，我该如何使用它？我意识到det(A) = det(L)*det(U)，因此它们的对角化之间必须存在某种关系，但我不够聪明以看出来。 - Xodarap

Answer 5

在R中，你可以检查igraph包和函数arpack，它是ARPACK库的接口。

- Marek