`dotv`: float32 rounding error is 10x of NumPy

[crusaderky]

I found a use case where the rounding error in float32 is ~1.3e-3, compared to ~1.3e-4 of NumPy.

    import numpy as np
    from blis import dotv

    A64 = np.full((95, ), fill_value=40.09845, dtype=np.float64)
    A64[0] = 0.0
    B64 = np.full((95, ), fill_value=2.0, dtype=np.float64)
    B64[1] = -98.0
    B64[2] = -86.0
    A32 = A64.astype(np.float32)
    B32 = B64.astype(np.float32)

    np_f64 = A64.dot(B64)  # np.float64(0.0)
    blis_f64 = dotv(A64, B64)  # 0.0
    np_f32 = A32.dot(B32)  # np.float32(-0.00012207031)
    blis_f32 = dotv(A32, B32)  # -0.001220703125

[rgommers]

Is that the built-in dot implementation in numpy, or is that going the route that's calling into BLAS? I can't quite remember the details, may be that this is small enough to not call BLAS. If it does, you're probably comparing BLIS and OpenBLAS.

[crusaderky]

It's the default for numpy-2.3.5-cp314-cp314t-manylinux_2_27_x86_64.manylinux_2_28_x86_64.whl from pypi

[rgommers]

Okay I checked, that is calling the OpenBLAS implementation, not a numpy-internal one (that depends on types, strides, etc.).

Also for reference with Accelerate, that's even farther off:

>>> A32.dot(B32)
np.float32(-0.0017089844)

So the test case may be particularly challenging?