Re: [PATCH] RFC for patch to add C99 fma/fmaf/fmal builtins

[Geert Bosch <bosch at adacore dot com>]

On Oct 4, 2010, at 19:39, Joseph S. Myers wrote:
>> I seem to recall that double-rounding is ok as long as the intermediate
>> step contains > 2N bits.  Which implies that for IEEE single (23 bits)
>> we can generically implement it with IEEE double (53 bits > 46):
>> 
>>  result = (float)((double)op0 * (double)op1) + (double)op2);
> 
> No, that's not safe.  Say (double)op0 * (double)op1 has a value exactly 
> half way between two representable float values, and (double)op2 is much 
> smaller in magnitude, so that in round-to-nearest the above becomes 
> (float)((double)op0 * (double)op1) - but the correct float result depends 
> on the sign of the tiny op2.

An algorithm for the correctly rounded FMA using just regular floating
point instructions has been published  in the IEEE TRANSACTIONS ON 
COMPUTERS, VOL. 57, NO. 4, APRIL 2008. There is no real short-cut 
for the correctly rounded sum of three floating point numbers.
I think it might be possible to devise a slightly simpler
algorithm for single precision, but as you illustrated, it doesn't
become trivial.

Still, I think it would be good in general for GCC to have such arithmetic
functions in its libgcc. The reason is that they are still self-contained,
relatively simple pure functions of their arguments, and it just doesn't
make sense having GCC do the compile-time evaluation of this function
correctly, but not the run-time evaluation.

Right now, we have the problem with GNAT that we cannot use GCC builtin 
functions in many cases, because they may fall back to an unreliable libc.
This problem is not solvable, because the number of different OS-es with
different libc versions is huge. The only solution is to have all builtin
math functions be entirely in libgcc.

  -Geert