X-Git-Url: http://www.chiark.greenend.org.uk/ucgi/~mdw/git/catacomb-python/blobdiff_plain/bf8910bde6a430a64868ef70319f94e94114c82c..183e9cd31b1ac2f14b86c5de6ac2643b8a4364a2:/mp.c diff --git a/mp.c b/mp.c index 00b728d..b253dba 100644 --- a/mp.c +++ b/mp.c @@ -1,13 +1,11 @@ /* -*-c-*- - * - * $Id$ * * Multiprecision arithmetic * * (c) 2004 Straylight/Edgeware */ -/*----- Licensing notice --------------------------------------------------* +/*----- Licensing notice --------------------------------------------------* * * This file is part of the Python interface to Catacomb. * @@ -15,12 +13,12 @@ * it under the terms of the GNU General Public License as published by * the Free Software Foundation; either version 2 of the License, or * (at your option) any later version. - * + * * Catacomb/Python is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * GNU General Public License for more details. - * + * * You should have received a copy of the GNU General Public License * along with Catacomb/Python; if not, write to the Free Software Foundation, * Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. @@ -165,7 +163,7 @@ PyObject *gf_pywrap(mp *x) return ((PyObject *)z); } -int mp_tolong_checked(mp *x, long *l) +int mp_tolong_checked(mp *x, long *l, int must) { static mp *longmin = 0, *longmax = 0; int rc = -1; @@ -174,8 +172,10 @@ int mp_tolong_checked(mp *x, long *l) longmin = mp_fromlong(MP_NEW, LONG_MIN); longmax = mp_fromlong(MP_NEW, LONG_MAX); } - if (MP_CMP(x, <, longmin) || MP_CMP(x, >, longmax)) - VALERR("mp out of range for int"); + if (MP_CMP(x, <, longmin) || MP_CMP(x, >, longmax)) { + if (must) VALERR("mp out of range for int"); + else goto end; + } *l = mp_tolong(x); rc = 0; end: @@ -364,7 +364,7 @@ static PyObject *mp_pyid(PyObject *x) { RETURN_OBJ(x); } PyObject *z = 0; \ long n; \ if (pre##binop(x, y, &xx, &yy)) RETURN_NOTIMPL; \ - if (mp_tolong_checked(yy, &n)) goto end; \ + if (mp_tolong_checked(yy, &n, 1)) goto end; \ if (n < 0) \ z = pre##_pywrap(mp_##rname(MP_NEW, xx, -n)); \ else \ @@ -485,8 +485,8 @@ static int mp_pynonzerop(PyObject *x) { return !MP_ZEROP(MP_X(x)); } static PyObject *mp_pyint(PyObject *x) { long l; - if (mp_tolong_checked(MP_X(x), &l)) return (0); - return (PyInt_FromLong(l)); + if (!mp_tolong_checked(MP_X(x), &l, 0)) return (PyInt_FromLong(l)); + else return mp_topylong(MP_X(x)); } static PyObject *mp_pylong(PyObject *x) { return (mp_topylong(MP_X(x))); } @@ -498,7 +498,7 @@ static PyObject *mp_pyfloat(PyObject *x) return (PyFloat_FromDouble(f)); } -#define COERCE(pre, PRE) \ +#define COERCE(pre, PRE) \ static int pre##_pycoerce(PyObject **x, PyObject **y) \ { \ mp *z; \ @@ -557,8 +557,6 @@ static PyObject *mpmeth_jacobi(PyObject *me, PyObject *arg) PyObject *z = 0; if (!PyArg_ParseTuple(arg, "O&:jacobi", convmp, &y)) goto end; - if (MP_NEGP(MP_X(me)) || MP_EVENP(MP_X(me))) - VALERR("must be positive and odd"); z = PyInt_FromLong(mp_jacobi(y, MP_X(me))); end: if (y) MP_DROP(y); @@ -574,7 +572,7 @@ end: } BITOP(mp, setbit, 2c); BITOP(mp, clearbit, 2c); -BITOP(gf, setbit, ); +BITOP(gf, setbit, ); BITOP(gf, clearbit, ); #undef BITOP @@ -684,6 +682,19 @@ end: return (z); } +static PyObject *mpmeth_leastcongruent(PyObject *me, PyObject *arg) +{ + mp *z, *b, *m; + PyObject *rc = 0; + + if (!PyArg_ParseTuple(arg, "O&O&:leastcongruent", convmp, &b, convmp, &m)) + goto end; + z = mp_leastcongruent(MP_NEW, b, MP_X(me), m); + rc = mp_pywrap(z); +end: + return (rc); +} + #define STOREOP(name, c) \ static PyObject *mpmeth_##name(PyObject *me, \ PyObject *arg, PyObject *kw) \ @@ -715,7 +726,7 @@ STOREOP(storeb2c, 2c) { \ buf b; \ char *p; \ - int sz; \ + Py_ssize_t sz; \ PyObject *rc = 0; \ mp *x; \ \ @@ -794,6 +805,8 @@ static PyMethodDef mp_pymethods[] = { "X.gcdx(Y) -> (gcd(X, Y), U, V) with X U + Y V = gcd(X, Y)") METH (modinv, "X.modinv(Y) -> multiplicative inverse of Y mod X") METH (modsqrt, "X.modsqrt(Y) -> square root of Y mod X, if X prime") + METH (leastcongruent, + "X.leastcongruent(B, M) -> smallest Z >= B with Z == X (mod M)") KWMETH(primep, "X.primep(rng = rand) -> true/false if X is prime") KWMETH(tostring, "X.tostring(radix = 10) -> STR") KWMETH(storel, "X.storel(len = -1) -> little-endian bytes") @@ -852,7 +865,7 @@ static PyNumberMethods mp_pynumber = { static PyTypeObject mp_pytype_skel = { PyObject_HEAD_INIT(0) 0, /* Header */ - "catacomb.MP", /* @tp_name@ */ + "MP", /* @tp_name@ */ sizeof(mp_pyobj), /* @tp_basicsize@ */ 0, /* @tp_itemsize@ */ @@ -879,13 +892,13 @@ static PyTypeObject mp_pytype_skel = { "Multiprecision integers, similar to `long' but more efficient and\n\ versatile. Support all the standard arithmetic operations.\n\ \n\ -Constructor mp(X, radix = R) attempts to convert X to an `mp'. If\n\ +Constructor mp(X, [radix = R]) attempts to convert X to an `mp'. If\n\ X is a string, it's read in radix-R form, or we look for a prefix\n\ if R = 0. Other acceptable things are ints and longs.\n\ \n\ Notes:\n\ \n\ - * Use `//' for division. MPs don't have `/' division.", + * Use `//' for integer division. `/' gives exact rational division.", 0, /* @tp_traverse@ */ 0, /* @tp_clear@ */ @@ -913,7 +926,7 @@ static PyObject *meth__MP_fromstring(PyObject *me, { int r = 0; char *p; - int len; + Py_ssize_t len; PyObject *z = 0; mp *zz; mptext_stringctx sc; @@ -925,17 +938,37 @@ static PyObject *meth__MP_fromstring(PyObject *me, if (!good_radix_p(r, 1)) VALERR("bad radix"); sc.buf = p; sc.lim = p + len; if ((zz = mp_read(MP_NEW, r, &mptext_stringops, &sc)) == 0) - SYNERR("bad integer"); - z = Py_BuildValue("(Ns#)", mp_pywrap(zz), sc.buf, (int)(sc.lim - sc.buf)); + VALERR("bad integer"); + z = Py_BuildValue("(Ns#)", mp_pywrap(zz), + sc.buf, (Py_ssize_t)(sc.lim - sc.buf)); end: return (z); } +static PyObject *meth__MP_factorial(PyObject *me, PyObject *arg) +{ + unsigned long i; + mp *x; + if (!PyArg_ParseTuple(arg, "OO&:factorial", &me, convulong, &i)) + return (0); + x = mp_factorial(i); + return mp_pywrap(x); +} + +static PyObject *meth__MP_fibonacci(PyObject *me, PyObject *arg) +{ + long i; + mp *x; + if (!PyArg_ParseTuple(arg, "Ol:fibonacci", &me, &i)) return (0); + x = mp_fibonacci(i); + return mp_pywrap(x); +} + #define LOADOP(pre, py, name) \ static PyObject *meth__##py##_##name(PyObject *me, PyObject *arg) \ { \ char *p; \ - int len; \ + Py_ssize_t len; \ if (!PyArg_ParseTuple(arg, "Os#:" #name, &me, &p, &len)) return (0); \ return (pre##_pywrap(mp_##name(MP_NEW, p, len))); \ } @@ -1045,7 +1078,7 @@ static PyMethodDef mpmul_pymethods[] = { static PyTypeObject *mpmul_pytype, mpmul_pytype_skel = { PyObject_HEAD_INIT(0) 0, /* Header */ - "catacomb.MPMul", /* @tp_name@ */ + "MPMul", /* @tp_name@ */ sizeof(mpmul_pyobj), /* @tp_basicsize@ */ 0, /* @tp_itemsize@ */ @@ -1068,7 +1101,7 @@ static PyTypeObject *mpmul_pytype, mpmul_pytype_skel = { Py_TPFLAGS_BASETYPE, /* @tp_doc@ */ -"An object for multiplying many small integers.", +"MPMul(N_0, N_1, ....): an object for multiplying many small integers.", 0, /* @tp_traverse@ */ 0, /* @tp_clear@ */ @@ -1193,7 +1226,7 @@ fail: static PyObject *mm_mexpr(PyObject *me, void *v, int n) { return mp_pywrap(mpmont_mexpr(MPMONT_PY(me), MP_NEW, v, n)); } - + static void mp_mexp_drop(void *p) { mp_expfactor *f = p; @@ -1233,7 +1266,7 @@ fail: static PyObject *mm_mexp(PyObject *me, void *v, int n) { return mp_pywrap(mpmont_mexp(MPMONT_PY(me), MP_NEW, v, n)); } - + static PyObject *mmmeth_mexp(PyObject *me, PyObject *arg) { return mexp_common(me, arg, sizeof(mp_expfactor), @@ -1294,17 +1327,17 @@ static PyGetSetDef mpmont_pygetset[] = { static PyMethodDef mpmont_pymethods[] = { #define METHNAME(name) mmmeth_##name - METH (int, "M.out(X) -> XR") + METH (int, "M.int(X) -> XR") METH (mul, "M.mul(XR, YR) -> ZR where Z = X Y") METH (expr, "M.expr(XR, N) -> ZR where Z = X^N mod M.m") METH (mexpr, "\ -B.mexp([(XR0, N0), (XR1, N1), ...]) = ZR where Z = X0^N0 X1^N1 mod B.m\n\ +M.mexpr([(XR0, N0), (XR1, N1), ...]) = ZR where Z = X0^N0 X1^N1 ... mod M.m\n\ \t(the list may be flattened if this more convenient.)") METH (reduce, "M.reduce(XR) -> X") METH (ext, "M.ext(XR) -> X") METH (exp, "M.exp(X, N) -> X^N mod M.m") METH (mexp, "\ -B.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 mod B.m\n\ +M.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 ... mod M.m\n\ \t(the list may be flattened if this more convenient.)") #undef METHNAME { 0 } @@ -1312,7 +1345,7 @@ B.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 mod B.m\n\ static PyTypeObject *mpmont_pytype, mpmont_pytype_skel = { PyObject_HEAD_INIT(0) 0, /* Header */ - "catacomb.MPMont", /* @tp_name@ */ + "MPMont", /* @tp_name@ */ sizeof(mpmont_pyobj), /* @tp_basicsize@ */ 0, /* @tp_itemsize@ */ @@ -1335,7 +1368,7 @@ static PyTypeObject *mpmont_pytype, mpmont_pytype_skel = { Py_TPFLAGS_BASETYPE, /* @tp_doc@ */ -"A Montgomery reduction context.", +"MPMont(N): a Montgomery reduction context.", 0, /* @tp_traverse@ */ 0, /* @tp_clear@ */ @@ -1386,7 +1419,7 @@ end: static PyObject *mb_mexp(PyObject *me, void *v, int n) { return mp_pywrap(mpbarrett_mexp(MPBARRETT_PY(me), MP_NEW, v, n)); } - + static PyObject *mbmeth_mexp(PyObject *me, PyObject *arg) { return mexp_common(me, arg, sizeof(mp_expfactor), @@ -1443,7 +1476,7 @@ static PyMethodDef mpbarrett_pymethods[] = { METH (reduce, "B.reduce(X) -> X mod B.m") METH (exp, "B.exp(X, N) -> X^N mod B.m") METH (mexp, "\ -B.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 mod B.m\n\ +B.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 ... mod B.m\n\ \t(the list may be flattened if this more convenient.)") #undef METHNAME { 0 } @@ -1451,7 +1484,7 @@ B.mexp([(X0, N0), (X1, N1), ...]) = X0^N0 X1^N1 mod B.m\n\ static PyTypeObject *mpbarrett_pytype, mpbarrett_pytype_skel = { PyObject_HEAD_INIT(0) 0, /* Header */ - "catacomb.MPBarrett", /* @tp_name@ */ + "MPBarrett", /* @tp_name@ */ sizeof(mpbarrett_pyobj), /* @tp_basicsize@ */ 0, /* @tp_itemsize@ */ @@ -1474,7 +1507,7 @@ static PyTypeObject *mpbarrett_pytype, mpbarrett_pytype_skel = { Py_TPFLAGS_BASETYPE, /* @tp_doc@ */ -"A Barrett reduction context.", +"MPBarrett(N): a Barrett reduction context.", 0, /* @tp_traverse@ */ 0, /* @tp_clear@ */ @@ -1579,7 +1612,7 @@ static PyMethodDef mpreduce_pymethods[] = { static PyTypeObject *mpreduce_pytype, mpreduce_pytype_skel = { PyObject_HEAD_INIT(0) 0, /* Header */ - "catacomb.MPReduce", /* @tp_name@ */ + "MPReduce", /* @tp_name@ */ sizeof(mpreduce_pyobj), /* @tp_basicsize@ */ 0, /* @tp_itemsize@ */ @@ -1602,7 +1635,7 @@ static PyTypeObject *mpreduce_pytype, mpreduce_pytype_skel = { Py_TPFLAGS_BASETYPE, /* @tp_doc@ */ -"A reduction context for reduction modulo primes of special form.", +"MPReduce(N): a reduction context for reduction modulo Solinas primes.", 0, /* @tp_traverse@ */ 0, /* @tp_clear@ */ @@ -1748,7 +1781,7 @@ static PyMethodDef mpcrt_pymethods[] = { static PyTypeObject *mpcrt_pytype, mpcrt_pytype_skel = { PyObject_HEAD_INIT(0) 0, /* Header */ - "catacomb.MPCRT", /* @tp_name@ */ + "MPCRT", /* @tp_name@ */ sizeof(mpcrt_pyobj), /* @tp_basicsize@ */ 0, /* @tp_itemsize@ */ @@ -1771,7 +1804,7 @@ static PyTypeObject *mpcrt_pytype, mpcrt_pytype_skel = { Py_TPFLAGS_BASETYPE, /* @tp_doc@ */ -"A context for the solution of Chinese Remainder Theorem problems.", +"MPCRT(SEQ): a context for solving Chinese Remainder Theorem problems.", 0, /* @tp_traverse@ */ 0, /* @tp_clear@ */ @@ -1970,7 +2003,7 @@ static PyMethodDef gf_pymethods[] = { METH (gcdx, "X.gcdx(Y) -> (gcd(X, Y), U, V) with X U + Y V = gcd(X, Y)") METH (modinv, "X.modinv(Y) -> multiplicative inverse of Y mod X") - METH (irreduciblep, "X.irreduciblep() -> true/false") + METH (irreduciblep, "X.irreduciblep() -> true/false") #undef METHNAME #define METHNAME(func) mpmeth_##func KWMETH(tostring, "X.tostring(radix = 10) -> STR") @@ -2030,7 +2063,7 @@ static PyNumberMethods gf_pynumber = { static PyTypeObject gf_pytype_skel = { PyObject_HEAD_INIT(0) 0, /* Header */ - "catacomb.GF", /* @tp_name@ */ + "GF", /* @tp_name@ */ sizeof(mp_pyobj), /* @tp_basicsize@ */ 0, /* @tp_itemsize@ */ @@ -2066,7 +2099,7 @@ but it's much easier to type than `p2' or `c2' or whatever.\n\ \n\ Notes:\n\ \n\ - * Use `//' for division. GFs don't have `/' division.", + * Use `//' for Euclidean division. `/' gives exact rational division.", 0, /* @tp_traverse@ */ 0, /* @tp_clear@ */ @@ -2094,7 +2127,7 @@ static PyObject *meth__GF_fromstring(PyObject *me, { int r = 0; char *p; - int len; + Py_ssize_t len; PyObject *z = 0; mp *zz; mptext_stringctx sc; @@ -2108,9 +2141,10 @@ static PyObject *meth__GF_fromstring(PyObject *me, if ((zz = mp_read(MP_NEW, r, &mptext_stringops, &sc)) == 0 || MP_NEGP(zz)) { if (zz) MP_DROP(zz); - SYNERR("bad binary polynomial"); + VALERR("bad binary polynomial"); } - z = Py_BuildValue("(Ns#)", gf_pywrap(zz), sc.buf, (int)(sc.lim - sc.buf)); + z = Py_BuildValue("(Ns#)", gf_pywrap(zz), + sc.buf, (Py_ssize_t)(sc.lim - sc.buf)); end: return (z); } @@ -2141,6 +2175,58 @@ end: return (rc); } +static PyObject *grmeth_trace(PyObject *me, PyObject *arg) +{ + PyObject *rc = 0; + mp *xx = 0; + + if (!PyArg_ParseTuple(arg, "O&:trace", convgf, &xx)) goto end; + rc = PyInt_FromLong(gfreduce_trace(GFREDUCE_PY(me), xx)); +end: + if (xx) MP_DROP(xx); + return (rc); +} + +static PyObject *grmeth_halftrace(PyObject *me, PyObject *arg) +{ + PyObject *rc = 0; + mp *xx = 0; + + if (!PyArg_ParseTuple(arg, "O&:halftrace", convgf, &xx)) goto end; + rc = gf_pywrap(gfreduce_halftrace(GFREDUCE_PY(me), MP_NEW, xx)); +end: + if (xx) MP_DROP(xx); + return (rc); +} + +static PyObject *grmeth_sqrt(PyObject *me, PyObject *arg) +{ + PyObject *rc = 0; + mp *xx = 0, *yy; + + if (!PyArg_ParseTuple(arg, "O&:sqrt", convgf, &xx)) goto end; + if ((yy = gfreduce_sqrt(GFREDUCE_PY(me), MP_NEW, xx)) == 0) + VALERR("no modular square root"); + rc = gf_pywrap(yy); +end: + if (xx) MP_DROP(xx); + return (rc); +} + +static PyObject *grmeth_quadsolve(PyObject *me, PyObject *arg) +{ + PyObject *rc = 0; + mp *xx = 0, *yy; + + if (!PyArg_ParseTuple(arg, "O&:quadsolve", convgf, &xx)) goto end; + if ((yy = gfreduce_quadsolve(GFREDUCE_PY(me), MP_NEW, xx)) == 0) + VALERR("no solution found"); + rc = gf_pywrap(yy); +end: + if (xx) MP_DROP(xx); + return (rc); +} + static PyObject *grmeth_reduce(PyObject *me, PyObject *arg) { PyObject *z = 0; @@ -2190,6 +2276,10 @@ static PyGetSetDef gfreduce_pygetset[] = { static PyMethodDef gfreduce_pymethods[] = { #define METHNAME(name) grmeth_##name METH (reduce, "R.reduce(X) -> X mod B.m") + METH (trace, "R.trace(X) -> Tr(X) = x + x^2 + ... + x^{2^{m - 1}}") + METH (halftrace, "R.halftrace(X) -> x + x^{2^2} + ... + x^{2^{m - 1}}") + METH (sqrt, "R.sqrt(X) -> Y where Y^2 = X mod R") + METH (quadsolve, "R.quadsolve(X) -> Y where Y^2 + Y = X mod R") METH (exp, "R.exp(X, N) -> X^N mod B.m") #undef METHNAME { 0 } @@ -2197,7 +2287,7 @@ static PyMethodDef gfreduce_pymethods[] = { static PyTypeObject *gfreduce_pytype, gfreduce_pytype_skel = { PyObject_HEAD_INIT(0) 0, /* Header */ - "catacomb.GFReduce", /* @tp_name@ */ + "GFReduce", /* @tp_name@ */ sizeof(gfreduce_pyobj), /* @tp_basicsize@ */ 0, /* @tp_itemsize@ */ @@ -2220,7 +2310,7 @@ static PyTypeObject *gfreduce_pytype, gfreduce_pytype_skel = { Py_TPFLAGS_BASETYPE, /* @tp_doc@ */ -"A reduction context for reduction modulo sparse irreducible polynomials.", +"GFReduce(N): a context for reduction modulo sparse polynomials.", 0, /* @tp_traverse@ */ 0, /* @tp_clear@ */ @@ -2331,7 +2421,7 @@ static PyMethodDef gfn_pymethods[] = { static PyTypeObject gfn_pytype_skel = { PyObject_HEAD_INIT(0) 0, /* Header */ - "catacomb.GFN", /* @tp_name@ */ + "GFN", /* @tp_name@ */ sizeof(gfn_pyobj), /* @tp_basicsize@ */ 0, /* @tp_itemsize@ */ @@ -2354,8 +2444,8 @@ static PyTypeObject gfn_pytype_skel = { Py_TPFLAGS_BASETYPE, /* @tp_doc@ */ -"An object for transforming elements of binary fields between polynomial\n\ -and normal basis representations.", +"GFN(P, BETA): an object for transforming elements of binary fields\n\ + between polynomial and normal basis representations.", 0, /* @tp_traverse@ */ 0, /* @tp_clear@ */ @@ -2382,18 +2472,22 @@ and normal basis representations.", static PyMethodDef methods[] = { #define METHNAME(func) meth_##func - KWMETH(_MP_fromstring, "\ + KWMETH(_MP_fromstring, "\ fromstring(STR, radix = 0) -> (X, REST)\n\ \n\ Parse STR as a large integer, according to radix. If radix is zero,\n\ read a prefix from STR to decide radix: allow `0' for octal, `0x' for hex\n\ or `R_' for other radix R.") - KWMETH(_GF_fromstring, "\ + KWMETH(_GF_fromstring, "\ fromstring(STR, radix = 0) -> (X, REST)\n\ \n\ Parse STR as a binary polynomial, according to radix. If radix is zero,\n\ read a prefix from STR to decide radix: allow `0' for octal, `0x' for hex\n\ or `R_' for other radix R.") + METH (_MP_factorial, "\ +factorial(I) -> I!: compute factorial") + METH (_MP_fibonacci, "\ +fibonacci(I) -> F(I): compute Fibonacci number") METH (_MP_loadl, "\ loadl(STR) -> X: read little-endian bytes") METH (_MP_loadb, "\