~ubuntu-branches/debian/sid/genius/sid

« back to all changes in this revision

Viewing changes to mpfr/mpfr-impl.h

  • Committer: Bazaar Package Importer
  • Author(s): Daniel Holbach
  • Date: 2006-08-21 12:57:45 UTC
  • Revision ID: james.westby@ubuntu.com-20060821125745-sl9ks8v7fq324bdf
Tags: upstream-0.7.6.1
ImportĀ upstreamĀ versionĀ 0.7.6.1

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/* Utilities for MPFR developers, not exported.
 
2
 
 
3
Copyright 1999, 2000, 2001, 2002, 2003, 2004, 2005
 
4
  Free Software Foundation, Inc.
 
5
 
 
6
This file is part of the MPFR Library.
 
7
 
 
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.
 
12
 
 
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.
 
17
 
 
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. */
 
22
 
 
23
#ifndef __MPFR_IMPL_H__
 
24
#define __MPFR_IMPL_H__
 
25
 
 
26
/* Include stdio.h iff we are debugging or we want to check */
 
27
#if defined(DEBUG) || defined(WANT_ASSERT)
 
28
# include <stdio.h>
 
29
#endif
 
30
 
 
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
 
36
#endif
 
37
 
 
38
/******************************************************
 
39
 ****************** Include files *********************
 
40
 ******************************************************/
 
41
 
 
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)
 
45
#if HAVE_CONFIG_H
 
46
#include "config.h"
 
47
#endif
 
48
#endif
 
49
 
 
50
#ifdef  MPFR_HAVE_GMP_IMPL /* Build with gmp internals*/
 
51
 
 
52
# ifndef __GMP_H__
 
53
#  include "gmp.h"
 
54
# endif
 
55
# ifndef __GMP_IMPL_H__
 
56
#  include "gmp-impl.h"
 
57
# endif
 
58
# ifdef MPFR_NEED_LONGLONG_H
 
59
#  include "longlong.h"
 
60
# endif
 
61
# ifndef __MPFR_H
 
62
#  include "mpfr.h"
 
63
# endif
 
64
 
 
65
#else /* Build without gmp internals */
 
66
 
 
67
# ifndef __GMP_H__
 
68
#  include "gmp.h"
 
69
# endif
 
70
# ifndef __MPFR_H
 
71
#  include "mpfr.h"
 
72
# endif
 
73
# ifndef __GMPFR_GMP_H__
 
74
#  include "mpfr-gmp.h"
 
75
# endif
 
76
# ifdef MPFR_NEED_LONGLONG_H
 
77
#  include "mpfr-longlong.h"
 
78
# endif
 
79
 
 
80
#endif
 
81
#undef MPFR_NEED_LONGLONG_H
 
82
 
 
83
 
 
84
/******************************************************
 
85
 ***************** Detection macros *******************
 
86
 ******************************************************/
 
87
 
 
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))
 
93
#else
 
94
# define __MPFR_STDC(version) 0
 
95
#endif
 
96
 
 
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))
 
100
#else
 
101
# define __MPFR_GNUC(a, i) 0
 
102
#endif
 
103
 
 
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))
 
107
#else
 
108
# define __MPFR_GLIBC(a, i) 0
 
109
#endif
 
110
 
 
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))
 
114
#else
 
115
# define __MPFR_GMP(a, b, c) 0
 
116
#endif
 
117
 
 
118
#if defined(__ICC)
 
119
# define __MPFR_ICC(a,b,c) (__ICC >= (a)*100+(b)*10+c)
 
120
#else
 
121
# define __MPFR_ICC(a,b,c) 0
 
122
#endif
 
123
 
 
124
 
 
125
 
 
126
/******************************************************
 
127
 ******************** Check GMP ***********************
 
128
 ******************************************************/
 
129
 
 
130
#if !__MPFR_GMP(4,1,0)
 
131
# error "GMP 4.1.0 or newer needed"
 
132
#endif
 
133
 
 
134
#if GMP_NAIL_BITS != 0
 
135
# error "MPFR doesn't support nonzero values of GMP_NAIL_BITS"
 
136
#endif
 
137
 
 
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"
 
140
#endif
 
141
 
 
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
 
152
#else
 
153
# error "Can't compute log2(BITS_PER_MP_LIMB)"
 
154
#endif
 
155
 
 
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))
 
159
#else
 
160
# define MPFR_NORETURN_ATTR
 
161
# define MPFR_CONST_ATTR
 
162
#endif
 
163
 
 
164
/******************************************************
 
165
 ************* Global Internal Variables **************
 
166
 ******************************************************/
 
167
 
 
168
#ifdef MPFR_USE_THREAD_SAFE
 
169
# if __MPFR_GNUC(3,3) || __MPFR_ICC(8,1,0)
 
