python - numpy.arrays 的 fsum，稳定求和

Question

我有许多具有小值的多维 numpy.array 我需要加上一点数字错误。对于 float，有 math.fsum (及其实现 here )，这一直对我很有帮助。 numpy.sum不够稳定。

如何获得 numpy.array 的稳定求和？

---

背景

这是给 the quadpy package 的.小值数组是函数在(许多)间隔的特定点乘以它们的权重的评估。这些的总和是所述函数在区间内的积分的近似值。

Answer 1

好吧，我已经实现了 accupy给出了一些稳定的求和算法。

这是 Kahan summation 的一个快速而肮脏的实现对于 numpy 数组。但是请注意，它对于病态总和不是很准确。

def kahan_sum(a, axis=0):
    '''Kahan summation of the numpy array along an axis.
    '''
    s = numpy.zeros(a.shape[:axis] + a.shape[axis+1:])
    c = numpy.zeros(s.shape)
    for i in range(a.shape[axis]):
        # https://stackoverflow.com/a/42817610/353337
        y = a[(slice(None),) * axis + (i,)] - c
        t = s + y
        c = (t - s) - y
        s = t.copy()
    return s

它完成了工作，但速度很慢，因为它是在第 axis 维度上进行 Python 循环。