Numpy slogdet 计算错误答案

【问题标题】：Numpy slogdet computation errorNumpy slogdet 计算错误
【发布时间】：2017-11-10 11:58:50
【问题描述】：

在计算 Vanermonde 矩阵的对数行列式时，numpy 的 slogdet 与精确结果之间似乎存在重大差异。

我与确切的日志行列式进行比较，请参阅例如here 以获得证明。

看到这个的最少代码是：

A = np.power.outer(np.linspace(0,1,50),range(50))

print np.linalg.slogdet(A)[1]

s = 0
for v1 in np.linspace(0,1,50):
    for v2 in np.linspace(0,1,50):
        if v1>v2:
            s+= np.log(v1-v2)

print s

哪个产量：

-1191.88408998
-1706.99560647

我想知道是否有更准确的对数行列式实现，我可以在这种情况下使用，也可以在非 Vandermonde 矩阵情况下使用。

【问题讨论】：

对于 N=10，值匹配，np.log(np.det(A)) 也是如此。可能到 N=30 时，det(A) 太小了，以至于即使 slogdet 也无法应付。

标签： numpy determinants

【解决方案1】：

你可以像这样使用 sympy 和 mpmath：

import numpy as np
import sympy as smp
import mpmath as mp

mp.mp.dps = 50 

linspace1 = list(map(smp.mpmath.mpf,np.linspace(0,1,50)))
A = np.power.outer(list(map(float,linspace1)),range(50))

first_print = smp.mpmath.mpf(np.linalg.slogdet(A)[1])
print(first_print)

s = 0
linspace2 = list(map(smp.mpmath.mpf,np.linspace(0,1,50)))
linspace3 = list(map(smp.mpmath.mpf,np.linspace(0,1,50)))
for v1 in linspace1:
    for v2 in linspace2:
        if v1>v2:
            s+= mp.log(v1-v2)

print(s)

结果

first_print = -1178.272517342130186079884879291057586669921875

s = -1706.9956064674289001970168329846189154212781094939

【讨论】：

谢谢，但是这两个数字仍然相差大约 50%，有什么想法吗？
@j__ 的预期值是多少？
1706.995 是使用特定于范德蒙德矩阵的公式的实际值，但 1706.272 的 slogdet 大约有 50% 的折扣
@j__ 我再次更新了答案，让我知道这些值是否正确！
我不认为 OP 抱怨迭代解决方案；这是关闭的slogdet。它更快，但不能处理具有非常小的行列式的较大数组。而mpmath 对slogdet 没有任何作用。这是使用 float64 编译的 LAPACK 函数计算得出的。