chiark - git - mdw - catacomb-python/blob - mp.c

   1 /* -*-c-*-
   2  *
   3  * Multiprecision arithmetic
   4  *
   5  * (c) 2004 Straylight/Edgeware
   6  */
   7
   8 /*----- Licensing notice --------------------------------------------------*
   9  *
  10  * This file is part of the Python interface to Catacomb.
  11  *
  12  * Catacomb/Python is free software; you can redistribute it and/or modify
  13  * it under the terms of the GNU General Public License as published by
  14  * the Free Software Foundation; either version 2 of the License, or
  15  * (at your option) any later version.
  16  *
  17  * Catacomb/Python is distributed in the hope that it will be useful,
  18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  20  * GNU General Public License for more details.
  21  *
  22  * You should have received a copy of the GNU General Public License
  23  * along with Catacomb/Python; if not, write to the Free Software Foundation,
  24  * Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
  25  */
  26
  27 /*----- Header files ------------------------------------------------------*/
  28
  29 #include "catacomb-python.h"
  30
  31 /*----- General utilities -------------------------------------------------*/
  32
  33 PyTypeObject *mp_pytype = 0;
  34 PyTypeObject *gf_pytype = 0;
  35
  36 mp *mp_frompylong(PyObject *obj)
  37 {
  38   unsigned long bits;
  39   PyLongObject *l = (PyLongObject *)obj;
  40   int sz;
  41   size_t w;
  42   mpd r = 0;
  43   int b = 0;
  44   int i;
  45   mp *x;
  46   mpw *p;
  47
  48   sz = l->ob_size;
  49   if (sz < 0) sz = -sz;
  50   assert(MPW_BITS >= SHIFT);
  51   bits = (unsigned long)sz * SHIFT;
  52   w = (bits + MPW_BITS - 1)/MPW_BITS;
  53   x = mp_new(w, l->ob_size < 0 ? MP_NEG : 0);
  54   p = x->v;
  55   for (i = 0; i < sz; i++) {
  56     r |= (mpd)l->ob_digit[i] << b;
  57     b += SHIFT;
  58     while (b >= MPW_BITS) {
  59       *p++ = MPW(r);
  60       r >>= MPW_BITS;
  61       b -= MPW_BITS;
  62     }
  63   }
  64   while (r) {
  65     *p++ = MPW(r);
  66     r >>= MPW_BITS;
  67   }
  68   x->vl = p;
  69   MP_SHRINK(x);
  70   return (x);
  71 }
  72
  73 PyObject *mp_topylong(mp *x)
  74 {
  75   unsigned long bits = mp_bits(x);
  76   int sz = (bits + SHIFT - 1)/SHIFT;
  77   PyLongObject *l = _PyLong_New(sz);
  78   mpd r = 0;
  79   int b = 0;
  80   mpw *p = x->v;
  81   int i = 0;
  82
  83   assert(MPW_BITS >= SHIFT);
  84   while (i < sz && p < x->vl) {
  85     r |= (mpd)*p++ << b;
  86     b += MPW_BITS;
  87     while (i < sz && b >= SHIFT) {
  88       l->ob_digit[i++] = r & MASK;
  89       r >>= SHIFT;
  90       b -= SHIFT;
  91     }
  92   }
  93   while (i < sz && r) {
  94     l->ob_digit[i++] = r & MASK;
  95     r >>= SHIFT;
  96   }
  97   l->ob_size = (x->f & MP_NEG) ? -sz : sz;
  98   return ((PyObject *)l);
  99 }
 100
 101 mp *mp_frompyobject(PyObject *o, int radix)
 102 {
 103   mp *x;
 104
 105   if (PyString_Check(o)) {
 106     mptext_stringctx sc;
 107     mp *x;
 108     sc.buf = PyString_AS_STRING(o);
 109     sc.lim = sc.buf + PyString_GET_SIZE(o);
 110     x = mp_read(MP_NEW, radix, &mptext_stringops, &sc);
 111     if (!x) return (0);
 112     if (sc.buf < sc.lim) { MP_DROP(x); return (0); }
 113     return (x);
 114   }
 115   if ((x = tomp(o)) != 0)
 116     return (x);
 117   return (0);
 118 }
 119
 120 PyObject *mp_topystring(mp *x, int radix, const char *xpre,
 121                         const char *pre, const char *post)
 122 {
 123   int len = mptext_len(x, radix) + 1;
 124   mptext_stringctx sc;
 125   PyObject *o;
 126   size_t xprelen = xpre ? strlen(xpre) : 0;
 127   size_t prelen = pre ? strlen(pre) : 0;
 128   size_t postlen = post ? strlen(post) : 0;
 129   char *p;
 130   MP_COPY(x);
 131   o = PyString_FromStringAndSize(0, len + 1 + xprelen + prelen + postlen);
 132   p = PyString_AS_STRING(o);
 133   sc.buf = p;
 134   if (xpre) { memcpy(sc.buf, xpre, xprelen); sc.buf += xprelen; }
 135   if (MP_NEGP(x)) { *sc.buf++ = '-'; x = mp_neg(x, x); }
 136   if (pre) { memcpy(sc.buf, pre, prelen); sc.buf += prelen; }
 137   sc.lim = sc.buf + len;
 138   mp_write(x, radix, &mptext_stringops, &sc);
 139   if (post) { memcpy(sc.buf, post, postlen); sc.buf += postlen; }
 140   MP_DROP(x);
 141   _PyString_Resize(&o, sc.buf - p);
 142   return (o);
 143 }
 144
 145 static int good_radix_p(int r, int readp)
 146 {
 147   return ((r >= -255 && r <= -2) ||
 148           (readp && r == 0) ||
 149           (r >= 2 && r <= 62));
 150 }
 151
 152 PyObject *mp_pywrap(mp *x)
 153 {
 154   mp_pyobj *z = PyObject_New(mp_pyobj, mp_pytype);
 155   z->x = x;
 156   return ((PyObject *)z);
 157 }
 158
 159 PyObject *gf_pywrap(mp *x)
 160 {
 161   mp_pyobj *z = PyObject_New(mp_pyobj, gf_pytype);
 162   z->x = x;
 163   return ((PyObject *)z);
 164 }
 165
 166 int mp_tolong_checked(mp *x, long *l, int must)
 167 {
 168   static mp *longmin = 0, *longmax = 0;
 169   int rc = -1;
 170
 171   if (!longmax) {
 172     longmin = mp_fromlong(MP_NEW, LONG_MIN);
 173     longmax = mp_fromlong(MP_NEW, LONG_MAX);
 174   }
 175   if (MP_CMP(x, <, longmin) || MP_CMP(x, >, longmax)) {
 176     if (must) VALERR("mp out of range for int");
 177     else goto end;
 178   }
 179   *l = mp_tolong(x);
 180   rc = 0;
 181 end:
 182   return (rc);
 183 }
 184
 185 /*----- Python interface --------------------------------------------------*/
 186
 187 static void mp_pydealloc(PyObject *o)
 188 {
 189   MP_DROP(MP_X(o));
 190   FREEOBJ(o);
 191 }
 192
 193 static PyObject *mp_pyrepr(PyObject *o)
 194   { return mp_topystring(MP_X(o), 10, "MP(", 0, "L)"); }
 195
 196 static PyObject *mp_pystr(PyObject *o)
 197   { return mp_topystring(MP_X(o), 10, 0, 0, 0); }
 198
 199 mp *tomp(PyObject *o)
 200 {
 201   PyObject *l;
 202   mp *x;
 203
 204   if (!o)
 205     return (0);
 206   else if (MP_PYCHECK(o) || GF_PYCHECK(o))
 207     return (MP_COPY(MP_X(o)));
 208   else if (FE_PYCHECK(o))
 209     return (F_OUT(FE_F(o), MP_NEW, FE_X(o)));
 210   else if (PFILT_PYCHECK(o))
 211     return (MP_COPY(PFILT_F(o)->m));
 212   else if (ECPT_PYCHECK(o)) {
 213     ec p = EC_INIT;
 214     getecptout(&p, o);
 215     x = MP_COPY(p.x);
 216     EC_DESTROY(&p);
 217     return (x);
 218   } else if (GE_PYCHECK(o)) {
 219     if ((x = G_TOINT(GE_G(o), MP_NEW, GE_X(o))) == 0)
 220       return (0);
 221     return (x);
 222   } else if (PyInt_Check(o))
 223     return (mp_fromlong(MP_NEW, PyInt_AS_LONG(o)));
 224   else if ((l = PyNumber_Long(o)) != 0) {
 225     x = mp_frompylong(l);
 226     Py_DECREF(l);
 227     return (x);
 228   } else {
 229     PyErr_Clear();
 230     return (0);
 231   }
 232 }
 233
 234 mp *getmp(PyObject *o)
 235 {
 236   mp *x = 0;
 237   if (!o) return (0);
 238   if ((x = tomp(o)) == 0) {
 239     PyErr_Format(PyExc_TypeError, "can't convert %.100s to mp",
 240                  o->ob_type->tp_name);
 241   }
 242   return (x);
 243 }
 244
 245 int convmp(PyObject *o, void *p)
 246 {
 247   mp *x;
 248   if ((x = getmp(o)) == 0) return (0);
 249   *(mp **)p = x;
 250   return (1);
 251 }
 252
 253 mp *getgf(PyObject *o)
 254 {
 255   mp *x = 0;
 256   if (!o) return (0);
 257   if ((x = tomp(o)) == 0) {
 258     PyErr_Format(PyExc_TypeError, "can't convert %.100s to gf",
 259                  o->ob_type->tp_name);
 260   }
 261   return (x);
 262 }
 263
 264 int convgf(PyObject *o, void *p)
 265 {
 266   mp *x;
 267   if ((x = getgf(o)) == 0) return (0);
 268   *(mp **)p = x;
 269   return (1);
 270 }
 271
 272 static mp *implicitmp(PyObject *o)
 273 {
 274   if (!o ||
 275       GF_PYCHECK(o) ||
 276       ECPT_PYCHECK(o) ||
 277       FE_PYCHECK(o) ||
 278       GE_PYCHECK(o))
 279     return (0);
 280   return (tomp(o));
 281 }
 282
 283 static mp *implicitgf(PyObject *o)
 284 {
 285   if (!o ||
 286       MP_PYCHECK(o) ||
 287       ECPT_PYCHECK(o) ||
 288       FE_PYCHECK(o) ||
 289       GE_PYCHECK(o))
 290     return (0);
 291   return (tomp(o));
 292 }
 293
 294 static int mpbinop(PyObject *x, PyObject *y, mp **xx, mp **yy)
 295 {
 296   if ((*xx = implicitmp(x)) == 0)
 297     return (-1);
 298   if ((*yy = implicitmp(y)) == 0) {
 299     MP_DROP(*xx);
 300     return (-1);
 301   }
 302   return (0);
 303 }
 304
 305 static int gfbinop(PyObject *x, PyObject *y, mp **xx, mp **yy)
 306 {
 307   if ((*xx = implicitgf(x)) == 0)
 308     return (-1);
 309   if ((*yy = implicitgf(y)) == 0) {
 310     MP_DROP(*xx);
 311     return (-1);
 312   }
 313   return (0);
 314 }
 315
 316 #define gf_and mp_and
 317 #define gf_or mp_or
 318 #define gf_xor mp_xor
 319 #define BINOP(pre, name)                                                \
 320   static PyObject *pre##_py##name(PyObject *x, PyObject *y) {           \
 321     mp *xx, *yy, *zz;                                                   \
 322     if (pre##binop(x, y, &xx, &yy)) RETURN_NOTIMPL;                     \
 323     zz = pre##_##name(MP_NEW, xx, yy);                                  \
 324     MP_DROP(xx); MP_DROP(yy);                                           \
 325     return (pre##_pywrap(zz));                                          \
 326   }
 327 BINOP(mp, add)
 328 BINOP(mp, sub)
 329 BINOP(mp, mul)
 330 BINOP(mp, and2c)
 331 BINOP(mp, or2c)
 332 BINOP(mp, xor2c)
 333 BINOP(gf, add)
 334 BINOP(gf, sub)
 335 BINOP(gf, mul)
 336 BINOP(gf, and)
 337 BINOP(gf, or)
 338 BINOP(gf, xor)
 339 #undef BINOP
 340
 341 static mp *mp_abs(mp *d, mp *x)
 342 {
 343   if (MP_NEGP(x))
 344     return (mp_neg(d, x));
 345   MP_COPY(x);
 346   if (d) MP_DROP(d);
 347   return (x);
 348 }
 349
 350 #define UNOP(pre, name)                                                 \
 351   static PyObject *pre##_py##name(PyObject *x)                          \
 352     { return mp_pywrap(pre##_##name(MP_NEW, MP_X(x))); }
 353 UNOP(mp, neg)
 354 UNOP(mp, abs)
 355 UNOP(mp, not2c)
 356 #undef UNOP
 357
 358 static PyObject *mp_pyid(PyObject *x) { RETURN_OBJ(x); }
 359
 360 #define gf_lsr mp_lsr
 361 #define SHIFTOP(pre, name, rname)                                       \
 362   static PyObject *pre##_py##name(PyObject *x, PyObject *y) {           \
 363     mp *xx, *yy;                                                        \
 364     PyObject *z = 0;                                                    \
 365     long n;                                                             \
 366     if (pre##binop(x, y, &xx, &yy)) RETURN_NOTIMPL;                     \
 367     if (mp_tolong_checked(yy, &n, 1)) goto end;                         \
 368     if (n < 0)                                                          \
 369       z = pre##_pywrap(mp_##rname(MP_NEW, xx, -n));                     \
 370     else                                                                \
 371       z = pre##_pywrap(mp_##name(MP_NEW, xx, n));                       \
 372   end:                                                                  \
 373     MP_DROP(xx); MP_DROP(yy);                                           \
 374     return (z);                                                         \
 375   }
 376 SHIFTOP(mp, lsl2c, lsr2c)
 377 SHIFTOP(mp, lsr2c, lsl2c)
 378 SHIFTOP(gf, lsl, lsr)
 379 SHIFTOP(gf, lsr, lsl)
 380 #undef SHIFTOP
 381
 382 #define DIVOP(pre, name, qq, rr, gather)                                \
 383   static PyObject *pre##_py##name(PyObject *x, PyObject *y) {           \
 384     mp *xx, *yy;                                                        \
 385     PyObject *z = 0;                                                    \
 386     INIT_##qq(q) INIT_##rr(r)                                           \
 387     if (pre##binop(x, y, &xx, &yy)) RETURN_NOTIMPL;                     \
 388     if (MP_ZEROP(yy))                                                   \
 389       ZDIVERR("division by zero");                                      \
 390     pre##_div(ARG_##qq(q), ARG_##rr(r), xx, yy);                        \
 391     z = gather;                                                         \
 392   end:                                                                  \
 393     MP_DROP(xx); MP_DROP(yy);                                           \
 394     return (z);                                                         \
 395   }
 396 #define INIT_YES(p) mp *p = MP_NEW;
 397 #define INIT_NO(p)
 398 #define ARG_YES(p) &p
 399 #define ARG_NO(p) 0
 400 DIVOP(mp, divmod, YES, YES,
 401       Py_BuildValue("(NN)", mp_pywrap(q), mp_pywrap(r)))
 402 DIVOP(mp, div, YES, NO, mp_pywrap(q))
 403 DIVOP(mp, mod, NO, YES, mp_pywrap(r))
 404 DIVOP(gf, divmod, YES, YES,
 405       Py_BuildValue("(NN)", gf_pywrap(q), gf_pywrap(r)))
 406 DIVOP(gf, div, YES, NO, gf_pywrap(q))
 407 DIVOP(gf, mod, NO, YES, gf_pywrap(r))
 408 #undef INIT_YES
 409 #undef INIT_NO
 410 #undef ARG_YES
 411 #undef ARG_NO
 412 #undef DIVOP
 413
 414 static mp *mp_modinv_checked(mp *d, mp *x, mp *p)
 415 {
 416   mp *g = MP_NEW;
 417   mp_gcd(&g, 0, &d, p, x);
 418   if (!MP_EQ(g, MP_ONE)) {
 419     MP_DROP(g); MP_DROP(d);
 420     PyErr_SetString(PyExc_ZeroDivisionError, "no modular inverse");
 421     return (0);
 422   }
 423   MP_DROP(g); return (d);
 424 }
 425
 426 static PyObject *mp_pyexp(PyObject *x, PyObject *y, PyObject *z)
 427 {
 428   mp *xx = 0, *yy = 0, *zz = 0;
 429   mp *r = 0;
 430   PyObject *rc = 0;
 431
 432   if ((xx = implicitmp(x)) == 0 || (yy = implicitmp(y)) == 0 ||
 433       (z && z != Py_None && (zz = tomp(z)) == 0)) {
 434     mp_drop(xx); mp_drop(yy); mp_drop(zz);
 435     RETURN_NOTIMPL;
 436   }
 437   if (!z || z == Py_None) {
 438     if (MP_NEGP(yy)) VALERR("negative exponent");
 439     r = mp_exp(MP_NEW, xx, yy);
 440   } else {
 441     if (!MP_POSP(zz)) VALERR("modulus must be positive");
 442     if (MP_NEGP(yy)) {
 443       if ((xx = mp_modinv_checked(xx, xx, zz)) == 0) goto end;
 444       yy = mp_neg(yy, yy);
 445     }
 446     if (MP_ODDP(zz)) {
 447       mpmont mm;
 448       mpmont_create(&mm, zz);
 449       r = mpmont_exp(&mm, MP_NEW, xx, yy);
 450       mpmont_destroy(&mm);
 451     } else {
 452       mpbarrett mb;
 453       mpbarrett_create(&mb, zz);
 454       r = mpbarrett_exp(&mb, MP_NEW, xx, yy);
 455       mpbarrett_destroy(&mb);
 456     }
 457   }
 458   rc = mp_pywrap(r);
 459 end:
 460   mp_drop(xx); mp_drop(yy); mp_drop(zz);
 461   return (rc);
 462 }
 463
 464 static mp *gf_modinv_checked(mp *d, mp *x, mp *p)
 465 {
 466   mp *g = MP_NEW;
 467   gf_gcd(&g, 0, &d, p, x);
 468   if (!MP_EQ(g, MP_ONE)) {
 469     MP_DROP(g); MP_DROP(d);
 470     PyErr_SetString(PyExc_ZeroDivisionError, "no modular inverse");
 471     return (0);
 472   }
 473   MP_DROP(g); return (d);
 474 }
 475
 476 #define BASEOP(name, radix, pre)                                        \
 477   static PyObject *mp_py##name(PyObject *x)                             \
 478     { return mp_topystring(MP_X(x), radix, 0, pre, 0); }
 479 BASEOP(oct, 8, "0");
 480 BASEOP(hex, 16, "0x");
 481 #undef BASEOP
 482
 483 static int mp_pynonzerop(PyObject *x) { return !MP_ZEROP(MP_X(x)); }
 484
 485 static PyObject *mp_pyint(PyObject *x)
 486 {
 487   long l;
 488   if (!mp_tolong_checked(MP_X(x), &l, 0)) return (PyInt_FromLong(l));
 489   else return mp_topylong(MP_X(x));
 490 }
 491 static PyObject *mp_pylong(PyObject *x)
 492   { return (mp_topylong(MP_X(x))); }
 493 static PyObject *mp_pyfloat(PyObject *x)
 494 {
 495   PyObject *l = mp_topylong(MP_X(x));
 496   double f = PyLong_AsDouble(l);
 497   Py_DECREF(l);
 498   return (PyFloat_FromDouble(f));
 499 }
 500
 501 #define COERCE(pre, PRE)                                                \
 502   static int pre##_pycoerce(PyObject **x, PyObject **y)                 \
 503   {                                                                     \
 504     mp *z;                                                              \
 505                                                                         \
 506     if (PRE##_PYCHECK(*y)) {                                            \
 507       Py_INCREF(*x); Py_INCREF(*y);                                     \
 508       return (0);                                                       \
 509     }                                                                   \
 510     if ((z = tomp(*y)) != 0) {                                          \
 511       Py_INCREF(*x);                                                    \
 512       *y = pre##_pywrap(z);                                             \
 513       return (0);                                                       \
 514     }                                                                   \
 515     return (1);                                                         \
 516   }
 517 COERCE(mp, MP)
 518 COERCE(gf, GF)
 519 #undef COERCE
 520
 521 static int mp_pycompare(PyObject *x, PyObject *y)
 522   { return mp_cmp(MP_X(x), MP_X(y)); }
 523
 524 static PyObject *mp_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
 525 {
 526   PyObject *x;
 527   mp *z;
 528   mp_pyobj *zz = 0;
 529   int radix = 0;
 530   char *kwlist[] = { "x", "radix", 0 };
 531
 532   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O|i:new", kwlist, &x, &radix))
 533     goto end;
 534   if (MP_PYCHECK(x)) RETURN_OBJ(x);
 535   if (!good_radix_p(radix, 1)) VALERR("bad radix");
 536   if ((z = mp_frompyobject(x, radix)) == 0) {
 537     PyErr_Format(PyExc_TypeError, "can't convert %.100s to mp",
 538                  x->ob_type->tp_name);
 539     goto end;
 540   }
 541   zz = (mp_pyobj *)ty->tp_alloc(ty, 0);
 542   zz->x = z;
 543 end:
 544   return ((PyObject *)zz);
 545 }
 546
 547 long mphash(mp *x)
 548 {
 549   PyObject *l = mp_topylong(x);
 550   long h = PyObject_Hash(l);
 551   Py_DECREF(l); return (h);
 552 }
 553
 554 static long mp_pyhash(PyObject *me) { return (mphash(MP_X(me))); }
 555
 556 static PyObject *mpmeth_jacobi(PyObject *me, PyObject *arg)
 557 {
 558   mp *y = 0;
 559   PyObject *z = 0;
 560
 561   if (!PyArg_ParseTuple(arg, "O&:jacobi", convmp, &y)) goto end;
 562   z = PyInt_FromLong(mp_jacobi(y, MP_X(me)));
 563 end:
 564   if (y) MP_DROP(y);
 565   return (z);
 566 }
 567
 568 #define BITOP(pre, name, c)                                             \
 569   static PyObject *pre##meth_##name(PyObject *me, PyObject *arg)        \
 570   {                                                                     \
 571     unsigned long i;                                                    \
 572     if (!PyArg_ParseTuple(arg, "O&:" #name, convulong, &i)) return (0); \
 573     return (pre##_pywrap(mp_##name##c(MP_NEW, MP_X(me), i)));           \
 574   }
 575 BITOP(mp, setbit, 2c);
 576 BITOP(mp, clearbit, 2c);
 577 BITOP(gf, setbit, );
 578 BITOP(gf, clearbit, );
 579 #undef BITOP
 580
 581 static PyObject *mpmeth_testbit(PyObject *me, PyObject *arg)
 582 {
 583   unsigned long i;
 584   if (!PyArg_ParseTuple(arg, "O&:testbit", convulong, &i)) return (0);
 585   return (getbool(mp_testbit2c(MP_X(me), i)));
 586 }
 587
 588 static PyObject *gfmeth_testbit(PyObject *me, PyObject *arg)
 589 {
 590   unsigned long i;
 591   if (!PyArg_ParseTuple(arg, "O&:testbit", convulong, &i)) return (0);
 592   return (getbool(mp_testbit(MP_X(me), i)));
 593 }
 594
 595 static PyObject *mpmeth_odd(PyObject *me, PyObject *arg)
 596 {
 597   mp *t;
 598   size_t s;
 599
 600   if (!PyArg_ParseTuple(arg, ":odd")) return (0);
 601   t = mp_odd(MP_NEW, MP_X(me), &s);
 602   return (Py_BuildValue("(lN)", (long)s, mp_pywrap(t)));
 603 }
 604
 605 static PyObject *mpmeth_sqr(PyObject *me, PyObject *arg)
 606 {
 607   if (!PyArg_ParseTuple(arg, ":sqr")) return (0);
 608   return (mp_pywrap(mp_sqr(MP_NEW, MP_X(me))));
 609 }
 610
 611 static PyObject *mpmeth_sqrt(PyObject *me, PyObject *arg)
 612 {
 613   if (!PyArg_ParseTuple(arg, ":sqrt")) return (0);
 614   if (MP_NEGP(MP_X(me))) VALERR("negative root");
 615   return (mp_pywrap(mp_sqrt(MP_NEW, MP_X(me))));
 616 end:
 617   return (0);
 618 }
 619
 620 static PyObject *mpmeth_gcd(PyObject *me, PyObject *arg)
 621 {
 622   mp *y = 0, *zz = MP_NEW;
 623   PyObject *z = 0;
 624
 625   if (!PyArg_ParseTuple(arg, "O&:gcd", convmp, &y)) goto end;
 626   mp_gcd(&zz, 0, 0, MP_X(me), y);
 627   z = mp_pywrap(zz);
 628 end:
 629   if (y) MP_DROP(y);
 630   return (z);
 631 }
 632
 633 static PyObject *mpmeth_gcdx(PyObject *me, PyObject *arg)
 634 {
 635   PyObject *z = 0;
 636   mp *yy = 0, *zz = MP_NEW, *uu = MP_NEW, *vv = MP_NEW;
 637
 638   if (!PyArg_ParseTuple(arg, "O&:gcdx", convmp, &yy)) goto end;
 639   mp_gcd(&zz, &uu, &vv, MP_X(me), yy);
 640   z = Py_BuildValue("(NNN)",
 641                     mp_pywrap(zz), mp_pywrap(uu), mp_pywrap(vv));
 642 end:
 643   if (yy) MP_DROP(yy);
 644   return (z);
 645 }
 646
 647 static PyObject *mpmeth_modinv(PyObject *me, PyObject *arg)
 648 {
 649   PyObject *z = 0;
 650   mp *yy = 0, *zz = MP_NEW;
 651
 652   if (!PyArg_ParseTuple(arg, "O&:modinv", convmp, &yy) ||
 653       (zz = mp_modinv_checked(MP_NEW, yy, MP_X(me))) == 0)
 654     goto end;
 655   z = mp_pywrap(zz);
 656 end:
 657   if (yy) MP_DROP(yy);
 658   return (z);
 659 }
 660
 661 static PyObject *mpmeth_tostring(PyObject *me, PyObject *arg, PyObject *kw)
 662 {
 663   int radix = 10;
 664   char *kwlist[] = { "radix", 0 };
 665   if (!PyArg_ParseTupleAndKeywords(arg, kw, "|i:tostring", kwlist, &radix))
 666     goto end;
 667   if (!good_radix_p(radix, 0)) VALERR("bad radix");
 668   return (mp_topystring(MP_X(me), radix, 0, 0, 0));
 669 end:
 670   return (0);
 671 }
 672
 673 static PyObject *mpmeth_modsqrt(PyObject *me, PyObject *arg)
 674 {
 675   PyObject *z = 0;
 676   mp *yy = 0, *zz = MP_NEW;
 677
 678   if (!PyArg_ParseTuple(arg, "O&:modsqrt", convmp, &yy)) goto end;
 679   if ((zz = mp_modsqrt(MP_NEW, yy, MP_X(me))) == 0)
 680     VALERR("no modular square root");
 681   z = mp_pywrap(zz);
 682 end:
 683   if (yy) MP_DROP(yy);
 684   return (z);
 685 }
 686
 687 static PyObject *mpmeth_leastcongruent(PyObject *me, PyObject *arg)
 688 {
 689   mp *z, *b, *m;
 690   PyObject *rc = 0;
 691
 692   if (!PyArg_ParseTuple(arg, "O&O&:leastcongruent", convmp, &b, convmp, &m))
 693     goto end;
 694   z = mp_leastcongruent(MP_NEW, b, MP_X(me), m);
 695   rc = mp_pywrap(z);
 696 end:
 697   return (rc);
 698 }
 699
 700 #define STOREOP(name, c)                                                \
 701   static PyObject *mpmeth_##name(PyObject *me,                          \
 702                                  PyObject *arg, PyObject *kw)           \
 703   {                                                                     \
 704     long len = -1;                                                      \
 705     char *kwlist[] = { "len", 0 };                                      \
 706     PyObject *rc = 0;                                                   \
 707                                                                         \
 708     if (!PyArg_ParseTupleAndKeywords(arg, kw, "|l:" #name,              \
 709                                     kwlist, &len))                      \
 710       goto end;                                                         \
 711     if (len < 0) {                                                      \
 712       len = mp_octets##c(MP_X(me));                                     \
 713       if (!len) len = 1;                                                \
 714     }                                                                   \
 715     rc = bytestring_pywrap(0, len);                                     \
 716     mp_##name(MP_X(me), PyString_AS_STRING(rc), len);                   \
 717   end:                                                                  \
 718     return (rc);                                                        \
 719   }
 720 STOREOP(storel, )
 721 STOREOP(storeb, )
 722 STOREOP(storel2c, 2c)
 723 STOREOP(storeb2c, 2c)
 724 #undef STOREOP
 725
 726 #define BUFOP(ty, pyty)                                                 \
 727   static PyObject *meth__##pyty##_frombuf(PyObject *me, PyObject *arg)  \
 728   {                                                                     \
 729     buf b;                                                              \
 730     char *p;                                                            \
 731     Py_ssize_t sz;                                                      \
 732     PyObject *rc = 0;                                                   \
 733     mp *x;                                                              \
 734                                                                         \
 735     if (!PyArg_ParseTuple(arg, "Os#:frombuf", &me, &p, &sz)) goto end;  \
 736     buf_init(&b, p, sz);                                                \
 737     if ((x = buf_getmp(&b)) == 0) VALERR("malformed data");             \
 738     rc = Py_BuildValue("(NN)", ty##_pywrap(x),                          \
 739                        bytestring_pywrapbuf(&b));                       \
 740   end:                                                                  \
 741     return (rc);                                                        \
 742   }
 743 BUFOP(mp, MP)
 744 BUFOP(gf, GF)
 745 #undef BUFOP
 746
 747 static PyObject *mpmeth_tobuf(PyObject *me, PyObject *arg)
 748 {
 749   buf b;
 750   PyObject *rc;
 751   mp *x;
 752   size_t n;
 753
 754   if (!PyArg_ParseTuple(arg, ":tobuf")) return (0);
 755   x = MP_X(me);
 756   n = mp_octets(x) + 3;
 757   rc = bytestring_pywrap(0, n);
 758   buf_init(&b, PyString_AS_STRING(rc), n);
 759   buf_putmp(&b, x);
 760   assert(BOK(&b));
 761   _PyString_Resize(&rc, BLEN(&b));
 762   return (rc);
 763 }
 764
 765 static PyObject *mpmeth_primep(PyObject *me, PyObject *arg, PyObject *kw)
 766 {
 767   grand *r = &rand_global;
 768   char *kwlist[] = { "rng", 0 };
 769   PyObject *rc = 0;
 770
 771   if (!PyArg_ParseTupleAndKeywords(arg, kw, "|O&", kwlist, convgrand, &r))
 772     goto end;
 773   rc = getbool(pgen_primep(MP_X(me), r));
 774 end:
 775   return (rc);
 776 }
 777
 778 static PyObject *mpget_nbits(PyObject *me, void *hunoz)
 779   { return (PyInt_FromLong(mp_bits(MP_X(me)))); }
 780
 781 static PyObject *mpget_noctets(PyObject *me, void *hunoz)
 782   { return (PyInt_FromLong(mp_octets(MP_X(me)))); }
 783
 784 static PyObject *mpget_noctets2c(PyObject *me, void *hunoz)
 785   { return (PyInt_FromLong(mp_octets2c(MP_X(me)))); }
 786
 787 static PyGetSetDef mp_pygetset[] = {
 788 #define GETSETNAME(op, func) mp##op##_##func
 789   GET   (nbits,         "X.nbits -> bit length of X")
 790   GET   (noctets,       "X.noctets -> octet length of X")
 791   GET   (noctets2c,     "X.noctets2c -> two's complement octet length of X")
 792 #undef GETSETNAME
 793   { 0 }
 794 };
 795
 796 static PyMethodDef mp_pymethods[] = {
 797 #define METHNAME(func) mpmeth_##func
 798   METH  (jacobi,        "X.jacobi(Y) -> Jacobi symbol (Y|X) (NB inversion!)")
 799   METH  (setbit,        "X.setbit(N) -> X with bit N set")
 800   METH  (clearbit,      "X.clearbit(N) -> X with bit N clear")
 801   METH  (testbit,       "X.testbit(N) -> true/false if bit N set/clear in X")
 802   METH  (odd,           "X.odd() -> S, T where X = 2^S T with T odd")
 803   METH  (sqr,           "X.sqr() -> X^2")
 804   METH  (sqrt,          "X.sqrt() -> largest integer <= sqrt(X)")
 805   METH  (gcd,           "X.gcd(Y) -> gcd(X, Y)")
 806   METH  (gcdx,
 807          "X.gcdx(Y) -> (gcd(X, Y), U, V) with X U + Y V = gcd(X, Y)")
 808   METH  (modinv,        "X.modinv(Y) -> multiplicative inverse of Y mod X")
 809   METH  (modsqrt,       "X.modsqrt(Y) -> square root of Y mod X, if X prime")
 810   METH  (leastcongruent,
 811          "X.leastcongruent(B, M) -> smallest Z >= B with Z == X (mod M)")
 812   KWMETH(primep,        "X.primep([rng = rand]) -> true/false if X is prime")
 813   KWMETH(tostring,      "X.tostring([radix = 10]) -> STR")
 814   KWMETH(storel,        "X.storel([len = -1]) -> little-endian bytes")
 815   KWMETH(storeb,        "X.storeb([len = -1]) -> big-endian bytes")
 816   KWMETH(storel2c,
 817          "X.storel2c([len = -1]) -> little-endian bytes, two's complement")
 818   KWMETH(storeb2c,
 819          "X.storeb2c([len = -1]) -> big-endian bytes, two's complement")
 820   METH  (tobuf,         "X.tobuf() -> buffer format")
 821 #undef METHNAME
 822   { 0 }
 823 };
 824
 825 static PyNumberMethods mp_pynumber = {
 826   mp_pyadd,                             /* @nb_add@ */
 827   mp_pysub,                             /* @nb_subtract@ */
 828   mp_pymul,                             /* @nb_multiply@ */
 829   0,                                    /* @nb_divide@ */
 830   mp_pymod,                             /* @nb_remainder@ */
 831   mp_pydivmod,                          /* @nb_divmod@ */
 832   mp_pyexp,                             /* @nb_power@ */
 833   mp_pyneg,                             /* @nb_negative@ */
 834   mp_pyid,                              /* @nb_positive@ */
 835   mp_pyabs,                             /* @nb_absolute@ */
 836   mp_pynonzerop,                        /* @nb_nonzero@ */
 837   mp_pynot2c,                           /* @nb_invert@ */
 838   mp_pylsl2c,                           /* @nb_lshift@ */
 839   mp_pylsr2c,                           /* @nb_rshift@ */
 840   mp_pyand2c,                           /* @nb_and@ */
 841   mp_pyxor2c,                           /* @nb_xor@ */
 842   mp_pyor2c,                            /* @nb_or@ */
 843   mp_pycoerce,                          /* @nb_coerce@ */
 844   mp_pyint,                             /* @nb_int@ */
 845   mp_pylong,                            /* @nb_long@ */
 846   mp_pyfloat,                           /* @nb_float@ */
 847   mp_pyoct,                             /* @nb_oct@ */
 848   mp_pyhex,                             /* @nb_hex@ */
 849
 850   0,                                    /* @nb_inplace_add@ */
 851   0,                                    /* @nb_inplace_subtract@ */
 852   0,                                    /* @nb_inplace_multiply@ */
 853   0,                                    /* @nb_inplace_divide@ */
 854   0,                                    /* @nb_inplace_remainder@ */
 855   0,                                    /* @nb_inplace_power@ */
 856   0,                                    /* @nb_inplace_lshift@ */
 857   0,                                    /* @nb_inplace_rshift@ */
 858   0,                                    /* @nb_inplace_and@ */
 859   0,                                    /* @nb_inplace_xor@ */
 860   0,                                    /* @nb_inplace_or@ */
 861
 862   mp_pydiv,                             /* @nb_floor_divide@ */
 863   0,                                    /* @nb_true_divide@ */
 864   0,                                    /* @nb_inplace_floor_divide@ */
 865   0,                                    /* @nb_inplace_true_divide@ */
 866 };
 867
 868 static PyTypeObject mp_pytype_skel = {
 869   PyObject_HEAD_INIT(0) 0,              /* Header */
 870   "MP",                                 /* @tp_name@ */
 871   sizeof(mp_pyobj),                     /* @tp_basicsize@ */
 872   0,                                    /* @tp_itemsize@ */
 873
 874   mp_pydealloc,                         /* @tp_dealloc@ */
 875   0,                                    /* @tp_print@ */
 876   0,                                    /* @tp_getattr@ */
 877   0,                                    /* @tp_setattr@ */
 878   mp_pycompare,                         /* @tp_compare@ */
 879   mp_pyrepr,                            /* @tp_repr@ */
 880   &mp_pynumber,                         /* @tp_as_number@ */
 881   0,                                    /* @tp_as_sequence@ */
 882   0,                                    /* @tp_as_mapping@ */
 883   mp_pyhash,                            /* @tp_hash@ */
 884   0,                                    /* @tp_call@ */
 885   mp_pystr,                             /* @tp_str@ */
 886   0,                                    /* @tp_getattro@ */
 887   0,                                    /* @tp_setattro@ */
 888   0,                                    /* @tp_as_buffer@ */
 889   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
 890     Py_TPFLAGS_CHECKTYPES |
 891     Py_TPFLAGS_BASETYPE,
 892
 893   /* @tp_doc@ */
 894 "Multiprecision integers, similar to `long' but more efficient and\n\
 895 versatile.  Support all the standard arithmetic operations, with\n\
 896 implicit conversions from `PrimeFilter', and other objects which\n\
 897 convert to `long'.\n\
 898 \n\
 899 Constructor MP(X, [radix = R]) attempts to convert X to an `MP'.  If\n\
 900 X is a string, it's read in radix-R form, or we look for a prefix\n\
 901 if R = 0.  Other acceptable things are field elements, elliptic curve\n\
 902 points, group elements, Python `int' and `long' objects, and anything\n\
 903 with an integer conversion.\n\
 904 \n\
 905 Notes:\n\
 906 \n\
 907   * Use `//' for division.  MPs don't have `/' division.",
 908
 909   0,                                    /* @tp_traverse@ */
 910   0,                                    /* @tp_clear@ */
 911   0,                                    /* @tp_richcompare@ */
 912   0,                                    /* @tp_weaklistoffset@ */
 913   0,                                    /* @tp_iter@ */
 914   0,                                    /* @tp_iternext@ */
 915   mp_pymethods,                         /* @tp_methods@ */
 916   0,                                    /* @tp_members@ */
 917   mp_pygetset,                          /* @tp_getset@ */
 918   0,                                    /* @tp_base@ */
 919   0,                                    /* @tp_dict@ */
 920   0,                                    /* @tp_descr_get@ */
 921   0,                                    /* @tp_descr_set@ */
 922   0,                                    /* @tp_dictoffset@ */
 923   0,                                    /* @tp_init@ */
 924   PyType_GenericAlloc,                  /* @tp_alloc@ */
 925   mp_pynew,                             /* @tp_new@ */
 926   0,                                    /* @tp_free@ */
 927   0                                     /* @tp_is_gc@ */
 928 };
 929
 930 static PyObject *meth__MP_fromstring(PyObject *me,
 931                                      PyObject *arg, PyObject *kw)
 932 {
 933   int r = 0;
 934   char *p;
 935   Py_ssize_t len;
 936   PyObject *z = 0;
 937   mp *zz;
 938   mptext_stringctx sc;
 939   char *kwlist[] = { "class", "x", "radix", 0 };
 940
 941   if (!PyArg_ParseTupleAndKeywords(arg, kw, "Os#|i:fromstring",
 942                                    kwlist, &me, &p, &len, &r))
 943     goto end;
 944   if (!good_radix_p(r, 1)) VALERR("bad radix");
 945   sc.buf = p; sc.lim = p + len;
 946   if ((zz = mp_read(MP_NEW, r, &mptext_stringops, &sc)) == 0)
 947     VALERR("bad integer");
 948   z = Py_BuildValue("(Ns#)", mp_pywrap(zz),
 949                     sc.buf, (Py_ssize_t)(sc.lim - sc.buf));
 950 end:
 951   return (z);
 952 }
 953
 954 static PyObject *meth__MP_factorial(PyObject *me, PyObject *arg)
 955 {
 956   unsigned long i;
 957   mp *x;
 958   if (!PyArg_ParseTuple(arg, "OO&:factorial", &me, convulong, &i))
 959     return (0);
 960   x = mp_factorial(i);
 961   return mp_pywrap(x);
 962 }
 963
 964 static PyObject *meth__MP_fibonacci(PyObject *me, PyObject *arg)
 965 {
 966   long i;
 967   mp *x;
 968   if (!PyArg_ParseTuple(arg, "Ol:fibonacci", &me, &i)) return (0);
 969   x = mp_fibonacci(i);
 970   return mp_pywrap(x);
 971 }
 972
 973 #define LOADOP(pre, py, name)                                           \
 974   static PyObject *meth__##py##_##name(PyObject *me, PyObject *arg)     \
 975   {                                                                     \
 976     char *p;                                                            \
 977     Py_ssize_t len;                                                     \
 978     if (!PyArg_ParseTuple(arg, "Os#:" #name, &me, &p, &len)) return (0); \
 979     return (pre##_pywrap(mp_##name(MP_NEW, p, len)));                   \
 980   }
 981 LOADOP(mp, MP, loadl)
 982 LOADOP(mp, MP, loadb)
 983 LOADOP(mp, MP, loadl2c)
 984 LOADOP(mp, MP, loadb2c)
 985 LOADOP(gf, GF, loadl)
 986 LOADOP(gf, GF, loadb)
 987 #undef LOADOP
 988
 989 /*----- Products of small integers ----------------------------------------*/
 990
 991 typedef struct mpmul_pyobj {
 992   PyObject_HEAD
 993   int livep;
 994   mpmul mm;
 995 } mpmul_pyobj;
 996
 997 #define MPMUL_LIVEP(o) (((mpmul_pyobj *)(o))->livep)
 998 #define MPMUL_PY(o) (&((mpmul_pyobj *)(o))->mm)
 999
1000 static void mpmul_pydealloc(PyObject *me)
1001 {
1002   if (MPMUL_LIVEP(me))
1003     mp_drop(mpmul_done(MPMUL_PY(me)));
1004   FREEOBJ(me);
1005 }
1006
1007 static PyObject *mmmeth_factor(PyObject *me, PyObject *arg)
1008 {
1009   PyObject *q, *i;
1010   mp *x;
1011
1012   if (!MPMUL_LIVEP(me)) VALERR("MPMul object invalid");
1013   if (PyTuple_Size(arg) != 1)
1014     i = PyObject_GetIter(arg);
1015   else {
1016     if ((q = PyTuple_GetItem(arg, 0)) == 0) goto end;
1017     if ((i = PyObject_GetIter(q)) == 0) {
1018       PyErr_Clear(); /* that's ok */
1019       i = PyObject_GetIter(arg);
1020     }
1021   }
1022   if (!i) goto end;
1023   while ((q = PyIter_Next(i)) != 0) {
1024     x = getmp(q); Py_DECREF(q); if (!x) {
1025       Py_DECREF(i);
1026       goto end;
1027     }
1028     mpmul_add(MPMUL_PY(me), x);
1029     MP_DROP(x);
1030   }
1031   Py_DECREF(i);
1032   RETURN_ME;
1033 end:
1034   return (0);
1035 }
1036
1037 static PyObject *mmmeth_done(PyObject *me, PyObject *arg)
1038 {
1039   mp *x;
1040
1041   if (!PyArg_ParseTuple(arg, ":done")) goto end;
1042   if (!MPMUL_LIVEP(me)) VALERR("MPMul object invalid");
1043   x = mpmul_done(MPMUL_PY(me));
1044   MPMUL_LIVEP(me) = 0;
1045   return (mp_pywrap(x));
1046 end:
1047   return (0);
1048 }
1049
1050 static PyObject *mpmul_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1051 {
1052   mpmul_pyobj *mm;
1053
1054   if (kw) TYERR("keyword arguments not allowed here");
1055   mm = (mpmul_pyobj *)ty->tp_alloc(ty, 0);
1056   mpmul_init(&mm->mm);
1057   mm->livep = 1;
1058   if (mmmeth_factor((PyObject *)mm, arg) == 0) {
1059     Py_DECREF(mm);
1060     goto end;
1061   }
1062   return ((PyObject *)mm);
1063 end:
1064   return (0);
1065 }
1066
1067 static PyObject *mmget_livep(PyObject *me, void *hunoz)
1068   { return (getbool(MPMUL_LIVEP(me))); }
1069
1070 static PyGetSetDef mpmul_pygetset[] = {
1071 #define GETSETNAME(op, name) mm##op##_##name
1072   GET   (livep,                 "MM.livep -> flag: object still valid?")
1073 #undef GETSETNAME
1074   { 0 }
1075 };
1076
1077 static PyMethodDef mpmul_pymethods[] = {
1078 #define METHNAME(name) mmmeth_##name
1079   METH  (factor,                "MM.factor(ITERABLE) or MM.factor(I, ...)")
1080   METH  (done,                  "MM.done() -> PRODUCT")
1081 #undef METHNAME
1082   { 0 }
1083 };
1084
1085 static PyTypeObject *mpmul_pytype, mpmul_pytype_skel = {
1086   PyObject_HEAD_INIT(0) 0,              /* Header */
1087   "MPMul",                              /* @tp_name@ */
1088   sizeof(mpmul_pyobj),                  /* @tp_basicsize@ */
1089   0,                                    /* @tp_itemsize@ */
1090
1091   mpmul_pydealloc,                      /* @tp_dealloc@ */
1092   0,                                    /* @tp_print@ */
1093   0,                                    /* @tp_getattr@ */
1094   0,                                    /* @tp_setattr@ */
1095   0,                                    /* @tp_compare@ */
1096   0,                                    /* @tp_repr@ */
1097   0,                                    /* @tp_as_number@ */
1098   0,                                    /* @tp_as_sequence@ */
1099   0,                                    /* @tp_as_mapping@ */
1100   0,                                    /* @tp_hash@ */
1101   0,                                    /* @tp_call@ */
1102   0,                                    /* @tp_str@ */
1103   0,                                    /* @tp_getattro@ */
1104   0,                                    /* @tp_setattro@ */
1105   0,                                    /* @tp_as_buffer@ */
1106   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1107     Py_TPFLAGS_BASETYPE,
1108
1109   /* @tp_doc@ */
1110 "An object for multiplying many small integers.",
1111
1112   0,                                    /* @tp_traverse@ */
1113   0,                                    /* @tp_clear@ */
1114   0,                                    /* @tp_richcompare@ */
1115   0,                                    /* @tp_weaklistoffset@ */
1116   0,                                    /* @tp_iter@ */
1117   0,                                    /* @tp_iternext@ */
1118   mpmul_pymethods,                      /* @tp_methods@ */
1119   0,                                    /* @tp_members@ */
1120   mpmul_pygetset,                       /* @tp_getset@ */
1121   0,                                    /* @tp_base@ */
1122   0,                                    /* @tp_dict@ */
1123   0,                                    /* @tp_descr_get@ */
1124   0,                                    /* @tp_descr_set@ */
1125   0,                                    /* @tp_dictoffset@ */
1126   0,                                    /* @tp_init@ */
1127   PyType_GenericAlloc,                  /* @tp_alloc@ */
1128   mpmul_pynew,                          /* @tp_new@ */
1129   0,                                    /* @tp_free@ */
1130   0                                     /* @tp_is_gc@ */
1131 };
1132
1133 /*----- Montgomery reduction ----------------------------------------------*/
1134
1135 typedef struct mpmont_pyobj {
1136   PyObject_HEAD
1137   mpmont mm;
1138 } mpmont_pyobj;
1139
1140 #define MPMONT_PY(o) (&((mpmont_pyobj *)(o))->mm)
1141
1142 static PyObject *mmmeth_int(PyObject *me, PyObject *arg)
1143 {
1144   PyObject *z = 0;
1145   mp *yy = 0;
1146   mpmont *mm = MPMONT_PY(me);
1147
1148   if (!PyArg_ParseTuple(arg, "O&:in", convmp, &yy))
1149     goto end;
1150   mp_div(0, &yy, yy, mm->m);
1151   z = mp_pywrap(mpmont_mul(mm, MP_NEW, yy, mm->r2));
1152 end:
1153   if (yy) MP_DROP(yy);
1154   return (z);
1155 }
1156
1157 static PyObject *mmmeth_mul(PyObject *me, PyObject *arg)
1158 {
1159   PyObject *rc = 0;
1160   mp *yy = 0, *zz = 0;
1161
1162   if (!PyArg_ParseTuple(arg, "O&O&:mul", convmp, &yy, convmp, &zz))
1163     goto end;
1164   rc = mp_pywrap(mpmont_mul(MPMONT_PY(me), MP_NEW, yy, zz));
1165 end:
1166   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1167   return (rc);
1168 }
1169
1170 static PyObject *mmmeth_exp(PyObject *me, PyObject *arg)
1171 {
1172   PyObject *rc = 0;
1173   mp *yy = 0, *zz = 0;
1174
1175   if (!PyArg_ParseTuple(arg, "O&O&:exp", convmp, &yy, convmp, &zz))
1176     goto end;
1177   if (MP_NEGP(zz)) {
1178     if ((yy = mp_modinv_checked(yy, yy, MPMONT_PY(me)->m)) == 0) goto end;
1179     zz = mp_neg(zz, zz);
1180   }
1181   rc = mp_pywrap(mpmont_exp(MPMONT_PY(me), MP_NEW, yy, zz));
1182 end:
1183   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1184   return (rc);
1185 }
1186
1187 static PyObject *mmmeth_expr(PyObject *me, PyObject *arg)
1188 {
1189   PyObject *rc = 0;
1190   mp *yy = 0, *zz = 0;
1191
1192   if (!PyArg_ParseTuple(arg, "O&O&:expr", convmp, &yy, convmp, &zz))
1193     goto end;
1194   if (MP_NEGP(zz)) {
1195     yy = mpmont_reduce(MPMONT_PY(me), yy, yy);
1196     if ((yy = mp_modinv_checked(yy, yy, MPMONT_PY(me)->m)) == 0) goto end;
1197     yy = mpmont_mul(MPMONT_PY(me), yy, yy, MPMONT_PY(me)->r2);
1198     zz = mp_neg(zz, zz);
1199   }
1200   rc = mp_pywrap(mpmont_expr(MPMONT_PY(me), MP_NEW, yy, zz));
1201 end:
1202   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1203   return (rc);
1204 }
1205
1206 static PyObject *mm_mexpr_id(PyObject *me)
1207   { return mp_pywrap(MP_COPY(MPMONT_PY(me)->r)); }
1208
1209 static int mm_mexpr_fill(void *p, PyObject *me, PyObject *x, PyObject *y)
1210 {
1211   mp *xx = 0, *yy = 0;
1212   mp_expfactor *f = p;
1213   mpmont *mm = MPMONT_PY(me);
1214
1215   if ((xx = getmp(x)) == 0 || (yy = getmp(y)) == 0)
1216     goto fail;
1217   if (MP_NEGP(yy)) {
1218     xx = mpmont_reduce(mm, xx, xx);
1219     if ((xx = mp_modinv_checked(xx, xx, yy)) == 0)
1220       goto fail;
1221     xx = mpmont_mul(mm, xx, xx, mm->r2);
1222     yy = mp_neg(yy, yy);
1223   }
1224   f->base = xx;
1225   f->exp = yy;
1226   return (0);
1227
1228 fail:
1229   mp_drop(xx); mp_drop(yy);
1230   return (-1);
1231 }
1232
1233 static PyObject *mm_mexpr(PyObject *me, void *v, int n)
1234   { return mp_pywrap(mpmont_mexpr(MPMONT_PY(me), MP_NEW, v, n)); }
1235
1236 static void mp_mexp_drop(void *p)
1237 {
1238   mp_expfactor *f = p;
1239   mp_drop(f->base);
1240   mp_drop(f->exp);
1241 }
1242
1243 static PyObject *mmmeth_mexpr(PyObject *me, PyObject *arg)
1244 {
1245   return mexp_common(me, arg, sizeof(mp_expfactor),
1246                      mm_mexpr_id, mm_mexpr_fill, mm_mexpr, mp_mexp_drop);
1247 }
1248
1249 static PyObject *mp_mexp_id(PyObject *me)
1250   { return mp_pywrap(MP_ONE); }
1251
1252 static int mp_mexp_fill(void *p, PyObject *me, PyObject *x, PyObject *y)
1253 {
1254   mp *xx = 0, *yy = 0;
1255   mp_expfactor *f = p;
1256
1257   if ((xx = getmp(x)) == 0 || (yy = getmp(y)) == 0)
1258     goto fail;
1259   if (MP_NEGP(yy)) {
1260     if ((xx = mp_modinv_checked(xx, xx, yy)) == 0)
1261       goto fail;
1262     yy = mp_neg(yy, yy);
1263   }
1264   f->base = xx;
1265   f->exp = yy;
1266   return (0);
1267
1268 fail:
1269   mp_drop(xx); mp_drop(yy);
1270   return (-1);
1271 }
1272
1273 static PyObject *mm_mexp(PyObject *me, void *v, int n)
1274   { return mp_pywrap(mpmont_mexp(MPMONT_PY(me), MP_NEW, v, n)); }
1275
1276 static PyObject *mmmeth_mexp(PyObject *me, PyObject *arg)
1277 {
1278   return mexp_common(me, arg, sizeof(mp_expfactor),
1279                      mp_mexp_id, mp_mexp_fill, mm_mexp, mp_mexp_drop);
1280 }
1281
1282 #define mmmeth_ext mmmeth_reduce
1283 static PyObject *mmmeth_reduce(PyObject *me, PyObject *arg)
1284 {
1285   PyObject *z = 0;
1286   mp *yy = 0;
1287
1288   if (!PyArg_ParseTuple(arg, "O&", convmp, &yy)) goto end;
1289   z = mp_pywrap(mpmont_reduce(MPMONT_PY(me), MP_NEW, yy));
1290 end:
1291   return (z);
1292 }
1293
1294 static void mpmont_pydealloc(PyObject *me)
1295 {
1296   mpmont_destroy(MPMONT_PY(me));
1297   FREEOBJ(me);
1298 }
1299
1300 static PyObject *mpmont_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1301 {
1302   mpmont_pyobj *mm = 0;
1303   char *kwlist[] = { "m", 0 };
1304   mp *xx = 0;
1305
1306   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convmp, &xx))
1307     goto end;
1308   if (!MP_POSP(xx) || !MP_ODDP(xx)) VALERR("m must be positive and odd");
1309   mm = (mpmont_pyobj *)ty->tp_alloc(ty, 0);
1310   mpmont_create(&mm->mm, xx);
1311 end:
1312   if (xx) MP_DROP(xx);
1313   return ((PyObject *)mm);
1314 }
1315
1316 static PyObject *mmget_m(PyObject *me, void *hunoz)
1317   { return (mp_pywrap(MP_COPY(MPMONT_PY(me)->m))); }
1318
1319 static PyObject *mmget_r(PyObject *me, void *hunoz)
1320   { return (mp_pywrap(MP_COPY(MPMONT_PY(me)->r))); }
1321
1322 static PyObject *mmget_r2(PyObject *me, void *hunoz)
1323   { return (mp_pywrap(MP_COPY(MPMONT_PY(me)->r2))); }
1324
1325 static PyGetSetDef mpmont_pygetset[] = {
1326 #define GETSETNAME(op, name) mm##op##_##name
1327   GET   (m,             "M.m -> modulus for reduction")
1328   GET   (r,             "M.r -> multiplicative identity")
1329   GET   (r2,            "M.r2 -> M.r^2, Montgomerization factor")
1330 #undef GETSETNAME
1331   { 0 }
1332 };
1333
1334 static PyMethodDef mpmont_pymethods[] = {
1335 #define METHNAME(name) mmmeth_##name
1336   METH  (int,           "M.int(X) -> XR")
1337   METH  (mul,           "M.mul(XR, YR) -> ZR where Z = X Y")
1338   METH  (expr,          "M.expr(XR, N) -> ZR where Z = X^N mod M.m")
1339   METH  (mexpr,         "\
1340 M.mexpr([(XR0, N0), (XR1, N1), ...]) = ZR where Z = X0^N0 X1^N1 ... mod M.m\n\
1341 \t(the list may be flattened if this more convenient.)")
1342   METH  (reduce,        "M.reduce(XR) -> X")
1343   METH  (ext,           "M.ext(XR) -> X")
1344   METH  (exp,           "M.exp(X, N) -> X^N mod M.m")
1345   METH  (mexp,          "\
1346 M.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 ... mod M.m\n\
1347 \t(the list may be flattened if this more convenient.)")
1348 #undef METHNAME
1349   { 0 }
1350 };
1351
1352 static PyTypeObject *mpmont_pytype, mpmont_pytype_skel = {
1353   PyObject_HEAD_INIT(0) 0,              /* Header */
1354   "MPMont",                             /* @tp_name@ */
1355   sizeof(mpmont_pyobj),                 /* @tp_basicsize@ */
1356   0,                                    /* @tp_itemsize@ */
1357
1358   mpmont_pydealloc,                     /* @tp_dealloc@ */
1359   0,                                    /* @tp_print@ */
1360   0,                                    /* @tp_getattr@ */
1361   0,                                    /* @tp_setattr@ */
1362   0,                                    /* @tp_compare@ */
1363   0,                                    /* @tp_repr@ */
1364   0,                                    /* @tp_as_number@ */
1365   0,                                    /* @tp_as_sequence@ */
1366   0,                                    /* @tp_as_mapping@ */
1367   0,                                    /* @tp_hash@ */
1368   0,                                    /* @tp_call@ */
1369   0,                                    /* @tp_str@ */
1370   0,                                    /* @tp_getattro@ */
1371   0,                                    /* @tp_setattro@ */
1372   0,                                    /* @tp_as_buffer@ */
1373   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1374     Py_TPFLAGS_BASETYPE,
1375
1376   /* @tp_doc@ */
1377 "A Montgomery reduction context.",
1378
1379   0,                                    /* @tp_traverse@ */
1380   0,                                    /* @tp_clear@ */
1381   0,                                    /* @tp_richcompare@ */
1382   0,                                    /* @tp_weaklistoffset@ */
1383   0,                                    /* @tp_iter@ */
1384   0,                                    /* @tp_iternext@ */
1385   mpmont_pymethods,                     /* @tp_methods@ */
1386   0,                                    /* @tp_members@ */
1387   mpmont_pygetset,                      /* @tp_getset@ */
1388   0,                                    /* @tp_base@ */
1389   0,                                    /* @tp_dict@ */
1390   0,                                    /* @tp_descr_get@ */
1391   0,                                    /* @tp_descr_set@ */
1392   0,                                    /* @tp_dictoffset@ */
1393   0,                                    /* @tp_init@ */
1394   PyType_GenericAlloc,                  /* @tp_alloc@ */
1395   mpmont_pynew,                         /* @tp_new@ */
1396   0,                                    /* @tp_free@ */
1397   0                                     /* @tp_is_gc@ */
1398 };
1399
1400 /*----- Barrett reduction -------------------------------------------------*/
1401
1402 typedef struct mpbarrett_pyobj {
1403   PyObject_HEAD
1404   mpbarrett mb;
1405 } mpbarrett_pyobj;
1406
1407 #define MPBARRETT_PY(o) (&((mpbarrett_pyobj *)(o))->mb)
1408
1409 static PyObject *mbmeth_exp(PyObject *me, PyObject *arg)
1410 {
1411   PyObject *rc = 0;
1412   mp *yy = 0, *zz = 0;
1413
1414   if (!PyArg_ParseTuple(arg, "O&O&:exp", convmp, &yy, convmp, &zz))
1415     goto end;
1416   if (MP_NEGP(zz)) {
1417     if ((yy = mp_modinv_checked(yy, yy, MPBARRETT_PY(me)->m)) == 0) goto end;
1418     zz = mp_neg(zz, zz);
1419   }
1420   rc = mp_pywrap(mpbarrett_exp(MPBARRETT_PY(me), MP_NEW, yy, zz));
1421 end:
1422   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1423   return (rc);
1424 }
1425
1426 static PyObject *mb_mexp(PyObject *me, void *v, int n)
1427   { return mp_pywrap(mpbarrett_mexp(MPBARRETT_PY(me), MP_NEW, v, n)); }
1428
1429 static PyObject *mbmeth_mexp(PyObject *me, PyObject *arg)
1430 {
1431   return mexp_common(me, arg, sizeof(mp_expfactor),
1432                      mp_mexp_id, mp_mexp_fill, mb_mexp, mp_mexp_drop);
1433 }
1434
1435 static PyObject *mbmeth_reduce(PyObject *me, PyObject *arg)
1436 {
1437   PyObject *z = 0;
1438   mp *yy = 0;
1439
1440   if (!PyArg_ParseTuple(arg, "O&:reduce", convmp, &yy))
1441     goto end;
1442   z = mp_pywrap(mpbarrett_reduce(MPBARRETT_PY(me), MP_NEW, yy));
1443 end:
1444   return (z);
1445 }
1446
1447 static void mpbarrett_pydealloc(PyObject *me)
1448 {
1449   mpbarrett_destroy(MPBARRETT_PY(me));
1450   FREEOBJ(me);
1451 }
1452
1453 static PyObject *mpbarrett_pynew(PyTypeObject *ty,
1454                                  PyObject *arg, PyObject *kw)
1455 {
1456   mpbarrett_pyobj *mb = 0;
1457   char *kwlist[] = { "m", 0 };
1458   mp *xx = 0;
1459
1460   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convmp, &xx))
1461     goto end;
1462   if (!MP_POSP(xx)) VALERR("m must be positive");
1463   mb = (mpbarrett_pyobj *)ty->tp_alloc(ty, 0);
1464   mpbarrett_create(&mb->mb, xx);
1465 end:
1466   if (xx) MP_DROP(xx);
1467   return ((PyObject *)mb);
1468 }
1469
1470 static PyObject *mbget_m(PyObject *me, void *hunoz)
1471   { return (mp_pywrap(MP_COPY(MPBARRETT_PY(me)->m))); }
1472
1473 static PyGetSetDef mpbarrett_pygetset[] = {
1474 #define GETSETNAME(op, name) mb##op##_##name
1475   GET   (m,             "B.m -> modulus for reduction")
1476 #undef GETSETNAME
1477   { 0 }
1478 };
1479
1480 static PyMethodDef mpbarrett_pymethods[] = {
1481 #define METHNAME(name) mbmeth_##name
1482   METH  (reduce,        "B.reduce(X) -> X mod B.m")
1483   METH  (exp,           "B.exp(X, N) -> X^N mod B.m")
1484   METH  (mexp,          "\
1485 B.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 ... mod B.m\n\
1486 \t(the list may be flattened if this more convenient.)")
1487 #undef METHNAME
1488   { 0 }
1489 };
1490
1491 static PyTypeObject *mpbarrett_pytype, mpbarrett_pytype_skel = {
1492   PyObject_HEAD_INIT(0) 0,              /* Header */
1493   "MPBarrett",                          /* @tp_name@ */
1494   sizeof(mpbarrett_pyobj),              /* @tp_basicsize@ */
1495   0,                                    /* @tp_itemsize@ */
1496
1497   mpbarrett_pydealloc,                  /* @tp_dealloc@ */
1498   0,                                    /* @tp_print@ */
1499   0,                                    /* @tp_getattr@ */
1500   0,                                    /* @tp_setattr@ */
1501   0,                                    /* @tp_compare@ */
1502   0,                                    /* @tp_repr@ */
1503   0,                                    /* @tp_as_number@ */
1504   0,                                    /* @tp_as_sequence@ */
1505   0,                                    /* @tp_as_mapping@ */
1506   0,                                    /* @tp_hash@ */
1507   0,                                    /* @tp_call@ */
1508   0,                                    /* @tp_str@ */
1509   0,                                    /* @tp_getattro@ */
1510   0,                                    /* @tp_setattro@ */
1511   0,                                    /* @tp_as_buffer@ */
1512   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1513     Py_TPFLAGS_BASETYPE,
1514
1515   /* @tp_doc@ */
1516 "A Barrett reduction context.",
1517
1518   0,                                    /* @tp_traverse@ */
1519   0,                                    /* @tp_clear@ */
1520   0,                                    /* @tp_richcompare@ */
1521   0,                                    /* @tp_weaklistoffset@ */
1522   0,                                    /* @tp_iter@ */
1523   0,                                    /* @tp_iternext@ */
1524   mpbarrett_pymethods,                  /* @tp_methods@ */
1525   0,                                    /* @tp_members@ */
1526   mpbarrett_pygetset,                   /* @tp_getset@ */
1527   0,                                    /* @tp_base@ */
1528   0,                                    /* @tp_dict@ */
1529   0,                                    /* @tp_descr_get@ */
1530   0,                                    /* @tp_descr_set@ */
1531   0,                                    /* @tp_dictoffset@ */
1532   0,                                    /* @tp_init@ */
1533   PyType_GenericAlloc,                  /* @tp_alloc@ */
1534   mpbarrett_pynew,                      /* @tp_new@ */
1535   0,                                    /* @tp_free@ */
1536   0                                     /* @tp_is_gc@ */
1537 };
1538
1539 /*----- Nice prime reduction ----------------------------------------------*/
1540
1541 typedef struct mpreduce_pyobj {
1542   PyObject_HEAD
1543   mpreduce mr;
1544 } mpreduce_pyobj;
1545
1546 #define MPREDUCE_PY(o) (&((mpreduce_pyobj *)(o))->mr)
1547
1548 static PyObject *mrmeth_exp(PyObject *me, PyObject *arg)
1549 {
1550   PyObject *rc = 0;
1551   mp *yy = 0, *zz = 0;
1552
1553   if (!PyArg_ParseTuple(arg, "O&O&:exp", convmp, &yy, convmp, &zz))
1554     goto end;
1555   if (MP_NEGP(zz)) {
1556     if ((yy = mp_modinv_checked(yy, yy, MPREDUCE_PY(me)->p)) == 0) goto end;
1557     zz = mp_neg(zz, zz);
1558   }
1559   rc = mp_pywrap(mpreduce_exp(MPREDUCE_PY(me), MP_NEW, yy, zz));
1560 end:
1561   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1562   return (rc);
1563 }
1564
1565 static PyObject *mrmeth_reduce(PyObject *me, PyObject *arg)
1566 {
1567   PyObject *z = 0;
1568   mp *yy = 0;
1569
1570   if (!PyArg_ParseTuple(arg, "O&:reduce", convmp, &yy)) goto end;
1571   z = mp_pywrap(mpreduce_do(MPREDUCE_PY(me), MP_NEW, yy));
1572 end:
1573   return (z);
1574 }
1575
1576 static void mpreduce_pydealloc(PyObject *me)
1577 {
1578   mpreduce_destroy(MPREDUCE_PY(me));
1579   FREEOBJ(me);
1580 }
1581
1582 static PyObject *mpreduce_pynew(PyTypeObject *ty,
1583                                  PyObject *arg, PyObject *kw)
1584 {
1585   mpreduce_pyobj *mr = 0;
1586   mpreduce r;
1587   char *kwlist[] = { "m", 0 };
1588   mp *xx = 0;
1589
1590   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convmp, &xx))
1591     goto end;
1592   if (!MP_POSP(xx)) VALERR("m must be positive");
1593   if (mpreduce_create(&r, xx)) VALERR("bad modulus (must be 2^k - ...)");
1594   mr = (mpreduce_pyobj *)ty->tp_alloc(ty, 0);
1595   mr->mr = r;
1596 end:
1597   if (xx) MP_DROP(xx);
1598   return ((PyObject *)mr);
1599 }
1600
1601 static PyObject *mrget_m(PyObject *me, void *hunoz)
1602   { return (mp_pywrap(MP_COPY(MPREDUCE_PY(me)->p))); }
1603
1604 static PyGetSetDef mpreduce_pygetset[] = {
1605 #define GETSETNAME(op, name) mr##op##_##name
1606   GET   (m,             "R.m -> modulus for reduction")
1607 #undef GETSETNAME
1608   { 0 }
1609 };
1610
1611 static PyMethodDef mpreduce_pymethods[] = {
1612 #define METHNAME(name) mrmeth_##name
1613   METH  (reduce,        "R.reduce(X) -> X mod B.m")
1614   METH  (exp,           "R.exp(X, N) -> X^N mod B.m")
1615 #undef METHNAME
1616   { 0 }
1617 };
1618
1619 static PyTypeObject *mpreduce_pytype, mpreduce_pytype_skel = {
1620   PyObject_HEAD_INIT(0) 0,              /* Header */
1621   "MPReduce",                           /* @tp_name@ */
1622   sizeof(mpreduce_pyobj),               /* @tp_basicsize@ */
1623   0,                                    /* @tp_itemsize@ */
1624
1625   mpreduce_pydealloc,                   /* @tp_dealloc@ */
1626   0,                                    /* @tp_print@ */
1627   0,                                    /* @tp_getattr@ */
1628   0,                                    /* @tp_setattr@ */
1629   0,                                    /* @tp_compare@ */
1630   0,                                    /* @tp_repr@ */
1631   0,                                    /* @tp_as_number@ */
1632   0,                                    /* @tp_as_sequence@ */
1633   0,                                    /* @tp_as_mapping@ */
1634   0,                                    /* @tp_hash@ */
1635   0,                                    /* @tp_call@ */
1636   0,                                    /* @tp_str@ */
1637   0,                                    /* @tp_getattro@ */
1638   0,                                    /* @tp_setattro@ */
1639   0,                                    /* @tp_as_buffer@ */
1640   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1641     Py_TPFLAGS_BASETYPE,
1642
1643   /* @tp_doc@ */
1644 "A reduction context for reduction modulo primes of special form.",
1645
1646   0,                                    /* @tp_traverse@ */
1647   0,                                    /* @tp_clear@ */
1648   0,                                    /* @tp_richcompare@ */
1649   0,                                    /* @tp_weaklistoffset@ */
1650   0,                                    /* @tp_iter@ */
1651   0,                                    /* @tp_iternext@ */
1652   mpreduce_pymethods,                   /* @tp_methods@ */
1653   0,                                    /* @tp_members@ */
1654   mpreduce_pygetset,                    /* @tp_getset@ */
1655   0,                                    /* @tp_base@ */
1656   0,                                    /* @tp_dict@ */
1657   0,                                    /* @tp_descr_get@ */
1658   0,                                    /* @tp_descr_set@ */
1659   0,                                    /* @tp_dictoffset@ */
1660   0,                                    /* @tp_init@ */
1661   PyType_GenericAlloc,                  /* @tp_alloc@ */
1662   mpreduce_pynew,                       /* @tp_new@ */
1663   0,                                    /* @tp_free@ */
1664   0                                     /* @tp_is_gc@ */
1665 };
1666
1667 /*----- Chinese Remainder Theorem solution --------------------------------*/
1668
1669 typedef struct mpcrt_pyobj {
1670   PyObject_HEAD
1671   mpcrt c;
1672 } mpcrt_pyobj;
1673
1674 #define MPCRT_PY(o) (&((mpcrt_pyobj *)(o))->c)
1675
1676 static PyObject *mcmeth_solve(PyObject *me, PyObject *arg)
1677 {
1678   mpcrt *c = MPCRT_PY(me);
1679   PyObject *q = 0, *x, *z = 0;
1680   mp *xx;
1681   mp **v = 0;
1682   int i = 0, n = c->k;
1683
1684   Py_INCREF(me);
1685   if (PyTuple_Size(arg) == n)
1686     q = arg;
1687   else if (!PyArg_ParseTuple(arg, "O:solve", &q))
1688     goto end;
1689   Py_INCREF(q);
1690   if (!PySequence_Check(q)) TYERR("want a sequence of residues");
1691   if (PySequence_Size(q) != n) VALERR("residue count mismatch");
1692   v = xmalloc(n * sizeof(*v));
1693   for (i = 0; i < n; i++) {
1694     if ((x = PySequence_GetItem(q, i)) == 0) goto end;
1695     xx = getmp(x); Py_DECREF(x); if (!xx) goto end;
1696     v[i] = xx;
1697   }
1698   z = mp_pywrap(mpcrt_solve(c, MP_NEW, v));
1699 end:
1700   if (v) {
1701     n = i;
1702     for (i = 0; i < n; i++)
1703       MP_DROP(v[i]);
1704     xfree(v);
1705   }
1706   Py_DECREF(me);
1707   Py_XDECREF(q);
1708   return (z);
1709 }
1710
1711 static void mpcrt_pydealloc(PyObject *me)
1712 {
1713   mpcrt *c = MPCRT_PY(me);
1714   mpcrt_destroy(c);
1715   xfree(c->v);
1716 }
1717
1718 static PyObject *mpcrt_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1719 {
1720   mpcrt_mod *v = 0;
1721   int n, i = 0, j;
1722   char *kwlist[] = { "mv", 0 };
1723   PyObject *q = 0, *x;
1724   mp *xx = MP_NEW, *y = MP_NEW, *g = MP_NEW;
1725   mpmul mm;
1726   mpcrt_pyobj *c = 0;
1727
1728   if (PyTuple_Size(arg) > 1)
1729     q = arg;
1730   else if (!PyArg_ParseTupleAndKeywords(arg, kw, "O:new", kwlist, &q))
1731     goto end;
1732   Py_INCREF(q);
1733   if (!PySequence_Check(q)) TYERR("want a sequence of moduli");
1734   n = PySequence_Size(q);
1735   if (PyErr_Occurred()) goto end;
1736   if (!n) VALERR("want at least one modulus");
1737   v = xmalloc(n * sizeof(*v));
1738   for (i = 0; i < n; i++) {
1739     if ((x = PySequence_GetItem(q, i)) == 0) goto end;
1740     xx = getmp(x); Py_DECREF(x); if (!xx) goto end;
1741     if (MP_CMP(xx, <=, MP_ZERO)) VALERR("moduli must be positive");
1742     v[i].m = xx; v[i].n = 0; v[i].ni = 0; v[i].nni = 0; xx = MP_NEW;
1743   }
1744   mpmul_init(&mm);
1745   for (j = 0; j < i; j++) mpmul_add(&mm, v[j].m);
1746   xx = mpmul_done(&mm);
1747   for (j = 0; j < i; j++) {
1748     mp_div(&y, 0, xx, v[j].m);
1749     mp_gcd(&g, 0, 0, y, v[j].m);
1750     if (!MP_EQ(g, MP_ONE)) VALERR("moduli must be pairwise coprime");
1751   }
1752
1753   c = (mpcrt_pyobj *)ty->tp_alloc(ty, 0);
1754   mpcrt_create(&c->c, v, n, 0);
1755   Py_DECREF(q);
1756   mp_drop(xx); mp_drop(y); mp_drop(g);
1757   return ((PyObject *)c);
1758
1759 end:
1760   if (v) {
1761     n = i;
1762     for (i = 0; i < n; i++)
1763       MP_DROP(v[i].m);
1764     xfree(v);
1765   }
1766   Py_XDECREF(q);
1767   mp_drop(xx); mp_drop(y); mp_drop(g);
1768   return (0);
1769 }
1770
1771 static PyObject *mcget_product(PyObject *me, void *hunoz)
1772   { return (mp_pywrap(MP_COPY(MPCRT_PY(me)->mb.m))); }
1773
1774 static PyObject *mcget_moduli(PyObject *me, void *hunoz)
1775 {
1776   int i;
1777   PyObject *q;
1778   mpcrt *c = MPCRT_PY(me);
1779
1780   if ((q = PyList_New(c->k)) == 0) return (0);
1781   for (i = 0; i < c->k; i++)
1782     PyList_SetItem(q, i, mp_pywrap(c->v[i].m));
1783   return (q);
1784 }
1785
1786 static PyGetSetDef mpcrt_pygetset[] = {
1787 #define GETSETNAME(op, name) mc##op##_##name
1788   GET   (product,       "C.product -> product of moduli")
1789   GET   (moduli,        "C.moduli -> list of individual moduli")
1790 #undef GETSETNAME
1791   { 0 }
1792 };
1793
1794 static PyMethodDef mpcrt_pymethods[] = {
1795 #define METHNAME(name) mcmeth_##name
1796   METH  (solve,         "C.solve([R0, R1]) -> X mod C.product")
1797 #undef METHNAME
1798   { 0 }
1799 };
1800
1801 static PyTypeObject *mpcrt_pytype, mpcrt_pytype_skel = {
1802   PyObject_HEAD_INIT(0) 0,              /* Header */
1803   "MPCRT",                              /* @tp_name@ */
1804   sizeof(mpcrt_pyobj),                  /* @tp_basicsize@ */
1805   0,                                    /* @tp_itemsize@ */
1806
1807   mpcrt_pydealloc,                      /* @tp_dealloc@ */
1808   0,                                    /* @tp_print@ */
1809   0,                                    /* @tp_getattr@ */
1810   0,                                    /* @tp_setattr@ */
1811   0,                                    /* @tp_compare@ */
1812   0,                                    /* @tp_repr@ */
1813   0,                                    /* @tp_as_number@ */
1814   0,                                    /* @tp_as_sequence@ */
1815   0,                                    /* @tp_as_mapping@ */
1816   0,                                    /* @tp_hash@ */
1817   0,                                    /* @tp_call@ */
1818   0,                                    /* @tp_str@ */
1819   0,                                    /* @tp_getattro@ */
1820   0,                                    /* @tp_setattro@ */
1821   0,                                    /* @tp_as_buffer@ */
1822   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1823     Py_TPFLAGS_BASETYPE,
1824
1825   /* @tp_doc@ */
1826 "A context for the solution of Chinese Remainder Theorem problems.",
1827
1828   0,                                    /* @tp_traverse@ */
1829   0,                                    /* @tp_clear@ */
1830   0,                                    /* @tp_richcompare@ */
1831   0,                                    /* @tp_weaklistoffset@ */
1832   0,                                    /* @tp_iter@ */
1833   0,                                    /* @tp_iternext@ */
1834   mpcrt_pymethods,                      /* @tp_methods@ */
1835   0,                                    /* @tp_members@ */
1836   mpcrt_pygetset,                       /* @tp_getset@ */
1837   0,                                    /* @tp_base@ */
1838   0,                                    /* @tp_dict@ */
1839   0,                                    /* @tp_descr_get@ */
1840   0,                                    /* @tp_descr_set@ */
1841   0,                                    /* @tp_dictoffset@ */
1842   0,                                    /* @tp_init@ */
1843   PyType_GenericAlloc,                  /* @tp_alloc@ */
1844   mpcrt_pynew,                          /* @tp_new@ */
1845   0,                                    /* @tp_free@ */
1846   0                                     /* @tp_is_gc@ */
1847 };
1848
1849 /*----- Binary polynomials ------------------------------------------------*/
1850
1851 static PyObject *gf_pyrepr(PyObject *o)
1852   { return mp_topystring(MP_X(o), 16, "GF(", "0x", "L)"); }
1853
1854 static PyObject *gf_pyrichcompare(PyObject *x, PyObject *y, int op)
1855 {
1856   mp *xx, *yy;
1857   int xl, yl;
1858   int b;
1859
1860   if (gfbinop(x, y, &xx, &yy)) RETURN_NOTIMPL;
1861   switch (op) {
1862     case Py_EQ: b = MP_EQ(xx, yy); break;
1863     case Py_NE: b = !MP_EQ(xx, yy); break;
1864     default:
1865       xl = mp_bits(xx);
1866       yl = mp_bits(yy);
1867       switch (op) {
1868         case Py_LT: b = xl < yl; break;
1869         case Py_LE: b = xl <= yl; break;
1870         case Py_GT: b = xl > yl; break;
1871         case Py_GE: b = xl >= yl; break;
1872         default: abort();
1873       }
1874       break;
1875   }
1876   MP_DROP(xx); MP_DROP(yy);
1877   return (getbool(b));
1878 }
1879
1880 static PyObject *gf_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1881 {
1882   PyObject *x;
1883   mp *z;
1884   mp_pyobj *zz = 0;
1885   int radix = 0;
1886   char *kwlist[] = { "x", "radix", 0 };
1887
1888   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O|i:gf", kwlist, &x, &radix))
1889     goto end;
1890   if (GF_PYCHECK(x)) RETURN_OBJ(x);
1891   if (!good_radix_p(radix, 1)) VALERR("radix out of range");
1892   if ((z = mp_frompyobject(x, radix)) == 0) {
1893     PyErr_Format(PyExc_TypeError, "can't convert %.100s to gf",
1894                  x->ob_type->tp_name);
1895     goto end;
1896   }
1897   if (MP_NEGP(z)) {
1898     MP_DROP(z);
1899     VALERR("gf cannot be negative");
1900   }
1901   zz = (mp_pyobj *)ty->tp_alloc(ty, 0);
1902   zz->x = z;
1903 end:
1904   return ((PyObject *)zz);
1905 }
1906
1907 static PyObject *gf_pyexp(PyObject *x, PyObject *y, PyObject *z)
1908 {
1909   mp *xx = 0, *yy = 0, *zz = 0;
1910   mp *r = 0;
1911   PyObject *rc = 0;
1912
1913   if ((xx = tomp(x)) == 0 || (yy = tomp(y)) == 0 ||
1914       (z && z != Py_None && (zz = tomp(z)) == 0)) {
1915     mp_drop(xx); mp_drop(yy); mp_drop(zz);
1916     RETURN_NOTIMPL;
1917   }
1918   if (!z || z == Py_None) {
1919     if (MP_NEGP(yy)) VALERR("negative exponent");
1920     r = gf_exp(MP_NEW, xx, yy);
1921   } else {
1922     gfreduce gr;
1923     if (MP_ZEROP(zz)) ZDIVERR("zero modulus");
1924     if (MP_NEGP(yy)) {
1925       if ((xx = gf_modinv_checked(xx, xx, zz)) == 0) goto end;
1926       yy = mp_neg(yy, yy);
1927     }
1928     gfreduce_create(&gr, zz);
1929     r = gfreduce_exp(&gr, MP_NEW, xx, yy);
1930     gfreduce_destroy(&gr);
1931   }
1932   rc = gf_pywrap(r);
1933 end:
1934   mp_drop(xx); mp_drop(yy); mp_drop(zz);
1935   return (rc);
1936 }
1937
1938 static PyObject *gfmeth_sqr(PyObject *me, PyObject *arg)
1939 {
1940   if (!PyArg_ParseTuple(arg, ":sqr")) return (0);
1941   return (gf_pywrap(gf_sqr(MP_NEW, MP_X(me))));
1942 }
1943
1944 static PyObject *gfmeth_gcd(PyObject *me, PyObject *arg)
1945 {
1946   PyObject *z = 0;
1947   mp *yy = 0, *zz = MP_NEW;
1948
1949   if (!PyArg_ParseTuple(arg, "O&:gcd", convgf, &yy)) goto end;
1950   gf_gcd(&zz, 0, 0, MP_X(me), yy);
1951   z = gf_pywrap(zz);
1952 end:
1953   if (yy) MP_DROP(yy);
1954   return (z);
1955 }
1956
1957 static PyObject *gfmeth_gcdx(PyObject *me, PyObject *arg)
1958 {
1959   PyObject *z = 0;
1960   mp *yy = 0, *zz = MP_NEW, *uu = MP_NEW, *vv = MP_NEW;
1961
1962   if (!PyArg_ParseTuple(arg, "O&:gcdx", convgf, &yy))
1963     goto end;
1964   gf_gcd(&zz, &uu, &vv, MP_X(me), yy);
1965   z = Py_BuildValue("(NNN)",
1966                     gf_pywrap(zz), gf_pywrap(uu), gf_pywrap(vv));
1967 end:
1968   if (yy) MP_DROP(yy);
1969   return (z);
1970 }
1971
1972 static PyObject *gfmeth_modinv(PyObject *me, PyObject *arg)
1973 {
1974   PyObject *z = 0;
1975   mp *yy = 0, *zz = MP_NEW;
1976
1977   if (!PyArg_ParseTuple(arg, "O&:modinv", convgf, &yy) ||
1978       (zz = gf_modinv_checked(MP_NEW, yy, MP_X(me))) == 0)
1979     goto end;
1980   z = gf_pywrap(zz);
1981 end:
1982   if (yy) MP_DROP(yy);
1983   return (z);
1984 }
1985
1986 static PyObject *gfmeth_irreduciblep(PyObject *me, PyObject *arg)
1987 {
1988   if (!PyArg_ParseTuple(arg, ":irreduciblep")) return (0);
1989   return getbool(gf_irreduciblep(MP_X(me)));
1990 }
1991
1992 static PyObject *gfget_degree(PyObject *me, void *hunoz)
1993   { return (PyInt_FromLong(mp_bits(MP_X(me)) - 1)); }
1994
1995 static PyGetSetDef gf_pygetset[] = {
1996 #define GETSETNAME(op, name) gf##op##_##name
1997   GET   (degree,        "X.degree -> polynomial degree of X")
1998 #undef GETSETNAME
1999 #define GETSETNAME(op, name) mp##op##_##name
2000   GET   (nbits,         "X.nbits -> bit length of X")
2001   GET   (noctets,       "X.noctets -> octet length of X")
2002 #undef GETSETNAME
2003   { 0 }
2004 };
2005
2006 static PyMethodDef gf_pymethods[] = {
2007 #define METHNAME(func) gfmeth_##func
2008   METH  (setbit,        "X.setbit(N) -> X with bit N set")
2009   METH  (clearbit,      "X.clearbit(N) -> X with bit N clear")
2010   METH  (testbit,       "X.testbit(N) -> true/false if bit N set/clear in X")
2011   METH  (sqr,           "X.sqr() -> X^2")
2012   METH  (gcd,           "X.gcd(Y) -> gcd(X, Y)")
2013   METH  (gcdx,
2014          "X.gcdx(Y) -> (gcd(X, Y), U, V) with X U + Y V = gcd(X, Y)")
2015   METH  (modinv,        "X.modinv(Y) -> multiplicative inverse of Y mod X")
2016   METH  (irreduciblep,  "X.irreduciblep() -> true/false")
2017 #undef METHNAME
2018 #define METHNAME(func) mpmeth_##func
2019   KWMETH(tostring,      "X.tostring([radix = 10]) -> STR")
2020   KWMETH(storel,        "X.storel([len = -1]) -> little-endian bytes")
2021   KWMETH(storeb,        "X.storeb([len = -1]) -> big-endian bytes")
2022   KWMETH(storel2c,
2023          "X.storel2c([len = -1]) -> little-endian bytes, two's complement")
2024   KWMETH(storeb2c,
2025          "X.storeb2c([len = -1]) -> big-endian bytes, two's complement")
2026   METH  (tobuf,         "X.tobuf() -> buffer format")
2027 #undef METHNAME
2028   { 0 }
2029 };
2030
2031 static PyNumberMethods gf_pynumber = {
2032   gf_pyadd,                             /* @nb_add@ */
2033   gf_pysub,                             /* @nb_subtract@ */
2034   gf_pymul,                             /* @nb_multiply@ */
2035   0,                                    /* @nb_divide@ */
2036   gf_pymod,                             /* @nb_remainder@ */
2037   gf_pydivmod,                          /* @nb_divmod@ */
2038   gf_pyexp,                             /* @nb_power@ */
2039   mp_pyid,                              /* @nb_negative@ */
2040   mp_pyid,                              /* @nb_positive@ */
2041   mp_pyid,                              /* @nb_absolute@ */
2042   mp_pynonzerop,                        /* @nb_nonzero@ */
2043   0 /* doesn't make any sense */,       /* @nb_invert@ */
2044   gf_pylsl,                             /* @nb_lshift@ */
2045   gf_pylsr,                             /* @nb_rshift@ */
2046   gf_pyand,                             /* @nb_and@ */
2047   gf_pyxor,                             /* @nb_xor@ */
2048   gf_pyor,                              /* @nb_or@ */
2049   gf_pycoerce,                          /* @nb_coerce@ */
2050   mp_pyint,                             /* @nb_int@ */
2051   mp_pylong,                            /* @nb_long@ */
2052   0 /* doesn't make any sense */,       /* @nb_float@ */
2053   mp_pyoct,                             /* @nb_oct@ */
2054   mp_pyhex,                             /* @nb_hex@ */
2055
2056   0,                                    /* @nb_inplace_add@ */
2057   0,                                    /* @nb_inplace_subtract@ */
2058   0,                                    /* @nb_inplace_multiply@ */
2059   0,                                    /* @nb_inplace_divide@ */
2060   0,                                    /* @nb_inplace_remainder@ */
2061   0,                                    /* @nb_inplace_power@ */
2062   0,                                    /* @nb_inplace_lshift@ */
2063   0,                                    /* @nb_inplace_rshift@ */
2064   0,                                    /* @nb_inplace_and@ */
2065   0,                                    /* @nb_inplace_xor@ */
2066   0,                                    /* @nb_inplace_or@ */
2067
2068   gf_pydiv,                             /* @nb_floor_divide@ */
2069   0,                                    /* @nb_true_divide@ */
2070   0,                                    /* @nb_inplace_floor_divide@ */
2071   0,                                    /* @nb_inplace_true_divide@ */
2072 };
2073
2074 static PyTypeObject gf_pytype_skel = {
2075   PyObject_HEAD_INIT(0) 0,              /* Header */
2076   "GF",                                 /* @tp_name@ */
2077   sizeof(mp_pyobj),                     /* @tp_basicsize@ */
2078   0,                                    /* @tp_itemsize@ */
2079
2080   mp_pydealloc,                         /* @tp_dealloc@ */
2081   0,                                    /* @tp_print@ */
2082   0,                                    /* @tp_getattr@ */
2083   0,                                    /* @tp_setattr@ */
2084   0,                                    /* @tp_compare@ */
2085   gf_pyrepr,                            /* @tp_repr@ */
2086   &gf_pynumber,                         /* @tp_as_number@ */
2087   0,                                    /* @tp_as_sequence@ */
2088   0,                                    /* @tp_as_mapping@ */
2089   mp_pyhash,                            /* @tp_hash@ */
2090   0,                                    /* @tp_call@ */
2091   mp_pyhex,                             /* @tp_str@ */
2092   0,                                    /* @tp_getattro@ */
2093   0,                                    /* @tp_setattro@ */
2094   0,                                    /* @tp_as_buffer@ */
2095   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
2096     Py_TPFLAGS_CHECKTYPES |
2097     Py_TPFLAGS_BASETYPE,
2098
2099   /* @tp_doc@ */
2100 "Binary polynomials.  Support almost all the standard arithmetic\n\
2101 operations.\n\
2102 \n\
2103 Constructor GF(X, [radix = R]) attempts to convert X to a `GF'.  If\n\
2104 X is a string, it's read in radix-R form, or we look for a prefix\n\
2105 if R = 0.  Other acceptable things are field elements, elliptic curve\n\
2106 points, group elements, Python `int' and `long' objects, and anything\n\
2107 with an integer conversion.\n\
2108 \n\
2109 The name is hopelessly wrong from a technical point of view, but\n\
2110 but it's much easier to type than `p2' or `c2' or whatever.\n\
2111 \n\
2112 Notes:\n\
2113 \n\
2114   * Use `//' for division.  GFs don't have `/' division.",
2115
2116   0,                                    /* @tp_traverse@ */
2117   0,                                    /* @tp_clear@ */
2118   gf_pyrichcompare,                     /* @tp_richcompare@ */
2119   0,                                    /* @tp_weaklistoffset@ */
2120   0,                                    /* @tp_iter@ */
2121   0,                                    /* @tp_iternext@ */
2122   gf_pymethods,                         /* @tp_methods@ */
2123   0,                                    /* @tp_members@ */
2124   gf_pygetset,                          /* @tp_getset@ */
2125   0,                                    /* @tp_base@ */
2126   0,                                    /* @tp_dict@ */
2127   0,                                    /* @tp_descr_get@ */
2128   0,                                    /* @tp_descr_set@ */
2129   0,                                    /* @tp_dictoffset@ */
2130   0,                                    /* @tp_init@ */
2131   PyType_GenericAlloc,                  /* @tp_alloc@ */
2132   gf_pynew,                             /* @tp_new@ */
2133   0,                                    /* @tp_free@ */
2134   0                                     /* @tp_is_gc@ */
2135 };
2136
2137 static PyObject *meth__GF_fromstring(PyObject *me,
2138                                     PyObject *arg, PyObject *kw)
2139 {
2140   int r = 0;
2141   char *p;
2142   Py_ssize_t len;
2143   PyObject *z = 0;
2144   mp *zz;
2145   mptext_stringctx sc;
2146   char *kwlist[] = { "class", "x", "radix", 0 };
2147
2148   if (!PyArg_ParseTupleAndKeywords(arg, kw, "Os#|i:fromstring",
2149                                    kwlist, &me, &p, &len, &r))
2150     goto end;
2151   if (!good_radix_p(r, 1)) VALERR("bad radix");
2152   sc.buf = p; sc.lim = p + len;
2153   if ((zz = mp_read(MP_NEW, r, &mptext_stringops, &sc)) == 0 ||
2154       MP_NEGP(zz)) {
2155     if (zz) MP_DROP(zz);
2156     VALERR("bad binary polynomial");
2157   }
2158   z = Py_BuildValue("(Ns#)", gf_pywrap(zz),
2159                     sc.buf, (Py_ssize_t)(sc.lim - sc.buf));
2160 end:
2161   return (z);
2162 }
2163
2164 /*----- Sparse poly reduction ---------------------------------------------*/
2165
2166 typedef struct gfreduce_pyobj {
2167   PyObject_HEAD
2168   gfreduce mr;
2169 } gfreduce_pyobj;
2170
2171 #define GFREDUCE_PY(o) (&((gfreduce_pyobj *)(o))->mr)
2172
2173 static PyObject *grmeth_exp(PyObject *me, PyObject *arg)
2174 {
2175   PyObject *rc = 0;
2176   mp *yy = 0, *zz = 0;
2177
2178   if (!PyArg_ParseTuple(arg, "O&O&:exp", convgf, &yy, convgf, &zz))
2179     goto end;
2180   if (MP_NEGP(zz)) {
2181     if ((yy = gf_modinv_checked(yy, yy, GFREDUCE_PY(me)->p)) == 0) goto end;
2182     zz = mp_neg(zz, zz);
2183   }
2184   rc = gf_pywrap(gfreduce_exp(GFREDUCE_PY(me), MP_NEW, yy, zz));
2185 end:
2186   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
2187   return (rc);
2188 }
2189
2190 static PyObject *grmeth_trace(PyObject *me, PyObject *arg)
2191 {
2192   PyObject *rc = 0;
2193   mp *xx = 0;
2194
2195   if (!PyArg_ParseTuple(arg, "O&:trace", convgf, &xx)) goto end;
2196   rc = PyInt_FromLong(gfreduce_trace(GFREDUCE_PY(me), xx));
2197 end:
2198   if (xx) MP_DROP(xx);
2199   return (rc);
2200 }
2201
2202 static PyObject *grmeth_halftrace(PyObject *me, PyObject *arg)
2203 {
2204   PyObject *rc = 0;
2205   mp *xx = 0;
2206
2207   if (!PyArg_ParseTuple(arg, "O&:halftrace", convgf, &xx)) goto end;
2208   rc = gf_pywrap(gfreduce_halftrace(GFREDUCE_PY(me), MP_NEW, xx));
2209 end:
2210   if (xx) MP_DROP(xx);
2211   return (rc);
2212 }
2213
2214 static PyObject *grmeth_sqrt(PyObject *me, PyObject *arg)
2215 {
2216   PyObject *rc = 0;
2217   mp *xx = 0, *yy;
2218
2219   if (!PyArg_ParseTuple(arg, "O&:sqrt", convgf, &xx)) goto end;
2220   if ((yy = gfreduce_sqrt(GFREDUCE_PY(me), MP_NEW, xx)) == 0)
2221     VALERR("no modular square root");
2222   rc = gf_pywrap(yy);
2223 end:
2224   if (xx) MP_DROP(xx);
2225   return (rc);
2226 }
2227
2228 static PyObject *grmeth_quadsolve(PyObject *me, PyObject *arg)
2229 {
2230   PyObject *rc = 0;
2231   mp *xx = 0, *yy;
2232
2233   if (!PyArg_ParseTuple(arg, "O&:quadsolve", convgf, &xx)) goto end;
2234   if ((yy = gfreduce_quadsolve(GFREDUCE_PY(me), MP_NEW, xx)) == 0)
2235     VALERR("no solution found");
2236   rc = gf_pywrap(yy);
2237 end:
2238   if (xx) MP_DROP(xx);
2239   return (rc);
2240 }
2241
2242 static PyObject *grmeth_reduce(PyObject *me, PyObject *arg)
2243 {
2244   PyObject *z = 0;
2245   mp *yy = 0;
2246
2247   if (!PyArg_ParseTuple(arg, "O&:reduce", convgf, &yy)) goto end;
2248   z = gf_pywrap(gfreduce_do(GFREDUCE_PY(me), MP_NEW, yy));
2249 end:
2250   return (z);
2251 }
2252
2253 static void gfreduce_pydealloc(PyObject *me)
2254 {
2255   gfreduce_destroy(GFREDUCE_PY(me));
2256   FREEOBJ(me);
2257 }
2258
2259 static PyObject *gfreduce_pynew(PyTypeObject *ty,
2260                                  PyObject *arg, PyObject *kw)
2261 {
2262   gfreduce_pyobj *mr = 0;
2263   gfreduce r;
2264   char *kwlist[] = { "m", 0 };
2265   mp *xx = 0;
2266
2267   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convgf, &xx))
2268     goto end;
2269   if (MP_ZEROP(xx)) ZDIVERR("modulus is zero!");
2270   gfreduce_create(&r, xx);
2271   mr = (gfreduce_pyobj *)ty->tp_alloc(ty, 0);
2272   mr->mr = r;
2273 end:
2274   if (xx) MP_DROP(xx);
2275   return ((PyObject *)mr);
2276 }
2277
2278 static PyObject *grget_m(PyObject *me, void *hunoz)
2279   { return (gf_pywrap(MP_COPY(GFREDUCE_PY(me)->p))); }
2280
2281 static PyGetSetDef gfreduce_pygetset[] = {
2282 #define GETSETNAME(op, name) gr##op##_##name
2283   GET   (m,             "R.m -> reduction polynomial")
2284 #undef GETSETNAME
2285   { 0 }
2286 };
2287
2288 static PyMethodDef gfreduce_pymethods[] = {
2289 #define METHNAME(name) grmeth_##name
2290   METH  (reduce,        "R.reduce(X) -> X mod B.m")
2291   METH  (trace,        "R.trace(X) -> Tr(X) = x + x^2 + ... + x^{2^{m - 1}}")
2292   METH  (halftrace,   "R.halftrace(X) -> x + x^{2^2} + ... + x^{2^{m - 1}}")
2293   METH  (sqrt,          "R.sqrt(X) -> Y where Y^2 = X mod R")
2294   METH  (quadsolve,     "R.quadsolve(X) -> Y where Y^2 + Y = X mod R")
2295   METH  (exp,           "R.exp(X, N) -> X^N mod B.m")
2296 #undef METHNAME
2297   { 0 }
2298 };
2299
2300 static PyTypeObject *gfreduce_pytype, gfreduce_pytype_skel = {
2301   PyObject_HEAD_INIT(0) 0,              /* Header */
2302   "GFReduce",                           /* @tp_name@ */
2303   sizeof(gfreduce_pyobj),               /* @tp_basicsize@ */
2304   0,                                    /* @tp_itemsize@ */
2305
2306   gfreduce_pydealloc,                   /* @tp_dealloc@ */
2307   0,                                    /* @tp_print@ */
2308   0,                                    /* @tp_getattr@ */
2309   0,                                    /* @tp_setattr@ */
2310   0,                                    /* @tp_compare@ */
2311   0,                                    /* @tp_repr@ */
2312   0,                                    /* @tp_as_number@ */
2313   0,                                    /* @tp_as_sequence@ */
2314   0,                                    /* @tp_as_mapping@ */
2315   0,                                    /* @tp_hash@ */
2316   0,                                    /* @tp_call@ */
2317   0,                                    /* @tp_str@ */
2318   0,                                    /* @tp_getattro@ */
2319   0,                                    /* @tp_setattro@ */
2320   0,                                    /* @tp_as_buffer@ */
2321   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
2322     Py_TPFLAGS_BASETYPE,
2323
2324   /* @tp_doc@ */
2325 "A reduction context for reduction modulo sparse irreducible polynomials.",
2326
2327   0,                                    /* @tp_traverse@ */
2328   0,                                    /* @tp_clear@ */
2329   0,                                    /* @tp_richcompare@ */
2330   0,                                    /* @tp_weaklistoffset@ */
2331   0,                                    /* @tp_iter@ */
2332   0,                                    /* @tp_iternext@ */
2333   gfreduce_pymethods,                   /* @tp_methods@ */
2334   0,                                    /* @tp_members@ */
2335   gfreduce_pygetset,                    /* @tp_getset@ */
2336   0,                                    /* @tp_base@ */
2337   0,                                    /* @tp_dict@ */
2338   0,                                    /* @tp_descr_get@ */
2339   0,                                    /* @tp_descr_set@ */
2340   0,                                    /* @tp_dictoffset@ */
2341   0,                                    /* @tp_init@ */
2342   PyType_GenericAlloc,                  /* @tp_alloc@ */
2343   gfreduce_pynew,                       /* @tp_new@ */
2344   0,                                    /* @tp_free@ */
2345   0                                     /* @tp_is_gc@ */
2346 };
2347
2348 /*----- Normal/poly transformation ----------------------------------------*/
2349
2350 typedef struct gfn_pyobj {
2351   PyObject_HEAD
2352   mp *p;
2353   gfn ntop, pton;
2354 } gfn_pyobj;
2355
2356 static PyTypeObject *gfn_pytype, gfn_pytype_skel;
2357
2358 #define GFN_P(o) (((gfn_pyobj *)(o))->p)
2359 #define GFN_PTON(o) (&((gfn_pyobj *)(o))->pton)
2360 #define GFN_NTOP(o) (&((gfn_pyobj *)(o))->ntop)
2361
2362 static PyObject *gfn_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
2363 {
2364   mp *p = 0, *beta = 0;
2365   gfn_pyobj *gg = 0;
2366   char *kwlist[] = { "p", "beta", 0 };
2367
2368   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&O&:new", kwlist,
2369                                    convgf, &p, convgf, &beta))
2370     goto end;
2371   gg = PyObject_New(gfn_pyobj, ty);
2372   if (gfn_create(p, beta, &gg->ntop, &gg->pton)) {
2373     FREEOBJ(gg);
2374     gg = 0;
2375     VALERR("can't invert transformation matrix");
2376   }
2377   gg->p = MP_COPY(p);
2378 end:
2379   mp_drop(p);
2380   mp_drop(beta);
2381   return ((PyObject *)gg);
2382 }
2383
2384 static PyObject *gfnget_p(PyObject *me, void *hunoz)
2385   { return (gf_pywrap(MP_COPY(GFN_P(me)))); }
2386
2387 static PyObject *gfnget_beta(PyObject *me, void *hunoz)
2388 {
2389   gfn *n = GFN_NTOP(me);
2390   mp *x = n->r[n->n - 1];
2391   return (gf_pywrap(MP_COPY(x)));
2392 }
2393
2394 #define XFORMOP(name, NAME)                                             \
2395   static PyObject *gfnmeth_##name(PyObject *me, PyObject *arg)          \
2396   {                                                                     \
2397     mp *xx = 0;                                                         \
2398     mp *z = 0;                                                          \
2399                                                                         \
2400     if (!PyArg_ParseTuple(arg, "O&:" #name, convgf, &xx)) goto end;     \
2401     z = gfn_transform(GFN_##NAME(me), MP_NEW, xx);                      \
2402   end:                                                                  \
2403     mp_drop(xx);                                                        \
2404     if (!z) return (0);                                                 \
2405     return (mp_pywrap(z));                                              \
2406   }
2407 XFORMOP(pton, PTON)
2408 XFORMOP(ntop, NTOP)
2409 #undef XFORMOP
2410
2411 static void gfn_pydealloc(PyObject *me)
2412 {
2413   gfn_destroy(GFN_PTON(me));
2414   gfn_destroy(GFN_NTOP(me));
2415   FREEOBJ(me);
2416 }
2417
2418 static PyGetSetDef gfn_pygetset[] = {
2419 #define GETSETNAME(op, name) gfn##op##_##name
2420   GET   (p,             "X.p -> polynomial basis, as polynomial")
2421   GET   (beta,          "X.beta -> normal basis element, in poly form")
2422 #undef GETSETNAME
2423   { 0 }
2424 };
2425
2426 static PyMethodDef gfn_pymethods[] = {
2427 #define METHNAME(name) gfnmeth_##name
2428   METH  (pton,          "X.pton(A) -> normal-basis representation of A")
2429   METH  (ntop,          "X.ntop(A) -> polynomial-basis representation of A")
2430 #undef METHNAME
2431   { 0 }
2432 };
2433
2434 static PyTypeObject gfn_pytype_skel = {
2435   PyObject_HEAD_INIT(0) 0,              /* Header */
2436   "GFN",                                /* @tp_name@ */
2437   sizeof(gfn_pyobj),                    /* @tp_basicsize@ */
2438   0,                                    /* @tp_itemsize@ */
2439
2440   gfn_pydealloc,                        /* @tp_dealloc@ */
2441   0,                                    /* @tp_print@ */
2442   0,                                    /* @tp_getattr@ */
2443   0,                                    /* @tp_setattr@ */
2444   0,                                    /* @tp_compare@ */
2445   0,                                    /* @tp_repr@ */
2446   0,                                    /* @tp_as_number@ */
2447   0,                                    /* @tp_as_sequence@ */
2448   0,                                    /* @tp_as_mapping@ */
2449   0,                                    /* @tp_hash@ */
2450   0,                                    /* @tp_call@ */
2451   0,                                    /* @tp_str@ */
2452   0,                                    /* @tp_getattro@ */
2453   0,                                    /* @tp_setattro@ */
2454   0,                                    /* @tp_as_buffer@ */
2455   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
2456     Py_TPFLAGS_BASETYPE,
2457
2458   /* @tp_doc@ */
2459 "An object for transforming elements of binary fields between polynomial\n\
2460 and normal basis representations.",
2461
2462   0,                                    /* @tp_traverse@ */
2463   0,                                    /* @tp_clear@ */
2464   0,                                    /* @tp_richcompare@ */
2465   0,                                    /* @tp_weaklistoffset@ */
2466   0,                                    /* @tp_iter@ */
2467   0,                                    /* @tp_iternext@ */
2468   gfn_pymethods,                        /* @tp_methods@ */
2469   0,                                    /* @tp_members@ */
2470   gfn_pygetset,                         /* @tp_getset@ */
2471   0,                                    /* @tp_base@ */
2472   0,                                    /* @tp_dict@ */
2473   0,                                    /* @tp_descr_get@ */
2474   0,                                    /* @tp_descr_set@ */
2475   0,                                    /* @tp_dictoffset@ */
2476   0,                                    /* @tp_init@ */
2477   PyType_GenericAlloc,                  /* @tp_alloc@ */
2478   gfn_pynew,                            /* @tp_new@ */
2479   0,                                    /* @tp_free@ */
2480   0                                     /* @tp_is_gc@ */
2481 };
2482
2483 /*----- Glue --------------------------------------------------------------*/
2484
2485 static PyMethodDef methods[] = {
2486 #define METHNAME(func) meth_##func
2487   KWMETH(_MP_fromstring,        "\
2488 fromstring(STR, [radix = 0]) -> (X, REST)\n\
2489 \n\
2490 Parse STR as a large integer, according to radix.  If radix is zero,\n\
2491 read a prefix from STR to decide radix: allow `0' for octal, `0x' for hex\n\
2492 or `R_' for other radix R.")
2493   KWMETH(_GF_fromstring,        "\
2494 fromstring(STR, [radix = 0]) -> (X, REST)\n\
2495 \n\
2496 Parse STR as a binary polynomial, according to radix.  If radix is zero,\n\
2497 read a prefix from STR to decide radix: allow `0' for octal, `0x' for hex\n\
2498 or `R_' for other radix R.")
2499   METH  (_MP_factorial,         "\
2500 factorial(I) -> I!: compute factorial")
2501   METH  (_MP_fibonacci,         "\
2502 fibonacci(I) -> F(I): compute Fibonacci number")
2503   METH  (_MP_loadl,             "\
2504 loadl(STR) -> X: read little-endian bytes")
2505   METH  (_MP_loadb,             "\
2506 loadb(STR) -> X: read big-endian bytes")
2507   METH  (_MP_loadl2c,           "\
2508 loadl2c(STR) -> X: read little-endian bytes, two's complement")
2509   METH  (_MP_loadb2c,           "\
2510 loadb2c(STR) -> X: read big-endian bytes, two's complement")
2511   METH  (_MP_frombuf,           "\
2512 frombuf(STR) -> (X, REST): read buffer format")
2513   METH  (_GF_loadl,             "\
2514 loadl(STR) -> X: read little-endian bytes")
2515   METH  (_GF_loadb,             "\
2516 loadb(STR) -> X: read big-endian bytes")
2517   METH  (_GF_frombuf,           "\
2518 frombuf(STR) -> (X, REST): read buffer format")
2519 #undef METHNAME
2520   { 0 }
2521 };
2522
2523 void mp_pyinit(void)
2524 {
2525   INITTYPE(mp, root);
2526   INITTYPE(gf, root);
2527   INITTYPE(mpmul, root);
2528   INITTYPE(mpmont, root);
2529   INITTYPE(mpbarrett, root);
2530   INITTYPE(mpreduce, root);
2531   INITTYPE(mpcrt, root);
2532   INITTYPE(gfreduce, root);
2533   INITTYPE(gfn, root);
2534   addmethods(methods);
2535 }
2536
2537 void mp_pyinsert(PyObject *mod)
2538 {
2539   INSERT("MP", mp_pytype);
2540   INSERT("MPMul", mpmul_pytype);
2541   INSERT("MPMont", mpmont_pytype);
2542   INSERT("MPBarrett", mpbarrett_pytype);
2543   INSERT("MPReduce", mpreduce_pytype);
2544   INSERT("MPCRT", mpcrt_pytype);
2545   INSERT("GF", gf_pytype);
2546   INSERT("GFReduce", gfreduce_pytype);
2547   INSERT("GFN", gfn_pytype);
2548 }
2549
2550 /*----- That's all, folks -------------------------------------------------*/