tgamma: raise FE_OVERFLOW (not FE_DIVBYZERO) on overflow (#325)

openlibm signals math errors via FP exception flags, not errno
(math_errhandling == MATH_ERREXCEPT), and tgamma already raised the right
flags everywhere except large-argument overflow: the x > 171.63 path
returned `x / zero`, which yields +Inf but raises FE_DIVBYZERO. That flag
is reserved for the pole at x == 0; a finite argument whose true gamma
exceeds DBL_MAX is an overflow and must raise FE_OVERFLOW (as the function
comment already documents, and as glibc does).

Return `x * 0x1p1023` instead, which overflows to +Inf and raises
FE_OVERFLOW. The pole (FE_DIVBYZERO) and negative-integer domain error
(FE_INVALID) were already correct and are unchanged.

This keeps the library consistent with its own math_errhandling rather
than introducing errno, which no other openlibm function sets.

Adds test/regression/test-325.c, which verifies the exception for each
error class against the result value; it self-checks that the
environment tracks FP flags first (qemu-user emulates them unreliably on
some arches, cf. #347) and skips otherwise.

Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>

Author: ViralBShah

bsdsrc/b_tgamma.c

@@ -131,7 +131,14 @@ tgamma(x)
 
 	if (isgreaterequal(x, 6)) {
 		if(x > 171.63)
-			return (x / zero);
+			/*
+			 * Overflow: the true value exceeds DBL_MAX.  Return +Inf
+			 * and raise FE_OVERFLOW (per math_errhandling ==
+			 * MATH_ERREXCEPT and the function comment above).  The
+			 * old `x / zero` returned +Inf but wrongly raised
+			 * FE_DIVBYZERO, which is reserved for the pole at x == 0.
+			 */
+			return (x * 0x1p1023);
 		u = large_gam(x);
 		return(__exp__D(u.a, u.b));
 	} else if (isgreaterequal(x, 1.0 + LEFT + x0))

test/regression/test-325.c

@@ -0,0 +1,65 @@
+/*
+ * Regression test for issue #325: tgamma() error signaling.
+ * https://github.com/JuliaMath/openlibm/issues/325
+ *
+ * openlibm reports math errors via FP exception flags, not errno
+ * (math_errhandling == MATH_ERREXCEPT).  tgamma already raised FE_DIVBYZERO at
+ * the pole (x == +-0) and FE_INVALID for negative-integer (domain) arguments,
+ * but overflow for large x went through `x / zero`, which returns +Inf while
+ * wrongly raising FE_DIVBYZERO instead of FE_OVERFLOW.  This pins the correct,
+ * glibc-matching exception for every error class.
+ *
+ * FP exception flags are emulated unreliably by qemu-user on some arches (see
+ * #347), so the test first checks that the environment tracks a plain
+ * divide-by-zero and skips if it does not.
+ */
+#include <math.h>
+#include <openlibm_fenv.h>
+
+#include "regress-util.h"
+
+static int
+fenv_tracks_flags(void)
+{
+	volatile double z = 0.0, r;
+
+	feclearexcept(FE_ALL_EXCEPT);
+	r = 1.0 / z;			/* a known divide-by-zero */
+	(void)r;
+	return fetestexcept(FE_DIVBYZERO) != 0;
+}
+
+static int
+raises(double x, int excepts)
+{
+	feclearexcept(FE_ALL_EXCEPT);
+	(void)tgamma(x);
+	return fetestexcept(FE_DIVBYZERO | FE_INVALID | FE_OVERFLOW) == excepts;
+}
+
+int
+main(void)
+{
+	if (!fenv_tracks_flags())
+		return REGRESS_SKIP;
+
+	/* Pole at x == +-0: +-Inf, FE_DIVBYZERO. */
+	CHECK(isinf(tgamma(0.0)) == 1);
+	CHECK(raises(0.0, FE_DIVBYZERO));
+	CHECK(isinf(tgamma(-0.0)) == -1);
+	CHECK(raises(-0.0, FE_DIVBYZERO));
+
+	/* Negative integer (domain error): NaN, FE_INVALID. */
+	CHECK(isnan(tgamma(-2.0)));
+	CHECK(raises(-2.0, FE_INVALID));
+
+	/* Overflow for large finite x: +Inf, FE_OVERFLOW (was FE_DIVBYZERO). */
+	CHECK(isinf(tgamma(200.0)) == 1);
+	CHECK(raises(200.0, FE_OVERFLOW));
+
+	/* In-range result raises none of these. */
+	CHECK(tgamma(5.0) == 24.0);
+	CHECK(raises(5.0, 0));
+
+	return 0;
+}