More precise mean_var for non-Dask dense arrays?

[flying-sheep]

In #127, I changed how mean_var works just for Dask:

fast-array-utils/src/fast_array_utils/stats/_mean_var.py

Line 40 in 940727f

mean_sq = mean(power(x, 2, dtype=np.float64), axis=axis) if isinstance(x, types.DaskArray) else mean(power(x, 2), axis=axis, dtype=np.float64)

We already have a custom numba implementation for sparse-in-core which operates on float64s, but doesn’t convert everything at once.

We could make the change also affect dense-in-core. Should we

• leave things as they are and keep power for float32 dense arrays as float32? It takes less space and is faster, but less precise.
• do exponentiation in 64 bit everywhere? this would allocate an exponentiated array of the input size in 64 bit, but I guess if people’s machines can handle a temporary 2n (prod(shape)×32×2), they can handle a temporary 3n (prod(shape)×32 + prod(shape)×64)?
• implement a custom numba implementation for dense-in-core as well?