170
#  define MPFR_THREAD_ATTR __thread
 
171
# else
 
172
#  error "Can't build MPFR as thread safe"
 
173
# endif
 
174
#else
 
175
# define MPFR_THREAD_ATTR
 
176
#endif
 
177
 
 
178
/* Cache struct */
 
179
struct __gmpfr_cache_s {
 
180
  mpfr_t x;
 
181
  int inexact;
 
182
  int (*func)(mpfr_ptr, mpfr_rnd_t);
 
183
};
 
184
typedef struct __gmpfr_cache_s mpfr_cache_t[1];
 
185
 
 
186
#if defined (__cplusplus)
 
187
extern "C" {
 
188
#endif
 
189
 
 
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;
 
199
 
 
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;
 
203
 
 
204
 
 
205
#if defined (__cplusplus)
 
206
 }
 
207
#endif
 
208
 
 
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
 
216
 
 
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)
 
222
 
 
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))
 
245
 
 
246
 
 
247
/******************************************************
 
248
 ******************** Assertions **********************
 
249
 ******************************************************/
 
250
 
 
251
/* Compile with -DWANT_ASSERT to check all assert statements */
 
252
 
 
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. */
 
256
 
 
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) )))
 
260
 
 
261
/* MPFR_ASSERTD(expr): assertions that should be checked when testing */
 
262
#ifdef WANT_ASSERT
 
263
# define MPFR_EXP_CHECK 1
 
264
# define MPFR_ASSERTD(expr)  MPFR_ASSERTN (expr)
 
265
#else
 
266
# define MPFR_ASSERTD(expr)  ((void) 0)
 
267
#endif
 
268
 
 
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)
 
278
 
 
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;}
 
282
 
 
283
 
 
284
 
 
285
/******************************************************
 
286
 ****************** double macros *********************
 
287
 ******************************************************/
 
288
 
 
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 */
 
293
 
 
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
 
304
# else
 
305
#  define MPFR_DOUBLE_SPEC 0
 
306
# endif
 
307
#endif
 
308
 
 
309
/* Debug non IEEE floats */
 
310
#ifdef XDEBUG
 
311
# undef _GMP_IEEE_FLOATS
 
312
#endif
 
313
#ifndef _GMP_IEEE_FLOATS
 
314
# define _GMP_IEEE_FLOATS 0
 
315
#endif
 
316
 
 
317
#ifndef IEEE_DBL_MANT_DIG
 
318
#define IEEE_DBL_MANT_DIG 53
 
319
#endif
 
320
#define MPFR_LIMBS_PER_DOUBLE ((IEEE_DBL_MANT_DIG-1)/BITS_PER_MP_LIMB+1)
 
321
 
 
322
/* Visual C++ doesn't support +1.0/.00, -1.0/0.0 and 0.0/0.0
 
323
   at compile time. */
 
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)
 
329
#else
 
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)
 
333
#endif
 
334
 
 
335
/* for x of type ieee_double_extract */
 
336
#if _GMP_IEEE_FLOATS
 
337
typedef union ieee_double_extract Ieee_double_extract;
 
338
 
 
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)))
 
346
#else
 
347
# define DOUBLE_ISINF(x) ((x) > DBL_MAX || (x) < -DBL_MAX)
 
348
# if MPFR_NANISNAN
 
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)))
 
352
# else
 
353
#  define DOUBLE_ISNAN(x) ((x) != (x))
 
354
# endif
 
355
#endif
 
356
 
 
357
 
 
358
 
 
359
/******************************************************
 
360
 *************** Long double macros *******************
 
361
 ******************************************************/
 
362
 
 
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.
 
375
*/
 
376
#if LDBL_MANT_DIG && FLT_RADIX == 2
 
377
# define MPFR_LDBL_MANT_DIG LDBL_MANT_DIG
 
378
#else
 
379
# define MPFR_LDBL_MANT_DIG \
 
380
  (sizeof(long double)*BITS_PER_MP_LIMB/sizeof(mp_limb_t))
 
381
#endif
 
382
#define MPFR_LIMBS_PER_LONG_DOUBLE \
 
383
  ((sizeof(long double)-1)/sizeof(mp_limb_t)+1)
 
384
 
 
385
/* LONGDOUBLE_NAN_ACTION executes the code "action" if x is a NaN. */
 
386
 
 
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)                       \
 
393
  do {                                                          \
 
394
    union {                                                     \
 
395
      long double    ld;                                        \
 
396
      struct {                                                  \
 
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;                                 \
 
403
      } s;                                                      \
 
404
    } u;                                                        \
 
405
    u.ld = (x);                                                 \
 
406
    if (u.s.exp == 0x7FFFL                                      \
 
407
        && (u.s.man0 | u.s.man1 | u.s.man2 | u.s.man3) != 0)    \
 
408
      { action; }                                               \
 
409
  } while (0)
 
