libc/math: Migrate to upstream musl, take 2

Migrate more libc/math files to use upstream musl. These ones have only trivial
type changes in osv.git:

diff -ur libc/math/acos.c musl/src/math/acos.c
--- libc/math/acos.c    2013-10-25 12:11:23.991751105 +0300
+++ musl/src/math/acos.c        2014-08-04 17:07:25.449850385 +0300
@@ -51,7 +51,7 @@

 static double R(double z)
 {
-       double p, q;
+       double_t p, q;
        p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
        q = 1.0+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
        return p/q;
diff -ur libc/math/acosf.c musl/src/math/acosf.c
--- libc/math/acosf.c   2013-10-25 12:11:24.009751042 +0300
+++ musl/src/math/acosf.c       2014-08-04 17:07:25.449850385 +0300
@@ -25,7 +25,7 @@

 static float R(float z)
 {
-       float p, q;
+       float_t p, q;
        p = z*(pS0+z*(pS1+z*pS2));
        q = 1.0f+z*qS1;
        return p/q;
diff -ur libc/math/asin.c musl/src/math/asin.c
--- libc/math/asin.c    2013-10-25 12:11:23.999751077 +0300
+++ musl/src/math/asin.c        2014-08-04 17:07:25.450850382 +0300
@@ -58,7 +58,7 @@

 static double R(double z)
 {
-       double p, q;
+       double_t p, q;
        p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
        q = 1.0+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
        return p/q;
diff -ur libc/math/asinf.c musl/src/math/asinf.c
--- libc/math/asinf.c   2013-10-25 12:11:24.025750986 +0300
+++ musl/src/math/asinf.c       2014-08-04 17:07:25.450850382 +0300
@@ -26,7 +26,7 @@

 static float R(float z)
 {
-       float p, q;
+       float_t p, q;
        p = z*(pS0+z*(pS1+z*pS2));
        q = 1.0f+z*qS1;
        return p/q;
diff -ur libc/math/atan.c musl/src/math/atan.c
--- libc/math/atan.c    2013-10-25 12:11:24.006751053 +0300
+++ musl/src/math/atan.c        2014-08-04 17:07:25.450850382 +0300
@@ -62,7 +62,7 @@

 double atan(double x)
 {
-       double w,s1,s2,z;
+       double_t w,s1,s2,z;
        uint32_t ix,sign;
        int id;
diff -ur libc/math/atanf.c musl/src/math/atanf.c
--- libc/math/atanf.c   2013-10-25 12:11:24.028750976 +0300
+++ musl/src/math/atanf.c       2014-08-04 17:07:25.450850382 +0300
@@ -40,7 +40,7 @@

 float atanf(float x)
 {
-       float w,s1,s2,z;
+       float_t w,s1,s2,z;
        uint32_t ix,sign;
        int id;
diff -ur libc/math/__cos.c musl/src/math/__cos.c
--- libc/math/__cos.c   2014-08-04 15:40:10.209333688 +0300
+++ musl/src/math/__cos.c       2014-08-04 17:07:25.448850388 +0300
@@ -60,7 +60,7 @@

 double __cos(double x, double y)
 {
-       double hz,z,r,w;
+       double_t hz,z,r,w;

        z  = x*x;
        w  = z*z;
diff -ur libc/math/__cosdf.c musl/src/math/__cosdf.c
--- libc/math/__cosdf.c 2014-08-04 15:40:10.209333688 +0300
+++ musl/src/math/__cosdf.c     2014-08-04 17:07:25.448850388 +0300
@@ -25,7 +25,7 @@

 float __cosdf(double x)
 {
-       double r, w, z;
+       double_t r, w, z;

        /* Try to optimize for parallel evaluation as in __tandf.c. */
        z = x*x;
diff -ur libc/math/erf.c musl/src/math/erf.c
--- libc/math/erf.c     2013-10-25 12:11:24.026750982 +0300
+++ musl/src/math/erf.c 2014-08-04 17:07:25.451850379 +0300
@@ -176,7 +176,7 @@

 static double erfc1(double x)
 {
-       double s,P,Q;
+       double_t s,P,Q;

        s = fabs(x) - 1;
        P = pa0+s*(pa1+s*(pa2+s*(pa3+s*(pa4+s*(pa5+s*pa6)))));
@@ -186,7 +186,8 @@

 static double erfc2(uint32_t ix, double x)
 {
-       double s,z,R,S;
+       double_t s,R,S;
+       double z;

        if (ix < 0x3ff40000)  /* |x| < 1.25 */
                return erfc1(x);
diff -ur libc/math/erff.c musl/src/math/erff.c
--- libc/math/erff.c    2013-10-25 12:11:24.003751062 +0300
+++ musl/src/math/erff.c        2014-08-04 17:07:25.451850379 +0300
@@ -86,7 +86,7 @@

 static float erfc1(float x)
 {
-       float s,P,Q;
+       float_t s,P,Q;

        s = fabsf(x) - 1;
        P = pa0+s*(pa1+s*(pa2+s*(pa3+s*(pa4+s*(pa5+s*pa6)))));
@@ -96,7 +96,8 @@

 static float erfc2(uint32_t ix, float x)
 {
-       float s,z,R,S;
+       float_t s,R,S;
+       float z;

        if (ix < 0x3fa00000)  /* |x| < 1.25 */
                return erfc1(x);
diff -ur libc/math/j0.c musl/src/math/j0.c
--- libc/math/j0.c      2013-10-25 12:11:24.017751014 +0300
+++ musl/src/math/j0.c  2014-08-04 17:07:25.456850364 +0300
@@ -263,7 +263,7 @@
 static double pzero(double x)
 {
        const double *p,*q;
-       double z,r,s;
+       double_t z,r,s;
        uint32_t ix;

        GET_HIGH_WORD(ix, x);
@@ -359,7 +359,7 @@
 static double qzero(double x)
 {
        const double *p,*q;
-       double s,r,z;
+       double_t s,r,z;
        uint32_t ix;

        GET_HIGH_WORD(ix, x);
diff -ur libc/math/j0f.c musl/src/math/j0f.c
--- libc/math/j0f.c     2013-10-25 12:11:23.996751087 +0300
+++ musl/src/math/j0f.c 2014-08-04 17:07:25.456850364 +0300
@@ -201,7 +201,7 @@
 static float pzerof(float x)
 {
        const float *p,*q;
-       float z,r,s;
+       float_t z,r,s;
        uint32_t ix;

        GET_FLOAT_WORD(ix, x);
@@ -297,7 +297,7 @@
 static float qzerof(float x)
 {
        const float *p,*q;
-       float s,r,z;
+       float_t s,r,z;
        uint32_t ix;

        GET_FLOAT_WORD(ix, x);
diff -ur libc/math/j1.c musl/src/math/j1.c
--- libc/math/j1.c      2013-10-25 12:11:23.989751112 +0300
+++ musl/src/math/j1.c  2014-08-04 17:07:25.456850364 +0300
@@ -250,7 +250,7 @@
 static double pone(double x)
 {
        const double *p,*q;
-       double z,r,s;
+       double_t z,r,s;
        uint32_t ix;

        GET_HIGH_WORD(ix, x);
@@ -346,7 +346,7 @@
 static double qone(double x)
 {
        const double *p,*q;
-       double  s,r,z;
+       double_t s,r,z;
        uint32_t ix;

        GET_HIGH_WORD(ix, x);
diff -ur libc/math/j1f.c musl/src/math/j1f.c
--- libc/math/j1f.c     2013-10-25 12:11:24.000751074 +0300
+++ musl/src/math/j1f.c 2014-08-04 17:07:25.456850364 +0300
@@ -198,7 +198,7 @@
 static float ponef(float x)
 {
        const float *p,*q;
-       float z,r,s;
+       float_t z,r,s;
        uint32_t ix;

        GET_FLOAT_WORD(ix, x);
@@ -294,7 +294,7 @@
 static float qonef(float x)
 {
        const float *p,*q;
-       float s,r,z;
+       float_t s,r,z;
        uint32_t ix;

        GET_FLOAT_WORD(ix, x);
diff -ur libc/math/__sin.c musl/src/math/__sin.c
--- libc/math/__sin.c   2014-08-04 15:40:10.211333681 +0300
+++ musl/src/math/__sin.c       2014-08-04 17:07:25.449850385 +0300
@@ -51,7 +51,7 @@

 double __sin(double x, double y, int iy)
 {
-       double z,r,v,w;
+       double_t z,r,v,w;

        z = x*x;
        w = z*z;
diff -ur libc/math/__sindf.c musl/src/math/__sindf.c
--- libc/math/__sindf.c 2014-08-04 15:40:10.211333681 +0300
+++ musl/src/math/__sindf.c     2014-08-04 17:07:25.449850385 +0300
@@ -25,7 +25,7 @@

 float __sindf(double x)
 {
-       double r, s, w, z;
+       double_t r, s, w, z;

        /* Try to optimize for parallel evaluation as in __tandf.c. */
        z = x*x;
diff -ur libc/math/tgamma.c musl/src/math/tgamma.c
--- libc/math/tgamma.c  2013-10-25 12:11:24.008751046 +0300
+++ musl/src/math/tgamma.c      2014-08-04 17:07:25.461850349 +0300
@@ -89,7 +89,7 @@
 /* S(x) rational function for positive x */
 static double S(double x)
 {
-       double num = 0, den = 0;
+       double_t num = 0, den = 0;
        int i;

        /* to avoid overflow handle large x differently */

Reviewed-by: Nadav Har'El <nyh@cloudius-systems.com>
Signed-off-by: Pekka Enberg <penberg@cloudius-systems.com>
libc/math: Migrate to upstream musl, take 3

Migrate more files with trivial diff to osv.git. Note that the files are
not even compiled in to OSv image so the change itself is anyway safe:

diff -ur libc/math/fmal.c musl/src/math/fmal.c
--- libc/math/fmal.c    2014-08-04 15:40:10.214333672 +0300
+++ musl/src/math/fmal.c        2014-08-04 17:07:25.453850373 +0300
@@ -162,6 +162,7 @@
  */
 long double fmal(long double x, long double y, long double z)
 {
+       #pragma STDC FENV_ACCESS ON
        long double xs, ys, zs, adj;
        struct dd xy, r;
        int oround;
@@ -261,6 +262,8 @@
                /*
                 * There is no need to worry about double rounding in directed
                 * rounding modes.
+                * TODO: underflow is not raised correctly, example in downward rounding:
+                * fmal(0x1.0000000001p-16000L, 0x1.0000000001p-400L, -0x1p-16440L)
                 */
                fesetround(oround);
                adj = r.lo + xy.lo;
diff -ur libc/math/llrintl.c musl/src/math/llrintl.c
--- libc/math/llrintl.c 2014-08-04 15:40:10.215333669 +0300
+++ musl/src/math/llrintl.c     2014-08-04 17:07:25.457850361 +0300
@@ -18,6 +18,7 @@
 */
 long long llrintl(long double x)
 {
+       #pragma STDC FENV_ACCESS ON
        int e;

        e = fetestexcept(FE_INEXACT);
diff -ur libc/math/lrintl.c musl/src/math/lrintl.c
--- libc/math/lrintl.c  2014-08-04 15:40:10.216333665 +0300
+++ musl/src/math/lrintl.c      2014-08-04 17:07:25.458850358 +0300
@@ -18,6 +18,7 @@
 */
 long lrintl(long double x)
 {
+       #pragma STDC FENV_ACCESS ON
        int e;

        e = fetestexcept(FE_INEXACT);
diff -ur libc/math/nearbyint.c musl/src/math/nearbyint.c
--- libc/math/nearbyint.c       2014-08-04 15:40:10.217333662 +0300
+++ musl/src/math/nearbyint.c   2014-08-04 17:07:25.459850355 +0300
@@ -6,6 +6,7 @@
 double nearbyint(double x)
 {
 #ifdef FE_INEXACT
+       #pragma STDC FENV_ACCESS ON
        int e;

        e = fetestexcept(FE_INEXACT);
diff -ur libc/math/nearbyintf.c musl/src/math/nearbyintf.c
--- libc/math/nearbyintf.c      2014-08-04 15:40:10.217333662 +0300
+++ musl/src/math/nearbyintf.c  2014-08-04 17:07:25.459850355 +0300
@@ -4,6 +4,7 @@
 float nearbyintf(float x)
 {
 #ifdef FE_INEXACT
+       #pragma STDC FENV_ACCESS ON
        int e;

        e = fetestexcept(FE_INEXACT);
diff -ur libc/math/nearbyintl.c musl/src/math/nearbyintl.c
--- libc/math/nearbyintl.c      2014-08-04 15:40:10.217333662 +0300
+++ musl/src/math/nearbyintl.c  2014-08-04 17:07:25.459850355 +0300
@@ -11,6 +11,7 @@
 long double nearbyintl(long double x)
 {
 #ifdef FE_INEXACT
+       #pragma STDC FENV_ACCESS ON
        int e;

        e = fetestexcept(FE_INEXACT);

Reviewed-by: Nadav Har'El <nyh@cloudius-systems.com>
Signed-off-by: Pekka Enberg <penberg@cloudius-systems.com>