Fix lgammaf spurious underflows (bug 18220).

The flt-32 implementation of lgammaf produces spurious underflow
exceptions for some large arguments, because of calculations involving
x^-2 multiplied by small constants.  This patch fixes this by
adjusting the threshold for a simpler computation to 2**26 (the error
in the simpler computation is on the order of 0.5 * log (x), for a
result on the order of x * log (x)).

Tested for x86_64 and x86.

	[BZ #18220]
	* sysdeps/ieee754/flt-32/e_lgammaf_r.c (__ieee754_lgammaf_r): Use
	2**26 not 2**58 as threshold for returning x * (log (x) - 1).
	* math/auto-libm-test-in: Add another test of lgamma.
	* math/auto-libm-test-out: Regenerated.

2 files changed, 26 insertions, 0 deletions

diff --git a/math/auto-libm-test-in b/math/auto-libm-test-in
index c4bfe7439c..2a88403aa7 100644
--- a/math/auto-libm-test-in
+++ b/math/auto-libm-test-in
@@ -1757,6 +1757,7 @@ lgamma 0.5
 lgamma -0.5
 lgamma 0.7
 lgamma 1.2
+lgamma 0x3.8p56
 lgamma 0x1p-5
 lgamma -0x1p-5
 lgamma 0x1p-10
diff --git a/math/auto-libm-test-out b/math/auto-libm-test-out
index d717ea4431..0aa7b75c3b 100644
--- a/math/auto-libm-test-out
+++ b/math/auto-libm-test-out
@@ -139612,6 +139612,31 @@ lgamma 1.2
 = lgamma tonearest ldbl-128ibm 0x1.33333333333333333333333333p+0L : -0x1.5db138c7d70c72cb0b57c9e83d8p-4L 1 : inexact-ok
 = lgamma towardzero ldbl-128ibm 0x1.33333333333333333333333333p+0L : -0x1.5db138c7d70c72cb0b57c9e83dp-4L 1 : inexact-ok
 = lgamma upward ldbl-128ibm 0x1.33333333333333333333333333p+0L : -0x1.5db138c7d70c72cb0b57c9e83dp-4L 1 : inexact-ok
+lgamma 0x3.8p56
+= lgamma downward flt-32 0x3.8p+56f : 0x8.8bdd4p+60f 1 : inexact-ok
+= lgamma tonearest flt-32 0x3.8p+56f : 0x8.8bdd4p+60f 1 : inexact-ok
+= lgamma towardzero flt-32 0x3.8p+56f : 0x8.8bdd4p+60f 1 : inexact-ok
+= lgamma upward flt-32 0x3.8p+56f : 0x8.8bdd5p+60f 1 : inexact-ok
+= lgamma downward dbl-64 0x3.8p+56 : 0x8.8bdd41bf4484p+60 1 : inexact-ok
+= lgamma tonearest dbl-64 0x3.8p+56 : 0x8.8bdd41bf44848p+60 1 : inexact-ok
+= lgamma towardzero dbl-64 0x3.8p+56 : 0x8.8bdd41bf4484p+60 1 : inexact-ok
+= lgamma upward dbl-64 0x3.8p+56 : 0x8.8bdd41bf44848p+60 1 : inexact-ok
+= lgamma downward ldbl-96-intel 0x3.8p+56L : 0x8.8bdd41bf4484605p+60L 1 : inexact-ok
+= lgamma tonearest ldbl-96-intel 0x3.8p+56L : 0x8.8bdd41bf4484605p+60L 1 : inexact-ok
+= lgamma towardzero ldbl-96-intel 0x3.8p+56L : 0x8.8bdd41bf4484605p+60L 1 : inexact-ok
+= lgamma upward ldbl-96-intel 0x3.8p+56L : 0x8.8bdd41bf4484606p+60L 1 : inexact-ok
+= lgamma downward ldbl-96-m68k 0x3.8p+56L : 0x8.8bdd41bf4484605p+60L 1 : inexact-ok
+= lgamma tonearest ldbl-96-m68k 0x3.8p+56L : 0x8.8bdd41bf4484605p+60L 1 : inexact-ok
+= lgamma towardzero ldbl-96-m68k 0x3.8p+56L : 0x8.8bdd41bf4484605p+60L 1 : inexact-ok
+= lgamma upward ldbl-96-m68k 0x3.8p+56L : 0x8.8bdd41bf4484606p+60L 1 : inexact-ok
+= lgamma downward ldbl-128 0x3.8p+56L : 0x8.8bdd41bf44846050819264e2d568p+60L 1 : inexact-ok
+= lgamma tonearest ldbl-128 0x3.8p+56L : 0x8.8bdd41bf44846050819264e2d57p+60L 1 : inexact-ok
+= lgamma towardzero ldbl-128 0x3.8p+56L : 0x8.8bdd41bf44846050819264e2d568p+60L 1 : inexact-ok
+= lgamma upward ldbl-128 0x3.8p+56L : 0x8.8bdd41bf44846050819264e2d57p+60L 1 : inexact-ok
+= lgamma downward ldbl-128ibm 0x3.8p+56L : 0x8.8bdd41bf44846050819264e2d4p+60L 1 : inexact-ok
+= lgamma tonearest ldbl-128ibm 0x3.8p+56L : 0x8.8bdd41bf44846050819264e2d4p+60L 1 : inexact-ok
+= lgamma towardzero ldbl-128ibm 0x3.8p+56L : 0x8.8bdd41bf44846050819264e2d4p+60L 1 : inexact-ok
+= lgamma upward ldbl-128ibm 0x3.8p+56L : 0x8.8bdd41bf44846050819264e2d8p+60L 1 : inexact-ok
 lgamma 0x1p-5
 = lgamma downward flt-32 0x8p-8f : 0x3.72d02cp+0f 1 : inexact-ok
 = lgamma tonearest flt-32 0x8p-8f : 0x3.72d03p+0f 1 : inexact-ok