From 1693fd03d92f1c2779b3a79ee8a157366cc5fed6 Mon Sep 17 00:00:00 2001
From: Bruce Evans <bde@FreeBSD.org>
Date: Tue, 17 Jul 2018 12:01:59 +0000
Subject: [PATCH] Minor cleanups to csqrt*(), mostly in comments.

Remove the STDC CX_LIMITED_RANGE pragma and its verbose comment.  We still
don't have any C99 compilers (that support fenv pragmas), and if we did
then there are thousands of other places in libm that would need to use
them more than here.

The other cleanups are smaller.
---
 lib/msun/src/s_csqrt.c  | 11 +----------
 lib/msun/src/s_csqrtf.c | 18 ++++++------------
 lib/msun/src/s_csqrtl.c | 21 +++++++--------------
 3 files changed, 14 insertions(+), 36 deletions(-)

diff --git a/lib/msun/src/s_csqrt.c b/lib/msun/src/s_csqrt.c
index cde8792c36f6..ca4924fc2bdc 100644
--- a/lib/msun/src/s_csqrt.c
+++ b/lib/msun/src/s_csqrt.c
@@ -35,16 +35,7 @@ __FBSDID("$FreeBSD$");
 
 #include "math_private.h"
 
-/*
- * gcc doesn't implement complex multiplication or division correctly,
- * so we need to handle infinities specially. We turn on this pragma to
- * notify conforming c99 compilers that the fast-but-incorrect code that
- * gcc generates is acceptable, since the special cases have already been
- * handled.
- */
-#pragma	STDC CX_LIMITED_RANGE	ON
-
-/* We risk spurious overflow for components >= DBL_MAX / (1 + sqrt(2)). */
+/* For avoiding overflow for components >= DBL_MAX / (1 + sqrt(2)). */
 #define	THRESH	0x1.a827999fcef32p+1022
 
 double complex
diff --git a/lib/msun/src/s_csqrtf.c b/lib/msun/src/s_csqrtf.c
index 0f36aa0ffb51..6836a3903c25 100644
--- a/lib/msun/src/s_csqrtf.c
+++ b/lib/msun/src/s_csqrtf.c
@@ -34,20 +34,14 @@ __FBSDID("$FreeBSD$");
 
 #include "math_private.h"
 
-/*
- * gcc doesn't implement complex multiplication or division correctly,
- * so we need to handle infinities specially. We turn on this pragma to
- * notify conforming c99 compilers that the fast-but-incorrect code that
- * gcc generates is acceptable, since the special cases have already been
- * handled.
- */
-#pragma	STDC CX_LIMITED_RANGE	ON
-
 float complex
 csqrtf(float complex z)
 {
-	float a = crealf(z), b = cimagf(z);
 	double t;
+	float a, b;
+
+	a = creal(z);
+	b = cimag(z);
 
 	/* Handle special cases. */
 	if (z == 0)
@@ -82,9 +76,9 @@ csqrtf(float complex z)
 	 */
 	if (a >= 0) {
 		t = sqrt((a + hypot(a, b)) * 0.5);
-		return (CMPLXF(t, b / (2.0 * t)));
+		return (CMPLXF(t, b / (2 * t)));
 	} else {
 		t = sqrt((-a + hypot(a, b)) * 0.5);
-		return (CMPLXF(fabsf(b) / (2.0 * t), copysignf(t, b)));
+		return (CMPLXF(fabsf(b) / (2 * t), copysignf(t, b)));
 	}
 }
diff --git a/lib/msun/src/s_csqrtl.c b/lib/msun/src/s_csqrtl.c
index 8aefbeab1a0b..5ab9406558af 100644
--- a/lib/msun/src/s_csqrtl.c
+++ b/lib/msun/src/s_csqrtl.c
@@ -36,25 +36,18 @@ __FBSDID("$FreeBSD$");
 #include "math_private.h"
 
 /*
- * gcc doesn't implement complex multiplication or division correctly,
- * so we need to handle infinities specially. We turn on this pragma to
- * notify conforming c99 compilers that the fast-but-incorrect code that
- * gcc generates is acceptable, since the special cases have already been
- * handled.
- */
-#pragma	STDC CX_LIMITED_RANGE	ON
-
-/*
- * We risk spurious overflow for components >= LDBL_MAX / (1 + sqrt(2)).
- * Rather than determining the fully precise value at which we might
- * overflow, just use a threshold of approximately LDBL_MAX / 4.
+ * THRESH is now calculated portably (up to 113-bit precision).  However,
+ * the denormal threshold is hard-coded for a 15-bit exponent with the usual
+ * bias.  s_logl.c and e_hypotl have less hard-coding but end up requiring
+ * the same for the exponent and more for the mantissa.
  */
 #if LDBL_MAX_EXP != 0x4000
 #error "Unsupported long double format"
-#else
-#define	THRESH	0x1p16382L
 #endif
 
+/* For avoiding overflow for components >= LDBL_MAX / (1 + sqrt(2)). */
+#define	THRESH	(LDBL_MAX / 2.414213562373095048801688724209698L)
+
 long double complex
 csqrtl(long double complex z)
 {