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

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