/* -*-c-*-
*
- * $Id: mpreduce.c,v 1.1 2004/03/27 00:04:46 mdw Exp $
+ * $Id$
*
* Efficient reduction modulo nice primes
*
* MA 02111-1307, USA.
*/
-/*----- Revision history --------------------------------------------------*
- *
- * $Log: mpreduce.c,v $
- * Revision 1.1 2004/03/27 00:04:46 mdw
- * Implement efficient reduction for pleasant-looking primes.
- *
- */
-
/*----- Header files ------------------------------------------------------*/
#include <mLib/darray.h>
* Arguments: @gfreduce *r@ = structure to fill in
* @mp *x@ = an integer
*
- * Returns: ---
+ * Returns: Zero if successful; nonzero on failure.
*
* Use: Initializes a context structure for reduction.
*/
-void mpreduce_create(mpreduce *r, mp *p)
+int mpreduce_create(mpreduce *r, mp *p)
{
mpscan sc;
enum { Z = 0, Z1 = 2, X = 4, X0 = 6 };
/* --- Fill in the easy stuff --- */
- assert(MP_ISPOS(p));
+ if (!MP_POSP(p))
+ return (-1);
d = mp_bits(p);
r->lim = d/MPW_BITS;
r->s = d%MPW_BITS;
}
}
if ((DA(&iv)[DA_LEN(&iv) - 1].op & ~1u) == MPRI_SUB) {
- fprintf(stderr,
- "mpreduce can't handle moduli of the form 2^m + x\n");
- abort();
+ mp_drop(r->p);
+ DA_DESTROY(&iv);
+ return (-1);
}
#undef INSTR
r->iv[i + r->in].argy = b;
}
}
- DA_DESTROY(&iv);
+ DA_DESTROY(&iv);
+ return (0);
}
/* --- @mpreduce_destroy@ --- *
/* --- If source is negative, divide --- */
- if (MP_ISNEG(x)) {
+ if (MP_NEGP(x)) {
mp_div(0, &d, x, r->p);
return (d);
}
mp *spare = (e->f & MP_BURN) ? MP_NEWSEC : MP_NEW;
MP_SHRINK(e);
- if (!MP_LEN(e))
+ MP_COPY(a);
+ if (MP_ZEROP(e))
;
- else if (MP_LEN(e) < EXP_THRESH)
- EXP_SIMPLE(x, a, e);
- else
- EXP_WINDOW(x, a, e);
+ else {
+ if (MP_NEGP(e))
+ a = mp_modinv(a, a, mr->p);
+ if (MP_LEN(e) < EXP_THRESH)
+ EXP_SIMPLE(x, a, e);
+ else
+ EXP_WINDOW(x, a, e);
+ }
+ mp_drop(a);
mp_drop(d);
mp_drop(spare);
return (x);