Neumaier summation is an improvement of Kahan summation for accurately summing arrays of floats.

import numba as nb

@nb.njit
def neumaier_sum(arr):
    s = arr[0]
    c = 0.0
    for i in range(1, len(arr)):
        t = s + arr[i]
        if abs(s) >= abs(arr[i]):
            c += (s - t) + arr[i]
        else:
            c += (arr[i] - t) + s
        s = t
    return s + c

This works well but it is at least four times slower than it would be were you to add fastmath=True. Unfortunately, fastmath allows rebracketing the sum (associativity) which has the effect of ruining its accuracy so we can’t do that.

Is there any way to speed up this code but have it output exactly what it does now?