From FreeBSD's log:

From FreeBSD's log:

Log:
Fixed about 50 million errors of infinity ulps and about 3 million errors
of between 1.0 and 1.8509 ulps for lgammaf(x) with x between -2**-21 and
-2**-70.

As usual, the cutoff for tiny args was not correctly translated to
float precision.  It was 2**-70 but 2**-21 works.  Not as usual, having
a too-small threshold was worse than a pessimization.  It was just a
pessimization for (positive) args between 2**-70 and 2**-21, but for
the first ~50 million (negative) args below -2**-70, the general code
overflowed and gave a result of infinity instead of correct (finite)
results near 70*log(2).  For the remaining ~361 million negative args
above -2**21, the general code gave almost-acceptable errors (lgamma[f]()
is not very accurate in general) but the pessimization was larger than
for misclassified tiny positive args.

Now the max error for lgammaf(x) with |x| < 2**-21 is 0.7885 ulps, and
speed and accuracy are almost the same for positive and negative args
in this range.  The maximum error overall is still infinity ulps.

A cutoff of 2**-70 is probably wastefully small for the double precision
case.  Smaller cutoffs can be used to reduce the max error to nearly
0.5 ulps for tiny args, but this is useless since the general algrorithm
for nearly-tiny args is not nearly that accurate -- it has a max error of
about 1 ulp.

Obtained-from:  FreeBSD