410
#endif
 
411
 
 
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
 
414
   x!=x. */
 
415
#ifndef LONGDOUBLE_NAN_ACTION
 
416
# define LONGDOUBLE_NAN_ACTION(x, action)               \
 
417
  do {                                                  \
 
418
    volatile long double __x = LONGDOUBLE_VOLATILE (x); \
 
419
    if ((x) != __x)                                     \
 
420
      { action; }                                       \
 
421
  } while (0)
 
422
# define WANT_LONGDOUBLE_VOLATILE 1
 
423
#endif
 
424
 
 
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
 
429
# ifdef 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
 
433
# else
 
434
#  define LONGDOUBLE_VOLATILE(x)  (x)
 
435
# endif
 
436
#endif
 
437
 
 
438
/* Some special case for IEEE_EXT Litle Endian */
 
439
#if HAVE_LDOUBLE_IEEE_EXT_LITTLE
 
440
 
 
441
typedef union {
 
442
  long double    ld;
 
443
  struct {
 
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;
 
449
  } s;
 
450
} mpfr_long_double_t;
 
451
 
 
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)
 
456
 
 
457
#endif
 
458
 
 
459
/******************************************************
 
460
 **************** mpfr_t properties *******************
 
461
 ******************************************************/
 
462
 
 
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)
 
467
 
 
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))
 
474
#else
 
475
# error "Invalid MPFR Prec format"
 
476
#endif
 
477
 
 
478
 
 
479
 
 
480
/******************************************************
 
481
 ***************** exponent limits ********************
 
482
 ******************************************************/
 
483
 
 
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)
 
489
#else
 
490
typedef unsigned long int  mpfr_uexp_t;
 
491
# define MPFR_EXP_MAX (LONG_MAX)
 
492
# define MPFR_EXP_MIN (LONG_MIN)
 
493
#endif
 
494
#ifndef mp_exp_unsigned_t
 
495
# define mp_exp_unsigned_t mpfr_uexp_t
 
496
#endif
 
497
 
 
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))
 
501
 
 
502
/* Definition of the intervals of the exponent limits */
 
503
#undef MPFR_EMIN_MIN
 
504
#undef MPFR_EMIN_MAX
 
505
#undef MPFR_EMAX_MIN
 
506
#undef MPFR_EMAX_MAX
 
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)
 
511
 
 
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. */
 
518
 
 
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))
 
523
#else
 
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)
 
527
#endif
 
528
 
 
529
 
 
530
 
 
531
/******************************************************
 
532
 ********** Singular Values (NAN, INF, ZERO) **********
 
533
 ******************************************************/
 
534
 
 
535
/*
 
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.
 
545
 */
 
546
 
 
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)
 
551
 
 
552
#define MPFR_CLEAR_FLAGS(x)
 
553
 
 
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)
 
561
 
 
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))
 
565
 
 
566
#define MPFR_ARE_SINGULAR(x,y) \
 
567
  (MPFR_UNLIKELY(MPFR_IS_SINGULAR(x)) || MPFR_UNLIKELY(MPFR_IS_SINGULAR(y)))
 
568
 
 
569
 
 
570
 
 
571
/******************************************************
 
572
 ********************* Sign Macros ********************
 
573
 ******************************************************/
 
574
 
 
575
#define MPFR_SIGN_POS (1)
 
576
#define MPFR_SIGN_NEG (-1)
 
577
 
 
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)
 
580
 
 
581
#define MPFR_IS_NEG(x) (MPFR_SIGN(x) < 0)
 
582
#define MPFR_IS_POS(x) (MPFR_SIGN(x) > 0)
 
583
 
 
584
#define MPFR_SET_POS(x) (MPFR_SIGN(x) = MPFR_SIGN_POS)
 
585
#define MPFR_SET_NEG(x) (MPFR_SIGN(x) = MPFR_SIGN_NEG)
 
586
 
 
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))
 
600
 
 
601
 
 
602
 
 
603
/******************************************************
 
604
 ***************** Ternary Value Macros ***************
 
605
 ******************************************************/
 
606
 
 
607
/* Special inexact value */
 
608
#define MPFR_EVEN_INEX 2
 
609
 
 
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
 
616
 
 
617
#define MPFR_SET_ERANGE() (__gmpfr_flags |= MPFR_FLAGS_ERANGE)
 
618
 
 
619
 
 
620
 
 
621
/******************************************************
 
622
 ************** Rounding mode macros  *****************
 
623
 ******************************************************/
 
624
 
 
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)
 
630
 
 
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))
 
635
 
 
636
/* Invert a rounding mode */
 
