diff -Naurd mpfr-2.2.1-p3/VERSION mpfr-2.2.1-p4/VERSION --- mpfr-2.2.1-p3/VERSION 2007-02-15 23:40:49.000000000 +0000 +++ mpfr-2.2.1-p4/VERSION 2007-02-15 23:43:29.000000000 +0000 @@ -1 +1 @@ -2.2.1-p3 +2.2.1-p4 diff -Naurd mpfr-2.2.1-p3/mpfr.h mpfr-2.2.1-p4/mpfr.h --- mpfr-2.2.1-p3/mpfr.h 2007-02-15 23:40:49.000000000 +0000 +++ mpfr-2.2.1-p4/mpfr.h 2007-02-15 23:43:29.000000000 +0000 @@ -26,7 +26,7 @@ #define MPFR_VERSION_MAJOR 2 #define MPFR_VERSION_MINOR 2 #define MPFR_VERSION_PATCHLEVEL 1 -#define MPFR_VERSION_STRING "2.2.1-p3" +#define MPFR_VERSION_STRING "2.2.1-p4" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-2.2.1-p3/pow_si.c mpfr-2.2.1-p4/pow_si.c --- mpfr-2.2.1-p3/pow_si.c 2006-11-29 09:49:47.000000000 +0000 +++ mpfr-2.2.1-p4/pow_si.c 2007-02-15 23:43:17.000000000 +0000 @@ -67,24 +67,20 @@ /* detect exact powers: x^(-n) is exact iff x is a power of 2 */ if (mpfr_cmp_si_2exp (x, MPFR_SIGN(x), MPFR_EXP(x) - 1) == 0) { - mp_exp_t expx = MPFR_EXP (x); /* warning: x and y may be the same - variable */ - mpfr_set_si (y, (n % 2) ? MPFR_INT_SIGN(x) : 1, rnd); - expx --; + mp_exp_t expx = MPFR_EXP (x) - 1, expy; MPFR_ASSERTD (n < 0); - /* Warning n*expx may overflow! - Some systems abort with LONG_MIN / 1 or LONG_MIN/-1*/ - if (n != -1 && expx > 0 && -expx < MPFR_EXP_MIN / (-n)) - MPFR_EXP (y) = MPFR_EMIN_MIN - 1; /* Underflow */ - else if (n != -1 && expx < 0 && -expx > MPFR_EXP_MAX / (-n)) - MPFR_EXP (y) = MPFR_EMAX_MAX + 1; /* Overflow */ - else - MPFR_EXP (y) += n * expx; - return mpfr_check_range (y, 0, rnd); + /* Warning: n * expx may overflow! + Some systems (apparently alpha-freebsd) abort with + LONG_MIN / 1, and LONG_MIN / -1 is undefined. */ + expy = + n != -1 && expx > 0 && expx > (__gmpfr_emin - 1) / n ? + MPFR_EMIN_MIN - 2 /* Underflow */ : + n != -1 && expx < 0 && expx < (__gmpfr_emax - 1) / n ? + MPFR_EMAX_MAX /* Overflow */ : n * expx; + return mpfr_set_si_2exp (y, n % 2 ? MPFR_INT_SIGN (x) : 1, + expy, rnd); } - n = -n; - /* General case */ { /* Declaration of the intermediary variable */ @@ -94,9 +90,12 @@ mp_prec_t Nt; /* working precision */ mp_exp_t err; /* error */ int inexact; + unsigned long abs_n; MPFR_SAVE_EXPO_DECL (expo); MPFR_ZIV_DECL (loop); + abs_n = - (unsigned long) n; + /* compute the precision of intermediary variable */ /* the optimal number of bits : see algorithms.tex */ Nt = Ny + 3 + MPFR_INT_CEIL_LOG2 (Ny); @@ -109,17 +108,17 @@ MPFR_ZIV_INIT (loop, Nt); for (;;) { - /* compute 1/(x^n) n>0*/ - mpfr_pow_ui (t, x, (unsigned long int) n, GMP_RNDN); + /* compute 1/(x^n), with n > 0 */ + mpfr_pow_ui (t, x, abs_n, GMP_RNDN); mpfr_ui_div (t, 1, t, GMP_RNDN); - /* FIXME: old code improved, but I think this is still incorrect. */ + /* FIXME: old code improved, but I think this is still incorrect. */ if (MPFR_UNLIKELY (MPFR_IS_ZERO (t))) { MPFR_ZIV_FREE (loop); mpfr_clear (t); MPFR_SAVE_EXPO_FREE (expo); return mpfr_underflow (y, rnd == GMP_RNDN ? GMP_RNDZ : rnd, - (unsigned) n & 1 ? MPFR_SIGN (x) : + abs_n & 1 ? MPFR_SIGN (x) : MPFR_SIGN_POS); } if (MPFR_UNLIKELY (MPFR_IS_INF (t))) @@ -127,8 +126,7 @@ MPFR_ZIV_FREE (loop); mpfr_clear (t); MPFR_SAVE_EXPO_FREE (expo); - return mpfr_overflow (y, rnd, - (unsigned) n & 1 ? MPFR_SIGN (x) : + return mpfr_overflow (y, rnd, abs_n & 1 ? MPFR_SIGN (x) : MPFR_SIGN_POS); } /* error estimate -- see pow function in algorithms.ps */ diff -Naurd mpfr-2.2.1-p3/tests/tpow.c mpfr-2.2.1-p4/tests/tpow.c --- mpfr-2.2.1-p3/tests/tpow.c 2006-11-29 09:49:47.000000000 +0000 +++ mpfr-2.2.1-p4/tests/tpow.c 2007-02-15 23:43:17.000000000 +0000 @@ -23,6 +23,7 @@ #include #include #include +#include #include "mpfr-test.h" @@ -231,6 +232,55 @@ mpfr_pow_si (x, x, -2, GMP_RNDN); MPFR_ASSERTN(mpfr_inf_p (x) && mpfr_sgn (x) > 0); + mpfr_set_si (x, 2, GMP_RNDN); + mpfr_pow_si (x, x, LONG_MAX, GMP_RNDN); /* 2^LONG_MAX */ + if (LONG_MAX > mpfr_get_emax () - 1) /* LONG_MAX + 1 > emax */ + { + MPFR_ASSERTN (mpfr_inf_p (x)); + } + else + { + MPFR_ASSERTN (mpfr_cmp_si_2exp (x, 1, LONG_MAX)); + } + + mpfr_set_si (x, 2, GMP_RNDN); + mpfr_pow_si (x, x, LONG_MIN, GMP_RNDN); /* 2^LONG_MIN */ + if (LONG_MIN + 1 < mpfr_get_emin ()) + { + MPFR_ASSERTN (mpfr_zero_p (x)); + } + else + { + MPFR_ASSERTN (mpfr_cmp_si_2exp (x, 1, LONG_MIN)); + } + + mpfr_set_si (x, 2, GMP_RNDN); + mpfr_pow_si (x, x, LONG_MIN + 1, GMP_RNDN); /* 2^(LONG_MIN+1) */ + if (mpfr_nan_p (x)) + { + printf ("Error in pow_si(2, LONG_MIN+1): got NaN\n"); + exit (1); + } + if (LONG_MIN + 2 < mpfr_get_emin ()) + { + MPFR_ASSERTN (mpfr_zero_p (x)); + } + else + { + MPFR_ASSERTN (mpfr_cmp_si_2exp (x, 1, LONG_MIN + 1)); + } + + mpfr_set_si_2exp (x, 1, -1, GMP_RNDN); /* 0.5 */ + mpfr_pow_si (x, x, LONG_MIN, GMP_RNDN); /* 2^(-LONG_MIN) */ + if (LONG_MIN < 1 - mpfr_get_emax ()) /* 1 - LONG_MIN > emax */ + { + MPFR_ASSERTN (mpfr_inf_p (x)); + } + else + { + MPFR_ASSERTN (mpfr_cmp_si_2exp (x, 2, - (LONG_MIN + 1))); + } + mpfr_clear (x); } @@ -270,6 +320,35 @@ } static void +pow_si_long_min (void) +{ + mpfr_t x, y, z; + int inex; + + mpfr_inits2 (sizeof(long) * CHAR_BIT + 32, x, y, z, (void *) 0); + mpfr_set_si_2exp (x, 3, -1, GMP_RNDN); /* 1.5 */ + + inex = mpfr_set_si (y, LONG_MIN, GMP_RNDN); + MPFR_ASSERTN (inex == 0); + mpfr_nextbelow (y); + mpfr_pow (y, x, y, GMP_RNDD); + + inex = mpfr_set_si (z, LONG_MIN, GMP_RNDN); + MPFR_ASSERTN (inex == 0); + mpfr_nextabove (z); + mpfr_pow (z, x, z, GMP_RNDU); + + mpfr_pow_si (x, x, LONG_MIN, GMP_RNDN); /* 1.5^LONG_MIN */ + if (mpfr_cmp (x, y) < 0 || mpfr_cmp (x, z) > 0) + { + printf ("Error in pow_si_long_min\n"); + exit (1); + } + + mpfr_clears (x, y, z, (void *) 0); +} + +static void check_inexact (mp_prec_t p) { mpfr_t x, y, z, t; @@ -751,6 +830,7 @@ check_pow_ui (); check_pow_si (); check_special_pow_si (); + pow_si_long_min (); for (p = 2; p < 100; p++) check_inexact (p); underflows (); diff -Naurd mpfr-2.2.1-p3/version.c mpfr-2.2.1-p4/version.c --- mpfr-2.2.1-p3/version.c 2007-02-15 23:40:49.000000000 +0000 +++ mpfr-2.2.1-p4/version.c 2007-02-15 23:43:29.000000000 +0000 @@ -24,5 +24,5 @@ const char * mpfr_get_version (void) { - return "2.2.1-p3"; + return "2.2.1-p4"; }