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