blob: a3e903e4ac8918b0ed0d7f84fa2b21325a2d23d3 [file] [log] [blame]
Guido van Rossum5f59d601992-12-14 16:59:51 +00001/***********************************************************
2Copyright 1992 by Stichting Mathematisch Centrum, Amsterdam, The
3Netherlands.
4
5 All Rights Reserved
6
7Permission to use, copy, modify, and distribute this software and its
8documentation for any purpose and without fee is hereby granted,
9provided that the above copyright notice appear in all copies and that
10both that copyright notice and this permission notice appear in
11supporting documentation, and that the names of Stichting Mathematisch
12Centrum or CWI not be used in advertising or publicity pertaining to
13distribution of the software without specific, written prior permission.
14
15STICHTING MATHEMATISCH CENTRUM DISCLAIMS ALL WARRANTIES WITH REGARD TO
16THIS SOFTWARE, INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND
17FITNESS, IN NO EVENT SHALL STICHTING MATHEMATISCH CENTRUM BE LIABLE
18FOR ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
19WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
20ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT
21OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
22
23******************************************************************/
24/* MPZ module */
25
26/* This module provides an interface to an alternate Multi-Precision
27 library, GNU MP in this case */
28
29/* XXX note: everywhere where mpz_size is called,
30 sizeof (limb) == sizeof (long) has been assumed. */
31
32
33/* MPZ objects */
34
35#include "allobjects.h"
36#include "modsupport.h" /* For getargs() etc. */
37#include <assert.h>
38
39/*
40** These are the cpp-flags used in this file...
41**
42**
43** MPZ_MDIV_BUG works around the mpz_m{div,mod,...} routines.
44** This bug has been fixed in a later release of
45** GMP.
46**
47** MPZ_GET_STR_BUG mpz_get_str corrupts memory, seems to be fixed
48** in a later release
49**
50** MPZ_DEBUG generates a bunch of diagnostic messages
51**
52** MPZ_SPARE_MALLOC if set, results in extra code that tries to
53** minimize the creation of extra objects.
54**
55** MPZ_TEST_DIV extra diagnostic output on stderr, when division
56** routines are involved
57**
58** MPZ_LIB_DOES_CHECKING if set, assumes that mpz library doesn't call
59** alloca with arg < 0 (when casted to a signed
60** integral type).
61**
62** MPZ_CONVERSIONS_AS_METHODS if set, presents the conversions as
63** methods. e.g., `mpz(5).long() == 5L'
64** Later, Guido provided an interface to the
65** standard functions. So this flag has no been
66** cleared, and `long(mpz(5)) == 5L'
67**
68** MP_TEST_ALLOC If set, you would discover why MPZ_GET_STR_BUG
69** is needed
70**
71** MAKEDUMMYINT Must be set if dynamic linking will be used
72*/
73
74
75/*
76** IMHO, mpz_m{div,mod,divmod}() do the wrong things when the denominator < 0
77** I assume that this will be fixed in a future release
78*/
79/*#define MPZ_MDIV_BUG fixed the (for me) nexessary parts in libgmp.a */
80/*
81** IMO, mpz_get_str() assumes a bit too large target space, if he doesn't
82** allocate it himself
83*/
84#define MPZ_GET_STR_BUG
85
86#include "gmp.h"
87typedef struct {
88 OB_HEAD
89 MP_INT mpz; /* the actual number */
90} mpzobject;
91
92extern typeobject MPZtype; /* Really static, forward */
93
94#define is_mpzobject(v) ((v)->ob_type == &MPZtype)
95
96static const char initialiser_name[] = "mpz";
97
98/* #define MPZ_DEBUG */
99
100static mpzobject *
101newmpzobject()
102{
103 mpzobject *mpzp;
104
105
106#ifdef MPZ_DEBUG
107 fputs( "mpz_object() called...\n", stderr );
108#endif /* def MPZ_DEBUG */
109 mpzp = NEWOBJ(mpzobject, &MPZtype);
110 if (mpzp == NULL)
111 return NULL;
112
113 mpz_init(&mpzp->mpz); /* actual initialisation */
114 return mpzp;
115} /* newmpzobject() */
116
117#ifdef MPZ_GET_STR_BUG
118#include "gmp-impl.h"
119#include "longlong.h"
120#endif /* def MPZ_GET_STR_BUG */
121
122static object *
123mpz_format(objp, base, withname)
124 object *objp;
125 int base;
126 unsigned char withname;
127{
128 mpzobject *mpzp = (mpzobject *)objp;
129 stringobject *strobjp;
130 int i;
131 int cmpres;
132 int taglong;
133 char *cp;
134 char prefix[5], *tcp;
135
136
137 tcp = &prefix[0];
138
139 if (mpzp == NULL || !is_mpzobject(mpzp)) {
140 err_badcall();
141 return NULL;
142 }
143
144 assert(base >= 2 && base <= 36);
145
146 if (withname)
147 i = strlen(initialiser_name) + 2; /* e.g. 'mpz(' + ')' */
148 else
149 i = 0;
150
151 if ((cmpres = mpz_cmp_si(&mpzp->mpz, 0L)) == 0)
152 base = 10; /* '0' in every base, right */
153 else if (cmpres < 0) {
154 *tcp++ = '-';
155 i += 1; /* space to hold '-' */
156 }
157
158#ifdef MPZ_DEBUG
159 fprintf(stderr, "mpz_format: mpz_sizeinbase %d\n",
160 (int)mpz_sizeinbase(&mpzp->mpz, base));
161#endif /* def MPZ_DEBUG */
162#ifdef MPZ_GET_STR_BUG
163 i += ((size_t) abs(mpzp->mpz.size) * BITS_PER_MP_LIMB
164 * __mp_bases[base].chars_per_bit_exactly) + 1;
165#else /* def MPZ_GET_STR_BUG */
166 i += (int)mpz_sizeinbase(&mpzp->mpz, base);
167#endif /* def MPZ_GET_STR_BUG else */
168
169 if (base == 16) {
170 *tcp++ = '0';
171 *tcp++ = 'x';
172 i += 2; /* space to hold '0x' */
173 }
174 else if (base == 8) {
175 *tcp++ = '0';
176 i += 1; /* space to hold the extra '0' */
177 }
178 else if (base > 10) {
179 *tcp++ = '0' + base / 10;
180 *tcp++ = '0' + base % 10;
181 *tcp++ = '#';
182 i += 3; /* space to hold e.g. '12#' */
183 }
184 else if (base < 10) {
185 *tcp++ = '0' + base;
186 *tcp++ = '#';
187 i += 2; /* space to hold e.g. '6#' */
188 }
189
190 /*
191 ** the following code looks if we need a 'L' attached to the number
192 ** it will also attach an 'L' to the value -0x80000000
193 */
194 taglong = 0;
195 if (mpz_size(&mpzp->mpz) > 1
196 || (long)mpz_get_ui(&mpzp->mpz) < 0L) {
197 taglong = 1;
198 i += 1; /* space to hold 'L' */
199 }
200
201#ifdef MPZ_DEBUG
202 fprintf(stderr, "mpz_format: requesting string size %d\n", i);
203#endif /* def MPZ_DEBUG */
204 if ((strobjp = (stringobject *)newsizedstringobject((char *)0, i))
205 == NULL)
206 return NULL;
207
208 /* get the beginning of the string memory and start copying things */
209 cp = GETSTRINGVALUE(strobjp);
210 if (withname) {
211 strcpy(cp, initialiser_name);
212 cp += strlen(initialiser_name);
213 *cp++ = '('; /*')'*/
214 }
215
216 /* copy the already prepared prefix; e.g. sign and base indicator */
217 *tcp = '\0';
218 strcpy(cp, prefix);
219 cp += tcp - prefix;
220
221 /* since' we have the sign already, let the lib think it's a positive
222 number */
223 if (cmpres < 0)
224 mpz_neg(&mpzp->mpz,&mpzp->mpz); /* hack Hack HAck HACk HACK */
225 (void)mpz_get_str(cp, base, &mpzp->mpz);
226 if (cmpres < 0)
227 mpz_neg(&mpzp->mpz,&mpzp->mpz); /* hack Hack HAck HACk HACK */
228#ifdef MPZ_DEBUG
229 fprintf(stderr, "mpz_format: base (ultim) %d, mpz_get_str: %s\n",
230 base, cp);
231#endif /* def MPZ_DEBUG */
232 cp += strlen(cp);
233
234 if (taglong)
235 *cp++ = 'L';
236 if (withname)
237 *cp++ = /*'('*/ ')';
238
239 *cp = '\0';
240
241#ifdef MPZ_DEBUG
242 fprintf(stderr,
243 "mpz_format: cp (str end) 0x%x, begin 0x%x, diff %d, i %d\n",
244 cp, GETSTRINGVALUE(strobjp), cp - GETSTRINGVALUE(strobjp), i);
245#endif /* def MPZ_DEBUG */
246 assert(cp - GETSTRINGVALUE(strobjp) <= i);
247
248 if (cp - GETSTRINGVALUE(strobjp) != i) {
249 strobjp->ob_size -= i - (cp - GETSTRINGVALUE(strobjp));
250 }
251
252 return (object *)strobjp;
253} /* mpz_format() */
254
255/* MPZ methods */
256
257static void
258mpz_dealloc(mpzp)
259 mpzobject *mpzp;
260{
261#ifdef MPZ_DEBUG
262 fputs( "mpz_dealloc() called...\n", stderr );
263#endif /* def MPZ_DEBUG */
264 mpz_clear(&mpzp->mpz);
265 DEL(mpzp);
266} /* mpz_dealloc() */
267
268/* ARGSUSED */
269static int
270mpz_print(v, fp, flags)
271 object *v;
272 FILE *fp;
273 int flags; /* Not used but required by interface */
274{
275 stringobject *str
276 = (stringobject *)mpz_format(v, 10, (unsigned char)1);
277
278 if (str == NULL)
279 return -1;
280
281 fputs(GETSTRINGVALUE(str), fp);
282 DECREF(str);
283 return 0;
284} /* mpz_print() */
285
286
287/* pointers to frequently used values 0, 1 and -1 */
288static mpzobject *mpz_value_zero, *mpz_value_one, *mpz_value_mone;
289
290static int
291mpz_compare(a, b)
292 mpzobject *a, *b;
293{
294 int cmpres;
295
296
297 /* guido sez it's better to return -1, 0 or 1 */
298 return (cmpres = mpz_cmp( &a->mpz, &b->mpz )) == 0 ? 0
299 : cmpres > 0 ? 1 : -1;
300} /* mpz_compare() */
301
302static object *
303mpz_addition(a, b)
304 mpzobject *a;
305 mpzobject *b;
306{
307 mpzobject *z;
308
309
310#ifdef MPZ_SPARE_MALLOC
311 if (mpz_cmp_ui(&a->mpz, (unsigned long int)0) == 0) {
312 INCREF(b);
313 return (object *)b;
314 }
315
316 if (mpz_cmp_ui(&b->mpz, (unsigned long int)0) == 0) {
317 INCREF(a);
318 return (object *)a;
319 }
320#endif /* def MPZ_SPARE_MALLOC */
321
322 if ((z = newmpzobject()) == NULL)
323 return NULL;
324
325 mpz_add(&z->mpz, &a->mpz, &b->mpz);
326 return (object *)z;
327} /* mpz_addition() */
328
329static object *
330mpz_substract(a, b)
331 mpzobject *a;
332 mpzobject *b;
333{
334 mpzobject *z;
335
336
337#ifdef MPZ_SPARE_MALLOC
338 if (mpz_cmp_ui(&b->mpz, (unsigned long int)0) == 0) {
339 INCREF(a);
340 return (object *)a;
341 }
342#endif /* MPZ_SPARE_MALLOC */
343
344 if ((z = newmpzobject()) == NULL)
345 return NULL;
346
347 mpz_sub(&z->mpz, &a->mpz, &b->mpz);
348 return (object *)z;
349} /* mpz_substract() */
350
351static object *
352mpz_multiply(a, b)
353 mpzobject *a;
354 mpzobject *b;
355{
356#ifdef MPZ_SPARE_MALLOC
357 int cmpres;
358#endif /* def MPZ_SPARE_MALLOC */
359 mpzobject *z;
360
361
362#ifdef MPZ_SPARE_MALLOC
363 if ((cmpres = mpz_cmp_ui(&a->mpz, (unsigned long int)0)) == 0) {
364 INCREF(mpz_value_zero);
365 return (object *)mpz_value_zero;
366 }
367 if (cmpres > 0 && mpz_cmp_ui(&a->mpz, (unsigned long int)1) == 0) {
368 INCREF(b);
369 return (object *)b;
370 }
371
372 if ((cmpres = mpz_cmp_ui(&b->mpz, (unsigned long_int)0)) == 0) {
373 INCREF(mpz_value_zero);
374 return (object *)mpz_value_zero;
375 }
376 if (cmpres > 0 && mpz_cmp_ui(&b->mpz, (unsigned long int)1) == 0) {
377 INCREF(a);
378 return (object *)a;
379 }
380#endif /* MPZ_SPARE_MALLOC */
381
382 if ((z = newmpzobject()) == NULL)
383 return NULL;
384
385 mpz_mul( &z->mpz, &a->mpz, &b->mpz );
386 return (object *)z;
387
388} /* mpz_multiply() */
389
390static object *
391mpz_divide(a, b)
392 mpzobject *a;
393 mpzobject *b;
394{
395#ifdef MPZ_SPARE_MALLOC
396 int cmpres;
397#endif /* def MPZ_SPARE_MALLOC */
398 mpzobject *z;
399
400
401 if ((
402#ifdef MPZ_SPARE_MALLOC
403 cmpres =
404#endif /* def MPZ_SPARE_MALLOC */
405 mpz_cmp_ui(&b->mpz, (unsigned long int)0)) == 0) {
406 err_setstr(ZeroDivisionError, "mpz./ by zero");
407 return NULL;
408 }
409#ifdef MPZ_SPARE_MALLOC
410 if (cmpres > 0 && mpz_cmp_ui(&b->mpz(unsigned long int)1) == 0) {
411 INCREF(a);
412 return (object *)a;
413 }
414#endif /* def MPZ_SPARE_MALLOC */
415
416 if ((z = newmpzobject()) == NULL)
417 return NULL;
418
419#ifdef MPZ_TEST_DIV
420 fputs("mpz_divide: div result", stderr);
421 mpz_div(&z->mpz, &a->mpz, &b->mpz);
422 mpz_out_str(stderr, 10, &z->mpz);
423 putc('\n', stderr);
424#endif /* def MPZ_TEST_DIV */
425#ifdef MPZ_MDIV_BUG
426 if ((mpz_cmp_ui(&a->mpz, (unsigned long int)0) < 0)
427 != (mpz_cmp_ui(&b->mpz, (unsigned long int)0) < 0)) {
428 /*
429 ** numerator has other sign than denominator: we have
430 ** to look at the remainder for a correction, since mpz_mdiv
431 ** also calls mpz_divmod, I can as well do it myself
432 */
433 MP_INT tmpmpz;
434
435
436 mpz_init(&tmpmpz);
437 mpz_divmod(&z->mpz, &tmpmpz, &a->mpz, &b->mpz);
438
439 if (mpz_cmp_ui(&tmpmpz, (unsigned long int)0) != 0)
440 mpz_sub_ui(&z->mpz, &z->mpz, (unsigned long int)1);
441
442 mpz_clear(&tmpmpz);
443 }
444 else
445 mpz_div(&z->mpz, &a->mpz, &b->mpz);
446 /* the ``naive'' implementation does it right for operands
447 having the same sign */
448
449#else /* def MPZ_MDIV_BUG */
450 mpz_mdiv(&z->mpz, &a->mpz, &b->mpz);
451#endif /* def MPZ_MDIV_BUG else */
452#ifdef MPZ_TEST_DIV
453 fputs("mpz_divide: mdiv result", stderr);
454 mpz_out_str(stderr, 10, &z->mpz);
455 putc('\n', stderr);
456#endif /* def MPZ_TEST_DIV */
457 return (object *)z;
458
459} /* mpz_divide() */
460
461static object *
462mpz_remainder(a, b)
463 mpzobject *a;
464 mpzobject *b;
465{
466#ifdef MPZ_SPARE_MALLOC
467 int cmpres;
468#endif /* def MPZ_SPARE_MALLOC */
469 mpzobject *z;
470
471
472 if ((
473#ifdef MPZ_SPARE_MALLOC
474 cmpres =
475#endif /* def MPZ_SPARE_MALLOC */
476 mpz_cmp_ui(&b->mpz, (unsigned long int)0)) == 0) {
477 err_setstr(ZeroDivisionError, "mpz.% by zero");
478 return NULL;
479 }
480#ifdef MPZ_SPARE_MALLOC
481 if (cmpres > 0) {
482 if ((cmpres = mpz_cmp_ui(&b->mpz, (unsigned long int)2)) == 0) {
483 INCREF(mpz_value_one);
484 return (object *)mpz_value_one;
485 }
486 if (cmpres < 0) {
487 /* b must be 1 now */
488 INCREF(mpz_value_zero);
489 return (object *)mpz_value_zero;
490 }
491 }
492#endif /* def MPZ_SPARE_MALLOC */
493
494 if ((z = newmpzobject()) == NULL)
495 return NULL;
496
497#ifdef MPZ_TEST_DIV
498 fputs("mpz_remain: mod result", stderr);
499 mpz_mod(&z->mpz, &a->mpz, &b->mpz);
500 mpz_out_str(stderr, 10, &z->mpz);
501 putc('\n', stderr);
502#endif /* def MPZ_TEST_DIV */
503#ifdef MPZ_MDIV_BUG
504
505 /* the ``naive'' implementation does it right for operands
506 having the same sign */
507 mpz_mod(&z->mpz, &a->mpz, &b->mpz);
508
509 /* assumption: z, a and b all point to different locations */
510 if ((mpz_cmp_ui(&a->mpz, (unsigned long int)0) < 0)
511 != (mpz_cmp_ui(&b->mpz, (unsigned long int)0) < 0)
512 && mpz_cmp_ui(&z->mpz, (unsigned long int)0) != 0)
513 mpz_add(&z->mpz, &z->mpz, &b->mpz);
514 /*
515 ** numerator has other sign than denominator: we have
516 ** to look at the remainder for a correction, since mpz_mdiv
517 ** also calls mpz_divmod, I can as well do it myself
518 */
519#else /* def MPZ_MDIV_BUG */
520 mpz_mmod(&z->mpz, &a->mpz, &b->mpz);
521#endif /* def MPZ_MDIV_BUG else */
522#ifdef MPZ_TEST_DIV
523 fputs("mpz_remain: mmod result", stderr);
524 mpz_out_str(stderr, 10, &z->mpz);
525 putc('\n', stderr);
526#endif /* def MPZ_TEST_DIV */
527 return (object *)z;
528
529} /* mpz_remainder() */
530
531static object *
532mpz_div_and_mod(a, b)
533 mpzobject *a;
534 mpzobject *b;
535{
536 object *z = NULL;
537 mpzobject *x = NULL, *y = NULL;
538
539
540 if (mpz_cmp_ui(&b->mpz, (unsigned long int)0) == 0) {
541 err_setstr(ZeroDivisionError, "mpz.divmod by zero");
542 return NULL;
543 }
544
545 if ((z = newtupleobject(2)) == NULL
546 || (x = newmpzobject()) == NULL
547 || (y = newmpzobject()) == NULL) {
548 XDECREF(z);
549 XDECREF(x);
550 XDECREF(y);
551 return NULL;
552 }
553
554#ifdef MPZ_TEST_DIV
555 fputs("mpz_divmod: dm result", stderr);
556 mpz_divmod(&x->mpz, &y->mpz, &a->mpz, &b->mpz);
557 mpz_out_str(stderr, 10, &x->mpz);
558 putc('\n', stderr);
559 mpz_out_str(stderr, 10, &y->mpz);
560 putc('\n', stderr);
561#endif /* def MPZ_TEST_DIV */
562#ifdef MPZ_MDIV_BUG
563 mpz_divmod(&x->mpz, &y->mpz, &a->mpz, &b->mpz);
564 if ((mpz_cmp_ui(&a->mpz, (unsigned long int)0) < 0)
565 != (mpz_cmp_ui(&b->mpz, (unsigned long int)0) < 0)
566 && mpz_cmp_ui(&y->mpz, (unsigned long int)0) != 0) {
567 /*
568 ** numerator has other sign than denominator: we have
569 ** to look at the remainder for a correction.
570 */
571 mpz_add(&y->mpz, &y->mpz, &b->mpz);
572 mpz_sub_ui(&x->mpz, &x->mpz, (unsigned long int)1);
573 }
574#else /* def MPZ_MDIV_BUG */
575 mpz_mdivmod( &x->mpz, &y->mpz, &a->mpz, &b->mpz );
576#endif /* def MPZ_MDIV_BUG else */
577#ifdef MPZ_TEST_DIV
578 fputs("mpz_divmod: mdm result", stderr);
579 mpz_out_str(stderr, 10, &x->mpz);
580 putc('\n', stderr);
581 mpz_out_str(stderr, 10, &y->mpz);
582 putc('\n', stderr);
583#endif /* def MPZ_TEST_DIV */
584
585 (void)settupleitem(z, 0, (object *)x);
586 (void)settupleitem(z, 1, (object *)y);
587
588 return z;
589} /* mpz_div_and_mod() */
590
591static object *
592mpz_power(a, b)
593 mpzobject *a;
594 mpzobject *b;
595{
596 mpzobject *z;
597 int cmpres;
598 long int longtmp1, longtmp2;
599
600
601 if ((cmpres = mpz_cmp_ui(&b->mpz, (unsigned long int)0)) == 0) {
602 /* the gnu-mp lib sets pow(0,0) to 0, we to 1 */
603
604 INCREF(mpz_value_one);
605 return (object *)mpz_value_one;
606 }
607
608 if (cmpres < 0) {
609 err_setstr(ValueError, "mpz.pow to negative exponent");
610 return NULL;
611 }
612
613 if ((cmpres = mpz_cmp_ui(&a->mpz, (unsigned long int)0)) == 0) {
614 /* the base is 0 */
615
616 INCREF(mpz_value_zero);
617 return (object *)mpz_value_zero;
618 }
619 else if (cmpres > 0
620 && mpz_cmp_ui(&a->mpz, (unsigned long int)1) == 0) {
621 /* the base is 1 */
622
623 INCREF(mpz_value_one);
624 return (object *)mpz_value_one;
625 }
626 else if (cmpres < 0
627 && mpz_cmp_si(&a->mpz, (long int)-1) == 0) {
628
629 MP_INT tmpmpz;
630 /* the base is -1: pow(-1, any) == 1,-1 for even,uneven b */
631 /* XXX this code needs to be optimized: what's better?
632 mpz_mmod_ui or mpz_mod_2exp, I choose for the latter
633 for *un*obvious reasons */
634
635 /* is the exponent even? */
636 mpz_init(&tmpmpz);
637
638 /* look to the remainder after a division by (1 << 1) */
639 mpz_mod_2exp(&tmpmpz, &b->mpz, (unsigned long int)1);
640
641 if (mpz_cmp_ui(&tmpmpz, (unsigned int)0) == 0) {
642 mpz_clear(&tmpmpz);
643 INCREF(mpz_value_one);
644 return (object *)mpz_value_one;
645 }
646 mpz_clear(&tmpmpz);
647 INCREF(mpz_value_mone);
648 return (object *)mpz_value_mone;
649 }
650
651#ifdef MPZ_LIB_DOES_CHECKING
652 /* check if it's doable: sizeof(exp) > sizeof(long) &&
653 abs(base) > 1 ?? --> No Way */
654 if (mpz_size(&b->mpz) > 1)
655 return (object *)err_nomem();
656#else /* def MPZ_LIB_DOES_CHECKING */
657 /* wet finger method */
658 if (mpz_cmp_ui(&b->mpz, (unsigned long int)0x10000) >= 0) {
659 err_setstr(ValueError, "mpz.pow outrageous exponent");
660 return NULL;
661 }
662#endif /* def MPZ_LIB_DOES_CHECKING else */
663
664 if ((z = newmpzobject()) == NULL)
665 return NULL;
666
667 mpz_pow_ui(&z->mpz, &a->mpz, mpz_get_ui(&b->mpz));
668
669 return (object *)z;
670} /* mpz_power() */
671
672
673static object *
674mpz_negative(v)
675 mpzobject *v;
676{
677 mpzobject *z;
678
679
680#ifdef MPZ_SPARE_MALLOC
681 if (mpz_cmp_ui(&v->mpz, (unsigned long int)0) == 0) {
682 /* -0 == 0 */
683 INCREF(v);
684 return (object *)v;
685 }
686#endif /* def MPZ_SPARE_MALLOC */
687
688 if ((z = newmpzobject()) == NULL)
689 return NULL;
690
691 mpz_neg(&z->mpz, &v->mpz);
692 return (object *)z;
693} /* mpz_negative() */
694
695
696static object *
697mpz_positive(v)
698 mpzobject *v;
699{
700 INCREF(v);
701 return (object *)v;
702} /* mpz_positive() */
703
704
705static object *
706mpz_absolute(v)
707 mpzobject *v;
708{
709 mpzobject *z;
710
711
712 if (mpz_cmp_ui(&v->mpz, (unsigned long int)0) >= 0) {
713 INCREF(v);
714 return (object *)v;
715 }
716
717 if ((z = newmpzobject()) == NULL)
718 return NULL;
719
720 mpz_neg(&z->mpz, &v->mpz);
721 return (object *)z;
722} /* mpz_absolute() */
723
724static int
725mpz_nonzero(v)
726 mpzobject *v;
727{
728 return mpz_cmp_ui(&v->mpz, (unsigned long int)0) != 0;
729} /* mpz_nonzero() */
730
731static object *
732mpz_invert(v)
733 mpzobject *v;
734{
735 mpzobject *z;
736
737
738 /* I think mpz_com does exactly what needed */
739 if ((z = newmpzobject()) == NULL)
740 return NULL;
741
742 mpz_com(&z->mpz, &v->mpz);
743 return (object *)z;
744} /* mpz_invert() */
745
746static object *
747mpz_lshift(a, b)
748 mpzobject *a;
749 mpzobject *b;
750{
751 int cmpres;
752 mpzobject *z;
753
754
755 if ((cmpres = mpz_cmp_ui(&b->mpz, (unsigned long int)0)) == 0) {
756 /* a << 0 == a */
757 INCREF(a);
758 return (object *)a;
759 }
760
761 if (cmpres < 0) {
762 err_setstr(ValueError, "mpz.<< negative shift count");
763 return NULL;
764 }
765
766#ifdef MPZ_LIB_DOES_CHECKING
767 if (mpz_size(&b->mpz) > 1)
768 return (object *)err_nomem();
769#else /* def MPZ_LIB_DOES_CHECKING */
770 /* wet finger method */
771 if (mpz_cmp_ui(&b->mpz, (unsigned long int)0x10000) >= 0) {
772 err_setstr(ValueError, "mpz.<< outrageous shift count");
773 return NULL;
774 }
775#endif /* def MPZ_LIB_DOES_CHECKING else */
776
777 if ((z = newmpzobject()) == NULL)
778 return NULL;
779
780 mpz_mul_2exp(&z->mpz, &a->mpz, mpz_get_ui(&b->mpz));
781 return (object *)z;
782} /* mpz_lshift() */
783
784static object *
785mpz_rshift(a, b)
786 mpzobject *a;
787 mpzobject *b;
788{
789 int cmpres;
790 mpzobject *z;
791
792
793 if ((cmpres = mpz_cmp_ui(&b->mpz, (unsigned long int)0)) == 0) {
794 /* a >> 0 == a */
795 INCREF(a);
796 return (object *)a;
797 }
798
799 if (cmpres < 0) {
800 err_setstr(ValueError, "mpz.>> negative shift count");
801 return NULL;
802 }
803
804 if (mpz_size(&b->mpz) > 1)
805 return (object *)err_nomem();
806
807 if ((z = newmpzobject()) == NULL)
808 return NULL;
809
810 mpz_div_2exp(&z->mpz, &a->mpz, mpz_get_ui(&b->mpz));
811 return (object *)z;
812} /* mpz_rshift() */
813
814static object *
815mpz_andfunc(a, b)
816 mpzobject *a;
817 mpzobject *b;
818{
819 mpzobject *z;
820
821
822 if ((z = newmpzobject()) == NULL)
823 return NULL;
824
825 mpz_and(&z->mpz, &a->mpz, &b->mpz);
826 return (object *)z;
827} /* mpz_andfunc() */
828
829/* hack Hack HAck HACk HACK, XXX this code is dead slow */
830void
831mpz_xor(res, op1, op2)
832 MP_INT *res;
833 const MP_INT *op1;
834 const MP_INT *op2;
835{
836 MP_INT tmpmpz;
837
838 mpz_init(&tmpmpz);
839
840 mpz_and(res, op1, op2);
841 mpz_com(&tmpmpz, res);
842 mpz_ior(res, op1, op2);
843 mpz_and(res, res, &tmpmpz);
844
845 mpz_clear(&tmpmpz);
846} /* mpz_xor() HACK */
847
848static object *
849mpz_xorfunc(a, b)
850 mpzobject *a;
851 mpzobject *b;
852{
853 mpzobject *z;
854
855
856 if ((z = newmpzobject()) == NULL)
857 return NULL;
858
859 mpz_xor(&z->mpz, &a->mpz, &b->mpz);
860 return (object *)z;
861} /* mpz_xorfunc() */
862
863static object *
864mpz_orfunc(a, b)
865 mpzobject *a;
866 mpzobject *b;
867{
868 mpzobject *z;
869
870
871 if ((z = newmpzobject()) == NULL)
872 return NULL;
873
874 mpz_ior(&z->mpz, &a->mpz, &b->mpz);
875 return (object *)z;
876} /* mpz_orfunc() */
877
878/* MPZ initialisation */
879
880#include "longintrepr.h"
881
882static object *
883MPZ_mpz(self, args)
884 object *self;
885 object *args;
886{
887 mpzobject *mpzp;
888 object *objp;
889
890
891#ifdef MPZ_DEBUG
892 fputs("MPZ_mpz() called...\n", stderr);
893#endif /* def MPZ_DEBUG */
894
895 if (!getargs(args, "O", &objp))
896 return NULL;
897
898 /* at least we know it's some object */
899 /* note DON't DECREF args NEITHER objp */
900
901 if (is_intobject(objp)) {
902 long lval;
903
904 if (!getargs(objp, "l", &lval))
905 return NULL;
906
907 if (lval == (long)0) {
908 INCREF(mpz_value_zero);
909 mpzp = mpz_value_zero;
910 }
911 else if (lval == (long)1) {
912 INCREF(mpz_value_one);
913 mpzp = mpz_value_one;
914 }
915 else if ((mpzp = newmpzobject()) == NULL)
916 return NULL;
917 else mpz_set_si(&mpzp->mpz, lval);
918 }
919 else if (is_longobject(objp)) {
920 MP_INT mplongdigit;
921 int i;
922 unsigned char isnegative;
923
924
925 if ((mpzp = newmpzobject()) == NULL)
926 return NULL;
927
928 mpz_set_si(&mpzp->mpz, 0L);
929 mpz_init(&mplongdigit);
930
931 /* how we're gonna handle this? */
932 if (isnegative = ((i = ((longobject *)objp)->ob_size) < 0) )
933 i = -i;
934
935 while (i--) {
936 mpz_set_ui(&mplongdigit,
937 (unsigned long)
938 ((longobject *)objp)->ob_digit[i]);
939 mpz_mul_2exp(&mplongdigit,&mplongdigit,
940 (unsigned long int)i * SHIFT);
941 mpz_ior(&mpzp->mpz, &mpzp->mpz, &mplongdigit);
942 }
943
944 if (isnegative)
945 mpz_neg(&mpzp->mpz, &mpzp->mpz);
946
947 /* get rid of allocation for tmp variable */
948 mpz_clear(&mplongdigit);
949 }
950 else if (is_stringobject(objp)) {
951 char *cp;
952 int len;
953 MP_INT mplongdigit;
954
955 if (!getargs(objp, "s#", &cp, &len))
956 return NULL;
957
958 if ((mpzp = newmpzobject()) == NULL)
959 return NULL;
960
961 mpz_set_si(&mpzp->mpz, 0L);
962 mpz_init(&mplongdigit);
963
964 /* let's do it the same way as with the long conversion:
965 without thinking how it can be faster (-: :-) */
966
967 cp += len;
968 while (len--) {
969 mpz_set_ui(&mplongdigit, (unsigned long)*--cp );
970 mpz_mul_2exp(&mplongdigit,&mplongdigit,
971 (unsigned long int)len * 8);
972 mpz_ior(&mpzp->mpz, &mpzp->mpz, &mplongdigit);
973 }
974
975 /* get rid of allocation for tmp variable */
976 mpz_clear(&mplongdigit);
977 }
978 else if (is_mpzobject(objp)) {
979 INCREF(objp);
980 mpzp = (mpzobject *)objp;
981 }
982 else {
983 err_setstr(TypeError,
984 "mpz.mpz() expects integer, long, string or mpz object argument");
985 return NULL;
986 }
987
988
989#ifdef MPZ_DEBUG
990 fputs("MPZ_mpz: created mpz=", stderr);
991 mpz_out_str(stderr, 10, &mpzp->mpz);
992 putc('\n', stderr);
993#endif /* def MPZ_DEBUG */
994 return (object *)mpzp;
995} /* MPZ_mpz() */
996
997static mpzobject *
998mpz_mpzcoerce(z)
999 object *z;
1000{
1001 /* shortcut: 9 out of 10 times the type is already ok */
1002 if (is_mpzobject(z)) {
1003 INCREF(z);
1004 return (mpzobject *)z; /* coercion succeeded */
1005 }
1006
1007 /* what types do we accept?: intobjects and longobjects */
1008 if (is_intobject(z) || is_longobject(z))
1009 return (mpzobject *)MPZ_mpz((object *)NULL, z);
1010
1011 err_setstr(TypeError, "number coercion (to mpzobject) failed");
1012 return NULL;
1013} /* mpz_mpzcoerce() */
1014
1015static object *
1016MPZ_powm(self, args)
1017 object *self;
1018 object *args;
1019{
1020 object *base, *exp, *mod;
1021 mpzobject *mpzbase = NULL, *mpzexp = NULL, *mpzmod = NULL;
1022 mpzobject *z;
1023 int tstres;
1024
1025
1026 if (!getargs(args, "(OOO)", &base, &exp, &mod))
1027 return NULL;
1028
1029 if ((mpzbase = mpz_mpzcoerce(base)) == NULL
1030 || (mpzexp = mpz_mpzcoerce(exp)) == NULL
1031 || (mpzmod = mpz_mpzcoerce(mod)) == NULL
1032 || (z = newmpzobject()) == NULL) {
1033 XDECREF(mpzbase);
1034 XDECREF(mpzexp);
1035 XDECREF(mpzmod);
1036 return NULL;
1037 }
1038
1039 if ((tstres=mpz_cmp_ui(&mpzexp->mpz, (unsigned long int)0)) == 0) {
1040 INCREF(mpz_value_one);
1041 return (object *)mpz_value_one;
1042 }
1043
1044 if (tstres < 0) {
1045 MP_INT absexp;
1046 /* negative exp */
1047
1048 mpz_init_set(&absexp, &mpzexp->mpz);
1049 mpz_abs(&absexp, &absexp);
1050 mpz_powm(&z->mpz, &mpzbase->mpz, &absexp, &mpzmod->mpz);
1051
1052 mpz_divm(&z->mpz, &mpz_value_one->mpz, &z->mpz, &mpzmod->mpz);
1053
1054 mpz_clear(&absexp);
1055 }
1056 else {
1057 mpz_powm(&z->mpz, &mpzbase->mpz, &mpzexp->mpz, &mpzmod->mpz);
1058 }
1059
1060 DECREF(mpzbase);
1061 DECREF(mpzexp);
1062 DECREF(mpzmod);
1063
1064 return (object *)z;
1065} /* MPZ_powm() */
1066
1067
1068static object *
1069MPZ_gcd(self, args)
1070 object *self;
1071 object *args;
1072{
1073 object *op1, *op2;
1074 mpzobject *mpzop1 = NULL, *mpzop2 = NULL;
1075 mpzobject *z;
1076
1077
1078 if (!getargs(args, "(OO)", &op1, &op2))
1079 return NULL;
1080
1081 if ((mpzop1 = mpz_mpzcoerce(op1)) == NULL
1082 || (mpzop2 = mpz_mpzcoerce(op2)) == NULL
1083 || (z = newmpzobject()) == NULL) {
1084 XDECREF(mpzop1);
1085 XDECREF(mpzop2);
1086 return NULL;
1087 }
1088
1089 /* ok, we have three mpzobjects, and an initialised result holder */
1090 mpz_gcd(&z->mpz, &mpzop1->mpz, &mpzop2->mpz);
1091
1092 DECREF(mpzop1);
1093 DECREF(mpzop2);
1094
1095 return (object *)z;
1096} /* MPZ_gcd() */
1097
1098
1099static object *
1100MPZ_gcdext(self, args)
1101 object *self;
1102 object *args;
1103{
1104 object *op1, *op2, *z = NULL;
1105 mpzobject *mpzop1 = NULL, *mpzop2 = NULL;
1106 mpzobject *g = NULL, *s = NULL, *t = NULL;
1107
1108
1109 if (!getargs(args, "(OO)", &op1, &op2))
1110 return NULL;
1111
1112 if ((mpzop1 = mpz_mpzcoerce(op1)) == NULL
1113 || (mpzop2 = mpz_mpzcoerce(op2)) == NULL
1114 || (z = newtupleobject(3)) == NULL
1115 || (g = newmpzobject()) == NULL
1116 || (s = newmpzobject()) == NULL
1117 || (t = newmpzobject()) == NULL) {
1118 XDECREF(mpzop1);
1119 XDECREF(mpzop2);
1120 XDECREF(z);
1121 XDECREF(g);
1122 XDECREF(s);
1123 /*XDECREF(t);*/
1124 return NULL;
1125 }
1126
1127 mpz_gcdext(&g->mpz, &s->mpz, &t->mpz, &mpzop1->mpz, &mpzop2->mpz);
1128
1129 DECREF(mpzop1);
1130 DECREF(mpzop2);
1131
1132 (void)settupleitem(z, 0, (object *)g);
1133 (void)settupleitem(z, 1, (object *)s);
1134 (void)settupleitem(z, 2, (object *)t);
1135
1136 return (object *)z;
1137} /* MPZ_gcdext() */
1138
1139
1140static object *
1141MPZ_sqrt(self, args)
1142 object *self;
1143 object *args;
1144{
1145 object *op;
1146 mpzobject *mpzop = NULL;
1147 mpzobject *z;
1148
1149
1150 if (!getargs(args, "O", &op))
1151 return NULL;
1152
1153 if ((mpzop = mpz_mpzcoerce(op)) == NULL
1154 || (z = newmpzobject()) == NULL) {
1155 XDECREF(mpzop);
1156 return NULL;
1157 }
1158
1159 mpz_sqrt(&z->mpz, &mpzop->mpz);
1160
1161 DECREF(mpzop);
1162
1163 return (object *)z;
1164} /* MPZ_sqrt() */
1165
1166
1167static object *
1168MPZ_sqrtrem(self, args)
1169 object *self;
1170 object *args;
1171{
1172 object *op, *z = NULL;
1173 mpzobject *mpzop = NULL;
1174 mpzobject *root = NULL, *rem = NULL;
1175
1176
1177 if (!getargs(args, "O", &op))
1178 return NULL;
1179
1180 if ((mpzop = mpz_mpzcoerce(op)) == NULL
1181 || (z = newtupleobject(2)) == NULL
1182 || (root = newmpzobject()) == NULL
1183 || (rem = newmpzobject()) == NULL) {
1184 XDECREF(mpzop);
1185 XDECREF(z);
1186 XDECREF(root);
1187 /*XDECREF(rem);*/
1188 return NULL;
1189 }
1190
1191 mpz_sqrtrem(&root->mpz, &rem->mpz, &mpzop->mpz);
1192
1193 DECREF(mpzop);
1194
1195 (void)settupleitem(z, 0, (object *)root);
1196 (void)settupleitem(z, 1, (object *)rem);
1197
1198 return (object *)z;
1199} /* MPZ_sqrtrem() */
1200
1201
1202void
1203#if __STDC__
1204mpz_divm(MP_INT *res, const MP_INT *num, const MP_INT *den, const MP_INT *mod)
1205#else
1206mpz_divm(res, num, den, mod)
1207 MP_INT *res;
1208 const MP_INT *num;
1209 const MP_INT *den;
1210 const MP_INT *mod;
1211#endif
1212{
1213 MP_INT s0, s1, q, r, x, d0, d1;
1214
1215 mpz_init_set(&s0, num);
1216 mpz_init_set_ui(&s1, 0);
1217 mpz_init(&q);
1218 mpz_init(&r);
1219 mpz_init(&x);
1220 mpz_init_set(&d0, den);
1221 mpz_init_set(&d1, mod);
1222
1223 while (d1.size != 0) {
1224 mpz_divmod(&q, &r, &d0, &d1);
1225 mpz_set(&d0, &d1);
1226 mpz_set(&d1, &r);
1227
1228 mpz_mul(&x, &s1, &q);
1229 mpz_sub(&x, &s0, &x);
1230 mpz_set(&s0, &s1);
1231 mpz_set(&s1, &x);
1232 }
1233
1234 if (d0.size != 1 || d0.d[0] != 1)
1235 res->size = 0; /* trouble: the gcd != 1; set s to zero */
1236 else {
1237#ifdef MPZ_MDIV_BUG
1238 /* watch out here! first check the signs, and then perform
1239 the mpz_mod() since mod could point to res */
1240 if ((s0.size < 0) != (mod->size < 0)) {
1241 mpz_mod(res, &s0, mod);
1242
1243 if (res->size)
1244 mpz_add(res, res, mod);
1245 }
1246 else
1247 mpz_mod(res, &s0, mod);
1248
1249#else /* def MPZ_MDIV_BUG */
1250 mpz_mmod(res, &s0, mod);
1251#endif /* def MPZ_MDIV_BUG else */
1252 }
1253
1254 mpz_clear(&s0);
1255 mpz_clear(&s1);
1256 mpz_clear(&q);
1257 mpz_clear(&r);
1258 mpz_clear(&x);
1259 mpz_clear(&d0);
1260 mpz_clear(&d1);
1261} /* mpz_divm() */
1262
1263
1264static object *
1265MPZ_divm(self, args)
1266 object *self;
1267 object *args;
1268{
1269 object *num, *den, *mod;
1270 mpzobject *mpznum, *mpzden, *mpzmod = NULL;
1271 mpzobject *z = NULL;
1272
1273
1274 if (!getargs(args, "(OOO)", &num, &den, &mod))
1275 return NULL;
1276
1277 if ((mpznum = mpz_mpzcoerce(num)) == NULL
1278 || (mpzden = mpz_mpzcoerce(den)) == NULL
1279 || (mpzmod = mpz_mpzcoerce(mod)) == NULL
1280 || (z = newmpzobject()) == NULL ) {
1281 XDECREF(mpznum);
1282 XDECREF(mpzden);
1283 XDECREF(mpzmod);
1284 return NULL;
1285 }
1286
1287 mpz_divm(&z->mpz, &mpznum->mpz, &mpzden->mpz, &mpzmod->mpz);
1288
1289 DECREF(mpznum);
1290 DECREF(mpzden);
1291 DECREF(mpzmod);
1292
1293 if (mpz_cmp_ui(&z->mpz, (unsigned long int)0) == 0) {
1294 DECREF(z);
1295 err_setstr(ValueError, "gcd(den, mod) != 1 or num == 0");
1296 return NULL;
1297 }
1298
1299 return (object *)z;
1300} /* MPZ_divm() */
1301
1302
1303/* MPZ methods-as-attributes */
1304#ifdef MPZ_CONVERSIONS_AS_METHODS
1305static object *
1306mpz_int(self, args)
1307 mpzobject *self;
1308 object *args;
1309#else /* def MPZ_CONVERSIONS_AS_METHODS */
1310static object *
1311mpz_int(self)
1312 mpzobject *self;
1313#endif /* def MPZ_CONVERSIONS_AS_METHODS else */
1314{
1315 long sli;
1316
1317
1318#ifdef MPZ_CONVERSIONS_AS_METHODS
1319 if (!getnoarg(args))
1320 return NULL;
1321#endif /* def MPZ_CONVERSIONS_AS_METHODS */
1322
1323 if (mpz_size(&self->mpz) > 1
1324 || (sli = (long)mpz_get_ui(&self->mpz)) < (long)0 ) {
1325 err_setstr(ValueError, "mpz.int() arg too long to convert");
1326 return NULL;
1327 }
1328
1329 if (mpz_cmp_ui(&self->mpz, (unsigned long)0) < 0)
1330 sli = -sli;
1331
1332 return newintobject(sli);
1333} /* mpz_int() */
1334
1335static object *
1336#ifdef MPZ_CONVERSIONS_AS_METHODS
1337mpz_long(self, args)
1338 mpzobject *self;
1339 object *args;
1340#else /* def MPZ_CONVERSIONS_AS_METHODS */
1341mpz_long(self)
1342 mpzobject *self;
1343#endif /* def MPZ_CONVERSIONS_AS_METHODS else */
1344{
1345 int i, isnegative;
1346 unsigned long int uli;
1347 longobject *longobjp;
1348 int ldcount;
1349 int bitpointer, newbitpointer;
1350 MP_INT mpzscratch;
1351
1352
1353#ifdef MPZ_CONVERSIONS_AS_METHODS
1354 if (!getnoarg(args))
1355 return NULL;
1356#endif /* def MPZ_CONVERSIONS_AS_METHODS */
1357
1358 /* determine length of python-long to be allocated */
1359 if ((longobjp = alloclongobject(i = (int)
1360 ((mpz_size(&self->mpz) * BITS_PER_MP_LIMB
1361 + SHIFT - 1) /
1362 SHIFT))) == NULL)
1363 return NULL;
1364
1365 /* determine sign, and copy self to scratch var */
1366 mpz_init_set(&mpzscratch, &self->mpz);
1367 if (isnegative = (mpz_cmp_ui(&self->mpz, (unsigned long int)0) < 0))
1368 mpz_neg(&mpzscratch, &mpzscratch);
1369
1370 /* let those bits come, let those bits go,
1371 e.g. dismantle mpzscratch, build longobject */
1372
1373 bitpointer = 0; /* the number of valid bits in stock */
1374 newbitpointer = 0;
1375 ldcount = 0; /* the python-long limb counter */
1376 uli = (unsigned long int)0;
1377 while (i--) {
1378 longobjp->ob_digit[ldcount] = uli & MASK;
1379
1380 /* check if we've had enough bits for this digit */
1381 if (bitpointer < SHIFT) {
1382 uli = mpz_get_ui(&mpzscratch);
1383 longobjp->ob_digit[ldcount] |=
1384 (uli << bitpointer) & MASK;
1385 uli >>= SHIFT-bitpointer;
1386 bitpointer += BITS_PER_MP_LIMB;
1387 mpz_div_2exp(&mpzscratch, &mpzscratch,
1388 BITS_PER_MP_LIMB);
1389 }
1390 else
1391 uli >>= SHIFT;
1392 bitpointer -= SHIFT;
1393 ldcount++;
1394 }
1395
1396 assert(mpz_cmp_ui(&mpzscratch, (unsigned long int)0) == 0);
1397 mpz_clear(&mpzscratch);
1398 assert(ldcount <= longobjp->ob_size);
1399
1400 /* long_normalize() is file-static */
1401 /* longobjp = long_normalize(longobjp); */
1402 while (ldcount > 0 && longobjp->ob_digit[ldcount-1] == 0)
1403 ldcount--;
1404 longobjp->ob_size = ldcount;
1405
1406
1407 if (isnegative)
1408 longobjp->ob_size = -longobjp->ob_size;
1409
1410 return (object *)longobjp;
1411
1412} /* mpz_long() */
1413
1414
1415/* I would have avoided pow() anyways, so ... */
1416static const double multiplier = 256.0 * 256.0 * 256.0 * 256.0;
1417
1418#ifdef MPZ_CONVERSIONS_AS_METHODS
1419static object *
1420mpz_float(self, args)
1421 mpzobject *self;
1422 object *args;
1423#else /* def MPZ_CONVERSIONS_AS_METHODS */
1424static object *
1425mpz_float(self)
1426 mpzobject *self;
1427#endif /* def MPZ_CONVERSIONS_AS_METHODS else */
1428{
1429 int i, isnegative;
1430 double x;
1431 double mulstate;
1432 MP_INT mpzscratch;
1433
1434
1435#ifdef MPZ_CONVERSIONS_AS_METHODS
1436 if (!getnoarg(args))
1437 return NULL;
1438#endif /* def MPZ_CONVERSIONS_AS_METHODS */
1439
1440 i = (int)mpz_size(&self->mpz);
1441
1442 /* determine sign, and copy abs(self) to scratch var */
1443 if (isnegative = (mpz_cmp_ui(&self->mpz, (unsigned long int)0) < 0)) {
1444 mpz_init(&mpzscratch);
1445 mpz_neg(&mpzscratch, &self->mpz);
1446 }
1447 else
1448 mpz_init_set(&mpzscratch, &self->mpz);
1449
1450 /* let those bits come, let those bits go,
1451 e.g. dismantle mpzscratch, build floatobject */
1452
1453 x = 0.0;
1454 mulstate = 1.0;
1455 while (i--) {
1456 x += mulstate * mpz_get_ui(&mpzscratch);
1457 mulstate *= multiplier;
1458 mpz_div_2exp(&mpzscratch, &mpzscratch, BITS_PER_MP_LIMB);
1459 }
1460
1461 assert(mpz_cmp_ui(&mpzscratch, (unsigned long int)0) == 0);
1462 mpz_clear(&mpzscratch);
1463
1464 if (isnegative)
1465 x = -x;
1466
1467 return newfloatobject(x);
1468
1469} /* mpz_float() */
1470
1471#ifdef MPZ_CONVERSIONS_AS_METHODS
1472static object *
1473mpz_hex(self, args)
1474 mpzobject *self;
1475 object *args;
1476#else /* def MPZ_CONVERSIONS_AS_METHODS */
1477static object *
1478mpz_hex(self)
1479 mpzobject *self;
1480#endif /* def MPZ_CONVERSIONS_AS_METHODS else */
1481{
1482#ifdef MPZ_CONVERSIONS_AS_METHODS
1483 if (!getnoarg(args))
1484 return NULL;
1485#endif /* def MPZ_CONVERSIONS_AS_METHODS */
1486
1487 return mpz_format(self, 16, (unsigned char)1);
1488} /* mpz_hex() */
1489
1490#ifdef MPZ_CONVERSIONS_AS_METHODS
1491static object *
1492mpz_oct(self, args)
1493 mpzobject *self;
1494 object *args;
1495#else /* def MPZ_CONVERSIONS_AS_METHODS */
1496static object *
1497mpz_oct(self)
1498 mpzobject *self;
1499#endif /* def MPZ_CONVERSIONS_AS_METHODS else */
1500{
1501#ifdef MPZ_CONVERSIONS_AS_METHODS
1502 if (!getnoarg(args))
1503 return NULL;
1504#endif /* def MPZ_CONVERSIONS_AS_METHODS */
1505
1506 return mpz_format(self, 8, (unsigned char)1);
1507} /* mpz_oct() */
1508
1509static object *
1510mpz_binary(self, args)
1511 mpzobject *self;
1512 object *args;
1513{
1514 int size;
1515 stringobject *strobjp;
1516 char *cp;
1517 MP_INT mp;
1518 unsigned long ldigit;
1519
1520 if (!getnoarg(args))
1521 return NULL;
1522
1523 if (mpz_cmp_ui(&self->mpz, (unsigned long int)0) < 0) {
1524 err_setstr(ValueError, "mpz.binary() arg must be >= 0");
1525 return NULL;
1526 }
1527
1528 mpz_init_set(&mp, &self->mpz);
1529 size = (int)mpz_size(&mp);
1530
1531 if ((strobjp = (stringobject *)
1532 newsizedstringobject((char *)0,
1533 size * sizeof (unsigned long int))) == NULL)
1534 return NULL;
1535
1536 /* get the beginning of the string memory and start copying things */
1537 cp = GETSTRINGVALUE(strobjp);
1538
1539 /* this has been programmed using a (fairly) decent lib-i/f it could
1540 be must faster if we looked into the GMP lib */
1541 while (size--) {
1542 ldigit = mpz_get_ui(&mp);
1543 mpz_div_2exp(&mp, &mp, BITS_PER_MP_LIMB);
1544 *cp++ = (unsigned char)(ldigit & 0xFF);
1545 *cp++ = (unsigned char)((ldigit >>= 8) & 0xFF);
1546 *cp++ = (unsigned char)((ldigit >>= 8) & 0xFF);
1547 *cp++ = (unsigned char)((ldigit >>= 8) & 0xFF);
1548 }
1549
1550 while (strobjp->ob_size && !*--cp)
1551 strobjp->ob_size--;
1552
1553 return (object *)strobjp;
1554} /* mpz_binary() */
1555
1556
1557static struct methodlist mpz_methods[] = {
1558#ifdef MPZ_CONVERSIONS_AS_METHODS
1559 {"int", mpz_int},
1560 {"long", mpz_long},
1561 {"float", mpz_float},
1562 {"hex", mpz_hex},
1563 {"oct", mpz_oct},
1564#endif /* def MPZ_CONVERSIONS_AS_METHODS */
1565 {"binary", mpz_binary},
1566 {NULL, NULL} /* sentinel */
1567};
1568
1569static object *
1570mpz_getattr(self, name)
1571 mpzobject *self;
1572 char *name;
1573{
1574 return findmethod(mpz_methods, (object *)self, name);
1575} /* mpz_getattr() */
1576
1577
1578static int
1579mpz_coerce(pv, pw)
1580 object **pv;
1581 object **pw;
1582{
1583 object *z;
1584
1585#ifdef MPZ_DEBUG
1586 fputs("mpz_coerce() called...\n", stderr);
1587#endif /* def MPZ_DEBUG */
1588
1589 assert(is_mpzobject(*pv));
1590
1591 /* always convert other arg to mpz value, except for floats */
1592 if (!is_floatobject(*pw)) {
1593 if ((z = (object *)mpz_mpzcoerce(*pw)) == NULL)
1594 return -1; /* -1: an error always has been set */
1595
1596 INCREF(*pv);
1597 *pw = z;
1598 }
1599 else {
1600 if ((z = mpz_float(*pv, NULL)) == NULL)
1601 return -1;
1602
1603 INCREF(*pw);
1604 *pv = z;
1605 }
1606 return 0; /* coercion succeeded */
1607
1608} /* mpz_coerce() */
1609
1610
1611static object *
1612mpz_repr(v)
1613 object *v;
1614{
1615 return mpz_format(v, 10, (unsigned char)1);
1616} /* mpz_repr() */
1617
1618
1619
1620#define UF (object* (*) FPROTO((object *))) /* Unary function */
1621#define BF (object* (*) FPROTO((object *, object *))) /* Binary function */
1622#define IF (int (*) FPROTO((object *))) /* Int function */
1623
1624static number_methods mpz_as_number = {
1625 BF mpz_addition, /*nb_add*/
1626 BF mpz_substract, /*nb_subtract*/
1627 BF mpz_multiply, /*nb_multiply*/
1628 BF mpz_divide, /*nb_divide*/
1629 BF mpz_remainder, /*nb_remainder*/
1630 BF mpz_div_and_mod, /*nb_divmod*/
1631 BF mpz_power, /*nb_power*/
1632 UF mpz_negative, /*nb_negative*/
1633 UF mpz_positive, /*tp_positive*/
1634 UF mpz_absolute, /*tp_absolute*/
1635 IF mpz_nonzero, /*tp_nonzero*/
1636 UF mpz_invert, /*nb_invert*/
1637 BF mpz_lshift, /*nb_lshift*/
1638 BF mpz_rshift, /*nb_rshift*/
1639 BF mpz_andfunc, /*nb_and*/
1640 BF mpz_xorfunc, /*nb_xor*/
1641 BF mpz_orfunc, /*nb_or*/
1642 (int (*) FPROTO((object **, object **)))
1643 mpz_coerce, /*nb_coerce*/
1644#ifndef MPZ_CONVERSIONS_AS_METHODS
1645 UF mpz_int, /*nb_int*/
1646 UF mpz_long, /*nb_long*/
1647 UF mpz_float, /*nb_float*/
1648 UF mpz_oct, /*nb_oct*/
1649 UF mpz_hex, /*nb_hex*/
1650#endif /* ndef MPZ_CONVERSIONS_AS_METHODS */
1651};
1652
1653static typeobject MPZtype = {
1654 OB_HEAD_INIT(&Typetype)
1655 0, /*ob_size*/
1656 "mpz", /*tp_name*/
1657 sizeof(mpzobject), /*tp_size*/
1658 0, /*tp_itemsize*/
1659 /* methods */
1660 mpz_dealloc, /*tp_dealloc*/
1661 mpz_print, /*tp_print*/
1662 mpz_getattr, /*tp_getattr*/
1663 0, /*tp_setattr*/
1664 mpz_compare, /*tp_compare*/
1665 mpz_repr, /*tp_repr*/
1666 &mpz_as_number, /*tp_as_number*/
1667};
1668
1669/* List of functions exported by this module */
1670
1671static struct methodlist mpz_functions[] = {
1672#if 0
1673 {initialiser_name, MPZ_mpz},
1674#else /* 0 */
1675 /* until guido ``fixes'' struct methodlist */
1676 {(char *)initialiser_name, MPZ_mpz},
1677#endif /* 0 else */
1678 {"powm", MPZ_powm},
1679 {"gcd", MPZ_gcd},
1680 {"gcdext", MPZ_gcdext},
1681 {"sqrt", MPZ_sqrt},
1682 {"sqrtrem", MPZ_sqrtrem},
1683 {"divm", MPZ_divm},
1684 {NULL, NULL} /* Sentinel */
1685};
1686
1687
1688/* #define MP_TEST_ALLOC */
1689
1690#ifdef MP_TEST_ALLOC
1691#define MP_TEST_SIZE 4
1692static const char mp_test_magic[MP_TEST_SIZE] = {'\xAA','\xAA','\xAA','\xAA'};
1693static mp_test_error( location )
1694 int *location;
1695{
1696 /* assumptions: *alloc returns address dividable by 4,
1697 mpz_* routines allocate in chunks dividable by four */
1698 fprintf(stderr, "MP_TEST_ERROR: location holds 0x%08d\n", *location );
1699 fatal("MP_TEST_ERROR");
1700} /* static mp_test_error() */
1701#define MP_EXTRA_ALLOC(size) ((size) + MP_TEST_SIZE)
1702#define MP_SET_TEST(basep,size) (void)memcpy( ((char *)(basep))+(size), mp_test_magic, MP_TEST_SIZE)
1703#define MP_DO_TEST(basep,size) if ( !memcmp( ((char *)(basep))+(size), mp_test_magic, MP_TEST_SIZE ) ) \
1704 ; \
1705 else \
1706 mp_test_error((int *)((char *)(basep) + size))
1707#else /* def MP_TEST_ALLOC */
1708#define MP_EXTRA_ALLOC(size) (size)
1709#define MP_SET_TEST(basep,size)
1710#define MP_DO_TEST(basep,size)
1711#endif /* def MP_TEST_ALLOC else */
1712
1713void *mp_allocate( alloc_size )
1714 size_t alloc_size;
1715{
1716 void *res;
1717
1718#ifdef MPZ_DEBUG
1719 fprintf(stderr, "mp_allocate : size %ld\n",
1720 alloc_size);
1721#endif /* def MPZ_DEBUG */
1722
1723 if ( (res = malloc(MP_EXTRA_ALLOC(alloc_size))) == NULL )
1724 fatal("mp_allocate failure");
1725
1726#ifdef MPZ_DEBUG
1727 fprintf(stderr, "mp_allocate : address 0x%08x\n", res);
1728#endif /* def MPZ_DEBUG */
1729
1730 MP_SET_TEST(res,alloc_size);
1731
1732 return res;
1733} /* mp_allocate() */
1734
1735
1736void *mp_reallocate( ptr, old_size, new_size )
1737 void *ptr;
1738 size_t old_size;
1739 size_t new_size;
1740{
1741 void *res;
1742
1743#ifdef MPZ_DEBUG
1744 fprintf(stderr, "mp_reallocate: old address 0x%08x, old size %ld\n",
1745 ptr, old_size);
1746#endif /* def MPZ_DEBUG */
1747
1748 MP_DO_TEST(ptr, old_size);
1749
1750 if ( (res = realloc(ptr, MP_EXTRA_ALLOC(new_size))) == NULL )
1751 fatal("mp_reallocate failure");
1752
1753#ifdef MPZ_DEBUG
1754 fprintf(stderr, "mp_reallocate: new address 0x%08x, new size %ld\n",
1755 res, new_size);
1756#endif /* def MPZ_DEBUG */
1757
1758 MP_SET_TEST(res, new_size);
1759
1760 return res;
1761} /* mp_reallocate() */
1762
1763
1764void mp_free( ptr, size )
1765 void *ptr;
1766 size_t size;
1767{
1768
1769#ifdef MPZ_DEBUG
1770 fprintf(stderr, "mp_free : old address 0x%08x, old size %ld\n",
1771 ptr, size);
1772#endif /* def MPZ_DEBUG */
1773
1774 MP_DO_TEST(ptr, size);
1775 free(ptr);
1776} /* mp_free() */
1777
1778
1779
1780/* Initialize this module. */
1781
1782void
1783initmpz()
1784{
1785#ifdef MPZ_DEBUG
1786 fputs( "initmpz() called...\n", stderr );
1787#endif /* def MPZ_DEBUG */
1788
1789 mp_set_memory_functions( mp_allocate, mp_reallocate, mp_free );
1790 (void)initmodule("mpz", mpz_functions);
1791
1792 /* create some frequently used constants */
1793 if ((mpz_value_zero = newmpzobject()) == NULL)
1794 fatal("initmpz: can't initialize mpz contstants");
1795 mpz_set_ui(&mpz_value_zero->mpz, (unsigned long int)0);
1796
1797 if ((mpz_value_one = newmpzobject()) == NULL)
1798 fatal("initmpz: can't initialize mpz contstants");
1799 mpz_set_ui(&mpz_value_one->mpz, (unsigned long int)1);
1800
1801 if ((mpz_value_mone = newmpzobject()) == NULL)
1802 fatal("initmpz: can't initialize mpz contstants");
1803 mpz_set_si(&mpz_value_mone->mpz, (long)-1);
1804
1805} /* initmpz() */
1806#ifdef MAKEDUMMYINT
1807int _mpz_dummy_int; /* XXX otherwise, we're .bss-less (DYNLOAD->Jack?) */
1808#endif /* def MAKEDUMMYINT */