Pretty sure that’s not possible. More accurate for some inputs will be less accurate for others. There’s a very tricky tension in float optimization that the most predictable operation structure is a fully skewed op tree, as in naive left-to-right summation, but this is the slowest and least accurate order of operations. Using a more balanced tree is faster and more accurate (great), but unfortunately which tree shape is fastest depends very much on hardware-specific factors like SIMD width (less great). And no tree shape is universally guaranteed to be fully accurate, although a full binary tree tends to have the best accuracy, but has bad base case performance, so the actual shape that tends to get used in high performance kernels is SIMD-width parallel in a loop up to some fixed size like 256 elements, then pairwise recursive reduction above that. The recursive reduction can also be threaded. Anyway, there’s no silver bullet here.

Kahan summation is an example (also described in the top level article) of one such “gotcha”. It involves adding a term that - if floats were algebraic in this sense - would always be zero, so ffast-math often deletes it, but this actually completely removes the accuracy improvement of the algorithm

Under EForth with FP done in software:

2 f 1 f 1 f f+ f- f. 0.000 ok

PFE, I think reusing the GLIBC math library:

2e0 1e0 1e0 f+ f- f. 0.000000 ok