637
#define MPFR_INVERT_RND(rnd) ((rnd == GMP_RNDU) ? GMP_RNDD : \
 
638
                             ((rnd == GMP_RNDD) ? GMP_RNDU : rnd))
 
639
 
 
640
/* Transform RNDU and RNDD to RNDA or RNDZ */
 
641
#define MPFR_UPDATE_RND_MODE(rnd, test)                            \
 
642
  do {                                                             \
 
643
    if (MPFR_UNLIKELY(MPFR_IS_RNDUTEST_OR_RNDDNOTTEST(rnd, test))) \
 
644
      rnd = GMP_RNDZ;                                              \
 
645
  } while (0)
 
646
 
 
647
 
 
648
 
 
649
/******************************************************
 
650
 ******************* Limb Macros **********************
 
651
 ******************************************************/
 
652
 
 
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
 
658
#else
 
659
# error "Neither GMP_LIMB_HIGHBIT nor MP_LIMB_T_HIGHBIT defined in GMP"
 
660
#endif
 
661
 
 
662
/* Mask to get the Most Significent Bit of a limb */
 
663
#define MPFR_LIMB_MSB(l) ((l)&MPFR_LIMB_HIGHBIT)
 
664
 
 
665
/* Definition of MPFR_LIMB_ONE & MPFR_LIMB_ZERO*/
 
666
#ifdef CNST_LIMB
 
667
# define MPFR_LIMB_ONE  CNST_LIMB(1)
 
668
# define MPFR_LIMB_ZERO CNST_LIMB(0)
 
669
#else
 
670
# define MPFR_LIMB_ONE  ((mp_limb_t) 1L)
 
671
# define MPFR_LIMB_ZERO ((mp_limb_t) 0L)
 
672
#endif
 
673
 
 
674
/* Mask for the low 's' bits of a limb */
 
675
#define MPFR_LIMB_MASK(s) ((MPFR_LIMB_ONE<<(s))-MPFR_LIMB_ONE)
 
676
 
 
677
 
 
678
 
 
679
/******************************************************
 
680
 ********************** Memory ************************
 
681
 ******************************************************/
 
682
 
 
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))
 
695
 
 
696
/* Temporary memory gestion */
 
697
#ifndef TMP_SALLOC
 
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
 
703
#else
 
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
 
708
#endif
 
709
 
 
710
/* This code is experimental: don't use it */
 
711
#ifdef MPFR_USE_OWN_MPFR_TMP_ALLOC
 
712
extern unsigned char *mpfr_stack;
 
713
#undef MPFR_TMP_DECL
 
714
#undef MPFR_TMP_MARK
 
715
#undef MPFR_TMP_ALLOC
 
716
#undef MPFR_TMP_FREE
 
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))
 
721
#endif
 
722
 
 
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),                                              \
 
729
   MPFR_SET_POS(x),                                                  \
 
730
   MPFR_SET_INVALID_EXP(x))
 
731
 
 
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))
 
735
 
 
736
#define MPFR_TMP_INIT_ABS(d, s)                                      \
 
737
 ( MPFR_PREC(d) = MPFR_PREC(s),                                      \
 
738
   MPFR_MANT(d) = MPFR_MANT(s),                                      \
 
739
   MPFR_SET_POS(d),                                                  \
 
740
   MPFR_EXP(d)  = MPFR_EXP(s))
 
741
 
 
742
 
 
743
 
 
744
/******************************************************
 
745
 *****************  Cache macros **********************
 
746
 ******************************************************/
 
747
 
 
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)
 
752
 
 
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}}
 
756
 
 
757
 
 
758
 
 
759
/******************************************************
 
760
 *******************  Threshold ***********************
 
761
 ******************************************************/
 
762
 
 
763
#include "mparam.h"
 
764
 
 
765
/******************************************************
 
766
 *****************  Useful macros *********************
 
767
 ******************************************************/
 
768
 
 
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))
 
774
#else
 
775
# define MPFR_LIKELY(x) (x)
 
776
# define MPFR_UNLIKELY(x) (x)
 
777
#endif
 
778
 
 
779
/* Ceil log 2: If GCC, uses a GCC extension, otherwise calls a function */
 
780
/* Warning:
 
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); }))
 
790
#else
 
791
# define MPFR_INT_CEIL_LOG2(x) (__gmpfr_int_ceil_log2(x))
 
792
#endif
 
793
 
 
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)                  \
 
799
 do {                                                                 \
 
800
  (c) = (a) + (b);                                                    \
 
801
  if ((c) < (a)) ACTION_IF_OVERFLOW;                                  \
 
802
 } while (0)
 
803
 
 
804
#define MPFR_SADD_OVERFLOW(c,a,b,STYPE,UTYPE,MIN,MAX,ACTION_IF_POS_OVERFLOW,ACTION_IF_NEG_OVERFLOW) \
 
805
  do {                                                                \
 
806
  if ((a) >= 0 && (b) >= 0) {                                         \
 
807
         UTYPE uc,ua,ub;                                              \
 
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) {                                    \
 
813
         UTYPE uc,ua,ub;                                              \
 
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);                                             \
 
821
 } while (0)
 
822
 
 
823
 
 
824
/* Set a number to 1 (Fast) - It doesn't check if 1 is in the exponent range */
 
