5

从一个模拟问题,我想在 MATLAB 中计算 1000x1000 的复杂方阵。由于这些值是指贝塞尔函数的值,因此矩阵一点也不稀疏。

由于我对行列式相对于某些参数(在我的情况下搜索的特征函数的能量)的变化感兴趣,因此我通过首先搜索所研究范围的重新缩放因子然后计算行列式来克服目前的问题,

result(k) = det(pre_factor*Matrix{k});

现在这是一个非常尴尬的解决方案,仅适用于最大 500x500 的矩阵尺寸。

有人知道这个问题的一个很好的解决方案吗?与 Mathematica 的接口原则上可能有效,但我对可行性表示怀疑。先感谢您

罗伯特

编辑:我没有找到计算问题的方便解决方案,因为这需要更改为更高的精度。相反,我用那个

ln det M = trace ln M

也就是说,当我对 k 进行推导时

A = trace(inv(M(k))*dM/dk)

所以我至少有关于k的行列式的对数的变化。从问题的物理背景,我可以得出对 A 的约束,这最终给了我一个对我的问题有效的解决方法。不幸的是,我不知道这种解决方法是否可以推广。

4

4 回答 4

5

您应该意识到,当您将矩阵乘以常数 k 时,您会将矩阵的行列式缩放 k^n,其中 n 是矩阵的维数。因此,对于 n = 1000 和 k = 2,您将行列式缩放为

>> 2^1000
ans =
     1.07150860718627e+301

这当然是一个巨大的数字,所以您可能会认为它会失败,因为在双精度中,MATLAB 只会表示与 realmax 一样大的浮点数。

>> realmax
ans =
     1.79769313486232e+308

没有必要重新计算该行列式的所有工作,而不是计算这样一个巨大矩阵的行列式无论如何都是一个非常适合的问题。

于 2010-12-03T10:44:33.387 回答
4

如果速度不是问题,您可能希望使用det(e^A) = e^(tr A)并将A矩阵乘以一些缩放常数(因此A - I光谱半径小于一)。

编辑:在 MatLab 中,矩阵 ( logm) 的对数是通过三角化计算的。因此,您最好计算矩阵的特征值并将它们相乘(或者更好的是,添加它们的对数)。您没有指定您的矩阵是否对称:如果是,则查找特征值比不是更容易。

于 2010-12-03T10:49:55.743 回答
1

您说行列式的当前值约为 10^-300。

您是否试图将行列式设为某个值,例如 1?如果是这样,重新缩放很尴尬:您正在考虑的矩阵是病态的,并且考虑到机器的精度,您应该考虑输出行列式为零。换句话说,不可能得到可靠的逆。

我建议修改矩阵的列或行,而不是重新调整它。


我用 R 用随机矩阵(随机正常值)进行了一个小测试,似乎行列式显然应该是非零的。

> n=100
> M=matrix(rnorm(n**2),n,n)
> det(M)
[1] -1.977380e+77
> kappa(M)
[1] 2318.188
于 2010-12-03T11:28:15.763 回答
1

这不是严格意义上的 matlab 解决方案,但您可能需要考虑使用Mahout。它专为大规模线性代数而设计。(1000x1000 对于它习惯的比例来说没有问题。)

您将调用 java以向/从 Mahout 传递数据。

于 2010-12-03T18:00:53.120 回答