https://gcc.gnu.org/bugzilla/show_bug.cgi?id=121942
--- Comment #7 from Linas Vepstas <linasvepstas at gmail dot com> --- Hmm. OK. I discovered I can break the aliasing by avoiding the call to llabs. What I'd really like to see in the standard is some function that compares two floats up to N ULPS. The underlying issue is that numeric algos tend to chop off a few bits of precision here and there, depending on the algo. I can live with that, as long as its "just a few bits", and this is what motivated this routine. Searching the net for the topic "comparing floats for equality" shows an ocean of spilled ink on this topic, going back decades. Having a strong, reliable "equality up to ULPS" in the standard library would go a long ways towards avoiding this kind of bug. It would also help avoid the naive "compare to epsilon" code commonly seen. For example, the flawed `if (1.0e-12 < fabs(x-y))` compare, or the slightly better but still flawed `if (1.0e-12 < fabs(1.0- x/y))` compare.
