diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2011-10-05 21:39:57.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2011-10-05 21:39:57.000000000 +0000 @@ -0,0 +1 @@ +mpfr_unlikely diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2011-10-03 08:17:15.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2011-10-05 21:39:57.000000000 +0000 @@ -1 +1 @@ -3.1.0 +3.1.0-p1 diff -Naurd mpfr-3.1.0-a/src/mpfr-impl.h mpfr-3.1.0-b/src/mpfr-impl.h --- mpfr-3.1.0-a/src/mpfr-impl.h 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr-impl.h 2011-10-05 21:39:57.000000000 +0000 @@ -988,10 +988,11 @@ ******************************************************/ /* Theses macros help the compiler to determine if a test is - * likely or unlikely. */ + likely or unlikely. The !! is necessary in case x is larger + than a long. */ #if __MPFR_GNUC(3,0) || __MPFR_ICC(8,1,0) # define MPFR_LIKELY(x) (__builtin_expect(!!(x),1)) -# define MPFR_UNLIKELY(x) (__builtin_expect((x),0)) +# define MPFR_UNLIKELY(x) (__builtin_expect(!!(x),0)) #else # define MPFR_LIKELY(x) (x) # define MPFR_UNLIKELY(x) (x) diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2011-10-05 21:39:57.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0" +#define MPFR_VERSION_STRING "3.1.0-p1" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2011-10-05 21:39:57.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0"; + return "3.1.0-p1"; } diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2011-10-14 10:43:32.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2011-10-14 10:43:32.000000000 +0000 @@ -0,0 +1 @@ +lib-search-path diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2011-10-05 21:39:57.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2011-10-14 10:43:32.000000000 +0000 @@ -1 +1 @@ -3.1.0-p1 +3.1.0-p2 diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2011-10-05 21:39:57.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2011-10-14 10:43:32.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0-p1" +#define MPFR_VERSION_STRING "3.1.0-p2" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2011-10-05 21:39:57.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2011-10-14 10:43:32.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0-p1"; + return "3.1.0-p2"; } diff -Naurd mpfr-3.1.0-a/tests/Makefile.am mpfr-3.1.0-b/tests/Makefile.am --- mpfr-3.1.0-a/tests/Makefile.am 2011-10-03 08:17:14.000000000 +0000 +++ mpfr-3.1.0-b/tests/Makefile.am 2011-10-03 08:17:14.000000000 +0000 @@ -65,8 +65,24 @@ TESTS = $(check_PROGRAMS) TESTS_ENVIRONMENT = MPFR_QUIET=1 $(VALGRIND) -# Option to prevent libtool from generating wrapper scripts for the tests. +# The -no-install option prevents libtool from generating wrapper scripts +# for the tests. # This is useful to easily run the test scripts under valgrind or gdb. # See discussion http://thread.gmane.org/gmane.comp.lib.gnulib.bugs/28033 # http://article.gmane.org/gmane.comp.lib.gnulib.bugs/28140 in particular. -AM_LDFLAGS = -no-install +# +# The -L$(top_builddir)/src/.libs option is necessary for some platforms, +# such as HP-UX, when --with-gmp or --with-gmp-lib is used and an old MPFR +# library is already installed in the corresponding lib directory: its +# purpose is to make sure that the local .libs comes first in the library +# search path (otherwise the tests are linked against the old MPFR library +# by the LINK command -- see the generated Makefile). See: +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00042.html +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00043.html +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00044.html +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00066.html +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00065.html +# and +# http://debbugs.gnu.org/cgi/bugreport.cgi?bug=9728 +# +AM_LDFLAGS = -no-install -L$(top_builddir)/src/.libs diff -Naurd mpfr-3.1.0-a/tests/Makefile.in mpfr-3.1.0-b/tests/Makefile.in --- mpfr-3.1.0-a/tests/Makefile.in 2011-10-03 08:17:35.000000000 +0000 +++ mpfr-3.1.0-b/tests/Makefile.in 2011-10-03 08:17:35.000000000 +0000 @@ -1124,11 +1124,27 @@ TESTS = $(check_PROGRAMS) TESTS_ENVIRONMENT = MPFR_QUIET=1 $(VALGRIND) -# Option to prevent libtool from generating wrapper scripts for the tests. +# The -no-install option prevents libtool from generating wrapper scripts +# for the tests. # This is useful to easily run the test scripts under valgrind or gdb. # See discussion http://thread.gmane.org/gmane.comp.lib.gnulib.bugs/28033 # http://article.gmane.org/gmane.comp.lib.gnulib.bugs/28140 in particular. -AM_LDFLAGS = -no-install +# +# The -L$(top_builddir)/src/.libs option is necessary for some platforms, +# such as HP-UX, when --with-gmp or --with-gmp-lib is used and an old MPFR +# library is already installed in the corresponding lib directory: its +# purpose is to make sure that the local .libs comes first in the library +# search path (otherwise the tests are linked against the old MPFR library +# by the LINK command -- see the generated Makefile). See: +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00042.html +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00043.html +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00044.html +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00066.html +# http://websympa.loria.fr/wwsympa/arc/mpfr/2011-10/msg00065.html +# and +# http://debbugs.gnu.org/cgi/bugreport.cgi?bug=9728 +# +AM_LDFLAGS = -no-install -L$(top_builddir)/src/.libs all: all-am .SUFFIXES: diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2011-11-03 15:15:11.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2011-11-03 15:15:11.000000000 +0000 @@ -0,0 +1 @@ +vasprintf diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2011-10-14 10:43:32.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2011-11-03 15:15:11.000000000 +0000 @@ -1 +1 @@ -3.1.0-p2 +3.1.0-p3 diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2011-10-14 10:43:32.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2011-11-03 15:15:11.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0-p2" +#define MPFR_VERSION_STRING "3.1.0-p3" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-3.1.0-a/src/vasprintf.c mpfr-3.1.0-b/src/vasprintf.c --- mpfr-3.1.0-a/src/vasprintf.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/vasprintf.c 2011-11-03 15:15:11.000000000 +0000 @@ -1178,7 +1178,7 @@ mpfr_exp_t exp; char * str; const int spec_g = (spec.spec == 'g' || spec.spec == 'G'); - const int keep_trailing_zeros = spec_g && spec.alt; + const int keep_trailing_zeros = !spec_g || spec.alt; /* WARNING: an empty precision field is forbidden (it means precision = 6 and it should have been changed to 6 before the function call) */ @@ -1356,7 +1356,7 @@ else /* 1 <= |p| */ { - size_t nsd; /* Number of significant digits */ + size_t str_len; /* Determine the position of the most significant decimal digit. */ exp = floor_log10 (p); @@ -1365,12 +1365,10 @@ /* P is too large to print all its integral part digits */ return -1; - np->ip_size = exp + 1; - - nsd = spec.prec + np->ip_size; if (dec_info == NULL) - { - str = mpfr_get_str (NULL, &exp, 10, nsd, p, spec.rnd_mode); + { /* this case occurs with mpfr_printf ("%.0RUf", x) with x=9.5 */ + str = + mpfr_get_str (NULL, &exp, 10, spec.prec+exp+1, p, spec.rnd_mode); register_string (np->sl, str); } else @@ -1379,81 +1377,60 @@ str = dec_info->str; } np->ip_ptr = MPFR_IS_NEG (p) ? ++str : str; /* skip sign */ + str_len = strlen (str); + + /* integral part */ + if (exp > str_len) + /* mpfr_get_str gives no trailing zero when p is rounded up to the next + power of 10 (p integer, so no fractional part) */ + { + np->ip_trailing_zeros = exp - str_len; + np->ip_size = str_len; + } + else + np->ip_size = exp; if (spec.group) /* thousands separator in integral part */ np->thousands_sep = MPFR_THOUSANDS_SEPARATOR; - if (nsd == 0 || (spec_g && !spec.alt)) - /* compute how much non-zero digits in integral and fractional - parts */ + /* fractional part */ + str += np->ip_size; + str_len -= np->ip_size; + if (!keep_trailing_zeros) + /* remove trailing zeros, if any */ { - size_t str_len; - str_len = strlen (str); /* note: the sign has been skipped */ - - if (exp > str_len) - /* mpfr_get_str doesn't give the trailing zeros when p is a - multiple of 10 (p integer, so no fractional part) */ - { - np->ip_trailing_zeros = exp - str_len; - np->ip_size = str_len; - if (spec.alt) - np->point = MPFR_DECIMAL_POINT; - } - else - /* str may contain some digits which are in fractional part */ + char *ptr = str + str_len - 1; /* pointer to the last digit of + str */ + while ((*ptr == '0') && (str_len != 0)) { - char *ptr; - - ptr = str + str_len - 1; /* points to the end of str */ - str_len -= np->ip_size; /* number of digits in fractional - part */ - - if (!keep_trailing_zeros) - /* remove trailing zeros, if any */ - { - while ((*ptr == '0') && (str_len != 0)) - { - --ptr; - --str_len; - } - } - - if (str_len > INT_MAX) - /* too many digits in fractional part */ - return -1; - - if (str_len != 0) - /* some digits in fractional part */ - { - np->point = MPFR_DECIMAL_POINT; - np->fp_ptr = str + np->ip_size; - np->fp_size = str_len; - } + --ptr; + --str_len; } } - else - /* spec.prec digits in fractional part */ + + if (str_len > 0) + /* some nonzero digits in fractional part */ { - if (np->ip_size == exp - 1) - /* the absolute value of the number has been rounded up to a power - of ten. - Insert an additional zero in integral part and put the rest of - them in fractional part. */ - np->ip_trailing_zeros = 1; + if (str_len > INT_MAX) + /* too many digits in fractional part */ + return -1; - if (spec.prec != 0) - { - MPFR_ASSERTD (np->ip_size + np->ip_trailing_zeros == exp); - MPFR_ASSERTD (np->ip_size + spec.prec == nsd); + np->point = MPFR_DECIMAL_POINT; + np->fp_ptr = str; + np->fp_size = str_len; + } - np->point = MPFR_DECIMAL_POINT; - np->fp_ptr = str + np->ip_size; - np->fp_size = spec.prec; - } - else if (spec.alt) - np->point = MPFR_DECIMAL_POINT; + if (keep_trailing_zeros && str_len < spec.prec) + /* add missing trailing zeros */ + { + np->point = MPFR_DECIMAL_POINT; + np->fp_trailing_zeros = spec.prec - np->fp_size; } + + if (spec.alt) + /* add decimal point even if no digits follow it */ + np->point = MPFR_DECIMAL_POINT; } return 0; diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2011-10-14 10:43:32.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2011-11-03 15:15:11.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0-p2"; + return "3.1.0-p3"; } diff -Naurd mpfr-3.1.0-a/tests/tsprintf.c mpfr-3.1.0-b/tests/tsprintf.c --- mpfr-3.1.0-a/tests/tsprintf.c 2011-10-03 08:17:14.000000000 +0000 +++ mpfr-3.1.0-b/tests/tsprintf.c 2011-11-03 15:15:11.000000000 +0000 @@ -475,6 +475,18 @@ check_sprintf ("-1.", "%- #0.1RG", x); /* precision zero */ + mpfr_set_d (x, 9.5, MPFR_RNDN); + check_sprintf ("9", "%.0RDf", x); + check_sprintf ("10", "%.0RUf", x); + + mpfr_set_d (x, 19.5, MPFR_RNDN); + check_sprintf ("19", "%.0RDf", x); + check_sprintf ("20", "%.0RUf", x); + + mpfr_set_d (x, 99.5, MPFR_RNDN); + check_sprintf ("99", "%.0RDf", x); + check_sprintf ("100", "%.0RUf", x); + mpfr_set_d (x, -9.5, MPFR_RNDN); check_sprintf ("-10", "%.0RDf", x); check_sprintf ("-10", "%.0RYf", x); @@ -1078,6 +1090,23 @@ mpfr_clear (x); } +static void +bug20111102 (void) +{ + mpfr_t t; + char s[100]; + + mpfr_init2 (t, 84); + mpfr_set_str (t, "999.99999999999999999999", 10, MPFR_RNDN); + mpfr_sprintf (s, "%.20RNg", t); + if (strcmp (s, "1000") != 0) + { + printf ("Error in bug20111102, expected 1000, got %s\n", s); + exit (1); + } + mpfr_clear (t); +} + /* In particular, the following test makes sure that the rounding * for %Ra and %Rb is not done on the MPFR number itself (as it * would overflow). Note: it has been reported on comp.std.c that @@ -1161,6 +1190,7 @@ locale = setlocale (LC_ALL, "C"); #endif + bug20111102 (); native_types (); hexadecimal (); binary (); diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2011-11-28 12:22:52.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2011-11-28 12:22:52.000000000 +0000 @@ -0,0 +1 @@ +gmp41compat diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2011-11-03 15:15:11.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2011-11-28 12:22:52.000000000 +0000 @@ -1 +1 @@ -3.1.0-p3 +3.1.0-p4 diff -Naurd mpfr-3.1.0-a/doc/mpfr.info mpfr-3.1.0-b/doc/mpfr.info --- mpfr-3.1.0-a/doc/mpfr.info 2011-10-03 09:43:04.000000000 +0000 +++ mpfr-3.1.0-b/doc/mpfr.info 2011-11-28 12:22:52.000000000 +0000 @@ -2994,11 +2994,12 @@ * `mpfr_urandom' and `mpfr_urandomb' changed in MPFR 3.1. Their behavior no longer depends on the platform (assuming this is also - true for GMP's random generator). As a consequence, the returned - values can be different between MPFR 3.1 and previous MPFR - versions. Note: as the reproducibility of these functions was not - specified before MPFR 3.1, the MPFR 3.1 behavior is _not_ regarded - as backward incompatible with previous versions. + true for GMP's random generator, which is not the case between GMP + 4.1 and 4.2 if `gmp_randinit_default' is used). As a consequence, + the returned values can be different between MPFR 3.1 and previous + MPFR versions. Note: as the reproducibility of these functions + was not specified before MPFR 3.1, the MPFR 3.1 behavior is _not_ + regarded as backward incompatible with previous versions.  @@ -4239,13 +4240,13 @@ Node: Type and Macro Changes129308 Node: Added Functions132029 Node: Changed Functions134972 -Node: Removed Functions139167 -Node: Other Changes139579 -Node: Contributors141108 -Node: References143574 -Node: GNU Free Documentation License145315 -Node: Concept Index167758 -Node: Function and Type Index173677 +Node: Removed Functions139253 +Node: Other Changes139665 +Node: Contributors141194 +Node: References143660 +Node: GNU Free Documentation License145401 +Node: Concept Index167844 +Node: Function and Type Index173763  End Tag Table diff -Naurd mpfr-3.1.0-a/doc/mpfr.texi mpfr-3.1.0-b/doc/mpfr.texi --- mpfr-3.1.0-a/doc/mpfr.texi 2011-10-03 08:17:14.000000000 +0000 +++ mpfr-3.1.0-b/doc/mpfr.texi 2011-11-28 12:22:52.000000000 +0000 @@ -3466,8 +3466,9 @@ a lack of specification. @item @code{mpfr_urandom} and @code{mpfr_urandomb} changed in MPFR 3.1. -Their behavior no longer depends on the platform (assuming this is also -true for GMP's random generator). As a consequence, the returned values +Their behavior no longer depends on the platform (assuming this is also true +for GMP's random generator, which is not the case between GMP 4.1 and 4.2 if +@code{gmp_randinit_default} is used). As a consequence, the returned values can be different between MPFR 3.1 and previous MPFR versions. Note: as the reproducibility of these functions was not specified before MPFR 3.1, the MPFR 3.1 behavior is @emph{not} regarded as diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2011-11-03 15:15:11.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2011-11-28 12:22:52.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0-p3" +#define MPFR_VERSION_STRING "3.1.0-p4" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2011-11-03 15:15:11.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2011-11-28 12:22:52.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0-p3"; + return "3.1.0-p4"; } diff -Naurd mpfr-3.1.0-a/tests/trandom.c mpfr-3.1.0-b/tests/trandom.c --- mpfr-3.1.0-a/tests/trandom.c 2011-10-03 08:17:14.000000000 +0000 +++ mpfr-3.1.0-b/tests/trandom.c 2011-11-28 12:22:52.000000000 +0000 @@ -114,21 +114,29 @@ mpfr_t x; gmp_randstate_t s; +#if __MPFR_GMP(4,2,0) +# define C1 "0.895943" +# define C2 "0.848824" +#else +# define C1 "0.479652" +# define C2 "0.648529" +#endif + gmp_randinit_default (s); gmp_randseed_ui (s, 42); mpfr_init2 (x, 17); mpfr_urandomb (x, s); - if (mpfr_cmp_str1 (x, "0.895943") != 0) + if (mpfr_cmp_str1 (x, C1) != 0) { - printf ("Error in bug20100914, expected 0.895943, got "); + printf ("Error in bug20100914, expected " C1 ", got "); mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN); printf ("\n"); exit (1); } mpfr_urandomb (x, s); - if (mpfr_cmp_str1 (x, "0.848824") != 0) + if (mpfr_cmp_str1 (x, C2) != 0) { - printf ("Error in bug20100914, expected 0.848824, got "); + printf ("Error in bug20100914, expected " C2 ", got "); mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN); printf ("\n"); exit (1); diff -Naurd mpfr-3.1.0-a/tests/turandom.c mpfr-3.1.0-b/tests/turandom.c --- mpfr-3.1.0-a/tests/turandom.c 2011-10-03 08:17:14.000000000 +0000 +++ mpfr-3.1.0-b/tests/turandom.c 2011-11-28 12:22:52.000000000 +0000 @@ -160,23 +160,29 @@ mpfr_t x; gmp_randstate_t s; +#if __MPFR_GMP(4,2,0) +# define C1 "0.8488312" +# define C2 "0.8156509" +#else +# define C1 "0.6485367" +# define C2 "0.9362717" +#endif + gmp_randinit_default (s); gmp_randseed_ui (s, 42); mpfr_init2 (x, 17); mpfr_urandom (x, s, MPFR_RNDN); - /* the following values are obtained on a 32-bit computer, we should get - the same values on a 64-bit computer */ - if (mpfr_cmp_str1 (x, "0.8488312") != 0) + if (mpfr_cmp_str1 (x, C1) != 0) { - printf ("Error in bug20100914, expected 0.8488312, got "); + printf ("Error in bug20100914, expected " C1 ", got "); mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN); printf ("\n"); exit (1); } mpfr_urandom (x, s, MPFR_RNDN); - if (mpfr_cmp_str1 (x, "0.8156509") != 0) + if (mpfr_cmp_str1 (x, C2) != 0) { - printf ("Error in bug20100914, expected 0.8156509, got "); + printf ("Error in bug20100914, expected " C2 ", got "); mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN); printf ("\n"); exit (1); diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2012-02-24 12:44:49.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2012-02-24 12:44:49.000000000 +0000 @@ -0,0 +1 @@ +logging-freeze diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2011-11-28 12:22:52.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2012-02-24 12:44:49.000000000 +0000 @@ -1 +1 @@ -3.1.0-p4 +3.1.0-p5 diff -Naurd mpfr-3.1.0-a/src/add_d.c mpfr-3.1.0-b/src/add_d.c --- mpfr-3.1.0-a/src/add_d.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/add_d.c 2012-02-24 12:44:49.000000000 +0000 @@ -34,7 +34,7 @@ (("b[%Pu]=%.*Rg c=%.20g rnd=%d", mpfr_get_prec(b), mpfr_log_prec, b, c, rnd_mode), ("a[%Pu]=%.*Rg inexact=%d", - mpfr_get_prec (a), mpfr_get_prec, a, inexact)); + mpfr_get_prec (a), mpfr_log_prec, a, inexact)); MPFR_SAVE_EXPO_MARK (expo); diff -Naurd mpfr-3.1.0-a/src/add_ui.c mpfr-3.1.0-b/src/add_ui.c --- mpfr-3.1.0-a/src/add_ui.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/add_ui.c 2012-02-24 12:44:49.000000000 +0000 @@ -29,7 +29,7 @@ MPFR_LOG_FUNC (("x[%Pu]=%.*Rg u=%d rnd=%d", mpfr_get_prec(x), mpfr_log_prec, x, u, rnd_mode), - ("y[%Pu]=%.*Rg", mpfr_get_prec (y), mpfr_get_prec, y)); + ("y[%Pu]=%.*Rg", mpfr_get_prec (y), mpfr_log_prec, y)); if (MPFR_LIKELY(u != 0) ) /* if u=0, do nothing */ { diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2011-11-28 12:22:52.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2012-02-24 12:44:49.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0-p4" +#define MPFR_VERSION_STRING "3.1.0-p5" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-3.1.0-a/src/mul_d.c mpfr-3.1.0-b/src/mul_d.c --- mpfr-3.1.0-a/src/mul_d.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/mul_d.c 2012-02-24 12:44:49.000000000 +0000 @@ -34,7 +34,7 @@ (("b[%Pu]=%.*Rg c=%.20g rnd=%d", mpfr_get_prec(b), mpfr_log_prec, b, c, rnd_mode), ("a[%Pu]=%.*Rg inexact=%d", - mpfr_get_prec (a), mpfr_get_prec, a, inexact)); + mpfr_get_prec (a), mpfr_log_prec, a, inexact)); MPFR_SAVE_EXPO_MARK (expo); diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2011-11-28 12:22:52.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2012-02-24 12:44:49.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0-p4"; + return "3.1.0-p5"; } diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2012-02-24 13:50:05.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2012-02-24 13:50:05.000000000 +0000 @@ -0,0 +1 @@ +logging-varfmt diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2012-02-24 12:44:49.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2012-02-24 13:50:05.000000000 +0000 @@ -1 +1 @@ -3.1.0-p5 +3.1.0-p6 diff -Naurd mpfr-3.1.0-a/src/mpfr-impl.h mpfr-3.1.0-b/src/mpfr-impl.h --- mpfr-3.1.0-a/src/mpfr-impl.h 2011-10-05 21:39:57.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr-impl.h 2012-02-24 13:50:05.000000000 +0000 @@ -1592,7 +1592,7 @@ do \ if ((MPFR_LOG_INTERNAL_F & mpfr_log_type) && \ (mpfr_log_current <= mpfr_log_level)) \ - LOG_PRINT ("%s.%d:%s[%#Pu]=%.*Rf\n", __func__, __LINE__, \ + LOG_PRINT ("%s.%d:%s[%#Pu]=%.*Rg\n", __func__, __LINE__, \ #x, mpfr_get_prec (x), mpfr_log_prec, x); \ while (0) diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2012-02-24 12:44:49.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2012-02-24 13:50:05.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0-p5" +#define MPFR_VERSION_STRING "3.1.0-p6" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2012-02-24 12:44:49.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2012-02-24 13:50:05.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0-p5"; + return "3.1.0-p6"; } diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2012-03-08 15:17:03.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2012-03-08 15:17:03.000000000 +0000 @@ -0,0 +1 @@ +large-prec diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2012-02-24 13:50:05.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2012-03-08 15:17:03.000000000 +0000 @@ -1 +1 @@ -3.1.0-p6 +3.1.0-p7 diff -Naurd mpfr-3.1.0-a/src/add1.c mpfr-3.1.0-b/src/add1.c --- mpfr-3.1.0-a/src/add1.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/add1.c 2012-03-08 15:17:03.000000000 +0000 @@ -44,12 +44,12 @@ bq = MPFR_PREC(b); cq = MPFR_PREC(c); - an = (aq-1)/GMP_NUMB_BITS+1; /* number of limbs of a */ + an = MPFR_PREC2LIMBS (aq); /* number of limbs of a */ aq2 = (mpfr_prec_t) an * GMP_NUMB_BITS; sh = aq2 - aq; /* non-significant bits in low limb */ - bn = (bq-1)/GMP_NUMB_BITS+1; /* number of limbs of b */ - cn = (cq-1)/GMP_NUMB_BITS+1; /* number of limbs of c */ + bn = MPFR_PREC2LIMBS (bq); /* number of limbs of b */ + cn = MPFR_PREC2LIMBS (cq); /* number of limbs of c */ ap = MPFR_MANT(a); bp = MPFR_MANT(b); @@ -124,7 +124,7 @@ dif = aq2 - diff_exp; /* dif is the number of bits of c which overlap with a' */ - difn = (dif-1)/GMP_NUMB_BITS + 1; + difn = MPFR_PREC2LIMBS (dif); /* only the highest difn limbs from c have to be considered */ if (MPFR_UNLIKELY(difn > cn)) { diff -Naurd mpfr-3.1.0-a/src/add1sp.c mpfr-3.1.0-b/src/add1sp.c --- mpfr-3.1.0-a/src/add1sp.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/add1sp.c 2012-03-08 15:17:03.000000000 +0000 @@ -107,7 +107,7 @@ /* Read prec and num of limbs */ p = MPFR_PREC(b); - n = (p+GMP_NUMB_BITS-1)/GMP_NUMB_BITS; + n = MPFR_PREC2LIMBS (p); MPFR_UNSIGNED_MINUS_MODULO(sh, p); bx = MPFR_GET_EXP(b); d = (mpfr_uexp_t) (bx - MPFR_GET_EXP(c)); diff -Naurd mpfr-3.1.0-a/src/agm.c mpfr-3.1.0-b/src/agm.c --- mpfr-3.1.0-a/src/agm.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/agm.c 2012-03-08 15:17:03.000000000 +0000 @@ -91,7 +91,7 @@ q = MPFR_PREC(r); p = q + MPFR_INT_CEIL_LOG2(q) + 15; MPFR_ASSERTD (p >= 7); /* see algorithms.tex */ - s = (p - 1) / GMP_NUMB_BITS + 1; + s = MPFR_PREC2LIMBS (p); /* b (op2) and a (op1) are the 2 operands but we want b >= a */ compare = mpfr_cmp (op1, op2); @@ -285,7 +285,7 @@ /* Next iteration */ MPFR_ZIV_NEXT (loop, p); - s = (p - 1) / GMP_NUMB_BITS + 1; + s = MPFR_PREC2LIMBS (p); } MPFR_ZIV_FREE (loop); diff -Naurd mpfr-3.1.0-a/src/eq.c mpfr-3.1.0-b/src/eq.c --- mpfr-3.1.0-a/src/eq.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/eq.c 2012-03-08 15:17:03.000000000 +0000 @@ -56,8 +56,8 @@ if (uexp != vexp) return 0; /* no bit agree */ - usize = (MPFR_PREC(u) - 1) / GMP_NUMB_BITS + 1; - vsize = (MPFR_PREC(v) - 1) / GMP_NUMB_BITS + 1; + usize = MPFR_LIMB_SIZE (u); + vsize = MPFR_LIMB_SIZE (v); if (vsize > usize) /* exchange u and v */ { diff -Naurd mpfr-3.1.0-a/src/exp.c mpfr-3.1.0-b/src/exp.c --- mpfr-3.1.0-a/src/exp.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/exp.c 2012-03-08 15:17:03.000000000 +0000 @@ -133,7 +133,7 @@ mp_size_t yn; int sh; - yn = 1 + (MPFR_PREC(y) - 1) / GMP_NUMB_BITS; + yn = MPFR_LIMB_SIZE (y); sh = (mpfr_prec_t) yn * GMP_NUMB_BITS - MPFR_PREC(y); MPFR_MANT(y)[0] += MPFR_LIMB_ONE << sh; inexact = 1; diff -Naurd mpfr-3.1.0-a/src/get_d.c mpfr-3.1.0-b/src/get_d.c --- mpfr-3.1.0-a/src/get_d.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/get_d.c 2012-03-08 15:17:03.000000000 +0000 @@ -100,7 +100,7 @@ nbits += (1021 + e); MPFR_ASSERTD (nbits >= 1); } - np = (nbits + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS; + np = MPFR_PREC2LIMBS (nbits); MPFR_ASSERTD ( np <= MPFR_LIMBS_PER_DOUBLE ); carry = mpfr_round_raw_4 (tp, MPFR_MANT(src), MPFR_PREC(src), negative, nbits, rnd_mode); diff -Naurd mpfr-3.1.0-a/src/get_flt.c mpfr-3.1.0-b/src/get_flt.c --- mpfr-3.1.0-a/src/get_flt.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/get_flt.c 2012-03-08 15:17:03.000000000 +0000 @@ -92,7 +92,7 @@ nbits += (125 + e); MPFR_ASSERTD (nbits >= 1); } - np = (nbits + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS; + np = MPFR_PREC2LIMBS (nbits); MPFR_ASSERTD(np <= MPFR_LIMBS_PER_FLT); carry = mpfr_round_raw_4 (tp, MPFR_MANT(src), MPFR_PREC(src), negative, nbits, rnd_mode); diff -Naurd mpfr-3.1.0-a/src/get_str.c mpfr-3.1.0-b/src/get_str.c --- mpfr-3.1.0-a/src/get_str.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/get_str.c 2012-03-08 15:17:03.000000000 +0000 @@ -2351,7 +2351,7 @@ /* the first digit will contain only r bits */ prec = (m - 1) * pow2 + r; /* total number of bits */ - n = (prec - 1) / GMP_NUMB_BITS + 1; + n = MPFR_PREC2LIMBS (prec); MPFR_TMP_MARK (marker); x1 = MPFR_TMP_LIMBS_ALLOC (n + 1); @@ -2417,12 +2417,12 @@ exact = 1; /* number of limbs */ - n = 1 + (prec - 1) / GMP_NUMB_BITS; + n = MPFR_PREC2LIMBS (prec); /* a will contain the approximation of the mantissa */ a = MPFR_TMP_LIMBS_ALLOC (n); - nx = 1 + (MPFR_PREC(x) - 1) / GMP_NUMB_BITS; + nx = MPFR_LIMB_SIZE (x); if ((mpfr_exp_t) m == g) /* final exponent is 0, no multiplication or division to perform */ diff -Naurd mpfr-3.1.0-a/src/init2.c mpfr-3.1.0-b/src/init2.c --- mpfr-3.1.0-a/src/init2.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/init2.c 2012-03-08 15:17:03.000000000 +0000 @@ -51,7 +51,7 @@ which both have an odd mantissa */ MPFR_ASSERTN(p >= MPFR_PREC_MIN && p <= MPFR_PREC_MAX); - xsize = (mp_size_t) ((p - 1) / GMP_NUMB_BITS) + 1; + xsize = MPFR_PREC2LIMBS (p); tmp = (mpfr_limb_ptr) (*__gmp_allocate_func)(MPFR_MALLOC_SIZE(xsize)); MPFR_PREC(x) = p; /* Set prec */ diff -Naurd mpfr-3.1.0-a/src/lngamma.c mpfr-3.1.0-b/src/lngamma.c --- mpfr-3.1.0-a/src/lngamma.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/lngamma.c 2012-03-08 15:17:03.000000000 +0000 @@ -67,7 +67,7 @@ /* Now, the unit bit is represented. */ - prec = ((prec - 1) / GMP_NUMB_BITS + 1) * GMP_NUMB_BITS - expo; + prec = MPFR_PREC2LIMBS (prec) * GMP_NUMB_BITS - expo; /* number of represented fractional bits (including the trailing 0's) */ x0 = *(MPFR_MANT (x) + prec / GMP_NUMB_BITS); diff -Naurd mpfr-3.1.0-a/src/mpfr-impl.h mpfr-3.1.0-b/src/mpfr-impl.h --- mpfr-3.1.0-a/src/mpfr-impl.h 2012-02-24 13:50:05.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr-impl.h 2012-03-09 12:06:26.000000000 +0000 @@ -646,10 +646,24 @@ **************** mpfr_t properties ******************* ******************************************************/ +/* In the following macro, p is usually a mpfr_prec_t, but this macro + works with other integer types (without integer overflow). Checking + that p >= 1 in debug mode is useful here because this macro can be + used on a computed precision (in particular, this formula does not + work for a degenerate case p = 0, and could give different results + on different platforms). But let us not use an assertion checking + in the MPFR_LAST_LIMB() and MPFR_LIMB_SIZE() macros below to avoid + too much expansion for assertions (in practice, this should be a + problem just when testing MPFR with the --enable-assert configure + option and the -ansi -pedantic-errors gcc compiler flags). */ +#define MPFR_PREC2LIMBS(p) \ + (MPFR_ASSERTD ((p) >= 1), ((p) - 1) / GMP_NUMB_BITS + 1) + #define MPFR_PREC(x) ((x)->_mpfr_prec) #define MPFR_EXP(x) ((x)->_mpfr_exp) #define MPFR_MANT(x) ((x)->_mpfr_d) -#define MPFR_LIMB_SIZE(x) ((MPFR_PREC((x))-1)/GMP_NUMB_BITS+1) +#define MPFR_LAST_LIMB(x) ((MPFR_PREC (x) - 1) / GMP_NUMB_BITS) +#define MPFR_LIMB_SIZE(x) (MPFR_LAST_LIMB (x) + 1) /****************************************************** @@ -749,7 +763,8 @@ #define MPFR_IS_FP(x) (!MPFR_IS_NAN(x) && !MPFR_IS_INF(x)) #define MPFR_IS_SINGULAR(x) (MPFR_EXP(x) <= MPFR_EXP_INF) #define MPFR_IS_PURE_FP(x) (!MPFR_IS_SINGULAR(x) && \ - (MPFR_ASSERTD (MPFR_MANT(x)[MPFR_LIMB_SIZE(x)-1] & MPFR_LIMB_HIGHBIT), 1)) + (MPFR_ASSERTD ((MPFR_MANT(x)[MPFR_LAST_LIMB(x)] \ + & MPFR_LIMB_HIGHBIT) != 0), 1)) #define MPFR_ARE_SINGULAR(x,y) \ (MPFR_UNLIKELY(MPFR_IS_SINGULAR(x)) || MPFR_UNLIKELY(MPFR_IS_SINGULAR(y))) @@ -1061,7 +1076,7 @@ /* Set a number to 1 (Fast) - It doesn't check if 1 is in the exponent range */ #define MPFR_SET_ONE(x) \ do { \ - mp_size_t _size = MPFR_LIMB_SIZE(x) - 1; \ + mp_size_t _size = MPFR_LAST_LIMB(x); \ MPFR_SET_POS(x); \ MPFR_EXP(x) = 1; \ MPN_ZERO ( MPFR_MANT(x), _size); \ @@ -1213,8 +1228,8 @@ _destp = MPFR_MANT (dest); \ if (MPFR_UNLIKELY (_destprec >= _srcprec)) \ { \ - _srcs = (_srcprec + GMP_NUMB_BITS-1)/GMP_NUMB_BITS; \ - _dests = (_destprec + GMP_NUMB_BITS-1)/GMP_NUMB_BITS - _srcs; \ + _srcs = MPFR_PREC2LIMBS (_srcprec); \ + _dests = MPFR_PREC2LIMBS (_destprec) - _srcs; \ MPN_COPY (_destp + _dests, srcp, _srcs); \ MPN_ZERO (_destp, _dests); \ inexact = 0; \ @@ -1227,8 +1242,8 @@ mp_limb_t _rb, _sb, _ulp; \ \ /* Compute Position and shift */ \ - _srcs = (_srcprec + GMP_NUMB_BITS-1)/GMP_NUMB_BITS; \ - _dests = (_destprec + GMP_NUMB_BITS-1)/GMP_NUMB_BITS; \ + _srcs = MPFR_PREC2LIMBS (_srcprec); \ + _dests = MPFR_PREC2LIMBS (_destprec); \ MPFR_UNSIGNED_MINUS_MODULO (_sh, _destprec); \ _sp = (srcp) + _srcs - _dests; \ \ @@ -1372,7 +1387,7 @@ if (MPFR_LIKELY (MPFR_PREC (dest) == MPFR_PREC (src))) \ { \ MPN_COPY (MPFR_MANT (dest), MPFR_MANT (src), \ - (MPFR_PREC (src) + GMP_NUMB_BITS-1)/GMP_NUMB_BITS); \ + MPFR_LIMB_SIZE (src)); \ inexact = 0; \ } \ else \ @@ -1682,7 +1697,7 @@ MPFR_ASSERTD (_prec >= MPFR_PREC_MIN); \ if (MPFR_UNLIKELY (_prec > MPFR_PREC_MAX)) \ mpfr_abort_prec_max (); \ - _size = (mpfr_prec_t) (_prec + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS; \ + _size = MPFR_PREC2LIMBS (_prec); \ if (MPFR_UNLIKELY (_size * (num) > MPFR_GROUP_STATIC_SIZE)) \ { \ (g).alloc = (num) * _size * sizeof (mp_limb_t); \ @@ -1733,7 +1748,7 @@ MPFR_ASSERTD (_prec >= MPFR_PREC_MIN); \ if (MPFR_UNLIKELY (_prec > MPFR_PREC_MAX)) \ mpfr_abort_prec_max (); \ - _size = (mpfr_prec_t) (_prec + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS; \ + _size = MPFR_PREC2LIMBS (_prec); \ (g).alloc = (num) * _size * sizeof (mp_limb_t); \ if (MPFR_LIKELY (_oalloc == 0)) \ (g).mant = (mp_limb_t *) (*__gmp_allocate_func) ((g).alloc); \ @@ -1886,7 +1901,7 @@ MPFR_NORETURN_ATTR; __MPFR_DECLSPEC void mpfr_rand_raw _MPFR_PROTO((mpfr_limb_ptr, gmp_randstate_t, - unsigned long)); + mpfr_prec_t)); __MPFR_DECLSPEC mpz_t* mpfr_bernoulli_internal _MPFR_PROTO((mpz_t*, unsigned long)); diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2012-02-24 13:50:05.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2012-03-08 15:17:03.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0-p6" +#define MPFR_VERSION_STRING "3.1.0-p7" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-3.1.0-a/src/mul.c mpfr-3.1.0-b/src/mul.c --- mpfr-3.1.0-a/src/mul.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/mul.c 2012-03-08 15:17:03.000000000 +0000 @@ -93,15 +93,15 @@ ax = MPFR_GET_EXP (b) + MPFR_GET_EXP (c); - bq = MPFR_PREC(b); - cq = MPFR_PREC(c); + bq = MPFR_PREC (b); + cq = MPFR_PREC (c); - MPFR_ASSERTD(bq+cq > bq); /* PREC_MAX is /2 so no integer overflow */ + MPFR_ASSERTN ((mpfr_uprec_t) bq + cq <= MPFR_PREC_MAX); - bn = (bq+GMP_NUMB_BITS-1)/GMP_NUMB_BITS; /* number of limbs of b */ - cn = (cq+GMP_NUMB_BITS-1)/GMP_NUMB_BITS; /* number of limbs of c */ + bn = MPFR_PREC2LIMBS (bq); /* number of limbs of b */ + cn = MPFR_PREC2LIMBS (cq); /* number of limbs of c */ k = bn + cn; /* effective nb of limbs used by b*c (= tn or tn+1) below */ - tn = (bq + cq + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS; + tn = MPFR_PREC2LIMBS (bq + cq); /* <= k, thus no int overflow */ MPFR_ASSERTD(tn <= k); @@ -292,12 +292,12 @@ bq = MPFR_PREC (b); cq = MPFR_PREC (c); - MPFR_ASSERTD (bq+cq > bq); /* PREC_MAX is /2 so no integer overflow */ + MPFR_ASSERTN ((mpfr_uprec_t) bq + cq <= MPFR_PREC_MAX); - bn = (bq+GMP_NUMB_BITS-1)/GMP_NUMB_BITS; /* number of limbs of b */ - cn = (cq+GMP_NUMB_BITS-1)/GMP_NUMB_BITS; /* number of limbs of c */ + bn = MPFR_PREC2LIMBS (bq); /* number of limbs of b */ + cn = MPFR_PREC2LIMBS (cq); /* number of limbs of c */ k = bn + cn; /* effective nb of limbs used by b*c (= tn or tn+1) below */ - tn = (bq + cq + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS; + tn = MPFR_PREC2LIMBS (bq + cq); MPFR_ASSERTD (tn <= k); /* tn <= k, thus no int overflow */ /* Check for no size_t overflow*/ diff -Naurd mpfr-3.1.0-a/src/pow.c mpfr-3.1.0-b/src/pow.c --- mpfr-3.1.0-a/src/pow.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/pow.c 2012-03-08 15:17:03.000000000 +0000 @@ -136,7 +136,7 @@ (b) all the 'z' bits are zero */ - prec = ((prec - 1) / GMP_NUMB_BITS + 1) * GMP_NUMB_BITS - expo; + prec = MPFR_PREC2LIMBS (prec) * GMP_NUMB_BITS - expo; /* number of z+0 bits */ yn = prec / GMP_NUMB_BITS; diff -Naurd mpfr-3.1.0-a/src/print_raw.c mpfr-3.1.0-b/src/print_raw.c --- mpfr-3.1.0-a/src/print_raw.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/print_raw.c 2012-03-08 15:17:03.000000000 +0000 @@ -84,7 +84,7 @@ int i; mpfr_prec_t count = 0; char c; - mp_size_t n = (r - 1) / GMP_NUMB_BITS + 1; + mp_size_t n = MPFR_PREC2LIMBS (r); printf("%s ", str); for(n-- ; n>=0 ; n--) @@ -109,7 +109,7 @@ int i; mpfr_prec_t count = 0; char c; - mp_size_t n = (r - 1) / GMP_NUMB_BITS + 1; + mp_size_t n = MPFR_PREC2LIMBS (r); for(n-- ; n>=0 ; n--) { diff -Naurd mpfr-3.1.0-a/src/round_prec.c mpfr-3.1.0-b/src/round_prec.c --- mpfr-3.1.0-a/src/round_prec.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/round_prec.c 2012-03-08 15:17:03.000000000 +0000 @@ -55,12 +55,12 @@ MPFR_ASSERTN(prec >= MPFR_PREC_MIN && prec <= MPFR_PREC_MAX); - nw = 1 + (prec - 1) / GMP_NUMB_BITS; /* needed allocated limbs */ + nw = MPFR_PREC2LIMBS (prec); /* needed allocated limbs */ /* check if x has enough allocated space for the significand */ /* Get the number of limbs from the precision. (Compatible with all allocation methods) */ - ow = (MPFR_PREC (x) + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS; + ow = MPFR_LIMB_SIZE (x); if (nw > ow) { /* FIXME: Variable can't be created using custom allocation, diff -Naurd mpfr-3.1.0-a/src/round_raw_generic.c mpfr-3.1.0-b/src/round_raw_generic.c --- mpfr-3.1.0-a/src/round_raw_generic.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/round_raw_generic.c 2012-03-08 15:17:03.000000000 +0000 @@ -80,7 +80,7 @@ (xprec <= yprec || MPFR_IS_LIKE_RNDZ (rnd_mode, neg))) return 0; - xsize = (xprec-1)/GMP_NUMB_BITS + 1; + xsize = MPFR_PREC2LIMBS (xprec); nw = yprec / GMP_NUMB_BITS; rw = yprec & (GMP_NUMB_BITS - 1); diff -Naurd mpfr-3.1.0-a/src/set.c mpfr-3.1.0-b/src/set.c --- mpfr-3.1.0-a/src/set.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/set.c 2012-03-08 15:17:03.000000000 +0000 @@ -48,8 +48,7 @@ /* Same precision and b is not singular: * just copy the mantissa, and set the exponent and the sign * The result is exact. */ - MPN_COPY (MPFR_MANT (a), MPFR_MANT (b), - (MPFR_PREC (b) + GMP_NUMB_BITS-1)/GMP_NUMB_BITS); + MPN_COPY (MPFR_MANT (a), MPFR_MANT (b), MPFR_LIMB_SIZE (b)); MPFR_RET (0); } else diff -Naurd mpfr-3.1.0-a/src/set_f.c mpfr-3.1.0-b/src/set_f.c --- mpfr-3.1.0-a/src/set_f.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/set_f.c 2012-03-08 15:17:03.000000000 +0000 @@ -43,7 +43,7 @@ if (SIZ(x) * MPFR_FROM_SIGN_TO_INT(MPFR_SIGN(y)) < 0) MPFR_CHANGE_SIGN (y); - sy = 1 + (MPFR_PREC(y) - 1) / GMP_NUMB_BITS; + sy = MPFR_LIMB_SIZE (y); my = MPFR_MANT(y); mx = PTR(x); diff -Naurd mpfr-3.1.0-a/src/set_prec.c mpfr-3.1.0-b/src/set_prec.c --- mpfr-3.1.0-a/src/set_prec.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/set_prec.c 2012-03-08 15:17:03.000000000 +0000 @@ -32,7 +32,7 @@ MPFR_ASSERTN (p >= MPFR_PREC_MIN && p <= MPFR_PREC_MAX); /* Calculate the new number of limbs */ - xsize = (p - 1) / GMP_NUMB_BITS + 1; + xsize = MPFR_PREC2LIMBS (p); /* Realloc only if the new size is greater than the old */ xoldsize = MPFR_GET_ALLOC_SIZE (x); diff -Naurd mpfr-3.1.0-a/src/setmax.c mpfr-3.1.0-b/src/setmax.c --- mpfr-3.1.0-a/src/setmax.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/setmax.c 2012-03-08 15:17:03.000000000 +0000 @@ -32,7 +32,7 @@ mp_limb_t *xp; MPFR_SET_EXP (x, e); - xn = 1 + (MPFR_PREC(x) - 1) / GMP_NUMB_BITS; + xn = MPFR_LIMB_SIZE (x); sh = (mpfr_prec_t) xn * GMP_NUMB_BITS - MPFR_PREC(x); xp = MPFR_MANT(x); xp[0] = MP_LIMB_T_MAX << sh; diff -Naurd mpfr-3.1.0-a/src/sqr.c mpfr-3.1.0-b/src/sqr.c --- mpfr-3.1.0-a/src/sqr.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/sqr.c 2012-03-08 15:17:03.000000000 +0000 @@ -56,11 +56,11 @@ ax = 2 * MPFR_GET_EXP (b); bq = MPFR_PREC(b); - MPFR_ASSERTD (2 * bq > bq); /* PREC_MAX is /2 so no integer overflow */ + MPFR_ASSERTN (2 * (mpfr_uprec_t) bq <= MPFR_PREC_MAX); - bn = MPFR_LIMB_SIZE(b); /* number of limbs of b */ - tn = 1 + (2 * bq - 1) / GMP_NUMB_BITS; /* number of limbs of square, - 2*bn or 2*bn-1 */ + bn = MPFR_LIMB_SIZE (b); /* number of limbs of b */ + tn = MPFR_PREC2LIMBS (2 * bq); /* number of limbs of square, + 2*bn or 2*bn-1 */ if (MPFR_UNLIKELY(bn > MPFR_SQR_THRESHOLD)) return mpfr_mul (a, b, b, rnd_mode); diff -Naurd mpfr-3.1.0-a/src/stack_interface.c mpfr-3.1.0-b/src/stack_interface.c --- mpfr-3.1.0-a/src/stack_interface.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/stack_interface.c 2012-03-08 15:17:03.000000000 +0000 @@ -26,7 +26,7 @@ size_t mpfr_custom_get_size (mpfr_prec_t prec) { - return (prec + GMP_NUMB_BITS -1) / GMP_NUMB_BITS * BYTES_PER_MP_LIMB; + return MPFR_PREC2LIMBS (prec) * BYTES_PER_MP_LIMB; } #undef mpfr_custom_init diff -Naurd mpfr-3.1.0-a/src/strtofr.c mpfr-3.1.0-b/src/strtofr.c --- mpfr-3.1.0-a/src/strtofr.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/strtofr.c 2012-03-08 15:17:03.000000000 +0000 @@ -467,7 +467,7 @@ /* Set y to the value of the ~prec most significant bits of pstr->mant (as long as we guarantee correct rounding, we don't need to get exactly prec bits). */ - ysize = (prec - 1) / GMP_NUMB_BITS + 1; + ysize = MPFR_PREC2LIMBS (prec); /* prec bits corresponds to ysize limbs */ ysize_bits = ysize * GMP_NUMB_BITS; /* and to ysize_bits >= prec > MPFR_PREC (x) bits */ diff -Naurd mpfr-3.1.0-a/src/sub1sp.c mpfr-3.1.0-b/src/sub1sp.c --- mpfr-3.1.0-a/src/sub1sp.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/sub1sp.c 2012-03-08 15:17:03.000000000 +0000 @@ -155,8 +155,8 @@ MPFR_ASSERTD(MPFR_IS_PURE_FP(c)); /* Read prec and num of limbs */ - p = MPFR_PREC(b); - n = (p-1)/GMP_NUMB_BITS+1; + p = MPFR_PREC (b); + n = MPFR_PREC2LIMBS (p); /* Fast cmp of |b| and |c|*/ bx = MPFR_GET_EXP (b); diff -Naurd mpfr-3.1.0-a/src/urandomb.c mpfr-3.1.0-b/src/urandomb.c --- mpfr-3.1.0-a/src/urandomb.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/urandomb.c 2012-03-08 15:17:03.000000000 +0000 @@ -31,13 +31,20 @@ a sufficient number of limbs */ void mpfr_rand_raw (mpfr_limb_ptr mp, gmp_randstate_t rstate, - unsigned long int nbits) + mpfr_prec_t nbits) { mpz_t z; + MPFR_ASSERTN (nbits >= 1); /* To be sure to avoid the potential allocation of mpz_urandomb */ - ALLOC(z) = SIZ(z) = ((nbits - 1) / GMP_NUMB_BITS) + 1; + ALLOC(z) = SIZ(z) = MPFR_PREC2LIMBS (nbits); PTR(z) = mp; +#if __MPFR_GMP(5,0,0) + /* Check for integer overflow (unless mp_bitcnt_t is signed, + but according to the GMP manual, this shouldn't happen). + Note: mp_bitcnt_t has been introduced in GMP 5.0.0. */ + MPFR_ASSERTN ((mp_bitcnt_t) -1 < 0 || nbits <= (mp_bitcnt_t) -1); +#endif mpz_urandomb (z, rstate, nbits); } diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2012-02-24 13:50:05.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2012-03-08 15:17:03.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0-p6"; + return "3.1.0-p7"; } diff -Naurd mpfr-3.1.0-a/tests/tinits.c mpfr-3.1.0-b/tests/tinits.c --- mpfr-3.1.0-a/tests/tinits.c 2011-10-03 08:17:14.000000000 +0000 +++ mpfr-3.1.0-b/tests/tinits.c 2012-03-08 15:17:03.000000000 +0000 @@ -1,4 +1,4 @@ -/* Test file for mpfr_inits, mpfr_inits2 and mpfr_clears. +/* Test file for mpfr_init2, mpfr_inits, mpfr_inits2 and mpfr_clears. Copyright 2003, 2006, 2007, 2008, 2009, 2010, 2011 Free Software Foundation, Inc. Contributed by the Arenaire and Caramel projects, INRIA. @@ -20,18 +20,43 @@ http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */ +#include + #include "mpfr-test.h" int main (void) { mpfr_t a, b, c; + long large_prec; tests_start_mpfr (); + mpfr_inits (a, b, c, (mpfr_ptr) 0); mpfr_clears (a, b, c, (mpfr_ptr) 0); mpfr_inits2 (200, a, b, c, (mpfr_ptr) 0); mpfr_clears (a, b, c, (mpfr_ptr) 0); + + /* test for precision 2^31-1, see + https://gforge.inria.fr/tracker/index.php?func=detail&aid=13918 */ + large_prec = 2147483647; + if (getenv ("MPFR_CHECK_LARGEMEM") != NULL) + { + /* We assume that the precision won't be increased internally. */ + if (large_prec > MPFR_PREC_MAX) + large_prec = MPFR_PREC_MAX; + mpfr_inits2 (large_prec, a, b, (mpfr_ptr) 0); + mpfr_set_ui (a, 17, MPFR_RNDN); + mpfr_set (b, a, MPFR_RNDN); + if (mpfr_get_ui (a, MPFR_RNDN) != 17) + { + printf ("Error in mpfr_init2 with precision 2^31-1\n"); + exit (1); + } + mpfr_clears (a, b, (mpfr_ptr) 0); + } + tests_end_mpfr (); + return 0; } diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2012-03-12 11:59:47.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2012-03-12 11:59:47.000000000 +0000 @@ -0,0 +1 @@ +__gmp_const diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2012-03-08 15:17:03.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2012-03-12 11:59:47.000000000 +0000 @@ -1 +1 @@ -3.1.0-p7 +3.1.0-p8 diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2012-03-08 15:17:03.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2012-03-12 11:59:47.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0-p7" +#define MPFR_VERSION_STRING "3.1.0-p8" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) @@ -39,6 +39,18 @@ # include #endif +/* GMP's internal __gmp_const macro has been removed on 2012-03-04: + http://gmplib.org:8000/gmp/rev/d287cfaf6732 + const is standard and now assumed to be available. If the __gmp_const + definition is no longer present in GMP, this probably means that GMP + assumes that const is available; thus let's define it to const. + Note: this is a temporary fix that can be backported to previous MPFR + versions. In the future, __gmp_const should be replaced by const like + in GMP. */ +#ifndef __gmp_const +# define __gmp_const const +#endif + /* Avoid some problems with macro expansion if the user defines macros with the same name as keywords. By convention, identifiers and macro names starting with mpfr_ are reserved by MPFR. */ diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2012-03-08 15:17:03.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2012-03-12 11:59:47.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0-p7"; + return "3.1.0-p8"; } diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2012-04-27 01:13:15.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2012-04-27 01:13:15.000000000 +0000 @@ -0,0 +1 @@ +gamma-underflow diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2012-03-12 11:59:47.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2012-04-27 01:13:15.000000000 +0000 @@ -1 +1 @@ -3.1.0-p8 +3.1.0-p9 diff -Naurd mpfr-3.1.0-a/src/gamma.c mpfr-3.1.0-b/src/gamma.c --- mpfr-3.1.0-a/src/gamma.c 2011-10-03 08:17:09.000000000 +0000 +++ mpfr-3.1.0-b/src/gamma.c 2012-04-27 01:13:15.000000000 +0000 @@ -296,7 +296,7 @@ /* we want an upper bound for x * [log(2-x)-1]. since x < 0, we need a lower bound on log(2-x) */ mpfr_ui_sub (xp, 2, x, MPFR_RNDD); - mpfr_log2 (xp, xp, MPFR_RNDD); + mpfr_log (xp, xp, MPFR_RNDD); mpfr_sub_ui (xp, xp, 1, MPFR_RNDD); mpfr_mul (xp, xp, x, MPFR_RNDU); diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2012-03-12 11:59:47.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2012-04-27 01:13:15.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0-p8" +#define MPFR_VERSION_STRING "3.1.0-p9" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2012-03-12 11:59:47.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2012-04-27 01:13:15.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0-p8"; + return "3.1.0-p9"; } diff -Naurd mpfr-3.1.0-a/tests/tgamma.c mpfr-3.1.0-b/tests/tgamma.c --- mpfr-3.1.0-a/tests/tgamma.c 2011-10-03 08:17:14.000000000 +0000 +++ mpfr-3.1.0-b/tests/tgamma.c 2012-04-27 01:13:15.000000000 +0000 @@ -478,6 +478,36 @@ mpfr_clear (x); } +/* bug found by Giridhar Tammana */ +static void +test20120426 (void) +{ + mpfr_t xa, xb; + int i; + mpfr_exp_t emin; + + mpfr_init2 (xa, 53); + mpfr_init2 (xb, 53); + mpfr_set_d (xb, -168.5, MPFR_RNDN); + emin = mpfr_get_emin (); + mpfr_set_emin (-1073); + i = mpfr_gamma (xa, xb, MPFR_RNDN); + i = mpfr_subnormalize (xa, i, MPFR_RNDN); /* new ternary value */ + mpfr_set_str (xb, "-9.5737343987585366746184749943e-304", 10, MPFR_RNDN); + if (!((i > 0) && (mpfr_cmp (xa, xb) == 0))) + { + printf ("Error in test20120426, i=%d\n", i); + printf ("expected "); + mpfr_print_binary (xb); putchar ('\n'); + printf ("got "); + mpfr_print_binary (xa); putchar ('\n'); + exit (1); + } + mpfr_set_emin (emin); + mpfr_clear (xa); + mpfr_clear (xb); +} + static void exprange (void) { @@ -821,6 +851,7 @@ gamma_integer (); test20071231 (); test20100709 (); + test20120426 (); data_check ("data/gamma", mpfr_gamma, "mpfr_gamma"); diff -Naurd mpfr-3.1.0-a/PATCHES mpfr-3.1.0-b/PATCHES --- mpfr-3.1.0-a/PATCHES 2012-05-07 18:52:45.000000000 +0000 +++ mpfr-3.1.0-b/PATCHES 2012-05-07 18:52:45.000000000 +0000 @@ -0,0 +1 @@ +gamma-overunderflow diff -Naurd mpfr-3.1.0-a/VERSION mpfr-3.1.0-b/VERSION --- mpfr-3.1.0-a/VERSION 2012-04-27 01:13:15.000000000 +0000 +++ mpfr-3.1.0-b/VERSION 2012-05-07 18:52:45.000000000 +0000 @@ -1 +1 @@ -3.1.0-p9 +3.1.0-p10 diff -Naurd mpfr-3.1.0-a/src/gamma.c mpfr-3.1.0-b/src/gamma.c --- mpfr-3.1.0-a/src/gamma.c 2012-04-27 01:13:15.000000000 +0000 +++ mpfr-3.1.0-b/src/gamma.c 2012-05-07 18:52:45.000000000 +0000 @@ -100,7 +100,8 @@ mpfr_t xp, GammaTrial, tmp, tmp2; mpz_t fact; mpfr_prec_t realprec; - int compared, inex, is_integer; + int compared, is_integer; + int inex = 0; /* 0 means: result gamma not set yet */ MPFR_GROUP_DECL (group); MPFR_SAVE_EXPO_DECL (expo); MPFR_ZIV_DECL (loop); @@ -377,6 +378,15 @@ mpfr_mul (GammaTrial, tmp2, xp, MPFR_RNDN); /* Pi*(2-x), error (1+u)^2 */ err_g = MPFR_GET_EXP(GammaTrial); mpfr_sin (GammaTrial, GammaTrial, MPFR_RNDN); /* sin(Pi*(2-x)) */ + /* If tmp is +Inf, we compute exp(lngamma(x)). */ + if (mpfr_inf_p (tmp)) + { + inex = mpfr_explgamma (gamma, x, &expo, tmp, tmp2, rnd_mode); + if (inex) + goto end; + else + goto ziv_next; + } err_g = err_g + 1 - MPFR_GET_EXP(GammaTrial); /* let g0 the true value of Pi*(2-x), g the computed value. We have g = g0 + h with |h| <= |(1+u^2)-1|*g. @@ -411,11 +421,16 @@ if (MPFR_LIKELY (MPFR_CAN_ROUND (GammaTrial, realprec - err_g, MPFR_PREC(gamma), rnd_mode))) break; + + ziv_next: MPFR_ZIV_NEXT (loop, realprec); } + + end: MPFR_ZIV_FREE (loop); - inex = mpfr_set (gamma, GammaTrial, rnd_mode); + if (inex == 0) + inex = mpfr_set (gamma, GammaTrial, rnd_mode); MPFR_GROUP_CLEAR (group); mpz_clear (fact); diff -Naurd mpfr-3.1.0-a/src/lngamma.c mpfr-3.1.0-b/src/lngamma.c --- mpfr-3.1.0-a/src/lngamma.c 2012-03-08 15:17:03.000000000 +0000 +++ mpfr-3.1.0-b/src/lngamma.c 2012-05-07 18:52:45.000000000 +0000 @@ -49,9 +49,72 @@ mpfr_set_ui_2exp (s, 9, -1, MPFR_RNDN); /* 4.5 */ } -#ifndef IS_GAMMA +#ifdef IS_GAMMA + +/* This function is called in case of intermediate overflow/underflow. + The s1 and s2 arguments are temporary MPFR numbers, having the + working precision. If the result could be determined, then the + flags are updated via pexpo, y is set to the result, and the + (non-zero) ternary value is returned. Otherwise 0 is returned + in order to perform the next Ziv iteration. */ static int -unit_bit (mpfr_srcptr (x)) +mpfr_explgamma (mpfr_ptr y, mpfr_srcptr x, mpfr_save_expo_t *pexpo, + mpfr_ptr s1, mpfr_ptr s2, mpfr_rnd_t rnd) +{ + mpfr_t t1, t2; + int inex1, inex2, sign; + MPFR_BLOCK_DECL (flags1); + MPFR_BLOCK_DECL (flags2); + MPFR_GROUP_DECL (group); + + MPFR_BLOCK (flags1, inex1 = mpfr_lgamma (s1, &sign, x, MPFR_RNDD)); + MPFR_ASSERTN (inex1 != 0); + /* s1 = RNDD(lngamma(x)), inexact */ + if (MPFR_UNLIKELY (MPFR_OVERFLOW (flags1))) + { + if (MPFR_SIGN (s1) > 0) + { + MPFR_SAVE_EXPO_UPDATE_FLAGS (*pexpo, MPFR_FLAGS_OVERFLOW); + return mpfr_overflow (y, rnd, sign); + } + else + { + MPFR_SAVE_EXPO_UPDATE_FLAGS (*pexpo, MPFR_FLAGS_UNDERFLOW); + return mpfr_underflow (y, rnd == MPFR_RNDN ? MPFR_RNDZ : rnd, sign); + } + } + + mpfr_set (s2, s1, MPFR_RNDN); /* exact */ + mpfr_nextabove (s2); /* v = RNDU(lngamma(z0)) */ + + if (sign < 0) + rnd = MPFR_INVERT_RND (rnd); /* since the result with be negated */ + MPFR_GROUP_INIT_2 (group, MPFR_PREC (y), t1, t2); + MPFR_BLOCK (flags1, inex1 = mpfr_exp (t1, s1, rnd)); + MPFR_BLOCK (flags2, inex2 = mpfr_exp (t2, s2, rnd)); + /* t1 is the rounding with mode 'rnd' of a lower bound on |Gamma(x)|, + t2 is the rounding with mode 'rnd' of an upper bound, thus if both + are equal, so is the wanted result. If t1 and t2 differ or the flags + differ, at some point of Ziv's loop they should agree. */ + if (mpfr_equal_p (t1, t2) && flags1 == flags2) + { + MPFR_ASSERTN ((inex1 > 0 && inex2 > 0) || (inex1 < 0 && inex2 < 0)); + mpfr_set4 (y, t1, MPFR_RNDN, sign); /* exact */ + if (sign < 0) + inex1 = - inex1; + MPFR_SAVE_EXPO_UPDATE_FLAGS (*pexpo, flags1); + } + else + inex1 = 0; /* couldn't determine the result */ + MPFR_GROUP_CLEAR (group); + + return inex1; +} + +#else + +static int +unit_bit (mpfr_srcptr x) { mpfr_exp_t expo; mpfr_prec_t prec; @@ -75,6 +138,7 @@ return (x0 >> (prec % GMP_NUMB_BITS)) & 1; } + #endif /* lngamma(x) = log(gamma(x)). @@ -99,12 +163,14 @@ mpfr_t s, t, u, v, z; unsigned long m, k, maxm; mpz_t *INITIALIZED(B); /* variable B declared as initialized */ - int inexact, compared; + int compared; + int inexact = 0; /* 0 means: result y not set yet */ mpfr_exp_t err_s, err_t; unsigned long Bm = 0; /* number of allocated B[] */ unsigned long oldBm; double d; MPFR_SAVE_EXPO_DECL (expo); + MPFR_ZIV_DECL (loop); compared = mpfr_cmp_ui (z0, 1); @@ -122,7 +188,7 @@ if (MPFR_EXP(z0) <= - (mpfr_exp_t) MPFR_PREC(y)) { mpfr_t l, h, g; - int ok, inex2; + int ok, inex1, inex2; mpfr_prec_t prec = MPFR_PREC(y) + 14; MPFR_ZIV_DECL (loop); @@ -157,14 +223,14 @@ mpfr_sub (h, h, g, MPFR_RNDD); mpfr_mul (g, z0, z0, MPFR_RNDU); mpfr_add (h, h, g, MPFR_RNDU); - inexact = mpfr_prec_round (l, MPFR_PREC(y), rnd); + inex1 = mpfr_prec_round (l, MPFR_PREC(y), rnd); inex2 = mpfr_prec_round (h, MPFR_PREC(y), rnd); /* Caution: we not only need l = h, but both inexact flags should agree. Indeed, one of the inexact flags might be zero. In that case if we assume lngamma(z0) cannot be exact, the other flag should be correct. We are conservative here and request that both inexact flags agree. */ - ok = SAME_SIGN (inexact, inex2) && mpfr_cmp (l, h) == 0; + ok = SAME_SIGN (inex1, inex2) && mpfr_cmp (l, h) == 0; if (ok) mpfr_set (y, h, rnd); /* exact */ mpfr_clear (l); @@ -172,8 +238,9 @@ mpfr_clear (g); if (ok) { + MPFR_ZIV_FREE (loop); MPFR_SAVE_EXPO_FREE (expo); - return mpfr_check_range (y, inexact, rnd); + return mpfr_check_range (y, inex1, rnd); } /* since we have log|gamma(x)| = - log|x| - gamma*x + O(x^2), if x ~ 2^(-n), then we have a n-bit approximation, thus @@ -205,9 +272,10 @@ thus lngamma(x) = log(Pi*(x-1)/sin(Pi*(2-x))) - lngamma(2-x) */ w = precy + MPFR_INT_CEIL_LOG2 (precy); + w += MPFR_INT_CEIL_LOG2 (w) + 14; + MPFR_ZIV_INIT (loop, w); while (1) { - w += MPFR_INT_CEIL_LOG2 (w) + 14; MPFR_ASSERTD(w >= 3); mpfr_set_prec (s, w); mpfr_set_prec (t, w); @@ -288,7 +356,9 @@ + (rnd == MPFR_RNDN))) goto end; } + MPFR_ZIV_NEXT (loop, w); } + MPFR_ZIV_FREE (loop); } /* now z0 > 1 */ @@ -298,10 +368,10 @@ /* since k is O(w), the value of log(z0*...*(z0+k-1)) is about w*log(w), so there is a cancellation of ~log(w) in the argument reconstruction */ w = precy + MPFR_INT_CEIL_LOG2 (precy); - - do + w += MPFR_INT_CEIL_LOG2 (w) + 13; + MPFR_ZIV_INIT (loop, w); + while (1) { - w += MPFR_INT_CEIL_LOG2 (w) + 13; MPFR_ASSERTD (w >= 3); /* argument reduction: we compute gamma(z0 + k), where the series @@ -441,6 +511,15 @@ #ifdef IS_GAMMA err_s = MPFR_GET_EXP(s); mpfr_exp (s, s, MPFR_RNDN); + /* If s is +Inf, we compute exp(lngamma(z0)). */ + if (mpfr_inf_p (s)) + { + inexact = mpfr_explgamma (y, z0, &expo, s, t, rnd); + if (inexact) + goto end0; + else + goto ziv_next; + } /* before the exponential, we have s = s0 + h where |h| <= (2m+48)*ulp(s), thus exp(s0) = exp(s) * exp(-h). For |h| <= 1/4, we have |exp(h)-1| <= 1.2*|h| thus @@ -480,16 +559,26 @@ err_s = (err_t == err_s) ? 1 + err_s : ((err_t > err_s) ? err_t : err_s); err_s += 1 - MPFR_GET_EXP(s); #endif + if (MPFR_LIKELY (MPFR_CAN_ROUND (s, w - err_s, precy, rnd))) + break; +#ifdef IS_GAMMA + ziv_next: +#endif + MPFR_ZIV_NEXT (loop, w); } - while (MPFR_UNLIKELY (!MPFR_CAN_ROUND (s, w - err_s, precy, rnd))); +#ifdef IS_GAMMA + end0: +#endif oldBm = Bm; while (Bm--) mpz_clear (B[Bm]); (*__gmp_free_func) (B, oldBm * sizeof (mpz_t)); end: - inexact = mpfr_set (y, s, rnd); + if (inexact == 0) + inexact = mpfr_set (y, s, rnd); + MPFR_ZIV_FREE (loop); mpfr_clear (s); mpfr_clear (t); diff -Naurd mpfr-3.1.0-a/src/mpfr.h mpfr-3.1.0-b/src/mpfr.h --- mpfr-3.1.0-a/src/mpfr.h 2012-04-27 01:13:15.000000000 +0000 +++ mpfr-3.1.0-b/src/mpfr.h 2012-05-07 18:52:45.000000000 +0000 @@ -27,7 +27,7 @@ #define MPFR_VERSION_MAJOR 3 #define MPFR_VERSION_MINOR 1 #define MPFR_VERSION_PATCHLEVEL 0 -#define MPFR_VERSION_STRING "3.1.0-p9" +#define MPFR_VERSION_STRING "3.1.0-p10" /* Macros dealing with MPFR VERSION */ #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c)) diff -Naurd mpfr-3.1.0-a/src/version.c mpfr-3.1.0-b/src/version.c --- mpfr-3.1.0-a/src/version.c 2012-04-27 01:13:15.000000000 +0000 +++ mpfr-3.1.0-b/src/version.c 2012-05-07 18:52:45.000000000 +0000 @@ -25,5 +25,5 @@ const char * mpfr_get_version (void) { - return "3.1.0-p9"; + return "3.1.0-p10"; } diff -Naurd mpfr-3.1.0-a/tests/tgamma.c mpfr-3.1.0-b/tests/tgamma.c --- mpfr-3.1.0-a/tests/tgamma.c 2012-04-27 01:13:15.000000000 +0000 +++ mpfr-3.1.0-b/tests/tgamma.c 2012-05-07 18:52:45.000000000 +0000 @@ -838,6 +838,175 @@ exit (1); } +/* Test mpfr_gamma in precision p1 by comparing it with exp(lgamma(x)) + computing with a working precision p2. Assume that x is not an + integer <= 2. */ +static void +exp_lgamma (mpfr_t x, mpfr_prec_t p1, mpfr_prec_t p2) +{ + mpfr_t yd, yu, zd, zu; + int inexd, inexu, sign; + int underflow = -1, overflow = -1; /* -1: we don't know */ + int got_underflow, got_overflow; + + if (mpfr_integer_p (x) && mpfr_cmp_si (x, 2) <= 0) + { + printf ("Warning! x is an integer <= 2 in exp_lgamma: "); + mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN); putchar ('\n'); + return; + } + mpfr_inits2 (p2, yd, yu, (mpfr_ptr) 0); + inexd = mpfr_lgamma (yd, &sign, x, MPFR_RNDD); + mpfr_set (yu, yd, MPFR_RNDN); /* exact */ + if (inexd) + mpfr_nextabove (yu); + mpfr_clear_flags (); + mpfr_exp (yd, yd, MPFR_RNDD); + if (! mpfr_underflow_p ()) + underflow = 0; + if (mpfr_overflow_p ()) + overflow = 1; + mpfr_clear_flags (); + mpfr_exp (yu, yu, MPFR_RNDU); + if (mpfr_underflow_p ()) + underflow = 1; + if (! mpfr_overflow_p ()) + overflow = 0; + if (sign < 0) + { + mpfr_neg (yd, yd, MPFR_RNDN); /* exact */ + mpfr_neg (yu, yu, MPFR_RNDN); /* exact */ + mpfr_swap (yd, yu); + } + /* yd < Gamma(x) < yu (strict inequalities since x != 1 and x != 2) */ + mpfr_inits2 (p1, zd, zu, (mpfr_ptr) 0); + mpfr_clear_flags (); + inexd = mpfr_gamma (zd, x, MPFR_RNDD); /* zd <= Gamma(x) < yu */ + got_underflow = underflow == -1 ? -1 : !! mpfr_underflow_p (); + got_overflow = overflow == -1 ? -1 : !! mpfr_overflow_p (); + if (! mpfr_less_p (zd, yu) || inexd > 0 || + got_underflow != underflow || + got_overflow != overflow) + { + printf ("Error in exp_lgamma on x = "); + mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); putchar ('\n'); + printf ("yu = "); + mpfr_dump (yu); + printf ("zd = "); + mpfr_dump (zd); + printf ("got inexd = %d, expected <= 0\n", inexd); + printf ("got underflow = %d, expected %d\n", got_underflow, underflow); + printf ("got overflow = %d, expected %d\n", got_overflow, overflow); + exit (1); + } + mpfr_clear_flags (); + inexu = mpfr_gamma (zu, x, MPFR_RNDU); /* zu >= Gamma(x) > yd */ + got_underflow = underflow == -1 ? -1 : !! mpfr_underflow_p (); + got_overflow = overflow == -1 ? -1 : !! mpfr_overflow_p (); + if (! mpfr_greater_p (zu, yd) || inexu < 0 || + got_underflow != underflow || + got_overflow != overflow) + { + printf ("Error in exp_lgamma on x = "); + mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); putchar ('\n'); + printf ("yd = "); + mpfr_dump (yd); + printf ("zu = "); + mpfr_dump (zu); + printf ("got inexu = %d, expected >= 0\n", inexu); + printf ("got underflow = %d, expected %d\n", got_underflow, underflow); + printf ("got overflow = %d, expected %d\n", got_overflow, overflow); + exit (1); + } + if (mpfr_equal_p (zd, zu)) + { + if (inexd != 0 || inexu != 0) + { + printf ("Error in exp_lgamma on x = "); + mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); putchar ('\n'); + printf ("zd = zu, thus exact, but inexd = %d and inexu = %d\n", + inexd, inexu); + exit (1); + } + MPFR_ASSERTN (got_underflow == 0); + MPFR_ASSERTN (got_overflow == 0); + } + else if (inexd == 0 || inexu == 0) + { + printf ("Error in exp_lgamma on x = "); + mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); putchar ('\n'); + printf ("zd != zu, thus inexact, but inexd = %d and inexu = %d\n", + inexd, inexu); + exit (1); + } + mpfr_clears (yd, yu, zd, zu, (mpfr_ptr) 0); +} + +static void +exp_lgamma_tests (void) +{ + mpfr_t x; + mpfr_exp_t emin, emax; + int i; + + emin = mpfr_get_emin (); + emax = mpfr_get_emax (); + set_emin (MPFR_EMIN_MIN); + set_emax (MPFR_EMAX_MAX); + + mpfr_init2 (x, 96); + for (i = 3; i <= 8; i++) + { + mpfr_set_ui (x, i, MPFR_RNDN); + exp_lgamma (x, 53, 64); + mpfr_nextbelow (x); + exp_lgamma (x, 53, 64); + mpfr_nextabove (x); + mpfr_nextabove (x); + exp_lgamma (x, 53, 64); + } + mpfr_set_str (x, "1.7", 10, MPFR_RNDN); + exp_lgamma (x, 53, 64); + mpfr_set_str (x, "-4.6308260837372266e+07", 10, MPFR_RNDN); + exp_lgamma (x, 53, 64); + mpfr_set_str (x, "-90.6308260837372266e+15", 10, MPFR_RNDN); + exp_lgamma (x, 53, 64); + /* The following test gives a large positive result < +Inf */ + mpfr_set_str (x, "1.2b13fc45a92dea1@14", 16, MPFR_RNDN); + exp_lgamma (x, 53, 64); + /* Idem for a large negative result > -Inf */ + mpfr_set_str (x, "-1.2b13fc45a92de81@14", 16, MPFR_RNDN); + exp_lgamma (x, 53, 64); + /* The following two tests trigger an endless loop in r8186 + on 64-bit machines (64-bit exponent). The second one (due + to undetected overflow) is a direct consequence of the + first one, due to the call of Gamma(2-x) if x < 1. */ + mpfr_set_str (x, "1.2b13fc45a92dec8@14", 16, MPFR_RNDN); + exp_lgamma (x, 53, 64); + mpfr_set_str (x, "-1.2b13fc45a92dea8@14", 16, MPFR_RNDN); + exp_lgamma (x, 53, 64); + /* Similar tests (overflow threshold) for 32-bit machines. */ + mpfr_set_str (x, "2ab68d8.657542f855111c61", 16, MPFR_RNDN); + exp_lgamma (x, 12, 64); + mpfr_set_str (x, "-2ab68d6.657542f855111c61", 16, MPFR_RNDN); + exp_lgamma (x, 12, 64); + /* The following test is an overflow on 32-bit and 64-bit machines. + Revision r8189 fails on 64-bit machines as the flag is unset. */ + mpfr_set_str (x, "1.2b13fc45a92ded8@14", 16, MPFR_RNDN); + exp_lgamma (x, 53, 64); + /* On the following tests, with r8196, one gets an underflow on + 32-bit machines, while a normal result is expected (see FIXME + in gamma.c:382). */ + mpfr_set_str (x, "-2ab68d6.657542f855111c6104", 16, MPFR_RNDN); + exp_lgamma (x, 12, 64); /* failure on 32-bit machines */ + mpfr_set_str (x, "-12b13fc45a92deb.1c6c5bc964", 16, MPFR_RNDN); + exp_lgamma (x, 12, 64); /* failure on 64-bit machines */ + mpfr_clear (x); + + set_emin (emin); + set_emax (emax); +} + int main (int argc, char *argv[]) { @@ -852,6 +1021,7 @@ test20071231 (); test20100709 (); test20120426 (); + exp_lgamma_tests (); data_check ("data/gamma", mpfr_gamma, "mpfr_gamma");