825
#define MPFR_SET_ONE(x)                                               \
 
826
do {                                                                  \
 
827
  mp_size_t _size = MPFR_LIMB_SIZE(x) - 1;                            \
 
828
  MPFR_SET_POS(x);                                                    \
 
829
  MPFR_EXP(x) = 1;                                                    \
 
830
  MPN_ZERO ( MPFR_MANT(x), _size);                                    \
 
831
  MPFR_MANT(x)[_size] = MPFR_LIMB_HIGHBIT;                            \
 
832
} while (0)
 
833
 
 
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)                              \
 
838
  do                                                                  \
 
839
    {                                                                 \
 
840
      if ((UINT_MAX % BITS_PER_MP_LIMB) == (BITS_PER_MP_LIMB-1))      \
 
841
        (s) = (mpfr_prec_t) (-(a)) % BITS_PER_MP_LIMB;                \
 
842
      else                                                            \
 
843
        {                                                             \
 
844
          (s) = (a) % BITS_PER_MP_LIMB;                               \
 
845
          if (s != 0)                                                 \
 
846
            (s) = BITS_PER_MP_LIMB - (s);                             \
 
847
        }                                                             \
 
848
      MPFR_ASSERTD ((s) >= 0 && (s) < BITS_PER_MP_LIMB);              \
 
849
    }                                                                 \
 
850
  while (0)
 
851
 
 
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 */
 
855
#ifdef DEBUG
 
856
# include <stdio.h>
 
857
# define MPFR_TRACE(x) x
 
858
#else
 
859
# define MPFR_TRACE(x) (void) 0
 
860
#endif
 
861
#define MPFR_DUMP(x) ( printf(#x"="), mpfr_dump(x) )
 
862
 
 
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)
 
866
 
 
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))
 
870
 
 
871
#define mpfr_get_d1(x) mpfr_get_d(x,__gmpfr_default_rounding_mode)
 
872
 
 
873
/* Store in r the size in bits of the mpz_t z */
 
874
#define MPFR_MPZ_SIZEINBASE2(r, z)              \
 
875
  do {                                          \
 
876
   int _cnt;                                    \
 
877
   mp_size_t _size;                             \
 
878
   MPFR_ASSERTD (mpz_sgn (z) != 0);             \
 
879
   _size = ABSIZ(z);                            \
 
880
   count_leading_zeros (_cnt, PTR(z)[_size-1]); \
 
881
   (r) = _size * BITS_PER_MP_LIMB - _cnt;       \
 
882
  } while (0)
 
883
 
 
884
/* Needs <locale.h> */
 
885
#define MPFR_DECIMAL_POINT ((unsigned char) localeconv()->decimal_point[0])
 
886
 
 
887
/******************************************************
 
888
 **************  Save exponent macros  ****************
 
889
 ******************************************************/
 
890
 
 
891
/* See README.dev for details on how to use the macros.
 
892
   They are used to set the exponent range to the maximum
 
893
   temporarily */
 
894
 
 
895
typedef struct {
 
896
  unsigned int saved_flags;
 
897
  mp_exp_t saved_emin;
 
898
  mp_exp_t saved_emax;
 
899
} mpfr_save_expo_t;
 
900
 
 
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)
 
914
 
 
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))
 
919
 
 
920
 
 
921
/******************************************************
 
922
 *****************  Inline Rounding *******************
 
923
 ******************************************************/
 
924
 
 
925
/*
 
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).
 
928
 */
 
929
 
 
930
/*
 
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.
 
936
 */
 
937
#define MPFR_RNDRAW_GEN(inexact, dest, srcp, sprec, rnd, sign,              \
 
938
                        MIDDLE_HANDLER, OVERFLOW_HANDLER)                   \
 
939
  do {                                                                      \
 
940
    mp_size_t dests, srcs;                                                  \
 
941
    mp_limb_t *destp;                                                       \
 
942
    mp_prec_t destprec, srcprec;                                            \
 
943
                                                                            \
 
944
    /* Check Trivial Case when Dest Mantissa has more bits than source */   \
 
945
    srcprec = sprec;                                                        \
 
946
    destprec = MPFR_PREC (dest);                                            \
 
947
    destp = MPFR_MANT (dest);                                               \
 
948
    if (MPFR_UNLIKELY (destprec >= srcprec))                                \
 
949
      {                                                                     \
 
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);                                            \
 
954
        inexact = 0;                                                        \
 
955
      }                                                                     \
 
