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