1
/* Utilities for MPFR developers, not exported.
3
Copyright 1999, 2000, 2001, 2002, 2003, 2004, 2005
4
Free Software Foundation, Inc.
6
This file is part of the MPFR Library.
8
The MPFR Library is free software; you can redistribute it and/or modify
9
it under the terms of the GNU Lesser General Public License as published by
10
the Free Software Foundation; either version 2.1 of the License, or (at your
11
option) any later version.
13
The MPFR Library is distributed in the hope that it will be useful, but
14
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15
or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16
License for more details.
18
You should have received a copy of the GNU Lesser General Public License
19
along with the MPFR Library; see the file COPYING.LIB. If not, write to
20
the Free Software Foundation, Inc., 51 Franklin Place, Fifth Floor, Boston,
21
MA 02110-1301, USA. */
23
#ifndef __MPFR_IMPL_H__
24
#define __MPFR_IMPL_H__
26
/* Include stdio.h iff we are debugging or we want to check */
27
#if defined(DEBUG) || defined(WANT_ASSERT)
31
/* Check if we are inside a build of MPFR or inside the test suite.
32
This is needed in mpfr.h to export or import the functions.
33
It matters only for Windows DLL */
34
#ifndef __MPFR_TEST_H__
35
# define __MPFR_WITHIN_MPFR 1
38
/******************************************************
39
****************** Include files *********************
40
******************************************************/
42
/* Include 'config.h' before using ANY configure macros if needed
43
NOTE: It isn't MPFR 'config.h', but GMP's one! */
44
#if defined(HAVE_CONFIG_H)
50
#ifdef MPFR_HAVE_GMP_IMPL /* Build with gmp internals*/
55
# ifndef __GMP_IMPL_H__
56
# include "gmp-impl.h"
58
# ifdef MPFR_NEED_LONGLONG_H
59
# include "longlong.h"
65
#else /* Build without gmp internals */
73
# ifndef __GMPFR_GMP_H__
74
# include "mpfr-gmp.h"
76
# ifdef MPFR_NEED_LONGLONG_H
77
# include "mpfr-longlong.h"
81
#undef MPFR_NEED_LONGLONG_H
84
/******************************************************
85
***************** Detection macros *******************
86
******************************************************/
88
/* Macros to detect STDC, GCC, GLIBC, GMP and ICC version */
89
#if defined(__STDC_VERSION__)
90
# define __MPFR_STDC(version) (__STDC_VERSION__>=(version))
91
#elif defined(__STDC__)
92
# define __MPFR_STDC(version) (0 == (version))
94
# define __MPFR_STDC(version) 0
97
#if defined(__GNUC__) && defined(__GNUC_MINOR__) && !defined(__ICC)
98
# define __MPFR_GNUC(a, i) \
99
(MPFR_VERSION_NUM(__GNUC__,__GNUC_MINOR__,0)>=MPFR_VERSION_NUM(a,i,0))
101
# define __MPFR_GNUC(a, i) 0
104
#if defined(__GLIBC__) && defined(__GLIBC_MINOR__)
105
# define __MPFR_GLIBC(a, i) \
106
(MPFR_VERSION_NUM(__GLIBC__,__GLIBC_MINOR__,0)>=MPFR_VERSION_NUM(a,i,0))
108
# define __MPFR_GLIBC(a, i) 0
111
#if defined(__GNU_MP_VERSION)&&defined(__GNU_MP_VERSION_MINOR)&&defined(__GNU_MP_VERSION_PATCHLEVEL)
112
# define __MPFR_GMP(a, b, c) \
113
(MPFR_VERSION_NUM(__GNU_MP_VERSION,__GNU_MP_VERSION_MINOR,__GNU_MP_VERSION_PATCHLEVEL) >= MPFR_VERSION_NUM(a,b,c))
115
# define __MPFR_GMP(a, b, c) 0
119
# define __MPFR_ICC(a,b,c) (__ICC >= (a)*100+(b)*10+c)
121
# define __MPFR_ICC(a,b,c) 0
126
/******************************************************
127
******************** Check GMP ***********************
128
******************************************************/
130
#if !__MPFR_GMP(4,1,0)
131
# error "GMP 4.1.0 or newer needed"
134
#if GMP_NAIL_BITS != 0
135
# error "MPFR doesn't support nonzero values of GMP_NAIL_BITS"
138
#if (BITS_PER_MP_LIMB<32) || (BITS_PER_MP_LIMB & (BITS_PER_MP_LIMB - 1))
139
# error "BITS_PER_MP_LIMB must be a power of 2, and >= 32"
142
#if BITS_PER_MP_LIMB == 16
143
# define MPFR_LOG2_BITS_PER_MP_LIMB 4
144
#elif BITS_PER_MP_LIMB == 32
145
# define MPFR_LOG2_BITS_PER_MP_LIMB 5
146
#elif BITS_PER_MP_LIMB == 64
147
# define MPFR_LOG2_BITS_PER_MP_LIMB 6
148
#elif BITS_PER_MP_LIMB == 128
149
# define MPFR_LOG2_BITS_PER_MP_LIMB 7
150
#elif BITS_PER_MP_LIMB == 256
151
# define MPFR_LOG2_BITS_PER_MP_LIMB 8
153
# error "Can't compute log2(BITS_PER_MP_LIMB)"
156
#if __MPFR_GNUC(3,0) || __MPFR_ICC(8,1,0)
157
# define MPFR_NORETURN_ATTR __attribute__ ((noreturn))
158
# define MPFR_CONST_ATTR __attribute__ ((const))
160
# define MPFR_NORETURN_ATTR
161
# define MPFR_CONST_ATTR
164
/******************************************************
165
************* Global Internal Variables **************
166
******************************************************/
168
#ifdef MPFR_USE_THREAD_SAFE
169
# if __MPFR_GNUC(3,3) || __MPFR_ICC(8,1,0)
170
# define MPFR_THREAD_ATTR __thread
172
# error "Can't build MPFR as thread safe"
175
# define MPFR_THREAD_ATTR
179
struct __gmpfr_cache_s {
182
int (*func)(mpfr_ptr, mpfr_rnd_t);
184
typedef struct __gmpfr_cache_s mpfr_cache_t[1];
186
#if defined (__cplusplus)
190
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR unsigned int __gmpfr_flags;
191
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR mp_exp_t __gmpfr_emin;
192
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR mp_exp_t __gmpfr_emax;
193
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR mp_prec_t __gmpfr_default_fp_bit_precision;
194
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_rnd_t __gmpfr_default_rounding_mode;
195
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_pi;
196
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_log2;
197
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_euler;
198
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_catalan;
200
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR const mpfr_t __gmpfr_one;
201
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR const mpfr_t __gmpfr_two;
202
__MPFR_DECLSPEC extern MPFR_THREAD_ATTR const mpfr_t __gmpfr_four;
205
#if defined (__cplusplus)
209
/* Flags of __gmpfr_flags */
210
#define MPFR_FLAGS_UNDERFLOW 1
211
#define MPFR_FLAGS_OVERFLOW 2
212
#define MPFR_FLAGS_NAN 4
213
#define MPFR_FLAGS_INEXACT 8
214
#define MPFR_FLAGS_ERANGE 16
215
#define MPFR_FLAGS_ALL 31
217
/* Replace some commun functions for direct access to the global vars */
218
#define mpfr_get_emin() (__gmpfr_emin + 0)
219
#define mpfr_get_emax() (__gmpfr_emax + 0)
220
#define mpfr_get_default_rounding_mode() (__gmpfr_default_rounding_mode + 0)
221
#define mpfr_get_default_prec() (__gmpfr_default_fp_bit_precision + 0)
223
#define mpfr_clear_flags() \
224
((void) (__gmpfr_flags = 0))
225
#define mpfr_clear_underflow() \
226
((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_UNDERFLOW))
227
#define mpfr_clear_overflow() \
228
((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_OVERFLOW))
229
#define mpfr_clear_nanflag() \
230
((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_NAN))
231
#define mpfr_clear_inexflag() \
232
((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_INEXACT))
233
#define mpfr_clear_erangeflag() \
234
((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_ERANGE))
235
#define mpfr_underflow_p() \
236
((int) (__gmpfr_flags & MPFR_FLAGS_UNDERFLOW))
237
#define mpfr_overflow_p() \
238
((int) (__gmpfr_flags & MPFR_FLAGS_OVERFLOW))
239
#define mpfr_nanflag_p() \
240
((int) (__gmpfr_flags & MPFR_FLAGS_NAN))
241
#define mpfr_inexflag_p() \
242
((int) (__gmpfr_flags & MPFR_FLAGS_INEXACT))
243
#define mpfr_erangeflag_p() \
244
((int) (__gmpfr_flags & MPFR_FLAGS_ERANGE))
247
/******************************************************
248
******************** Assertions **********************
249
******************************************************/
251
/* Compile with -DWANT_ASSERT to check all assert statements */
253
/* Note: do not use GMP macros ASSERT_ALWAYS and ASSERT as they are not
254
expressions, and as a consequence, they cannot be used in a for(),
255
with a comma operator and so on. */
257
/* MPFR_ASSERTN(expr): assertions that should always be checked */
258
#define MPFR_ASSERTN(expr) \
259
((void) ((MPFR_UNLIKELY(expr)) || MPFR_UNLIKELY( (ASSERT_FAIL(expr),0) )))
261
/* MPFR_ASSERTD(expr): assertions that should be checked when testing */
263
# define MPFR_EXP_CHECK 1
264
# define MPFR_ASSERTD(expr) MPFR_ASSERTN (expr)
266
# define MPFR_ASSERTD(expr) ((void) 0)
269
/* Check if the args are correct */
270
/* Can't be used since TMP variables are not correct */
271
#define MPFR_CHECK1(x,r) \
272
MPFR_ASSERTD(mpfr_check(x) && GMP_RNDN <= r && r <= GMP_RNDD)
273
#define MPFR_CHECK2(x,y,r) \
274
MPFR_ASSERTD(mpfr_check(x) && mpfr_check(y) && GMP_RNDN <= r && r <= GMP_RNDD)
275
#define MPFR_CHECK3(x,y,z,r) \
276
MPFR_ASSERTD(mpfr_check(x) && mpfr_check(y) && mpfr_check(z) && \
277
GMP_RNDN <= r && r <= GMP_RNDD)
279
/* Code to deal with impossible
280
WARNING: It doesn't use do { } while (0) for Insure++*/
281
#define MPFR_RET_NEVER_GO_HERE() {MPFR_ASSERTN(0); return 0;}
285
/******************************************************
286
****************** double macros *********************
287
******************************************************/
289
/* Definition of constants */
290
#define LOG2 0.69314718055994528622 /* log(2) rounded to zero on 53 bits */
291
#define ALPHA 4.3191365662914471407 /* a+2 = a*log(a), rounded to +infinity */
292
#define EXPM1 0.36787944117144227851 /* exp(-1), rounded to zero */
294
/* MPFR_DOUBLE_SPEC = 1 if the C type 'double' corresponds to IEEE-754
295
double precision, 0 if it doesn't, and undefined if one doesn't know.
296
On all the tested machines, MPFR_DOUBLE_SPEC = 1. To have this macro
297
defined here, #include <float.h> is needed. If need be, other values
298
could be defined for other specs (once they are known). */
299
#if !defined(MPFR_DOUBLE_SPEC) && defined(FLT_RADIX) && \
300
defined(DBL_MANT_DIG) && defined(DBL_MIN_EXP) && defined(DBL_MAX_EXP)
301
# if FLT_RADIX == 2 && DBL_MANT_DIG == 53 && \
302
DBL_MIN_EXP == -1021 && DBL_MAX_EXP == 1024
303
# define MPFR_DOUBLE_SPEC 1
305
# define MPFR_DOUBLE_SPEC 0
309
/* Debug non IEEE floats */
311
# undef _GMP_IEEE_FLOATS
313
#ifndef _GMP_IEEE_FLOATS
314
# define _GMP_IEEE_FLOATS 0
317
#ifndef IEEE_DBL_MANT_DIG
318
#define IEEE_DBL_MANT_DIG 53
320
#define MPFR_LIMBS_PER_DOUBLE ((IEEE_DBL_MANT_DIG-1)/BITS_PER_MP_LIMB+1)
322
/* Visual C++ doesn't support +1.0/.00, -1.0/0.0 and 0.0/0.0
324
#if defined(_MSC_VER) && defined(_WIN32) && (_MSC_VER >= 1200)
325
static double double_zero = 0.0;
326
# define DBL_NAN (double_zero/double_zero)
327
# define DBL_POS_INF ((double) 1.0/double_zero)
328
# define DBL_NEG_INF ((double)-1.0/double_zero)
330
# define DBL_POS_INF ((double) 1.0/0.0)
331
# define DBL_NEG_INF ((double)-1.0/0.0)
332
# define DBL_NAN ((double) 0.0/0.0)
335
/* for x of type ieee_double_extract */
337
typedef union ieee_double_extract Ieee_double_extract;
339
# define DOUBLE_ISNANorINF(x) (((Ieee_double_extract *)&(x))->s.exp == 0x7ff)
340
# define DOUBLE_ISINF(x) (DOUBLE_ISNANorINF(x) && \
341
(((Ieee_double_extract *)&(x))->s.manl == 0) && \
342
(((Ieee_double_extract *)&(x))->s.manh == 0))
343
# define DOUBLE_ISNAN(x) (DOUBLE_ISNANorINF(x) && \
344
((((Ieee_double_extract *)&(x))->s.manl != 0) || \
345
(((Ieee_double_extract *)&(x))->s.manh != 0)))
347
# define DOUBLE_ISINF(x) ((x) > DBL_MAX || (x) < -DBL_MAX)
349
/* Avoid MIPSpro / IRIX64 (incorrect) optimizations.
350
The + must not be replaced by a ||. */
351
# define DOUBLE_ISNAN(x) (!(((x) >= 0.0) + ((x) <= 0.0)))
353
# define DOUBLE_ISNAN(x) ((x) != (x))
359
/******************************************************
360
*************** Long double macros *******************
361
******************************************************/
363
/* We try to get the exact value of the precision of long double
364
(provided by the implementation) in order to provide correct
365
rounding in this case (not guaranteed if the C implementation
366
does not have an adequate long double arithmetic). Note that
367
it may be lower than the precision of some numbers that can
368
be represented in a long double; e.g. on FreeBSD/x86, it is
369
53 because the processor is configured to round in double
370
precision (even when using the long double type -- this is a
371
limitation of the x87 arithmetic), and on Mac OS X, it is 106
372
because the implementation is a double-double arithmetic.
373
Otherwise (e.g. in base 10), we get an upper bound of the
374
precision, and correct rounding isn't currently provided.
376
#if LDBL_MANT_DIG && FLT_RADIX == 2
377
# define MPFR_LDBL_MANT_DIG LDBL_MANT_DIG
379
# define MPFR_LDBL_MANT_DIG \
380
(sizeof(long double)*BITS_PER_MP_LIMB/sizeof(mp_limb_t))
382
#define MPFR_LIMBS_PER_LONG_DOUBLE \
383
((sizeof(long double)-1)/sizeof(mp_limb_t)+1)
385
/* LONGDOUBLE_NAN_ACTION executes the code "action" if x is a NaN. */
387
/* On hppa2.0n-hp-hpux10 with the unbundled HP cc, the test x!=x on a NaN
388
has been seen false, meaning NaNs are not detected. This seemed to
389
happen only after other comparisons, not sure what's really going on. In
390
any case we can pick apart the bytes to identify a NaN. */
391
#ifdef HAVE_LDOUBLE_IEEE_QUAD_BIG
392
# define LONGDOUBLE_NAN_ACTION(x, action) \
397
unsigned int sign : 1; \
398
unsigned int exp : 15; \
399
unsigned int man3 : 16; \
400
unsigned int man2 : 32; \
401
unsigned int man1 : 32; \
402
unsigned int man0 : 32; \
406
if (u.s.exp == 0x7FFFL \
407
&& (u.s.man0 | u.s.man1 | u.s.man2 | u.s.man3) != 0) \
412
/* Under IEEE rules, NaN is not equal to anything, including itself.
413
"volatile" here stops "cc" on mips64-sgi-irix6.5 from optimizing away
415
#ifndef LONGDOUBLE_NAN_ACTION
416
# define LONGDOUBLE_NAN_ACTION(x, action) \
418
volatile long double __x = LONGDOUBLE_VOLATILE (x); \
422
# define WANT_LONGDOUBLE_VOLATILE 1
425
/* If we don't have a proper "volatile" then volatile is #defined to empty,
426
in this case call through an external function to stop the compiler
427
optimizing anything. */
428
#ifdef WANT_LONGDOUBLE_VOLATILE
430
__MPFR_DECLSPEC long double __gmpfr_longdouble_volatile _MPFR_PROTO ((long double)) MPFR_CONST_ATTR;
431
# define LONGDOUBLE_VOLATILE(x) (__gmpfr_longdouble_volatile (x))
432
# define WANT_GMPFR_LONGDOUBLE_VOLATILE 1
434
# define LONGDOUBLE_VOLATILE(x) (x)
438
/* Some special case for IEEE_EXT Litle Endian */
439
#if HAVE_LDOUBLE_IEEE_EXT_LITTLE
444
unsigned int manl : 32;
445
unsigned int manh : 32;
446
unsigned int expl : 8 ;
447
unsigned int exph : 7;
448
unsigned int sign : 1;
450
} mpfr_long_double_t;
452
/* #undef MPFR_LDBL_MANT_DIG */
453
#undef MPFR_LIMBS_PER_LONG_DOUBLE
454
/* #define MPFR_LDBL_MANT_DIG 64 */
455
#define MPFR_LIMBS_PER_LONG_DOUBLE ((64-1)/BITS_PER_MP_LIMB+1)
459
/******************************************************
460
**************** mpfr_t properties *******************
461
******************************************************/
463
#define MPFR_PREC(x) ((x)->_mpfr_prec)
464
#define MPFR_EXP(x) ((x)->_mpfr_exp)
465
#define MPFR_MANT(x) ((x)->_mpfr_d)
466
#define MPFR_LIMB_SIZE(x) ((MPFR_PREC((x))-1)/BITS_PER_MP_LIMB+1)
468
#if _MPFR_PREC_FORMAT == 1
469
# define MPFR_INTPREC_MAX (USHRT_MAX & ~(unsigned int) (BITS_PER_MP_LIMB - 1))
470
#elif _MPFR_PREC_FORMAT == 2
471
# define MPFR_INTPREC_MAX (UINT_MAX & ~(unsigned int) (BITS_PER_MP_LIMB - 1))
472
#elif _MPFR_PREC_FORMAT == 3
473
# define MPFR_INTPREC_MAX (ULONG_MAX & ~(unsigned long) (BITS_PER_MP_LIMB - 1))
475
# error "Invalid MPFR Prec format"
480
/******************************************************
481
***************** exponent limits ********************
482
******************************************************/
484
/* Defined limits and unsigned type of exponent */
485
#if __GMP_MP_SIZE_T_INT == 1
486
typedef unsigned int mpfr_uexp_t;
487
# define MPFR_EXP_MAX (INT_MAX)
488
# define MPFR_EXP_MIN (INT_MIN)
490
typedef unsigned long int mpfr_uexp_t;
491
# define MPFR_EXP_MAX (LONG_MAX)
492
# define MPFR_EXP_MIN (LONG_MIN)
494
#ifndef mp_exp_unsigned_t
495
# define mp_exp_unsigned_t mpfr_uexp_t
498
/* Invalid exponent value (to track bugs...) */
499
#define MPFR_EXP_INVALID \
500
((mp_exp_t) 1 << (BITS_PER_MP_LIMB*sizeof(mp_exp_t)/sizeof(mp_limb_t)-2))
502
/* Definition of the intervals of the exponent limits */
507
#define MPFR_EMIN_MIN (1-MPFR_EXP_INVALID)
508
#define MPFR_EMIN_MAX (MPFR_EXP_INVALID-1)
509
#define MPFR_EMAX_MIN (1-MPFR_EXP_INVALID)
510
#define MPFR_EMAX_MAX (MPFR_EXP_INVALID-1)
512
/* Use MPFR_GET_EXP and MPFR_SET_EXP instead of MPFR_EXP directly,
513
unless when the exponent may be out-of-range, for instance when
514
setting the exponent before calling mpfr_check_range.
515
MPFR_EXP_CHECK is defined when WANT_ASSERT is defined, but if you
516
don't use WANT_ASSERT (for speed reasons), you can still define
517
MPFR_EXP_CHECK by setting -DMPFR_EXP_CHECK in $CFLAGS. */
519
#ifdef MPFR_EXP_CHECK
520
# define MPFR_GET_EXP(x) (mpfr_get_exp) (x)
521
# define MPFR_SET_EXP(x, exp) MPFR_ASSERTN (!mpfr_set_exp ((x), (exp)))
522
# define MPFR_SET_INVALID_EXP(x) ((void) (MPFR_EXP (x) = MPFR_EXP_INVALID))
524
# define MPFR_GET_EXP(x) MPFR_EXP (x)
525
# define MPFR_SET_EXP(x, exp) ((void) (MPFR_EXP (x) = (exp)))
526
# define MPFR_SET_INVALID_EXP(x) ((void) 0)
531
/******************************************************
532
********** Singular Values (NAN, INF, ZERO) **********
533
******************************************************/
536
* Clear flags macros are still defined and should be still used
537
* since the functions must not assume the internal format.
538
* How to deal with special values ?
539
* 1. Check if is a special value (Zero, Nan, Inf) wiht MPFR_IS_SINGULAR
540
* 2. Deal with the special value with MPFR_IS_NAN, MPFR_IS_INF, etc
541
* 3. Else clear the flags of the dest (it must be done after since src
542
* may be also the dest!)
543
* MPFR_SET_INF, MPFR_SET_NAN, MPFR_SET_ZERO must clear by
544
* themselves the other flags.
547
/* Enum special value of exponent.*/
548
# define MPFR_EXP_ZERO (MPFR_EXP_MIN+1)
549
# define MPFR_EXP_NAN (MPFR_EXP_MIN+2)
550
# define MPFR_EXP_INF (MPFR_EXP_MIN+3)
552
#define MPFR_CLEAR_FLAGS(x)
554
#define MPFR_IS_NAN(x) (MPFR_EXP(x) == MPFR_EXP_NAN)
555
#define MPFR_SET_NAN(x) (MPFR_EXP(x) = MPFR_EXP_NAN)
556
#define MPFR_IS_INF(x) (MPFR_EXP(x) == MPFR_EXP_INF)
557
#define MPFR_SET_INF(x) (MPFR_EXP(x) = MPFR_EXP_INF)
558
#define MPFR_IS_ZERO(x) (MPFR_EXP(x) == MPFR_EXP_ZERO)
559
#define MPFR_SET_ZERO(x) (MPFR_EXP(x) = MPFR_EXP_ZERO)
560
#define MPFR_NOTZERO(x) (MPFR_EXP(x) != MPFR_EXP_ZERO)
562
#define MPFR_IS_FP(x) (!MPFR_IS_NAN(x) && !MPFR_IS_INF(x))
563
#define MPFR_IS_SINGULAR(x) (MPFR_EXP(x) <= MPFR_EXP_INF)
564
#define MPFR_IS_PURE_FP(x) (!MPFR_IS_SINGULAR(x))
566
#define MPFR_ARE_SINGULAR(x,y) \
567
(MPFR_UNLIKELY(MPFR_IS_SINGULAR(x)) || MPFR_UNLIKELY(MPFR_IS_SINGULAR(y)))
571
/******************************************************
572
********************* Sign Macros ********************
573
******************************************************/
575
#define MPFR_SIGN_POS (1)
576
#define MPFR_SIGN_NEG (-1)
578
#define MPFR_IS_STRICTPOS(x) (MPFR_NOTZERO((x)) && MPFR_SIGN(x) > 0)
579
#define MPFR_IS_STRICTNEG(x) (MPFR_NOTZERO((x)) && MPFR_SIGN(x) < 0)
581
#define MPFR_IS_NEG(x) (MPFR_SIGN(x) < 0)
582
#define MPFR_IS_POS(x) (MPFR_SIGN(x) > 0)
584
#define MPFR_SET_POS(x) (MPFR_SIGN(x) = MPFR_SIGN_POS)
585
#define MPFR_SET_NEG(x) (MPFR_SIGN(x) = MPFR_SIGN_NEG)
587
#define MPFR_CHANGE_SIGN(x) (MPFR_SIGN(x) = -MPFR_SIGN(x))
588
#define MPFR_SET_SAME_SIGN(x, y) (MPFR_SIGN(x) = MPFR_SIGN(y))
589
#define MPFR_SET_OPPOSITE_SIGN(x, y) (MPFR_SIGN(x) = -MPFR_SIGN(y))
590
#define MPFR_ASSERT_SIGN(s) \
591
(MPFR_ASSERTD((s) == MPFR_SIGN_POS || (s) == MPFR_SIGN_NEG))
592
#define MPFR_SET_SIGN(x, s) \
593
(MPFR_ASSERT_SIGN(s), MPFR_SIGN(x) = s)
594
#define MPFR_IS_POS_SIGN(s1) (s1 > 0)
595
#define MPFR_IS_NEG_SIGN(s1) (s1 < 0)
596
#define MPFR_MULT_SIGN(s1, s2) ((s1) * (s2))
597
/* Transform a sign to 1 or -1 */
598
#define MPFR_FROM_SIGN_TO_INT(s) (s)
599
#define MPFR_INT_SIGN(x) MPFR_FROM_SIGN_TO_INT(MPFR_SIGN(x))
603
/******************************************************
604
***************** Ternary Value Macros ***************
605
******************************************************/
607
/* Special inexact value */
608
#define MPFR_EVEN_INEX 2
610
/* When returning the ternary inexact value, ALWAYS use one of the
611
following two macros, unless the flag comes from another function
612
returning the ternary inexact value */
613
#define MPFR_RET(I) return \
614
(I) ? ((__gmpfr_flags |= MPFR_FLAGS_INEXACT), (I)) : 0
615
#define MPFR_RET_NAN return (__gmpfr_flags |= MPFR_FLAGS_NAN), 0
617
#define MPFR_SET_ERANGE() (__gmpfr_flags |= MPFR_FLAGS_ERANGE)
621
/******************************************************
622
************** Rounding mode macros *****************
623
******************************************************/
625
/* We want to test this :
626
* (rnd == GMP_RNDU && test) || (rnd == RNDD && !test)
627
* ie it transforms RNDU or RNDD to Away or Zero according to the sign */
628
#define MPFR_IS_RNDUTEST_OR_RNDDNOTTEST(rnd, test) \
629
(((rnd) + (test)) == GMP_RNDD)
631
/* We want to test if rnd = Zero, or Away.
632
'test' is true iff negative. */
633
#define MPFR_IS_LIKE_RNDZ(rnd, test) \
634
((rnd==GMP_RNDZ) || MPFR_IS_RNDUTEST_OR_RNDDNOTTEST (rnd, test))
636
/* Invert a rounding mode */
637
#define MPFR_INVERT_RND(rnd) ((rnd == GMP_RNDU) ? GMP_RNDD : \
638
((rnd == GMP_RNDD) ? GMP_RNDU : rnd))
640
/* Transform RNDU and RNDD to RNDA or RNDZ */
641
#define MPFR_UPDATE_RND_MODE(rnd, test) \
643
if (MPFR_UNLIKELY(MPFR_IS_RNDUTEST_OR_RNDDNOTTEST(rnd, test))) \
649
/******************************************************
650
******************* Limb Macros **********************
651
******************************************************/
653
/* Definition of MPFR_LIMB_HIGHBIT */
654
#if defined(GMP_LIMB_HIGHBIT)
655
# define MPFR_LIMB_HIGHBIT GMP_LIMB_HIGHBIT
656
#elif defined(MP_LIMB_T_HIGHBIT)
657
# define MPFR_LIMB_HIGHBIT MP_LIMB_T_HIGHBIT
659
# error "Neither GMP_LIMB_HIGHBIT nor MP_LIMB_T_HIGHBIT defined in GMP"
662
/* Mask to get the Most Significent Bit of a limb */
663
#define MPFR_LIMB_MSB(l) ((l)&MPFR_LIMB_HIGHBIT)
665
/* Definition of MPFR_LIMB_ONE & MPFR_LIMB_ZERO*/
667
# define MPFR_LIMB_ONE CNST_LIMB(1)
668
# define MPFR_LIMB_ZERO CNST_LIMB(0)
670
# define MPFR_LIMB_ONE ((mp_limb_t) 1L)
671
# define MPFR_LIMB_ZERO ((mp_limb_t) 0L)
674
/* Mask for the low 's' bits of a limb */
675
#define MPFR_LIMB_MASK(s) ((MPFR_LIMB_ONE<<(s))-MPFR_LIMB_ONE)
679
/******************************************************
680
********************** Memory ************************
681
******************************************************/
683
/* Heap Memory gestion */
684
typedef union { mp_size_t s; mp_limb_t l; } mpfr_size_limb_t;
685
#define MPFR_GET_ALLOC_SIZE(x) \
686
( ((mp_size_t*) MPFR_MANT(x))[-1] + 0)
687
#define MPFR_SET_ALLOC_SIZE(x, n) \
688
( ((mp_size_t*) MPFR_MANT(x))[-1] = n)
689
#define MPFR_MALLOC_SIZE(s) \
690
( sizeof(mpfr_size_limb_t) + BYTES_PER_MP_LIMB * ((size_t) s) )
691
#define MPFR_SET_MANT_PTR(x,p) \
692
(MPFR_MANT(x) = (mp_limb_t*) ((mpfr_size_limb_t*) p + 1))
693
#define MPFR_GET_REAL_PTR(x) \
694
((mp_limb_t*) ((mpfr_size_limb_t*) MPFR_MANT(x) - 1))
696
/* Temporary memory gestion */
698
/* GMP 4.1.x or below or internals */
699
#define MPFR_TMP_DECL TMP_DECL
700
#define MPFR_TMP_MARK TMP_MARK
701
#define MPFR_TMP_ALLOC TMP_ALLOC
702
#define MPFR_TMP_FREE TMP_FREE
704
#define MPFR_TMP_DECL(x) TMP_DECL
705
#define MPFR_TMP_MARK(x) TMP_MARK
706
#define MPFR_TMP_ALLOC(s) TMP_SALLOC(s)
707
#define MPFR_TMP_FREE(x) TMP_FREE
710
/* This code is experimental: don't use it */
711
#ifdef MPFR_USE_OWN_MPFR_TMP_ALLOC
712
extern unsigned char *mpfr_stack;
715
#undef MPFR_TMP_ALLOC
717
#define MPFR_TMP_DECL(_x) unsigned char *(_x)
718
#define MPFR_TMP_MARK(_x) ((_x) = mpfr_stack)
719
#define MPFR_TMP_ALLOC(_s) (mpfr_stack += (_s), mpfr_stack - (_s))
720
#define MPFR_TMP_FREE(_x) (mpfr_stack = (_x))
723
/* temporary allocate 1 limb at xp, and initialize mpfr variable x */
724
/* The temporary var doesn't have any size field, but it doesn't matter
725
* since only functions dealing with the Heap care about it */
726
#define MPFR_TMP_INIT1(xp, x, p) \
727
( MPFR_PREC(x) = (p), \
728
MPFR_MANT(x) = (xp), \
730
MPFR_SET_INVALID_EXP(x))
732
#define MPFR_TMP_INIT(xp, x, p, s) \
733
(xp = (mp_ptr) MPFR_TMP_ALLOC(BYTES_PER_MP_LIMB * ((size_t) s)), \
734
MPFR_TMP_INIT1(xp, x, p))
736
#define MPFR_TMP_INIT_ABS(d, s) \
737
( MPFR_PREC(d) = MPFR_PREC(s), \
738
MPFR_MANT(d) = MPFR_MANT(s), \
740
MPFR_EXP(d) = MPFR_EXP(s))
744
/******************************************************
745
***************** Cache macros **********************
746
******************************************************/
748
#define mpfr_const_pi(_d,_r) mpfr_cache(_d, __gmpfr_cache_const_pi,_r)
749
#define mpfr_const_log2(_d,_r) mpfr_cache(_d, __gmpfr_cache_const_log2, _r)
750
#define mpfr_const_euler(_d,_r) mpfr_cache(_d, __gmpfr_cache_const_euler, _r)
751
#define mpfr_const_catalan(_d,_r) mpfr_cache(_d,__gmpfr_cache_const_catalan,_r)
753
#define MPFR_DECL_INIT_CACHE(_cache,_func) \
754
mpfr_cache_t MPFR_THREAD_ATTR _cache = \
755
{{{{0,MPFR_SIGN_POS,0,(mp_limb_t*)0}},0,_func}}
759
/******************************************************
760
******************* Threshold ***********************
761
******************************************************/
765
/******************************************************
766
***************** Useful macros *********************
767
******************************************************/
769
/* Theses macros help the compiler to determine if a test is
770
* likely or unlikely. */
771
#if __MPFR_GNUC(3,0) || __MPFR_ICC(8,1,0)
772
# define MPFR_LIKELY(x) (__builtin_expect(!!(x),1))
773
# define MPFR_UNLIKELY(x) (__builtin_expect((x),0))
775
# define MPFR_LIKELY(x) (x)
776
# define MPFR_UNLIKELY(x) (x)
779
/* Ceil log 2: If GCC, uses a GCC extension, otherwise calls a function */
781
* Needs to define MPFR_NEED_LONGLONG.
782
* Computes ceil(log2(x)) only for x integer (unsigned long)
783
* Undefined if x is 0 */
784
#if __MPFR_GNUC(2,95) || __MPFR_ICC(8,1,0)
785
# define MPFR_INT_CEIL_LOG2(x) \
786
(__extension__ ({int _b; mp_limb_t _limb = (x); \
787
MPFR_ASSERTN (_limb == (x)); \
788
count_leading_zeros (_b, _limb); \
789
(BITS_PER_MP_LIMB - _b); }))
791
# define MPFR_INT_CEIL_LOG2(x) (__gmpfr_int_ceil_log2(x))
794
/* Add two integers with overflow handling */
795
/* Example: MPFR_SADD_OVERFLOW (c, a, b, long, unsigned long,
796
* LONG_MIN, LONG_MAX,
797
* goto overflow, goto underflow); */
798
#define MPFR_UADD_OVERFLOW(c,a,b,ACTION_IF_OVERFLOW) \
801
if ((c) < (a)) ACTION_IF_OVERFLOW; \
804
#define MPFR_SADD_OVERFLOW(c,a,b,STYPE,UTYPE,MIN,MAX,ACTION_IF_POS_OVERFLOW,ACTION_IF_NEG_OVERFLOW) \
806
if ((a) >= 0 && (b) >= 0) { \
808
ua = (UTYPE) a; ub = (UTYPE) b; \
809
MPFR_UADD_OVERFLOW (uc, ua, ub, ACTION_IF_POS_OVERFLOW); \
810
if (uc > (UTYPE)(MAX)) ACTION_IF_POS_OVERFLOW; \
811
else (c) = (STYPE) uc; \
812
} else if ((a) < 0 && (b) < 0) { \
814
ua = -(UTYPE) a; ub = -(UTYPE) b; \
815
MPFR_UADD_OVERFLOW (uc, ua, ub, ACTION_IF_NEG_OVERFLOW); \
816
if (uc >= -(UTYPE)(MIN) || uc > (UTYPE)(MAX)) { \
817
if (uc == -(UTYPE)(MIN)) (c) = (MIN); \
818
else ACTION_IF_NEG_OVERFLOW; } \
819
else (c) = -(STYPE) uc; \
820
} else (c) = (a) + (b); \
824
/* Set a number to 1 (Fast) - It doesn't check if 1 is in the exponent range */
825
#define MPFR_SET_ONE(x) \
827
mp_size_t _size = MPFR_LIMB_SIZE(x) - 1; \
830
MPN_ZERO ( MPFR_MANT(x), _size); \
831
MPFR_MANT(x)[_size] = MPFR_LIMB_HIGHBIT; \
834
/* Compute s = (-a) % BITS_PER_MP_LIMB
835
* a is unsigned! Check if it works,
836
* otherwise tries another way to compute it */
837
#define MPFR_UNSIGNED_MINUS_MODULO(s, a) \
840
if ((UINT_MAX % BITS_PER_MP_LIMB) == (BITS_PER_MP_LIMB-1)) \
841
(s) = (mpfr_prec_t) (-(a)) % BITS_PER_MP_LIMB; \
844
(s) = (a) % BITS_PER_MP_LIMB; \
846
(s) = BITS_PER_MP_LIMB - (s); \
848
MPFR_ASSERTD ((s) >= 0 && (s) < BITS_PER_MP_LIMB); \
852
/* Use it only for debug reasons */
853
/* MPFR_TRACE (operation) : execute operation iff DEBUG flag is set */
854
/* MPFR_DUMP (x) : print x (a mpfr_t) on stdout */
857
# define MPFR_TRACE(x) x
859
# define MPFR_TRACE(x) (void) 0
861
#define MPFR_DUMP(x) ( printf(#x"="), mpfr_dump(x) )
863
/* Test if X (positive) is a power of 2 */
864
#define IS_POW2(X) (((X) & ((X) - 1)) == 0)
865
#define NOT_POW2(X) (((X) & ((X) - 1)) != 0)
867
/* Safe absolute value (to avoid possible integer overflow) */
868
/* type is the target (unsigned) type */
869
#define SAFE_ABS(type,x) ((x) >= 0 ? (type)(x) : -(type)(x))
871
#define mpfr_get_d1(x) mpfr_get_d(x,__gmpfr_default_rounding_mode)
873
/* Store in r the size in bits of the mpz_t z */
874
#define MPFR_MPZ_SIZEINBASE2(r, z) \
878
MPFR_ASSERTD (mpz_sgn (z) != 0); \
880
count_leading_zeros (_cnt, PTR(z)[_size-1]); \
881
(r) = _size * BITS_PER_MP_LIMB - _cnt; \
884
/* Needs <locale.h> */
885
#define MPFR_DECIMAL_POINT ((unsigned char) localeconv()->decimal_point[0])
887
/******************************************************
888
************** Save exponent macros ****************
889
******************************************************/
891
/* See README.dev for details on how to use the macros.
892
They are used to set the exponent range to the maximum
896
unsigned int saved_flags;
901
#define MPFR_SAVE_EXPO_DECL(x) mpfr_save_expo_t x
902
#define MPFR_SAVE_EXPO_MARK(x) \
903
((x).saved_flags = __gmpfr_flags, \
904
(x).saved_emin = __gmpfr_emin, \
905
(x).saved_emax = __gmpfr_emax, \
906
__gmpfr_emin = MPFR_EMIN_MIN, \
907
__gmpfr_emax = MPFR_EMAX_MAX)
908
#define MPFR_SAVE_EXPO_FREE(x) \
909
(__gmpfr_flags = (x).saved_flags, \
910
__gmpfr_emin = (x).saved_emin, \
911
__gmpfr_emax = (x).saved_emax)
912
#define MPFR_SAVE_EXPO_UPDATE_FLAGS(x, flags) \
913
(x).saved_flags |= (flags)
915
/* Speed up final checking */
916
#define mpfr_check_range(x,t,r) \
917
(MPFR_LIKELY (MPFR_EXP (x) >= __gmpfr_emin && MPFR_EXP (x) <= __gmpfr_emax) \
918
? (t) : mpfr_check_range(x,t,r))
921
/******************************************************
922
***************** Inline Rounding *******************
923
******************************************************/
926
* Note: due to the labels, one cannot use a macro MPFR_RNDRAW* more than
927
* once in a function (otherwise these labels would not be unique).
931
* Round mantissa (srcp, sprec) to mpfr_t dest using rounding mode rnd
932
* assuming dest's sign is sign.
933
* In rounding to nearest mode, execute MIDDLE_HANDLER when the value
934
* is the middle of two consecutive numbers in dest precision.
935
* Execute OVERFLOW_HANDLER in case of overflow when rounding.
937
#define MPFR_RNDRAW_GEN(inexact, dest, srcp, sprec, rnd, sign, \
938
MIDDLE_HANDLER, OVERFLOW_HANDLER) \
940
mp_size_t dests, srcs; \
942
mp_prec_t destprec, srcprec; \
944
/* Check Trivial Case when Dest Mantissa has more bits than source */ \
946
destprec = MPFR_PREC (dest); \
947
destp = MPFR_MANT (dest); \
948
if (MPFR_UNLIKELY (destprec >= srcprec)) \
950
srcs = (srcprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB; \
951
dests = (destprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB - srcs; \
952
MPN_COPY (destp + dests, srcp, srcs); \
953
MPN_ZERO (destp, dests); \
958
/* Non trivial case: rounding needed */ \
961
mp_limb_t rb, sb, ulp; \
963
/* Compute Position and shift */ \
964
srcs = (srcprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB; \
965
dests = (destprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB; \
966
MPFR_UNSIGNED_MINUS_MODULO (sh, destprec); \
967
sp = srcp + srcs - dests; \
969
/* General case when prec % BITS_PER_MP_LIMB != 0 */ \
970
if (MPFR_LIKELY (sh != 0)) \
973
/* Compute Rounding Bit and Sticky Bit */ \
974
mask = MPFR_LIMB_ONE << (sh-1); \
976
sb = sp[0] & (mask-1); \
977
if (MPFR_UNLIKELY (sb == 0)) \
978
{ /* TODO: Improve it */ \
981
for (tmp = sp, n = srcs - dests ; n != 0 && sb == 0 ; n--) \
988
MPFR_ASSERTD (dests < srcs); \
989
/* Compute Rounding Bit and Sticky Bit */ \
990
rb = sp[-1] & MPFR_LIMB_HIGHBIT; \
991
sb = sp[-1] & (MPFR_LIMB_HIGHBIT-1); \
992
if (MPFR_UNLIKELY (sb == 0)) \
996
for (tmp = sp-1, n = srcs - dests-1 ; n!=0 && sb==0 ; n--) \
999
ulp = MPFR_LIMB_ONE; \
1002
if (MPFR_LIKELY (rnd == GMP_RNDN)) \
1007
inexact = MPFR_LIKELY ((sb | rb) != 0) ? -sign : 0; \
1009
MPN_COPY (destp, sp, dests); \
1010
destp[0] &= ~(ulp-1); \
1012
else if (MPFR_UNLIKELY (sb == 0)) \
1013
{ /* Middle of two consecutive representable numbers */ \
1021
if (MPFR_UNLIKELY (mpn_add_1 (destp, sp, dests, ulp))) \
1023
destp[dests-1] = MPFR_LIMB_HIGHBIT; \
1026
destp[0] &= ~(ulp-1); \
1030
{ /* Directed rounding mode */ \
1031
if (MPFR_LIKELY (MPFR_IS_LIKE_RNDZ (rnd, \
1032
MPFR_IS_NEG_SIGN (sign)))) \
1034
else if (MPFR_UNLIKELY ((sb | rb) == 0)) \
1046
* Round mantissa (srcp, sprec) to mpfr_t dest using rounding mode rnd
1047
* assuming dest's sign is sign.
1048
* Execute OVERFLOW_HANDLER in case of overflow when rounding.
1050
#define MPFR_RNDRAW(inexact, dest, srcp, sprec, rnd, sign, OVERFLOW_HANDLER) \
1051
MPFR_RNDRAW_GEN (inexact, dest, srcp, sprec, rnd, sign, \
1052
if ((sp[0] & ulp) == 0) \
1062
* Round mantissa (srcp, sprec) to mpfr_t dest using rounding mode rnd
1063
* assuming dest's sign is sign.
1064
* Execute OVERFLOW_HANDLER in case of overflow when rounding.
1065
* Set inexact to +/- MPFR_EVEN_INEX in case of even rounding.
1067
#define MPFR_RNDRAW_EVEN(inexact, dest, srcp, sprec, rnd, sign, \
1069
MPFR_RNDRAW_GEN (inexact, dest, srcp, sprec, rnd, sign, \
1070
if ((sp[0] & ulp) == 0) \
1072
inexact = -MPFR_EVEN_INEX * sign; \
1077
inexact = MPFR_EVEN_INEX * sign; \
1078
goto addoneulp_doit; \
1082
/* Return TRUE if b is non singular and we can round it to precision 'prec'
1083
with rounding mode 'rnd', and with error at most 'error' */
1084
#define MPFR_CAN_ROUND(b,err,prec,rnd) \
1085
(!MPFR_IS_SINGULAR (b) && mpfr_round_p (MPFR_MANT (b), MPFR_LIMB_SIZE (b), \
1086
(err), (prec) + ((rnd)==GMP_RNDN)))
1088
/* Assuming that the function as a taylor expansion which looks like:
1089
y=o(f(x)) = o(x + g(x)) with |g(x)| <=2^(EXP(x)-err)
1090
we can quickly set y to x if x is small (ie err > prec(y)+1) in most
1091
cases. It assumes that f(x) is not representable exactly as a FP number.
1092
x must not be a singular value (NAN, INF or ZERO).
1094
y is the destination (a mpfr_t), x the value to set (a mpfr_t),
1095
err the error term (a mp_exp_t), dir (an int) is the direction of
1096
the commited error (if dir = 0, it rounds towards 0, if dir=1,
1097
it rounds away from 0), rnd the rounding mode.
1099
It returns from the function a ternary value in case of success.
1100
If you want to free something, you must fill the "extra" field
1101
in consequences, otherwise put nothing in it.
1103
The test is less restrictive thant necessary, but the function
1104
will finish the check itself.
1106
#define MPFR_FAST_COMPUTE_IF_SMALL_INPUT(y,x,err,dir,rnd,extra) \
1108
mp_exp_t _err = (err); \
1109
if (MPFR_UNLIKELY (_err > 0 \
1110
&& (mpfr_uexp_t) _err > MPFR_PREC (y) + 1)) \
1112
int _inexact = mpfr_round_near_x ((y),(x),(err),(dir),(rnd)); \
1113
if (_inexact != 0) \
1121
/******************************************************
1122
*************** Ziv Loop Macro *********************
1123
******************************************************/
1125
#ifndef MPFR_USE_LOGGING
1127
#define MPFR_ZIV_DECL(_x) mp_prec_t _x
1128
#define MPFR_ZIV_INIT(_x, _p) (_x) = BITS_PER_MP_LIMB
1129
#define MPFR_ZIV_NEXT(_x, _p) ((_p) += (_x), (_x) = (_p)/2)
1130
#define MPFR_ZIV_FREE(x)
1134
#define MPFR_ZIV_DECL(_x) \
1136
int _x ## _cpt = 1; \
1137
static unsigned long _x ## _loop = 0, _x ## _bad = 0; \
1138
static const char *_x ## _fname = __func__; \
1139
auto void __attribute__ ((destructor)) x ## _f (void); \
1140
void __attribute__ ((destructor)) x ## _f (void) { \
1141
if (_x ## _loop != 0 && MPFR_LOG_STAT_F&mpfr_log_type) \
1142
fprintf (mpfr_log_file, \
1143
"%s: Ziv failed %2.2f%% (%lu bad cases / %lu calls)\n", _x ## _fname, \
1144
(double) 100.0 * _x ## _bad / _x ## _loop, _x ## _bad, _x ## _loop ); }
1146
#define MPFR_ZIV_INIT(_x, _p) ((_x) = BITS_PER_MP_LIMB, _x ## _loop ++); \
1147
if (MPFR_LOG_BADCASE_F&mpfr_log_type && mpfr_log_current<=mpfr_log_level) \
1148
fprintf (mpfr_log_file, "%s:ZIV 1st prec=%lu\n", __func__, \
1149
(unsigned long) (_p))
1151
#define MPFR_ZIV_NEXT(_x, _p) \
1152
((_p)+=(_x),(_x)=(_p)/2, _x ## _bad += (_x ## _cpt == 1), _x ## _cpt ++); \
1153
if (MPFR_LOG_BADCASE_F&mpfr_log_type && mpfr_log_current<=mpfr_log_level) \
1154
fprintf (mpfr_log_file, "%s:ZIV new prec=%lu\n", __func__, \
1155
(unsigned long) (_p))
1157
#define MPFR_ZIV_FREE(_x) \
1158
if (MPFR_LOG_BADCASE_F&mpfr_log_type && _x##_cpt>1 \
1159
&& mpfr_log_current<=mpfr_log_level) \
1160
fprintf (mpfr_log_file, "%s:ZIV %d loops\n", __func__, _x ## _cpt)
1165
/******************************************************
1166
*************** Logging Macros *********************
1167
******************************************************/
1169
/* The different kind of LOG */
1170
#define MPFR_LOG_INPUT_F 1
1171
#define MPFR_LOG_OUTPUT_F 2
1172
#define MPFR_LOG_INTERNAL_F 4
1173
#define MPFR_LOG_TIME_F 8
1174
#define MPFR_LOG_BADCASE_F 16
1175
#define MPFR_LOG_MSG_F 32
1176
#define MPFR_LOG_STAT_F 64
1178
#ifdef MPFR_USE_LOGGING
1182
/* Check if we can support this feature */
1183
# ifdef MPFR_USE_THREAD_SAFE
1184
# error "Enable either `Logging' or `thread-safe', not both"
1186
# if !__MPFR_GNUC(3,0)
1187
# error "Logging not supported (GCC >= 3.0)"
1190
#if defined (__cplusplus)
1194
__MPFR_DECLSPEC extern FILE *mpfr_log_file;
1195
__MPFR_DECLSPEC extern int mpfr_log_type;
1196
__MPFR_DECLSPEC extern int mpfr_log_level;
1197
__MPFR_DECLSPEC extern int mpfr_log_current;
1198
__MPFR_DECLSPEC extern int mpfr_log_base;
1199
__MPFR_DECLSPEC extern mp_prec_t mpfr_log_prec;
1201
#if defined (__cplusplus)
1205
#define MPFR_LOG_VAR(x) \
1206
if((MPFR_LOG_INTERNAL_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level))\
1207
fprintf (mpfr_log_file, "%s.%d:%s[%#R]=%R\n", __func__,__LINE__, #x, x, x);
1209
#define MPFR_LOG_MSG2(format, ...) \
1210
if ((MPFR_LOG_MSG_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1211
fprintf (mpfr_log_file, "%s.%d:"format, __func__, __LINE__, __VA_ARGS__);
1212
#define MPFR_LOG_MSG(x) MPFR_LOG_MSG2 x
1214
#define MPFR_LOG_BEGIN2(format, ...) \
1215
mpfr_log_current ++; \
1216
if ((MPFR_LOG_INPUT_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1217
fprintf (mpfr_log_file, "%s:IN "format"\n",__func__,__VA_ARGS__); \
1218
if ((MPFR_LOG_TIME_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1219
__gmpfr_log_time = mpfr_get_cputime ();
1220
#define MPFR_LOG_BEGIN(x) \
1221
int __gmpfr_log_time = 0; \
1224
#define MPFR_LOG_END2(format, ...) \
1225
if ((MPFR_LOG_TIME_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1226
fprintf (mpfr_log_file, "%s:TIM %dms\n", __mpfr_log_fname, \
1227
mpfr_get_cputime () - __gmpfr_log_time); \
1228
if ((MPFR_LOG_OUTPUT_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1229
fprintf (mpfr_log_file, "%s:OUT "format"\n",__mpfr_log_fname,__VA_ARGS__);\
1230
mpfr_log_current --;
1231
#define MPFR_LOG_END(x) \
1232
static const char *__mpfr_log_fname = __func__; \
1235
#define MPFR_LOG_FUNC(begin,end) \
1236
static const char *__mpfr_log_fname = __func__; \
1237
auto void __mpfr_log_cleanup (int *time); \
1238
void __mpfr_log_cleanup (int *time) { \
1239
int __gmpfr_log_time = *time; \
1240
MPFR_LOG_END2 end; } \
1241
int __gmpfr_log_time __attribute__ ((cleanup (__mpfr_log_cleanup))); \
1242
__gmpfr_log_time = 0; \
1243
MPFR_LOG_BEGIN2 begin
1245
#else /* MPFR_USE_LOGGING */
1247
/* Define void macro for logging */
1249
#define MPFR_LOG_VAR(x)
1250
#define MPFR_LOG_BEGIN(x)
1251
#define MPFR_LOG_END(x)
1252
#define MPFR_LOG_MSG(x)
1253
#define MPFR_LOG_FUNC(x,y)
1255
#endif /* MPFR_USE_LOGGING */
1258
/**************************************************************
1259
************ Group Initialize Functions Macros *************
1260
**************************************************************/
1262
#ifndef MPFR_GROUP_STATIC_SIZE
1263
# define MPFR_GROUP_STATIC_SIZE 16
1266
struct mpfr_group_t {
1269
mp_limb_t tab[MPFR_GROUP_STATIC_SIZE];
1272
#define MPFR_GROUP_DECL(g) struct mpfr_group_t g
1273
#define MPFR_GROUP_CLEAR(g) do { \
1274
if (MPFR_UNLIKELY ((g).alloc != 0)) { \
1275
MPFR_ASSERTD ((g).mant != (g).tab); \
1276
(*__gmp_free_func)((g).mant, (g).alloc); \
1279
#define MPFR_GROUP_INIT_TEMPLATE(g, prec, num, handler) do { \
1280
mp_prec_t _prec = (prec); \
1282
MPFR_ASSERTD (_prec >= MPFR_PREC_MIN); \
1283
if (MPFR_UNLIKELY (_prec > MPFR_PREC_MAX)) mpfr_abort_prec_max(); \
1284
_size = (mp_prec_t)(_prec +BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB;\
1285
if (MPFR_UNLIKELY (_size*(num) > MPFR_GROUP_STATIC_SIZE)) { \
1286
(g).alloc = (num)*_size*sizeof (mp_limb_t); \
1287
(g).mant = (*__gmp_allocate_func) ((g).alloc); \
1290
(g).mant = (g).tab; \
1294
#define MPFR_GROUP_TINIT(g, n, x) MPFR_TMP_INIT1 ((g).mant+_size*(n), x, _prec)
1296
#define MPFR_GROUP_INIT_1(g, prec, x) \
1297
MPFR_GROUP_INIT_TEMPLATE(g, prec, 1, MPFR_GROUP_TINIT(g, 0, x))
1298
#define MPFR_GROUP_INIT_2(g, prec, x, y) \
1299
MPFR_GROUP_INIT_TEMPLATE(g, prec, 2, \
1300
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y))
1301
#define MPFR_GROUP_INIT_3(g, prec, x, y, z) \
1302
MPFR_GROUP_INIT_TEMPLATE(g, prec, 3, \
1303
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1304
MPFR_GROUP_TINIT(g, 2, z))
1305
#define MPFR_GROUP_INIT_4(g, prec, x, y, z, t) \
1306
MPFR_GROUP_INIT_TEMPLATE(g, prec, 4, \
1307
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1308
MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t))
1309
#define MPFR_GROUP_INIT_5(g, prec, x, y, z, t, a) \
1310
MPFR_GROUP_INIT_TEMPLATE(g, prec, 5, \
1311
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1312
MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \
1313
MPFR_GROUP_TINIT(g, 4, a))
1314
#define MPFR_GROUP_INIT_6(g, prec, x, y, z, t, a, b) \
1315
MPFR_GROUP_INIT_TEMPLATE(g, prec, 6, \
1316
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1317
MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \
1318
MPFR_GROUP_TINIT(g, 4, a);MPFR_GROUP_TINIT(g, 5, b))
1320
#define MPFR_GROUP_REPREC_TEMPLATE(g, prec, num, handler) do { \
1321
mp_prec_t _prec = (prec); \
1322
size_t _oalloc = (g).alloc; \
1324
MPFR_ASSERTD (_prec >= MPFR_PREC_MIN); \
1325
if (MPFR_UNLIKELY (_prec > MPFR_PREC_MAX)) mpfr_abort_prec_max(); \
1326
_size = (mp_prec_t)(_prec +BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB;\
1327
(g).alloc = (num)*_size*sizeof(mp_limb_t); \
1328
if (MPFR_LIKELY (_oalloc == 0)) \
1329
(g).mant = (*__gmp_allocate_func) ((g).alloc); \
1331
(g).mant = (*__gmp_reallocate_func)((g).mant,_oalloc,(g).alloc);\
1336
#define MPFR_GROUP_REPREC_1(g, prec, x) \
1337
MPFR_GROUP_REPREC_TEMPLATE(g, prec, 1, MPFR_GROUP_TINIT(g, 0, x))
1338
#define MPFR_GROUP_REPREC_2(g, prec, x, y) \
1339
MPFR_GROUP_REPREC_TEMPLATE(g, prec, 2, \
1340
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y))
1341
#define MPFR_GROUP_REPREC_3(g, prec, x, y, z) \
1342
MPFR_GROUP_REPREC_TEMPLATE(g, prec, 3, \
1343
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1344
MPFR_GROUP_TINIT(g, 2, z))
1345
#define MPFR_GROUP_REPREC_4(g, prec, x, y, z, t) \
1346
MPFR_GROUP_REPREC_TEMPLATE(g, prec, 4, \
1347
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1348
MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t))
1349
#define MPFR_GROUP_REPREC_5(g, prec, x, y, z, t, a) \
1350
MPFR_GROUP_REPREC_TEMPLATE(g, prec, 5, \
1351
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1352
MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \
1353
MPFR_GROUP_TINIT(g, 4, a))
1354
#define MPFR_GROUP_REPREC_6(g, prec, x, y, z, t, a, b) \
1355
MPFR_GROUP_REPREC_TEMPLATE(g, prec, 6, \
1356
MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1357
MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \
1358
MPFR_GROUP_TINIT(g, 4, a);MPFR_GROUP_TINIT(g, 5, b))
1361
/******************************************************
1362
*************** Internal Functions *****************
1363
******************************************************/
1365
#if defined (__cplusplus)
1369
__MPFR_DECLSPEC int mpfr_underflow _MPFR_PROTO ((mpfr_ptr, mp_rnd_t, int));
1370
__MPFR_DECLSPEC int mpfr_overflow _MPFR_PROTO ((mpfr_ptr, mp_rnd_t, int));
1372
__MPFR_DECLSPEC int mpfr_add1 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1373
mpfr_srcptr, mp_rnd_t));
1374
__MPFR_DECLSPEC int mpfr_sub1 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1375
mpfr_srcptr, mp_rnd_t));
1376
__MPFR_DECLSPEC int mpfr_add1sp _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1377
mpfr_srcptr, mp_rnd_t));
1378
__MPFR_DECLSPEC int mpfr_sub1sp _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1379
mpfr_srcptr, mp_rnd_t));
1380
__MPFR_DECLSPEC int mpfr_can_round_raw _MPFR_PROTO ((const mp_limb_t *,
1381
mp_size_t, int, mp_exp_t, mp_rnd_t, mp_rnd_t, mp_prec_t));
1383
__MPFR_DECLSPEC int mpfr_cmp2 _MPFR_PROTO ((mpfr_srcptr, mpfr_srcptr,
1386
__MPFR_DECLSPEC long __gmpfr_ceil_log2 _MPFR_PROTO ((double));
1387
__MPFR_DECLSPEC long __gmpfr_floor_log2 _MPFR_PROTO ((double));
1388
__MPFR_DECLSPEC double __gmpfr_ceil_exp2 _MPFR_PROTO ((double));
1389
__MPFR_DECLSPEC unsigned long __gmpfr_isqrt _MPFR_PROTO ((unsigned long));
1390
__MPFR_DECLSPEC unsigned long __gmpfr_cuberoot _MPFR_PROTO ((unsigned long));
1391
__MPFR_DECLSPEC int __gmpfr_int_ceil_log2 _MPFR_PROTO ((unsigned long));
1393
__MPFR_DECLSPEC int mpfr_exp_2 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,mp_rnd_t));
1394
__MPFR_DECLSPEC int mpfr_exp_3 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,mp_rnd_t));
1395
__MPFR_DECLSPEC int mpfr_powerof2_raw _MPFR_PROTO ((mpfr_srcptr));
1397
__MPFR_DECLSPEC void mpfr_setmax _MPFR_PROTO ((mpfr_ptr, mp_exp_t));
1398
__MPFR_DECLSPEC void mpfr_setmin _MPFR_PROTO ((mpfr_ptr, mp_exp_t));
1400
__MPFR_DECLSPEC long mpfr_mpn_exp _MPFR_PROTO ((mp_limb_t *, mp_exp_t *, int,
1403
#ifdef _MPFR_H_HAVE_FILE
1404
__MPFR_DECLSPEC void mpfr_fprint_binary _MPFR_PROTO ((FILE *, mpfr_srcptr));
1406
__MPFR_DECLSPEC void mpfr_print_binary _MPFR_PROTO ((mpfr_srcptr));
1407
__MPFR_DECLSPEC void mpfr_print_mant_binary _MPFR_PROTO ((const char*,
1408
const mp_limb_t*, mp_prec_t));
1409
__MPFR_DECLSPEC void mpfr_set_str_binary _MPFR_PROTO((mpfr_ptr, const char*));
1411
__MPFR_DECLSPEC int mpfr_round_raw _MPFR_PROTO ((mp_limb_t *,
1412
const mp_limb_t *, mp_prec_t, int, mp_prec_t, mp_rnd_t, int *));
1413
__MPFR_DECLSPEC int mpfr_round_raw_2 _MPFR_PROTO ((const mp_limb_t *,
1414
mp_prec_t, int, mp_prec_t, mp_rnd_t));
1415
__MPFR_DECLSPEC int mpfr_round_raw_3 _MPFR_PROTO ((const mp_limb_t *,
1416
mp_prec_t, int, mp_prec_t, mp_rnd_t, int *));
1417
__MPFR_DECLSPEC int mpfr_round_raw_4 _MPFR_PROTO ((mp_limb_t *,
1418
const mp_limb_t *, mp_prec_t, int, mp_prec_t, mp_rnd_t));
1420
#define mpfr_round_raw2(xp, xn, neg, r, prec) \
1421
mpfr_round_raw_2((xp),(xn)*BITS_PER_MP_LIMB,(neg),(prec),(r))
1423
__MPFR_DECLSPEC int mpfr_check _MPFR_PROTO ((mpfr_srcptr));
1425
__MPFR_DECLSPEC int mpfr_sum_sort _MPFR_PROTO ((mpfr_srcptr *const,
1426
unsigned long, mpfr_srcptr *));
1428
__MPFR_DECLSPEC int mpfr_get_cputime _MPFR_PROTO ((void));
1430
__MPFR_DECLSPEC void mpfr_nexttozero _MPFR_PROTO ((mpfr_ptr));
1431
__MPFR_DECLSPEC void mpfr_nexttoinf _MPFR_PROTO ((mpfr_ptr));
1433
__MPFR_DECLSPEC int mpfr_const_pi_internal _MPFR_PROTO ((mpfr_ptr,mp_rnd_t));
1434
__MPFR_DECLSPEC int mpfr_const_log2_internal _MPFR_PROTO((mpfr_ptr,mp_rnd_t));
1435
__MPFR_DECLSPEC int mpfr_const_euler_internal _MPFR_PROTO((mpfr_ptr, mp_rnd_t));
1436
__MPFR_DECLSPEC int mpfr_const_catalan_internal _MPFR_PROTO((mpfr_ptr, mp_rnd_t));
1438
__MPFR_DECLSPEC void mpfr_init_cache _MPFR_PROTO ((mpfr_cache_t,
1439
int(*)(mpfr_ptr,mpfr_rnd_t)));
1440
__MPFR_DECLSPEC void mpfr_clear_cache _MPFR_PROTO ((mpfr_cache_t));
1441
__MPFR_DECLSPEC int mpfr_cache _MPFR_PROTO ((mpfr_ptr, mpfr_cache_t,
1444
__MPFR_DECLSPEC void mpfr_mulhigh_n _MPFR_PROTO ((mp_ptr, mp_srcptr,
1445
mp_srcptr, mp_size_t));
1447
__MPFR_DECLSPEC int mpfr_round_p _MPFR_PROTO ((mp_limb_t *, mp_size_t,
1448
mp_exp_t, mp_prec_t));
1450
__MPFR_DECLSPEC void mpfr_dump_mant _MPFR_PROTO ((const mp_limb_t *,
1451
mp_prec_t, mp_prec_t,
1454
__MPFR_DECLSPEC int mpfr_round_near_x _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1455
mp_exp_t, int, mp_rnd_t));
1456
__MPFR_DECLSPEC void mpfr_abort_prec_max _MPFR_PROTO ((void))
1459
#if defined (__cplusplus)