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

   1 /* -*-c-*-
   2  *
   3  * $Id$
   4  *
   5  * Multiprecision arithmetic
   6  *
   7  * (c) 2004 Straylight/Edgeware
   8  */
   9
  10 /*----- Licensing notice --------------------------------------------------*
  11  *
  12  * This file is part of the Python interface to Catacomb.
  13  *
  14  * Catacomb/Python is free software; you can redistribute it and/or modify
  15  * it under the terms of the GNU General Public License as published by
  16  * the Free Software Foundation; either version 2 of the License, or
  17  * (at your option) any later version.
  18  *
  19  * Catacomb/Python is distributed in the hope that it will be useful,
  20  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  21  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  22  * GNU General Public License for more details.
  23  *
  24  * You should have received a copy of the GNU General Public License
  25  * along with Catacomb/Python; if not, write to the Free Software Foundation,
  26  * Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
  27  */
  28
  29 /*----- Header files ------------------------------------------------------*/
  30
  31 #include "catacomb-python.h"
  32
  33 /*----- General utilities -------------------------------------------------*/
  34
  35 PyTypeObject *mp_pytype = 0;
  36 PyTypeObject *gf_pytype = 0;
  37
  38 mp *mp_frompylong(PyObject *obj)
  39 {
  40   unsigned long bits;
  41   PyLongObject *l = (PyLongObject *)obj;
  42   int sz;
  43   size_t w;
  44   mpd r = 0;
  45   int b = 0;
  46   int i;
  47   mp *x;
  48   mpw *p;
  49
  50   sz = l->ob_size;
  51   if (sz < 0) sz = -sz;
  52   assert(MPW_BITS >= SHIFT);
  53   bits = (unsigned long)sz * SHIFT;
  54   w = (bits + MPW_BITS - 1)/MPW_BITS;
  55   x = mp_new(w, l->ob_size < 0 ? MP_NEG : 0);
  56   p = x->v;
  57   for (i = 0; i < sz; i++) {
  58     r |= (mpd)l->ob_digit[i] << b;
  59     b += SHIFT;
  60     while (b >= MPW_BITS) {
  61       *p++ = MPW(r);
  62       r >>= MPW_BITS;
  63       b -= MPW_BITS;
  64     }
  65   }
  66   while (r) {
  67     *p++ = MPW(r);
  68     r >>= MPW_BITS;
  69   }
  70   x->vl = p;
  71   MP_SHRINK(x);
  72   return (x);
  73 }
  74
  75 PyObject *mp_topylong(mp *x)
  76 {
  77   unsigned long bits = mp_bits(x);
  78   int sz = (bits + SHIFT - 1)/SHIFT;
  79   PyLongObject *l = _PyLong_New(sz);
  80   mpd r = 0;
  81   int b = 0;
  82   mpw *p = x->v;
  83   int i = 0;
  84
  85   assert(MPW_BITS >= SHIFT);
  86   while (i < sz && p < x->vl) {
  87     r |= (mpd)*p++ << b;
  88     b += MPW_BITS;
  89     while (i < sz && b >= SHIFT) {
  90       l->ob_digit[i++] = r & MASK;
  91       r >>= SHIFT;
  92       b -= SHIFT;
  93     }
  94   }
  95   while (i < sz && r) {
  96     l->ob_digit[i++] = r & MASK;
  97     r >>= SHIFT;
  98   }
  99   l->ob_size = (x->f & MP_NEG) ? -sz : sz;
 100   return ((PyObject *)l);
 101 }
 102
 103 mp *mp_frompyobject(PyObject *o, int radix)
 104 {
 105   mp *x;
 106
 107   if (PyString_Check(o)) {
 108     mptext_stringctx sc;
 109     mp *x;
 110     sc.buf = PyString_AS_STRING(o);
 111     sc.lim = sc.buf + PyString_GET_SIZE(o);
 112     x = mp_read(MP_NEW, radix, &mptext_stringops, &sc);
 113     if (!x) return (0);
 114     if (sc.buf < sc.lim) { MP_DROP(x); return (0); }
 115     return (x);
 116   }
 117   if ((x = tomp(o)) != 0)
 118     return (x);
 119   return (0);
 120 }
 121
 122 PyObject *mp_topystring(mp *x, int radix, const char *xpre,
 123                         const char *pre, const char *post)
 124 {
 125   int len = mptext_len(x, radix) + 1;
 126   mptext_stringctx sc;
 127   PyObject *o;
 128   size_t xprelen = xpre ? strlen(xpre) : 0;
 129   size_t prelen = pre ? strlen(pre) : 0;
 130   size_t postlen = post ? strlen(post) : 0;
 131   char *p;
 132   MP_COPY(x);
 133   o = PyString_FromStringAndSize(0, len + 1 + xprelen + prelen + postlen);
 134   p = PyString_AS_STRING(o);
 135   sc.buf = p;
 136   if (xpre) { memcpy(sc.buf, xpre, xprelen); sc.buf += xprelen; }
 137   if (MP_NEGP(x)) { *sc.buf++ = '-'; x = mp_neg(x, x); }
 138   if (pre) { memcpy(sc.buf, pre, prelen); sc.buf += prelen; }
 139   sc.lim = sc.buf + len;
 140   mp_write(x, radix, &mptext_stringops, &sc);
 141   if (post) { memcpy(sc.buf, post, postlen); sc.buf += postlen; }
 142   MP_DROP(x);
 143   _PyString_Resize(&o, sc.buf - p);
 144   return (o);
 145 }
 146
 147 static int good_radix_p(int r, int readp)
 148 {
 149   return ((r >= -255 && r <= -2) ||
 150           (readp && r == 0) ||
 151           (r >= 2 && r <= 62));
 152 }
 153
 154 PyObject *mp_pywrap(mp *x)
 155 {
 156   mp_pyobj *z = PyObject_New(mp_pyobj, mp_pytype);
 157   z->x = x;
 158   return ((PyObject *)z);
 159 }
 160
 161 PyObject *gf_pywrap(mp *x)
 162 {
 163   mp_pyobj *z = PyObject_New(mp_pyobj, gf_pytype);
 164   z->x = x;
 165   return ((PyObject *)z);
 166 }
 167
 168 int mp_tolong_checked(mp *x, long *l)
 169 {
 170   static mp *longmin = 0, *longmax = 0;
 171   int rc = -1;
 172
 173   if (!longmax) {
 174     longmin = mp_fromlong(MP_NEW, LONG_MIN);
 175     longmax = mp_fromlong(MP_NEW, LONG_MAX);
 176   }
 177   if (MP_CMP(x, <, longmin) || MP_CMP(x, >, longmax))
 178     VALERR("mp out of range for int");
 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)) 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)) return (0);
 489   return (PyInt_FromLong(l));
 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 static long mp_pyhash(PyObject *me)
 548 {
 549   long h;
 550   PyObject *l = mp_topylong(MP_X(me)); h = PyObject_Hash(l);
 551   Py_DECREF(l); return (h);
 552 }
 553
 554 static PyObject *mpmeth_jacobi(PyObject *me, PyObject *arg)
 555 {
 556   mp *y = 0;
 557   PyObject *z = 0;
 558
 559   if (!PyArg_ParseTuple(arg, "O&:jacobi", convmp, &y)) goto end;
 560   if (MP_NEGP(MP_X(me)) || MP_EVENP(MP_X(me)))
 561     VALERR("must be positive and odd");
 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 #define STOREOP(name, c)                                                \
 688   static PyObject *mpmeth_##name(PyObject *me,                          \
 689                                  PyObject *arg, PyObject *kw)           \
 690   {                                                                     \
 691     long len = -1;                                                      \
 692     char *kwlist[] = { "len", 0 };                                      \
 693     PyObject *rc = 0;                                                   \
 694                                                                         \
 695     if (!PyArg_ParseTupleAndKeywords(arg, kw, "|l:" #name,              \
 696                                     kwlist, &len))                      \
 697       goto end;                                                         \
 698     if (len < 0) {                                                      \
 699       len = mp_octets##c(MP_X(me));                                     \
 700       if (!len) len = 1;                                                \
 701     }                                                                   \
 702     rc = bytestring_pywrap(0, len);                                     \
 703     mp_##name(MP_X(me), PyString_AS_STRING(rc), len);                   \
 704   end:                                                                  \
 705     return (rc);                                                        \
 706   }
 707 STOREOP(storel, )
 708 STOREOP(storeb, )
 709 STOREOP(storel2c, 2c)
 710 STOREOP(storeb2c, 2c)
 711 #undef STOREOP
 712
 713 #define BUFOP(ty, pyty)                                                 \
 714   static PyObject *meth__##pyty##_frombuf(PyObject *me, PyObject *arg)  \
 715   {                                                                     \
 716     buf b;                                                              \
 717     char *p;                                                            \
 718     int sz;                                                             \
 719     PyObject *rc = 0;                                                   \
 720     mp *x;                                                              \
 721                                                                         \
 722     if (!PyArg_ParseTuple(arg, "Os#:frombuf", &me, &p, &sz)) goto end;  \
 723     buf_init(&b, p, sz);                                                \
 724     if ((x = buf_getmp(&b)) == 0) VALERR("malformed data");             \
 725     rc = Py_BuildValue("(NN)", ty##_pywrap(x),                          \
 726                        bytestring_pywrapbuf(&b));                       \
 727   end:                                                                  \
 728     return (rc);                                                        \
 729   }
 730 BUFOP(mp, MP)
 731 BUFOP(gf, GF)
 732 #undef BUFOP
 733
 734 static PyObject *mpmeth_tobuf(PyObject *me, PyObject *arg)
 735 {
 736   buf b;
 737   PyObject *rc;
 738   mp *x;
 739   size_t n;
 740
 741   if (!PyArg_ParseTuple(arg, ":tobuf")) return (0);
 742   x = MP_X(me);
 743   n = mp_octets(x) + 3;
 744   rc = bytestring_pywrap(0, n);
 745   buf_init(&b, PyString_AS_STRING(rc), n);
 746   buf_putmp(&b, x);
 747   assert(BOK(&b));
 748   _PyString_Resize(&rc, BLEN(&b));
 749   return (rc);
 750 }
 751
 752 static PyObject *mpmeth_primep(PyObject *me, PyObject *arg, PyObject *kw)
 753 {
 754   grand *r = &rand_global;
 755   char *kwlist[] = { "rng", 0 };
 756   PyObject *rc = 0;
 757
 758   if (!PyArg_ParseTupleAndKeywords(arg, kw, "|O&", kwlist, convgrand, &r))
 759     goto end;
 760   rc = getbool(pgen_primep(MP_X(me), r));
 761 end:
 762   return (rc);
 763 }
 764
 765 static PyObject *mpget_nbits(PyObject *me, void *hunoz)
 766   { return (PyInt_FromLong(mp_bits(MP_X(me)))); }
 767
 768 static PyObject *mpget_noctets(PyObject *me, void *hunoz)
 769   { return (PyInt_FromLong(mp_octets(MP_X(me)))); }
 770
 771 static PyObject *mpget_noctets2c(PyObject *me, void *hunoz)
 772   { return (PyInt_FromLong(mp_octets2c(MP_X(me)))); }
 773
 774 static PyGetSetDef mp_pygetset[] = {
 775 #define GETSETNAME(op, func) mp##op##_##func
 776   GET   (nbits,         "X.nbits -> bit length of X")
 777   GET   (noctets,       "X.noctets -> octet length of X")
 778   GET   (noctets2c,     "X.noctets2c -> two's complement octet length of X")
 779 #undef GETSETNAME
 780   { 0 }
 781 };
 782
 783 static PyMethodDef mp_pymethods[] = {
 784 #define METHNAME(func) mpmeth_##func
 785   METH  (jacobi,        "X.jacobi(Y) -> Jacobi symbol (Y/X) (NB inversion!)")
 786   METH  (setbit,        "X.setbit(N) -> X with bit N set")
 787   METH  (clearbit,      "X.clearbit(N) -> X with bit N clear")
 788   METH  (testbit,       "X.testbit(N) -> true/false if bit N set/clear in X")
 789   METH  (odd,           "X.odd() -> S, T where X = 2^S T with T odd")
 790   METH  (sqr,           "X.sqr() -> X^2")
 791   METH  (sqrt,          "X.sqrt() -> largest integer <= sqrt(X)")
 792   METH  (gcd,           "X.gcd(Y) -> gcd(X, Y)")
 793   METH  (gcdx,
 794          "X.gcdx(Y) -> (gcd(X, Y), U, V) with X U + Y V = gcd(X, Y)")
 795   METH  (modinv,        "X.modinv(Y) -> multiplicative inverse of Y mod X")
 796   METH  (modsqrt,       "X.modsqrt(Y) -> square root of Y mod X, if X prime")
 797   KWMETH(primep,        "X.primep(rng = rand) -> true/false if X is prime")
 798   KWMETH(tostring,      "X.tostring(radix = 10) -> STR")
 799   KWMETH(storel,        "X.storel(len = -1) -> little-endian bytes")
 800   KWMETH(storeb,        "X.storeb(len = -1) -> big-endian bytes")
 801   KWMETH(storel2c,
 802          "X.storel2c(len = -1) -> little-endian bytes, two's complement")
 803   KWMETH(storeb2c,
 804          "X.storeb2c(len = -1) -> big-endian bytes, two's complement")
 805   METH  (tobuf,         "X.tobuf() -> buffer format")
 806 #undef METHNAME
 807   { 0 }
 808 };
 809
 810 static PyNumberMethods mp_pynumber = {
 811   mp_pyadd,                             /* @nb_add@ */
 812   mp_pysub,                             /* @nb_subtract@ */
 813   mp_pymul,                             /* @nb_multiply@ */
 814   0,                                    /* @nb_divide@ */
 815   mp_pymod,                             /* @nb_remainder@ */
 816   mp_pydivmod,                          /* @nb_divmod@ */
 817   mp_pyexp,                             /* @nb_power@ */
 818   mp_pyneg,                             /* @nb_negative@ */
 819   mp_pyid,                              /* @nb_positive@ */
 820   mp_pyabs,                             /* @nb_absolute@ */
 821   mp_pynonzerop,                        /* @nb_nonzero@ */
 822   mp_pynot2c,                           /* @nb_invert@ */
 823   mp_pylsl2c,                           /* @nb_lshift@ */
 824   mp_pylsr2c,                           /* @nb_rshift@ */
 825   mp_pyand2c,                           /* @nb_and@ */
 826   mp_pyxor2c,                           /* @nb_xor@ */
 827   mp_pyor2c,                            /* @nb_or@ */
 828   mp_pycoerce,                          /* @nb_coerce@ */
 829   mp_pyint,                             /* @nb_int@ */
 830   mp_pylong,                            /* @nb_long@ */
 831   mp_pyfloat,                           /* @nb_float@ */
 832   mp_pyoct,                             /* @nb_oct@ */
 833   mp_pyhex,                             /* @nb_hex@ */
 834
 835   0,                                    /* @nb_inplace_add@ */
 836   0,                                    /* @nb_inplace_subtract@ */
 837   0,                                    /* @nb_inplace_multiply@ */
 838   0,                                    /* @nb_inplace_divide@ */
 839   0,                                    /* @nb_inplace_remainder@ */
 840   0,                                    /* @nb_inplace_power@ */
 841   0,                                    /* @nb_inplace_lshift@ */
 842   0,                                    /* @nb_inplace_rshift@ */
 843   0,                                    /* @nb_inplace_and@ */
 844   0,                                    /* @nb_inplace_xor@ */
 845   0,                                    /* @nb_inplace_or@ */
 846
 847   mp_pydiv,                             /* @nb_floor_divide@ */
 848   0,                                    /* @nb_true_divide@ */
 849   0,                                    /* @nb_inplace_floor_divide@ */
 850   0,                                    /* @nb_inplace_true_divide@ */
 851 };
 852
 853 static PyTypeObject mp_pytype_skel = {
 854   PyObject_HEAD_INIT(0) 0,              /* Header */
 855   "catacomb.MP",                        /* @tp_name@ */
 856   sizeof(mp_pyobj),                     /* @tp_basicsize@ */
 857   0,                                    /* @tp_itemsize@ */
 858
 859   mp_pydealloc,                         /* @tp_dealloc@ */
 860   0,                                    /* @tp_print@ */
 861   0,                                    /* @tp_getattr@ */
 862   0,                                    /* @tp_setattr@ */
 863   mp_pycompare,                         /* @tp_compare@ */
 864   mp_pyrepr,                            /* @tp_repr@ */
 865   &mp_pynumber,                         /* @tp_as_number@ */
 866   0,                                    /* @tp_as_sequence@ */
 867   0,                                    /* @tp_as_mapping@ */
 868   mp_pyhash,                            /* @tp_hash@ */
 869   0,                                    /* @tp_call@ */
 870   mp_pystr,                             /* @tp_str@ */
 871   0,                                    /* @tp_getattro@ */
 872   0,                                    /* @tp_setattro@ */
 873   0,                                    /* @tp_as_buffer@ */
 874   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
 875     Py_TPFLAGS_CHECKTYPES |
 876     Py_TPFLAGS_BASETYPE,
 877
 878   /* @tp_doc@ */
 879 "Multiprecision integers, similar to `long' but more efficient and\n\
 880 versatile.  Support all the standard arithmetic operations.\n\
 881 \n\
 882 Constructor mp(X, radix = R) attempts to convert X to an `mp'.  If\n\
 883 X is a string, it's read in radix-R form, or we look for a prefix\n\
 884 if R = 0.  Other acceptable things are ints and longs.\n\
 885 \n\
 886 Notes:\n\
 887 \n\
 888   * Use `//' for division.  MPs don't have `/' division.",
 889
 890   0,                                    /* @tp_traverse@ */
 891   0,                                    /* @tp_clear@ */
 892   0,                                    /* @tp_richcompare@ */
 893   0,                                    /* @tp_weaklistoffset@ */
 894   0,                                    /* @tp_iter@ */
 895   0,                                    /* @tp_iternext@ */
 896   mp_pymethods,                         /* @tp_methods@ */
 897   0,                                    /* @tp_members@ */
 898   mp_pygetset,                          /* @tp_getset@ */
 899   0,                                    /* @tp_base@ */
 900   0,                                    /* @tp_dict@ */
 901   0,                                    /* @tp_descr_get@ */
 902   0,                                    /* @tp_descr_set@ */
 903   0,                                    /* @tp_dictoffset@ */
 904   0,                                    /* @tp_init@ */
 905   PyType_GenericAlloc,                  /* @tp_alloc@ */
 906   mp_pynew,                             /* @tp_new@ */
 907   0,                                    /* @tp_free@ */
 908   0                                     /* @tp_is_gc@ */
 909 };
 910
 911 static PyObject *meth__MP_fromstring(PyObject *me,
 912                                      PyObject *arg, PyObject *kw)
 913 {
 914   int r = 0;
 915   char *p;
 916   int len;
 917   PyObject *z = 0;
 918   mp *zz;
 919   mptext_stringctx sc;
 920   char *kwlist[] = { "class", "x", "radix", 0 };
 921
 922   if (!PyArg_ParseTupleAndKeywords(arg, kw, "Os#|i:fromstring",
 923                                    kwlist, &me, &p, &len, &r))
 924     goto end;
 925   if (!good_radix_p(r, 1)) VALERR("bad radix");
 926   sc.buf = p; sc.lim = p + len;
 927   if ((zz = mp_read(MP_NEW, r, &mptext_stringops, &sc)) == 0)
 928     SYNERR("bad integer");
 929   z = Py_BuildValue("(Ns#)", mp_pywrap(zz), sc.buf, (int)(sc.lim - sc.buf));
 930 end:
 931   return (z);
 932 }
 933
 934 #define LOADOP(pre, py, name)                                           \
 935   static PyObject *meth__##py##_##name(PyObject *me, PyObject *arg)     \
 936   {                                                                     \
 937     char *p;                                                            \
 938     int len;                                                            \
 939     if (!PyArg_ParseTuple(arg, "Os#:" #name, &me, &p, &len)) return (0); \
 940     return (pre##_pywrap(mp_##name(MP_NEW, p, len)));                   \
 941   }
 942 LOADOP(mp, MP, loadl)
 943 LOADOP(mp, MP, loadb)
 944 LOADOP(mp, MP, loadl2c)
 945 LOADOP(mp, MP, loadb2c)
 946 LOADOP(gf, GF, loadl)
 947 LOADOP(gf, GF, loadb)
 948 #undef LOADOP
 949
 950 /*----- Products of small integers ----------------------------------------*/
 951
 952 typedef struct mpmul_pyobj {
 953   PyObject_HEAD
 954   int livep;
 955   mpmul mm;
 956 } mpmul_pyobj;
 957
 958 #define MPMUL_LIVEP(o) (((mpmul_pyobj *)(o))->livep)
 959 #define MPMUL_PY(o) (&((mpmul_pyobj *)(o))->mm)
 960
 961 static void mpmul_pydealloc(PyObject *me)
 962 {
 963   if (MPMUL_LIVEP(me))
 964     mp_drop(mpmul_done(MPMUL_PY(me)));
 965   FREEOBJ(me);
 966 }
 967
 968 static PyObject *mmmeth_factor(PyObject *me, PyObject *arg)
 969 {
 970   PyObject *q, *i;
 971   mp *x;
 972
 973   if (!MPMUL_LIVEP(me)) VALERR("MPMul object invalid");
 974   if (PyTuple_Size(arg) != 1)
 975     i = PyObject_GetIter(arg);
 976   else {
 977     if ((q = PyTuple_GetItem(arg, 0)) == 0) goto end;
 978     if ((i = PyObject_GetIter(q)) == 0) {
 979       PyErr_Clear(); /* that's ok */
 980       i = PyObject_GetIter(arg);
 981     }
 982   }
 983   if (!i) goto end;
 984   while ((q = PyIter_Next(i)) != 0) {
 985     x = getmp(q); Py_DECREF(q); if (!x) {
 986       Py_DECREF(i);
 987       goto end;
 988     }
 989     mpmul_add(MPMUL_PY(me), x);
 990     MP_DROP(x);
 991   }
 992   Py_DECREF(i);
 993   RETURN_ME;
 994 end:
 995   return (0);
 996 }
 997
 998 static PyObject *mmmeth_done(PyObject *me, PyObject *arg)
 999 {
1000   mp *x;
1001
1002   if (!PyArg_ParseTuple(arg, ":done")) goto end;
1003   if (!MPMUL_LIVEP(me)) VALERR("MPMul object invalid");
1004   x = mpmul_done(MPMUL_PY(me));
1005   MPMUL_LIVEP(me) = 0;
1006   return (mp_pywrap(x));
1007 end:
1008   return (0);
1009 }
1010
1011 static PyObject *mpmul_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1012 {
1013   mpmul_pyobj *mm;
1014
1015   if (kw) TYERR("keyword arguments not allowed here");
1016   mm = (mpmul_pyobj *)ty->tp_alloc(ty, 0);
1017   mpmul_init(&mm->mm);
1018   mm->livep = 1;
1019   if (mmmeth_factor((PyObject *)mm, arg) == 0) {
1020     Py_DECREF(mm);
1021     goto end;
1022   }
1023   return ((PyObject *)mm);
1024 end:
1025   return (0);
1026 }
1027
1028 static PyObject *mmget_livep(PyObject *me, void *hunoz)
1029   { return (getbool(MPMUL_LIVEP(me))); }
1030
1031 static PyGetSetDef mpmul_pygetset[] = {
1032 #define GETSETNAME(op, name) mm##op##_##name
1033   GET   (livep,                 "MM.livep -> flag: object still valid?")
1034 #undef GETSETNAME
1035   { 0 }
1036 };
1037
1038 static PyMethodDef mpmul_pymethods[] = {
1039 #define METHNAME(name) mmmeth_##name
1040   METH  (factor,                "MM.factor(ITERABLE) or MM.factor(I, ...)")
1041   METH  (done,                  "MM.done() -> PRODUCT")
1042 #undef METHNAME
1043   { 0 }
1044 };
1045
1046 static PyTypeObject *mpmul_pytype, mpmul_pytype_skel = {
1047   PyObject_HEAD_INIT(0) 0,              /* Header */
1048   "catacomb.MPMul",                     /* @tp_name@ */
1049   sizeof(mpmul_pyobj),                  /* @tp_basicsize@ */
1050   0,                                    /* @tp_itemsize@ */
1051
1052   mpmul_pydealloc,                      /* @tp_dealloc@ */
1053   0,                                    /* @tp_print@ */
1054   0,                                    /* @tp_getattr@ */
1055   0,                                    /* @tp_setattr@ */
1056   0,                                    /* @tp_compare@ */
1057   0,                                    /* @tp_repr@ */
1058   0,                                    /* @tp_as_number@ */
1059   0,                                    /* @tp_as_sequence@ */
1060   0,                                    /* @tp_as_mapping@ */
1061   0,                                    /* @tp_hash@ */
1062   0,                                    /* @tp_call@ */
1063   0,                                    /* @tp_str@ */
1064   0,                                    /* @tp_getattro@ */
1065   0,                                    /* @tp_setattro@ */
1066   0,                                    /* @tp_as_buffer@ */
1067   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1068     Py_TPFLAGS_BASETYPE,
1069
1070   /* @tp_doc@ */
1071 "An object for multiplying many small integers.",
1072
1073   0,                                    /* @tp_traverse@ */
1074   0,                                    /* @tp_clear@ */
1075   0,                                    /* @tp_richcompare@ */
1076   0,                                    /* @tp_weaklistoffset@ */
1077   0,                                    /* @tp_iter@ */
1078   0,                                    /* @tp_iternext@ */
1079   mpmul_pymethods,                      /* @tp_methods@ */
1080   0,                                    /* @tp_members@ */
1081   mpmul_pygetset,                       /* @tp_getset@ */
1082   0,                                    /* @tp_base@ */
1083   0,                                    /* @tp_dict@ */
1084   0,                                    /* @tp_descr_get@ */
1085   0,                                    /* @tp_descr_set@ */
1086   0,                                    /* @tp_dictoffset@ */
1087   0,                                    /* @tp_init@ */
1088   PyType_GenericAlloc,                  /* @tp_alloc@ */
1089   mpmul_pynew,                          /* @tp_new@ */
1090   0,                                    /* @tp_free@ */
1091   0                                     /* @tp_is_gc@ */
1092 };
1093
1094 /*----- Montgomery reduction ----------------------------------------------*/
1095
1096 typedef struct mpmont_pyobj {
1097   PyObject_HEAD
1098   mpmont mm;
1099 } mpmont_pyobj;
1100
1101 #define MPMONT_PY(o) (&((mpmont_pyobj *)(o))->mm)
1102
1103 static PyObject *mmmeth_int(PyObject *me, PyObject *arg)
1104 {
1105   PyObject *z = 0;
1106   mp *yy = 0;
1107   mpmont *mm = MPMONT_PY(me);
1108
1109   if (!PyArg_ParseTuple(arg, "O&:in", convmp, &yy))
1110     goto end;
1111   mp_div(0, &yy, yy, mm->m);
1112   z = mp_pywrap(mpmont_mul(mm, MP_NEW, yy, mm->r2));
1113 end:
1114   if (yy) MP_DROP(yy);
1115   return (z);
1116 }
1117
1118 static PyObject *mmmeth_mul(PyObject *me, PyObject *arg)
1119 {
1120   PyObject *rc = 0;
1121   mp *yy = 0, *zz = 0;
1122
1123   if (!PyArg_ParseTuple(arg, "O&O&:mul", convmp, &yy, convmp, &zz))
1124     goto end;
1125   rc = mp_pywrap(mpmont_mul(MPMONT_PY(me), MP_NEW, yy, zz));
1126 end:
1127   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1128   return (rc);
1129 }
1130
1131 static PyObject *mmmeth_exp(PyObject *me, PyObject *arg)
1132 {
1133   PyObject *rc = 0;
1134   mp *yy = 0, *zz = 0;
1135
1136   if (!PyArg_ParseTuple(arg, "O&O&:exp", convmp, &yy, convmp, &zz))
1137     goto end;
1138   if (MP_NEGP(zz)) {
1139     if ((yy = mp_modinv_checked(yy, yy, MPMONT_PY(me)->m)) == 0) goto end;
1140     zz = mp_neg(zz, zz);
1141   }
1142   rc = mp_pywrap(mpmont_exp(MPMONT_PY(me), MP_NEW, yy, zz));
1143 end:
1144   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1145   return (rc);
1146 }
1147
1148 static PyObject *mmmeth_expr(PyObject *me, PyObject *arg)
1149 {
1150   PyObject *rc = 0;
1151   mp *yy = 0, *zz = 0;
1152
1153   if (!PyArg_ParseTuple(arg, "O&O&:expr", convmp, &yy, convmp, &zz))
1154     goto end;
1155   if (MP_NEGP(zz)) {
1156     yy = mpmont_reduce(MPMONT_PY(me), yy, yy);
1157     if ((yy = mp_modinv_checked(yy, yy, MPMONT_PY(me)->m)) == 0) goto end;
1158     yy = mpmont_mul(MPMONT_PY(me), yy, yy, MPMONT_PY(me)->r2);
1159     zz = mp_neg(zz, zz);
1160   }
1161   rc = mp_pywrap(mpmont_expr(MPMONT_PY(me), MP_NEW, yy, zz));
1162 end:
1163   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1164   return (rc);
1165 }
1166
1167 static PyObject *mm_mexpr_id(PyObject *me)
1168   { return mp_pywrap(MP_COPY(MPMONT_PY(me)->r)); }
1169
1170 static int mm_mexpr_fill(void *p, PyObject *me, PyObject *x, PyObject *y)
1171 {
1172   mp *xx = 0, *yy = 0;
1173   mp_expfactor *f = p;
1174   mpmont *mm = MPMONT_PY(me);
1175
1176   if ((xx = getmp(x)) == 0 || (yy = getmp(y)) == 0)
1177     goto fail;
1178   if (MP_NEGP(yy)) {
1179     xx = mpmont_reduce(mm, xx, xx);
1180     if ((xx = mp_modinv_checked(xx, xx, yy)) == 0)
1181       goto fail;
1182     xx = mpmont_mul(mm, xx, xx, mm->r2);
1183     yy = mp_neg(yy, yy);
1184   }
1185   f->base = xx;
1186   f->exp = yy;
1187   return (0);
1188
1189 fail:
1190   mp_drop(xx); mp_drop(yy);
1191   return (-1);
1192 }
1193
1194 static PyObject *mm_mexpr(PyObject *me, void *v, int n)
1195   { return mp_pywrap(mpmont_mexpr(MPMONT_PY(me), MP_NEW, v, n)); }
1196
1197 static void mp_mexp_drop(void *p)
1198 {
1199   mp_expfactor *f = p;
1200   mp_drop(f->base);
1201   mp_drop(f->exp);
1202 }
1203
1204 static PyObject *mmmeth_mexpr(PyObject *me, PyObject *arg)
1205 {
1206   return mexp_common(me, arg, sizeof(mp_expfactor),
1207                      mm_mexpr_id, mm_mexpr_fill, mm_mexpr, mp_mexp_drop);
1208 }
1209
1210 static PyObject *mp_mexp_id(PyObject *me)
1211   { return mp_pywrap(MP_ONE); }
1212
1213 static int mp_mexp_fill(void *p, PyObject *me, PyObject *x, PyObject *y)
1214 {
1215   mp *xx = 0, *yy = 0;
1216   mp_expfactor *f = p;
1217
1218   if ((xx = getmp(x)) == 0 || (yy = getmp(y)) == 0)
1219     goto fail;
1220   if (MP_NEGP(yy)) {
1221     if ((xx = mp_modinv_checked(xx, xx, yy)) == 0)
1222       goto fail;
1223     yy = mp_neg(yy, yy);
1224   }
1225   f->base = xx;
1226   f->exp = yy;
1227   return (0);
1228
1229 fail:
1230   mp_drop(xx); mp_drop(yy);
1231   return (-1);
1232 }
1233
1234 static PyObject *mm_mexp(PyObject *me, void *v, int n)
1235   { return mp_pywrap(mpmont_mexp(MPMONT_PY(me), MP_NEW, v, n)); }
1236
1237 static PyObject *mmmeth_mexp(PyObject *me, PyObject *arg)
1238 {
1239   return mexp_common(me, arg, sizeof(mp_expfactor),
1240                      mp_mexp_id, mp_mexp_fill, mm_mexp, mp_mexp_drop);
1241 }
1242
1243 #define mmmeth_ext mmmeth_reduce
1244 static PyObject *mmmeth_reduce(PyObject *me, PyObject *arg)
1245 {
1246   PyObject *z = 0;
1247   mp *yy = 0;
1248
1249   if (!PyArg_ParseTuple(arg, "O&", convmp, &yy)) goto end;
1250   z = mp_pywrap(mpmont_reduce(MPMONT_PY(me), MP_NEW, yy));
1251 end:
1252   return (z);
1253 }
1254
1255 static void mpmont_pydealloc(PyObject *me)
1256 {
1257   mpmont_destroy(MPMONT_PY(me));
1258   FREEOBJ(me);
1259 }
1260
1261 static PyObject *mpmont_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1262 {
1263   mpmont_pyobj *mm = 0;
1264   char *kwlist[] = { "m", 0 };
1265   mp *xx = 0;
1266
1267   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convmp, &xx))
1268     goto end;
1269   if (!MP_POSP(xx) || !MP_ODDP(xx)) VALERR("m must be positive and odd");
1270   mm = (mpmont_pyobj *)ty->tp_alloc(ty, 0);
1271   mpmont_create(&mm->mm, xx);
1272 end:
1273   if (xx) MP_DROP(xx);
1274   return ((PyObject *)mm);
1275 }
1276
1277 static PyObject *mmget_m(PyObject *me, void *hunoz)
1278   { return (mp_pywrap(MP_COPY(MPMONT_PY(me)->m))); }
1279
1280 static PyObject *mmget_r(PyObject *me, void *hunoz)
1281   { return (mp_pywrap(MP_COPY(MPMONT_PY(me)->r))); }
1282
1283 static PyObject *mmget_r2(PyObject *me, void *hunoz)
1284   { return (mp_pywrap(MP_COPY(MPMONT_PY(me)->r2))); }
1285
1286 static PyGetSetDef mpmont_pygetset[] = {
1287 #define GETSETNAME(op, name) mm##op##_##name
1288   GET   (m,             "M.m -> modulus for reduction")
1289   GET   (r,             "M.r -> multiplicative identity")
1290   GET   (r2,            "M.r2 -> M.r^2, Montgomerization factor")
1291 #undef GETSETNAME
1292   { 0 }
1293 };
1294
1295 static PyMethodDef mpmont_pymethods[] = {
1296 #define METHNAME(name) mmmeth_##name
1297   METH  (int,           "M.out(X) -> XR")
1298   METH  (mul,           "M.mul(XR, YR) -> ZR where Z = X Y")
1299   METH  (expr,          "M.expr(XR, N) -> ZR where Z = X^N mod M.m")
1300   METH  (mexpr,         "\
1301 B.mexp([(XR0, N0), (XR1, N1), ...]) = ZR where Z = X0^N0 X1^N1 mod B.m\n\
1302 \t(the list may be flattened if this more convenient.)")
1303   METH  (reduce,        "M.reduce(XR) -> X")
1304   METH  (ext,           "M.ext(XR) -> X")
1305   METH  (exp,           "M.exp(X, N) -> X^N mod M.m")
1306   METH  (mexp,          "\
1307 B.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 mod B.m\n\
1308 \t(the list may be flattened if this more convenient.)")
1309 #undef METHNAME
1310   { 0 }
1311 };
1312
1313 static PyTypeObject *mpmont_pytype, mpmont_pytype_skel = {
1314   PyObject_HEAD_INIT(0) 0,              /* Header */
1315   "catacomb.MPMont",                    /* @tp_name@ */
1316   sizeof(mpmont_pyobj),                 /* @tp_basicsize@ */
1317   0,                                    /* @tp_itemsize@ */
1318
1319   mpmont_pydealloc,                     /* @tp_dealloc@ */
1320   0,                                    /* @tp_print@ */
1321   0,                                    /* @tp_getattr@ */
1322   0,                                    /* @tp_setattr@ */
1323   0,                                    /* @tp_compare@ */
1324   0,                                    /* @tp_repr@ */
1325   0,                                    /* @tp_as_number@ */
1326   0,                                    /* @tp_as_sequence@ */
1327   0,                                    /* @tp_as_mapping@ */
1328   0,                                    /* @tp_hash@ */
1329   0,                                    /* @tp_call@ */
1330   0,                                    /* @tp_str@ */
1331   0,                                    /* @tp_getattro@ */
1332   0,                                    /* @tp_setattro@ */
1333   0,                                    /* @tp_as_buffer@ */
1334   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1335     Py_TPFLAGS_BASETYPE,
1336
1337   /* @tp_doc@ */
1338 "A Montgomery reduction context.",
1339
1340   0,                                    /* @tp_traverse@ */
1341   0,                                    /* @tp_clear@ */
1342   0,                                    /* @tp_richcompare@ */
1343   0,                                    /* @tp_weaklistoffset@ */
1344   0,                                    /* @tp_iter@ */
1345   0,                                    /* @tp_iternext@ */
1346   mpmont_pymethods,                     /* @tp_methods@ */
1347   0,                                    /* @tp_members@ */
1348   mpmont_pygetset,                      /* @tp_getset@ */
1349   0,                                    /* @tp_base@ */
1350   0,                                    /* @tp_dict@ */
1351   0,                                    /* @tp_descr_get@ */
1352   0,                                    /* @tp_descr_set@ */
1353   0,                                    /* @tp_dictoffset@ */
1354   0,                                    /* @tp_init@ */
1355   PyType_GenericAlloc,                  /* @tp_alloc@ */
1356   mpmont_pynew,                         /* @tp_new@ */
1357   0,                                    /* @tp_free@ */
1358   0                                     /* @tp_is_gc@ */
1359 };
1360
1361 /*----- Barrett reduction -------------------------------------------------*/
1362
1363 typedef struct mpbarrett_pyobj {
1364   PyObject_HEAD
1365   mpbarrett mb;
1366 } mpbarrett_pyobj;
1367
1368 #define MPBARRETT_PY(o) (&((mpbarrett_pyobj *)(o))->mb)
1369
1370 static PyObject *mbmeth_exp(PyObject *me, PyObject *arg)
1371 {
1372   PyObject *rc = 0;
1373   mp *yy = 0, *zz = 0;
1374
1375   if (!PyArg_ParseTuple(arg, "O&O&:exp", convmp, &yy, convmp, &zz))
1376     goto end;
1377   if (MP_NEGP(zz)) {
1378     if ((yy = mp_modinv_checked(yy, yy, MPBARRETT_PY(me)->m)) == 0) goto end;
1379     zz = mp_neg(zz, zz);
1380   }
1381   rc = mp_pywrap(mpbarrett_exp(MPBARRETT_PY(me), MP_NEW, yy, zz));
1382 end:
1383   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1384   return (rc);
1385 }
1386
1387 static PyObject *mb_mexp(PyObject *me, void *v, int n)
1388   { return mp_pywrap(mpbarrett_mexp(MPBARRETT_PY(me), MP_NEW, v, n)); }
1389
1390 static PyObject *mbmeth_mexp(PyObject *me, PyObject *arg)
1391 {
1392   return mexp_common(me, arg, sizeof(mp_expfactor),
1393                      mp_mexp_id, mp_mexp_fill, mb_mexp, mp_mexp_drop);
1394 }
1395
1396 static PyObject *mbmeth_reduce(PyObject *me, PyObject *arg)
1397 {
1398   PyObject *z = 0;
1399   mp *yy = 0;
1400
1401   if (!PyArg_ParseTuple(arg, "O&:reduce", convmp, &yy))
1402     goto end;
1403   z = mp_pywrap(mpbarrett_reduce(MPBARRETT_PY(me), MP_NEW, yy));
1404 end:
1405   return (z);
1406 }
1407
1408 static void mpbarrett_pydealloc(PyObject *me)
1409 {
1410   mpbarrett_destroy(MPBARRETT_PY(me));
1411   FREEOBJ(me);
1412 }
1413
1414 static PyObject *mpbarrett_pynew(PyTypeObject *ty,
1415                                  PyObject *arg, PyObject *kw)
1416 {
1417   mpbarrett_pyobj *mb = 0;
1418   char *kwlist[] = { "m", 0 };
1419   mp *xx = 0;
1420
1421   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convmp, &xx))
1422     goto end;
1423   if (!MP_POSP(xx)) VALERR("m must be positive");
1424   mb = (mpbarrett_pyobj *)ty->tp_alloc(ty, 0);
1425   mpbarrett_create(&mb->mb, xx);
1426 end:
1427   if (xx) MP_DROP(xx);
1428   return ((PyObject *)mb);
1429 }
1430
1431 static PyObject *mbget_m(PyObject *me, void *hunoz)
1432   { return (mp_pywrap(MP_COPY(MPBARRETT_PY(me)->m))); }
1433
1434 static PyGetSetDef mpbarrett_pygetset[] = {
1435 #define GETSETNAME(op, name) mb##op##_##name
1436   GET   (m,             "B.m -> modulus for reduction")
1437 #undef GETSETNAME
1438   { 0 }
1439 };
1440
1441 static PyMethodDef mpbarrett_pymethods[] = {
1442 #define METHNAME(name) mbmeth_##name
1443   METH  (reduce,        "B.reduce(X) -> X mod B.m")
1444   METH  (exp,           "B.exp(X, N) -> X^N mod B.m")
1445   METH  (mexp,          "\
1446 B.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 mod B.m\n\
1447 \t(the list may be flattened if this more convenient.)")
1448 #undef METHNAME
1449   { 0 }
1450 };
1451
1452 static PyTypeObject *mpbarrett_pytype, mpbarrett_pytype_skel = {
1453   PyObject_HEAD_INIT(0) 0,              /* Header */
1454   "catacomb.MPBarrett",                 /* @tp_name@ */
1455   sizeof(mpbarrett_pyobj),              /* @tp_basicsize@ */
1456   0,                                    /* @tp_itemsize@ */
1457
1458   mpbarrett_pydealloc,                  /* @tp_dealloc@ */
1459   0,                                    /* @tp_print@ */
1460   0,                                    /* @tp_getattr@ */
1461   0,                                    /* @tp_setattr@ */
1462   0,                                    /* @tp_compare@ */
1463   0,                                    /* @tp_repr@ */
1464   0,                                    /* @tp_as_number@ */
1465   0,                                    /* @tp_as_sequence@ */
1466   0,                                    /* @tp_as_mapping@ */
1467   0,                                    /* @tp_hash@ */
1468   0,                                    /* @tp_call@ */
1469   0,                                    /* @tp_str@ */
1470   0,                                    /* @tp_getattro@ */
1471   0,                                    /* @tp_setattro@ */
1472   0,                                    /* @tp_as_buffer@ */
1473   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1474     Py_TPFLAGS_BASETYPE,
1475
1476   /* @tp_doc@ */
1477 "A Barrett reduction context.",
1478
1479   0,                                    /* @tp_traverse@ */
1480   0,                                    /* @tp_clear@ */
1481   0,                                    /* @tp_richcompare@ */
1482   0,                                    /* @tp_weaklistoffset@ */
1483   0,                                    /* @tp_iter@ */
1484   0,                                    /* @tp_iternext@ */
1485   mpbarrett_pymethods,                  /* @tp_methods@ */
1486   0,                                    /* @tp_members@ */
1487   mpbarrett_pygetset,                   /* @tp_getset@ */
1488   0,                                    /* @tp_base@ */
1489   0,                                    /* @tp_dict@ */
1490   0,                                    /* @tp_descr_get@ */
1491   0,                                    /* @tp_descr_set@ */
1492   0,                                    /* @tp_dictoffset@ */
1493   0,                                    /* @tp_init@ */
1494   PyType_GenericAlloc,                  /* @tp_alloc@ */
1495   mpbarrett_pynew,                      /* @tp_new@ */
1496   0,                                    /* @tp_free@ */
1497   0                                     /* @tp_is_gc@ */
1498 };
1499
1500 /*----- Nice prime reduction ----------------------------------------------*/
1501
1502 typedef struct mpreduce_pyobj {
1503   PyObject_HEAD
1504   mpreduce mr;
1505 } mpreduce_pyobj;
1506
1507 #define MPREDUCE_PY(o) (&((mpreduce_pyobj *)(o))->mr)
1508
1509 static PyObject *mrmeth_exp(PyObject *me, PyObject *arg)
1510 {
1511   PyObject *rc = 0;
1512   mp *yy = 0, *zz = 0;
1513
1514   if (!PyArg_ParseTuple(arg, "O&O&:exp", convmp, &yy, convmp, &zz))
1515     goto end;
1516   if (MP_NEGP(zz)) {
1517     if ((yy = mp_modinv_checked(yy, yy, MPREDUCE_PY(me)->p)) == 0) goto end;
1518     zz = mp_neg(zz, zz);
1519   }
1520   rc = mp_pywrap(mpreduce_exp(MPREDUCE_PY(me), MP_NEW, yy, zz));
1521 end:
1522   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1523   return (rc);
1524 }
1525
1526 static PyObject *mrmeth_reduce(PyObject *me, PyObject *arg)
1527 {
1528   PyObject *z = 0;
1529   mp *yy = 0;
1530
1531   if (!PyArg_ParseTuple(arg, "O&:reduce", convmp, &yy)) goto end;
1532   z = mp_pywrap(mpreduce_do(MPREDUCE_PY(me), MP_NEW, yy));
1533 end:
1534   return (z);
1535 }
1536
1537 static void mpreduce_pydealloc(PyObject *me)
1538 {
1539   mpreduce_destroy(MPREDUCE_PY(me));
1540   FREEOBJ(me);
1541 }
1542
1543 static PyObject *mpreduce_pynew(PyTypeObject *ty,
1544                                  PyObject *arg, PyObject *kw)
1545 {
1546   mpreduce_pyobj *mr = 0;
1547   mpreduce r;
1548   char *kwlist[] = { "m", 0 };
1549   mp *xx = 0;
1550
1551   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convmp, &xx))
1552     goto end;
1553   if (!MP_POSP(xx)) VALERR("m must be positive");
1554   if (mpreduce_create(&r, xx)) VALERR("bad modulus (must be 2^k - ...)");
1555   mr = (mpreduce_pyobj *)ty->tp_alloc(ty, 0);
1556   mr->mr = r;
1557 end:
1558   if (xx) MP_DROP(xx);
1559   return ((PyObject *)mr);
1560 }
1561
1562 static PyObject *mrget_m(PyObject *me, void *hunoz)
1563   { return (mp_pywrap(MP_COPY(MPREDUCE_PY(me)->p))); }
1564
1565 static PyGetSetDef mpreduce_pygetset[] = {
1566 #define GETSETNAME(op, name) mr##op##_##name
1567   GET   (m,             "R.m -> modulus for reduction")
1568 #undef GETSETNAME
1569   { 0 }
1570 };
1571
1572 static PyMethodDef mpreduce_pymethods[] = {
1573 #define METHNAME(name) mrmeth_##name
1574   METH  (reduce,        "R.reduce(X) -> X mod B.m")
1575   METH  (exp,           "R.exp(X, N) -> X^N mod B.m")
1576 #undef METHNAME
1577   { 0 }
1578 };
1579
1580 static PyTypeObject *mpreduce_pytype, mpreduce_pytype_skel = {
1581   PyObject_HEAD_INIT(0) 0,              /* Header */
1582   "catacomb.MPReduce",                  /* @tp_name@ */
1583   sizeof(mpreduce_pyobj),               /* @tp_basicsize@ */
1584   0,                                    /* @tp_itemsize@ */
1585
1586   mpreduce_pydealloc,                   /* @tp_dealloc@ */
1587   0,                                    /* @tp_print@ */
1588   0,                                    /* @tp_getattr@ */
1589   0,                                    /* @tp_setattr@ */
1590   0,                                    /* @tp_compare@ */
1591   0,                                    /* @tp_repr@ */
1592   0,                                    /* @tp_as_number@ */
1593   0,                                    /* @tp_as_sequence@ */
1594   0,                                    /* @tp_as_mapping@ */
1595   0,                                    /* @tp_hash@ */
1596   0,                                    /* @tp_call@ */
1597   0,                                    /* @tp_str@ */
1598   0,                                    /* @tp_getattro@ */
1599   0,                                    /* @tp_setattro@ */
1600   0,                                    /* @tp_as_buffer@ */
1601   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1602     Py_TPFLAGS_BASETYPE,
1603
1604   /* @tp_doc@ */
1605 "A reduction context for reduction modulo primes of special form.",
1606
1607   0,                                    /* @tp_traverse@ */
1608   0,                                    /* @tp_clear@ */
1609   0,                                    /* @tp_richcompare@ */
1610   0,                                    /* @tp_weaklistoffset@ */
1611   0,                                    /* @tp_iter@ */
1612   0,                                    /* @tp_iternext@ */
1613   mpreduce_pymethods,                   /* @tp_methods@ */
1614   0,                                    /* @tp_members@ */
1615   mpreduce_pygetset,                    /* @tp_getset@ */
1616   0,                                    /* @tp_base@ */
1617   0,                                    /* @tp_dict@ */
1618   0,                                    /* @tp_descr_get@ */
1619   0,                                    /* @tp_descr_set@ */
1620   0,                                    /* @tp_dictoffset@ */
1621   0,                                    /* @tp_init@ */
1622   PyType_GenericAlloc,                  /* @tp_alloc@ */
1623   mpreduce_pynew,                       /* @tp_new@ */
1624   0,                                    /* @tp_free@ */
1625   0                                     /* @tp_is_gc@ */
1626 };
1627
1628 /*----- Chinese Remainder Theorem solution --------------------------------*/
1629
1630 typedef struct mpcrt_pyobj {
1631   PyObject_HEAD
1632   mpcrt c;
1633 } mpcrt_pyobj;
1634
1635 #define MPCRT_PY(o) (&((mpcrt_pyobj *)(o))->c)
1636
1637 static PyObject *mcmeth_solve(PyObject *me, PyObject *arg)
1638 {
1639   mpcrt *c = MPCRT_PY(me);
1640   PyObject *q = 0, *x, *z = 0;
1641   mp *xx;
1642   mp **v = 0;
1643   int i = 0, n = c->k;
1644
1645   Py_INCREF(me);
1646   if (PyTuple_Size(arg) == n)
1647     q = arg;
1648   else if (!PyArg_ParseTuple(arg, "O:solve", &q))
1649     goto end;
1650   Py_INCREF(q);
1651   if (!PySequence_Check(q)) TYERR("want a sequence of residues");
1652   if (PySequence_Size(q) != n) VALERR("residue count mismatch");
1653   v = xmalloc(n * sizeof(*v));
1654   for (i = 0; i < n; i++) {
1655     if ((x = PySequence_GetItem(q, i)) == 0) goto end;
1656     xx = getmp(x); Py_DECREF(x); if (!xx) goto end;
1657     v[i] = xx;
1658   }
1659   z = mp_pywrap(mpcrt_solve(c, MP_NEW, v));
1660 end:
1661   if (v) {
1662     n = i;
1663     for (i = 0; i < n; i++)
1664       MP_DROP(v[i]);
1665     xfree(v);
1666   }
1667   Py_DECREF(me);
1668   Py_XDECREF(q);
1669   return (z);
1670 }
1671
1672 static void mpcrt_pydealloc(PyObject *me)
1673 {
1674   mpcrt *c = MPCRT_PY(me);
1675   mpcrt_destroy(c);
1676   xfree(c->v);
1677 }
1678
1679 static PyObject *mpcrt_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1680 {
1681   mpcrt_mod *v = 0;
1682   int n, i = 0;
1683   char *kwlist[] = { "mv", 0 };
1684   PyObject *q = 0, *x;
1685   mp *xx;
1686   mpcrt_pyobj *c = 0;
1687
1688   if (PyTuple_Size(arg) > 1)
1689     q = arg;
1690   else if (!PyArg_ParseTupleAndKeywords(arg, kw, "O:new", kwlist, &q))
1691     goto end;
1692   Py_INCREF(q);
1693   if (!PySequence_Check(q)) TYERR("want a sequence of moduli");
1694   n = PySequence_Size(q);
1695   if (PyErr_Occurred()) goto end;
1696   if (!n) VALERR("want at least one modulus");
1697   v = xmalloc(n * sizeof(*v));
1698   for (i = 0; i < n; i++) {
1699     if ((x = PySequence_GetItem(q, i)) == 0) goto end;
1700     xx = getmp(x); Py_DECREF(x); if (!xx) goto end;
1701     v[i].m = xx; v[i].n = 0; v[i].ni = 0; v[i].nni = 0;
1702   }
1703   c = (mpcrt_pyobj *)ty->tp_alloc(ty, 0);
1704   mpcrt_create(&c->c, v, n, 0);
1705   Py_DECREF(q);
1706   return ((PyObject *)c);
1707
1708 end:
1709   if (v) {
1710     n = i;
1711     for (i = 0; i < n; i++)
1712       MP_DROP(v[i].m);
1713     xfree(v);
1714   }
1715   Py_XDECREF(q);
1716   return (0);
1717 }
1718
1719 static PyObject *mcget_product(PyObject *me, void *hunoz)
1720   { return (mp_pywrap(MP_COPY(MPCRT_PY(me)->mb.m))); }
1721
1722 static PyObject *mcget_moduli(PyObject *me, void *hunoz)
1723 {
1724   int i;
1725   PyObject *q;
1726   mpcrt *c = MPCRT_PY(me);
1727
1728   if ((q = PyList_New(c->k)) == 0) return (0);
1729   for (i = 0; i < c->k; i++)
1730     PyList_SetItem(q, i, mp_pywrap(c->v[i].m));
1731   return (q);
1732 }
1733
1734 static PyGetSetDef mpcrt_pygetset[] = {
1735 #define GETSETNAME(op, name) mc##op##_##name
1736   GET   (product,       "C.product -> product of moduli")
1737   GET   (moduli,        "C.moduli -> list of individual moduli")
1738 #undef GETSETNAME
1739   { 0 }
1740 };
1741
1742 static PyMethodDef mpcrt_pymethods[] = {
1743 #define METHNAME(name) mcmeth_##name
1744   METH  (solve,         "C.solve([R0, R1]) -> X mod C.product")
1745 #undef METHNAME
1746   { 0 }
1747 };
1748
1749 static PyTypeObject *mpcrt_pytype, mpcrt_pytype_skel = {
1750   PyObject_HEAD_INIT(0) 0,              /* Header */
1751   "catacomb.MPCRT",                     /* @tp_name@ */
1752   sizeof(mpcrt_pyobj),                  /* @tp_basicsize@ */
1753   0,                                    /* @tp_itemsize@ */
1754
1755   mpcrt_pydealloc,                      /* @tp_dealloc@ */
1756   0,                                    /* @tp_print@ */
1757   0,                                    /* @tp_getattr@ */
1758   0,                                    /* @tp_setattr@ */
1759   0,                                    /* @tp_compare@ */
1760   0,                                    /* @tp_repr@ */
1761   0,                                    /* @tp_as_number@ */
1762   0,                                    /* @tp_as_sequence@ */
1763   0,                                    /* @tp_as_mapping@ */
1764   0,                                    /* @tp_hash@ */
1765   0,                                    /* @tp_call@ */
1766   0,                                    /* @tp_str@ */
1767   0,                                    /* @tp_getattro@ */
1768   0,                                    /* @tp_setattro@ */
1769   0,                                    /* @tp_as_buffer@ */
1770   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1771     Py_TPFLAGS_BASETYPE,
1772
1773   /* @tp_doc@ */
1774 "A context for the solution of Chinese Remainder Theorem problems.",
1775
1776   0,                                    /* @tp_traverse@ */
1777   0,                                    /* @tp_clear@ */
1778   0,                                    /* @tp_richcompare@ */
1779   0,                                    /* @tp_weaklistoffset@ */
1780   0,                                    /* @tp_iter@ */
1781   0,                                    /* @tp_iternext@ */
1782   mpcrt_pymethods,                      /* @tp_methods@ */
1783   0,                                    /* @tp_members@ */
1784   mpcrt_pygetset,                       /* @tp_getset@ */
1785   0,                                    /* @tp_base@ */
1786   0,                                    /* @tp_dict@ */
1787   0,                                    /* @tp_descr_get@ */
1788   0,                                    /* @tp_descr_set@ */
1789   0,                                    /* @tp_dictoffset@ */
1790   0,                                    /* @tp_init@ */
1791   PyType_GenericAlloc,                  /* @tp_alloc@ */
1792   mpcrt_pynew,                          /* @tp_new@ */
1793   0,                                    /* @tp_free@ */
1794   0                                     /* @tp_is_gc@ */
1795 };
1796
1797 /*----- Binary polynomials ------------------------------------------------*/
1798
1799 static PyObject *gf_pyrepr(PyObject *o)
1800   { return mp_topystring(MP_X(o), 16, "GF(", "0x", "L)"); }
1801
1802 static PyObject *gf_pyrichcompare(PyObject *x, PyObject *y, int op)
1803 {
1804   mp *xx, *yy;
1805   int xl, yl;
1806   int b;
1807
1808   if (gfbinop(x, y, &xx, &yy)) RETURN_NOTIMPL;
1809   switch (op) {
1810     case Py_EQ: b = MP_EQ(xx, yy); break;
1811     case Py_NE: b = !MP_EQ(xx, yy); break;
1812     default:
1813       xl = mp_bits(xx);
1814       yl = mp_bits(yy);
1815       switch (op) {
1816         case Py_LT: b = xl < yl; break;
1817         case Py_LE: b = xl <= yl; break;
1818         case Py_GT: b = xl > yl; break;
1819         case Py_GE: b = xl >= yl; break;
1820         default: abort();
1821       }
1822       break;
1823   }
1824   MP_DROP(xx); MP_DROP(yy);
1825   return (getbool(b));
1826 }
1827
1828 static PyObject *gf_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1829 {
1830   PyObject *x;
1831   mp *z;
1832   mp_pyobj *zz = 0;
1833   int radix = 0;
1834   char *kwlist[] = { "x", "radix", 0 };
1835
1836   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O|i:gf", kwlist, &x, &radix))
1837     goto end;
1838   if (GF_PYCHECK(x)) RETURN_OBJ(x);
1839   if (!good_radix_p(radix, 1)) VALERR("radix out of range");
1840   if ((z = mp_frompyobject(x, radix)) == 0) {
1841     PyErr_Format(PyExc_TypeError, "can't convert %.100s to gf",
1842                  x->ob_type->tp_name);
1843     goto end;
1844   }
1845   if (MP_NEGP(z)) {
1846     MP_DROP(z);
1847     VALERR("gf cannot be negative");
1848   }
1849   zz = (mp_pyobj *)ty->tp_alloc(ty, 0);
1850   zz->x = z;
1851 end:
1852   return ((PyObject *)zz);
1853 }
1854
1855 static long gf_pyhash(PyObject *me)
1856 {
1857   long i = mp_tolong(MP_X(me));
1858   i ^= 0xc7ecd67c; /* random perturbance */
1859   if (i == -1)
1860     i = -2;
1861   return (i);
1862 }
1863
1864 static PyObject *gf_pyexp(PyObject *x, PyObject *y, PyObject *z)
1865 {
1866   mp *xx = 0, *yy = 0, *zz = 0;
1867   mp *r = 0;
1868   PyObject *rc = 0;
1869
1870   if ((xx = tomp(x)) == 0 || (yy = tomp(y)) == 0 ||
1871       (z && z != Py_None && (zz = tomp(z)) == 0)) {
1872     mp_drop(xx); mp_drop(yy); mp_drop(zz);
1873     RETURN_NOTIMPL;
1874   }
1875   if (!z || z == Py_None) {
1876     if (MP_NEGP(yy)) VALERR("negative exponent");
1877     r = gf_exp(MP_NEW, xx, yy);
1878   } else {
1879     gfreduce gr;
1880     if (MP_ZEROP(zz)) ZDIVERR("zero modulus");
1881     if (MP_NEGP(yy)) {
1882       if ((xx = gf_modinv_checked(xx, xx, zz)) == 0) goto end;
1883       yy = mp_neg(yy, yy);
1884     }
1885     gfreduce_create(&gr, zz);
1886     r = gfreduce_exp(&gr, MP_NEW, xx, yy);
1887     gfreduce_destroy(&gr);
1888   }
1889   rc = gf_pywrap(r);
1890 end:
1891   mp_drop(xx); mp_drop(yy); mp_drop(zz);
1892   return (rc);
1893 }
1894
1895 static PyObject *gfmeth_sqr(PyObject *me, PyObject *arg)
1896 {
1897   if (!PyArg_ParseTuple(arg, ":sqr")) return (0);
1898   return (gf_pywrap(gf_sqr(MP_NEW, MP_X(me))));
1899 }
1900
1901 static PyObject *gfmeth_gcd(PyObject *me, PyObject *arg)
1902 {
1903   PyObject *z = 0;
1904   mp *yy = 0, *zz = MP_NEW;
1905
1906   if (!PyArg_ParseTuple(arg, "O&:gcd", convgf, &yy)) goto end;
1907   gf_gcd(&zz, 0, 0, MP_X(me), yy);
1908   z = gf_pywrap(zz);
1909 end:
1910   if (yy) MP_DROP(yy);
1911   return (z);
1912 }
1913
1914 static PyObject *gfmeth_gcdx(PyObject *me, PyObject *arg)
1915 {
1916   PyObject *z = 0;
1917   mp *yy = 0, *zz = MP_NEW, *uu = MP_NEW, *vv = MP_NEW;
1918
1919   if (!PyArg_ParseTuple(arg, "O&:gcdx", convgf, &yy))
1920     goto end;
1921   gf_gcd(&zz, &uu, &vv, MP_X(me), yy);
1922   z = Py_BuildValue("(NNN)",
1923                     gf_pywrap(zz), gf_pywrap(uu), gf_pywrap(vv));
1924 end:
1925   if (yy) MP_DROP(yy);
1926   return (z);
1927 }
1928
1929 static PyObject *gfmeth_modinv(PyObject *me, PyObject *arg)
1930 {
1931   PyObject *z = 0;
1932   mp *yy = 0, *zz = MP_NEW;
1933
1934   if (!PyArg_ParseTuple(arg, "O&:modinv", convgf, &yy) ||
1935       (zz = gf_modinv_checked(MP_NEW, yy, MP_X(me))) == 0)
1936     goto end;
1937   z = gf_pywrap(zz);
1938 end:
1939   if (yy) MP_DROP(yy);
1940   return (z);
1941 }
1942
1943 static PyObject *gfmeth_irreduciblep(PyObject *me, PyObject *arg)
1944 {
1945   if (!PyArg_ParseTuple(arg, ":irreduciblep")) return (0);
1946   return getbool(gf_irreduciblep(MP_X(me)));
1947 }
1948
1949 static PyObject *gfget_degree(PyObject *me, void *hunoz)
1950   { return (PyInt_FromLong(mp_bits(MP_X(me)) - 1)); }
1951
1952 static PyGetSetDef gf_pygetset[] = {
1953 #define GETSETNAME(op, name) gf##op##_##name
1954   GET   (degree,        "X.degree -> polynomial degree of X")
1955 #undef GETSETNAME
1956 #define GETSETNAME(op, name) mp##op##_##name
1957   GET   (nbits,         "X.nbits -> bit length of X")
1958   GET   (noctets,       "X.noctets -> octet length of X")
1959 #undef GETSETNAME
1960   { 0 }
1961 };
1962
1963 static PyMethodDef gf_pymethods[] = {
1964 #define METHNAME(func) gfmeth_##func
1965   METH  (setbit,        "X.setbit(N) -> X with bit N set")
1966   METH  (clearbit,      "X.clearbit(N) -> X with bit N clear")
1967   METH  (testbit,       "X.testbit(N) -> true/false if bit N set/clear in X")
1968   METH  (sqr,           "X.sqr() -> X^2")
1969   METH  (gcd,           "X.gcd(Y) -> gcd(X, Y)")
1970   METH  (gcdx,
1971          "X.gcdx(Y) -> (gcd(X, Y), U, V) with X U + Y V = gcd(X, Y)")
1972   METH  (modinv,        "X.modinv(Y) -> multiplicative inverse of Y mod X")
1973   METH  (irreduciblep,  "X.irreduciblep() -> true/false")
1974 #undef METHNAME
1975 #define METHNAME(func) mpmeth_##func
1976   KWMETH(tostring,      "X.tostring(radix = 10) -> STR")
1977   KWMETH(storel,        "X.storel(len = -1) -> little-endian bytes")
1978   KWMETH(storeb,        "X.storeb(len = -1) -> big-endian bytes")
1979   KWMETH(storel2c,
1980          "X.storel2c(len = -1) -> little-endian bytes, two's complement")
1981   KWMETH(storeb2c,
1982          "X.storeb2c(len = -1) -> big-endian bytes, two's complement")
1983   METH  (tobuf,         "X.tobuf() -> buffer format")
1984 #undef METHNAME
1985   { 0 }
1986 };
1987
1988 static PyNumberMethods gf_pynumber = {
1989   gf_pyadd,                             /* @nb_add@ */
1990   gf_pysub,                             /* @nb_subtract@ */
1991   gf_pymul,                             /* @nb_multiply@ */
1992   0,                                    /* @nb_divide@ */
1993   gf_pymod,                             /* @nb_remainder@ */
1994   gf_pydivmod,                          /* @nb_divmod@ */
1995   gf_pyexp,                             /* @nb_power@ */
1996   mp_pyid,                              /* @nb_negative@ */
1997   mp_pyid,                              /* @nb_positive@ */
1998   mp_pyid,                              /* @nb_absolute@ */
1999   mp_pynonzerop,                        /* @nb_nonzero@ */
2000   0 /* doesn't make any sense */,       /* @nb_invert@ */
2001   gf_pylsl,                             /* @nb_lshift@ */
2002   gf_pylsr,                             /* @nb_rshift@ */
2003   gf_pyand,                             /* @nb_and@ */
2004   gf_pyxor,                             /* @nb_xor@ */
2005   gf_pyor,                              /* @nb_or@ */
2006   gf_pycoerce,                          /* @nb_coerce@ */
2007   mp_pyint,                             /* @nb_int@ */
2008   mp_pylong,                            /* @nb_long@ */
2009   0 /* doesn't make any sense */,       /* @nb_float@ */
2010   mp_pyoct,                             /* @nb_oct@ */
2011   mp_pyhex,                             /* @nb_hex@ */
2012
2013   0,                                    /* @nb_inplace_add@ */
2014   0,                                    /* @nb_inplace_subtract@ */
2015   0,                                    /* @nb_inplace_multiply@ */
2016   0,                                    /* @nb_inplace_divide@ */
2017   0,                                    /* @nb_inplace_remainder@ */
2018   0,                                    /* @nb_inplace_power@ */
2019   0,                                    /* @nb_inplace_lshift@ */
2020   0,                                    /* @nb_inplace_rshift@ */
2021   0,                                    /* @nb_inplace_and@ */
2022   0,                                    /* @nb_inplace_xor@ */
2023   0,                                    /* @nb_inplace_or@ */
2024
2025   gf_pydiv,                             /* @nb_floor_divide@ */
2026   0,                                    /* @nb_true_divide@ */
2027   0,                                    /* @nb_inplace_floor_divide@ */
2028   0,                                    /* @nb_inplace_true_divide@ */
2029 };
2030
2031 static PyTypeObject gf_pytype_skel = {
2032   PyObject_HEAD_INIT(0) 0,              /* Header */
2033   "catacomb.GF",                        /* @tp_name@ */
2034   sizeof(mp_pyobj),                     /* @tp_basicsize@ */
2035   0,                                    /* @tp_itemsize@ */
2036
2037   mp_pydealloc,                         /* @tp_dealloc@ */
2038   0,                                    /* @tp_print@ */
2039   0,                                    /* @tp_getattr@ */
2040   0,                                    /* @tp_setattr@ */
2041   0,                                    /* @tp_compare@ */
2042   gf_pyrepr,                            /* @tp_repr@ */
2043   &gf_pynumber,                         /* @tp_as_number@ */
2044   0,                                    /* @tp_as_sequence@ */
2045   0,                                    /* @tp_as_mapping@ */
2046   gf_pyhash,                            /* @tp_hash@ */
2047   0,                                    /* @tp_call@ */
2048   mp_pyhex,                             /* @tp_str@ */
2049   0,                                    /* @tp_getattro@ */
2050   0,                                    /* @tp_setattro@ */
2051   0,                                    /* @tp_as_buffer@ */
2052   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
2053     Py_TPFLAGS_CHECKTYPES |
2054     Py_TPFLAGS_BASETYPE,
2055
2056   /* @tp_doc@ */
2057 "Binary polynomials.  Support almost all the standard arithmetic\n\
2058 operations.\n\
2059 \n\
2060 Constructor gf(X, radix = R) attempts to convert X to a `gf'.  If\n\
2061 X is a string, it's read in radix-R form, or we look for a prefix\n\
2062 if R = 0.  Other acceptable things are ints and longs.\n\
2063 \n\
2064 The name is hopelessly wrong from a technical point of view, but\n\
2065 but it's much easier to type than `p2' or `c2' or whatever.\n\
2066 \n\
2067 Notes:\n\
2068 \n\
2069   * Use `//' for division.  GFs don't have `/' division.",
2070
2071   0,                                    /* @tp_traverse@ */
2072   0,                                    /* @tp_clear@ */
2073   gf_pyrichcompare,                     /* @tp_richcompare@ */
2074   0,                                    /* @tp_weaklistoffset@ */
2075   0,                                    /* @tp_iter@ */
2076   0,                                    /* @tp_iternext@ */
2077   gf_pymethods,                         /* @tp_methods@ */
2078   0,                                    /* @tp_members@ */
2079   gf_pygetset,                          /* @tp_getset@ */
2080   0,                                    /* @tp_base@ */
2081   0,                                    /* @tp_dict@ */
2082   0,                                    /* @tp_descr_get@ */
2083   0,                                    /* @tp_descr_set@ */
2084   0,                                    /* @tp_dictoffset@ */
2085   0,                                    /* @tp_init@ */
2086   PyType_GenericAlloc,                  /* @tp_alloc@ */
2087   gf_pynew,                             /* @tp_new@ */
2088   0,                                    /* @tp_free@ */
2089   0                                     /* @tp_is_gc@ */
2090 };
2091
2092 static PyObject *meth__GF_fromstring(PyObject *me,
2093                                     PyObject *arg, PyObject *kw)
2094 {
2095   int r = 0;
2096   char *p;
2097   int len;
2098   PyObject *z = 0;
2099   mp *zz;
2100   mptext_stringctx sc;
2101   char *kwlist[] = { "class", "x", "radix", 0 };
2102
2103   if (!PyArg_ParseTupleAndKeywords(arg, kw, "Os#|i:fromstring",
2104                                    kwlist, &me, &p, &len, &r))
2105     goto end;
2106   if (!good_radix_p(r, 1)) VALERR("bad radix");
2107   sc.buf = p; sc.lim = p + len;
2108   if ((zz = mp_read(MP_NEW, r, &mptext_stringops, &sc)) == 0 ||
2109       MP_NEGP(zz)) {
2110     if (zz) MP_DROP(zz);
2111     SYNERR("bad binary polynomial");
2112   }
2113   z = Py_BuildValue("(Ns#)", gf_pywrap(zz), sc.buf, (int)(sc.lim - sc.buf));
2114 end:
2115   return (z);
2116 }
2117
2118 /*----- Sparse poly reduction ---------------------------------------------*/
2119
2120 typedef struct gfreduce_pyobj {
2121   PyObject_HEAD
2122   gfreduce mr;
2123 } gfreduce_pyobj;
2124
2125 #define GFREDUCE_PY(o) (&((gfreduce_pyobj *)(o))->mr)
2126
2127 static PyObject *grmeth_exp(PyObject *me, PyObject *arg)
2128 {
2129   PyObject *rc = 0;
2130   mp *yy = 0, *zz = 0;
2131
2132   if (!PyArg_ParseTuple(arg, "O&O&:exp", convgf, &yy, convgf, &zz))
2133     goto end;
2134   if (MP_NEGP(zz)) {
2135     if ((yy = gf_modinv_checked(yy, yy, GFREDUCE_PY(me)->p)) == 0) goto end;
2136     zz = mp_neg(zz, zz);
2137   }
2138   rc = gf_pywrap(gfreduce_exp(GFREDUCE_PY(me), MP_NEW, yy, zz));
2139 end:
2140   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
2141   return (rc);
2142 }
2143
2144 static PyObject *grmeth_reduce(PyObject *me, PyObject *arg)
2145 {
2146   PyObject *z = 0;
2147   mp *yy = 0;
2148
2149   if (!PyArg_ParseTuple(arg, "O&:reduce", convgf, &yy)) goto end;
2150   z = gf_pywrap(gfreduce_do(GFREDUCE_PY(me), MP_NEW, yy));
2151 end:
2152   return (z);
2153 }
2154
2155 static void gfreduce_pydealloc(PyObject *me)
2156 {
2157   gfreduce_destroy(GFREDUCE_PY(me));
2158   FREEOBJ(me);
2159 }
2160
2161 static PyObject *gfreduce_pynew(PyTypeObject *ty,
2162                                  PyObject *arg, PyObject *kw)
2163 {
2164   gfreduce_pyobj *mr = 0;
2165   gfreduce r;
2166   char *kwlist[] = { "m", 0 };
2167   mp *xx = 0;
2168
2169   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convgf, &xx))
2170     goto end;
2171   if (MP_ZEROP(xx)) ZDIVERR("modulus is zero!");
2172   gfreduce_create(&r, xx);
2173   mr = (gfreduce_pyobj *)ty->tp_alloc(ty, 0);
2174   mr->mr = r;
2175 end:
2176   if (xx) MP_DROP(xx);
2177   return ((PyObject *)mr);
2178 }
2179
2180 static PyObject *grget_m(PyObject *me, void *hunoz)
2181   { return (gf_pywrap(MP_COPY(GFREDUCE_PY(me)->p))); }
2182
2183 static PyGetSetDef gfreduce_pygetset[] = {
2184 #define GETSETNAME(op, name) gr##op##_##name
2185   GET   (m,             "R.m -> reduction polynomial")
2186 #undef GETSETNAME
2187   { 0 }
2188 };
2189
2190 static PyMethodDef gfreduce_pymethods[] = {
2191 #define METHNAME(name) grmeth_##name
2192   METH  (reduce,        "R.reduce(X) -> X mod B.m")
2193   METH  (exp,           "R.exp(X, N) -> X^N mod B.m")
2194 #undef METHNAME
2195   { 0 }
2196 };
2197
2198 static PyTypeObject *gfreduce_pytype, gfreduce_pytype_skel = {
2199   PyObject_HEAD_INIT(0) 0,              /* Header */
2200   "catacomb.GFReduce",                  /* @tp_name@ */
2201   sizeof(gfreduce_pyobj),               /* @tp_basicsize@ */
2202   0,                                    /* @tp_itemsize@ */
2203
2204   gfreduce_pydealloc,                   /* @tp_dealloc@ */
2205   0,                                    /* @tp_print@ */
2206   0,                                    /* @tp_getattr@ */
2207   0,                                    /* @tp_setattr@ */
2208   0,                                    /* @tp_compare@ */
2209   0,                                    /* @tp_repr@ */
2210   0,                                    /* @tp_as_number@ */
2211   0,                                    /* @tp_as_sequence@ */
2212   0,                                    /* @tp_as_mapping@ */
2213   0,                                    /* @tp_hash@ */
2214   0,                                    /* @tp_call@ */
2215   0,                                    /* @tp_str@ */
2216   0,                                    /* @tp_getattro@ */
2217   0,                                    /* @tp_setattro@ */
2218   0,                                    /* @tp_as_buffer@ */
2219   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
2220     Py_TPFLAGS_BASETYPE,
2221
2222   /* @tp_doc@ */
2223 "A reduction context for reduction modulo sparse irreducible polynomials.",
2224
2225   0,                                    /* @tp_traverse@ */
2226   0,                                    /* @tp_clear@ */
2227   0,                                    /* @tp_richcompare@ */
2228   0,                                    /* @tp_weaklistoffset@ */
2229   0,                                    /* @tp_iter@ */
2230   0,                                    /* @tp_iternext@ */
2231   gfreduce_pymethods,                   /* @tp_methods@ */
2232   0,                                    /* @tp_members@ */
2233   gfreduce_pygetset,                    /* @tp_getset@ */
2234   0,                                    /* @tp_base@ */
2235   0,                                    /* @tp_dict@ */
2236   0,                                    /* @tp_descr_get@ */
2237   0,                                    /* @tp_descr_set@ */
2238   0,                                    /* @tp_dictoffset@ */
2239   0,                                    /* @tp_init@ */
2240   PyType_GenericAlloc,                  /* @tp_alloc@ */
2241   gfreduce_pynew,                       /* @tp_new@ */
2242   0,                                    /* @tp_free@ */
2243   0                                     /* @tp_is_gc@ */
2244 };
2245
2246 /*----- Normal/poly transformation ----------------------------------------*/
2247
2248 typedef struct gfn_pyobj {
2249   PyObject_HEAD
2250   mp *p;
2251   gfn ntop, pton;
2252 } gfn_pyobj;
2253
2254 static PyTypeObject *gfn_pytype, gfn_pytype_skel;
2255
2256 #define GFN_P(o) (((gfn_pyobj *)(o))->p)
2257 #define GFN_PTON(o) (&((gfn_pyobj *)(o))->pton)
2258 #define GFN_NTOP(o) (&((gfn_pyobj *)(o))->ntop)
2259
2260 static PyObject *gfn_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
2261 {
2262   mp *p = 0, *beta = 0;
2263   gfn_pyobj *gg = 0;
2264   char *kwlist[] = { "p", "beta", 0 };
2265
2266   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&O&:new", kwlist,
2267                                    convgf, &p, convgf, &beta))
2268     goto end;
2269   gg = PyObject_New(gfn_pyobj, ty);
2270   if (gfn_create(p, beta, &gg->ntop, &gg->pton)) {
2271     FREEOBJ(gg);
2272     gg = 0;
2273     VALERR("can't invert transformation matrix");
2274   }
2275   gg->p = MP_COPY(p);
2276 end:
2277   mp_drop(p);
2278   mp_drop(beta);
2279   return ((PyObject *)gg);
2280 }
2281
2282 static PyObject *gfnget_p(PyObject *me, void *hunoz)
2283   { return (gf_pywrap(MP_COPY(GFN_P(me)))); }
2284
2285 static PyObject *gfnget_beta(PyObject *me, void *hunoz)
2286 {
2287   gfn *n = GFN_NTOP(me);
2288   mp *x = n->r[n->n - 1];
2289   return (gf_pywrap(MP_COPY(x)));
2290 }
2291
2292 #define XFORMOP(name, NAME)                                             \
2293   static PyObject *gfnmeth_##name(PyObject *me, PyObject *arg)          \
2294   {                                                                     \
2295     mp *xx = 0;                                                         \
2296     mp *z = 0;                                                          \
2297                                                                         \
2298     if (!PyArg_ParseTuple(arg, "O&:" #name, convgf, &xx)) goto end;     \
2299     z = gfn_transform(GFN_##NAME(me), MP_NEW, xx);                      \
2300   end:                                                                  \
2301     mp_drop(xx);                                                        \
2302     if (!z) return (0);                                                 \
2303     return (mp_pywrap(z));                                              \
2304   }
2305 XFORMOP(pton, PTON)
2306 XFORMOP(ntop, NTOP)
2307 #undef XFORMOP
2308
2309 static void gfn_pydealloc(PyObject *me)
2310 {
2311   gfn_destroy(GFN_PTON(me));
2312   gfn_destroy(GFN_NTOP(me));
2313   FREEOBJ(me);
2314 }
2315
2316 static PyGetSetDef gfn_pygetset[] = {
2317 #define GETSETNAME(op, name) gfn##op##_##name
2318   GET   (p,             "X.p -> polynomial basis, as polynomial")
2319   GET   (beta,          "X.beta -> normal basis element, in poly form")
2320 #undef GETSETNAME
2321   { 0 }
2322 };
2323
2324 static PyMethodDef gfn_pymethods[] = {
2325 #define METHNAME(name) gfnmeth_##name
2326   METH  (pton,          "X.pton(A) -> normal-basis representation of A")
2327   METH  (ntop,          "X.ntop(A) -> polynomial-basis representation of A")
2328 #undef METHNAME
2329   { 0 }
2330 };
2331
2332 static PyTypeObject gfn_pytype_skel = {
2333   PyObject_HEAD_INIT(0) 0,              /* Header */
2334   "catacomb.GFN",                       /* @tp_name@ */
2335   sizeof(gfn_pyobj),                    /* @tp_basicsize@ */
2336   0,                                    /* @tp_itemsize@ */
2337
2338   gfn_pydealloc,                        /* @tp_dealloc@ */
2339   0,                                    /* @tp_print@ */
2340   0,                                    /* @tp_getattr@ */
2341   0,                                    /* @tp_setattr@ */
2342   0,                                    /* @tp_compare@ */
2343   0,                                    /* @tp_repr@ */
2344   0,                                    /* @tp_as_number@ */
2345   0,                                    /* @tp_as_sequence@ */
2346   0,                                    /* @tp_as_mapping@ */
2347   0,                                    /* @tp_hash@ */
2348   0,                                    /* @tp_call@ */
2349   0,                                    /* @tp_str@ */
2350   0,                                    /* @tp_getattro@ */
2351   0,                                    /* @tp_setattro@ */
2352   0,                                    /* @tp_as_buffer@ */
2353   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
2354     Py_TPFLAGS_BASETYPE,
2355
2356   /* @tp_doc@ */
2357 "An object for transforming elements of binary fields between polynomial\n\
2358 and normal basis representations.",
2359
2360   0,                                    /* @tp_traverse@ */
2361   0,                                    /* @tp_clear@ */
2362   0,                                    /* @tp_richcompare@ */
2363   0,                                    /* @tp_weaklistoffset@ */
2364   0,                                    /* @tp_iter@ */
2365   0,                                    /* @tp_iternext@ */
2366   gfn_pymethods,                        /* @tp_methods@ */
2367   0,                                    /* @tp_members@ */
2368   gfn_pygetset,                         /* @tp_getset@ */
2369   0,                                    /* @tp_base@ */
2370   0,                                    /* @tp_dict@ */
2371   0,                                    /* @tp_descr_get@ */
2372   0,                                    /* @tp_descr_set@ */
2373   0,                                    /* @tp_dictoffset@ */
2374   0,                                    /* @tp_init@ */
2375   PyType_GenericAlloc,                  /* @tp_alloc@ */
2376   gfn_pynew,                            /* @tp_new@ */
2377   0,                                    /* @tp_free@ */
2378   0                                     /* @tp_is_gc@ */
2379 };
2380
2381 /*----- Glue --------------------------------------------------------------*/
2382
2383 static PyMethodDef methods[] = {
2384 #define METHNAME(func) meth_##func
2385   KWMETH(_MP_fromstring,        "\
2386 fromstring(STR, radix = 0) -> (X, REST)\n\
2387 \n\
2388 Parse STR as a large integer, according to radix.  If radix is zero,\n\
2389 read a prefix from STR to decide radix: allow `0' for octal, `0x' for hex\n\
2390 or `R_' for other radix R.")
2391   KWMETH(_GF_fromstring,        "\
2392 fromstring(STR, radix = 0) -> (X, REST)\n\
2393 \n\
2394 Parse STR as a binary polynomial, according to radix.  If radix is zero,\n\
2395 read a prefix from STR to decide radix: allow `0' for octal, `0x' for hex\n\
2396 or `R_' for other radix R.")
2397   METH  (_MP_loadl,             "\
2398 loadl(STR) -> X: read little-endian bytes")
2399   METH  (_MP_loadb,             "\
2400 loadb(STR) -> X: read big-endian bytes")
2401   METH  (_MP_loadl2c,           "\
2402 loadl2c(STR) -> X: read little-endian bytes, two's complement")
2403   METH  (_MP_loadb2c,           "\
2404 loadb2c(STR) -> X: read big-endian bytes, two's complement")
2405   METH  (_MP_frombuf,           "\
2406 frombuf(STR) -> (X, REST): read buffer format")
2407   METH  (_GF_loadl,             "\
2408 loadl(STR) -> X: read little-endian bytes")
2409   METH  (_GF_loadb,             "\
2410 loadb(STR) -> X: read big-endian bytes")
2411   METH  (_GF_frombuf,           "\
2412 frombuf(STR) -> (X, REST): read buffer format")
2413 #undef METHNAME
2414   { 0 }
2415 };
2416
2417 void mp_pyinit(void)
2418 {
2419   INITTYPE(mp, root);
2420   INITTYPE(gf, root);
2421   INITTYPE(mpmul, root);
2422   INITTYPE(mpmont, root);
2423   INITTYPE(mpbarrett, root);
2424   INITTYPE(mpreduce, root);
2425   INITTYPE(mpcrt, root);
2426   INITTYPE(gfreduce, root);
2427   INITTYPE(gfn, root);
2428   addmethods(methods);
2429 }
2430
2431 void mp_pyinsert(PyObject *mod)
2432 {
2433   INSERT("MP", mp_pytype);
2434   INSERT("MPMul", mpmul_pytype);
2435   INSERT("MPMont", mpmont_pytype);
2436   INSERT("MPBarrett", mpbarrett_pytype);
2437   INSERT("MPReduce", mpreduce_pytype);
2438   INSERT("MPCRT", mpcrt_pytype);
2439   INSERT("GF", gf_pytype);
2440   INSERT("GFReduce", gfreduce_pytype);
2441   INSERT("GFN", gfn_pytype);
2442 }
2443
2444 /*----- That's all, folks -------------------------------------------------*/