| /* Utilities for MPFR developers, not exported. |
| |
| Copyright 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 Free Software Foundation, Inc. |
| Contributed by the Arenaire and Cacao projects, INRIA. |
| |
| This file is part of the GNU MPFR Library. |
| |
| The GNU MPFR Library is free software; you can redistribute it and/or modify |
| it under the terms of the GNU Lesser General Public License as published by |
| the Free Software Foundation; either version 2.1 of the License, or (at your |
| option) any later version. |
| |
| The GNU MPFR Library is distributed in the hope that it will be useful, but |
| WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY |
| or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public |
| License for more details. |
| |
| You should have received a copy of the GNU Lesser General Public License |
| along with the GNU MPFR Library; see the file COPYING.LIB. If not, write to |
| the Free Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, |
| MA 02110-1301, USA. */ |
| |
| #ifndef __MPFR_IMPL_H__ |
| #define __MPFR_IMPL_H__ |
| |
| /* Let's include some standard headers unconditionally as they are |
| already needed by several source files or when some options are |
| enabled/disabled, and it is easy to forget them (some configure |
| options may hide the error). |
| Note: If some source file must not have such a header included |
| (which is very unlikely and probably means something broken in |
| this source file), we should do that with some macro (that would |
| also force to disable incompatible features). */ |
| #if defined (__cplusplus) |
| #include <cstdio> |
| #else |
| #include <stdio.h> |
| #endif |
| #include <limits.h> |
| |
| /* Check if we are inside a build of MPFR or inside the test suite. |
| This is needed in mpfr.h to export or import the functions. |
| It matters only for Windows DLL */ |
| #ifndef __MPFR_TEST_H__ |
| # define __MPFR_WITHIN_MPFR 1 |
| #endif |
| |
| /****************************************************** |
| ****************** Include files ********************* |
| ******************************************************/ |
| |
| /* Include 'config.h' before using ANY configure macros if needed |
| NOTE: It isn't MPFR 'config.h', but GMP's one! */ |
| #ifdef HAVE_CONFIG_H |
| # include "config.h" |
| #endif |
| |
| #ifdef MPFR_HAVE_GMP_IMPL /* Build with gmp internals*/ |
| |
| # ifndef __GMP_H__ |
| # include "gmp.h" |
| # endif |
| # ifndef __GMP_IMPL_H__ |
| # include "gmp-impl.h" |
| # endif |
| # ifdef MPFR_NEED_LONGLONG_H |
| # include "longlong.h" |
| # endif |
| # ifndef __MPFR_H |
| # include "mpfr.h" |
| # endif |
| |
| #else /* Build without gmp internals */ |
| |
| # ifndef __GMP_H__ |
| # include "gmp.h" |
| # endif |
| # ifndef __MPFR_H |
| # include "mpfr.h" |
| # endif |
| # ifndef __GMPFR_GMP_H__ |
| # include "mpfr-gmp.h" |
| # endif |
| # ifdef MPFR_NEED_LONGLONG_H |
| # include "mpfr-longlong.h" |
| # endif |
| |
| #endif |
| #undef MPFR_NEED_LONGLONG_H |
| |
| /* For the definition of MPFR_THREAD_ATTR. GCC/ICC detection macros are |
| no longer used, as they sometimes gave incorrect information about |
| the support of thread-local variables. A configure check is now done. |
| If the use of detection macros is needed in the future, this could be |
| moved below (after the detection macros are defined). */ |
| #include "mpfr-thread.h" |
| |
| |
| /****************************************************** |
| ***************** Detection macros ******************* |
| ******************************************************/ |
| |
| /* Macros to detect STDC, GCC, GLIBC, GMP and ICC version */ |
| #if defined(__STDC_VERSION__) |
| # define __MPFR_STDC(version) (__STDC_VERSION__>=(version)) |
| #elif defined(__STDC__) |
| # define __MPFR_STDC(version) (0 == (version)) |
| #else |
| # define __MPFR_STDC(version) 0 |
| #endif |
| |
| #if defined(__ICC) |
| # define __MPFR_ICC(a,b,c) (__ICC >= (a)*100+(b)*10+(c)) |
| #elif defined(__INTEL_COMPILER) |
| # define __MPFR_ICC(a,b,c) (__INTEL_COMPILER >= (a)*100+(b)*10+(c)) |
| #else |
| # define __MPFR_ICC(a,b,c) 0 |
| #endif |
| |
| #if defined(__GNUC__) && defined(__GNUC_MINOR__) && ! __MPFR_ICC(0,0,0) |
| # define __MPFR_GNUC(a,i) \ |
| (MPFR_VERSION_NUM(__GNUC__,__GNUC_MINOR__,0) >= MPFR_VERSION_NUM(a,i,0)) |
| #else |
| # define __MPFR_GNUC(a,i) 0 |
| #endif |
| |
| #if defined(__GLIBC__) && defined(__GLIBC_MINOR__) |
| # define __MPFR_GLIBC(a,i) \ |
| (MPFR_VERSION_NUM(__GLIBC__,__GLIBC_MINOR__,0) >= MPFR_VERSION_NUM(a,i,0)) |
| #else |
| # define __MPFR_GLIBC(a,i) 0 |
| #endif |
| |
| #if defined(__GNU_MP_VERSION) && \ |
| defined(__GNU_MP_VERSION_MINOR) && \ |
| defined(__GNU_MP_VERSION_PATCHLEVEL) |
| # define __MPFR_GMP(a,b,c) \ |
| (MPFR_VERSION_NUM(__GNU_MP_VERSION,__GNU_MP_VERSION_MINOR,__GNU_MP_VERSION_PATCHLEVEL) >= MPFR_VERSION_NUM(a,b,c)) |
| #else |
| # define __MPFR_GMP(a,b,c) 0 |
| #endif |
| |
| |
| |
| /****************************************************** |
| ******************** Check GMP *********************** |
| ******************************************************/ |
| |
| #if !__MPFR_GMP(4,1,0) |
| # error "GMP 4.1.0 or newer needed" |
| #endif |
| |
| #if GMP_NAIL_BITS != 0 |
| # error "MPFR doesn't support nonzero values of GMP_NAIL_BITS" |
| #endif |
| |
| #if (BITS_PER_MP_LIMB<32) || (BITS_PER_MP_LIMB & (BITS_PER_MP_LIMB - 1)) |
| # error "BITS_PER_MP_LIMB must be a power of 2, and >= 32" |
| #endif |
| |
| #if BITS_PER_MP_LIMB == 16 |
| # define MPFR_LOG2_BITS_PER_MP_LIMB 4 |
| #elif BITS_PER_MP_LIMB == 32 |
| # define MPFR_LOG2_BITS_PER_MP_LIMB 5 |
| #elif BITS_PER_MP_LIMB == 64 |
| # define MPFR_LOG2_BITS_PER_MP_LIMB 6 |
| #elif BITS_PER_MP_LIMB == 128 |
| # define MPFR_LOG2_BITS_PER_MP_LIMB 7 |
| #elif BITS_PER_MP_LIMB == 256 |
| # define MPFR_LOG2_BITS_PER_MP_LIMB 8 |
| #else |
| # error "Can't compute log2(BITS_PER_MP_LIMB)" |
| #endif |
| |
| #if __MPFR_GNUC(3,0) || __MPFR_ICC(8,1,0) |
| # define MPFR_NORETURN_ATTR __attribute__ ((noreturn)) |
| # define MPFR_CONST_ATTR __attribute__ ((const)) |
| #else |
| # define MPFR_NORETURN_ATTR |
| # define MPFR_CONST_ATTR |
| #endif |
| |
| /****************************************************** |
| ************* Global Internal Variables ************** |
| ******************************************************/ |
| |
| /* Cache struct */ |
| struct __gmpfr_cache_s { |
| mpfr_t x; |
| int inexact; |
| int (*func)(mpfr_ptr, mpfr_rnd_t); |
| }; |
| typedef struct __gmpfr_cache_s mpfr_cache_t[1]; |
| |
| #if defined (__cplusplus) |
| extern "C" { |
| #endif |
| |
| __MPFR_DECLSPEC extern MPFR_THREAD_ATTR unsigned int __gmpfr_flags; |
| __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mp_exp_t __gmpfr_emin; |
| __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mp_exp_t __gmpfr_emax; |
| __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mp_prec_t __gmpfr_default_fp_bit_precision; |
| __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_rnd_t __gmpfr_default_rounding_mode; |
| __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_pi; |
| __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_log2; |
| __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_euler; |
| __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_catalan; |
| |
| #define BASE_MAX 36 |
| __MPFR_DECLSPEC extern const __mpfr_struct __gmpfr_l2b[BASE_MAX-1][2]; |
| |
| /* Note: do not use the following values when they can be outside the |
| current exponent range, e.g. when the exponent range has not been |
| extended yet; under such a condition, they can be used only in |
| mpfr_cmpabs. */ |
| __MPFR_DECLSPEC extern const mpfr_t __gmpfr_one; |
| __MPFR_DECLSPEC extern const mpfr_t __gmpfr_two; |
| __MPFR_DECLSPEC extern const mpfr_t __gmpfr_four; |
| |
| |
| #if defined (__cplusplus) |
| } |
| #endif |
| |
| /* Flags of __gmpfr_flags */ |
| #define MPFR_FLAGS_UNDERFLOW 1 |
| #define MPFR_FLAGS_OVERFLOW 2 |
| #define MPFR_FLAGS_NAN 4 |
| #define MPFR_FLAGS_INEXACT 8 |
| #define MPFR_FLAGS_ERANGE 16 |
| #define MPFR_FLAGS_ALL 31 |
| |
| /* Replace some commun functions for direct access to the global vars */ |
| #define mpfr_get_emin() (__gmpfr_emin + 0) |
| #define mpfr_get_emax() (__gmpfr_emax + 0) |
| #define mpfr_get_default_rounding_mode() (__gmpfr_default_rounding_mode + 0) |
| #define mpfr_get_default_prec() (__gmpfr_default_fp_bit_precision + 0) |
| |
| #define mpfr_clear_flags() \ |
| ((void) (__gmpfr_flags = 0)) |
| #define mpfr_clear_underflow() \ |
| ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_UNDERFLOW)) |
| #define mpfr_clear_overflow() \ |
| ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_OVERFLOW)) |
| #define mpfr_clear_nanflag() \ |
| ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_NAN)) |
| #define mpfr_clear_inexflag() \ |
| ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_INEXACT)) |
| #define mpfr_clear_erangeflag() \ |
| ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_ERANGE)) |
| #define mpfr_underflow_p() \ |
| ((int) (__gmpfr_flags & MPFR_FLAGS_UNDERFLOW)) |
| #define mpfr_overflow_p() \ |
| ((int) (__gmpfr_flags & MPFR_FLAGS_OVERFLOW)) |
| #define mpfr_nanflag_p() \ |
| ((int) (__gmpfr_flags & MPFR_FLAGS_NAN)) |
| #define mpfr_inexflag_p() \ |
| ((int) (__gmpfr_flags & MPFR_FLAGS_INEXACT)) |
| #define mpfr_erangeflag_p() \ |
| ((int) (__gmpfr_flags & MPFR_FLAGS_ERANGE)) |
| |
| /* Testing an exception flag correctly is tricky. There are mainly two |
| pitfalls: First, one needs to remember to clear the corresponding |
| flag, in case it was set before the function call or during some |
| intermediate computations (in practice, one can clear all the flags). |
| Secondly, one needs to test the flag early enough, i.e. before it |
| can be modified by another function. Moreover, it is quite difficult |
| (if not impossible) to reliably check problems with "make check". To |
| avoid these pitfalls, it is recommended to use the following macros. |
| Other use of the exception-flag predicate functions/macros will be |
| detected by mpfrlint. |
| Note: _op can be either a statement or an expression. |
| MPFR_BLOCK_EXCEP should be used only inside a block; it is useful to |
| detect some exception in order to exit the block as soon as possible. */ |
| #define MPFR_BLOCK_DECL(_flags) unsigned int _flags |
| #define MPFR_BLOCK(_flags,_op) \ |
| do \ |
| { \ |
| mpfr_clear_flags (); \ |
| _op; \ |
| (_flags) = __gmpfr_flags; \ |
| } \ |
| while (0) |
| #define MPFR_BLOCK_TEST(_flags,_f) MPFR_UNLIKELY ((_flags) & (_f)) |
| #define MPFR_BLOCK_EXCEP (__gmpfr_flags & (MPFR_FLAGS_UNDERFLOW | \ |
| MPFR_FLAGS_OVERFLOW | \ |
| MPFR_FLAGS_NAN)) |
| /* Let's use a MPFR_ prefix, because e.g. OVERFLOW is defined by glibc's |
| math.h, though this is not a reserved identifier! */ |
| #define MPFR_UNDERFLOW(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_UNDERFLOW) |
| #define MPFR_OVERFLOW(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_OVERFLOW) |
| #define MPFR_NANFLAG(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_NAN) |
| #define MPFR_INEXFLAG(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_INEXACT) |
| #define MPFR_ERANGEFLAG(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_ERANGE) |
| |
| |
| /****************************************************** |
| ******************** Assertions ********************** |
| ******************************************************/ |
| |
| /* Compile with -DWANT_ASSERT to check all assert statements */ |
| |
| /* Note: do not use GMP macros ASSERT_ALWAYS and ASSERT as they are not |
| expressions, and as a consequence, they cannot be used in a for(), |
| with a comma operator and so on. */ |
| |
| /* MPFR_ASSERTN(expr): assertions that should always be checked */ |
| #define MPFR_ASSERTN(expr) \ |
| ((void) ((MPFR_UNLIKELY(expr)) || MPFR_UNLIKELY( (ASSERT_FAIL(expr),0) ))) |
| |
| /* MPFR_ASSERTD(expr): assertions that should be checked when testing */ |
| #ifdef WANT_ASSERT |
| # define MPFR_EXP_CHECK 1 |
| # define MPFR_ASSERTD(expr) MPFR_ASSERTN (expr) |
| #else |
| # define MPFR_ASSERTD(expr) ((void) 0) |
| #endif |
| |
| /* Code to deal with impossible |
| WARNING: It doesn't use do { } while (0) for Insure++*/ |
| #define MPFR_RET_NEVER_GO_HERE() {MPFR_ASSERTN(0); return 0;} |
| |
| |
| /****************************************************** |
| ******************** Warnings ************************ |
| ******************************************************/ |
| |
| /* MPFR_WARNING is no longer useful, but let's keep the macro in case |
| it needs to be used again in the future. */ |
| |
| #ifdef MPFR_USE_WARNINGS |
| # include <stdlib.h> |
| # define MPFR_WARNING(W) \ |
| do \ |
| { \ |
| char *q = getenv ("MPFR_QUIET"); \ |
| if (q == NULL || *q == 0) \ |
| fprintf (stderr, "MPFR: %s\n", W); \ |
| } \ |
| while (0) |
| #else |
| # define MPFR_WARNING(W) ((void) 0) |
| #endif |
| |
| |
| /****************************************************** |
| ****************** double macros ********************* |
| ******************************************************/ |
| |
| /* Definition of constants */ |
| #define LOG2 0.69314718055994528622 /* log(2) rounded to zero on 53 bits */ |
| #define ALPHA 4.3191365662914471407 /* a+2 = a*log(a), rounded to +infinity */ |
| #define EXPM1 0.36787944117144227851 /* exp(-1), rounded to zero */ |
| |
| /* MPFR_DOUBLE_SPEC = 1 if the C type 'double' corresponds to IEEE-754 |
| double precision, 0 if it doesn't, and undefined if one doesn't know. |
| On all the tested machines, MPFR_DOUBLE_SPEC = 1. To have this macro |
| defined here, #include <float.h> is needed. If need be, other values |
| could be defined for other specs (once they are known). */ |
| #if !defined(MPFR_DOUBLE_SPEC) && defined(FLT_RADIX) && \ |
| defined(DBL_MANT_DIG) && defined(DBL_MIN_EXP) && defined(DBL_MAX_EXP) |
| # if FLT_RADIX == 2 && DBL_MANT_DIG == 53 && \ |
| DBL_MIN_EXP == -1021 && DBL_MAX_EXP == 1024 |
| # define MPFR_DOUBLE_SPEC 1 |
| # else |
| # define MPFR_DOUBLE_SPEC 0 |
| # endif |
| #endif |
| |
| /* Debug non IEEE floats */ |
| #ifdef XDEBUG |
| # undef _GMP_IEEE_FLOATS |
| #endif |
| #ifndef _GMP_IEEE_FLOATS |
| # define _GMP_IEEE_FLOATS 0 |
| #endif |
| |
| #ifndef IEEE_DBL_MANT_DIG |
| #define IEEE_DBL_MANT_DIG 53 |
| #endif |
| #define MPFR_LIMBS_PER_DOUBLE ((IEEE_DBL_MANT_DIG-1)/BITS_PER_MP_LIMB+1) |
| |
| /* Visual C++ doesn't support +1.0/.00, -1.0/0.0 and 0.0/0.0 |
| at compile time. */ |
| #if defined(_MSC_VER) && defined(_WIN32) && (_MSC_VER >= 1200) |
| static double double_zero = 0.0; |
| # define DBL_NAN (double_zero/double_zero) |
| # define DBL_POS_INF ((double) 1.0/double_zero) |
| # define DBL_NEG_INF ((double)-1.0/double_zero) |
| # define DBL_NEG_ZERO (-double_zero) |
| #else |
| # define DBL_POS_INF ((double) 1.0/0.0) |
| # define DBL_NEG_INF ((double)-1.0/0.0) |
| # define DBL_NAN ((double) 0.0/0.0) |
| # define DBL_NEG_ZERO (-0.0) |
| #endif |
| |
| /* Note: the argument x must be a lvalue of type double. */ |
| #if _GMP_IEEE_FLOATS |
| typedef union ieee_double_extract Ieee_double_extract; |
| |
| # define DOUBLE_ISNANorINF(x) (((Ieee_double_extract *)&(x))->s.exp == 0x7ff) |
| # define DOUBLE_ISINF(x) (DOUBLE_ISNANorINF(x) && \ |
| (((Ieee_double_extract *)&(x))->s.manl == 0) && \ |
| (((Ieee_double_extract *)&(x))->s.manh == 0)) |
| # define DOUBLE_ISNAN(x) (DOUBLE_ISNANorINF(x) && \ |
| ((((Ieee_double_extract *)&(x))->s.manl != 0) || \ |
| (((Ieee_double_extract *)&(x))->s.manh != 0))) |
| #else |
| /* Below, the &(x) == &(x) || &(x) != &(x) allows to make sure that x |
| is a lvalue without (probably) any warning from the compiler. The |
| &(x) != &(x) is needed to avoid a failure under Mac OS X 10.4.11 |
| (with Xcode 2.4.1, i.e. the latest one). */ |
| # define LVALUE(x) (&(x) == &(x) || &(x) != &(x)) |
| # define DOUBLE_ISINF(x) (LVALUE(x) && ((x) > DBL_MAX || (x) < -DBL_MAX)) |
| # ifdef MPFR_NANISNAN |
| /* Avoid MIPSpro / IRIX64 / gcc -ffast-math (incorrect) optimizations. |
| The + must not be replaced by a ||. With gcc -ffast-math, NaN is |
| regarded as a positive number or something like that; the second |
| test catches this case. */ |
| # define DOUBLE_ISNAN(x) \ |
| (LVALUE(x) && !((((x) >= 0.0) + ((x) <= 0.0)) && -(x)*(x) <= 0.0)) |
| # else |
| # define DOUBLE_ISNAN(x) (LVALUE(x) && (x) != (x)) |
| # endif |
| #endif |
| |
| /****************************************************** |
| *************** Long double macros ******************* |
| ******************************************************/ |
| |
| /* We try to get the exact value of the precision of long double |
| (provided by the implementation) in order to provide correct |
| rounding in this case (not guaranteed if the C implementation |
| does not have an adequate long double arithmetic). Note that |
| it may be lower than the precision of some numbers that can |
| be represented in a long double; e.g. on FreeBSD/x86, it is |
| 53 because the processor is configured to round in double |
| precision (even when using the long double type -- this is a |
| limitation of the x87 arithmetic), and on Mac OS X, it is 106 |
| because the implementation is a double-double arithmetic. |
| Otherwise (e.g. in base 10), we get an upper bound of the |
| precision, and correct rounding isn't currently provided. |
| */ |
| #if defined(LDBL_MANT_DIG) && FLT_RADIX == 2 |
| # define MPFR_LDBL_MANT_DIG LDBL_MANT_DIG |
| #else |
| # define MPFR_LDBL_MANT_DIG \ |
| (sizeof(long double)*BITS_PER_MP_LIMB/sizeof(mp_limb_t)) |
| #endif |
| #define MPFR_LIMBS_PER_LONG_DOUBLE \ |
| ((sizeof(long double)-1)/sizeof(mp_limb_t)+1) |
| |
| /* LONGDOUBLE_NAN_ACTION executes the code "action" if x is a NaN. */ |
| |
| /* On hppa2.0n-hp-hpux10 with the unbundled HP cc, the test x!=x on a NaN |
| has been seen false, meaning NaNs are not detected. This seemed to |
| happen only after other comparisons, not sure what's really going on. In |
| any case we can pick apart the bytes to identify a NaN. */ |
| #ifdef HAVE_LDOUBLE_IEEE_QUAD_BIG |
| # define LONGDOUBLE_NAN_ACTION(x, action) \ |
| do { \ |
| union { \ |
| long double ld; \ |
| struct { \ |
| unsigned int sign : 1; \ |
| unsigned int exp : 15; \ |
| unsigned int man3 : 16; \ |
| unsigned int man2 : 32; \ |
| unsigned int man1 : 32; \ |
| unsigned int man0 : 32; \ |
| } s; \ |
| } u; \ |
| u.ld = (x); \ |
| if (u.s.exp == 0x7FFFL \ |
| && (u.s.man0 | u.s.man1 | u.s.man2 | u.s.man3) != 0) \ |
| { action; } \ |
| } while (0) |
| #endif |
| |
| /* Under IEEE rules, NaN is not equal to anything, including itself. |
| "volatile" here stops "cc" on mips64-sgi-irix6.5 from optimizing away |
| x!=x. */ |
| #ifndef LONGDOUBLE_NAN_ACTION |
| # define LONGDOUBLE_NAN_ACTION(x, action) \ |
| do { \ |
| volatile long double __x = LONGDOUBLE_VOLATILE (x); \ |
| if ((x) != __x) \ |
| { action; } \ |
| } while (0) |
| # define WANT_LONGDOUBLE_VOLATILE 1 |
| #endif |
| |
| /* If we don't have a proper "volatile" then volatile is #defined to empty, |
| in this case call through an external function to stop the compiler |
| optimizing anything. */ |
| #ifdef WANT_LONGDOUBLE_VOLATILE |
| # ifdef volatile |
| __MPFR_DECLSPEC long double __gmpfr_longdouble_volatile _MPFR_PROTO ((long double)) MPFR_CONST_ATTR; |
| # define LONGDOUBLE_VOLATILE(x) (__gmpfr_longdouble_volatile (x)) |
| # define WANT_GMPFR_LONGDOUBLE_VOLATILE 1 |
| # else |
| # define LONGDOUBLE_VOLATILE(x) (x) |
| # endif |
| #endif |
| |
| /* Some special case for IEEE_EXT Litle Endian */ |
| #if HAVE_LDOUBLE_IEEE_EXT_LITTLE |
| |
| typedef union { |
| long double ld; |
| struct { |
| unsigned int manl : 32; |
| unsigned int manh : 32; |
| unsigned int expl : 8 ; |
| unsigned int exph : 7; |
| unsigned int sign : 1; |
| } s; |
| } mpfr_long_double_t; |
| |
| /* #undef MPFR_LDBL_MANT_DIG */ |
| #undef MPFR_LIMBS_PER_LONG_DOUBLE |
| /* #define MPFR_LDBL_MANT_DIG 64 */ |
| #define MPFR_LIMBS_PER_LONG_DOUBLE ((64-1)/BITS_PER_MP_LIMB+1) |
| |
| #endif |
| |
| /****************************************************** |
| *************** _Decimal64 support ******************* |
| ******************************************************/ |
| |
| #ifdef MPFR_WANT_DECIMAL_FLOATS |
| /* to cast between binary64 and decimal64 */ |
| union ieee_double_decimal64 { double d; _Decimal64 d64; }; |
| #endif |
| |
| /****************************************************** |
| **************** mpfr_t properties ******************* |
| ******************************************************/ |
| |
| #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)/BITS_PER_MP_LIMB+1) |
| |
| #if _MPFR_PREC_FORMAT == 1 |
| # define MPFR_INTPREC_MAX (USHRT_MAX & ~(unsigned int) (BITS_PER_MP_LIMB - 1)) |
| #elif _MPFR_PREC_FORMAT == 2 |
| # define MPFR_INTPREC_MAX (UINT_MAX & ~(unsigned int) (BITS_PER_MP_LIMB - 1)) |
| #elif _MPFR_PREC_FORMAT == 3 |
| # define MPFR_INTPREC_MAX (ULONG_MAX & ~(unsigned long) (BITS_PER_MP_LIMB - 1)) |
| #else |
| # error "Invalid MPFR Prec format" |
| #endif |
| |
| |
| |
| /****************************************************** |
| ***************** exponent limits ******************** |
| ******************************************************/ |
| |
| /* Define limits and unsigned type of exponent. The following definitions |
| * depend on mp_exp_t; if this type changes in GMP, these definitions will |
| * need to be modified (alternatively, a mpfr_exp_t type could be defined). |
| */ |
| #if __GMP_MP_SIZE_T_INT == 1 |
| typedef unsigned int mpfr_uexp_t; |
| # define MPFR_EXP_MAX (INT_MAX) |
| # define MPFR_EXP_MIN (INT_MIN) |
| #else |
| typedef unsigned long int mpfr_uexp_t; |
| # define MPFR_EXP_MAX (LONG_MAX) |
| # define MPFR_EXP_MIN (LONG_MIN) |
| #endif |
| #ifndef mp_exp_unsigned_t |
| # define mp_exp_unsigned_t mpfr_uexp_t |
| #endif |
| |
| #if MPFR_EXP_MIN >= LONG_MIN && MPFR_EXP_MAX <= LONG_MAX |
| typedef long int mpfr_eexp_t; |
| # define mpfr_get_exp_t(x,r) mpfr_get_si((x),(r)) |
| # define mpfr_set_exp_t(x,e,r) mpfr_set_si((x),(e),(r)) |
| #elif defined (_MPFR_H_HAVE_INTMAX_T) |
| typedef intmax_t mpfr_eexp_t; |
| # define mpfr_get_exp_t(x,r) mpfr_get_sj((x),(r)) |
| # define mpfr_set_exp_t(x,e,r) mpfr_set_sj((x),(e),(r)) |
| #else |
| # error "Cannot define mpfr_get_exp_t and mpfr_set_exp_t" |
| #endif |
| |
| /* Invalid exponent value (to track bugs...) */ |
| #define MPFR_EXP_INVALID \ |
| ((mp_exp_t) 1 << (BITS_PER_MP_LIMB*sizeof(mp_exp_t)/sizeof(mp_limb_t)-2)) |
| |
| /* Definition of the exponent limits for MPFR numbers. |
| * These limits are chosen so that if e is such an exponent, then 2e-1 and |
| * 2e+1 are representable. This is useful for intermediate computations, |
| * in particular the multiplication. |
| */ |
| #undef MPFR_EMIN_MIN |
| #undef MPFR_EMIN_MAX |
| #undef MPFR_EMAX_MIN |
| #undef MPFR_EMAX_MAX |
| #define MPFR_EMIN_MIN (1-MPFR_EXP_INVALID) |
| #define MPFR_EMIN_MAX (MPFR_EXP_INVALID-1) |
| #define MPFR_EMAX_MIN (1-MPFR_EXP_INVALID) |
| #define MPFR_EMAX_MAX (MPFR_EXP_INVALID-1) |
| |
| /* Use MPFR_GET_EXP and MPFR_SET_EXP instead of MPFR_EXP directly, |
| unless when the exponent may be out-of-range, for instance when |
| setting the exponent before calling mpfr_check_range. |
| MPFR_EXP_CHECK is defined when WANT_ASSERT is defined, but if you |
| don't use WANT_ASSERT (for speed reasons), you can still define |
| MPFR_EXP_CHECK by setting -DMPFR_EXP_CHECK in $CFLAGS. */ |
| |
| #ifdef MPFR_EXP_CHECK |
| # define MPFR_GET_EXP(x) (mpfr_get_exp) (x) |
| # define MPFR_SET_EXP(x, exp) MPFR_ASSERTN (!mpfr_set_exp ((x), (exp))) |
| # define MPFR_SET_INVALID_EXP(x) ((void) (MPFR_EXP (x) = MPFR_EXP_INVALID)) |
| #else |
| # define MPFR_GET_EXP(x) MPFR_EXP (x) |
| # define MPFR_SET_EXP(x, exp) ((void) (MPFR_EXP (x) = (exp))) |
| # define MPFR_SET_INVALID_EXP(x) ((void) 0) |
| #endif |
| |
| |
| |
| /****************************************************** |
| ********** Singular Values (NAN, INF, ZERO) ********** |
| ******************************************************/ |
| |
| /* |
| * Clear flags macros are still defined and should be still used |
| * since the functions must not assume the internal format. |
| * How to deal with special values ? |
| * 1. Check if is a special value (Zero, Nan, Inf) wiht MPFR_IS_SINGULAR |
| * 2. Deal with the special value with MPFR_IS_NAN, MPFR_IS_INF, etc |
| * 3. Else clear the flags of the dest (it must be done after since src |
| * may be also the dest!) |
| * MPFR_SET_INF, MPFR_SET_NAN, MPFR_SET_ZERO must clear by |
| * themselves the other flags. |
| */ |
| |
| /* Enum special value of exponent.*/ |
| # define MPFR_EXP_ZERO (MPFR_EXP_MIN+1) |
| # define MPFR_EXP_NAN (MPFR_EXP_MIN+2) |
| # define MPFR_EXP_INF (MPFR_EXP_MIN+3) |
| |
| #define MPFR_CLEAR_FLAGS(x) |
| |
| #define MPFR_IS_NAN(x) (MPFR_EXP(x) == MPFR_EXP_NAN) |
| #define MPFR_SET_NAN(x) (MPFR_EXP(x) = MPFR_EXP_NAN) |
| #define MPFR_IS_INF(x) (MPFR_EXP(x) == MPFR_EXP_INF) |
| #define MPFR_SET_INF(x) (MPFR_EXP(x) = MPFR_EXP_INF) |
| #define MPFR_IS_ZERO(x) (MPFR_EXP(x) == MPFR_EXP_ZERO) |
| #define MPFR_SET_ZERO(x) (MPFR_EXP(x) = MPFR_EXP_ZERO) |
| #define MPFR_NOTZERO(x) (MPFR_EXP(x) != MPFR_EXP_ZERO) |
| |
| #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)) |
| |
| #define MPFR_ARE_SINGULAR(x,y) \ |
| (MPFR_UNLIKELY(MPFR_IS_SINGULAR(x)) || MPFR_UNLIKELY(MPFR_IS_SINGULAR(y))) |
| |
| |
| |
| /****************************************************** |
| ********************* Sign Macros ******************** |
| ******************************************************/ |
| |
| #define MPFR_SIGN_POS (1) |
| #define MPFR_SIGN_NEG (-1) |
| |
| #define MPFR_IS_STRICTPOS(x) (MPFR_NOTZERO((x)) && MPFR_SIGN(x) > 0) |
| #define MPFR_IS_STRICTNEG(x) (MPFR_NOTZERO((x)) && MPFR_SIGN(x) < 0) |
| |
| #define MPFR_IS_NEG(x) (MPFR_SIGN(x) < 0) |
| #define MPFR_IS_POS(x) (MPFR_SIGN(x) > 0) |
| |
| #define MPFR_SET_POS(x) (MPFR_SIGN(x) = MPFR_SIGN_POS) |
| #define MPFR_SET_NEG(x) (MPFR_SIGN(x) = MPFR_SIGN_NEG) |
| |
| #define MPFR_CHANGE_SIGN(x) (MPFR_SIGN(x) = -MPFR_SIGN(x)) |
| #define MPFR_SET_SAME_SIGN(x, y) (MPFR_SIGN(x) = MPFR_SIGN(y)) |
| #define MPFR_SET_OPPOSITE_SIGN(x, y) (MPFR_SIGN(x) = -MPFR_SIGN(y)) |
| #define MPFR_ASSERT_SIGN(s) \ |
| (MPFR_ASSERTD((s) == MPFR_SIGN_POS || (s) == MPFR_SIGN_NEG)) |
| #define MPFR_SET_SIGN(x, s) \ |
| (MPFR_ASSERT_SIGN(s), MPFR_SIGN(x) = s) |
| #define MPFR_IS_POS_SIGN(s1) (s1 > 0) |
| #define MPFR_IS_NEG_SIGN(s1) (s1 < 0) |
| #define MPFR_MULT_SIGN(s1, s2) ((s1) * (s2)) |
| /* Transform a sign to 1 or -1 */ |
| #define MPFR_FROM_SIGN_TO_INT(s) (s) |
| #define MPFR_INT_SIGN(x) MPFR_FROM_SIGN_TO_INT(MPFR_SIGN(x)) |
| |
| |
| |
| /****************************************************** |
| ***************** Ternary Value Macros *************** |
| ******************************************************/ |
| |
| /* Special inexact value */ |
| #define MPFR_EVEN_INEX 2 |
| |
| /* When returning the ternary inexact value, ALWAYS use one of the |
| following two macros, unless the flag comes from another function |
| returning the ternary inexact value */ |
| #define MPFR_RET(I) return \ |
| (I) ? ((__gmpfr_flags |= MPFR_FLAGS_INEXACT), (I)) : 0 |
| #define MPFR_RET_NAN return (__gmpfr_flags |= MPFR_FLAGS_NAN), 0 |
| |
| #define MPFR_SET_ERANGE() (__gmpfr_flags |= MPFR_FLAGS_ERANGE) |
| |
| #define SIGN(I) ((I) < 0 ? -1 : (I) > 0) |
| #define SAME_SIGN(I1,I2) (SIGN (I1) == SIGN (I2)) |
| |
| |
| |
| /****************************************************** |
| ************** Rounding mode macros ***************** |
| ******************************************************/ |
| |
| /* We want to test this : |
| * (rnd == GMP_RNDU && test) || (rnd == RNDD && !test) |
| * ie it transforms RNDU or RNDD to Away or Zero according to the sign */ |
| #define MPFR_IS_RNDUTEST_OR_RNDDNOTTEST(rnd, test) \ |
| (((rnd) + (test)) == GMP_RNDD) |
| |
| /* We want to test if rnd = Zero, or Away. |
| 'test' is true iff negative. */ |
| #define MPFR_IS_LIKE_RNDZ(rnd, test) \ |
| ((rnd==GMP_RNDZ) || MPFR_IS_RNDUTEST_OR_RNDDNOTTEST (rnd, test)) |
| |
| /* Invert a rounding mode */ |
| #define MPFR_INVERT_RND(rnd) ((rnd == GMP_RNDU) ? GMP_RNDD : \ |
| ((rnd == GMP_RNDD) ? GMP_RNDU : rnd)) |
| |
| /* Transform RNDU and RNDD to RNDA or RNDZ */ |
| #define MPFR_UPDATE_RND_MODE(rnd, test) \ |
| do { \ |
| if (MPFR_UNLIKELY(MPFR_IS_RNDUTEST_OR_RNDDNOTTEST(rnd, test))) \ |
| rnd = GMP_RNDZ; \ |
| } while (0) |
| |
| |
| |
| /****************************************************** |
| ******************* Limb Macros ********************** |
| ******************************************************/ |
| |
| /* Definition of MPFR_LIMB_HIGHBIT */ |
| #if defined(GMP_LIMB_HIGHBIT) |
| # define MPFR_LIMB_HIGHBIT GMP_LIMB_HIGHBIT |
| #elif defined(MP_LIMB_T_HIGHBIT) |
| # define MPFR_LIMB_HIGHBIT MP_LIMB_T_HIGHBIT |
| #else |
| # error "Neither GMP_LIMB_HIGHBIT nor MP_LIMB_T_HIGHBIT defined in GMP" |
| #endif |
| |
| /* Mask to get the Most Significant Bit of a limb */ |
| #define MPFR_LIMB_MSB(l) ((l)&MPFR_LIMB_HIGHBIT) |
| |
| /* Definition of MPFR_LIMB_ONE & MPFR_LIMB_ZERO*/ |
| #ifdef CNST_LIMB |
| # define MPFR_LIMB_ONE CNST_LIMB(1) |
| # define MPFR_LIMB_ZERO CNST_LIMB(0) |
| #else |
| # define MPFR_LIMB_ONE ((mp_limb_t) 1L) |
| # define MPFR_LIMB_ZERO ((mp_limb_t) 0L) |
| #endif |
| |
| /* Mask for the low 's' bits of a limb */ |
| #define MPFR_LIMB_MASK(s) ((MPFR_LIMB_ONE<<(s))-MPFR_LIMB_ONE) |
| |
| |
| |
| /****************************************************** |
| ********************** Memory ************************ |
| ******************************************************/ |
| |
| /* Heap Memory gestion */ |
| typedef union { mp_size_t s; mp_limb_t l; } mpfr_size_limb_t; |
| #define MPFR_GET_ALLOC_SIZE(x) \ |
| ( ((mp_size_t*) MPFR_MANT(x))[-1] + 0) |
| #define MPFR_SET_ALLOC_SIZE(x, n) \ |
| ( ((mp_size_t*) MPFR_MANT(x))[-1] = n) |
| #define MPFR_MALLOC_SIZE(s) \ |
| ( sizeof(mpfr_size_limb_t) + BYTES_PER_MP_LIMB * ((size_t) s) ) |
| #define MPFR_SET_MANT_PTR(x,p) \ |
| (MPFR_MANT(x) = (mp_limb_t*) ((mpfr_size_limb_t*) p + 1)) |
| #define MPFR_GET_REAL_PTR(x) \ |
| ((mp_limb_t*) ((mpfr_size_limb_t*) MPFR_MANT(x) - 1)) |
| |
| /* Temporary memory gestion */ |
| #ifndef TMP_SALLOC |
| /* GMP 4.1.x or below or internals */ |
| #define MPFR_TMP_DECL TMP_DECL |
| #define MPFR_TMP_MARK TMP_MARK |
| #define MPFR_TMP_ALLOC TMP_ALLOC |
| #define MPFR_TMP_FREE TMP_FREE |
| #else |
| #define MPFR_TMP_DECL(x) TMP_DECL |
| #define MPFR_TMP_MARK(x) TMP_MARK |
| #define MPFR_TMP_ALLOC(s) TMP_ALLOC(s) |
| #define MPFR_TMP_FREE(x) TMP_FREE |
| #endif |
| |
| /* This code is experimental: don't use it */ |
| #ifdef MPFR_USE_OWN_MPFR_TMP_ALLOC |
| extern unsigned char *mpfr_stack; |
| #undef MPFR_TMP_DECL |
| #undef MPFR_TMP_MARK |
| #undef MPFR_TMP_ALLOC |
| #undef MPFR_TMP_FREE |
| #define MPFR_TMP_DECL(_x) unsigned char *(_x) |
| #define MPFR_TMP_MARK(_x) ((_x) = mpfr_stack) |
| #define MPFR_TMP_ALLOC(_s) (mpfr_stack += (_s), mpfr_stack - (_s)) |
| #define MPFR_TMP_FREE(_x) (mpfr_stack = (_x)) |
| #endif |
| |
| /* temporary allocate 1 limb at xp, and initialize mpfr variable x */ |
| /* The temporary var doesn't have any size field, but it doesn't matter |
| * since only functions dealing with the Heap care about it */ |
| #define MPFR_TMP_INIT1(xp, x, p) \ |
| ( MPFR_PREC(x) = (p), \ |
| MPFR_MANT(x) = (xp), \ |
| MPFR_SET_POS(x), \ |
| MPFR_SET_INVALID_EXP(x)) |
| |
| #define MPFR_TMP_INIT(xp, x, p, s) \ |
| (xp = (mp_ptr) MPFR_TMP_ALLOC(BYTES_PER_MP_LIMB * ((size_t) s)), \ |
| MPFR_TMP_INIT1(xp, x, p)) |
| |
| #define MPFR_TMP_INIT_ABS(d, s) \ |
| ( MPFR_PREC(d) = MPFR_PREC(s), \ |
| MPFR_MANT(d) = MPFR_MANT(s), \ |
| MPFR_SET_POS(d), \ |
| MPFR_EXP(d) = MPFR_EXP(s)) |
| |
| |
| |
| /****************************************************** |
| ***************** Cache macros ********************** |
| ******************************************************/ |
| |
| #define mpfr_const_pi(_d,_r) mpfr_cache(_d, __gmpfr_cache_const_pi,_r) |
| #define mpfr_const_log2(_d,_r) mpfr_cache(_d, __gmpfr_cache_const_log2, _r) |
| #define mpfr_const_euler(_d,_r) mpfr_cache(_d, __gmpfr_cache_const_euler, _r) |
| #define mpfr_const_catalan(_d,_r) mpfr_cache(_d,__gmpfr_cache_const_catalan,_r) |
| |
| #define MPFR_DECL_INIT_CACHE(_cache,_func) \ |
| mpfr_cache_t MPFR_THREAD_ATTR _cache = \ |
| {{{{0,MPFR_SIGN_POS,0,(mp_limb_t*)0}},0,_func}} |
| |
| |
| |
| /****************************************************** |
| ******************* Threshold *********************** |
| ******************************************************/ |
| |
| #include "mparam.h" |
| |
| /****************************************************** |
| ***************** Useful macros ********************* |
| ******************************************************/ |
| |
| /* Theses macros help the compiler to determine if a test is |
| * likely or unlikely. */ |
| #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)) |
| #else |
| # define MPFR_LIKELY(x) (x) |
| # define MPFR_UNLIKELY(x) (x) |
| #endif |
| |
| /* Declare that some variable is initialized before being used (without a |
| dummy initialization) in order to avoid some compiler warnings. Use the |
| VAR = VAR trick (see http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36296) |
| only with gcc as this is undefined behavior, and we don't know what |
| other compilers do (they may also be smarter). This trick could be |
| disabled with future gcc versions. */ |
| #if defined(__GNUC__) |
| # define INITIALIZED(VAR) VAR = VAR |
| #else |
| # define INITIALIZED(VAR) VAR |
| #endif |
| |
| /* Ceil log 2: If GCC, uses a GCC extension, otherwise calls a function */ |
| /* Warning: |
| * Needs to define MPFR_NEED_LONGLONG. |
| * Computes ceil(log2(x)) only for x integer (unsigned long) |
| * Undefined if x is 0 */ |
| #if __MPFR_GNUC(2,95) || __MPFR_ICC(8,1,0) |
| # define MPFR_INT_CEIL_LOG2(x) \ |
| (MPFR_UNLIKELY ((x) == 1) ? 0 : \ |
| __extension__ ({ int _b; mp_limb_t _limb; \ |
| MPFR_ASSERTN ((x) > 1); \ |
| _limb = (x) - 1; \ |
| MPFR_ASSERTN (_limb == (x) - 1); \ |
| count_leading_zeros (_b, _limb); \ |
| (BITS_PER_MP_LIMB - _b); })) |
| #else |
| # define MPFR_INT_CEIL_LOG2(x) (__gmpfr_int_ceil_log2(x)) |
| #endif |
| |
| /* Add two integers with overflow handling */ |
| /* Example: MPFR_SADD_OVERFLOW (c, a, b, long, unsigned long, |
| * LONG_MIN, LONG_MAX, |
| * goto overflow, goto underflow); */ |
| #define MPFR_UADD_OVERFLOW(c,a,b,ACTION_IF_OVERFLOW) \ |
| do { \ |
| (c) = (a) + (b); \ |
| if ((c) < (a)) ACTION_IF_OVERFLOW; \ |
| } while (0) |
| |
| #define MPFR_SADD_OVERFLOW(c,a,b,STYPE,UTYPE,MIN,MAX,ACTION_IF_POS_OVERFLOW,ACTION_IF_NEG_OVERFLOW) \ |
| do { \ |
| if ((a) >= 0 && (b) >= 0) { \ |
| UTYPE uc,ua,ub; \ |
| ua = (UTYPE) a; ub = (UTYPE) b; \ |
| MPFR_UADD_OVERFLOW (uc, ua, ub, ACTION_IF_POS_OVERFLOW); \ |
| if (uc > (UTYPE)(MAX)) ACTION_IF_POS_OVERFLOW; \ |
| else (c) = (STYPE) uc; \ |
| } else if ((a) < 0 && (b) < 0) { \ |
| UTYPE uc,ua,ub; \ |
| ua = -(UTYPE) a; ub = -(UTYPE) b; \ |
| MPFR_UADD_OVERFLOW (uc, ua, ub, ACTION_IF_NEG_OVERFLOW); \ |
| if (uc >= -(UTYPE)(MIN) || uc > (UTYPE)(MAX)) { \ |
| if (uc == -(UTYPE)(MIN)) (c) = (MIN); \ |
| else ACTION_IF_NEG_OVERFLOW; } \ |
| else (c) = -(STYPE) uc; \ |
| } else (c) = (a) + (b); \ |
| } while (0) |
| |
| |
| /* 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; \ |
| MPFR_SET_POS(x); \ |
| MPFR_EXP(x) = 1; \ |
| MPN_ZERO ( MPFR_MANT(x), _size); \ |
| MPFR_MANT(x)[_size] = MPFR_LIMB_HIGHBIT; \ |
| } while (0) |
| |
| /* Compute s = (-a) % BITS_PER_MP_LIMB |
| * a is unsigned! Check if it works, |
| * otherwise tries another way to compute it */ |
| #define MPFR_UNSIGNED_MINUS_MODULO(s, a) \ |
| do \ |
| { \ |
| if (IS_POW2 (BITS_PER_MP_LIMB)) \ |
| (s) = (-(a)) % BITS_PER_MP_LIMB; \ |
| else \ |
| { \ |
| (s) = (a) % BITS_PER_MP_LIMB; \ |
| if ((s) != 0) \ |
| (s) = BITS_PER_MP_LIMB - (s); \ |
| } \ |
| MPFR_ASSERTD ((s) >= 0 && (s) < BITS_PER_MP_LIMB); \ |
| } \ |
| while (0) |
| |
| /* Use it only for debug reasons */ |
| /* MPFR_TRACE (operation) : execute operation iff DEBUG flag is set */ |
| /* MPFR_DUMP (x) : print x (a mpfr_t) on stdout */ |
| #ifdef DEBUG |
| # define MPFR_TRACE(x) x |
| #else |
| # define MPFR_TRACE(x) (void) 0 |
| #endif |
| #define MPFR_DUMP(x) ( printf(#x"="), mpfr_dump(x) ) |
| |
| /* Test if X (positive) is a power of 2 */ |
| #define IS_POW2(X) (((X) & ((X) - 1)) == 0) |
| #define NOT_POW2(X) (((X) & ((X) - 1)) != 0) |
| |
| /* Safe absolute value (to avoid possible integer overflow) */ |
| /* type is the target (unsigned) type */ |
| #define SAFE_ABS(type,x) ((x) >= 0 ? (type)(x) : -(type)(x)) |
| |
| #define mpfr_get_d1(x) mpfr_get_d(x,__gmpfr_default_rounding_mode) |
| |
| /* Store in r the size in bits of the mpz_t z */ |
| #define MPFR_MPZ_SIZEINBASE2(r, z) \ |
| do { \ |
| int _cnt; \ |
| mp_size_t _size; \ |
| MPFR_ASSERTD (mpz_sgn (z) != 0); \ |
| _size = ABSIZ(z); \ |
| count_leading_zeros (_cnt, PTR(z)[_size-1]); \ |
| (r) = _size * BITS_PER_MP_LIMB - _cnt; \ |
| } while (0) |
| |
| /* Needs <locale.h> */ |
| #ifdef HAVE_LOCALE_H |
| #include <locale.h> |
| /* Warning! In case of signed char, the value of MPFR_DECIMAL_POINT may |
| be negative (the ISO C99 does not seem to forbid negative values). */ |
| #define MPFR_DECIMAL_POINT (localeconv()->decimal_point[0]) |
| #define MPFR_THOUSANDS_SEPARATOR (localeconv()->thousands_sep[0]) |
| #else |
| #define MPFR_DECIMAL_POINT ((char) '.') |
| #define MPFR_THOUSANDS_SEPARATOR ('\0') |
| #endif |
| |
| |
| /* Set y to s*significand(x)*2^e, for example MPFR_ALIAS(y,x,1,MPFR_EXP(x)) |
| sets y to |x|, and MPFR_ALIAS(y,x,MPFR_SIGN(x),0) sets y to x*2^f such |
| that 1/2 <= |y| < 1. Does not check y is in the valid exponent range. |
| WARNING! x and y share the same mantissa. So, some operations are |
| not valid if x has been provided via an argument, e.g., trying to |
| modify the mantissa of y, even temporarily, or calling mpfr_clear on y. |
| */ |
| #define MPFR_ALIAS(y,x,s,e) \ |
| do \ |
| { \ |
| MPFR_PREC(y) = MPFR_PREC(x); \ |
| MPFR_SIGN(y) = (s); \ |
| MPFR_EXP(y) = (e); \ |
| MPFR_MANT(y) = MPFR_MANT(x); \ |
| } while (0) |
| |
| |
| /****************************************************** |
| ************** Save exponent macros **************** |
| ******************************************************/ |
| |
| /* See README.dev for details on how to use the macros. |
| They are used to set the exponent range to the maximum |
| temporarily */ |
| |
| typedef struct { |
| unsigned int saved_flags; |
| mp_exp_t saved_emin; |
| mp_exp_t saved_emax; |
| } mpfr_save_expo_t; |
| |
| #define MPFR_SAVE_EXPO_DECL(x) mpfr_save_expo_t x |
| #define MPFR_SAVE_EXPO_MARK(x) \ |
| ((x).saved_flags = __gmpfr_flags, \ |
| (x).saved_emin = __gmpfr_emin, \ |
| (x).saved_emax = __gmpfr_emax, \ |
| __gmpfr_emin = MPFR_EMIN_MIN, \ |
| __gmpfr_emax = MPFR_EMAX_MAX) |
| #define MPFR_SAVE_EXPO_FREE(x) \ |
| (__gmpfr_flags = (x).saved_flags, \ |
| __gmpfr_emin = (x).saved_emin, \ |
| __gmpfr_emax = (x).saved_emax) |
| #define MPFR_SAVE_EXPO_UPDATE_FLAGS(x, flags) \ |
| (x).saved_flags |= (flags) |
| |
| /* Speed up final checking */ |
| #define mpfr_check_range(x,t,r) \ |
| (MPFR_LIKELY (MPFR_EXP (x) >= __gmpfr_emin && MPFR_EXP (x) <= __gmpfr_emax) \ |
| ? ((t) ? (__gmpfr_flags |= MPFR_FLAGS_INEXACT, (t)) : 0) \ |
| : mpfr_check_range(x,t,r)) |
| |
| |
| /****************************************************** |
| ***************** Inline Rounding ******************* |
| ******************************************************/ |
| |
| /* |
| * Note: due to the labels, one cannot use a macro MPFR_RNDRAW* more than |
| * once in a function (otherwise these labels would not be unique). |
| */ |
| |
| /* |
| * Round mantissa (srcp, sprec) to mpfr_t dest using rounding mode rnd |
| * assuming dest's sign is sign. |
| * In rounding to nearest mode, execute MIDDLE_HANDLER when the value |
| * is the middle of two consecutive numbers in dest precision. |
| * Execute OVERFLOW_HANDLER in case of overflow when rounding. |
| */ |
| #define MPFR_RNDRAW_GEN(inexact, dest, srcp, sprec, rnd, sign, \ |
| MIDDLE_HANDLER, OVERFLOW_HANDLER) \ |
| do { \ |
| mp_size_t _dests, _srcs; \ |
| mp_limb_t *_destp; \ |
| mp_prec_t _destprec, _srcprec; \ |
| \ |
| /* Check Trivial Case when Dest Mantissa has more bits than source */ \ |
| _srcprec = sprec; \ |
| _destprec = MPFR_PREC (dest); \ |
| _destp = MPFR_MANT (dest); \ |
| if (MPFR_UNLIKELY (_destprec >= _srcprec)) \ |
| { \ |
| _srcs = (_srcprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB; \ |
| _dests = (_destprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB - _srcs; \ |
| MPN_COPY (_destp + _dests, srcp, _srcs); \ |
| MPN_ZERO (_destp, _dests); \ |
| inexact = 0; \ |
| } \ |
| else \ |
| { \ |
| /* Non trivial case: rounding needed */ \ |
| mp_prec_t _sh; \ |
| mp_limb_t *_sp; \ |
| mp_limb_t _rb, _sb, _ulp; \ |
| \ |
| /* Compute Position and shift */ \ |
| _srcs = (_srcprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB; \ |
| _dests = (_destprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB; \ |
| MPFR_UNSIGNED_MINUS_MODULO (_sh, _destprec); \ |
| _sp = srcp + _srcs - _dests; \ |
| \ |
| /* General case when prec % BITS_PER_MP_LIMB != 0 */ \ |
| if (MPFR_LIKELY (_sh != 0)) \ |
| { \ |
| mp_limb_t _mask; \ |
| /* Compute Rounding Bit and Sticky Bit */ \ |
| _mask = MPFR_LIMB_ONE << (_sh - 1); \ |
| _rb = _sp[0] & _mask; \ |
| _sb = _sp[0] & (_mask - 1); \ |
| if (MPFR_UNLIKELY (_sb == 0)) \ |
| { /* TODO: Improve it */ \ |
| mp_limb_t *_tmp; \ |
| mp_size_t _n; \ |
| for (_tmp = _sp, _n = _srcs - _dests ; \ |
| _n != 0 && _sb == 0 ; _n--) \ |
| _sb = *--_tmp; \ |
| } \ |
| _ulp = 2 * _mask; \ |
| } \ |
| else /* _sh == 0 */ \ |
| { \ |
| MPFR_ASSERTD (_dests < _srcs); \ |
| /* Compute Rounding Bit and Sticky Bit */ \ |
| _rb = _sp[-1] & MPFR_LIMB_HIGHBIT; \ |
| _sb = _sp[-1] & (MPFR_LIMB_HIGHBIT-1); \ |
| if (MPFR_UNLIKELY (_sb == 0)) \ |
| { \ |
| mp_limb_t *_tmp; \ |
| mp_size_t _n; \ |
| for (_tmp = _sp - 1, _n = _srcs - _dests - 1 ; \ |
| _n != 0 && _sb == 0 ; _n--) \ |
| _sb = *--_tmp; \ |
| } \ |
| _ulp = MPFR_LIMB_ONE; \ |
| } \ |
| /* Rounding */ \ |
| if (MPFR_LIKELY (rnd == GMP_RNDN)) \ |
| { \ |
| if (_rb == 0) \ |
| { \ |
| trunc: \ |
| inexact = MPFR_LIKELY ((_sb | _rb) != 0) ? -sign : 0; \ |
| trunc_doit: \ |
| MPN_COPY (_destp, _sp, _dests); \ |
| _destp[0] &= ~(_ulp - 1); \ |
| } \ |
| else if (MPFR_UNLIKELY (_sb == 0)) \ |
| { /* Middle of two consecutive representable numbers */ \ |
| MIDDLE_HANDLER; \ |
| } \ |
| else \ |
| { \ |
| if (0) \ |
| goto addoneulp_doit; /* dummy code to avoid warning */ \ |
| addoneulp: \ |
| inexact = sign; \ |
| addoneulp_doit: \ |
| if (MPFR_UNLIKELY (mpn_add_1 (_destp, _sp, _dests, _ulp))) \ |
| { \ |
| _destp[_dests - 1] = MPFR_LIMB_HIGHBIT; \ |
| OVERFLOW_HANDLER; \ |
| } \ |
| _destp[0] &= ~(_ulp - 1); \ |
| } \ |
| } \ |
| else \ |
| { /* Directed rounding mode */ \ |
| if (MPFR_LIKELY (MPFR_IS_LIKE_RNDZ (rnd, \ |
| MPFR_IS_NEG_SIGN (sign)))) \ |
| goto trunc; \ |
| else if (MPFR_UNLIKELY ((_sb | _rb) == 0)) \ |
| { \ |
| inexact = 0; \ |
| goto trunc_doit; \ |
| } \ |
| else \ |
| goto addoneulp; \ |
| } \ |
| } \ |
| } while (0) |
| |
| /* |
| * Round mantissa (srcp, sprec) to mpfr_t dest using rounding mode rnd |
| * assuming dest's sign is sign. |
| * Execute OVERFLOW_HANDLER in case of overflow when rounding. |
| */ |
| #define MPFR_RNDRAW(inexact, dest, srcp, sprec, rnd, sign, OVERFLOW_HANDLER) \ |
| MPFR_RNDRAW_GEN (inexact, dest, srcp, sprec, rnd, sign, \ |
| if ((_sp[0] & _ulp) == 0) \ |
| { \ |
| inexact = -sign; \ |
| goto trunc_doit; \ |
| } \ |
| else \ |
| goto addoneulp; \ |
| , OVERFLOW_HANDLER) |
| |
| /* |
| * Round mantissa (srcp, sprec) to mpfr_t dest using rounding mode rnd |
| * assuming dest's sign is sign. |
| * Execute OVERFLOW_HANDLER in case of overflow when rounding. |
| * Set inexact to +/- MPFR_EVEN_INEX in case of even rounding. |
| */ |
| #define MPFR_RNDRAW_EVEN(inexact, dest, srcp, sprec, rnd, sign, \ |
| OVERFLOW_HANDLER) \ |
| MPFR_RNDRAW_GEN (inexact, dest, srcp, sprec, rnd, sign, \ |
| if ((_sp[0] & _ulp) == 0) \ |
| { \ |
| inexact = -MPFR_EVEN_INEX * sign; \ |
| goto trunc_doit; \ |
| } \ |
| else \ |
| { \ |
| inexact = MPFR_EVEN_INEX * sign; \ |
| goto addoneulp_doit; \ |
| } \ |
| , OVERFLOW_HANDLER) |
| |
| /* Return TRUE if b is non singular and we can round it to precision 'prec' |
| and determine the ternary value, with rounding mode 'rnd', and with |
| error at most 'error' */ |
| #define MPFR_CAN_ROUND(b,err,prec,rnd) \ |
| (!MPFR_IS_SINGULAR (b) && mpfr_round_p (MPFR_MANT (b), MPFR_LIMB_SIZE (b), \ |
| (err), (prec) + ((rnd)==GMP_RNDN))) |
| |
| /* TODO: fix this description (see round_near_x.c). */ |
| /* Assuming that the function has a Taylor expansion which looks like: |
| y=o(f(x)) = o(v + g(x)) with |g(x)| <= 2^(EXP(v)-err) |
| we can quickly set y to v if x is small (ie err > prec(y)+1) in most |
| cases. It assumes that f(x) is not representable exactly as a FP number. |
| v must not be a singular value (NAN, INF or ZERO); usual values are |
| v=1 or v=x. |
| |
| y is the destination (a mpfr_t), v the value to set (a mpfr_t), |
| err1+err2 with err2 <= 3 the error term (mp_exp_t's), dir (an int) is |
| the direction of the committed error (if dir = 0, it rounds towards 0, |
| if dir=1, it rounds away from 0), rnd the rounding mode. |
| |
| It returns from the function a ternary value in case of success. |
| If you want to free something, you must fill the "extra" field |
| in consequences, otherwise put nothing in it. |
| |
| The test is less restrictive than necessary, but the function |
| will finish the check itself. |
| |
| Note: err1 + err2 is allowed to overflow as mp_exp_t, but it must give |
| its real value as mpfr_uexp_t. |
| */ |
| #define MPFR_FAST_COMPUTE_IF_SMALL_INPUT(y,v,err1,err2,dir,rnd,extra) \ |
| do { \ |
| mpfr_ptr _y = (y); \ |
| mp_exp_t _err1 = (err1); \ |
| mp_exp_t _err2 = (err2); \ |
| if (_err1 > 0) \ |
| { \ |
| mpfr_uexp_t _err = (mpfr_uexp_t) _err1 + _err2; \ |
| if (MPFR_UNLIKELY (_err > MPFR_PREC (_y) + 1)) \ |
| { \ |
| int _inexact = mpfr_round_near_x (_y,(v),_err,(dir),(rnd)); \ |
| if (_inexact != 0) \ |
| { \ |
| extra; \ |
| return _inexact; \ |
| } \ |
| } \ |
| } \ |
| } while (0) |
| |
| /* Variant, to be called somewhere after MPFR_SAVE_EXPO_MARK. This variant |
| is needed when there are some computations before or when some non-zero |
| real constant is used, such as __gmpfr_one for mpfr_cos. */ |
| #define MPFR_SMALL_INPUT_AFTER_SAVE_EXPO(y,v,err1,err2,dir,rnd,expo,extra) \ |
| do { \ |
| mpfr_ptr _y = (y); \ |
| mp_exp_t _err1 = (err1); \ |
| mp_exp_t _err2 = (err2); \ |
| if (_err1 > 0) \ |
| { \ |
| mpfr_uexp_t _err = (mpfr_uexp_t) _err1 + _err2; \ |
| if (MPFR_UNLIKELY (_err > MPFR_PREC (_y) + 1)) \ |
| { \ |
| int _inexact; \ |
| mpfr_clear_flags (); \ |
| _inexact = mpfr_round_near_x (_y,(v),_err,(dir),(rnd)); \ |
| if (_inexact != 0) \ |
| { \ |
| extra; \ |
| MPFR_SAVE_EXPO_UPDATE_FLAGS (expo, __gmpfr_flags); \ |
| MPFR_SAVE_EXPO_FREE (expo); \ |
| return mpfr_check_range (_y, _inexact, (rnd)); \ |
| } \ |
| } \ |
| } \ |
| } while (0) |
| |
| /****************************************************** |
| *************** Ziv Loop Macro ********************* |
| ******************************************************/ |
| |
| #ifndef MPFR_USE_LOGGING |
| |
| #define MPFR_ZIV_DECL(_x) mp_prec_t _x |
| #define MPFR_ZIV_INIT(_x, _p) (_x) = BITS_PER_MP_LIMB |
| #define MPFR_ZIV_NEXT(_x, _p) ((_p) += (_x), (_x) = (_p)/2) |
| #define MPFR_ZIV_FREE(x) |
| |
| #else |
| |
| /* The following test on glibc is there mainly for Darwin (Mac OS X), to |
| obtain a better error message. The real test should have been a test |
| concerning nested functions in gcc, which are disabled by default on |
| Darwin; but it is not possible to do that without a configure test. */ |
| # if defined (__cplusplus) || !(__MPFR_GNUC(3,0) && __MPFR_GLIBC(2,0)) |
| # error "Logging not supported (needs gcc >= 3.0 and GNU C Library >= 2.0)." |
| # endif |
| |
| /* Use LOGGING */ |
| #define MPFR_ZIV_DECL(_x) \ |
| mp_prec_t _x; \ |
| int _x ## _cpt = 1; \ |
| static unsigned long _x ## _loop = 0, _x ## _bad = 0; \ |
| static const char *_x ## _fname = __func__; \ |
| auto void __attribute__ ((destructor)) x ## _f (void); \ |
| void __attribute__ ((destructor)) x ## _f (void) { \ |
| if (_x ## _loop != 0 && MPFR_LOG_STAT_F&mpfr_log_type) \ |
| fprintf (mpfr_log_file, \ |
| "%s: Ziv failed %2.2f%% (%lu bad cases / %lu calls)\n", _x ## _fname, \ |
| (double) 100.0 * _x ## _bad / _x ## _loop, _x ## _bad, _x ## _loop ); } |
| |
| #define MPFR_ZIV_INIT(_x, _p) ((_x) = BITS_PER_MP_LIMB, _x ## _loop ++); \ |
| if (MPFR_LOG_BADCASE_F&mpfr_log_type && mpfr_log_current<=mpfr_log_level) \ |
| fprintf (mpfr_log_file, "%s:ZIV 1st prec=%lu\n", __func__, \ |
| (unsigned long) (_p)) |
| |
| #define MPFR_ZIV_NEXT(_x, _p) \ |
| ((_p)+=(_x),(_x)=(_p)/2, _x ## _bad += (_x ## _cpt == 1), _x ## _cpt ++); \ |
| if (MPFR_LOG_BADCASE_F&mpfr_log_type && mpfr_log_current<=mpfr_log_level) \ |
| fprintf (mpfr_log_file, "%s:ZIV new prec=%lu\n", __func__, \ |
| (unsigned long) (_p)) |
| |
| #define MPFR_ZIV_FREE(_x) \ |
| if (MPFR_LOG_BADCASE_F&mpfr_log_type && _x##_cpt>1 \ |
| && mpfr_log_current<=mpfr_log_level) \ |
| fprintf (mpfr_log_file, "%s:ZIV %d loops\n", __func__, _x ## _cpt) |
| |
| #endif |
| |
| |
| /****************************************************** |
| *************** Logging Macros ********************* |
| ******************************************************/ |
| |
| /* The different kind of LOG */ |
| #define MPFR_LOG_INPUT_F 1 |
| #define MPFR_LOG_OUTPUT_F 2 |
| #define MPFR_LOG_INTERNAL_F 4 |
| #define MPFR_LOG_TIME_F 8 |
| #define MPFR_LOG_BADCASE_F 16 |
| #define MPFR_LOG_MSG_F 32 |
| #define MPFR_LOG_STAT_F 64 |
| |
| #ifdef MPFR_USE_LOGGING |
| |
| /* Check if we can support this feature */ |
| # ifdef MPFR_USE_THREAD_SAFE |
| # error "Enable either `Logging' or `thread-safe', not both" |
| # endif |
| # if !__MPFR_GNUC(3,0) |
| # error "Logging not supported (GCC >= 3.0)" |
| # endif |
| |
| #if defined (__cplusplus) |
| extern "C" { |
| #endif |
| |
| __MPFR_DECLSPEC extern FILE *mpfr_log_file; |
| __MPFR_DECLSPEC extern int mpfr_log_type; |
| __MPFR_DECLSPEC extern int mpfr_log_level; |
| __MPFR_DECLSPEC extern int mpfr_log_current; |
| __MPFR_DECLSPEC extern int mpfr_log_base; |
| __MPFR_DECLSPEC extern mp_prec_t mpfr_log_prec; |
| |
| #if defined (__cplusplus) |
| } |
| #endif |
| |
| #define MPFR_LOG_VAR(x) \ |
| if((MPFR_LOG_INTERNAL_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level))\ |
| fprintf (mpfr_log_file, "%s.%d:%s[%#R]=%R\n", __func__,__LINE__, #x, x, x); |
| |
| #define MPFR_LOG_MSG2(format, ...) \ |
| if ((MPFR_LOG_MSG_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \ |
| fprintf (mpfr_log_file, "%s.%d: "format, __func__, __LINE__, __VA_ARGS__); |
| #define MPFR_LOG_MSG(x) MPFR_LOG_MSG2 x |
| |
| #define MPFR_LOG_BEGIN2(format, ...) \ |
| mpfr_log_current ++; \ |
| if ((MPFR_LOG_INPUT_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \ |
| fprintf (mpfr_log_file, "%s:IN "format"\n",__func__,__VA_ARGS__); \ |
| if ((MPFR_LOG_TIME_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \ |
| __gmpfr_log_time = mpfr_get_cputime (); |
| #define MPFR_LOG_BEGIN(x) \ |
| int __gmpfr_log_time = 0; \ |
| MPFR_LOG_BEGIN2 x |
| |
| #define MPFR_LOG_END2(format, ...) \ |
| if ((MPFR_LOG_TIME_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \ |
| fprintf (mpfr_log_file, "%s:TIM %dms\n", __mpfr_log_fname, \ |
| mpfr_get_cputime () - __gmpfr_log_time); \ |
| if ((MPFR_LOG_OUTPUT_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \ |
| fprintf (mpfr_log_file, "%s:OUT "format"\n",__mpfr_log_fname,__VA_ARGS__);\ |
| mpfr_log_current --; |
| #define MPFR_LOG_END(x) \ |
| static const char *__mpfr_log_fname = __func__; \ |
| MPFR_LOG_END2 x |
| |
| #define MPFR_LOG_FUNC(begin,end) \ |
| static const char *__mpfr_log_fname = __func__; \ |
| auto void __mpfr_log_cleanup (int *time); \ |
| void __mpfr_log_cleanup (int *time) { \ |
| int __gmpfr_log_time = *time; \ |
| MPFR_LOG_END2 end; } \ |
| int __gmpfr_log_time __attribute__ ((cleanup (__mpfr_log_cleanup))); \ |
| __gmpfr_log_time = 0; \ |
| MPFR_LOG_BEGIN2 begin |
| |
| #else /* MPFR_USE_LOGGING */ |
| |
| /* Define void macro for logging */ |
| |
| #define MPFR_LOG_VAR(x) |
| #define MPFR_LOG_BEGIN(x) |
| #define MPFR_LOG_END(x) |
| #define MPFR_LOG_MSG(x) |
| #define MPFR_LOG_FUNC(x,y) |
| |
| #endif /* MPFR_USE_LOGGING */ |
| |
| |
| /************************************************************** |
| ************ Group Initialize Functions Macros ************* |
| **************************************************************/ |
| |
| #ifndef MPFR_GROUP_STATIC_SIZE |
| # define MPFR_GROUP_STATIC_SIZE 16 |
| #endif |
| |
| struct mpfr_group_t { |
| size_t alloc; |
| mp_limb_t *mant; |
| mp_limb_t tab[MPFR_GROUP_STATIC_SIZE]; |
| }; |
| |
| #define MPFR_GROUP_DECL(g) struct mpfr_group_t g |
| #define MPFR_GROUP_CLEAR(g) do { \ |
| if (MPFR_UNLIKELY ((g).alloc != 0)) { \ |
| MPFR_ASSERTD ((g).mant != (g).tab); \ |
| (*__gmp_free_func) ((g).mant, (g).alloc); \ |
| }} while (0) |
| |
| #define MPFR_GROUP_INIT_TEMPLATE(g, prec, num, handler) do { \ |
| mp_prec_t _prec = (prec); \ |
| mp_size_t _size; \ |
| MPFR_ASSERTD (_prec >= MPFR_PREC_MIN); \ |
| if (MPFR_UNLIKELY (_prec > MPFR_PREC_MAX)) \ |
| mpfr_abort_prec_max (); \ |
| _size = (mp_prec_t) (_prec + BITS_PER_MP_LIMB - 1) / BITS_PER_MP_LIMB; \ |
| if (MPFR_UNLIKELY (_size * (num) > MPFR_GROUP_STATIC_SIZE)) \ |
| { \ |
| (g).alloc = (num) * _size * sizeof (mp_limb_t); \ |
| (g).mant = (mp_limb_t *) (*__gmp_allocate_func) ((g).alloc); \ |
| } \ |
| else \ |
| { \ |
| (g).alloc = 0; \ |
| (g).mant = (g).tab; \ |
| } \ |
| handler; \ |
| } while (0) |
| #define MPFR_GROUP_TINIT(g, n, x) \ |
| MPFR_TMP_INIT1 ((g).mant + _size * (n), x, _prec) |
| |
| #define MPFR_GROUP_INIT_1(g, prec, x) \ |
| MPFR_GROUP_INIT_TEMPLATE(g, prec, 1, MPFR_GROUP_TINIT(g, 0, x)) |
| #define MPFR_GROUP_INIT_2(g, prec, x, y) \ |
| MPFR_GROUP_INIT_TEMPLATE(g, prec, 2, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y)) |
| #define MPFR_GROUP_INIT_3(g, prec, x, y, z) \ |
| MPFR_GROUP_INIT_TEMPLATE(g, prec, 3, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \ |
| MPFR_GROUP_TINIT(g, 2, z)) |
| #define MPFR_GROUP_INIT_4(g, prec, x, y, z, t) \ |
| MPFR_GROUP_INIT_TEMPLATE(g, prec, 4, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \ |
| MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t)) |
| #define MPFR_GROUP_INIT_5(g, prec, x, y, z, t, a) \ |
| MPFR_GROUP_INIT_TEMPLATE(g, prec, 5, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \ |
| MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \ |
| MPFR_GROUP_TINIT(g, 4, a)) |
| #define MPFR_GROUP_INIT_6(g, prec, x, y, z, t, a, b) \ |
| MPFR_GROUP_INIT_TEMPLATE(g, prec, 6, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \ |
| MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \ |
| MPFR_GROUP_TINIT(g, 4, a);MPFR_GROUP_TINIT(g, 5, b)) |
| |
| #define MPFR_GROUP_REPREC_TEMPLATE(g, prec, num, handler) do { \ |
| mp_prec_t _prec = (prec); \ |
| size_t _oalloc = (g).alloc; \ |
| mp_size_t _size; \ |
| MPFR_ASSERTD (_prec >= MPFR_PREC_MIN); \ |
| if (MPFR_UNLIKELY (_prec > MPFR_PREC_MAX)) \ |
| mpfr_abort_prec_max (); \ |
| _size = (mp_prec_t) (_prec + BITS_PER_MP_LIMB - 1) / BITS_PER_MP_LIMB; \ |
| (g).alloc = (num) * _size * sizeof (mp_limb_t); \ |
| if (MPFR_LIKELY (_oalloc == 0)) \ |
| (g).mant = (mp_limb_t *) (*__gmp_allocate_func) ((g).alloc); \ |
| else \ |
| (g).mant = (mp_limb_t *) \ |
| (*__gmp_reallocate_func) ((g).mant, _oalloc, (g).alloc); \ |
| handler; \ |
| } while (0) |
| |
| #define MPFR_GROUP_REPREC_1(g, prec, x) \ |
| MPFR_GROUP_REPREC_TEMPLATE(g, prec, 1, MPFR_GROUP_TINIT(g, 0, x)) |
| #define MPFR_GROUP_REPREC_2(g, prec, x, y) \ |
| MPFR_GROUP_REPREC_TEMPLATE(g, prec, 2, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y)) |
| #define MPFR_GROUP_REPREC_3(g, prec, x, y, z) \ |
| MPFR_GROUP_REPREC_TEMPLATE(g, prec, 3, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \ |
| MPFR_GROUP_TINIT(g, 2, z)) |
| #define MPFR_GROUP_REPREC_4(g, prec, x, y, z, t) \ |
| MPFR_GROUP_REPREC_TEMPLATE(g, prec, 4, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \ |
| MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t)) |
| #define MPFR_GROUP_REPREC_5(g, prec, x, y, z, t, a) \ |
| MPFR_GROUP_REPREC_TEMPLATE(g, prec, 5, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \ |
| MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \ |
| MPFR_GROUP_TINIT(g, 4, a)) |
| #define MPFR_GROUP_REPREC_6(g, prec, x, y, z, t, a, b) \ |
| MPFR_GROUP_REPREC_TEMPLATE(g, prec, 6, \ |
| MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \ |
| MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \ |
| MPFR_GROUP_TINIT(g, 4, a);MPFR_GROUP_TINIT(g, 5, b)) |
| |
| |
| /****************************************************** |
| *************** Internal Functions ***************** |
| ******************************************************/ |
| |
| #if defined (__cplusplus) |
| extern "C" { |
| #endif |
| |
| __MPFR_DECLSPEC int mpfr_underflow _MPFR_PROTO ((mpfr_ptr, mp_rnd_t, int)); |
| __MPFR_DECLSPEC int mpfr_overflow _MPFR_PROTO ((mpfr_ptr, mp_rnd_t, int)); |
| |
| __MPFR_DECLSPEC int mpfr_add1 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr, |
| mpfr_srcptr, mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_sub1 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr, |
| mpfr_srcptr, mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_add1sp _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr, |
| mpfr_srcptr, mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_sub1sp _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr, |
| mpfr_srcptr, mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_can_round_raw _MPFR_PROTO ((const mp_limb_t *, |
| mp_size_t, int, mp_exp_t, mp_rnd_t, mp_rnd_t, mp_prec_t)); |
| |
| __MPFR_DECLSPEC int mpfr_cmp2 _MPFR_PROTO ((mpfr_srcptr, mpfr_srcptr, |
| mp_prec_t *)); |
| |
| __MPFR_DECLSPEC long __gmpfr_ceil_log2 _MPFR_PROTO ((double)); |
| __MPFR_DECLSPEC long __gmpfr_floor_log2 _MPFR_PROTO ((double)); |
| __MPFR_DECLSPEC double __gmpfr_ceil_exp2 _MPFR_PROTO ((double)); |
| __MPFR_DECLSPEC unsigned long __gmpfr_isqrt _MPFR_PROTO ((unsigned long)); |
| __MPFR_DECLSPEC unsigned long __gmpfr_cuberoot _MPFR_PROTO ((unsigned long)); |
| __MPFR_DECLSPEC int __gmpfr_int_ceil_log2 _MPFR_PROTO ((unsigned long)); |
| |
| __MPFR_DECLSPEC int mpfr_exp_2 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_exp_3 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_powerof2_raw _MPFR_PROTO ((mpfr_srcptr)); |
| |
| __MPFR_DECLSPEC int mpfr_pow_general _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr, |
| mpfr_srcptr, mp_rnd_t, int, mpfr_save_expo_t *)); |
| |
| __MPFR_DECLSPEC void mpfr_setmax _MPFR_PROTO ((mpfr_ptr, mp_exp_t)); |
| __MPFR_DECLSPEC void mpfr_setmin _MPFR_PROTO ((mpfr_ptr, mp_exp_t)); |
| |
| __MPFR_DECLSPEC long mpfr_mpn_exp _MPFR_PROTO ((mp_limb_t *, mp_exp_t *, int, |
| mp_exp_t, size_t)); |
| |
| #ifdef _MPFR_H_HAVE_FILE |
| __MPFR_DECLSPEC void mpfr_fprint_binary _MPFR_PROTO ((FILE *, mpfr_srcptr)); |
| #endif |
| __MPFR_DECLSPEC void mpfr_print_binary _MPFR_PROTO ((mpfr_srcptr)); |
| __MPFR_DECLSPEC void mpfr_print_mant_binary _MPFR_PROTO ((const char*, |
| const mp_limb_t*, mp_prec_t)); |
| __MPFR_DECLSPEC void mpfr_set_str_binary _MPFR_PROTO((mpfr_ptr, const char*)); |
| |
| __MPFR_DECLSPEC int mpfr_round_raw _MPFR_PROTO ((mp_limb_t *, |
| const mp_limb_t *, mp_prec_t, int, mp_prec_t, mp_rnd_t, int *)); |
| __MPFR_DECLSPEC int mpfr_round_raw_2 _MPFR_PROTO ((const mp_limb_t *, |
| mp_prec_t, int, mp_prec_t, mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_round_raw_3 _MPFR_PROTO ((const mp_limb_t *, |
| mp_prec_t, int, mp_prec_t, mp_rnd_t, int *)); |
| __MPFR_DECLSPEC int mpfr_round_raw_4 _MPFR_PROTO ((mp_limb_t *, |
| const mp_limb_t *, mp_prec_t, int, mp_prec_t, mp_rnd_t)); |
| |
| #define mpfr_round_raw2(xp, xn, neg, r, prec) \ |
| mpfr_round_raw_2((xp),(xn)*BITS_PER_MP_LIMB,(neg),(prec),(r)) |
| |
| __MPFR_DECLSPEC int mpfr_check _MPFR_PROTO ((mpfr_srcptr)); |
| |
| __MPFR_DECLSPEC int mpfr_sum_sort _MPFR_PROTO ((mpfr_srcptr *const, |
| unsigned long, mpfr_srcptr *)); |
| |
| __MPFR_DECLSPEC int mpfr_get_cputime _MPFR_PROTO ((void)); |
| |
| __MPFR_DECLSPEC void mpfr_nexttozero _MPFR_PROTO ((mpfr_ptr)); |
| __MPFR_DECLSPEC void mpfr_nexttoinf _MPFR_PROTO ((mpfr_ptr)); |
| |
| __MPFR_DECLSPEC int mpfr_const_pi_internal _MPFR_PROTO ((mpfr_ptr,mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_const_log2_internal _MPFR_PROTO((mpfr_ptr,mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_const_euler_internal _MPFR_PROTO((mpfr_ptr, mp_rnd_t)); |
| __MPFR_DECLSPEC int mpfr_const_catalan_internal _MPFR_PROTO((mpfr_ptr, mp_rnd_t)); |
| |
| __MPFR_DECLSPEC void mpfr_init_cache _MPFR_PROTO ((mpfr_cache_t, |
| int(*)(mpfr_ptr,mpfr_rnd_t))); |
| __MPFR_DECLSPEC void mpfr_clear_cache _MPFR_PROTO ((mpfr_cache_t)); |
| __MPFR_DECLSPEC int mpfr_cache _MPFR_PROTO ((mpfr_ptr, mpfr_cache_t, |
| mpfr_rnd_t)); |
| |
| __MPFR_DECLSPEC void mpfr_mulhigh_n _MPFR_PROTO ((mp_ptr, mp_srcptr, |
| mp_srcptr, mp_size_t)); |
| __MPFR_DECLSPEC void mpfr_sqrhigh_n _MPFR_PROTO ((mp_ptr, mp_srcptr, mp_size_t)); |
| |
| __MPFR_DECLSPEC int mpfr_round_p _MPFR_PROTO ((mp_limb_t *, mp_size_t, |
| mp_exp_t, mp_prec_t)); |
| |
| __MPFR_DECLSPEC void mpfr_dump_mant _MPFR_PROTO ((const mp_limb_t *, |
| mp_prec_t, mp_prec_t, |
| mp_prec_t)); |
| |
| __MPFR_DECLSPEC int mpfr_round_near_x _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr, |
| mpfr_uexp_t, int, |
| mp_rnd_t)); |
| __MPFR_DECLSPEC void mpfr_abort_prec_max _MPFR_PROTO ((void)) |
| MPFR_NORETURN_ATTR; |
| |
| #if defined (__cplusplus) |
| } |
| #endif |
| |
| #endif |