?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
kernel - Fix hold race and catch bad vm_object terminations
[dragonfly.git] / contrib / mpfr / mpfr-impl.h
bloba4202d4cc95141bf25e3220b77a57fc0ab6def29
1 /* Utilities for MPFR developers, not exported.
2
3 Copyright 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 Free Software Foundation, Inc.
4 Contributed by the Arenaire and Cacao projects, INRIA.
5
6 This file is part of the GNU MPFR Library.
7
8 The GNU 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 GNU 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 GNU MPFR Library; see the file COPYING.LIB. If not, write to
20 the Free Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston,
21 MA 02110-1301, USA. */
22
23 #ifndef __MPFR_IMPL_H__
24 #define __MPFR_IMPL_H__
25
26 /* Let's include some standard headers unconditionally as they are
27 already needed by several source files or when some options are
28 enabled/disabled, and it is easy to forget them (some configure
29 options may hide the error).
30 Note: If some source file must not have such a header included
31 (which is very unlikely and probably means something broken in
32 this source file), we should do that with some macro (that would
33 also force to disable incompatible features). */
34 #if defined (__cplusplus)
35 #include <cstdio>
36 #include <cstring>
37 #else
38 #include <stdio.h>
39 #include <string.h>
40 #endif
41 #include <limits.h>
42
43 /* Check if we are inside a build of MPFR or inside the test suite.
44 This is needed in mpfr.h to export or import the functions.
45 It matters only for Windows DLL */
46 #ifndef __MPFR_TEST_H__
47 # define __MPFR_WITHIN_MPFR 1
48 #endif
49
50 /******************************************************
51 ****************** Include files *********************
52 ******************************************************/
53
54 /* Include 'config.h' before using ANY configure macros if needed
55 NOTE: It isn't MPFR 'config.h', but GMP's one! */
56 #ifdef HAVE_CONFIG_H
57 # include "config.h"
58 #endif
59
60 #ifdef MPFR_HAVE_GMP_IMPL /* Build with gmp internals*/
61
62 # ifndef __GMP_H__
63 # include "gmp.h"
64 # endif
65 # ifndef __GMP_IMPL_H__
66 # include "gmp-impl.h"
67 # endif
68 # ifndef BITS_PER_MP_LIMB
69 # define BITS_PER_MP_LIMB GMP_LIMB_BITS
70 # endif
71 #ifndef mpn_sqr_n
72 # define mpn_sqr_n mpn_sqr
73 #endif
74 # ifdef MPFR_NEED_LONGLONG_H
75 # include "longlong.h"
76 # endif
77 # ifndef __MPFR_H
78 # include "mpfr.h"
79 # endif
80
81 #else /* Build without gmp internals */
82
83 # ifndef __GMP_H__
84 # include "gmp.h"
85 # endif
86 # ifndef __MPFR_H
87 # include "mpfr.h"
88 # endif
89 # ifndef __GMPFR_GMP_H__
90 # include "mpfr-gmp.h"
91 # endif
92 # ifdef MPFR_NEED_LONGLONG_H
93 # include "mpfr-longlong.h"
94 # endif
95
96 #endif
97 #undef MPFR_NEED_LONGLONG_H
98
99 /* For the definition of MPFR_THREAD_ATTR. GCC/ICC detection macros are
100 no longer used, as they sometimes gave incorrect information about
101 the support of thread-local variables. A configure check is now done.
102 If the use of detection macros is needed in the future, this could be
103 moved below (after the detection macros are defined). */
104 #include "mpfr-thread.h"
105
106
107 /******************************************************
108 ***************** Detection macros *******************
109 ******************************************************/
110
111 /* Macros to detect STDC, GCC, GLIBC, GMP and ICC version */
112 #if defined(__STDC_VERSION__)
113 # define __MPFR_STDC(version) (__STDC_VERSION__>=(version))
114 #elif defined(__STDC__)
115 # define __MPFR_STDC(version) (0 == (version))
116 #else
117 # define __MPFR_STDC(version) 0
118 #endif
119
120 #if defined(__ICC)
121 # define __MPFR_ICC(a,b,c) (__ICC >= (a)*100+(b)*10+(c))
122 #elif defined(__INTEL_COMPILER)
123 # define __MPFR_ICC(a,b,c) (__INTEL_COMPILER >= (a)*100+(b)*10+(c))
124 #else
125 # define __MPFR_ICC(a,b,c) 0
126 #endif
127
128 #if defined(__GNUC__) && defined(__GNUC_MINOR__) && ! __MPFR_ICC(0,0,0)
129 # define __MPFR_GNUC(a,i) \
130 (MPFR_VERSION_NUM(__GNUC__,__GNUC_MINOR__,0) >= MPFR_VERSION_NUM(a,i,0))
131 #else
132 # define __MPFR_GNUC(a,i) 0
133 #endif
134
135 #if defined(__GLIBC__) && defined(__GLIBC_MINOR__)
136 # define __MPFR_GLIBC(a,i) \
137 (MPFR_VERSION_NUM(__GLIBC__,__GLIBC_MINOR__,0) >= MPFR_VERSION_NUM(a,i,0))
138 #else
139 # define __MPFR_GLIBC(a,i) 0
140 #endif
141
142 #if defined(__GNU_MP_VERSION) && \
143 defined(__GNU_MP_VERSION_MINOR) && \
144 defined(__GNU_MP_VERSION_PATCHLEVEL)
145 # define __MPFR_GMP(a,b,c) \
146 (MPFR_VERSION_NUM(__GNU_MP_VERSION,__GNU_MP_VERSION_MINOR,__GNU_MP_VERSION_PATCHLEVEL) >= MPFR_VERSION_NUM(a,b,c))
147 #else
148 # define __MPFR_GMP(a,b,c) 0
149 #endif
150
151
152
153 /******************************************************
154 ******************** Check GMP ***********************
155 ******************************************************/
156
157 #if !__MPFR_GMP(4,1,0)
158 # error "GMP 4.1.0 or newer needed"
159 #endif
160
161 #if GMP_NAIL_BITS != 0
162 # error "MPFR doesn't support nonzero values of GMP_NAIL_BITS"
163 #endif
164
165 #if (BITS_PER_MP_LIMB<32) || (BITS_PER_MP_LIMB & (BITS_PER_MP_LIMB - 1))
166 # error "BITS_PER_MP_LIMB must be a power of 2, and >= 32"
167 #endif
168
169 #if BITS_PER_MP_LIMB == 16
170 # define MPFR_LOG2_BITS_PER_MP_LIMB 4
171 #elif BITS_PER_MP_LIMB == 32
172 # define MPFR_LOG2_BITS_PER_MP_LIMB 5
173 #elif BITS_PER_MP_LIMB == 64
174 # define MPFR_LOG2_BITS_PER_MP_LIMB 6
175 #elif BITS_PER_MP_LIMB == 128
176 # define MPFR_LOG2_BITS_PER_MP_LIMB 7
177 #elif BITS_PER_MP_LIMB == 256
178 # define MPFR_LOG2_BITS_PER_MP_LIMB 8
179 #else
180 # error "Can't compute log2(BITS_PER_MP_LIMB)"
181 #endif
182
183 #if __MPFR_GNUC(3,0) || __MPFR_ICC(8,1,0)
184 # define MPFR_NORETURN_ATTR __attribute__ ((noreturn))
185 # define MPFR_CONST_ATTR __attribute__ ((const))
186 #else
187 # define MPFR_NORETURN_ATTR
188 # define MPFR_CONST_ATTR
189 #endif
190
191 /******************************************************
192 ************* Global Internal Variables **************
193 ******************************************************/
194
195 /* Cache struct */
196 struct __gmpfr_cache_s {
197 mpfr_t x;
198 int inexact;
199 int (*func)(mpfr_ptr, mpfr_rnd_t);
200 };
201 typedef struct __gmpfr_cache_s mpfr_cache_t[1];
202
203 #if defined (__cplusplus)
204 extern "C" {
205 #endif
206
207 __MPFR_DECLSPEC extern MPFR_THREAD_ATTR unsigned int __gmpfr_flags;
208 __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mp_exp_t __gmpfr_emin;
209 __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mp_exp_t __gmpfr_emax;
210 __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mp_prec_t __gmpfr_default_fp_bit_precision;
211 __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_rnd_t __gmpfr_default_rounding_mode;
212 __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_pi;
213 __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_log2;
214 __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_euler;
215 __MPFR_DECLSPEC extern MPFR_THREAD_ATTR mpfr_cache_t __gmpfr_cache_const_catalan;
216
217 #define BASE_MAX 36
218 __MPFR_DECLSPEC extern const __mpfr_struct __gmpfr_l2b[BASE_MAX-1][2];
219
220 /* Note: do not use the following values when they can be outside the
221 current exponent range, e.g. when the exponent range has not been
222 extended yet; under such a condition, they can be used only in
223 mpfr_cmpabs. */
224 __MPFR_DECLSPEC extern const mpfr_t __gmpfr_one;
225 __MPFR_DECLSPEC extern const mpfr_t __gmpfr_two;
226 __MPFR_DECLSPEC extern const mpfr_t __gmpfr_four;
227
228
229 #if defined (__cplusplus)
230 }
231 #endif
232
233 /* Flags of __gmpfr_flags */
234 #define MPFR_FLAGS_UNDERFLOW 1
235 #define MPFR_FLAGS_OVERFLOW 2
236 #define MPFR_FLAGS_NAN 4
237 #define MPFR_FLAGS_INEXACT 8
238 #define MPFR_FLAGS_ERANGE 16
239 #define MPFR_FLAGS_ALL 31
240
241 /* Replace some commun functions for direct access to the global vars */
242 #define mpfr_get_emin() (__gmpfr_emin + 0)
243 #define mpfr_get_emax() (__gmpfr_emax + 0)
244 #define mpfr_get_default_rounding_mode() (__gmpfr_default_rounding_mode + 0)
245 #define mpfr_get_default_prec() (__gmpfr_default_fp_bit_precision + 0)
246
247 #define mpfr_clear_flags() \
248 ((void) (__gmpfr_flags = 0))
249 #define mpfr_clear_underflow() \
250 ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_UNDERFLOW))
251 #define mpfr_clear_overflow() \
252 ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_OVERFLOW))
253 #define mpfr_clear_nanflag() \
254 ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_NAN))
255 #define mpfr_clear_inexflag() \
256 ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_INEXACT))
257 #define mpfr_clear_erangeflag() \
258 ((void) (__gmpfr_flags &= MPFR_FLAGS_ALL ^ MPFR_FLAGS_ERANGE))
259 #define mpfr_underflow_p() \
260 ((int) (__gmpfr_flags & MPFR_FLAGS_UNDERFLOW))
261 #define mpfr_overflow_p() \
262 ((int) (__gmpfr_flags & MPFR_FLAGS_OVERFLOW))
263 #define mpfr_nanflag_p() \
264 ((int) (__gmpfr_flags & MPFR_FLAGS_NAN))
265 #define mpfr_inexflag_p() \
266 ((int) (__gmpfr_flags & MPFR_FLAGS_INEXACT))
267 #define mpfr_erangeflag_p() \
268 ((int) (__gmpfr_flags & MPFR_FLAGS_ERANGE))
269
270 /* Testing an exception flag correctly is tricky. There are mainly two
271 pitfalls: First, one needs to remember to clear the corresponding
272 flag, in case it was set before the function call or during some
273 intermediate computations (in practice, one can clear all the flags).
274 Secondly, one needs to test the flag early enough, i.e. before it
275 can be modified by another function. Moreover, it is quite difficult
276 (if not impossible) to reliably check problems with "make check". To
277 avoid these pitfalls, it is recommended to use the following macros.
278 Other use of the exception-flag predicate functions/macros will be
279 detected by mpfrlint.
280 Note: _op can be either a statement or an expression.
281 MPFR_BLOCK_EXCEP should be used only inside a block; it is useful to
282 detect some exception in order to exit the block as soon as possible. */
283 #define MPFR_BLOCK_DECL(_flags) unsigned int _flags
284 #define MPFR_BLOCK(_flags,_op) \
285 do \
286 { \
287 mpfr_clear_flags (); \
288 _op; \
289 (_flags) = __gmpfr_flags; \
290 } \
291 while (0)
292 #define MPFR_BLOCK_TEST(_flags,_f) MPFR_UNLIKELY ((_flags) & (_f))
293 #define MPFR_BLOCK_EXCEP (__gmpfr_flags & (MPFR_FLAGS_UNDERFLOW | \
294 MPFR_FLAGS_OVERFLOW | \
295 MPFR_FLAGS_NAN))
296 /* Let's use a MPFR_ prefix, because e.g. OVERFLOW is defined by glibc's
297 math.h, though this is not a reserved identifier! */
298 #define MPFR_UNDERFLOW(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_UNDERFLOW)
299 #define MPFR_OVERFLOW(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_OVERFLOW)
300 #define MPFR_NANFLAG(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_NAN)
301 #define MPFR_INEXFLAG(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_INEXACT)
302 #define MPFR_ERANGEFLAG(_flags) MPFR_BLOCK_TEST (_flags, MPFR_FLAGS_ERANGE)
303
304
305 /******************************************************
306 ******************** Assertions **********************
307 ******************************************************/
308
309 /* Compile with -DWANT_ASSERT to check all assert statements */
310
311 /* Note: do not use GMP macros ASSERT_ALWAYS and ASSERT as they are not
312 expressions, and as a consequence, they cannot be used in a for(),
313 with a comma operator and so on. */
314
315 /* MPFR_ASSERTN(expr): assertions that should always be checked */
316 #define MPFR_ASSERTN(expr) \
317 ((void) ((MPFR_UNLIKELY(expr)) || MPFR_UNLIKELY( (ASSERT_FAIL(expr),0) )))
318
319 /* MPFR_ASSERTD(expr): assertions that should be checked when testing */
320 #ifdef WANT_ASSERT
321 # define MPFR_EXP_CHECK 1
322 # define MPFR_ASSERTD(expr) MPFR_ASSERTN (expr)
323 #else
324 # define MPFR_ASSERTD(expr) ((void) 0)
325 #endif
326
327 /* Code to deal with impossible
328 WARNING: It doesn't use do { } while (0) for Insure++*/
329 #define MPFR_RET_NEVER_GO_HERE() {MPFR_ASSERTN(0); return 0;}
330
331
332 /******************************************************
333 ******************** Warnings ************************
334 ******************************************************/
335
336 /* MPFR_WARNING is no longer useful, but let's keep the macro in case
337 it needs to be used again in the future. */
338
339 #ifdef MPFR_USE_WARNINGS
340 # include <stdlib.h>
341 # define MPFR_WARNING(W) \
342 do \
343 { \
344 char *q = getenv ("MPFR_QUIET"); \
345 if (q == NULL || *q == 0) \
346 fprintf (stderr, "MPFR: %s\n", W); \
347 } \
348 while (0)
349 #else
350 # define MPFR_WARNING(W) ((void) 0)
351 #endif
352
353
354 /******************************************************
355 ****************** double macros *********************
356 ******************************************************/
357
358 /* Definition of constants */
359 #define LOG2 0.69314718055994528622 /* log(2) rounded to zero on 53 bits */
360 #define ALPHA 4.3191365662914471407 /* a+2 = a*log(a), rounded to +infinity */
361 #define EXPM1 0.36787944117144227851 /* exp(-1), rounded to zero */
362
363 /* MPFR_DOUBLE_SPEC = 1 if the C type 'double' corresponds to IEEE-754
364 double precision, 0 if it doesn't, and undefined if one doesn't know.
365 On all the tested machines, MPFR_DOUBLE_SPEC = 1. To have this macro
366 defined here, #include <float.h> is needed. If need be, other values
367 could be defined for other specs (once they are known). */
368 #if !defined(MPFR_DOUBLE_SPEC) && defined(FLT_RADIX) && \
369 defined(DBL_MANT_DIG) && defined(DBL_MIN_EXP) && defined(DBL_MAX_EXP)
370 # if FLT_RADIX == 2 && DBL_MANT_DIG == 53 && \
371 DBL_MIN_EXP == -1021 && DBL_MAX_EXP == 1024
372 # define MPFR_DOUBLE_SPEC 1
373 # else
374 # define MPFR_DOUBLE_SPEC 0
375 # endif
376 #endif
377
378 /* Debug non IEEE floats */
379 #ifdef XDEBUG
380 # undef _GMP_IEEE_FLOATS
381 #endif
382 #ifndef _GMP_IEEE_FLOATS
383 # define _GMP_IEEE_FLOATS 0
384 #endif
385
386 #ifndef IEEE_DBL_MANT_DIG
387 #define IEEE_DBL_MANT_DIG 53
388 #endif
389 #define MPFR_LIMBS_PER_DOUBLE ((IEEE_DBL_MANT_DIG-1)/BITS_PER_MP_LIMB+1)
390
391 /* Visual C++ doesn't support +1.0/.00, -1.0/0.0 and 0.0/0.0
392 at compile time. */
393 #if defined(_MSC_VER) && defined(_WIN32) && (_MSC_VER >= 1200)
394 static double double_zero = 0.0;
395 # define DBL_NAN (double_zero/double_zero)
396 # define DBL_POS_INF ((double) 1.0/double_zero)
397 # define DBL_NEG_INF ((double)-1.0/double_zero)
398 # define DBL_NEG_ZERO (-double_zero)
399 #else
400 # define DBL_POS_INF ((double) 1.0/0.0)
401 # define DBL_NEG_INF ((double)-1.0/0.0)
402 # define DBL_NAN ((double) 0.0/0.0)
403 # define DBL_NEG_ZERO (-0.0)
404 #endif
405
406 /* Note: In the past, there was specific code for _GMP_IEEE_FLOATS, which
407 was based on NaN and Inf memory representations. This code was breaking
408 the aliasing rules (see ISO C99, 6.5#6 and 6.5#7 on the effective type)
409 and for this reason it did not behave correctly with GCC 4.5.0 20091119.
410 The code needed a memory transfer and was probably not better than the
411 macros below with a good compiler (a fix based on the NaN / Inf memory
412 representation would be even worse due to C limitations), and this code
413 could be selected only when MPFR was built with --with-gmp-build, thus
414 introducing a difference (bad for maintaining/testing MPFR); therefore
415 it has been removed. The old code required that the argument x be an
416 lvalue of type double. We still require that, in case one would need
417 to change the macros below, e.g. for some broken compiler. But the
418 LVALUE(x) condition could be removed if really necessary. */
419 /* Below, the &(x) == &(x) || &(x) != &(x) allows to make sure that x
420 is a lvalue without (probably) any warning from the compiler. The
421 &(x) != &(x) is needed to avoid a failure under Mac OS X 10.4.11
422 (with Xcode 2.4.1, i.e. the latest one). */
423 #define LVALUE(x) (&(x) == &(x) || &(x) != &(x))
424 #define DOUBLE_ISINF(x) (LVALUE(x) && ((x) > DBL_MAX || (x) < -DBL_MAX))
425 #ifdef MPFR_NANISNAN
426 /* Avoid MIPSpro / IRIX64 / gcc -ffast-math (incorrect) optimizations.
427 The + must not be replaced by a ||. With gcc -ffast-math, NaN is
428 regarded as a positive number or something like that; the second
429 test catches this case. */
430 # define DOUBLE_ISNAN(x) \
431 (LVALUE(x) && !((((x) >= 0.0) + ((x) <= 0.0)) && -(x)*(x) <= 0.0))
432 #else
433 # define DOUBLE_ISNAN(x) (LVALUE(x) && (x) != (x))
434 #endif
435
436 /******************************************************
437 *************** Long double macros *******************
438 ******************************************************/
439
440 /* We try to get the exact value of the precision of long double
441 (provided by the implementation) in order to provide correct
442 rounding in this case (not guaranteed if the C implementation
443 does not have an adequate long double arithmetic). Note that
444 it may be lower than the precision of some numbers that can
445 be represented in a long double; e.g. on FreeBSD/x86, it is
446 53 because the processor is configured to round in double
447 precision (even when using the long double type -- this is a
448 limitation of the x87 arithmetic), and on Mac OS X, it is 106
449 because the implementation is a double-double arithmetic.
450 Otherwise (e.g. in base 10), we get an upper bound of the
451 precision, and correct rounding isn't currently provided.
452 */
453 #if defined(LDBL_MANT_DIG) && FLT_RADIX == 2
454 # define MPFR_LDBL_MANT_DIG LDBL_MANT_DIG
455 #else
456 # define MPFR_LDBL_MANT_DIG \
457 (sizeof(long double)*BITS_PER_MP_LIMB/sizeof(mp_limb_t))
458 #endif
459 #define MPFR_LIMBS_PER_LONG_DOUBLE \
460 ((sizeof(long double)-1)/sizeof(mp_limb_t)+1)
461
462 /* LONGDOUBLE_NAN_ACTION executes the code "action" if x is a NaN. */
463
464 /* On hppa2.0n-hp-hpux10 with the unbundled HP cc, the test x!=x on a NaN
465 has been seen false, meaning NaNs are not detected. This seemed to
466 happen only after other comparisons, not sure what's really going on. In
467 any case we can pick apart the bytes to identify a NaN. */
468 #ifdef HAVE_LDOUBLE_IEEE_QUAD_BIG
469 # define LONGDOUBLE_NAN_ACTION(x, action) \
470 do { \
471 union { \
472 long double ld; \
473 struct { \
474 unsigned int sign : 1; \
475 unsigned int exp : 15; \
476 unsigned int man3 : 16; \
477 unsigned int man2 : 32; \
478 unsigned int man1 : 32; \
479 unsigned int man0 : 32; \
480 } s; \
481 } u; \
482 u.ld = (x); \
483 if (u.s.exp == 0x7FFFL \
484 && (u.s.man0 | u.s.man1 | u.s.man2 | u.s.man3) != 0) \
485 { action; } \
486 } while (0)
487 #endif
488
489 /* Under IEEE rules, NaN is not equal to anything, including itself.
490 "volatile" here stops "cc" on mips64-sgi-irix6.5 from optimizing away
491 x!=x. */
492 #ifndef LONGDOUBLE_NAN_ACTION
493 # define LONGDOUBLE_NAN_ACTION(x, action) \
494 do { \
495 volatile long double __x = LONGDOUBLE_VOLATILE (x); \
496 if ((x) != __x) \
497 { action; } \
498 } while (0)
499 # define WANT_LONGDOUBLE_VOLATILE 1
500 #endif
501
502 /* If we don't have a proper "volatile" then volatile is #defined to empty,
503 in this case call through an external function to stop the compiler
504 optimizing anything. */
505 #ifdef WANT_LONGDOUBLE_VOLATILE
506 # ifdef volatile
507 __MPFR_DECLSPEC long double __gmpfr_longdouble_volatile _MPFR_PROTO ((long double)) MPFR_CONST_ATTR;
508 # define LONGDOUBLE_VOLATILE(x) (__gmpfr_longdouble_volatile (x))
509 # define WANT_GMPFR_LONGDOUBLE_VOLATILE 1
510 # else
511 # define LONGDOUBLE_VOLATILE(x) (x)
512 # endif
513 #endif
514
515 /* Some special case for IEEE_EXT Litle Endian */
516 #if HAVE_LDOUBLE_IEEE_EXT_LITTLE
517
518 typedef union {
519 long double ld;
520 struct {
521 unsigned int manl : 32;
522 unsigned int manh : 32;
523 unsigned int expl : 8 ;
524 unsigned int exph : 7;
525 unsigned int sign : 1;
526 } s;
527 } mpfr_long_double_t;
528
529 /* #undef MPFR_LDBL_MANT_DIG */
530 #undef MPFR_LIMBS_PER_LONG_DOUBLE
531 /* #define MPFR_LDBL_MANT_DIG 64 */
532 #define MPFR_LIMBS_PER_LONG_DOUBLE ((64-1)/BITS_PER_MP_LIMB+1)
533
534 #endif
535
536 /******************************************************
537 *************** _Decimal64 support *******************
538 ******************************************************/
539
540 #ifdef MPFR_WANT_DECIMAL_FLOATS
541 /* to cast between binary64 and decimal64 */
542 union ieee_double_decimal64 { double d; _Decimal64 d64; };
543 #endif
544
545 /******************************************************
546 **************** mpfr_t properties *******************
547 ******************************************************/
548
549 #define MPFR_PREC(x) ((x)->_mpfr_prec)
550 #define MPFR_EXP(x) ((x)->_mpfr_exp)
551 #define MPFR_MANT(x) ((x)->_mpfr_d)
552 #define MPFR_LIMB_SIZE(x) ((MPFR_PREC((x))-1)/BITS_PER_MP_LIMB+1)
553
554 #if _MPFR_PREC_FORMAT == 1
555 # define MPFR_INTPREC_MAX (USHRT_MAX & ~(unsigned int) (BITS_PER_MP_LIMB - 1))
556 #elif _MPFR_PREC_FORMAT == 2
557 # define MPFR_INTPREC_MAX (UINT_MAX & ~(unsigned int) (BITS_PER_MP_LIMB - 1))
558 #elif _MPFR_PREC_FORMAT == 3
559 # define MPFR_INTPREC_MAX (ULONG_MAX & ~(unsigned long) (BITS_PER_MP_LIMB - 1))
560 #else
561 # error "Invalid MPFR Prec format"
562 #endif
563
564
565
566 /******************************************************
567 ***************** exponent limits ********************
568 ******************************************************/
569
570 /* Define limits and unsigned type of exponent. The following definitions
571 * depend on mp_exp_t; if this type changes in GMP, these definitions will
572 * need to be modified (alternatively, a mpfr_exp_t type could be defined).
573 */
574 #if __GMP_MP_SIZE_T_INT == 1
575 typedef unsigned int mpfr_uexp_t;
576 # define MPFR_EXP_MAX (INT_MAX)
577 # define MPFR_EXP_MIN (INT_MIN)
578 #else
579 typedef unsigned long int mpfr_uexp_t;
580 # define MPFR_EXP_MAX (LONG_MAX)
581 # define MPFR_EXP_MIN (LONG_MIN)
582 #endif
583 #ifndef mp_exp_unsigned_t
584 # define mp_exp_unsigned_t mpfr_uexp_t
585 #endif
586
587 #if MPFR_EXP_MIN >= LONG_MIN && MPFR_EXP_MAX <= LONG_MAX
588 typedef long int mpfr_eexp_t;
589 # define mpfr_get_exp_t(x,r) mpfr_get_si((x),(r))
590 # define mpfr_set_exp_t(x,e,r) mpfr_set_si((x),(e),(r))
591 #elif defined (_MPFR_H_HAVE_INTMAX_T)
592 typedef intmax_t mpfr_eexp_t;
593 # define mpfr_get_exp_t(x,r) mpfr_get_sj((x),(r))
594 # define mpfr_set_exp_t(x,e,r) mpfr_set_sj((x),(e),(r))
595 #else
596 # error "Cannot define mpfr_get_exp_t and mpfr_set_exp_t"
597 #endif
598
599 /* Invalid exponent value (to track bugs...) */
600 #define MPFR_EXP_INVALID \
601 ((mp_exp_t) 1 << (BITS_PER_MP_LIMB*sizeof(mp_exp_t)/sizeof(mp_limb_t)-2))
602
603 /* Definition of the exponent limits for MPFR numbers.
604 * These limits are chosen so that if e is such an exponent, then 2e-1 and
605 * 2e+1 are representable. This is useful for intermediate computations,
606 * in particular the multiplication.
607 */
608 #undef MPFR_EMIN_MIN
609 #undef MPFR_EMIN_MAX
610 #undef MPFR_EMAX_MIN
611 #undef MPFR_EMAX_MAX
612 #define MPFR_EMIN_MIN (1-MPFR_EXP_INVALID)
613 #define MPFR_EMIN_MAX (MPFR_EXP_INVALID-1)
614 #define MPFR_EMAX_MIN (1-MPFR_EXP_INVALID)
615 #define MPFR_EMAX_MAX (MPFR_EXP_INVALID-1)
616
617 /* Use MPFR_GET_EXP and MPFR_SET_EXP instead of MPFR_EXP directly,
618 unless when the exponent may be out-of-range, for instance when
619 setting the exponent before calling mpfr_check_range.
620 MPFR_EXP_CHECK is defined when WANT_ASSERT is defined, but if you
621 don't use WANT_ASSERT (for speed reasons), you can still define
622 MPFR_EXP_CHECK by setting -DMPFR_EXP_CHECK in $CFLAGS. */
623
624 #ifdef MPFR_EXP_CHECK
625 # define MPFR_GET_EXP(x) (mpfr_get_exp) (x)
626 # define MPFR_SET_EXP(x, exp) MPFR_ASSERTN (!mpfr_set_exp ((x), (exp)))
627 # define MPFR_SET_INVALID_EXP(x) ((void) (MPFR_EXP (x) = MPFR_EXP_INVALID))
628 #else
629 # define MPFR_GET_EXP(x) MPFR_EXP (x)
630 # define MPFR_SET_EXP(x, exp) ((void) (MPFR_EXP (x) = (exp)))
631 # define MPFR_SET_INVALID_EXP(x) ((void) 0)
632 #endif
633
634
635
636 /******************************************************
637 ********** Singular Values (NAN, INF, ZERO) **********
638 ******************************************************/
639
640 /*
641 * Clear flags macros are still defined and should be still used
642 * since the functions must not assume the internal format.
643 * How to deal with special values ?
644 * 1. Check if is a special value (Zero, Nan, Inf) wiht MPFR_IS_SINGULAR
645 * 2. Deal with the special value with MPFR_IS_NAN, MPFR_IS_INF, etc
646 * 3. Else clear the flags of the dest (it must be done after since src
647 * may be also the dest!)
648 * MPFR_SET_INF, MPFR_SET_NAN, MPFR_SET_ZERO must clear by
649 * themselves the other flags.
650 */
651
652 /* Enum special value of exponent.*/
653 # define MPFR_EXP_ZERO (MPFR_EXP_MIN+1)
654 # define MPFR_EXP_NAN (MPFR_EXP_MIN+2)
655 # define MPFR_EXP_INF (MPFR_EXP_MIN+3)
656
657 #define MPFR_CLEAR_FLAGS(x)
658
659 #define MPFR_IS_NAN(x) (MPFR_EXP(x) == MPFR_EXP_NAN)
660 #define MPFR_SET_NAN(x) (MPFR_EXP(x) = MPFR_EXP_NAN)
661 #define MPFR_IS_INF(x) (MPFR_EXP(x) == MPFR_EXP_INF)
662 #define MPFR_SET_INF(x) (MPFR_EXP(x) = MPFR_EXP_INF)
663 #define MPFR_IS_ZERO(x) (MPFR_EXP(x) == MPFR_EXP_ZERO)
664 #define MPFR_SET_ZERO(x) (MPFR_EXP(x) = MPFR_EXP_ZERO)
665 #define MPFR_NOTZERO(x) (MPFR_EXP(x) != MPFR_EXP_ZERO)
666
667 #define MPFR_IS_FP(x) (!MPFR_IS_NAN(x) && !MPFR_IS_INF(x))
668 #define MPFR_IS_SINGULAR(x) (MPFR_EXP(x) <= MPFR_EXP_INF)
669 #define MPFR_IS_PURE_FP(x) (!MPFR_IS_SINGULAR(x))
670
671 #define MPFR_ARE_SINGULAR(x,y) \
672 (MPFR_UNLIKELY(MPFR_IS_SINGULAR(x)) || MPFR_UNLIKELY(MPFR_IS_SINGULAR(y)))
673
674
675
676 /******************************************************
677 ********************* Sign Macros ********************
678 ******************************************************/
679
680 #define MPFR_SIGN_POS (1)
681 #define MPFR_SIGN_NEG (-1)
682
683 #define MPFR_IS_STRICTPOS(x) (MPFR_NOTZERO((x)) && MPFR_SIGN(x) > 0)
684 #define MPFR_IS_STRICTNEG(x) (MPFR_NOTZERO((x)) && MPFR_SIGN(x) < 0)
685
686 #define MPFR_IS_NEG(x) (MPFR_SIGN(x) < 0)
687 #define MPFR_IS_POS(x) (MPFR_SIGN(x) > 0)
688
689 #define MPFR_SET_POS(x) (MPFR_SIGN(x) = MPFR_SIGN_POS)
690 #define MPFR_SET_NEG(x) (MPFR_SIGN(x) = MPFR_SIGN_NEG)
691
692 #define MPFR_CHANGE_SIGN(x) (MPFR_SIGN(x) = -MPFR_SIGN(x))
693 #define MPFR_SET_SAME_SIGN(x, y) (MPFR_SIGN(x) = MPFR_SIGN(y))
694 #define MPFR_SET_OPPOSITE_SIGN(x, y) (MPFR_SIGN(x) = -MPFR_SIGN(y))
695 #define MPFR_ASSERT_SIGN(s) \
696 (MPFR_ASSERTD((s) == MPFR_SIGN_POS || (s) == MPFR_SIGN_NEG))
697 #define MPFR_SET_SIGN(x, s) \
698 (MPFR_ASSERT_SIGN(s), MPFR_SIGN(x) = s)
699 #define MPFR_IS_POS_SIGN(s1) (s1 > 0)
700 #define MPFR_IS_NEG_SIGN(s1) (s1 < 0)
701 #define MPFR_MULT_SIGN(s1, s2) ((s1) * (s2))
702 /* Transform a sign to 1 or -1 */
703 #define MPFR_FROM_SIGN_TO_INT(s) (s)
704 #define MPFR_INT_SIGN(x) MPFR_FROM_SIGN_TO_INT(MPFR_SIGN(x))
705
706
707
708 /******************************************************
709 ***************** Ternary Value Macros ***************
710 ******************************************************/
711
712 /* Special inexact value */
713 #define MPFR_EVEN_INEX 2
714
715 /* When returning the ternary inexact value, ALWAYS use one of the
716 following two macros, unless the flag comes from another function
717 returning the ternary inexact value */
718 #define MPFR_RET(I) return \
719 (I) ? ((__gmpfr_flags |= MPFR_FLAGS_INEXACT), (I)) : 0
720 #define MPFR_RET_NAN return (__gmpfr_flags |= MPFR_FLAGS_NAN), 0
721
722 #define MPFR_SET_ERANGE() (__gmpfr_flags |= MPFR_FLAGS_ERANGE)
723
724 #define SIGN(I) ((I) < 0 ? -1 : (I) > 0)
725 #define SAME_SIGN(I1,I2) (SIGN (I1) == SIGN (I2))
726
727
728
729 /******************************************************
730 ************** Rounding mode macros *****************
731 ******************************************************/
732
733 /* We want to test this :
734 * (rnd == GMP_RNDU && test) || (rnd == RNDD && !test)
735 * ie it transforms RNDU or RNDD to Away or Zero according to the sign */
736 #define MPFR_IS_RNDUTEST_OR_RNDDNOTTEST(rnd, test) \
737 (((rnd) + (test)) == GMP_RNDD)
738
739 /* We want to test if rnd = Zero, or Away.
740 'test' is true iff negative. */
741 #define MPFR_IS_LIKE_RNDZ(rnd, test) \
742 ((rnd==GMP_RNDZ) || MPFR_IS_RNDUTEST_OR_RNDDNOTTEST (rnd, test))
743
744 /* Invert a rounding mode */
745 #define MPFR_INVERT_RND(rnd) ((rnd == GMP_RNDU) ? GMP_RNDD : \
746 ((rnd == GMP_RNDD) ? GMP_RNDU : rnd))
747
748 /* Transform RNDU and RNDD to RNDA or RNDZ */
749 #define MPFR_UPDATE_RND_MODE(rnd, test) \
750 do { \
751 if (MPFR_UNLIKELY(MPFR_IS_RNDUTEST_OR_RNDDNOTTEST(rnd, test))) \
752 rnd = GMP_RNDZ; \
753 } while (0)
754
755
756
757 /******************************************************
758 ******************* Limb Macros **********************
759 ******************************************************/
760
761 /* Definition of MPFR_LIMB_HIGHBIT */
762 #if defined(GMP_LIMB_HIGHBIT)
763 # define MPFR_LIMB_HIGHBIT GMP_LIMB_HIGHBIT
764 #elif defined(MP_LIMB_T_HIGHBIT)
765 # define MPFR_LIMB_HIGHBIT MP_LIMB_T_HIGHBIT
766 #else
767 # error "Neither GMP_LIMB_HIGHBIT nor MP_LIMB_T_HIGHBIT defined in GMP"
768 #endif
769
770 /* Mask to get the Most Significant Bit of a limb */
771 #define MPFR_LIMB_MSB(l) ((l)&MPFR_LIMB_HIGHBIT)
772
773 /* Definition of MPFR_LIMB_ONE & MPFR_LIMB_ZERO*/
774 #ifdef CNST_LIMB
775 # define MPFR_LIMB_ONE CNST_LIMB(1)
776 # define MPFR_LIMB_ZERO CNST_LIMB(0)
777 #else
778 # define MPFR_LIMB_ONE ((mp_limb_t) 1L)
779 # define MPFR_LIMB_ZERO ((mp_limb_t) 0L)
780 #endif
781
782 /* Mask for the low 's' bits of a limb */
783 #define MPFR_LIMB_MASK(s) ((MPFR_LIMB_ONE<<(s))-MPFR_LIMB_ONE)
784
785
786
787 /******************************************************
788 ********************** Memory ************************
789 ******************************************************/
790
791 /* Heap Memory gestion */
792 typedef union { mp_size_t s; mp_limb_t l; } mpfr_size_limb_t;
793 #define MPFR_GET_ALLOC_SIZE(x) \
794 ( ((mp_size_t*) MPFR_MANT(x))[-1] + 0)
795 #define MPFR_SET_ALLOC_SIZE(x, n) \
796 ( ((mp_size_t*) MPFR_MANT(x))[-1] = n)
797 #define MPFR_MALLOC_SIZE(s) \
798 ( sizeof(mpfr_size_limb_t) + BYTES_PER_MP_LIMB * ((size_t) s) )
799 #define MPFR_SET_MANT_PTR(x,p) \
800 (MPFR_MANT(x) = (mp_limb_t*) ((mpfr_size_limb_t*) p + 1))
801 #define MPFR_GET_REAL_PTR(x) \
802 ((mp_limb_t*) ((mpfr_size_limb_t*) MPFR_MANT(x) - 1))
803
804 /* Temporary memory gestion */
805 #ifndef TMP_SALLOC
806 /* GMP 4.1.x or below or internals */
807 #define MPFR_TMP_DECL TMP_DECL
808 #define MPFR_TMP_MARK TMP_MARK
809 #define MPFR_TMP_ALLOC TMP_ALLOC
810 #define MPFR_TMP_FREE TMP_FREE
811 #else
812 #define MPFR_TMP_DECL(x) TMP_DECL
813 #define MPFR_TMP_MARK(x) TMP_MARK
814 #define MPFR_TMP_ALLOC(s) TMP_ALLOC(s)
815 #define MPFR_TMP_FREE(x) TMP_FREE
816 #endif
817
818 /* This code is experimental: don't use it */
819 #ifdef MPFR_USE_OWN_MPFR_TMP_ALLOC
820 extern unsigned char *mpfr_stack;
821 #undef MPFR_TMP_DECL
822 #undef MPFR_TMP_MARK
823 #undef MPFR_TMP_ALLOC
824 #undef MPFR_TMP_FREE
825 #define MPFR_TMP_DECL(_x) unsigned char *(_x)
826 #define MPFR_TMP_MARK(_x) ((_x) = mpfr_stack)
827 #define MPFR_TMP_ALLOC(_s) (mpfr_stack += (_s), mpfr_stack - (_s))
828 #define MPFR_TMP_FREE(_x) (mpfr_stack = (_x))
829 #endif
830
831 /* temporary allocate 1 limb at xp, and initialize mpfr variable x */
832 /* The temporary var doesn't have any size field, but it doesn't matter
833 * since only functions dealing with the Heap care about it */
834 #define MPFR_TMP_INIT1(xp, x, p) \
835 ( MPFR_PREC(x) = (p), \
836 MPFR_MANT(x) = (xp), \
837 MPFR_SET_POS(x), \
838 MPFR_SET_INVALID_EXP(x))
839
840 #define MPFR_TMP_INIT(xp, x, p, s) \
841 (xp = (mp_ptr) MPFR_TMP_ALLOC(BYTES_PER_MP_LIMB * ((size_t) s)), \
842 MPFR_TMP_INIT1(xp, x, p))
843
844 #define MPFR_TMP_INIT_ABS(d, s) \
845 ( MPFR_PREC(d) = MPFR_PREC(s), \
846 MPFR_MANT(d) = MPFR_MANT(s), \
847 MPFR_SET_POS(d), \
848 MPFR_EXP(d) = MPFR_EXP(s))
849
850
851
852 /******************************************************
853 ***************** Cache macros **********************
854 ******************************************************/
855
856 #define mpfr_const_pi(_d,_r) mpfr_cache(_d, __gmpfr_cache_const_pi,_r)
857 #define mpfr_const_log2(_d,_r) mpfr_cache(_d, __gmpfr_cache_const_log2, _r)
858 #define mpfr_const_euler(_d,_r) mpfr_cache(_d, __gmpfr_cache_const_euler, _r)
859 #define mpfr_const_catalan(_d,_r) mpfr_cache(_d,__gmpfr_cache_const_catalan,_r)
860
861 #define MPFR_DECL_INIT_CACHE(_cache,_func) \
862 mpfr_cache_t MPFR_THREAD_ATTR _cache = \
863 {{{{0,MPFR_SIGN_POS,0,(mp_limb_t*)0}},0,_func}}
864
865
866
867 /******************************************************
868 ******************* Threshold ***********************
869 ******************************************************/
870
871 #include "mparam.h"
872
873 /******************************************************
874 ***************** Useful macros *********************
875 ******************************************************/
876
877 /* Theses macros help the compiler to determine if a test is
878 * likely or unlikely. */
879 #if __MPFR_GNUC(3,0) || __MPFR_ICC(8,1,0)
880 # define MPFR_LIKELY(x) (__builtin_expect(!!(x),1))
881 # define MPFR_UNLIKELY(x) (__builtin_expect((x),0))
882 #else
883 # define MPFR_LIKELY(x) (x)
884 # define MPFR_UNLIKELY(x) (x)
885 #endif
886
887 /* Declare that some variable is initialized before being used (without a
888 dummy initialization) in order to avoid some compiler warnings. Use the
889 VAR = VAR trick (see http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36296)
890 only with gcc as this is undefined behavior, and we don't know what
891 other compilers do (they may also be smarter). This trick could be
892 disabled with future gcc versions. */
893 #if defined(__GNUC__)
894 # define INITIALIZED(VAR) VAR = VAR
895 #else
896 # define INITIALIZED(VAR) VAR
897 #endif
898
899 /* Ceil log 2: If GCC, uses a GCC extension, otherwise calls a function */
900 /* Warning:
901 * Needs to define MPFR_NEED_LONGLONG.
902 * Computes ceil(log2(x)) only for x integer (unsigned long)
903 * Undefined if x is 0 */
904 #if __MPFR_GNUC(2,95) || __MPFR_ICC(8,1,0)
905 # define MPFR_INT_CEIL_LOG2(x) \
906 (MPFR_UNLIKELY ((x) == 1) ? 0 : \
907 __extension__ ({ int _b; mp_limb_t _limb; \
908 MPFR_ASSERTN ((x) > 1); \
909 _limb = (x) - 1; \
910 MPFR_ASSERTN (_limb == (x) - 1); \
911 count_leading_zeros (_b, _limb); \
912 (BITS_PER_MP_LIMB - _b); }))
913 #else
914 # define MPFR_INT_CEIL_LOG2(x) (__gmpfr_int_ceil_log2(x))
915 #endif
916
917 /* Add two integers with overflow handling */
918 /* Example: MPFR_SADD_OVERFLOW (c, a, b, long, unsigned long,
919 * LONG_MIN, LONG_MAX,
920 * goto overflow, goto underflow); */
921 #define MPFR_UADD_OVERFLOW(c,a,b,ACTION_IF_OVERFLOW) \
922 do { \
923 (c) = (a) + (b); \
924 if ((c) < (a)) ACTION_IF_OVERFLOW; \
925 } while (0)
926
927 #define MPFR_SADD_OVERFLOW(c,a,b,STYPE,UTYPE,MIN,MAX,ACTION_IF_POS_OVERFLOW,ACTION_IF_NEG_OVERFLOW) \
928 do { \
929 if ((a) >= 0 && (b) >= 0) { \
930 UTYPE uc,ua,ub; \
931 ua = (UTYPE) a; ub = (UTYPE) b; \
932 MPFR_UADD_OVERFLOW (uc, ua, ub, ACTION_IF_POS_OVERFLOW); \
933 if (uc > (UTYPE)(MAX)) ACTION_IF_POS_OVERFLOW; \
934 else (c) = (STYPE) uc; \
935 } else if ((a) < 0 && (b) < 0) { \
936 UTYPE uc,ua,ub; \
937 ua = -(UTYPE) a; ub = -(UTYPE) b; \
938 MPFR_UADD_OVERFLOW (uc, ua, ub, ACTION_IF_NEG_OVERFLOW); \
939 if (uc >= -(UTYPE)(MIN) || uc > (UTYPE)(MAX)) { \
940 if (uc == -(UTYPE)(MIN)) (c) = (MIN); \
941 else ACTION_IF_NEG_OVERFLOW; } \
942 else (c) = -(STYPE) uc; \
943 } else (c) = (a) + (b); \
944 } while (0)
945
946
947 /* Set a number to 1 (Fast) - It doesn't check if 1 is in the exponent range */
948 #define MPFR_SET_ONE(x) \
949 do { \
950 mp_size_t _size = MPFR_LIMB_SIZE(x) - 1; \
951 MPFR_SET_POS(x); \
952 MPFR_EXP(x) = 1; \
953 MPN_ZERO ( MPFR_MANT(x), _size); \
954 MPFR_MANT(x)[_size] = MPFR_LIMB_HIGHBIT; \
955 } while (0)
956
957 /* Compute s = (-a) % BITS_PER_MP_LIMB
958 * a is unsigned! Check if it works,
959 * otherwise tries another way to compute it */
960 #define MPFR_UNSIGNED_MINUS_MODULO(s, a) \
961 do \
962 { \
963 if (IS_POW2 (BITS_PER_MP_LIMB)) \
964 (s) = (-(a)) % BITS_PER_MP_LIMB; \
965 else \
966 { \
967 (s) = (a) % BITS_PER_MP_LIMB; \
968 if ((s) != 0) \
969 (s) = BITS_PER_MP_LIMB - (s); \
970 } \
971 MPFR_ASSERTD ((s) >= 0 && (s) < BITS_PER_MP_LIMB); \
972 } \
973 while (0)
974
975 /* Use it only for debug reasons */
976 /* MPFR_TRACE (operation) : execute operation iff DEBUG flag is set */
977 /* MPFR_DUMP (x) : print x (a mpfr_t) on stdout */
978 #ifdef DEBUG
979 # define MPFR_TRACE(x) x
980 #else
981 # define MPFR_TRACE(x) (void) 0
982 #endif
983 #define MPFR_DUMP(x) ( printf(#x"="), mpfr_dump(x) )
984
985 /* Test if X (positive) is a power of 2 */
986 #define IS_POW2(X) (((X) & ((X) - 1)) == 0)
987 #define NOT_POW2(X) (((X) & ((X) - 1)) != 0)
988
989 /* Safe absolute value (to avoid possible integer overflow) */
990 /* type is the target (unsigned) type */
991 #define SAFE_ABS(type,x) ((x) >= 0 ? (type)(x) : -(type)(x))
992
993 #define mpfr_get_d1(x) mpfr_get_d(x,__gmpfr_default_rounding_mode)
994
995 /* Store in r the size in bits of the mpz_t z */
996 #define MPFR_MPZ_SIZEINBASE2(r, z) \
997 do { \
998 int _cnt; \
999 mp_size_t _size; \
1000 MPFR_ASSERTD (mpz_sgn (z) != 0); \
1001 _size = ABSIZ(z); \
1002 count_leading_zeros (_cnt, PTR(z)[_size-1]); \
1003 (r) = _size * BITS_PER_MP_LIMB - _cnt; \
1004 } while (0)
1005
1006 /* Needs <locale.h> */
1007 #ifdef HAVE_LOCALE_H
1008 #include <locale.h>
1009 /* Warning! In case of signed char, the value of MPFR_DECIMAL_POINT may
1010 be negative (the ISO C99 does not seem to forbid negative values). */
1011 #define MPFR_DECIMAL_POINT (localeconv()->decimal_point[0])
1012 #define MPFR_THOUSANDS_SEPARATOR (localeconv()->thousands_sep[0])
1013 #else
1014 #define MPFR_DECIMAL_POINT ((char) '.')
1015 #define MPFR_THOUSANDS_SEPARATOR ('\0')
1016 #endif
1017
1018
1019 /* Set y to s*significand(x)*2^e, for example MPFR_ALIAS(y,x,1,MPFR_EXP(x))
1020 sets y to |x|, and MPFR_ALIAS(y,x,MPFR_SIGN(x),0) sets y to x*2^f such
1021 that 1/2 <= |y| < 1. Does not check y is in the valid exponent range.
1022 WARNING! x and y share the same mantissa. So, some operations are
1023 not valid if x has been provided via an argument, e.g., trying to
1024 modify the mantissa of y, even temporarily, or calling mpfr_clear on y.
1025 */
1026 #define MPFR_ALIAS(y,x,s,e) \
1027 do \
1028 { \
1029 MPFR_PREC(y) = MPFR_PREC(x); \
1030 MPFR_SIGN(y) = (s); \
1031 MPFR_EXP(y) = (e); \
1032 MPFR_MANT(y) = MPFR_MANT(x); \
1033 } while (0)
1034
1035
1036 /******************************************************
1037 ************** Save exponent macros ****************
1038 ******************************************************/
1039
1040 /* See README.dev for details on how to use the macros.
1041 They are used to set the exponent range to the maximum
1042 temporarily */
1043
1044 typedef struct {
1045 unsigned int saved_flags;
1046 mp_exp_t saved_emin;
1047 mp_exp_t saved_emax;
1048 } mpfr_save_expo_t;
1049
1050 #define MPFR_SAVE_EXPO_DECL(x) mpfr_save_expo_t x
1051 #define MPFR_SAVE_EXPO_MARK(x) \
1052 ((x).saved_flags = __gmpfr_flags, \
1053 (x).saved_emin = __gmpfr_emin, \
1054 (x).saved_emax = __gmpfr_emax, \
1055 __gmpfr_emin = MPFR_EMIN_MIN, \
1056 __gmpfr_emax = MPFR_EMAX_MAX)
1057 #define MPFR_SAVE_EXPO_FREE(x) \
1058 (__gmpfr_flags = (x).saved_flags, \
1059 __gmpfr_emin = (x).saved_emin, \
1060 __gmpfr_emax = (x).saved_emax)
1061 #define MPFR_SAVE_EXPO_UPDATE_FLAGS(x, flags) \
1062 (x).saved_flags |= (flags)
1063
1064 /* Speed up final checking */
1065 #define mpfr_check_range(x,t,r) \
1066 (MPFR_LIKELY (MPFR_EXP (x) >= __gmpfr_emin && MPFR_EXP (x) <= __gmpfr_emax) \
1067 ? ((t) ? (__gmpfr_flags |= MPFR_FLAGS_INEXACT, (t)) : 0) \
1068 : mpfr_check_range(x,t,r))
1069
1070
1071 /******************************************************
1072 ***************** Inline Rounding *******************
1073 ******************************************************/
1074
1075 /*
1076 * Note: due to the labels, one cannot use a macro MPFR_RNDRAW* more than
1077 * once in a function (otherwise these labels would not be unique).
1078 */
1079
1080 /*
1081 * Round mantissa (srcp, sprec) to mpfr_t dest using rounding mode rnd
1082 * assuming dest's sign is sign.
1083 * In rounding to nearest mode, execute MIDDLE_HANDLER when the value
1084 * is the middle of two consecutive numbers in dest precision.
1085 * Execute OVERFLOW_HANDLER in case of overflow when rounding.
1086 */
1087 #define MPFR_RNDRAW_GEN(inexact, dest, srcp, sprec, rnd, sign, \
1088 MIDDLE_HANDLER, OVERFLOW_HANDLER) \
1089 do { \
1090 mp_size_t _dests, _srcs; \
1091 mp_limb_t *_destp; \
1092 mp_prec_t _destprec, _srcprec; \
1093 \
1094 /* Check Trivial Case when Dest Mantissa has more bits than source */ \
1095 _srcprec = sprec; \
1096 _destprec = MPFR_PREC (dest); \
1097 _destp = MPFR_MANT (dest); \
1098 if (MPFR_UNLIKELY (_destprec >= _srcprec)) \
1099 { \
1100 _srcs = (_srcprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB; \
1101 _dests = (_destprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB - _srcs; \
1102 MPN_COPY (_destp + _dests, srcp, _srcs); \
1103 MPN_ZERO (_destp, _dests); \
1104 inexact = 0; \
1105 } \
1106 else \
1107 { \
1108 /* Non trivial case: rounding needed */ \
1109 mp_prec_t _sh; \
1110 mp_limb_t *_sp; \
1111 mp_limb_t _rb, _sb, _ulp; \
1112 \
1113 /* Compute Position and shift */ \
1114 _srcs = (_srcprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB; \
1115 _dests = (_destprec + BITS_PER_MP_LIMB-1)/BITS_PER_MP_LIMB; \
1116 MPFR_UNSIGNED_MINUS_MODULO (_sh, _destprec); \
1117 _sp = srcp + _srcs - _dests; \
1118 \
1119 /* General case when prec % BITS_PER_MP_LIMB != 0 */ \
1120 if (MPFR_LIKELY (_sh != 0)) \
1121 { \
1122 mp_limb_t _mask; \
1123 /* Compute Rounding Bit and Sticky Bit */ \
1124 _mask = MPFR_LIMB_ONE << (_sh - 1); \
1125 _rb = _sp[0] & _mask; \
1126 _sb = _sp[0] & (_mask - 1); \
1127 if (MPFR_UNLIKELY (_sb == 0)) \
1128 { /* TODO: Improve it */ \
1129 mp_limb_t *_tmp; \
1130 mp_size_t _n; \
1131 for (_tmp = _sp, _n = _srcs - _dests ; \
1132 _n != 0 && _sb == 0 ; _n--) \
1133 _sb = *--_tmp; \
1134 } \
1135 _ulp = 2 * _mask; \
1136 } \
1137 else /* _sh == 0 */ \
1138 { \
1139 MPFR_ASSERTD (_dests < _srcs); \
1140 /* Compute Rounding Bit and Sticky Bit */ \
1141 _rb = _sp[-1] & MPFR_LIMB_HIGHBIT; \
1142 _sb = _sp[-1] & (MPFR_LIMB_HIGHBIT-1); \
1143 if (MPFR_UNLIKELY (_sb == 0)) \
1144 { \
1145 mp_limb_t *_tmp; \
1146 mp_size_t _n; \
1147 for (_tmp = _sp - 1, _n = _srcs - _dests - 1 ; \
1148 _n != 0 && _sb == 0 ; _n--) \
1149 _sb = *--_tmp; \
1150 } \
1151 _ulp = MPFR_LIMB_ONE; \
1152 } \
1153 /* Rounding */ \
1154 if (MPFR_LIKELY (rnd == GMP_RNDN)) \
1155 { \
1156 if (_rb == 0) \
1157 { \
1158 trunc: \
1159 inexact = MPFR_LIKELY ((_sb | _rb) != 0) ? -sign : 0; \
1160 trunc_doit: \
1161 MPN_COPY (_destp, _sp, _dests); \
1162 _destp[0] &= ~(_ulp - 1); \
1163 } \
1164 else if (MPFR_UNLIKELY (_sb == 0)) \
1165 { /* Middle of two consecutive representable numbers */ \
1166 MIDDLE_HANDLER; \
1167 } \
1168 else \
1169 { \
1170 if (0) \
1171 goto addoneulp_doit; /* dummy code to avoid warning */ \
1172 addoneulp: \
1173 inexact = sign; \
1174 addoneulp_doit: \
1175 if (MPFR_UNLIKELY (mpn_add_1 (_destp, _sp, _dests, _ulp))) \
1176 { \
1177 _destp[_dests - 1] = MPFR_LIMB_HIGHBIT; \
1178 OVERFLOW_HANDLER; \
1179 } \
1180 _destp[0] &= ~(_ulp - 1); \
1181 } \
1182 } \
1183 else \
1184 { /* Directed rounding mode */ \
1185 if (MPFR_LIKELY (MPFR_IS_LIKE_RNDZ (rnd, \
1186 MPFR_IS_NEG_SIGN (sign)))) \
1187 goto trunc; \
1188 else if (MPFR_UNLIKELY ((_sb | _rb) == 0)) \
1189 { \
1190 inexact = 0; \
1191 goto trunc_doit; \
1192 } \
1193 else \
1194 goto addoneulp; \
1195 } \
1196 } \
1197 } while (0)
1198
1199 /*
1200 * Round mantissa (srcp, sprec) to mpfr_t dest using rounding mode rnd
1201 * assuming dest's sign is sign.
1202 * Execute OVERFLOW_HANDLER in case of overflow when rounding.
1203 */
1204 #define MPFR_RNDRAW(inexact, dest, srcp, sprec, rnd, sign, OVERFLOW_HANDLER) \
1205 MPFR_RNDRAW_GEN (inexact, dest, srcp, sprec, rnd, sign, \
1206 if ((_sp[0] & _ulp) == 0) \
1207 { \
1208 inexact = -sign; \
1209 goto trunc_doit; \
1210 } \
1211 else \
1212 goto addoneulp; \
1213 , OVERFLOW_HANDLER)
1214
1215 /*
1216 * Round mantissa (srcp, sprec) to mpfr_t dest using rounding mode rnd
1217 * assuming dest's sign is sign.
1218 * Execute OVERFLOW_HANDLER in case of overflow when rounding.
1219 * Set inexact to +/- MPFR_EVEN_INEX in case of even rounding.
1220 */
1221 #define MPFR_RNDRAW_EVEN(inexact, dest, srcp, sprec, rnd, sign, \
1222 OVERFLOW_HANDLER) \
1223 MPFR_RNDRAW_GEN (inexact, dest, srcp, sprec, rnd, sign, \
1224 if ((_sp[0] & _ulp) == 0) \
1225 { \
1226 inexact = -MPFR_EVEN_INEX * sign; \
1227 goto trunc_doit; \
1228 } \
1229 else \
1230 { \
1231 inexact = MPFR_EVEN_INEX * sign; \
1232 goto addoneulp_doit; \
1233 } \
1234 , OVERFLOW_HANDLER)
1235
1236 /* Return TRUE if b is non singular and we can round it to precision 'prec'
1237 and determine the ternary value, with rounding mode 'rnd', and with
1238 error at most 'error' */
1239 #define MPFR_CAN_ROUND(b,err,prec,rnd) \
1240 (!MPFR_IS_SINGULAR (b) && mpfr_round_p (MPFR_MANT (b), MPFR_LIMB_SIZE (b), \
1241 (err), (prec) + ((rnd)==GMP_RNDN)))
1242
1243 /* TODO: fix this description (see round_near_x.c). */
1244 /* Assuming that the function has a Taylor expansion which looks like:
1245 y=o(f(x)) = o(v + g(x)) with |g(x)| <= 2^(EXP(v)-err)
1246 we can quickly set y to v if x is small (ie err > prec(y)+1) in most
1247 cases. It assumes that f(x) is not representable exactly as a FP number.
1248 v must not be a singular value (NAN, INF or ZERO); usual values are
1249 v=1 or v=x.
1250
1251 y is the destination (a mpfr_t), v the value to set (a mpfr_t),
1252 err1+err2 with err2 <= 3 the error term (mp_exp_t's), dir (an int) is
1253 the direction of the committed error (if dir = 0, it rounds toward 0,
1254 if dir=1, it rounds away from 0), rnd the rounding mode.
1255
1256 It returns from the function a ternary value in case of success.
1257 If you want to free something, you must fill the "extra" field
1258 in consequences, otherwise put nothing in it.
1259
1260 The test is less restrictive than necessary, but the function
1261 will finish the check itself.
1262
1263 Note: err1 + err2 is allowed to overflow as mp_exp_t, but it must give
1264 its real value as mpfr_uexp_t.
1265 */
1266 #define MPFR_FAST_COMPUTE_IF_SMALL_INPUT(y,v,err1,err2,dir,rnd,extra) \
1267 do { \
1268 mpfr_ptr _y = (y); \
1269 mp_exp_t _err1 = (err1); \
1270 mp_exp_t _err2 = (err2); \
1271 if (_err1 > 0) \
1272 { \
1273 mpfr_uexp_t _err = (mpfr_uexp_t) _err1 + _err2; \
1274 if (MPFR_UNLIKELY (_err > MPFR_PREC (_y) + 1)) \
1275 { \
1276 int _inexact = mpfr_round_near_x (_y,(v),_err,(dir),(rnd)); \
1277 if (_inexact != 0) \
1278 { \
1279 extra; \
1280 return _inexact; \
1281 } \
1282 } \
1283 } \
1284 } while (0)
1285
1286 /* Variant, to be called somewhere after MPFR_SAVE_EXPO_MARK. This variant
1287 is needed when there are some computations before or when some non-zero
1288 real constant is used, such as __gmpfr_one for mpfr_cos. */
1289 #define MPFR_SMALL_INPUT_AFTER_SAVE_EXPO(y,v,err1,err2,dir,rnd,expo,extra) \
1290 do { \
1291 mpfr_ptr _y = (y); \
1292 mp_exp_t _err1 = (err1); \
1293 mp_exp_t _err2 = (err2); \
1294 if (_err1 > 0) \
1295 { \
1296 mpfr_uexp_t _err = (mpfr_uexp_t) _err1 + _err2; \
1297 if (MPFR_UNLIKELY (_err > MPFR_PREC (_y) + 1)) \
1298 { \
1299 int _inexact; \
1300 mpfr_clear_flags (); \
1301 _inexact = mpfr_round_near_x (_y,(v),_err,(dir),(rnd)); \
1302 if (_inexact != 0) \
1303 { \
1304 extra; \
1305 MPFR_SAVE_EXPO_UPDATE_FLAGS (expo, __gmpfr_flags); \
1306 MPFR_SAVE_EXPO_FREE (expo); \
1307 return mpfr_check_range (_y, _inexact, (rnd)); \
1308 } \
1309 } \
1310 } \
1311 } while (0)
1312
1313 /******************************************************
1314 *************** Ziv Loop Macro *********************
1315 ******************************************************/
1316
1317 #ifndef MPFR_USE_LOGGING
1318
1319 #define MPFR_ZIV_DECL(_x) mp_prec_t _x
1320 #define MPFR_ZIV_INIT(_x, _p) (_x) = BITS_PER_MP_LIMB
1321 #define MPFR_ZIV_NEXT(_x, _p) ((_p) += (_x), (_x) = (_p)/2)
1322 #define MPFR_ZIV_FREE(x)
1323
1324 #else
1325
1326 /* The following test on glibc is there mainly for Darwin (Mac OS X), to
1327 obtain a better error message. The real test should have been a test
1328 concerning nested functions in gcc, which are disabled by default on
1329 Darwin; but it is not possible to do that without a configure test. */
1330 # if defined (__cplusplus) || !(__MPFR_GNUC(3,0) && __MPFR_GLIBC(2,0))
1331 # error "Logging not supported (needs gcc >= 3.0 and GNU C Library >= 2.0)."
1332 # endif
1333
1334 /* Use LOGGING */
1335 #define MPFR_ZIV_DECL(_x) \
1336 mp_prec_t _x; \
1337 int _x ## _cpt = 1; \
1338 static unsigned long _x ## _loop = 0, _x ## _bad = 0; \
1339 static const char *_x ## _fname = __func__; \
1340 auto void __attribute__ ((destructor)) x ## _f (void); \
1341 void __attribute__ ((destructor)) x ## _f (void) { \
1342 if (_x ## _loop != 0 && MPFR_LOG_STAT_F&mpfr_log_type) \
1343 fprintf (mpfr_log_file, \
1344 "%s: Ziv failed %2.2f%% (%lu bad cases / %lu calls)\n", _x ## _fname, \
1345 (double) 100.0 * _x ## _bad / _x ## _loop, _x ## _bad, _x ## _loop ); }
1346
1347 #define MPFR_ZIV_INIT(_x, _p) ((_x) = BITS_PER_MP_LIMB, _x ## _loop ++); \
1348 if (MPFR_LOG_BADCASE_F&mpfr_log_type && mpfr_log_current<=mpfr_log_level) \
1349 fprintf (mpfr_log_file, "%s:ZIV 1st prec=%lu\n", __func__, \
1350 (unsigned long) (_p))
1351
1352 #define MPFR_ZIV_NEXT(_x, _p) \
1353 ((_p)+=(_x),(_x)=(_p)/2, _x ## _bad += (_x ## _cpt == 1), _x ## _cpt ++); \
1354 if (MPFR_LOG_BADCASE_F&mpfr_log_type && mpfr_log_current<=mpfr_log_level) \
1355 fprintf (mpfr_log_file, "%s:ZIV new prec=%lu\n", __func__, \
1356 (unsigned long) (_p))
1357
1358 #define MPFR_ZIV_FREE(_x) \
1359 if (MPFR_LOG_BADCASE_F&mpfr_log_type && _x##_cpt>1 \
1360 && mpfr_log_current<=mpfr_log_level) \
1361 fprintf (mpfr_log_file, "%s:ZIV %d loops\n", __func__, _x ## _cpt)
1362
1363 #endif
1364
1365
1366 /******************************************************
1367 *************** Logging Macros *********************
1368 ******************************************************/
1369
1370 /* The different kind of LOG */
1371 #define MPFR_LOG_INPUT_F 1
1372 #define MPFR_LOG_OUTPUT_F 2
1373 #define MPFR_LOG_INTERNAL_F 4
1374 #define MPFR_LOG_TIME_F 8
1375 #define MPFR_LOG_BADCASE_F 16
1376 #define MPFR_LOG_MSG_F 32
1377 #define MPFR_LOG_STAT_F 64
1378
1379 #ifdef MPFR_USE_LOGGING
1380
1381 /* Check if we can support this feature */
1382 # ifdef MPFR_USE_THREAD_SAFE
1383 # error "Enable either `Logging' or `thread-safe', not both"
1384 # endif
1385 # if !__MPFR_GNUC(3,0)
1386 # error "Logging not supported (GCC >= 3.0)"
1387 # endif
1388
1389 #if defined (__cplusplus)
1390 extern "C" {
1391 #endif
1392
1393 __MPFR_DECLSPEC extern FILE *mpfr_log_file;
1394 __MPFR_DECLSPEC extern int mpfr_log_type;
1395 __MPFR_DECLSPEC extern int mpfr_log_level;
1396 __MPFR_DECLSPEC extern int mpfr_log_current;
1397 __MPFR_DECLSPEC extern int mpfr_log_base;
1398 __MPFR_DECLSPEC extern mp_prec_t mpfr_log_prec;
1399
1400 #if defined (__cplusplus)
1401 }
1402 #endif
1403
1404 #define MPFR_LOG_VAR(x) \
1405 if((MPFR_LOG_INTERNAL_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level))\
1406 fprintf (mpfr_log_file, "%s.%d:%s[%#R]=%R\n", __func__,__LINE__, #x, x, x);
1407
1408 #define MPFR_LOG_MSG2(format, ...) \
1409 if ((MPFR_LOG_MSG_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1410 fprintf (mpfr_log_file, "%s.%d: "format, __func__, __LINE__, __VA_ARGS__);
1411 #define MPFR_LOG_MSG(x) MPFR_LOG_MSG2 x
1412
1413 #define MPFR_LOG_BEGIN2(format, ...) \
1414 mpfr_log_current ++; \
1415 if ((MPFR_LOG_INPUT_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1416 fprintf (mpfr_log_file, "%s:IN "format"\n",__func__,__VA_ARGS__); \
1417 if ((MPFR_LOG_TIME_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1418 __gmpfr_log_time = mpfr_get_cputime ();
1419 #define MPFR_LOG_BEGIN(x) \
1420 int __gmpfr_log_time = 0; \
1421 MPFR_LOG_BEGIN2 x
1422
1423 #define MPFR_LOG_END2(format, ...) \
1424 if ((MPFR_LOG_TIME_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1425 fprintf (mpfr_log_file, "%s:TIM %dms\n", __mpfr_log_fname, \
1426 mpfr_get_cputime () - __gmpfr_log_time); \
1427 if ((MPFR_LOG_OUTPUT_F&mpfr_log_type)&&(mpfr_log_current<=mpfr_log_level)) \
1428 fprintf (mpfr_log_file, "%s:OUT "format"\n",__mpfr_log_fname,__VA_ARGS__);\
1429 mpfr_log_current --;
1430 #define MPFR_LOG_END(x) \
1431 static const char *__mpfr_log_fname = __func__; \
1432 MPFR_LOG_END2 x
1433
1434 #define MPFR_LOG_FUNC(begin,end) \
1435 static const char *__mpfr_log_fname = __func__; \
1436 auto void __mpfr_log_cleanup (int *time); \
1437 void __mpfr_log_cleanup (int *time) { \
1438 int __gmpfr_log_time = *time; \
1439 MPFR_LOG_END2 end; } \
1440 int __gmpfr_log_time __attribute__ ((cleanup (__mpfr_log_cleanup))); \
1441 __gmpfr_log_time = 0; \
1442 MPFR_LOG_BEGIN2 begin
1443
1444 #else /* MPFR_USE_LOGGING */
1445
1446 /* Define void macro for logging */
1447
1448 #define MPFR_LOG_VAR(x)
1449 #define MPFR_LOG_BEGIN(x)
1450 #define MPFR_LOG_END(x)
1451 #define MPFR_LOG_MSG(x)
1452 #define MPFR_LOG_FUNC(x,y)
1453
1454 #endif /* MPFR_USE_LOGGING */
1455
1456
1457 /**************************************************************
1458 ************ Group Initialize Functions Macros *************
1459 **************************************************************/
1460
1461 #ifndef MPFR_GROUP_STATIC_SIZE
1462 # define MPFR_GROUP_STATIC_SIZE 16
1463 #endif
1464
1465 struct mpfr_group_t {
1466 size_t alloc;
1467 mp_limb_t *mant;
1468 mp_limb_t tab[MPFR_GROUP_STATIC_SIZE];
1469 };
1470
1471 #define MPFR_GROUP_DECL(g) struct mpfr_group_t g
1472 #define MPFR_GROUP_CLEAR(g) do { \
1473 if (MPFR_UNLIKELY ((g).alloc != 0)) { \
1474 MPFR_ASSERTD ((g).mant != (g).tab); \
1475 (*__gmp_free_func) ((g).mant, (g).alloc); \
1476 }} while (0)
1477
1478 #define MPFR_GROUP_INIT_TEMPLATE(g, prec, num, handler) do { \
1479 mp_prec_t _prec = (prec); \
1480 mp_size_t _size; \
1481 MPFR_ASSERTD (_prec >= MPFR_PREC_MIN); \
1482 if (MPFR_UNLIKELY (_prec > MPFR_PREC_MAX)) \
1483 mpfr_abort_prec_max (); \
1484 _size = (mp_prec_t) (_prec + BITS_PER_MP_LIMB - 1) / BITS_PER_MP_LIMB; \
1485 if (MPFR_UNLIKELY (_size * (num) > MPFR_GROUP_STATIC_SIZE)) \
1486 { \
1487 (g).alloc = (num) * _size * sizeof (mp_limb_t); \
1488 (g).mant = (mp_limb_t *) (*__gmp_allocate_func) ((g).alloc); \
1489 } \
1490 else \
1491 { \
1492 (g).alloc = 0; \
1493 (g).mant = (g).tab; \
1494 } \
1495 handler; \
1496 } while (0)
1497 #define MPFR_GROUP_TINIT(g, n, x) \
1498 MPFR_TMP_INIT1 ((g).mant + _size * (n), x, _prec)
1499
1500 #define MPFR_GROUP_INIT_1(g, prec, x) \
1501 MPFR_GROUP_INIT_TEMPLATE(g, prec, 1, MPFR_GROUP_TINIT(g, 0, x))
1502 #define MPFR_GROUP_INIT_2(g, prec, x, y) \
1503 MPFR_GROUP_INIT_TEMPLATE(g, prec, 2, \
1504 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y))
1505 #define MPFR_GROUP_INIT_3(g, prec, x, y, z) \
1506 MPFR_GROUP_INIT_TEMPLATE(g, prec, 3, \
1507 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1508 MPFR_GROUP_TINIT(g, 2, z))
1509 #define MPFR_GROUP_INIT_4(g, prec, x, y, z, t) \
1510 MPFR_GROUP_INIT_TEMPLATE(g, prec, 4, \
1511 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1512 MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t))
1513 #define MPFR_GROUP_INIT_5(g, prec, x, y, z, t, a) \
1514 MPFR_GROUP_INIT_TEMPLATE(g, prec, 5, \
1515 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1516 MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \
1517 MPFR_GROUP_TINIT(g, 4, a))
1518 #define MPFR_GROUP_INIT_6(g, prec, x, y, z, t, a, b) \
1519 MPFR_GROUP_INIT_TEMPLATE(g, prec, 6, \
1520 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1521 MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \
1522 MPFR_GROUP_TINIT(g, 4, a);MPFR_GROUP_TINIT(g, 5, b))
1523
1524 #define MPFR_GROUP_REPREC_TEMPLATE(g, prec, num, handler) do { \
1525 mp_prec_t _prec = (prec); \
1526 size_t _oalloc = (g).alloc; \
1527 mp_size_t _size; \
1528 MPFR_ASSERTD (_prec >= MPFR_PREC_MIN); \
1529 if (MPFR_UNLIKELY (_prec > MPFR_PREC_MAX)) \
1530 mpfr_abort_prec_max (); \
1531 _size = (mp_prec_t) (_prec + BITS_PER_MP_LIMB - 1) / BITS_PER_MP_LIMB; \
1532 (g).alloc = (num) * _size * sizeof (mp_limb_t); \
1533 if (MPFR_LIKELY (_oalloc == 0)) \
1534 (g).mant = (mp_limb_t *) (*__gmp_allocate_func) ((g).alloc); \
1535 else \
1536 (g).mant = (mp_limb_t *) \
1537 (*__gmp_reallocate_func) ((g).mant, _oalloc, (g).alloc); \
1538 handler; \
1539 } while (0)
1540
1541 #define MPFR_GROUP_REPREC_1(g, prec, x) \
1542 MPFR_GROUP_REPREC_TEMPLATE(g, prec, 1, MPFR_GROUP_TINIT(g, 0, x))
1543 #define MPFR_GROUP_REPREC_2(g, prec, x, y) \
1544 MPFR_GROUP_REPREC_TEMPLATE(g, prec, 2, \
1545 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y))
1546 #define MPFR_GROUP_REPREC_3(g, prec, x, y, z) \
1547 MPFR_GROUP_REPREC_TEMPLATE(g, prec, 3, \
1548 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1549 MPFR_GROUP_TINIT(g, 2, z))
1550 #define MPFR_GROUP_REPREC_4(g, prec, x, y, z, t) \
1551 MPFR_GROUP_REPREC_TEMPLATE(g, prec, 4, \
1552 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1553 MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t))
1554 #define MPFR_GROUP_REPREC_5(g, prec, x, y, z, t, a) \
1555 MPFR_GROUP_REPREC_TEMPLATE(g, prec, 5, \
1556 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1557 MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \
1558 MPFR_GROUP_TINIT(g, 4, a))
1559 #define MPFR_GROUP_REPREC_6(g, prec, x, y, z, t, a, b) \
1560 MPFR_GROUP_REPREC_TEMPLATE(g, prec, 6, \
1561 MPFR_GROUP_TINIT(g, 0, x);MPFR_GROUP_TINIT(g, 1, y); \
1562 MPFR_GROUP_TINIT(g, 2, z);MPFR_GROUP_TINIT(g, 3, t); \
1563 MPFR_GROUP_TINIT(g, 4, a);MPFR_GROUP_TINIT(g, 5, b))
1564
1565
1566 /******************************************************
1567 *************** Internal Functions *****************
1568 ******************************************************/
1569
1570 #if defined (__cplusplus)
1571 extern "C" {
1572 #endif
1573
1574 __MPFR_DECLSPEC int mpfr_underflow _MPFR_PROTO ((mpfr_ptr, mp_rnd_t, int));
1575 __MPFR_DECLSPEC int mpfr_overflow _MPFR_PROTO ((mpfr_ptr, mp_rnd_t, int));
1576
1577 __MPFR_DECLSPEC int mpfr_add1 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1578 mpfr_srcptr, mp_rnd_t));
1579 __MPFR_DECLSPEC int mpfr_sub1 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1580 mpfr_srcptr, mp_rnd_t));
1581 __MPFR_DECLSPEC int mpfr_add1sp _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1582 mpfr_srcptr, mp_rnd_t));
1583 __MPFR_DECLSPEC int mpfr_sub1sp _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1584 mpfr_srcptr, mp_rnd_t));
1585 __MPFR_DECLSPEC int mpfr_can_round_raw _MPFR_PROTO ((const mp_limb_t *,
1586 mp_size_t, int, mp_exp_t, mp_rnd_t, mp_rnd_t, mp_prec_t));
1587
1588 __MPFR_DECLSPEC int mpfr_cmp2 _MPFR_PROTO ((mpfr_srcptr, mpfr_srcptr,
1589 mp_prec_t *));
1590
1591 __MPFR_DECLSPEC long __gmpfr_ceil_log2 _MPFR_PROTO ((double));
1592 __MPFR_DECLSPEC long __gmpfr_floor_log2 _MPFR_PROTO ((double));
1593 __MPFR_DECLSPEC double __gmpfr_ceil_exp2 _MPFR_PROTO ((double));
1594 __MPFR_DECLSPEC unsigned long __gmpfr_isqrt _MPFR_PROTO ((unsigned long));
1595 __MPFR_DECLSPEC unsigned long __gmpfr_cuberoot _MPFR_PROTO ((unsigned long));
1596 __MPFR_DECLSPEC int __gmpfr_int_ceil_log2 _MPFR_PROTO ((unsigned long));
1597
1598 __MPFR_DECLSPEC int mpfr_exp_2 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,mp_rnd_t));
1599 __MPFR_DECLSPEC int mpfr_exp_3 _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,mp_rnd_t));
1600 __MPFR_DECLSPEC int mpfr_powerof2_raw _MPFR_PROTO ((mpfr_srcptr));
1601
1602 __MPFR_DECLSPEC int mpfr_pow_general _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1603 mpfr_srcptr, mp_rnd_t, int, mpfr_save_expo_t *));
1604
1605 __MPFR_DECLSPEC void mpfr_setmax _MPFR_PROTO ((mpfr_ptr, mp_exp_t));
1606 __MPFR_DECLSPEC void mpfr_setmin _MPFR_PROTO ((mpfr_ptr, mp_exp_t));
1607
1608 __MPFR_DECLSPEC long mpfr_mpn_exp _MPFR_PROTO ((mp_limb_t *, mp_exp_t *, int,
1609 mp_exp_t, size_t));
1610
1611 #ifdef _MPFR_H_HAVE_FILE
1612 __MPFR_DECLSPEC void mpfr_fprint_binary _MPFR_PROTO ((FILE *, mpfr_srcptr));
1613 #endif
1614 __MPFR_DECLSPEC void mpfr_print_binary _MPFR_PROTO ((mpfr_srcptr));
1615 __MPFR_DECLSPEC void mpfr_print_mant_binary _MPFR_PROTO ((const char*,
1616 const mp_limb_t*, mp_prec_t));
1617 __MPFR_DECLSPEC void mpfr_set_str_binary _MPFR_PROTO((mpfr_ptr, const char*));
1618
1619 __MPFR_DECLSPEC int mpfr_round_raw _MPFR_PROTO ((mp_limb_t *,
1620 const mp_limb_t *, mp_prec_t, int, mp_prec_t, mp_rnd_t, int *));
1621 __MPFR_DECLSPEC int mpfr_round_raw_2 _MPFR_PROTO ((const mp_limb_t *,
1622 mp_prec_t, int, mp_prec_t, mp_rnd_t));
1623 __MPFR_DECLSPEC int mpfr_round_raw_3 _MPFR_PROTO ((const mp_limb_t *,
1624 mp_prec_t, int, mp_prec_t, mp_rnd_t, int *));
1625 __MPFR_DECLSPEC int mpfr_round_raw_4 _MPFR_PROTO ((mp_limb_t *,
1626 const mp_limb_t *, mp_prec_t, int, mp_prec_t, mp_rnd_t));
1627
1628 #define mpfr_round_raw2(xp, xn, neg, r, prec) \
1629 mpfr_round_raw_2((xp),(xn)*BITS_PER_MP_LIMB,(neg),(prec),(r))
1630
1631 __MPFR_DECLSPEC int mpfr_check _MPFR_PROTO ((mpfr_srcptr));
1632
1633 __MPFR_DECLSPEC int mpfr_sum_sort _MPFR_PROTO ((mpfr_srcptr *const,
1634 unsigned long, mpfr_srcptr *));
1635
1636 __MPFR_DECLSPEC int mpfr_get_cputime _MPFR_PROTO ((void));
1637
1638 __MPFR_DECLSPEC void mpfr_nexttozero _MPFR_PROTO ((mpfr_ptr));
1639 __MPFR_DECLSPEC void mpfr_nexttoinf _MPFR_PROTO ((mpfr_ptr));
1640
1641 __MPFR_DECLSPEC int mpfr_const_pi_internal _MPFR_PROTO ((mpfr_ptr,mp_rnd_t));
1642 __MPFR_DECLSPEC int mpfr_const_log2_internal _MPFR_PROTO((mpfr_ptr,mp_rnd_t));
1643 __MPFR_DECLSPEC int mpfr_const_euler_internal _MPFR_PROTO((mpfr_ptr, mp_rnd_t));
1644 __MPFR_DECLSPEC int mpfr_const_catalan_internal _MPFR_PROTO((mpfr_ptr, mp_rnd_t));
1645
1646 #if 0
1647 __MPFR_DECLSPEC void mpfr_init_cache _MPFR_PROTO ((mpfr_cache_t,
1648 int(*)(mpfr_ptr,mpfr_rnd_t)));
1649 #endif
1650 __MPFR_DECLSPEC void mpfr_clear_cache _MPFR_PROTO ((mpfr_cache_t));
1651 __MPFR_DECLSPEC int mpfr_cache _MPFR_PROTO ((mpfr_ptr, mpfr_cache_t,
1652 mpfr_rnd_t));
1653
1654 __MPFR_DECLSPEC void mpfr_mulhigh_n _MPFR_PROTO ((mp_ptr, mp_srcptr,
1655 mp_srcptr, mp_size_t));
1656 __MPFR_DECLSPEC void mpfr_sqrhigh_n _MPFR_PROTO ((mp_ptr, mp_srcptr, mp_size_t));
1657
1658 __MPFR_DECLSPEC int mpfr_round_p _MPFR_PROTO ((mp_limb_t *, mp_size_t,
1659 mp_exp_t, mp_prec_t));
1660
1661 __MPFR_DECLSPEC void mpfr_dump_mant _MPFR_PROTO ((const mp_limb_t *,
1662 mp_prec_t, mp_prec_t,
1663 mp_prec_t));
1664
1665 __MPFR_DECLSPEC int mpfr_round_near_x _MPFR_PROTO ((mpfr_ptr, mpfr_srcptr,
1666 mpfr_uexp_t, int,
1667 mp_rnd_t));
1668 __MPFR_DECLSPEC void mpfr_abort_prec_max _MPFR_PROTO ((void))
1669 MPFR_NORETURN_ATTR;
1670
1671 #if defined (__cplusplus)
1672 }
1673 #endif
1674
1675 #endif
DragonFly BSD