956
    else                                                                    \
 
957
      {                                                                     \
 
958
        /* Non trivial case: rounding needed */                             \
 
959
        mp_prec_t sh;                                                       \
 
960
        mp_limb_t *sp;                                                      \
 
961
        mp_limb_t rb, sb, ulp;                                              \
 
962
                                                                            \
 
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;                                           \
 
968
                                                                            \
 
969
        /* General case when prec % BITS_PER_MP_LIMB != 0 */                \
 
970
        if (MPFR_LIKELY (sh != 0))                                          \
 
971
          {                                                                 \
 
972
            mp_limb_t mask;                                                 \
 
973
            /* Compute Rounding Bit and Sticky Bit */                       \
 
974
            mask = MPFR_LIMB_ONE << (sh-1);                                 \
 
975
            rb = sp[0] & mask;                                              \
 
976
            sb = sp[0] & (mask-1);                                          \
 
977
            if (MPFR_UNLIKELY (sb == 0))                                    \
 
978
              { /* TODO: Improve it */                                      \
 
979
                mp_limb_t *tmp;                                             \
 
980
                mp_size_t n;                                                \
 
981
                for (tmp = sp, n = srcs - dests ; n != 0 && sb == 0 ; n--)  \
 
982
                  sb = *--tmp;                                              \
 
983
              }                                                             \
 
984
            ulp = 2*mask;                                                   \
 
985
          }                                                                 \
 
986
        else /* sh == 0 */                                                  \
 
987
          {                                                                 \
 
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))                                    \
 
993
              {                                                             \
 
994
                mp_limb_t *tmp;                                             \
 
995
                mp_size_t n;                                                \
 
996
                for (tmp = sp-1, n = srcs - dests-1 ; n!=0 && sb==0 ; n--)  \
 
997
                  sb = *--tmp;                                              \
 
998
              }                                                             \
 
999
            ulp = MPFR_LIMB_ONE;                                            \
 
1000
          }                                                                 \
 
1001
        /* Rounding */                                                      \
 
1002
        if (MPFR_LIKELY (rnd == GMP_RNDN))                                  \
 
1003
          {                                                                 \
 
1004
            if (rb == 0)                                                    \
 
1005
              {                                                             \
 
1006
              trunc:                                                        \
 
1007
                inexact = MPFR_LIKELY ((sb | rb) != 0) ? -sign : 0;         \
 
1008
              trunc_doit:                                                   \
 
1009
                MPN_COPY (destp, sp, dests);                                \
 
1010
                destp[0] &= ~(ulp-1);                                       \
 
1011
              }                                                             \
 
1012
            else if (MPFR_UNLIKELY (sb == 0))                               \
 
1013
              { /* Middle of two consecutive representable numbers */       \
 
1014
                MIDDLE_HANDLER;                                             \
 
1015
              }                                                             \
 
1016
            else                                                            \
 
1017
              {                                                             \
 
1018
              addoneulp:                                                    \
 
1019
                inexact = sign;                                             \
 
1020
              addoneulp_doit:                                               \
 
1021
                if (MPFR_UNLIKELY (mpn_add_1 (destp, sp, dests, ulp)))      \
 
1022
                  {                                                         \
 
1023
                    destp[dests-1] = MPFR_LIMB_HIGHBIT;                     \
 
1024
                    OVERFLOW_HANDLER;                                       \
 
1025
                  }                                                         \
 
1026
                destp[0] &= ~(ulp-1);                                       \
 
1027
              }                                                             \
 
1028
          }                                                                 \
 
1029
        else                                                                \
 
1030
          { /* Directed rounding mode */                                    \
 
1031
            if (MPFR_LIKELY (MPFR_IS_LIKE_RNDZ (rnd,                        \
 
1032
                                                MPFR_IS_NEG_SIGN (sign))))  \
 
1033
              goto trunc;                                                   \
 
1034
             else if (MPFR_UNLIKELY ((sb | rb) == 0))                       \
 
1035
               {                                                            \
 
1036
                 inexact = 0;                                               \
 
1037
                 goto trunc_doit;                                           \
 
1038
               }                                                            \
 
1039
             else                                                           \
 
1040
              goto addoneulp;                                               \
 
1041
          }                                                                 \
 
1042
      }                                                                     \
 
1043
  } while (0)
 
1044
 
 
1045
/*
 
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.
 
1049
 */
 
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)                                                 \
 
1053
       {                                                                     \
 
1054
         inexact = -sign;                                                    \
 
1055
         goto trunc_doit;                                                    \
 
1056
       }                                                                     \
 
