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     int 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   int 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), sc.buf, (int)(sc.lim - sc.buf));
 926 end:
 927   return (z);
 928 }
 929
 930 static PyObject *meth__MP_factorial(PyObject *me, PyObject *arg)
 931 {
 932   unsigned long i;
 933   mp *x;
 934   if (!PyArg_ParseTuple(arg, "OO&:factorial", &me, convulong, &i))
 935     return (0);
 936   x = mp_factorial(i);
 937   return mp_pywrap(x);
 938 }
 939
 940 static PyObject *meth__MP_fibonacci(PyObject *me, PyObject *arg)
 941 {
 942   long i;
 943   mp *x;
 944   if (!PyArg_ParseTuple(arg, "Ol:fibonacci", &me, &i)) return (0);
 945   x = mp_fibonacci(i);
 946   return mp_pywrap(x);
 947 }
 948
 949 #define LOADOP(pre, py, name)                                           \
 950   static PyObject *meth__##py##_##name(PyObject *me, PyObject *arg)     \
 951   {                                                                     \
 952     char *p;                                                            \
 953     int len;                                                            \
 954     if (!PyArg_ParseTuple(arg, "Os#:" #name, &me, &p, &len)) return (0); \
 955     return (pre##_pywrap(mp_##name(MP_NEW, p, len)));                   \
 956   }
 957 LOADOP(mp, MP, loadl)
 958 LOADOP(mp, MP, loadb)
 959 LOADOP(mp, MP, loadl2c)
 960 LOADOP(mp, MP, loadb2c)
 961 LOADOP(gf, GF, loadl)
 962 LOADOP(gf, GF, loadb)
 963 #undef LOADOP
 964
 965 /*----- Products of small integers ----------------------------------------*/
 966
 967 typedef struct mpmul_pyobj {
 968   PyObject_HEAD
 969   int livep;
 970   mpmul mm;
 971 } mpmul_pyobj;
 972
 973 #define MPMUL_LIVEP(o) (((mpmul_pyobj *)(o))->livep)
 974 #define MPMUL_PY(o) (&((mpmul_pyobj *)(o))->mm)
 975
 976 static void mpmul_pydealloc(PyObject *me)
 977 {
 978   if (MPMUL_LIVEP(me))
 979     mp_drop(mpmul_done(MPMUL_PY(me)));
 980   FREEOBJ(me);
 981 }
 982
 983 static PyObject *mmmeth_factor(PyObject *me, PyObject *arg)
 984 {
 985   PyObject *q, *i;
 986   mp *x;
 987
 988   if (!MPMUL_LIVEP(me)) VALERR("MPMul object invalid");
 989   if (PyTuple_Size(arg) != 1)
 990     i = PyObject_GetIter(arg);
 991   else {
 992     if ((q = PyTuple_GetItem(arg, 0)) == 0) goto end;
 993     if ((i = PyObject_GetIter(q)) == 0) {
 994       PyErr_Clear(); /* that's ok */
 995       i = PyObject_GetIter(arg);
 996     }
 997   }
 998   if (!i) goto end;
 999   while ((q = PyIter_Next(i)) != 0) {
1000     x = getmp(q); Py_DECREF(q); if (!x) {
1001       Py_DECREF(i);
1002       goto end;
1003     }
1004     mpmul_add(MPMUL_PY(me), x);
1005     MP_DROP(x);
1006   }
1007   Py_DECREF(i);
1008   RETURN_ME;
1009 end:
1010   return (0);
1011 }
1012
1013 static PyObject *mmmeth_done(PyObject *me, PyObject *arg)
1014 {
1015   mp *x;
1016
1017   if (!PyArg_ParseTuple(arg, ":done")) goto end;
1018   if (!MPMUL_LIVEP(me)) VALERR("MPMul object invalid");
1019   x = mpmul_done(MPMUL_PY(me));
1020   MPMUL_LIVEP(me) = 0;
1021   return (mp_pywrap(x));
1022 end:
1023   return (0);
1024 }
1025
1026 static PyObject *mpmul_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1027 {
1028   mpmul_pyobj *mm;
1029
1030   if (kw) TYERR("keyword arguments not allowed here");
1031   mm = (mpmul_pyobj *)ty->tp_alloc(ty, 0);
1032   mpmul_init(&mm->mm);
1033   mm->livep = 1;
1034   if (mmmeth_factor((PyObject *)mm, arg) == 0) {
1035     Py_DECREF(mm);
1036     goto end;
1037   }
1038   return ((PyObject *)mm);
1039 end:
1040   return (0);
1041 }
1042
1043 static PyObject *mmget_livep(PyObject *me, void *hunoz)
1044   { return (getbool(MPMUL_LIVEP(me))); }
1045
1046 static PyGetSetDef mpmul_pygetset[] = {
1047 #define GETSETNAME(op, name) mm##op##_##name
1048   GET   (livep,                 "MM.livep -> flag: object still valid?")
1049 #undef GETSETNAME
1050   { 0 }
1051 };
1052
1053 static PyMethodDef mpmul_pymethods[] = {
1054 #define METHNAME(name) mmmeth_##name
1055   METH  (factor,                "MM.factor(ITERABLE) or MM.factor(I, ...)")
1056   METH  (done,                  "MM.done() -> PRODUCT")
1057 #undef METHNAME
1058   { 0 }
1059 };
1060
1061 static PyTypeObject *mpmul_pytype, mpmul_pytype_skel = {
1062   PyObject_HEAD_INIT(0) 0,              /* Header */
1063   "MPMul",                              /* @tp_name@ */
1064   sizeof(mpmul_pyobj),                  /* @tp_basicsize@ */
1065   0,                                    /* @tp_itemsize@ */
1066
1067   mpmul_pydealloc,                      /* @tp_dealloc@ */
1068   0,                                    /* @tp_print@ */
1069   0,                                    /* @tp_getattr@ */
1070   0,                                    /* @tp_setattr@ */
1071   0,                                    /* @tp_compare@ */
1072   0,                                    /* @tp_repr@ */
1073   0,                                    /* @tp_as_number@ */
1074   0,                                    /* @tp_as_sequence@ */
1075   0,                                    /* @tp_as_mapping@ */
1076   0,                                    /* @tp_hash@ */
1077   0,                                    /* @tp_call@ */
1078   0,                                    /* @tp_str@ */
1079   0,                                    /* @tp_getattro@ */
1080   0,                                    /* @tp_setattro@ */
1081   0,                                    /* @tp_as_buffer@ */
1082   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1083     Py_TPFLAGS_BASETYPE,
1084
1085   /* @tp_doc@ */
1086 "An object for multiplying many small integers.",
1087
1088   0,                                    /* @tp_traverse@ */
1089   0,                                    /* @tp_clear@ */
1090   0,                                    /* @tp_richcompare@ */
1091   0,                                    /* @tp_weaklistoffset@ */
1092   0,                                    /* @tp_iter@ */
1093   0,                                    /* @tp_iternext@ */
1094   mpmul_pymethods,                      /* @tp_methods@ */
1095   0,                                    /* @tp_members@ */
1096   mpmul_pygetset,                       /* @tp_getset@ */
1097   0,                                    /* @tp_base@ */
1098   0,                                    /* @tp_dict@ */
1099   0,                                    /* @tp_descr_get@ */
1100   0,                                    /* @tp_descr_set@ */
1101   0,                                    /* @tp_dictoffset@ */
1102   0,                                    /* @tp_init@ */
1103   PyType_GenericAlloc,                  /* @tp_alloc@ */
1104   mpmul_pynew,                          /* @tp_new@ */
1105   0,                                    /* @tp_free@ */
1106   0                                     /* @tp_is_gc@ */
1107 };
1108
1109 /*----- Montgomery reduction ----------------------------------------------*/
1110
1111 typedef struct mpmont_pyobj {
1112   PyObject_HEAD
1113   mpmont mm;
1114 } mpmont_pyobj;
1115
1116 #define MPMONT_PY(o) (&((mpmont_pyobj *)(o))->mm)
1117
1118 static PyObject *mmmeth_int(PyObject *me, PyObject *arg)
1119 {
1120   PyObject *z = 0;
1121   mp *yy = 0;
1122   mpmont *mm = MPMONT_PY(me);
1123
1124   if (!PyArg_ParseTuple(arg, "O&:in", convmp, &yy))
1125     goto end;
1126   mp_div(0, &yy, yy, mm->m);
1127   z = mp_pywrap(mpmont_mul(mm, MP_NEW, yy, mm->r2));
1128 end:
1129   if (yy) MP_DROP(yy);
1130   return (z);
1131 }
1132
1133 static PyObject *mmmeth_mul(PyObject *me, PyObject *arg)
1134 {
1135   PyObject *rc = 0;
1136   mp *yy = 0, *zz = 0;
1137
1138   if (!PyArg_ParseTuple(arg, "O&O&:mul", convmp, &yy, convmp, &zz))
1139     goto end;
1140   rc = mp_pywrap(mpmont_mul(MPMONT_PY(me), MP_NEW, yy, zz));
1141 end:
1142   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1143   return (rc);
1144 }
1145
1146 static PyObject *mmmeth_exp(PyObject *me, PyObject *arg)
1147 {
1148   PyObject *rc = 0;
1149   mp *yy = 0, *zz = 0;
1150
1151   if (!PyArg_ParseTuple(arg, "O&O&:exp", convmp, &yy, convmp, &zz))
1152     goto end;
1153   if (MP_NEGP(zz)) {
1154     if ((yy = mp_modinv_checked(yy, yy, MPMONT_PY(me)->m)) == 0) goto end;
1155     zz = mp_neg(zz, zz);
1156   }
1157   rc = mp_pywrap(mpmont_exp(MPMONT_PY(me), MP_NEW, yy, zz));
1158 end:
1159   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1160   return (rc);
1161 }
1162
1163 static PyObject *mmmeth_expr(PyObject *me, PyObject *arg)
1164 {
1165   PyObject *rc = 0;
1166   mp *yy = 0, *zz = 0;
1167
1168   if (!PyArg_ParseTuple(arg, "O&O&:expr", convmp, &yy, convmp, &zz))
1169     goto end;
1170   if (MP_NEGP(zz)) {
1171     yy = mpmont_reduce(MPMONT_PY(me), yy, yy);
1172     if ((yy = mp_modinv_checked(yy, yy, MPMONT_PY(me)->m)) == 0) goto end;
1173     yy = mpmont_mul(MPMONT_PY(me), yy, yy, MPMONT_PY(me)->r2);
1174     zz = mp_neg(zz, zz);
1175   }
1176   rc = mp_pywrap(mpmont_expr(MPMONT_PY(me), MP_NEW, yy, zz));
1177 end:
1178   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1179   return (rc);
1180 }
1181
1182 static PyObject *mm_mexpr_id(PyObject *me)
1183   { return mp_pywrap(MP_COPY(MPMONT_PY(me)->r)); }
1184
1185 static int mm_mexpr_fill(void *p, PyObject *me, PyObject *x, PyObject *y)
1186 {
1187   mp *xx = 0, *yy = 0;
1188   mp_expfactor *f = p;
1189   mpmont *mm = MPMONT_PY(me);
1190
1191   if ((xx = getmp(x)) == 0 || (yy = getmp(y)) == 0)
1192     goto fail;
1193   if (MP_NEGP(yy)) {
1194     xx = mpmont_reduce(mm, xx, xx);
1195     if ((xx = mp_modinv_checked(xx, xx, yy)) == 0)
1196       goto fail;
1197     xx = mpmont_mul(mm, xx, xx, mm->r2);
1198     yy = mp_neg(yy, yy);
1199   }
1200   f->base = xx;
1201   f->exp = yy;
1202   return (0);
1203
1204 fail:
1205   mp_drop(xx); mp_drop(yy);
1206   return (-1);
1207 }
1208
1209 static PyObject *mm_mexpr(PyObject *me, void *v, int n)
1210   { return mp_pywrap(mpmont_mexpr(MPMONT_PY(me), MP_NEW, v, n)); }
1211
1212 static void mp_mexp_drop(void *p)
1213 {
1214   mp_expfactor *f = p;
1215   mp_drop(f->base);
1216   mp_drop(f->exp);
1217 }
1218
1219 static PyObject *mmmeth_mexpr(PyObject *me, PyObject *arg)
1220 {
1221   return mexp_common(me, arg, sizeof(mp_expfactor),
1222                      mm_mexpr_id, mm_mexpr_fill, mm_mexpr, mp_mexp_drop);
1223 }
1224
1225 static PyObject *mp_mexp_id(PyObject *me)
1226   { return mp_pywrap(MP_ONE); }
1227
1228 static int mp_mexp_fill(void *p, PyObject *me, PyObject *x, PyObject *y)
1229 {
1230   mp *xx = 0, *yy = 0;
1231   mp_expfactor *f = p;
1232
1233   if ((xx = getmp(x)) == 0 || (yy = getmp(y)) == 0)
1234     goto fail;
1235   if (MP_NEGP(yy)) {
1236     if ((xx = mp_modinv_checked(xx, xx, yy)) == 0)
1237       goto fail;
1238     yy = mp_neg(yy, yy);
1239   }
1240   f->base = xx;
1241   f->exp = yy;
1242   return (0);
1243
1244 fail:
1245   mp_drop(xx); mp_drop(yy);
1246   return (-1);
1247 }
1248
1249 static PyObject *mm_mexp(PyObject *me, void *v, int n)
1250   { return mp_pywrap(mpmont_mexp(MPMONT_PY(me), MP_NEW, v, n)); }
1251
1252 static PyObject *mmmeth_mexp(PyObject *me, PyObject *arg)
1253 {
1254   return mexp_common(me, arg, sizeof(mp_expfactor),
1255                      mp_mexp_id, mp_mexp_fill, mm_mexp, mp_mexp_drop);
1256 }
1257
1258 #define mmmeth_ext mmmeth_reduce
1259 static PyObject *mmmeth_reduce(PyObject *me, PyObject *arg)
1260 {
1261   PyObject *z = 0;
1262   mp *yy = 0;
1263
1264   if (!PyArg_ParseTuple(arg, "O&", convmp, &yy)) goto end;
1265   z = mp_pywrap(mpmont_reduce(MPMONT_PY(me), MP_NEW, yy));
1266 end:
1267   return (z);
1268 }
1269
1270 static void mpmont_pydealloc(PyObject *me)
1271 {
1272   mpmont_destroy(MPMONT_PY(me));
1273   FREEOBJ(me);
1274 }
1275
1276 static PyObject *mpmont_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1277 {
1278   mpmont_pyobj *mm = 0;
1279   char *kwlist[] = { "m", 0 };
1280   mp *xx = 0;
1281
1282   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convmp, &xx))
1283     goto end;
1284   if (!MP_POSP(xx) || !MP_ODDP(xx)) VALERR("m must be positive and odd");
1285   mm = (mpmont_pyobj *)ty->tp_alloc(ty, 0);
1286   mpmont_create(&mm->mm, xx);
1287 end:
1288   if (xx) MP_DROP(xx);
1289   return ((PyObject *)mm);
1290 }
1291
1292 static PyObject *mmget_m(PyObject *me, void *hunoz)
1293   { return (mp_pywrap(MP_COPY(MPMONT_PY(me)->m))); }
1294
1295 static PyObject *mmget_r(PyObject *me, void *hunoz)
1296   { return (mp_pywrap(MP_COPY(MPMONT_PY(me)->r))); }
1297
1298 static PyObject *mmget_r2(PyObject *me, void *hunoz)
1299   { return (mp_pywrap(MP_COPY(MPMONT_PY(me)->r2))); }
1300
1301 static PyGetSetDef mpmont_pygetset[] = {
1302 #define GETSETNAME(op, name) mm##op##_##name
1303   GET   (m,             "M.m -> modulus for reduction")
1304   GET   (r,             "M.r -> multiplicative identity")
1305   GET   (r2,            "M.r2 -> M.r^2, Montgomerization factor")
1306 #undef GETSETNAME
1307   { 0 }
1308 };
1309
1310 static PyMethodDef mpmont_pymethods[] = {
1311 #define METHNAME(name) mmmeth_##name
1312   METH  (int,           "M.int(X) -> XR")
1313   METH  (mul,           "M.mul(XR, YR) -> ZR where Z = X Y")
1314   METH  (expr,          "M.expr(XR, N) -> ZR where Z = X^N mod M.m")
1315   METH  (mexpr,         "\
1316 M.mexpr([(XR0, N0), (XR1, N1), ...]) = ZR where Z = X0^N0 X1^N1 ... mod M.m\n\
1317 \t(the list may be flattened if this more convenient.)")
1318   METH  (reduce,        "M.reduce(XR) -> X")
1319   METH  (ext,           "M.ext(XR) -> X")
1320   METH  (exp,           "M.exp(X, N) -> X^N mod M.m")
1321   METH  (mexp,          "\
1322 M.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 ... mod M.m\n\
1323 \t(the list may be flattened if this more convenient.)")
1324 #undef METHNAME
1325   { 0 }
1326 };
1327
1328 static PyTypeObject *mpmont_pytype, mpmont_pytype_skel = {
1329   PyObject_HEAD_INIT(0) 0,              /* Header */
1330   "MPMont",                             /* @tp_name@ */
1331   sizeof(mpmont_pyobj),                 /* @tp_basicsize@ */
1332   0,                                    /* @tp_itemsize@ */
1333
1334   mpmont_pydealloc,                     /* @tp_dealloc@ */
1335   0,                                    /* @tp_print@ */
1336   0,                                    /* @tp_getattr@ */
1337   0,                                    /* @tp_setattr@ */
1338   0,                                    /* @tp_compare@ */
1339   0,                                    /* @tp_repr@ */
1340   0,                                    /* @tp_as_number@ */
1341   0,                                    /* @tp_as_sequence@ */
1342   0,                                    /* @tp_as_mapping@ */
1343   0,                                    /* @tp_hash@ */
1344   0,                                    /* @tp_call@ */
1345   0,                                    /* @tp_str@ */
1346   0,                                    /* @tp_getattro@ */
1347   0,                                    /* @tp_setattro@ */
1348   0,                                    /* @tp_as_buffer@ */
1349   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1350     Py_TPFLAGS_BASETYPE,
1351
1352   /* @tp_doc@ */
1353 "A Montgomery reduction context.",
1354
1355   0,                                    /* @tp_traverse@ */
1356   0,                                    /* @tp_clear@ */
1357   0,                                    /* @tp_richcompare@ */
1358   0,                                    /* @tp_weaklistoffset@ */
1359   0,                                    /* @tp_iter@ */
1360   0,                                    /* @tp_iternext@ */
1361   mpmont_pymethods,                     /* @tp_methods@ */
1362   0,                                    /* @tp_members@ */
1363   mpmont_pygetset,                      /* @tp_getset@ */
1364   0,                                    /* @tp_base@ */
1365   0,                                    /* @tp_dict@ */
1366   0,                                    /* @tp_descr_get@ */
1367   0,                                    /* @tp_descr_set@ */
1368   0,                                    /* @tp_dictoffset@ */
1369   0,                                    /* @tp_init@ */
1370   PyType_GenericAlloc,                  /* @tp_alloc@ */
1371   mpmont_pynew,                         /* @tp_new@ */
1372   0,                                    /* @tp_free@ */
1373   0                                     /* @tp_is_gc@ */
1374 };
1375
1376 /*----- Barrett reduction -------------------------------------------------*/
1377
1378 typedef struct mpbarrett_pyobj {
1379   PyObject_HEAD
1380   mpbarrett mb;
1381 } mpbarrett_pyobj;
1382
1383 #define MPBARRETT_PY(o) (&((mpbarrett_pyobj *)(o))->mb)
1384
1385 static PyObject *mbmeth_exp(PyObject *me, PyObject *arg)
1386 {
1387   PyObject *rc = 0;
1388   mp *yy = 0, *zz = 0;
1389
1390   if (!PyArg_ParseTuple(arg, "O&O&:exp", convmp, &yy, convmp, &zz))
1391     goto end;
1392   if (MP_NEGP(zz)) {
1393     if ((yy = mp_modinv_checked(yy, yy, MPBARRETT_PY(me)->m)) == 0) goto end;
1394     zz = mp_neg(zz, zz);
1395   }
1396   rc = mp_pywrap(mpbarrett_exp(MPBARRETT_PY(me), MP_NEW, yy, zz));
1397 end:
1398   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1399   return (rc);
1400 }
1401
1402 static PyObject *mb_mexp(PyObject *me, void *v, int n)
1403   { return mp_pywrap(mpbarrett_mexp(MPBARRETT_PY(me), MP_NEW, v, n)); }
1404
1405 static PyObject *mbmeth_mexp(PyObject *me, PyObject *arg)
1406 {
1407   return mexp_common(me, arg, sizeof(mp_expfactor),
1408                      mp_mexp_id, mp_mexp_fill, mb_mexp, mp_mexp_drop);
1409 }
1410
1411 static PyObject *mbmeth_reduce(PyObject *me, PyObject *arg)
1412 {
1413   PyObject *z = 0;
1414   mp *yy = 0;
1415
1416   if (!PyArg_ParseTuple(arg, "O&:reduce", convmp, &yy))
1417     goto end;
1418   z = mp_pywrap(mpbarrett_reduce(MPBARRETT_PY(me), MP_NEW, yy));
1419 end:
1420   return (z);
1421 }
1422
1423 static void mpbarrett_pydealloc(PyObject *me)
1424 {
1425   mpbarrett_destroy(MPBARRETT_PY(me));
1426   FREEOBJ(me);
1427 }
1428
1429 static PyObject *mpbarrett_pynew(PyTypeObject *ty,
1430                                  PyObject *arg, PyObject *kw)
1431 {
1432   mpbarrett_pyobj *mb = 0;
1433   char *kwlist[] = { "m", 0 };
1434   mp *xx = 0;
1435
1436   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convmp, &xx))
1437     goto end;
1438   if (!MP_POSP(xx)) VALERR("m must be positive");
1439   mb = (mpbarrett_pyobj *)ty->tp_alloc(ty, 0);
1440   mpbarrett_create(&mb->mb, xx);
1441 end:
1442   if (xx) MP_DROP(xx);
1443   return ((PyObject *)mb);
1444 }
1445
1446 static PyObject *mbget_m(PyObject *me, void *hunoz)
1447   { return (mp_pywrap(MP_COPY(MPBARRETT_PY(me)->m))); }
1448
1449 static PyGetSetDef mpbarrett_pygetset[] = {
1450 #define GETSETNAME(op, name) mb##op##_##name
1451   GET   (m,             "B.m -> modulus for reduction")
1452 #undef GETSETNAME
1453   { 0 }
1454 };
1455
1456 static PyMethodDef mpbarrett_pymethods[] = {
1457 #define METHNAME(name) mbmeth_##name
1458   METH  (reduce,        "B.reduce(X) -> X mod B.m")
1459   METH  (exp,           "B.exp(X, N) -> X^N mod B.m")
1460   METH  (mexp,          "\
1461 B.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 ... mod B.m\n\
1462 \t(the list may be flattened if this more convenient.)")
1463 #undef METHNAME
1464   { 0 }
1465 };
1466
1467 static PyTypeObject *mpbarrett_pytype, mpbarrett_pytype_skel = {
1468   PyObject_HEAD_INIT(0) 0,              /* Header */
1469   "MPBarrett",                          /* @tp_name@ */
1470   sizeof(mpbarrett_pyobj),              /* @tp_basicsize@ */
1471   0,                                    /* @tp_itemsize@ */
1472
1473   mpbarrett_pydealloc,                  /* @tp_dealloc@ */
1474   0,                                    /* @tp_print@ */
1475   0,                                    /* @tp_getattr@ */
1476   0,                                    /* @tp_setattr@ */
1477   0,                                    /* @tp_compare@ */
1478   0,                                    /* @tp_repr@ */
1479   0,                                    /* @tp_as_number@ */
1480   0,                                    /* @tp_as_sequence@ */
1481   0,                                    /* @tp_as_mapping@ */
1482   0,                                    /* @tp_hash@ */
1483   0,                                    /* @tp_call@ */
1484   0,                                    /* @tp_str@ */
1485   0,                                    /* @tp_getattro@ */
1486   0,                                    /* @tp_setattro@ */
1487   0,                                    /* @tp_as_buffer@ */
1488   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1489     Py_TPFLAGS_BASETYPE,
1490
1491   /* @tp_doc@ */
1492 "A Barrett reduction context.",
1493
1494   0,                                    /* @tp_traverse@ */
1495   0,                                    /* @tp_clear@ */
1496   0,                                    /* @tp_richcompare@ */
1497   0,                                    /* @tp_weaklistoffset@ */
1498   0,                                    /* @tp_iter@ */
1499   0,                                    /* @tp_iternext@ */
1500   mpbarrett_pymethods,                  /* @tp_methods@ */
1501   0,                                    /* @tp_members@ */
1502   mpbarrett_pygetset,                   /* @tp_getset@ */
1503   0,                                    /* @tp_base@ */
1504   0,                                    /* @tp_dict@ */
1505   0,                                    /* @tp_descr_get@ */
1506   0,                                    /* @tp_descr_set@ */
1507   0,                                    /* @tp_dictoffset@ */
1508   0,                                    /* @tp_init@ */
1509   PyType_GenericAlloc,                  /* @tp_alloc@ */
1510   mpbarrett_pynew,                      /* @tp_new@ */
1511   0,                                    /* @tp_free@ */
1512   0                                     /* @tp_is_gc@ */
1513 };
1514
1515 /*----- Nice prime reduction ----------------------------------------------*/
1516
1517 typedef struct mpreduce_pyobj {
1518   PyObject_HEAD
1519   mpreduce mr;
1520 } mpreduce_pyobj;
1521
1522 #define MPREDUCE_PY(o) (&((mpreduce_pyobj *)(o))->mr)
1523
1524 static PyObject *mrmeth_exp(PyObject *me, PyObject *arg)
1525 {
1526   PyObject *rc = 0;
1527   mp *yy = 0, *zz = 0;
1528
1529   if (!PyArg_ParseTuple(arg, "O&O&:exp", convmp, &yy, convmp, &zz))
1530     goto end;
1531   if (MP_NEGP(zz)) {
1532     if ((yy = mp_modinv_checked(yy, yy, MPREDUCE_PY(me)->p)) == 0) goto end;
1533     zz = mp_neg(zz, zz);
1534   }
1535   rc = mp_pywrap(mpreduce_exp(MPREDUCE_PY(me), MP_NEW, yy, zz));
1536 end:
1537   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
1538   return (rc);
1539 }
1540
1541 static PyObject *mrmeth_reduce(PyObject *me, PyObject *arg)
1542 {
1543   PyObject *z = 0;
1544   mp *yy = 0;
1545
1546   if (!PyArg_ParseTuple(arg, "O&:reduce", convmp, &yy)) goto end;
1547   z = mp_pywrap(mpreduce_do(MPREDUCE_PY(me), MP_NEW, yy));
1548 end:
1549   return (z);
1550 }
1551
1552 static void mpreduce_pydealloc(PyObject *me)
1553 {
1554   mpreduce_destroy(MPREDUCE_PY(me));
1555   FREEOBJ(me);
1556 }
1557
1558 static PyObject *mpreduce_pynew(PyTypeObject *ty,
1559                                  PyObject *arg, PyObject *kw)
1560 {
1561   mpreduce_pyobj *mr = 0;
1562   mpreduce r;
1563   char *kwlist[] = { "m", 0 };
1564   mp *xx = 0;
1565
1566   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convmp, &xx))
1567     goto end;
1568   if (!MP_POSP(xx)) VALERR("m must be positive");
1569   if (mpreduce_create(&r, xx)) VALERR("bad modulus (must be 2^k - ...)");
1570   mr = (mpreduce_pyobj *)ty->tp_alloc(ty, 0);
1571   mr->mr = r;
1572 end:
1573   if (xx) MP_DROP(xx);
1574   return ((PyObject *)mr);
1575 }
1576
1577 static PyObject *mrget_m(PyObject *me, void *hunoz)
1578   { return (mp_pywrap(MP_COPY(MPREDUCE_PY(me)->p))); }
1579
1580 static PyGetSetDef mpreduce_pygetset[] = {
1581 #define GETSETNAME(op, name) mr##op##_##name
1582   GET   (m,             "R.m -> modulus for reduction")
1583 #undef GETSETNAME
1584   { 0 }
1585 };
1586
1587 static PyMethodDef mpreduce_pymethods[] = {
1588 #define METHNAME(name) mrmeth_##name
1589   METH  (reduce,        "R.reduce(X) -> X mod B.m")
1590   METH  (exp,           "R.exp(X, N) -> X^N mod B.m")
1591 #undef METHNAME
1592   { 0 }
1593 };
1594
1595 static PyTypeObject *mpreduce_pytype, mpreduce_pytype_skel = {
1596   PyObject_HEAD_INIT(0) 0,              /* Header */
1597   "MPReduce",                           /* @tp_name@ */
1598   sizeof(mpreduce_pyobj),               /* @tp_basicsize@ */
1599   0,                                    /* @tp_itemsize@ */
1600
1601   mpreduce_pydealloc,                   /* @tp_dealloc@ */
1602   0,                                    /* @tp_print@ */
1603   0,                                    /* @tp_getattr@ */
1604   0,                                    /* @tp_setattr@ */
1605   0,                                    /* @tp_compare@ */
1606   0,                                    /* @tp_repr@ */
1607   0,                                    /* @tp_as_number@ */
1608   0,                                    /* @tp_as_sequence@ */
1609   0,                                    /* @tp_as_mapping@ */
1610   0,                                    /* @tp_hash@ */
1611   0,                                    /* @tp_call@ */
1612   0,                                    /* @tp_str@ */
1613   0,                                    /* @tp_getattro@ */
1614   0,                                    /* @tp_setattro@ */
1615   0,                                    /* @tp_as_buffer@ */
1616   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1617     Py_TPFLAGS_BASETYPE,
1618
1619   /* @tp_doc@ */
1620 "A reduction context for reduction modulo primes of special form.",
1621
1622   0,                                    /* @tp_traverse@ */
1623   0,                                    /* @tp_clear@ */
1624   0,                                    /* @tp_richcompare@ */
1625   0,                                    /* @tp_weaklistoffset@ */
1626   0,                                    /* @tp_iter@ */
1627   0,                                    /* @tp_iternext@ */
1628   mpreduce_pymethods,                   /* @tp_methods@ */
1629   0,                                    /* @tp_members@ */
1630   mpreduce_pygetset,                    /* @tp_getset@ */
1631   0,                                    /* @tp_base@ */
1632   0,                                    /* @tp_dict@ */
1633   0,                                    /* @tp_descr_get@ */
1634   0,                                    /* @tp_descr_set@ */
1635   0,                                    /* @tp_dictoffset@ */
1636   0,                                    /* @tp_init@ */
1637   PyType_GenericAlloc,                  /* @tp_alloc@ */
1638   mpreduce_pynew,                       /* @tp_new@ */
1639   0,                                    /* @tp_free@ */
1640   0                                     /* @tp_is_gc@ */
1641 };
1642
1643 /*----- Chinese Remainder Theorem solution --------------------------------*/
1644
1645 typedef struct mpcrt_pyobj {
1646   PyObject_HEAD
1647   mpcrt c;
1648 } mpcrt_pyobj;
1649
1650 #define MPCRT_PY(o) (&((mpcrt_pyobj *)(o))->c)
1651
1652 static PyObject *mcmeth_solve(PyObject *me, PyObject *arg)
1653 {
1654   mpcrt *c = MPCRT_PY(me);
1655   PyObject *q = 0, *x, *z = 0;
1656   mp *xx;
1657   mp **v = 0;
1658   int i = 0, n = c->k;
1659
1660   Py_INCREF(me);
1661   if (PyTuple_Size(arg) == n)
1662     q = arg;
1663   else if (!PyArg_ParseTuple(arg, "O:solve", &q))
1664     goto end;
1665   Py_INCREF(q);
1666   if (!PySequence_Check(q)) TYERR("want a sequence of residues");
1667   if (PySequence_Size(q) != n) VALERR("residue count mismatch");
1668   v = xmalloc(n * sizeof(*v));
1669   for (i = 0; i < n; i++) {
1670     if ((x = PySequence_GetItem(q, i)) == 0) goto end;
1671     xx = getmp(x); Py_DECREF(x); if (!xx) goto end;
1672     v[i] = xx;
1673   }
1674   z = mp_pywrap(mpcrt_solve(c, MP_NEW, v));
1675 end:
1676   if (v) {
1677     n = i;
1678     for (i = 0; i < n; i++)
1679       MP_DROP(v[i]);
1680     xfree(v);
1681   }
1682   Py_DECREF(me);
1683   Py_XDECREF(q);
1684   return (z);
1685 }
1686
1687 static void mpcrt_pydealloc(PyObject *me)
1688 {
1689   mpcrt *c = MPCRT_PY(me);
1690   mpcrt_destroy(c);
1691   xfree(c->v);
1692 }
1693
1694 static PyObject *mpcrt_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1695 {
1696   mpcrt_mod *v = 0;
1697   int n, i = 0;
1698   char *kwlist[] = { "mv", 0 };
1699   PyObject *q = 0, *x;
1700   mp *xx;
1701   mpcrt_pyobj *c = 0;
1702
1703   if (PyTuple_Size(arg) > 1)
1704     q = arg;
1705   else if (!PyArg_ParseTupleAndKeywords(arg, kw, "O:new", kwlist, &q))
1706     goto end;
1707   Py_INCREF(q);
1708   if (!PySequence_Check(q)) TYERR("want a sequence of moduli");
1709   n = PySequence_Size(q);
1710   if (PyErr_Occurred()) goto end;
1711   if (!n) VALERR("want at least one modulus");
1712   v = xmalloc(n * sizeof(*v));
1713   for (i = 0; i < n; i++) {
1714     if ((x = PySequence_GetItem(q, i)) == 0) goto end;
1715     xx = getmp(x); Py_DECREF(x); if (!xx) goto end;
1716     v[i].m = xx; v[i].n = 0; v[i].ni = 0; v[i].nni = 0;
1717   }
1718   c = (mpcrt_pyobj *)ty->tp_alloc(ty, 0);
1719   mpcrt_create(&c->c, v, n, 0);
1720   Py_DECREF(q);
1721   return ((PyObject *)c);
1722
1723 end:
1724   if (v) {
1725     n = i;
1726     for (i = 0; i < n; i++)
1727       MP_DROP(v[i].m);
1728     xfree(v);
1729   }
1730   Py_XDECREF(q);
1731   return (0);
1732 }
1733
1734 static PyObject *mcget_product(PyObject *me, void *hunoz)
1735   { return (mp_pywrap(MP_COPY(MPCRT_PY(me)->mb.m))); }
1736
1737 static PyObject *mcget_moduli(PyObject *me, void *hunoz)
1738 {
1739   int i;
1740   PyObject *q;
1741   mpcrt *c = MPCRT_PY(me);
1742
1743   if ((q = PyList_New(c->k)) == 0) return (0);
1744   for (i = 0; i < c->k; i++)
1745     PyList_SetItem(q, i, mp_pywrap(c->v[i].m));
1746   return (q);
1747 }
1748
1749 static PyGetSetDef mpcrt_pygetset[] = {
1750 #define GETSETNAME(op, name) mc##op##_##name
1751   GET   (product,       "C.product -> product of moduli")
1752   GET   (moduli,        "C.moduli -> list of individual moduli")
1753 #undef GETSETNAME
1754   { 0 }
1755 };
1756
1757 static PyMethodDef mpcrt_pymethods[] = {
1758 #define METHNAME(name) mcmeth_##name
1759   METH  (solve,         "C.solve([R0, R1]) -> X mod C.product")
1760 #undef METHNAME
1761   { 0 }
1762 };
1763
1764 static PyTypeObject *mpcrt_pytype, mpcrt_pytype_skel = {
1765   PyObject_HEAD_INIT(0) 0,              /* Header */
1766   "MPCRT",                              /* @tp_name@ */
1767   sizeof(mpcrt_pyobj),                  /* @tp_basicsize@ */
1768   0,                                    /* @tp_itemsize@ */
1769
1770   mpcrt_pydealloc,                      /* @tp_dealloc@ */
1771   0,                                    /* @tp_print@ */
1772   0,                                    /* @tp_getattr@ */
1773   0,                                    /* @tp_setattr@ */
1774   0,                                    /* @tp_compare@ */
1775   0,                                    /* @tp_repr@ */
1776   0,                                    /* @tp_as_number@ */
1777   0,                                    /* @tp_as_sequence@ */
1778   0,                                    /* @tp_as_mapping@ */
1779   0,                                    /* @tp_hash@ */
1780   0,                                    /* @tp_call@ */
1781   0,                                    /* @tp_str@ */
1782   0,                                    /* @tp_getattro@ */
1783   0,                                    /* @tp_setattro@ */
1784   0,                                    /* @tp_as_buffer@ */
1785   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
1786     Py_TPFLAGS_BASETYPE,
1787
1788   /* @tp_doc@ */
1789 "A context for the solution of Chinese Remainder Theorem problems.",
1790
1791   0,                                    /* @tp_traverse@ */
1792   0,                                    /* @tp_clear@ */
1793   0,                                    /* @tp_richcompare@ */
1794   0,                                    /* @tp_weaklistoffset@ */
1795   0,                                    /* @tp_iter@ */
1796   0,                                    /* @tp_iternext@ */
1797   mpcrt_pymethods,                      /* @tp_methods@ */
1798   0,                                    /* @tp_members@ */
1799   mpcrt_pygetset,                       /* @tp_getset@ */
1800   0,                                    /* @tp_base@ */
1801   0,                                    /* @tp_dict@ */
1802   0,                                    /* @tp_descr_get@ */
1803   0,                                    /* @tp_descr_set@ */
1804   0,                                    /* @tp_dictoffset@ */
1805   0,                                    /* @tp_init@ */
1806   PyType_GenericAlloc,                  /* @tp_alloc@ */
1807   mpcrt_pynew,                          /* @tp_new@ */
1808   0,                                    /* @tp_free@ */
1809   0                                     /* @tp_is_gc@ */
1810 };
1811
1812 /*----- Binary polynomials ------------------------------------------------*/
1813
1814 static PyObject *gf_pyrepr(PyObject *o)
1815   { return mp_topystring(MP_X(o), 16, "GF(", "0x", "L)"); }
1816
1817 static PyObject *gf_pyrichcompare(PyObject *x, PyObject *y, int op)
1818 {
1819   mp *xx, *yy;
1820   int xl, yl;
1821   int b;
1822
1823   if (gfbinop(x, y, &xx, &yy)) RETURN_NOTIMPL;
1824   switch (op) {
1825     case Py_EQ: b = MP_EQ(xx, yy); break;
1826     case Py_NE: b = !MP_EQ(xx, yy); break;
1827     default:
1828       xl = mp_bits(xx);
1829       yl = mp_bits(yy);
1830       switch (op) {
1831         case Py_LT: b = xl < yl; break;
1832         case Py_LE: b = xl <= yl; break;
1833         case Py_GT: b = xl > yl; break;
1834         case Py_GE: b = xl >= yl; break;
1835         default: abort();
1836       }
1837       break;
1838   }
1839   MP_DROP(xx); MP_DROP(yy);
1840   return (getbool(b));
1841 }
1842
1843 static PyObject *gf_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
1844 {
1845   PyObject *x;
1846   mp *z;
1847   mp_pyobj *zz = 0;
1848   int radix = 0;
1849   char *kwlist[] = { "x", "radix", 0 };
1850
1851   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O|i:gf", kwlist, &x, &radix))
1852     goto end;
1853   if (GF_PYCHECK(x)) RETURN_OBJ(x);
1854   if (!good_radix_p(radix, 1)) VALERR("radix out of range");
1855   if ((z = mp_frompyobject(x, radix)) == 0) {
1856     PyErr_Format(PyExc_TypeError, "can't convert %.100s to gf",
1857                  x->ob_type->tp_name);
1858     goto end;
1859   }
1860   if (MP_NEGP(z)) {
1861     MP_DROP(z);
1862     VALERR("gf cannot be negative");
1863   }
1864   zz = (mp_pyobj *)ty->tp_alloc(ty, 0);
1865   zz->x = z;
1866 end:
1867   return ((PyObject *)zz);
1868 }
1869
1870 static long gf_pyhash(PyObject *me)
1871 {
1872   long i = mp_tolong(MP_X(me));
1873   i ^= 0xc7ecd67c; /* random perturbance */
1874   if (i == -1)
1875     i = -2;
1876   return (i);
1877 }
1878
1879 static PyObject *gf_pyexp(PyObject *x, PyObject *y, PyObject *z)
1880 {
1881   mp *xx = 0, *yy = 0, *zz = 0;
1882   mp *r = 0;
1883   PyObject *rc = 0;
1884
1885   if ((xx = tomp(x)) == 0 || (yy = tomp(y)) == 0 ||
1886       (z && z != Py_None && (zz = tomp(z)) == 0)) {
1887     mp_drop(xx); mp_drop(yy); mp_drop(zz);
1888     RETURN_NOTIMPL;
1889   }
1890   if (!z || z == Py_None) {
1891     if (MP_NEGP(yy)) VALERR("negative exponent");
1892     r = gf_exp(MP_NEW, xx, yy);
1893   } else {
1894     gfreduce gr;
1895     if (MP_ZEROP(zz)) ZDIVERR("zero modulus");
1896     if (MP_NEGP(yy)) {
1897       if ((xx = gf_modinv_checked(xx, xx, zz)) == 0) goto end;
1898       yy = mp_neg(yy, yy);
1899     }
1900     gfreduce_create(&gr, zz);
1901     r = gfreduce_exp(&gr, MP_NEW, xx, yy);
1902     gfreduce_destroy(&gr);
1903   }
1904   rc = gf_pywrap(r);
1905 end:
1906   mp_drop(xx); mp_drop(yy); mp_drop(zz);
1907   return (rc);
1908 }
1909
1910 static PyObject *gfmeth_sqr(PyObject *me, PyObject *arg)
1911 {
1912   if (!PyArg_ParseTuple(arg, ":sqr")) return (0);
1913   return (gf_pywrap(gf_sqr(MP_NEW, MP_X(me))));
1914 }
1915
1916 static PyObject *gfmeth_gcd(PyObject *me, PyObject *arg)
1917 {
1918   PyObject *z = 0;
1919   mp *yy = 0, *zz = MP_NEW;
1920
1921   if (!PyArg_ParseTuple(arg, "O&:gcd", convgf, &yy)) goto end;
1922   gf_gcd(&zz, 0, 0, MP_X(me), yy);
1923   z = gf_pywrap(zz);
1924 end:
1925   if (yy) MP_DROP(yy);
1926   return (z);
1927 }
1928
1929 static PyObject *gfmeth_gcdx(PyObject *me, PyObject *arg)
1930 {
1931   PyObject *z = 0;
1932   mp *yy = 0, *zz = MP_NEW, *uu = MP_NEW, *vv = MP_NEW;
1933
1934   if (!PyArg_ParseTuple(arg, "O&:gcdx", convgf, &yy))
1935     goto end;
1936   gf_gcd(&zz, &uu, &vv, MP_X(me), yy);
1937   z = Py_BuildValue("(NNN)",
1938                     gf_pywrap(zz), gf_pywrap(uu), gf_pywrap(vv));
1939 end:
1940   if (yy) MP_DROP(yy);
1941   return (z);
1942 }
1943
1944 static PyObject *gfmeth_modinv(PyObject *me, PyObject *arg)
1945 {
1946   PyObject *z = 0;
1947   mp *yy = 0, *zz = MP_NEW;
1948
1949   if (!PyArg_ParseTuple(arg, "O&:modinv", convgf, &yy) ||
1950       (zz = gf_modinv_checked(MP_NEW, yy, MP_X(me))) == 0)
1951     goto end;
1952   z = gf_pywrap(zz);
1953 end:
1954   if (yy) MP_DROP(yy);
1955   return (z);
1956 }
1957
1958 static PyObject *gfmeth_irreduciblep(PyObject *me, PyObject *arg)
1959 {
1960   if (!PyArg_ParseTuple(arg, ":irreduciblep")) return (0);
1961   return getbool(gf_irreduciblep(MP_X(me)));
1962 }
1963
1964 static PyObject *gfget_degree(PyObject *me, void *hunoz)
1965   { return (PyInt_FromLong(mp_bits(MP_X(me)) - 1)); }
1966
1967 static PyGetSetDef gf_pygetset[] = {
1968 #define GETSETNAME(op, name) gf##op##_##name
1969   GET   (degree,        "X.degree -> polynomial degree of X")
1970 #undef GETSETNAME
1971 #define GETSETNAME(op, name) mp##op##_##name
1972   GET   (nbits,         "X.nbits -> bit length of X")
1973   GET   (noctets,       "X.noctets -> octet length of X")
1974 #undef GETSETNAME
1975   { 0 }
1976 };
1977
1978 static PyMethodDef gf_pymethods[] = {
1979 #define METHNAME(func) gfmeth_##func
1980   METH  (setbit,        "X.setbit(N) -> X with bit N set")
1981   METH  (clearbit,      "X.clearbit(N) -> X with bit N clear")
1982   METH  (testbit,       "X.testbit(N) -> true/false if bit N set/clear in X")
1983   METH  (sqr,           "X.sqr() -> X^2")
1984   METH  (gcd,           "X.gcd(Y) -> gcd(X, Y)")
1985   METH  (gcdx,
1986          "X.gcdx(Y) -> (gcd(X, Y), U, V) with X U + Y V = gcd(X, Y)")
1987   METH  (modinv,        "X.modinv(Y) -> multiplicative inverse of Y mod X")
1988   METH  (irreduciblep,  "X.irreduciblep() -> true/false")
1989 #undef METHNAME
1990 #define METHNAME(func) mpmeth_##func
1991   KWMETH(tostring,      "X.tostring(radix = 10) -> STR")
1992   KWMETH(storel,        "X.storel(len = -1) -> little-endian bytes")
1993   KWMETH(storeb,        "X.storeb(len = -1) -> big-endian bytes")
1994   KWMETH(storel2c,
1995          "X.storel2c(len = -1) -> little-endian bytes, two's complement")
1996   KWMETH(storeb2c,
1997          "X.storeb2c(len = -1) -> big-endian bytes, two's complement")
1998   METH  (tobuf,         "X.tobuf() -> buffer format")
1999 #undef METHNAME
2000   { 0 }
2001 };
2002
2003 static PyNumberMethods gf_pynumber = {
2004   gf_pyadd,                             /* @nb_add@ */
2005   gf_pysub,                             /* @nb_subtract@ */
2006   gf_pymul,                             /* @nb_multiply@ */
2007   0,                                    /* @nb_divide@ */
2008   gf_pymod,                             /* @nb_remainder@ */
2009   gf_pydivmod,                          /* @nb_divmod@ */
2010   gf_pyexp,                             /* @nb_power@ */
2011   mp_pyid,                              /* @nb_negative@ */
2012   mp_pyid,                              /* @nb_positive@ */
2013   mp_pyid,                              /* @nb_absolute@ */
2014   mp_pynonzerop,                        /* @nb_nonzero@ */
2015   0 /* doesn't make any sense */,       /* @nb_invert@ */
2016   gf_pylsl,                             /* @nb_lshift@ */
2017   gf_pylsr,                             /* @nb_rshift@ */
2018   gf_pyand,                             /* @nb_and@ */
2019   gf_pyxor,                             /* @nb_xor@ */
2020   gf_pyor,                              /* @nb_or@ */
2021   gf_pycoerce,                          /* @nb_coerce@ */
2022   mp_pyint,                             /* @nb_int@ */
2023   mp_pylong,                            /* @nb_long@ */
2024   0 /* doesn't make any sense */,       /* @nb_float@ */
2025   mp_pyoct,                             /* @nb_oct@ */
2026   mp_pyhex,                             /* @nb_hex@ */
2027
2028   0,                                    /* @nb_inplace_add@ */
2029   0,                                    /* @nb_inplace_subtract@ */
2030   0,                                    /* @nb_inplace_multiply@ */
2031   0,                                    /* @nb_inplace_divide@ */
2032   0,                                    /* @nb_inplace_remainder@ */
2033   0,                                    /* @nb_inplace_power@ */
2034   0,                                    /* @nb_inplace_lshift@ */
2035   0,                                    /* @nb_inplace_rshift@ */
2036   0,                                    /* @nb_inplace_and@ */
2037   0,                                    /* @nb_inplace_xor@ */
2038   0,                                    /* @nb_inplace_or@ */
2039
2040   gf_pydiv,                             /* @nb_floor_divide@ */
2041   0,                                    /* @nb_true_divide@ */
2042   0,                                    /* @nb_inplace_floor_divide@ */
2043   0,                                    /* @nb_inplace_true_divide@ */
2044 };
2045
2046 static PyTypeObject gf_pytype_skel = {
2047   PyObject_HEAD_INIT(0) 0,              /* Header */
2048   "GF",                                 /* @tp_name@ */
2049   sizeof(mp_pyobj),                     /* @tp_basicsize@ */
2050   0,                                    /* @tp_itemsize@ */
2051
2052   mp_pydealloc,                         /* @tp_dealloc@ */
2053   0,                                    /* @tp_print@ */
2054   0,                                    /* @tp_getattr@ */
2055   0,                                    /* @tp_setattr@ */
2056   0,                                    /* @tp_compare@ */
2057   gf_pyrepr,                            /* @tp_repr@ */
2058   &gf_pynumber,                         /* @tp_as_number@ */
2059   0,                                    /* @tp_as_sequence@ */
2060   0,                                    /* @tp_as_mapping@ */
2061   gf_pyhash,                            /* @tp_hash@ */
2062   0,                                    /* @tp_call@ */
2063   mp_pyhex,                             /* @tp_str@ */
2064   0,                                    /* @tp_getattro@ */
2065   0,                                    /* @tp_setattro@ */
2066   0,                                    /* @tp_as_buffer@ */
2067   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
2068     Py_TPFLAGS_CHECKTYPES |
2069     Py_TPFLAGS_BASETYPE,
2070
2071   /* @tp_doc@ */
2072 "Binary polynomials.  Support almost all the standard arithmetic\n\
2073 operations.\n\
2074 \n\
2075 Constructor gf(X, radix = R) attempts to convert X to a `gf'.  If\n\
2076 X is a string, it's read in radix-R form, or we look for a prefix\n\
2077 if R = 0.  Other acceptable things are ints and longs.\n\
2078 \n\
2079 The name is hopelessly wrong from a technical point of view, but\n\
2080 but it's much easier to type than `p2' or `c2' or whatever.\n\
2081 \n\
2082 Notes:\n\
2083 \n\
2084   * Use `//' for division.  GFs don't have `/' division.",
2085
2086   0,                                    /* @tp_traverse@ */
2087   0,                                    /* @tp_clear@ */
2088   gf_pyrichcompare,                     /* @tp_richcompare@ */
2089   0,                                    /* @tp_weaklistoffset@ */
2090   0,                                    /* @tp_iter@ */
2091   0,                                    /* @tp_iternext@ */
2092   gf_pymethods,                         /* @tp_methods@ */
2093   0,                                    /* @tp_members@ */
2094   gf_pygetset,                          /* @tp_getset@ */
2095   0,                                    /* @tp_base@ */
2096   0,                                    /* @tp_dict@ */
2097   0,                                    /* @tp_descr_get@ */
2098   0,                                    /* @tp_descr_set@ */
2099   0,                                    /* @tp_dictoffset@ */
2100   0,                                    /* @tp_init@ */
2101   PyType_GenericAlloc,                  /* @tp_alloc@ */
2102   gf_pynew,                             /* @tp_new@ */
2103   0,                                    /* @tp_free@ */
2104   0                                     /* @tp_is_gc@ */
2105 };
2106
2107 static PyObject *meth__GF_fromstring(PyObject *me,
2108                                     PyObject *arg, PyObject *kw)
2109 {
2110   int r = 0;
2111   char *p;
2112   int len;
2113   PyObject *z = 0;
2114   mp *zz;
2115   mptext_stringctx sc;
2116   char *kwlist[] = { "class", "x", "radix", 0 };
2117
2118   if (!PyArg_ParseTupleAndKeywords(arg, kw, "Os#|i:fromstring",
2119                                    kwlist, &me, &p, &len, &r))
2120     goto end;
2121   if (!good_radix_p(r, 1)) VALERR("bad radix");
2122   sc.buf = p; sc.lim = p + len;
2123   if ((zz = mp_read(MP_NEW, r, &mptext_stringops, &sc)) == 0 ||
2124       MP_NEGP(zz)) {
2125     if (zz) MP_DROP(zz);
2126     VALERR("bad binary polynomial");
2127   }
2128   z = Py_BuildValue("(Ns#)", gf_pywrap(zz), sc.buf, (int)(sc.lim - sc.buf));
2129 end:
2130   return (z);
2131 }
2132
2133 /*----- Sparse poly reduction ---------------------------------------------*/
2134
2135 typedef struct gfreduce_pyobj {
2136   PyObject_HEAD
2137   gfreduce mr;
2138 } gfreduce_pyobj;
2139
2140 #define GFREDUCE_PY(o) (&((gfreduce_pyobj *)(o))->mr)
2141
2142 static PyObject *grmeth_exp(PyObject *me, PyObject *arg)
2143 {
2144   PyObject *rc = 0;
2145   mp *yy = 0, *zz = 0;
2146
2147   if (!PyArg_ParseTuple(arg, "O&O&:exp", convgf, &yy, convgf, &zz))
2148     goto end;
2149   if (MP_NEGP(zz)) {
2150     if ((yy = gf_modinv_checked(yy, yy, GFREDUCE_PY(me)->p)) == 0) goto end;
2151     zz = mp_neg(zz, zz);
2152   }
2153   rc = gf_pywrap(gfreduce_exp(GFREDUCE_PY(me), MP_NEW, yy, zz));
2154 end:
2155   if (yy) MP_DROP(yy); if (zz) MP_DROP(zz);
2156   return (rc);
2157 }
2158
2159 static PyObject *grmeth_trace(PyObject *me, PyObject *arg)
2160 {
2161   PyObject *rc = 0;
2162   mp *xx = 0;
2163
2164   if (!PyArg_ParseTuple(arg, "O&:trace", convgf, &xx)) goto end;
2165   rc = PyInt_FromLong(gfreduce_trace(GFREDUCE_PY(me), xx));
2166 end:
2167   if (xx) MP_DROP(xx);
2168   return (rc);
2169 }
2170
2171 static PyObject *grmeth_halftrace(PyObject *me, PyObject *arg)
2172 {
2173   PyObject *rc = 0;
2174   mp *xx = 0;
2175
2176   if (!PyArg_ParseTuple(arg, "O&:halftrace", convgf, &xx)) goto end;
2177   rc = gf_pywrap(gfreduce_halftrace(GFREDUCE_PY(me), MP_NEW, xx));
2178 end:
2179   if (xx) MP_DROP(xx);
2180   return (rc);
2181 }
2182
2183 static PyObject *grmeth_sqrt(PyObject *me, PyObject *arg)
2184 {
2185   PyObject *rc = 0;
2186   mp *xx = 0, *yy;
2187
2188   if (!PyArg_ParseTuple(arg, "O&:sqrt", convgf, &xx)) goto end;
2189   if ((yy = gfreduce_sqrt(GFREDUCE_PY(me), MP_NEW, xx)) == 0)
2190     VALERR("no modular square root");
2191   rc = gf_pywrap(yy);
2192 end:
2193   if (xx) MP_DROP(xx);
2194   return (rc);
2195 }
2196
2197 static PyObject *grmeth_quadsolve(PyObject *me, PyObject *arg)
2198 {
2199   PyObject *rc = 0;
2200   mp *xx = 0, *yy;
2201
2202   if (!PyArg_ParseTuple(arg, "O&:quadsolve", convgf, &xx)) goto end;
2203   if ((yy = gfreduce_quadsolve(GFREDUCE_PY(me), MP_NEW, xx)) == 0)
2204     VALERR("no solution found");
2205   rc = gf_pywrap(yy);
2206 end:
2207   if (xx) MP_DROP(xx);
2208   return (rc);
2209 }
2210
2211 static PyObject *grmeth_reduce(PyObject *me, PyObject *arg)
2212 {
2213   PyObject *z = 0;
2214   mp *yy = 0;
2215
2216   if (!PyArg_ParseTuple(arg, "O&:reduce", convgf, &yy)) goto end;
2217   z = gf_pywrap(gfreduce_do(GFREDUCE_PY(me), MP_NEW, yy));
2218 end:
2219   return (z);
2220 }
2221
2222 static void gfreduce_pydealloc(PyObject *me)
2223 {
2224   gfreduce_destroy(GFREDUCE_PY(me));
2225   FREEOBJ(me);
2226 }
2227
2228 static PyObject *gfreduce_pynew(PyTypeObject *ty,
2229                                  PyObject *arg, PyObject *kw)
2230 {
2231   gfreduce_pyobj *mr = 0;
2232   gfreduce r;
2233   char *kwlist[] = { "m", 0 };
2234   mp *xx = 0;
2235
2236   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&:new", kwlist, convgf, &xx))
2237     goto end;
2238   if (MP_ZEROP(xx)) ZDIVERR("modulus is zero!");
2239   gfreduce_create(&r, xx);
2240   mr = (gfreduce_pyobj *)ty->tp_alloc(ty, 0);
2241   mr->mr = r;
2242 end:
2243   if (xx) MP_DROP(xx);
2244   return ((PyObject *)mr);
2245 }
2246
2247 static PyObject *grget_m(PyObject *me, void *hunoz)
2248   { return (gf_pywrap(MP_COPY(GFREDUCE_PY(me)->p))); }
2249
2250 static PyGetSetDef gfreduce_pygetset[] = {
2251 #define GETSETNAME(op, name) gr##op##_##name
2252   GET   (m,             "R.m -> reduction polynomial")
2253 #undef GETSETNAME
2254   { 0 }
2255 };
2256
2257 static PyMethodDef gfreduce_pymethods[] = {
2258 #define METHNAME(name) grmeth_##name
2259   METH  (reduce,        "R.reduce(X) -> X mod B.m")
2260   METH  (trace,        "R.trace(X) -> Tr(X) = x + x^2 + ... + x^{2^{m - 1}}")
2261   METH  (halftrace,   "R.halftrace(X) -> x + x^{2^2} + ... + x^{2^{m - 1}}")
2262   METH  (sqrt,          "R.sqrt(X) -> Y where Y^2 = X mod R")
2263   METH  (quadsolve,     "R.quadsolve(X) -> Y where Y^2 + Y = X mod R")
2264   METH  (exp,           "R.exp(X, N) -> X^N mod B.m")
2265 #undef METHNAME
2266   { 0 }
2267 };
2268
2269 static PyTypeObject *gfreduce_pytype, gfreduce_pytype_skel = {
2270   PyObject_HEAD_INIT(0) 0,              /* Header */
2271   "GFReduce",                           /* @tp_name@ */
2272   sizeof(gfreduce_pyobj),               /* @tp_basicsize@ */
2273   0,                                    /* @tp_itemsize@ */
2274
2275   gfreduce_pydealloc,                   /* @tp_dealloc@ */
2276   0,                                    /* @tp_print@ */
2277   0,                                    /* @tp_getattr@ */
2278   0,                                    /* @tp_setattr@ */
2279   0,                                    /* @tp_compare@ */
2280   0,                                    /* @tp_repr@ */
2281   0,                                    /* @tp_as_number@ */
2282   0,                                    /* @tp_as_sequence@ */
2283   0,                                    /* @tp_as_mapping@ */
2284   0,                                    /* @tp_hash@ */
2285   0,                                    /* @tp_call@ */
2286   0,                                    /* @tp_str@ */
2287   0,                                    /* @tp_getattro@ */
2288   0,                                    /* @tp_setattro@ */
2289   0,                                    /* @tp_as_buffer@ */
2290   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
2291     Py_TPFLAGS_BASETYPE,
2292
2293   /* @tp_doc@ */
2294 "A reduction context for reduction modulo sparse irreducible polynomials.",
2295
2296   0,                                    /* @tp_traverse@ */
2297   0,                                    /* @tp_clear@ */
2298   0,                                    /* @tp_richcompare@ */
2299   0,                                    /* @tp_weaklistoffset@ */
2300   0,                                    /* @tp_iter@ */
2301   0,                                    /* @tp_iternext@ */
2302   gfreduce_pymethods,                   /* @tp_methods@ */
2303   0,                                    /* @tp_members@ */
2304   gfreduce_pygetset,                    /* @tp_getset@ */
2305   0,                                    /* @tp_base@ */
2306   0,                                    /* @tp_dict@ */
2307   0,                                    /* @tp_descr_get@ */
2308   0,                                    /* @tp_descr_set@ */
2309   0,                                    /* @tp_dictoffset@ */
2310   0,                                    /* @tp_init@ */
2311   PyType_GenericAlloc,                  /* @tp_alloc@ */
2312   gfreduce_pynew,                       /* @tp_new@ */
2313   0,                                    /* @tp_free@ */
2314   0                                     /* @tp_is_gc@ */
2315 };
2316
2317 /*----- Normal/poly transformation ----------------------------------------*/
2318
2319 typedef struct gfn_pyobj {
2320   PyObject_HEAD
2321   mp *p;
2322   gfn ntop, pton;
2323 } gfn_pyobj;
2324
2325 static PyTypeObject *gfn_pytype, gfn_pytype_skel;
2326
2327 #define GFN_P(o) (((gfn_pyobj *)(o))->p)
2328 #define GFN_PTON(o) (&((gfn_pyobj *)(o))->pton)
2329 #define GFN_NTOP(o) (&((gfn_pyobj *)(o))->ntop)
2330
2331 static PyObject *gfn_pynew(PyTypeObject *ty, PyObject *arg, PyObject *kw)
2332 {
2333   mp *p = 0, *beta = 0;
2334   gfn_pyobj *gg = 0;
2335   char *kwlist[] = { "p", "beta", 0 };
2336
2337   if (!PyArg_ParseTupleAndKeywords(arg, kw, "O&O&:new", kwlist,
2338                                    convgf, &p, convgf, &beta))
2339     goto end;
2340   gg = PyObject_New(gfn_pyobj, ty);
2341   if (gfn_create(p, beta, &gg->ntop, &gg->pton)) {
2342     FREEOBJ(gg);
2343     gg = 0;
2344     VALERR("can't invert transformation matrix");
2345   }
2346   gg->p = MP_COPY(p);
2347 end:
2348   mp_drop(p);
2349   mp_drop(beta);
2350   return ((PyObject *)gg);
2351 }
2352
2353 static PyObject *gfnget_p(PyObject *me, void *hunoz)
2354   { return (gf_pywrap(MP_COPY(GFN_P(me)))); }
2355
2356 static PyObject *gfnget_beta(PyObject *me, void *hunoz)
2357 {
2358   gfn *n = GFN_NTOP(me);
2359   mp *x = n->r[n->n - 1];
2360   return (gf_pywrap(MP_COPY(x)));
2361 }
2362
2363 #define XFORMOP(name, NAME)                                             \
2364   static PyObject *gfnmeth_##name(PyObject *me, PyObject *arg)          \
2365   {                                                                     \
2366     mp *xx = 0;                                                         \
2367     mp *z = 0;                                                          \
2368                                                                         \
2369     if (!PyArg_ParseTuple(arg, "O&:" #name, convgf, &xx)) goto end;     \
2370     z = gfn_transform(GFN_##NAME(me), MP_NEW, xx);                      \
2371   end:                                                                  \
2372     mp_drop(xx);                                                        \
2373     if (!z) return (0);                                                 \
2374     return (mp_pywrap(z));                                              \
2375   }
2376 XFORMOP(pton, PTON)
2377 XFORMOP(ntop, NTOP)
2378 #undef XFORMOP
2379
2380 static void gfn_pydealloc(PyObject *me)
2381 {
2382   gfn_destroy(GFN_PTON(me));
2383   gfn_destroy(GFN_NTOP(me));
2384   FREEOBJ(me);
2385 }
2386
2387 static PyGetSetDef gfn_pygetset[] = {
2388 #define GETSETNAME(op, name) gfn##op##_##name
2389   GET   (p,             "X.p -> polynomial basis, as polynomial")
2390   GET   (beta,          "X.beta -> normal basis element, in poly form")
2391 #undef GETSETNAME
2392   { 0 }
2393 };
2394
2395 static PyMethodDef gfn_pymethods[] = {
2396 #define METHNAME(name) gfnmeth_##name
2397   METH  (pton,          "X.pton(A) -> normal-basis representation of A")
2398   METH  (ntop,          "X.ntop(A) -> polynomial-basis representation of A")
2399 #undef METHNAME
2400   { 0 }
2401 };
2402
2403 static PyTypeObject gfn_pytype_skel = {
2404   PyObject_HEAD_INIT(0) 0,              /* Header */
2405   "GFN",                                /* @tp_name@ */
2406   sizeof(gfn_pyobj),                    /* @tp_basicsize@ */
2407   0,                                    /* @tp_itemsize@ */
2408
2409   gfn_pydealloc,                        /* @tp_dealloc@ */
2410   0,                                    /* @tp_print@ */
2411   0,                                    /* @tp_getattr@ */
2412   0,                                    /* @tp_setattr@ */
2413   0,                                    /* @tp_compare@ */
2414   0,                                    /* @tp_repr@ */
2415   0,                                    /* @tp_as_number@ */
2416   0,                                    /* @tp_as_sequence@ */
2417   0,                                    /* @tp_as_mapping@ */
2418   0,                                    /* @tp_hash@ */
2419   0,                                    /* @tp_call@ */
2420   0,                                    /* @tp_str@ */
2421   0,                                    /* @tp_getattro@ */
2422   0,                                    /* @tp_setattro@ */
2423   0,                                    /* @tp_as_buffer@ */
2424   Py_TPFLAGS_DEFAULT |                  /* @tp_flags@ */
2425     Py_TPFLAGS_BASETYPE,
2426
2427   /* @tp_doc@ */
2428 "An object for transforming elements of binary fields between polynomial\n\
2429 and normal basis representations.",
2430
2431   0,                                    /* @tp_traverse@ */
2432   0,                                    /* @tp_clear@ */
2433   0,                                    /* @tp_richcompare@ */
2434   0,                                    /* @tp_weaklistoffset@ */
2435   0,                                    /* @tp_iter@ */
2436   0,                                    /* @tp_iternext@ */
2437   gfn_pymethods,                        /* @tp_methods@ */
2438   0,                                    /* @tp_members@ */
2439   gfn_pygetset,                         /* @tp_getset@ */
2440   0,                                    /* @tp_base@ */
2441   0,                                    /* @tp_dict@ */
2442   0,                                    /* @tp_descr_get@ */
2443   0,                                    /* @tp_descr_set@ */
2444   0,                                    /* @tp_dictoffset@ */
2445   0,                                    /* @tp_init@ */
2446   PyType_GenericAlloc,                  /* @tp_alloc@ */
2447   gfn_pynew,                            /* @tp_new@ */
2448   0,                                    /* @tp_free@ */
2449   0                                     /* @tp_is_gc@ */
2450 };
2451
2452 /*----- Glue --------------------------------------------------------------*/
2453
2454 static PyMethodDef methods[] = {
2455 #define METHNAME(func) meth_##func
2456   KWMETH(_MP_fromstring,        "\
2457 fromstring(STR, radix = 0) -> (X, REST)\n\
2458 \n\
2459 Parse STR as a large integer, according to radix.  If radix is zero,\n\
2460 read a prefix from STR to decide radix: allow `0' for octal, `0x' for hex\n\
2461 or `R_' for other radix R.")
2462   KWMETH(_GF_fromstring,        "\
2463 fromstring(STR, radix = 0) -> (X, REST)\n\
2464 \n\
2465 Parse STR as a binary polynomial, according to radix.  If radix is zero,\n\
2466 read a prefix from STR to decide radix: allow `0' for octal, `0x' for hex\n\
2467 or `R_' for other radix R.")
2468   METH  (_MP_factorial,         "\
2469 factorial(I) -> I!: compute factorial")
2470   METH  (_MP_fibonacci,         "\
2471 fibonacci(I) -> F(I): compute Fibonacci number")
2472   METH  (_MP_loadl,             "\
2473 loadl(STR) -> X: read little-endian bytes")
2474   METH  (_MP_loadb,             "\
2475 loadb(STR) -> X: read big-endian bytes")
2476   METH  (_MP_loadl2c,           "\
2477 loadl2c(STR) -> X: read little-endian bytes, two's complement")
2478   METH  (_MP_loadb2c,           "\
2479 loadb2c(STR) -> X: read big-endian bytes, two's complement")
2480   METH  (_MP_frombuf,           "\
2481 frombuf(STR) -> (X, REST): read buffer format")
2482   METH  (_GF_loadl,             "\
2483 loadl(STR) -> X: read little-endian bytes")
2484   METH  (_GF_loadb,             "\
2485 loadb(STR) -> X: read big-endian bytes")
2486   METH  (_GF_frombuf,           "\
2487 frombuf(STR) -> (X, REST): read buffer format")
2488 #undef METHNAME
2489   { 0 }
2490 };
2491
2492 void mp_pyinit(void)
2493 {
2494   INITTYPE(mp, root);
2495   INITTYPE(gf, root);
2496   INITTYPE(mpmul, root);
2497   INITTYPE(mpmont, root);
2498   INITTYPE(mpbarrett, root);
2499   INITTYPE(mpreduce, root);
2500   INITTYPE(mpcrt, root);
2501   INITTYPE(gfreduce, root);
2502   INITTYPE(gfn, root);
2503   addmethods(methods);
2504 }
2505
2506 void mp_pyinsert(PyObject *mod)
2507 {
2508   INSERT("MP", mp_pytype);
2509   INSERT("MPMul", mpmul_pytype);
2510   INSERT("MPMont", mpmont_pytype);
2511   INSERT("MPBarrett", mpbarrett_pytype);
2512   INSERT("MPReduce", mpreduce_pytype);
2513   INSERT("MPCRT", mpcrt_pytype);
2514   INSERT("GF", gf_pytype);
2515   INSERT("GFReduce", gfreduce_pytype);
2516   INSERT("GFN", gfn_pytype);
2517 }
2518
2519 /*----- That's all, folks -------------------------------------------------*/