Update to FreeBSD libm r336665.

This reverts commit 253a8306316cedfd6fd3e3a169fbffe4cac04035 and moves
us forward to a revision that contains fixes for the problem with the
previous attempt.

This also makes sincos(3)/sincosf(3)/sincosl(3) available to `_BSD_SOURCE`
as well as `_GNU_SOURCE`.

The new FreeBSD libm code requires the FreeBSD `__CONCAT` macro, and all
our existing callers are FreeBSD too, so update that.

There's also an assumption that <complex.h> drags in <math.h> which isn't
true for us, so work around that with `-include` in the makefile. This
then causes clang to recognize a bug -- returning from a void function --
in our fake (LP32) sincosl(3), so fix that too.

Bug: http://b/111710419
Change-Id: I84703ad844f8afde6ec6b11604ab3c096ccb62c3
Test: ran tests
diff --git a/libm/upstream-freebsd/lib/msun/src/s_csqrtf.c b/libm/upstream-freebsd/lib/msun/src/s_csqrtf.c
index 195e3b7..16d3266 100644
--- a/libm/upstream-freebsd/lib/msun/src/s_csqrtf.c
+++ b/libm/upstream-freebsd/lib/msun/src/s_csqrtf.c
@@ -27,27 +27,21 @@
  */
 
 #include <sys/cdefs.h>
-__FBSDID("$FreeBSD: head/lib/msun/src/s_csqrtf.c 326219 2017-11-26 02:00:33Z pfg $");
+__FBSDID("$FreeBSD: head/lib/msun/src/s_csqrtf.c 336412 2018-07-17 12:01:59Z bde $");
 
 #include <complex.h>
 #include <math.h>
 
 #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)
@@ -56,7 +50,7 @@
 		return (CMPLXF(INFINITY, b));
 	if (isnan(a)) {
 		t = (b - b) / (b - b);	/* raise invalid if b is not a NaN */
-		return (CMPLXF(a, t));	/* return NaN + NaN i */
+		return (CMPLXF(a + 0.0L + t, a + 0.0L + t)); /* NaN + NaN i */
 	}
 	if (isinf(a)) {
 		/*
@@ -70,10 +64,10 @@
 		else
 			return (CMPLXF(a, copysignf(b - b, b)));
 	}
-	/*
-	 * The remaining special case (b is NaN) is handled just fine by
-	 * the normal code path below.
-	 */
+	if (isnan(b)) {
+		t = (a - a) / (a - a);	/* raise invalid */
+		return (CMPLXF(b + 0.0L + t, b + 0.0L + t)); /* NaN + NaN i */
+	}
 
 	/*
 	 * We compute t in double precision to avoid overflow and to
@@ -82,9 +76,9 @@
 	 */
 	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)));
 	}
 }