1057
     else                                                                    \
 
1058
       goto addoneulp;                                                       \
 
1059
     , OVERFLOW_HANDLER)
 
1060
 
 
1061
/*
 
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.
 
1066
 */
 
1067
#define MPFR_RNDRAW_EVEN(inexact, dest, srcp, sprec, rnd, sign, \
 
1068
                         OVERFLOW_HANDLER)                      \
 
1069
   MPFR_RNDRAW_GEN (inexact, dest, srcp, sprec, rnd, sign,      \
 
1070
     if ((sp[0] & ulp) == 0)                                    \
 
1071
       {                                                        \
 
1072
         inexact = -MPFR_EVEN_INEX * sign;                      \
 
1073
         goto trunc_doit;                                       \
 
1074
       }                                                        \
 
1075
     else                                                       \
 
1076
       {                                                        \
 
1077
         inexact = MPFR_EVEN_INEX * sign;                       \
 
1078
         goto addoneulp_doit;                                   \
 
1079
       }                                                        \
 
1080
     , OVERFLOW_HANDLER)
 
1081
 
 
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)))
 
1087
 
 
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).
 
1093
 
 
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.
 
1098
 
 
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.
 
1102
 
 
1103
   The test is less restrictive thant necessary, but the function
 
1104
   will finish the check itself.
 
1105
*/
 
1106
#define MPFR_FAST_COMPUTE_IF_SMALL_INPUT(y,x,err,dir,rnd,extra)     \
 
1107
  do {                                                              \
 
1108
   mp_exp_t _err = (err);                                           \
 
1109
   if (MPFR_UNLIKELY (_err > 0                                      \
 
1110
                      && (mpfr_uexp_t) _err > MPFR_PREC (y) + 1))   \
 
1111
    {                                                               \
 
1112
      int _inexact = mpfr_round_near_x ((y),(x),(err),(dir),(rnd)); \
 
1113
      if (_inexact != 0)                                            \
 
1114
        {                                                           \
 
1115
         extra;                                                     \
 
1116
         return _inexact;                                           \
 
1117
        }                                                           \
 
1118
    }                                                               \
 
1119
  } while (0)
 
1120
 
 
1121
/******************************************************
 
1122
 ***************  Ziv Loop Macro  *********************
 
1123
 ******************************************************/
 
1124
 
 
1125
#ifndef MPFR_USE_LOGGING
 
1126
 
 
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)
 
1131
 
 
1132
#else
 
1133
/* Use LOGGING */
 
1134
#define MPFR_ZIV_DECL(_x)                                     \
 
1135
  mp_prec_t _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 ); }
 
1145
 
 
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))
 
1150
 
 
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))
 
1156
 
 
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)
 
1161
 
 
1162
#endif
 
1163
 
 
1164
 
 
1165
/******************************************************
 
1166
 ***************  Logging Macros  *********************
 
1167
 ******************************************************/
 
1168
 
 
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
 
1177
 
 
1178
#ifdef MPFR_USE_LOGGING
 
1179
 
 
1180
# include <stdio.h>
 
1181
 
 
1182
/* Check if we can support this feature */
 
1183
# ifdef MPFR_USE_THREAD_SAFE
 
1184
#  error "Enable either `Logging' or `thread-safe', not both"
 
1185
# endif
 
1186
# if !__MPFR_GNUC(3,0)
 
1187
#  error "Logging not supported (GCC >= 3.0)"
 
1188
# endif
 
1189
 
 
1190
#if defined (__cplusplus)
 
1191
extern "C" {
 
1192
#endif
 
1193
 
 
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;
 
1200
 
 
1201
#if defined (__cplusplus)
 
1202
 }
 
1203
#endif
 
1204
 
 
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);
 
1208
 
 
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
 
1213
 
 
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;                                                  \
 
1222
  MPFR_LOG_BEGIN2 x
 
1223
 
 
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__;                           \
 
1233
  MPFR_LOG_END2 x
 
1234
 
 
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
 
1244
 
 
1245
#else /* MPFR_USE_LOGGING */
 
1246
 
 
1247
/* Define void macro for logging */
 
1248
 
 
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)
 
1254
 
 
1255
#endif /* MPFR_USE_LOGGING */
 
1256
 
 
1257
 
 
1258
/**************************************************************
 
1259
 ************  Group Initialize Functions Macros  *************
 
1260
 **************************************************************/
 
1261
 
 
1262
#ifndef MPFR_GROUP_STATIC_SIZE
 
1263
# define MPFR_GROUP_STATIC_SIZE 16
 
1264
#endif
 
