[PATCH] libm -- fix expl kernels

[Steve Kargl]

The following patch fixes the ld80 and ld128 kernels for expl(x).
I can't test ld128, so that patch is untested.

. Micro-optimization: use sincosl(x) instead of a call to cosl(x) and
  a call to sinl(x).  Argument reduction is done once not twice.
. Use a long double constant instead of an invalid double constant.
. Spell scale2 correctly


--- /usr/src/lib/msun/ld80/k_expl.h	2019-02-23 07:45:04.364785000 -0800
+++ ld80/k_expl.h	2019-12-26 08:15:29.366146000 -0800
@@ -277,7 +277,7 @@
 static inline long double complex
 __ldexp_cexpl(long double complex z, int expt)
 {
-	long double exp_x, hi, lo;
+	long double c, exp_x, hi, lo, s;
 	long double x, y, scale1, scale2;
 	int half_expt, k;
 
@@ -285,16 +285,17 @@
 	y = cimagl(z);
 	__k_expl(x, &hi, &lo, &k);
 
-	exp_x = (lo + hi) * 0x1p16382;
+	exp_x = (lo + hi) * 0x1p16382L;
 	expt += k - 16382;
 
 	scale1 = 1;
 	half_expt = expt / 2;
 	SET_LDBL_EXPSIGN(scale1, BIAS + half_expt);
 	scale2 = 1;
-	SET_LDBL_EXPSIGN(scale1, BIAS + expt - half_expt);
+	SET_LDBL_EXPSIGN(scale2, BIAS + expt - half_expt);
 
-	return (CMPLXL(cos(y) * exp_x * scale1 * scale2,
-	    sinl(y) * exp_x * scale1 * scale2));
+	sincosl(y, &s, &c);
+	return (CMPLXL(c * exp_x * scale1 * scale2,
+	    s * exp_x * scale1 * scale2));
 }
 #endif /* _COMPLEX_H */
--- /usr/src/lib/msun/ld128/k_expl.h	2020-03-07 09:48:15.719876000 -0800
+++ ld128/k_expl.h	2020-09-19 16:51:43.577195000 -0700
@@ -300,7 +300,7 @@
 static inline long double complex
 __ldexp_cexpl(long double complex z, int expt)
 {
-	long double exp_x, hi, lo;
+	long double c, exp_x, hi, lo, s;
 	long double x, y, scale1, scale2;
 	int half_expt, k;
 
@@ -308,16 +308,17 @@
 	y = cimagl(z);
 	__k_expl(x, &hi, &lo, &k);
 
-	exp_x = (lo + hi) * 0x1p16382;
+	exp_x = (lo + hi) * 0x1p16382L;
 	expt += k - 16382;
 
 	scale1 = 1;
 	half_expt = expt / 2;
 	SET_LDBL_EXPSIGN(scale1, BIAS + half_expt);
 	scale2 = 1;
-	SET_LDBL_EXPSIGN(scale1, BIAS + expt - half_expt);
+	SET_LDBL_EXPSIGN(scale2, BIAS + expt - half_expt);
 
-	return (CMPLXL(cos(y) * exp_x * scale1 * scale2,
-	    sinl(y) * exp_x * scale1 * scale2));
+	sincosl(y, &s, &c);
+	return (CMPLXL(c * exp_x * scale1 * scale2,
+	    s * exp_x * scale1 * scale2));
 }
 #endif /* _COMPLEX_H */
-- 
Steve