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

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