1265
 
 
1266
struct mpfr_group_t {
 
1267
  size_t     alloc;
 
1268
  mp_limb_t *mant;
 
1269
  mp_limb_t  tab[MPFR_GROUP_STATIC_SIZE];
 
1270
};
 
1271
 
 
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);                     \
 
1277
 }} while (0)
 
1278
 
 
1279
#define MPFR_GROUP_INIT_TEMPLATE(g, prec, num, handler) do {     \
 
1280
 mp_prec_t _prec = (prec);                                       \
 
1281
 mp_size_t _size;                                                \
 
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);                 \
 
1288
 } else {                                                        \
 
1289
  (g).alloc = 0;                                                 \
 
1290
  (g).mant = (g).tab;                                            \
 
1291
 }                                                               \
 
1292
 handler;                                                        \
 
1293
 } while (0)
 
1294
#define MPFR_GROUP_TINIT(g, n, x) MPFR_TMP_INIT1 ((g).mant+_size*(n), x, _prec)
 
1295
 
 
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))
 
1319
 
 
1320
#define MPFR_GROUP_REPREC_TEMPLATE(g, prec, num, handler) do {   \
 
1321
 mp_prec_t _prec = (prec);                                       \
 
1322
 size_t    _oalloc = (g).alloc;                                  \
 
1323
 mp_size_t _size;                                                \
 
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);                 \
 
1330
 else                                                            \
 
1331
  (g).mant = (*__gmp_reallocate_func)((g).mant,_oalloc,(g).alloc);\
 
1332
 handler;                                                        \
 
1333
} while (0)
 
1334
 
 
1335
 
 
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))
 
1359
 
 
1360
 
 
1361
/******************************************************
 
1362
 ***************  Internal Functions  *****************
 
1363
 ******************************************************/
 
1364
 
 
1365
#if defined (__cplusplus)
 
1366
extern "C" {
 
1367
#endif
 
1368
 
 
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));
 
1371
 
 
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));
 
1382
 
 
1383
__MPFR_DECLSPEC int mpfr_cmp2 _MPFR_PROTO ((mpfr_srcptr, mpfr_srcptr,
 
1384
                                            mp_prec_t *));
 
1385
 
 
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));
 
1392
 
 
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));
 
1396
 
 
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));
 
1399
 
 
1400
__MPFR_DECLSPEC long mpfr_mpn_exp _MPFR_PROTO ((mp_limb_t *, mp_exp_t *, int,
 
1401
                           mp_exp_t, size_t));
 
1402
 
 
1403
#ifdef _MPFR_H_HAVE_FILE
 
1404
__MPFR_DECLSPEC void mpfr_fprint_binary _MPFR_PROTO ((FILE *, mpfr_srcptr));
 
1405
#endif
 
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*));
 
1410
 
 
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));
 
1419
 
 
1420
#define mpfr_round_raw2(xp, xn, neg, r, prec) \
 
1421
  mpfr_round_raw_2((xp),(xn)*BITS_PER_MP_LIMB,(neg),(prec),(r))
 
1422
 
 
1423
__MPFR_DECLSPEC int mpfr_check _MPFR_PROTO ((mpfr_srcptr));
 
1424
 
 
1425
__MPFR_DECLSPEC int mpfr_sum_sort _MPFR_PROTO ((mpfr_srcptr *const,
 
1426
                                                unsigned long, mpfr_srcptr *));
 
1427
 
 
1428
__MPFR_DECLSPEC int mpfr_get_cputime _MPFR_PROTO ((void));
 
1429
 
 
1430
__MPFR_DECLSPEC void mpfr_nexttozero _MPFR_PROTO ((mpfr_ptr));
 
1431
__MPFR_DECLSPEC void mpfr_nexttoinf _MPFR_PROTO ((mpfr_ptr));
 
1432
 
 
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));
 
1437
 
 
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,
 
1442
                                              mpfr_rnd_t));
 
1443
 
 
1444
__MPFR_DECLSPEC void mpfr_mulhigh_n _MPFR_PROTO ((mp_ptr, mp_srcptr,
 
1445
                                                  mp_srcptr, mp_size_t));
 
1446
 
 
1447
__MPFR_DECLSPEC int mpfr_round_p _MPFR_PROTO ((mp_limb_t *, mp_size_t,
 
1448
                                               mp_exp_t, mp_prec_t));
 
1449
 
 
1450
__MPFR_DECLSPEC void mpfr_dump_mant _MPFR_PROTO ((const mp_limb_t *,
 
1451
                                                  mp_prec_t, mp_prec_t,
 
1452
                                                  mp_prec_t));
 
1453
 
 
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))
 
1457
       MPFR_NORETURN_ATTR;
 
1458
 
 
1459
#if defined (__cplusplus)
 
1460
}
 
1461
#endif
 
1462
 
 
1463
#endif