diff --git a/mkl_umath/src/mkl_umath_loops.c.src b/mkl_umath/src/mkl_umath_loops.c.src index 86a62c4..a35f9c6 100644 --- a/mkl_umath/src/mkl_umath_loops.c.src +++ b/mkl_umath/src/mkl_umath_loops.c.src @@ -1161,7 +1161,11 @@ pairwise_sum_@TYPE@(char *a, npy_intp n, npy_intp stride) { if (n < 8) { npy_intp i; - @type@ res = 0.; + /* + * Start with -0 to preserve -0 values. The reason is that summing + * only -0 should return -0, but `0 + -0 == 0` while `-0 + -0 == -0`. + */ + @type@ res = -0.0; for (i = 0; i < n; i++) { res += @trf@(*((@dtype@*)(a + i * stride))); @@ -2263,8 +2267,8 @@ pairwise_sum_@TYPE@(@ftype@ *rr, @ftype@ * ri, char * a, npy_intp n, if (n < 8) { npy_intp i; - *rr = 0.; - *ri = 0.; + *rr = -0.0; + *ri = -0.0; for (i = 0; i < n; i += 2) { *rr += *((@ftype@ *)(a + i * stride + 0)); *ri += *((@ftype@ *)(a + i * stride + sizeof(@ftype@)));