X-Git-Url: https://git.distorted.org.uk/u/mdw/catacomb/blobdiff_plain/79ba130cb5776f994f6a3f0f87159d8cbc5ff129..ef5f48103e83977bda6ef4d7d1aacbb66a629b10:/mpmont.c diff --git a/mpmont.c b/mpmont.c index 7522dba..d061934 100644 --- a/mpmont.c +++ b/mpmont.c @@ -1,6 +1,6 @@ /* -*-c-*- * - * $Id: mpmont.c,v 1.5 1999/11/22 13:58:40 mdw Exp $ + * $Id: mpmont.c,v 1.6 1999/12/10 23:18:39 mdw Exp $ * * Montgomery reduction * @@ -30,6 +30,9 @@ /*----- Revision history --------------------------------------------------* * * $Log: mpmont.c,v $ + * Revision 1.6 1999/12/10 23:18:39 mdw + * Change interface for suggested destinations. + * * Revision 1.5 1999/11/22 13:58:40 mdw * Add an option to disable Montgomery reduction, so that performance * comparisons can be done. @@ -77,6 +80,7 @@ * Returns: --- * * Use: Initializes a Montgomery reduction context ready for use. + * The argument @m@ must be a positive odd integer. */ #ifdef MPMONT_DISABLE @@ -93,6 +97,12 @@ void mpmont_create(mpmont *mm, mp *m) void mpmont_create(mpmont *mm, mp *m) { + /* --- Validate the arguments --- */ + + assert(((void)"Montgomery modulus must be positive", + (m->f & MP_NEG) == 0)); + assert(((void)"Montgomery modulus must be odd", m->v[0] & 1)); + /* --- Take a copy of the modulus --- */ mp_shrink(m); @@ -106,7 +116,7 @@ void mpmont_create(mpmont *mm, mp *m) { mpw av[2] = { 0, 1 }; mp a, b; - mp *i; + mp *i = MP_NEW; mpw mi; mp_build(&a, av, av + 2); @@ -158,14 +168,14 @@ void mpmont_destroy(mpmont *mm) * * Arguments: @mpmont *mm@ = pointer to Montgomery reduction context * @mp *d@ = destination - * @const mp *a@ = source, assumed positive + * @mp *a@ = source, assumed positive * * Returns: Result, %$a R^{-1} \bmod m$%. */ #ifdef MPMONT_DISABLE -mp *mpmont_reduce(mpmont *mm, mp *d, const mp *a) +mp *mpmont_reduce(mpmont *mm, mp *d, mp *a) { mp_div(0, &d, a, mm->m); return (d); @@ -173,10 +183,10 @@ mp *mpmont_reduce(mpmont *mm, mp *d, const mp *a) #else -mp *mpmont_reduce(mpmont *mm, mp *d, const mp *a) +mp *mpmont_reduce(mpmont *mm, mp *d, mp *a) { mpw *dv, *dvl; - const mpw *mv, *mvl; + mpw *mv, *mvl; size_t n; /* --- Initial conditioning of the arguments --- */ @@ -219,14 +229,14 @@ mp *mpmont_reduce(mpmont *mm, mp *d, const mp *a) * * Arguments: @mpmont *mm@ = pointer to Montgomery reduction context * @mp *d@ = destination - * @const mp *a, *b@ = sources, assumed positive + * @mp *a, *b@ = sources, assumed positive * * Returns: Result, %$a b R^{-1} \bmod m$%. */ #ifdef MPMONT_DISABLE -mp *mpmont_mul(mpmont *mm, mp *d, const mp *a, const mp *b) +mp *mpmont_mul(mpmont *mm, mp *d, mp *a, mp *b) { d = mp_mul(d, a, b); mp_div(0, &d, d, mm->m); @@ -235,59 +245,69 @@ mp *mpmont_mul(mpmont *mm, mp *d, const mp *a, const mp *b) #else -mp *mpmont_mul(mpmont *mm, mp *d, const mp *a, const mp *b) +mp *mpmont_mul(mpmont *mm, mp *d, mp *a, mp *b) { - mpw *dv, *dvl; - const mpw *av, *avl; - const mpw *bv, *bvl; - const mpw *mv, *mvl; - mpw y; - size_t n, i; - - /* --- Initial conditioning of the arguments --- */ + if (MP_LEN(a) > KARATSUBA_CUTOFF && MP_LEN(b) > KARATSUBA_CUTOFF) { + d = mp_mul(d, a, b); + d = mpmont_reduce(mm, d, d); + } else { + mpw *dv, *dvl; + mpw *av, *avl; + mpw *bv, *bvl; + mpw *mv, *mvl; + mpw y; + size_t n, i; + + /* --- Initial conditioning of the arguments --- */ + + if (MP_LEN(a) > MP_LEN(b)) { + mp *t = a; a = b; b = t; + } + n = MP_LEN(mm->m); - if (MP_LEN(a) > MP_LEN(b)) { - const mp *t = a; a = b; b = t; - } - n = MP_LEN(mm->m); - - MP_MODIFY(d, 2 * n + 1); - dv = d->v; dvl = d->vl; - MPX_ZERO(dv, dvl); - av = a->v; avl = a->vl; - bv = b->v; bvl = b->vl; - mv = mm->m->v; mvl = mm->m->vl; - y = *bv; + a = MP_COPY(a); + b = MP_COPY(b); + MP_MODIFY(d, 2 * n + 1); + dv = d->v; dvl = d->vl; + MPX_ZERO(dv, dvl); + av = a->v; avl = a->vl; + bv = b->v; bvl = b->vl; + mv = mm->m->v; mvl = mm->m->vl; + y = *bv; + + /* --- Montgomery multiplication phase --- */ + + i = 0; + while (i < n && av < avl) { + mpw x = *av++; + mpw u = MPW((*dv + x * y) * mm->mi); + MPX_UMLAN(dv, dvl, bv, bvl, x); + MPX_UMLAN(dv, dvl, mv, mvl, u); + dv++; + i++; + } - /* --- Montgomery multiplication phase --- */ + /* --- Simpler Montgomery reduction phase --- */ - i = 0; - while (i < n && av < avl) { - mpw x = *av++; - mpw u = MPW((*dv + x * y) * mm->mi); - MPX_UMLAN(dv, dvl, bv, bvl, x); - MPX_UMLAN(dv, dvl, mv, mvl, u); - dv++; - i++; - } + while (i < n) { + mpw u = MPW(*dv * mm->mi); + MPX_UMLAN(dv, dvl, mv, mvl, u); + dv++; + i++; + } - /* --- Simpler Montgomery reduction phase --- */ + /* --- Done --- */ - while (i < n) { - mpw u = MPW(*dv * mm->mi); - MPX_UMLAN(dv, dvl, mv, mvl, u); - dv++; - i++; + memmove(d->v, dv, MPWS(dvl - dv)); + d->vl -= dv - d->v; + MP_SHRINK(d); + d->f = (a->f | b->f) & MP_BURN; + if (MP_CMP(d, >=, mm->m)) + d = mp_sub(d, d, mm->m); + MP_DROP(a); + MP_DROP(b); } - /* --- Done --- */ - - memmove(d->v, dv, MPWS(dvl - dv)); - d->vl -= dv - d->v; - MP_SHRINK(d); - d->f = (a->f | b->f) & MP_BURN; - if (MP_CMP(d, >=, mm->m)) - d = mp_sub(d, d, mm->m); return (d); } @@ -296,17 +316,18 @@ mp *mpmont_mul(mpmont *mm, mp *d, const mp *a, const mp *b) /* --- @mpmont_expr@ --- * * * Arguments: @mpmont *mm@ = pointer to Montgomery reduction context - * @const mp *a@ = base - * @const mp *e@ = exponent + * @mp *d@ = fake destination + * @mp *a@ = base + * @mp *e@ = exponent * * Returns: Result, %$a^e R \bmod m$%. */ -mp *mpmont_expr(mpmont *mm, const mp *a, const mp *e) +mp *mpmont_expr(mpmont *mm, mp *d, mp *a, mp *e) { mpscan sc; mp *ar = mpmont_mul(mm, MP_NEW, a, mm->r2); - mp *d = MP_COPY(mm->r); + mp *x = MP_COPY(mm->r); mp *spare = MP_NEW; mp_scan(&sc, e); @@ -322,8 +343,8 @@ mp *mpmont_expr(mpmont *mm, const mp *a, const mp *e) spare = ar; ar = dd; sq--; } - dd = mpmont_mul(mm, spare, d, ar); - spare = d; d = dd; + dd = mpmont_mul(mm, spare, x, ar); + spare = x; x = dd; } sq++; if (!MP_STEP(&sc)) @@ -333,21 +354,24 @@ mp *mpmont_expr(mpmont *mm, const mp *a, const mp *e) MP_DROP(ar); if (spare != MP_NEW) MP_DROP(spare); - return (d); + if (d != MP_NEW) + MP_DROP(d); + return (x); } /* --- @mpmont_exp@ --- * * * Arguments: @mpmont *mm@ = pointer to Montgomery reduction context - * @const mp *a@ = base - * @const mp *e@ = exponent + * @mp *d@ = fake destination + * @mp *a@ = base + * @mp *e@ = exponent * * Returns: Result, %$a^e \bmod m$%. */ -mp *mpmont_exp(mpmont *mm, const mp *a, const mp *e) +mp *mpmont_exp(mpmont *mm, mp *d, mp *a, mp *e) { - mp *d = mpmont_expr(mm, a, e); + d = mpmont_expr(mm, d, a, e); d = mpmont_reduce(mm, d, d); return (d); } @@ -399,6 +423,7 @@ static int tcreate(dstr *v) MP_DROP(r); MP_DROP(r2); mpmont_destroy(&mm); + assert(mparena_count(MPARENA_GLOBAL) == 0); return (ok); } @@ -456,6 +481,7 @@ static int tmul(dstr *v) MP_DROP(b); MP_DROP(r); mpmont_destroy(&mm); + assert(mparena_count(MPARENA_GLOBAL) == 0); return ok; } @@ -471,7 +497,7 @@ static int texp(dstr *v) mpmont mm; mpmont_create(&mm, m); - mr = mpmont_exp(&mm, a, b); + mr = mpmont_exp(&mm, MP_NEW, a, b); if (MP_CMP(mr, !=, r)) { fputs("\n*** montgomery modexp failed", stderr); @@ -490,14 +516,15 @@ static int texp(dstr *v) MP_DROP(r); MP_DROP(mr); mpmont_destroy(&mm); + assert(mparena_count(MPARENA_GLOBAL) == 0); return ok; } static test_chunk tests[] = { - { "create", tcreate, { &type_mp, &type_mp, &type_mp, &type_mp } }, - { "mul", tmul, { &type_mp, &type_mp, &type_mp, &type_mp } }, - { "exp", texp, { &type_mp, &type_mp, &type_mp, &type_mp } }, + { "create", tcreate, { &type_mp, &type_mp, &type_mp, &type_mp, 0 } }, + { "mul", tmul, { &type_mp, &type_mp, &type_mp, &type_mp, 0 } }, + { "exp", texp, { &type_mp, &type_mp, &type_mp, &type_mp, 0 } }, { 0, 0, { 0 } }, };