[PATCH 23/25] math: Use erfc from CORE-MATH

[Paul Zimmermann]

       Hi DJ,

> > +  for (int i = 0; i < 22; i++)
> > +    if (x == exceptions[i][0])
> > +      return exceptions[i][1] + exceptions[i][2];
> 
> IMHO the need for this code implies that all your other algorithms are
> faulty, if you need to fix their errors...

not necessarily. Correct rounding is ensured using Ziv's algorithm:

1) [fast path] compute an approximation y of f(x) using say 16 guard bits,
   and an error bound e
2) [rounding test] if y-e and y+e round to the same value, return this value
3) [accurate path] using an approximation y of f(x) using a larger precision,
   and round it to the target precision

The correctness of the accurate path is validated by the knowledge of "hard-
to-round" inputs. However, if the accurate path returns a wrong result for some
hard-to-round input, it is sometimes more efficient to use it as an
exceptional value, than increasing the precision of the accurate path.

> The exception table could have been sorted, too.
> 
> > +  /* subnormal exceptions */
> > +  if (x == 0x1.a8f7bfbd15495p+4)
> > +    return fma (0x1p-1074, -0.25, 0x1.99ef5883f656cp-1024);
> 
> Could have been part of the exceptions table?

I'm not sure, since the output should be 0x1.99ef5883f656cp-1024 for rounding
to nearest and upward, and 0x1.99ef5883f6568p-1024 for roundings toward zero
and downward. For the other exceptional cases, we return round(yh+yl) where
yh+yl is a double-double approximation of erfc(x). Here since we are in the
subnormal case, I cannot find such yh and yl values.

Paul