math/gen: fix fma fenv flags on nan

[libc-test] / src / math / gen / mp.c

diff --git a/src/math/gen/mp.c b/src/math/gen/mp.c
index cb10c41..172b454 100644 (file)
--- a/src/math/gen/mp.c
+++ b/src/math/gen/mp.c
@@ -123,7 +123,7 @@ static void genf(struct t *p, mpfr_t my, int t, int r)
 
        t = adjust(mr, my, t, r);
        p->y = mpfr_get_flt(mr, r);
-       p->e = eflags(isnan(p->x) || isnan(p->x2));
+       p->e = eflags(isnan(p->x) || isnan(p->x2) || isnan(p->x3));
        i = eulpf(p->y);
        if (!isfinite(p->y)) {
                p->dy = 0;
@@ -193,7 +193,7 @@ static void gend(struct t *p, mpfr_t my, int t, int r)
 
        t = adjust(mr, my, t, r);
        p->y = mpfr_get_d(mr, r);
-       p->e = eflags(isnan(p->x) || isnan(p->x2));
+       p->e = eflags(isnan(p->x) || isnan(p->x2) || isnan(p->x3));
        i = eulp(p->y);
        if (!isfinite(p->y)) {
                p->dy = 0;
@@ -268,7 +268,7 @@ static void genl(struct t *p, mpfr_t my, int t, int r)
 
        t = adjust(mr, my, t, r);
        p->y = mpfr_get_ld(mr, r);
-       p->e = eflags(isnan(p->x) || isnan(p->x2));
+       p->e = eflags(isnan(p->x) || isnan(p->x2) || isnan(p->x3));
        i = eulpl(p->y);
        if (!isfinite(p->y)) {
                p->dy = 0;
@@ -355,6 +355,15 @@ static int wrap_remquo(mpfr_t my, const mpfr_t mx, const mpfr_t mx2, mpfr_rnd_t
 {
        return mpfr_remquo(my, &mpremquo_q, mx, mx2, r);
 }
+static int mpbessel_n;
+static int wrap_jn(mpfr_t my, const mpfr_t mx, mpfr_rnd_t r)
+{
+       return mpfr_jn(my, mpbessel_n, mx, r);
+}
+static int wrap_yn(mpfr_t my, const mpfr_t mx, mpfr_rnd_t r)
+{
+       return mpfr_yn(my, mpbessel_n, mx, r);
+}
 static int wrap_ceil(mpfr_t my, const mpfr_t mx, mpfr_rnd_t r)
 {
        return mpfr_ceil(my, mx);
@@ -473,30 +482,69 @@ int mplogb(struct t *t)
 {
        MPFR_DECL_INIT(mx, 53);
 
-       mpfr_set_d(mx, t->x, MPFR_RNDN);
-       t->y = mpfr_get_exp(mx) - 1;
        t->dy = 0;
        t->e = 0;
+       if (t->x == 0) {
+               t->y = -INFINITY;
+               t->e |= DIVBYZERO;
+               return 0;
+       }
+       if (isinf(t->x)) {
+               t->y = INFINITY;
+               return 0;
+       }
+       if (isnan(t->x)) {
+               t->y = t->x;
+               return 0;
+       }
+       mpfr_set_d(mx, t->x, MPFR_RNDN);
+       t->y = mpfr_get_exp(mx) - 1;
        return 0;
 }
 int mplogbf(struct t *t)
 {
        MPFR_DECL_INIT(mx, 24);
 
-       mpfr_set_flt(mx, t->x, MPFR_RNDN);
-       t->y = mpfr_get_exp(mx) - 1;
        t->dy = 0;
        t->e = 0;
+       if (t->x == 0) {
+               t->y = -INFINITY;
+               t->e |= DIVBYZERO;
+               return 0;
+       }
+       if (isinf(t->x)) {
+               t->y = INFINITY;
+               return 0;
+       }
+       if (isnan(t->x)) {
+               t->y = t->x;
+               return 0;
+       }
+       mpfr_set_flt(mx, t->x, MPFR_RNDN);
+       t->y = mpfr_get_exp(mx) - 1;
        return 0;
 }
 int mplogbl(struct t *t)
 {
        MPFR_DECL_INIT(mx, 64);
 
-       mpfr_set_ld(mx, t->x, MPFR_RNDN);
-       t->y = mpfr_get_exp(mx) - 1;
        t->dy = 0;
        t->e = 0;
+       if (t->x == 0) {
+               t->y = -INFINITY;
+               t->e |= DIVBYZERO;
+               return 0;
+       }
+       if (isinf(t->x)) {
+               t->y = INFINITY;
+               return 0;
+       }
+       if (isnan(t->x)) {
+               t->y = t->x;
+               return 0;
+       }
+       mpfr_set_ld(mx, t->x, MPFR_RNDN);
+       t->y = mpfr_get_exp(mx) - 1;
        return 0;
 }
 int mpnearbyint(struct t *t) { return mpd1(t, wrap_nearbyint) || (t->e&=~INEXACT, 0); }
@@ -571,6 +619,7 @@ int mptanl(struct t *t) { return mpl1(t, mpfr_tan); }
 int mptanh(struct t *t) { return mpd1(t, mpfr_tanh); }
 int mptanhf(struct t *t) { return mpf1(t, mpfr_tanh); }
 int mptanhl(struct t *t) { return mpl1(t, mpfr_tanh); }
+// TODO: tgamma(2) raises wrong flags
 int mptgamma(struct t *t) { return mpd1(t, mpfr_gamma); }
 int mptgammaf(struct t *t) { return mpf1(t, mpfr_gamma); }
 int mptgammal(struct t *t) { return mpl1(t, mpfr_gamma); }
@@ -739,6 +788,8 @@ int mpilogb(struct t *t)
        mpfr_set_d(mx, t->x, MPFR_RNDN);
        t->i = mpfr_get_exp(mx) - 1;
        t->e = 0;
+       if (isinf(t->x) || isnan(t->x) || t->x == 0)
+               t->e = INVALID;
        return 0;
 }
 int mpilogbf(struct t *t)
@@ -748,6 +799,8 @@ int mpilogbf(struct t *t)
        mpfr_set_flt(mx, t->x, MPFR_RNDN);
        t->i = mpfr_get_exp(mx) - 1;
        t->e = 0;
+       if (isinf(t->x) || isnan(t->x) || t->x == 0)
+               t->e = INVALID;
        return 0;
 }
 int mpilogbl(struct t *t)
@@ -757,6 +810,8 @@ int mpilogbl(struct t *t)
        mpfr_set_ld(mx, t->x, MPFR_RNDN);
        t->i = mpfr_get_exp(mx) - 1;
        t->e = 0;
+       if (isinf(t->x) || isnan(t->x) || t->x == 0)
+               t->e = INVALID;
        return 0;
 }
 
@@ -792,7 +847,7 @@ int mpmodf(struct t *t)
                return r;
        t->y2 = t->y;
        t->dy2 = t->dy;
-       e = t->e;
+       e = t->e & ~INEXACT;
        r = mpd1(t, mpfr_frac);
        t->e |= e;
        return r;
@@ -807,7 +862,7 @@ int mpmodff(struct t *t)
                return r;
        t->y2 = t->y;
        t->dy2 = t->dy;
-       e = t->e;
+       e = t->e & ~INEXACT;
        r = mpf1(t, mpfr_frac);
        t->e |= e;
        return r;
@@ -822,7 +877,7 @@ int mpmodfl(struct t *t)
                return r;
        t->y2 = t->y;
        t->dy2 = t->dy;
-       e = t->e;
+       e = t->e & ~INEXACT;
        r = mpl1(t, mpfr_frac);
        t->e |= e;
        return r;
@@ -961,3 +1016,10 @@ int mpfmal(struct t *t)
 #endif
 }
 
+int mpjn(struct t *t) { mpbessel_n = t->i; return mpd1(t, wrap_jn); }
+int mpjnf(struct t *t) { mpbessel_n = t->i; return mpf1(t, wrap_jn); }
+int mpjnl(struct t *t) { mpbessel_n = t->i; return mpl1(t, wrap_jn); }
+int mpyn(struct t *t) { mpbessel_n = t->i; return mpd1(t, wrap_yn); }
+int mpynf(struct t *t) { mpbessel_n = t->i; return mpf1(t, wrap_yn); }
+int mpynl(struct t *t) { mpbessel_n = t->i; return mpl1(t, wrap_yn); }
+