numpy (superb library - no slight intended here!) does not seem to support Kahan or butterfly summation in cumsum, which is critically important for accuracy in large arrays of numbers of similar magnitude.

Reproducing code example:

import numpy as np; print(np.cumsum(np.ones((2**28,), dtype=np.float32))[-1] - 2**28)

import numpy as np
print(np.cumsum(np.ones((2**28,), dtype=np.float32))[-1] - 2**28)
# Should be zero

Note that sum obviously uses butterfly summation:

print(np.sum(np.ones((2**28,), dtype=np.float32)) - 2**28)
# Is zero on any machine I can find