48d3af97610e6d00f6c9161d6c1efe470a8dfefa
[kopensolaris-gnu/glibc.git] / stdlib / gmp-impl.h
1 /* Include file for internal GNU MP types and definitions.
2
3 Copyright (C) 1991, 1993, 1994 Free Software Foundation, Inc.
4
5 This file is part of the GNU MP Library.
6
7 The GNU MP Library is free software; you can redistribute it and/or modify
8 it under the terms of the GNU Library General Public License as published by
9 the Free Software Foundation; either version 2 of the License, or (at your
10 option) any later version.
11
12 The GNU MP Library is distributed in the hope that it will be useful, but
13 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14 or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Library General Public
15 License for more details.
16
17 You should have received a copy of the GNU Library General Public License
18 along with the GNU MP Library; see the file COPYING.LIB.  If not, write to
19 the Free Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. */
20
21 #if ! defined (alloca)
22 #if defined (__GNUC__)
23 #define alloca __builtin_alloca
24 #endif
25 #endif
26
27 #if ! defined (alloca)
28 #if defined (__sparc__) || defined (sparc) || defined (__sgi)
29 #include <alloca.h>
30 #endif
31 #endif
32
33 #ifndef NULL
34 #define NULL 0L
35 #endif
36
37 #if ! defined (__GNUC__)
38 #define inline                  /* Empty */
39 void *alloca();
40 #endif
41
42 #define ABS(x) (x >= 0 ? x : -x)
43 #define MIN(l,o) ((l) < (o) ? (l) : (o))
44 #define MAX(h,i) ((h) > (i) ? (h) : (i))
45
46 #include "gmp-mparam.h"
47 /* #include "longlong.h" */
48
49 #ifdef __STDC__
50 void *malloc (size_t);
51 void *realloc (void *, size_t);
52 void free (void *);
53
54 extern void *   (*_mp_allocate_func) (size_t);
55 extern void *   (*_mp_reallocate_func) (void *, size_t, size_t);
56 extern void     (*_mp_free_func) (void *, size_t);
57
58 void *_mp_default_allocate (size_t);
59 void *_mp_default_reallocate (void *, size_t, size_t);
60 void _mp_default_free (void *, size_t);
61
62 #else
63
64 #define const                   /* Empty */
65 #define signed                  /* Empty */
66
67 void *malloc ();
68 void *realloc ();
69 void free ();
70
71 extern void *   (*_mp_allocate_func) ();
72 extern void *   (*_mp_reallocate_func) ();
73 extern void     (*_mp_free_func) ();
74
75 void *_mp_default_allocate ();
76 void *_mp_default_reallocate ();
77 void _mp_default_free ();
78 #endif
79
80 /* Copy NLIMBS *limbs* from SRC to DST.  */
81 #define MPN_COPY_INCR(DST, SRC, NLIMBS) \
82   do {                                                                  \
83     mp_size_t __i;                                                      \
84     for (__i = 0; __i < (NLIMBS); __i++)                                \
85       (DST)[__i] = (SRC)[__i];                                          \
86   } while (0)
87 #define MPN_COPY_DECR(DST, SRC, NLIMBS) \
88   do {                                                                  \
89     mp_size_t __i;                                                      \
90     for (__i = (NLIMBS) - 1; __i >= 0; __i--)                           \
91       (DST)[__i] = (SRC)[__i];                                          \
92   } while (0)
93 #define MPN_COPY MPN_COPY_INCR
94
95 /* Zero NLIMBS *limbs* AT DST.  */
96 #define MPN_ZERO(DST, NLIMBS) \
97   do {                                                                  \
98     mp_size_t __i;                                                      \
99     for (__i = 0; __i < (NLIMBS); __i++)                                \
100       (DST)[__i] = 0;                                                   \
101   } while (0)
102
103 #define MPN_NORMALIZE(DST, NLIMBS) \
104   do {                                                                  \
105     while (NLIMBS > 0)                                                  \
106       {                                                                 \
107         if ((DST)[(NLIMBS) - 1] != 0)                                   \
108           break;                                                        \
109         NLIMBS--;                                                       \
110       }                                                                 \
111   } while (0)
112 #define MPN_NORMALIZE_NOT_ZERO(DST, NLIMBS) \
113   do {                                                                  \
114     while (1)                                                           \
115       {                                                                 \
116         if ((DST)[(NLIMBS) - 1] != 0)                                   \
117           break;                                                        \
118         NLIMBS--;                                                       \
119       }                                                                 \
120   } while (0)
121
122 /*  Swap (mp_ptr, mp_size_t) (U, UL) with (V, VL)  */
123 #define MPN_SWAP(u, l, v, m) \
124   do {                                                                  \
125     { mp_ptr _; _ = (u), (u) = (v), (v) = _;}                           \
126     { mp_size_t _; _ = (l), (l) = (m), (m) = _;}                        \
127   } while (0)
128
129 /*  Return true iff the limb X has less bits than the limb Y.  */
130 #define MPN_LESS_BITS_LIMB(x,y) ((x) < (y) && (x) < ((x) ^ (y)))
131
132 /*  Return true iff (mp_ptr, mp_size_t) (U, UL) has less bits than (V, VL).  */
133 #define MPN_LESS_BITS(u, l, v, m) \
134   ((l) < (m)                                                            \
135    || ((l) == (m) && (l) != 0 && MPN_LESS_BITS_LIMB ((u)[(l - 1)], (v)[(l) - 1])))
136
137 /*  Return true iff (mp_ptr, mp_size_t) (U, UL) has more bits than (V, VL).  */
138 #define MPN_MORE_BITS(u, l, v, m) MPN_LESS_BITS (v, m, u, l)
139
140 /*  Perform twos complement on (mp_ptr, mp_size_t) (U, UL), 
141     putting result at (v, VL).  Precondition: U[0] != 0.  */
142 #define MPN_COMPL_INCR(u, v, l) \
143   do {                                                                  \
144     mp_size_t _ = 0;                                                    \
145     (u)[0] = -(v)[_];                                                   \
146     while (_++ < (l))                                                   \
147       (u)[_] = ~(v)[_];                                                 \
148   } while (0)
149 #define MPN_COMPL MPN_COMPL_INCR
150
151 /* Initialize the MP_INT X with space for NLIMBS limbs.
152    X should be a temporary variable, and it will be automatically
153    cleared out when the running function returns.
154    We use __x here to make it possible to accept both mpz_ptr and mpz_t
155    arguments.  */
156 #define MPZ_TMP_INIT(X, NLIMBS) \
157   do {                                                                  \
158     mpz_ptr __x = (X);                                                  \
159     __x->alloc = (NLIMBS);                                              \
160     __x->d = (mp_ptr) alloca ((NLIMBS) * BYTES_PER_MP_LIMB);            \
161   } while (0)
162
163 #define MPN_MUL_N_RECURSE(prodp, up, vp, size, tspace) \
164   do {                                                                  \
165     if ((size) < KARATSUBA_THRESHOLD)                                   \
166       ____mpn_mul_n_basecase (prodp, up, vp, size);                     \
167     else                                                                \
168       ____mpn_mul_n (prodp, up, vp, size, tspace);                      \
169   } while (0);
170 #define MPN_SQR_N_RECURSE(prodp, up, size, tspace) \
171   do {                                                                  \
172     if ((size) < KARATSUBA_THRESHOLD)                                   \
173       ____mpn_sqr_n_basecase (prodp, up, size);                         \
174     else                                                                \
175       ____mpn_sqr_n (prodp, up, size, tspace);                          \
176   } while (0);
177 #define assert(trueval) do {if (!(trueval)) abort ();} while (0)
178
179 /* Structure for conversion between internal binary format and
180    strings in base 2..36.  */
181 struct bases
182 {
183   /* Number of digits in the conversion base that always fits in
184      an mp_limb.  For example, for base 10 this is 10, since
185      2**32 = 4294967296 has ten digits.  */
186   int chars_per_limb;
187
188   /* log(2)/log(conversion_base) */
189   float chars_per_bit_exactly;
190
191   /* big_base is conversion_base**chars_per_limb, i.e. the biggest
192      number that fits a word, built by factors of conversion_base.
193      Exception: For 2, 4, 8, etc, big_base is log2(base), i.e. the
194      number of bits used to represent each digit in the base.  */
195   mp_limb big_base;
196
197   /* big_base_inverted is a BITS_PER_MP_LIMB bit approximation to
198      1/big_base, represented as a fixed-point number.  Instead of
199      dividing by big_base an application can choose to multiply
200      by big_base_inverted.  */
201   mp_limb big_base_inverted;
202 };
203
204 extern const struct bases __mp_bases[];
205 extern mp_size_t __gmp_default_fp_limb_precision;
206
207 /* Divide the two-limb number in (NH,,NL) by D, with DI being the largest
208    limb not larger than (2**(2*BITS_PER_MP_LIMB))/D - (2**BITS_PER_MP_LIMB).
209    If this would yield overflow, DI should be the largest possible number
210    (i.e., only ones).  For correct operation, the most significant bit of D
211    has to be set.  Put the quotient in Q and the remainder in R.  */
212 #define udiv_qrnnd_preinv(q, r, nh, nl, d, di) \
213   do {                                                                  \
214     mp_limb _q, _ql, _r;                                                \
215     mp_limb _xh, _xl;                                                   \
216     umul_ppmm (_q, _ql, (nh), (di));                                    \
217     _q += (nh);                 /* DI is 2**BITS_PER_MP_LIMB too small */\
218     umul_ppmm (_xh, _xl, _q, (d));                                      \
219     sub_ddmmss (_xh, _r, (nh), (nl), _xh, _xl);                         \
220     if (_xh != 0)                                                       \
221       {                                                                 \
222         sub_ddmmss (_xh, _r, _xh, _r, 0, (d));                          \
223         _q += 1;                                                        \
224         if (_xh != 0)                                                   \
225           {                                                             \
226             sub_ddmmss (_xh, _r, _xh, _r, 0, (d));                      \
227             _q += 1;                                                    \
228           }                                                             \
229       }                                                                 \
230     if (_r >= (d))                                                      \
231       {                                                                 \
232         _r -= (d);                                                      \
233         _q += 1;                                                        \
234       }                                                                 \
235     (r) = _r;                                                           \
236     (q) = _q;                                                           \
237   } while (0)
238 /* Like udiv_qrnnd_preinv, but for for any value D.  DNORM is D shifted left
239    so that its most significant bit is set.  LGUP is ceil(log2(D)).  */
240 #define udiv_qrnnd_preinv2gen(q, r, nh, nl, d, di, dnorm, lgup) \
241   do {                                                                  \
242     mp_limb n2, n10, n1, nadj, q1;                                      \
243     mp_limb _xh, _xl;                                                   \
244     n2 = ((nh) << (BITS_PER_MP_LIMB - (lgup))) + ((nl) >> 1 >> (l - 1));\
245     n10 = (nl) << (BITS_PER_MP_LIMB - (lgup));                          \
246     n1 = ((mp_limb_signed) n10 >> (BITS_PER_MP_LIMB - 1));              \
247     nadj = n10 + (n1 & (dnorm));                                        \
248     umul_ppmm (_xh, _xl, di, n2 - n1);                                  \
249     add_ssaaaa (_xh, _xl, _xh, _xl, 0, nadj);                           \
250     q1 = ~(n2 + _xh);                                                   \
251     umul_ppmm (_xh, _xl, q1, d);                                        \
252     add_ssaaaa (_xh, _xl, _xh, _xl, nh, nl);                            \
253     _xh -= (d);                                                         \
254     (r) = _xl + ((d) & _xh);                                            \
255     (q) = _xh - q1;                                                     \
256   } while (0)
257 /* Exactly like udiv_qrnnd_preinv, but branch-free.  It is not clear which
258    version to use.  */
259 #define udiv_qrnnd_preinv2norm(q, r, nh, nl, d, di) \
260   do {                                                                  \
261     mp_limb n2, n10, n1, nadj, q1;                                      \
262     mp_limb _xh, _xl;                                                   \
263     n2 = (nh);                                                          \
264     n10 = (nl);                                                         \
265     n1 = ((mp_limb_signed) n10 >> (BITS_PER_MP_LIMB - 1));              \
266     nadj = n10 + (n1 & (d));                                            \
267     umul_ppmm (_xh, _xl, di, n2 - n1);                                  \
268     add_ssaaaa (_xh, _xl, _xh, _xl, 0, nadj);                           \
269     q1 = ~(n2 + _xh);                                                   \
270     umul_ppmm (_xh, _xl, q1, d);                                        \
271     add_ssaaaa (_xh, _xl, _xh, _xl, nh, nl);                            \
272     _xh -= (d);                                                         \
273     (r) = _xl + ((d) & _xh);                                            \
274     (q) = _xh - q1;                                                     \
275   } while (0)
276
277 #if defined (__GNUC__)
278 /* Define stuff for longlong.h.  */
279 typedef unsigned int UQItype    __attribute__ ((mode (QI)));
280 typedef          int SItype     __attribute__ ((mode (SI)));
281 typedef unsigned int USItype    __attribute__ ((mode (SI)));
282 typedef          int DItype     __attribute__ ((mode (DI)));
283 typedef unsigned int UDItype    __attribute__ ((mode (DI)));
284 #else
285 typedef unsigned char UQItype;
286 typedef          long SItype;
287 typedef unsigned long USItype;
288 #endif
289
290 typedef mp_limb UWtype;
291 typedef unsigned int UHWtype;
292 #define W_TYPE_SIZE BITS_PER_MP_LIMB
293
294
295 #ifndef IEEE_DOUBLE_BIG_ENDIAN
296 #define IEEE_DOUBLE_BIG_ENDIAN 1
297 #endif
298
299 #if IEEE_DOUBLE_BIG_ENDIAN
300 union ieee_double_extract
301 {
302   struct
303     {
304       unsigned long sig:1;
305       unsigned long exp:11;
306       unsigned long manh:20;
307       unsigned long manl:32;
308     } s;
309   double d;
310 };
311 #else
312 union ieee_double_extract
313 {
314   struct
315     {
316       unsigned long manl:32;
317       unsigned long manh:20;
318       unsigned long exp:11;
319       unsigned long sig:1;
320     } s;
321   double d;
322 };
323 #endif