Fix incorrect hypotl(3) result with subnormal numbers

This adjusts the factor used to scale the subnormal numbers, so it becomes the right value after adjusting its exponent. Thanks to Steve Kargl for finding the most elegant fix. Also enable the hypot tests, and add a test case for this bug. PR: 253313 MFC after: 1 week
freebsd · Feb 10, 2021 · d3338f3 · d3338f3
1 parent a6dc68c
commit d3338f3
Show file tree

Hide file tree

Showing 3 changed files with 22 additions and 1 deletion.
diff --git a/contrib/netbsd-tests/lib/libm/t_hypot.c b/contrib/netbsd-tests/lib/libm/t_hypot.c
@@ -70,12 +70,32 @@ ATF_TC_BODY(pr50698, tc)
 	ATF_CHECK(!isnan(val));
 }
 
+ATF_TC(hypotl_near_underflow);
+ATF_TC_HEAD(hypotl_near_underflow, tc)
+{
+	atf_tc_set_md_var(tc, "descr", "Test hypotl near underflow");
+}
+
+ATF_TC_BODY(hypotl_near_underflow, tc)
+{
+	volatile long double a = 0x1.b2933cafa0bb7p-16383L;
+	volatile long double b = 0x1.fffffffffffffp-16351L;
+	volatile long double e = 0x1.fffffffffffffp-16351L;
+	volatile long double ulp = __LDBL_EPSILON__;
+
+	volatile long double val = hypotl(a, b);
+
+	ATF_CHECK(!isinf(val));
+	ATF_CHECK(fabsl(val - e) <= 2 * ulp);
+}
+
 ATF_TP_ADD_TCS(tp)
 {
 
 	ATF_TP_ADD_TC(tp, hypot_integer);
 	ATF_TP_ADD_TC(tp, hypotf_integer);
 	ATF_TP_ADD_TC(tp, pr50698);
+	ATF_TP_ADD_TC(tp, hypotl_near_underflow);
 
 	return atf_no_error();
 }
diff --git a/lib/msun/src/e_hypotl.c b/lib/msun/src/e_hypotl.c
@@ -82,7 +82,7 @@ hypotl(long double x, long double y)
 	        man_t manh, manl;
 		GET_LDBL_MAN(manh,manl,b);
 		if((manh|manl)==0) return a;
-		t1=0;
+		t1=1;
 		SET_HIGH_WORD(t1,ESW(MAX_EXP-2));	/* t1=2^(MAX_EXP-2) */
 		b *= t1;
 		a *= t1;

diff --git a/lib/msun/tests/Makefile b/lib/msun/tests/Makefile
@@ -31,6 +31,7 @@ NETBSD_ATF_TESTS_C+=	erf_test
 NETBSD_ATF_TESTS_C+=	exp_test
 NETBSD_ATF_TESTS_C+=	fmod_test
 NETBSD_ATF_TESTS_C+=	fe_round_test
+NETBSD_ATF_TESTS_C+=	hypot_test
 NETBSD_ATF_TESTS_C+=	infinity_test
 NETBSD_ATF_TESTS_C+=	ilogb_test
 NETBSD_ATF_TESTS_C+=	ldexp_test