FOO-merge-N Nth branch merge point
`ec' -- elliptic curve work
-- No merges
++ ec-merge-1 Closed.
-## -*-makefile-*-
+## -*-m4-*-
##
- ## $Id: Makefile.m4,v 1.66 2004/03/21 22:43:50 mdw Exp $
-## $Id: Makefile.m4,v 1.60.2.2 2004/03/21 22:39:46 mdw Exp $
++## $Id: Makefile.m4,v 1.67 2004/03/21 22:52:06 mdw Exp $
##
## Makefile for Catacomb
##
##----- Revision history ----------------------------------------------------
##
## $Log: Makefile.m4,v $
-## Revision 1.60.2.2 2004/03/21 22:39:46 mdw
-## Elliptic curves on binary fields work.
++## Revision 1.67 2004/03/21 22:52:06 mdw
++## Merge and close elliptic curve branch.
+ ##
-## Revision 1.60.2.1 2003/06/10 13:43:53 mdw
-## Simple (non-projective) curves over prime fields now seem to work.
++## Revision 1.60.2.2 2004/03/21 22:39:46 mdw
++## Elliptic curves on binary fields work.
++##
++## Revision 1.60.2.1 2003/06/10 13:43:53 mdw
++## Simple (non-projective) curves over prime fields now seem to work.
++##
+## Revision 1.66 2004/03/21 22:43:50 mdw
+## New hash variant SHA224.
+##
+## Revision 1.65 2003/11/29 23:39:36 mdw
+## Debianization.
+##
+## Revision 1.64 2003/11/10 22:18:30 mdw
+## Build fixes.
+##
+## Revision 1.63 2003/10/17 16:30:46 mdw
+## Report errors if key files don't exist!
+##
+## Revision 1.62 2003/10/12 15:02:09 mdw
+## Reliability fixes.
+##
+## Revision 1.61 2003/10/11 21:02:33 mdw
+## Import buf stuff from tripe.
##
## Revision 1.60 2003/05/16 01:12:37 mdw
## Ship `rc2-tab.h' and `skipjack-tab.h'.
exp.c mpcrt.c mpmul.c mprand.c \
mpbarrett.c mpbarrett-mexp.c mpbarrett-exp.h \
mpmont.c mpmont-mexp.c mpmont-exp.h \
- rho.c \
+ rho.c buf.c \
- GF_SOURCES PGEN_SOURCES')
+ GF_SOURCES PGEN_SOURCES EC_SOURCES')
define(`GF_SOURCES',
- `gfx.c gfx-kmul.c gfx-sqr.c')
+ `gfx.c gfx-kmul.c gfx-sqr.c gf-arith.c gf-gcd.c gfreduce.c')
+
+ define(`EC_SOURCES',
+ `field.c f-prime.c f-binpoly.c ec.c ec-prime.c ec-bin.c')
define(`PGEN_SOURCES',
`pfilt.c rabin.c \
--- /dev/null
- * $Id: ec2.cal,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
+ /* -*-apcalc-*-
+ *
++ * $Id: ec2.cal,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Testbed for elliptic curve arithmetic over binary fields
+ *
+ * (c) 2004 Straylight/Edgeware
+ */
+
+ /*----- Licensing notice --------------------------------------------------*
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ *
+ * Catacomb 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 Library General Public License for more details.
+ *
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+ /*----- Revision history --------------------------------------------------*
+ *
+ * $Log: ec2.cal,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.2.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ * Revision 1.1.4.2 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
+ * Revision 1.1.4.1 2003/06/10 13:43:53 mdw
+ * Simple (non-projective) curves over prime fields now seem to work.
+ *
+ * Revision 1.1 2000/10/08 16:01:37 mdw
+ * Prototypes of various bits of code.
+ *
+ */
+
+ /*----- Object types ------------------------------------------------------*/
+
+ obj ec2_curve { a, b, p };
+ obj ec2_pt { x, y, e };
+ obj ecpp_pt { x, y, z, e };
+
+ /*----- Main code ---------------------------------------------------------*/
+
+ define ec2_curve(a, b, p)
+ {
+ local obj ec2_curve e;
+ e.a = a;
+ e.b = b;
+ e.p = p;
+ return (e);
+ }
+
+ define ec2_pt(x, y, e)
+ {
+ local obj ec2_pt p;
+ p.x = x % e.p;
+ p.y = y % e.p;
+ p.e = e;
+ return (p);
+ }
+
+ define ec2_pt_print(a)
+ {
+ print "(" : a.x : ", " : a.y : ")" :;
+ }
+
+ define ec2_pt_add(a, b)
+ {
+ local e, alpha;
+ local obj ec2_pt d;
+
+ print "> ecadd: ", a, b;
+ if (a == 0)
+ d = b;
+ else if (b == 0)
+ d = a;
+ else if (!istype(a, b))
+ quit "bad type arguments to ec2_pt_add";
+ else if (a.e != b.e)
+ quit "points from different curves in ec2_pt_add";
+ else {
+ e = a.e;
+ if (a.x != b.x) {
+ alpha = ((a.y + b.y) * gf_inv(a.x + b.x, e.p)) % e.p;
+ d.x = (e.a + alpha^2 + alpha + a.x + b.x) % e.p;
+ } else if (a.y != b.y || a.x == gf(0))
+ return 0;
+ else {
+ alpha = a.x + a.y * gf_inv(a.x, e.p) % e.p;
+ d.x = (e.a + alpha^2 + alpha) % e.p;
+ }
+ d.y = ((a.x + d.x) * alpha + d.x + a.y) % e.p;
+ d.e = e;
+ }
+
+ print "< ecadd: ", d;
+ return (d);
+ }
+
+ define ec2_pt_dbl(a)
+ {
+ local e, alpha;
+ local obj ec2_pt d;
+ print "> ecdbl: ", a;
+ if (istype(a, 1))
+ return (0);
+ e = a.e;
+ alpha = a.x + a.y * gf_inv(a.x, e.p) % e.p;
+ d.x = (e.a + alpha^2 + alpha) % e.p;
+ d.y = ((a.x + d.x) * alpha + d.x + a.y) % e.p;
+ d.e = e;
+ print "< ecdbl: ", d;
+ return (d);
+ }
+
+ define ec2_pt_sub(a, b) { return ec2_pt_add(a, ec2_pt_neg(b)); }
+
+ define ec2_pt_neg(a)
+ {
+ local obj ec2_pt d;
+ d.x = a.x;
+ d.y = a.x + a.y;
+ d.e = a.e;
+ return (d);
+ }
+
+ define ec2_pt_check(a)
+ {
+ local e;
+
+ e = a.e;
+ if ((a.y^2 + a.x * a.y) % e.p != (a.x^3 + e.a * a.x^2 + e.b) % e.p)
+ quit "bad curve point";
+ }
+
+ define ec2_pt_mul(a, b)
+ {
+ local p, n;
+ local d;
+
+ if (istype(a, 1)) {
+ n = a;
+ p = b;
+ } else if (istype(b, 1)) {
+ n = b;
+ p = a;
+ } else
+ return (newerror("bad arguments to ec2_pt_mul"));
+
+ d = 0;
+ while (n) {
+ if (n & 1)
+ d += p;
+ n >>= 1;
+ p = ec2_pt_dbl(p);
+ }
+ return (d);
+ }
+
+ /*----- FIPS186-2 standard curves -----------------------------------------*/
+
+ b163 = ec2_curve(gf(1),gf(0x20a601907b8c953ca1481eb10512f78744a3205fd),
+ gf(0x800000000000000000000000000000000000000c9));
+ b163_r = 5846006549323611672814742442876390689256843201587;
+ b163_g = ec2_pt(0x3f0eba16286a2d57ea0991168d4994637e8343e36,
+ 0x0d51fbc6c71a0094fa2cdd545b11c5c0c797324f1, b163);
+
+ /*----- That's all, folks -------------------------------------------------*/
+
/* -*-apcalc-*-
*
- * $Id: ecp.cal,v 1.1 2000/10/08 16:01:37 mdw Exp $
- * $Id: ecp.cal,v 1.1.4.2 2004/03/20 00:13:31 mdw Exp $
++ * $Id: ecp.cal,v 1.2 2004/03/21 22:52:06 mdw Exp $
*
* Testbed for elliptic curve arithmetic over prime fields
*
/*----- Revision history --------------------------------------------------*
*
* $Log: ecp.cal,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.4.2 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
+ * Revision 1.1.4.1 2003/06/10 13:43:53 mdw
+ * Simple (non-projective) curves over prime fields now seem to work.
+ *
* Revision 1.1 2000/10/08 16:01:37 mdw
* Prototypes of various bits of code.
*
/* -*-apcalc-*-
*
- * $Id: gfx.cal,v 1.1 2000/10/08 16:01:37 mdw Exp $
- * $Id: gfx.cal,v 1.1.4.1 2004/03/21 22:39:46 mdw Exp $
++ * $Id: gfx.cal,v 1.2 2004/03/21 22:52:06 mdw Exp $
*
* Testbed for %$\gf{2}$% poltnomial arithmetic
*
/*----- Revision history --------------------------------------------------*
*
* $Log: gfx.cal,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.4.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
* Revision 1.1 2000/10/08 16:01:37 mdw
* Prototypes of various bits of code.
*
-dnl -*-fundamental-*-
+dnl -*-m4-*-
dnl
- dnl $Id: configure.in,v 1.26 2003/11/29 23:39:36 mdw Exp $
-dnl $Id: configure.in,v 1.24.2.1 2003/06/10 13:43:53 mdw Exp $
++dnl $Id: configure.in,v 1.27 2004/03/21 22:52:06 mdw Exp $
dnl
dnl Autoconfiguration for Catacomb
dnl
dnl ----- Revision history --------------------------------------------------
dnl
dnl $Log: configure.in,v $
-dnl Revision 1.24.2.1 2003/06/10 13:43:53 mdw
-dnl Simple (non-projective) curves over prime fields now seem to work.
++dnl Revision 1.27 2004/03/21 22:52:06 mdw
++dnl Merge and close elliptic curve branch.
++dnl
++dnl Revision 1.24.2.1 2003/06/10 13:43:53 mdw
++dnl Simple (non-projective) curves over prime fields now seem to work.
++dnl
+dnl Revision 1.26 2003/11/29 23:39:36 mdw
+dnl Debianization.
+dnl
+dnl Revision 1.25 2003/10/11 21:02:33 mdw
+dnl Import buf stuff from tripe.
dnl
dnl Revision 1.24 2003/05/16 00:30:28 mdw
dnl Version bump.
dnl --- Boring boilerplate ---
AC_INIT(blkc.h)
- mdw_INIT_LIB(catacomb, Catacomb, 2.0.1)
-mdw_INIT_LIB(catacomb, Catacomb, 2.1.0ec1)
++mdw_INIT_LIB(catacomb, Catacomb, 2.1.0)
AM_CONFIG_HEADER(config.h)
dnl --- Make sure I can compile and build libraries ---
--- /dev/null
++catacomb (2.1.0) experimental; urgency=low
++
++ * Added support for elliptic curves, on both prime and binary fields
++ (polynomial basis only). No actual crypto, but there's enough already
++ to do ECDH and stuff on well-known curves Testing is currently a bit
++ patchy.
++
++ -- Mark Wooding <mdw@nsict.org> Sun, 21 Mar 2004 22:47:56 +0000
++
+catacomb (2.0.1) experimental; urgency=low
+
+ * Debianization!
+ * (pixie): Don't report uninteresting errors when accepting connections.
+
+ -- Mark Wooding <mdw@nsict.org> Thu, 11 Dec 2003 10:47:59 +0000
--- /dev/null
- * $Id: ec-bin.c,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
+ /* -*-c-*-
+ *
++ * $Id: ec-bin.c,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Arithmetic for elliptic curves over binary fields
+ *
+ * (c) 2004 Straylight/Edgeware
+ */
+
+ /*----- Licensing notice --------------------------------------------------*
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ *
+ * Catacomb 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 Library General Public License for more details.
+ *
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+ /*----- Revision history --------------------------------------------------*
+ *
+ * $Log: ec-bin.c,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.2.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ */
+
+ /*----- Header files ------------------------------------------------------*/
+
+ #include <mLib/sub.h>
+
+ #include "ec.h"
+
+ /*----- Data structures ---------------------------------------------------*/
+
+ typedef struct ecctx {
+ ec_curve c;
+ mp *a, *b;
+ mp *bb;
+ } ecctx;
+
+ /*----- Main code ---------------------------------------------------------*/
+
+ static const ec_ops ec_binops, ec_binprojops;
+
+ static ec *ecneg(ec_curve *c, ec *d, const ec *p)
+ {
+ EC_COPY(d, p);
+ if (d->x)
+ d->y = F_ADD(c->f, d->y, d->y, d->x);
+ return (d);
+ }
+
+ static ec *ecprojneg(ec_curve *c, ec *d, const ec *p)
+ {
+ EC_COPY(d, p);
+ if (d->x) {
+ mp *t = F_MUL(c->f, MP_NEW, d->x, d->z);
+ d->y = F_ADD(c->f, d->y, d->y, t);
+ MP_DROP(t);
+ }
+ return (d);
+ }
+
+ static ec *ecfind(ec_curve *c, ec *d, mp *x)
+ {
+ /* write me */
+ return (0);
+ }
+
+ static ec *ecdbl(ec_curve *c, ec *d, const ec *a)
+ {
+ if (EC_ATINF(a) || F_ZEROP(c->f, a->x))
+ EC_SETINF(d);
+ else {
+ field *f = c->f;
+ ecctx *cc = (ecctx *)c;
+ mp *lambda;
+ mp *dx, *dy;
+
+ dx = F_INV(f, MP_NEW, a->x); /* %$x^{-1}$% */
+ dy = F_MUL(f, MP_NEW, dx, a->y); /* %$y/x$% */
+ lambda = F_ADD(f, dy, dy, a->x); /* %$\lambda = x + y/x$% */
+
+ dx = F_SQR(f, dx, lambda); /* %$\lambda^2$% */
+ dx = F_ADD(f, dx, dx, lambda); /* %$\lambda^2 + \lambda$% */
+ dx = F_ADD(f, dx, dx, cc->a); /* %$x' = a + \lambda^2 + \lambda$% */
+
+ dy = F_ADD(f, MP_NEW, a->x, dx); /* %$ x + x' $% */
+ dy = F_MUL(f, dy, dy, lambda); /* %$ (x + x') \lambda$% */
+ dy = F_ADD(f, dy, dy, a->y); /* %$ (x + x') \lambda + y$% */
+ dy = F_ADD(f, dy, dy, dx); /* %$ y' = (x + x') \lambda + y + x'$% */
+
+ EC_DESTROY(d);
+ d->x = dx;
+ d->y = dy;
+ d->z = 0;
+ MP_DROP(lambda);
+ }
+ return (d);
+ }
+
+ static ec *ecprojdbl(ec_curve *c, ec *d, const ec *a)
+ {
+ if (EC_ATINF(a) || F_ZEROP(c->f, a->x))
+ EC_SETINF(d);
+ else {
+ field *f = c->f;
+ ecctx *cc = (ecctx *)c;
+ mp *dx, *dy, *dz, *u, *v;
+
+ dy = F_SQR(f, MP_NEW, a->z); /* %$z^2$% */
+ dx = F_MUL(f, MP_NEW, dy, cc->bb); /* %$c z^2$% */
+ dx = F_ADD(f, dx, dx, a->x); /* %$x + c z^2$% */
+ dz = F_SQR(f, MP_NEW, dx); /* %$(x + c z^2)^2$% */
+ dx = F_SQR(f, dx, dz); /* %$x' = (x + c z^2)^4$% */
+
+ dz = F_MUL(f, dz, dy, a->x); /* %$z' = x z^2$% */
+
+ dy = F_SQR(f, dy, a->x); /* %$x^2$% */
+ u = F_MUL(f, MP_NEW, a->y, a->z); /* %$y z$% */
+ u = F_ADD(f, u, u, dz); /* %$z' + y z$% */
+ u = F_ADD(f, u, u, dy); /* %$u = z' + x^2 + y z$% */
+
+ v = F_SQR(f, MP_NEW, dy); /* %$x^4$% */
+ dy = F_MUL(f, dy, v, dz); /* %$x^4 z'$% */
+ v = F_MUL(f, v, u, dx); /* %$u x'$% */
+ dy = F_ADD(f, dy, dy, v); /* %$y' = x^4 z' + u x'$% */
+
+ EC_DESTROY(d);
+ d->x = dx;
+ d->y = dy;
+ d->z = dz;
+ MP_DROP(u);
+ MP_DROP(v);
+ assert(!(d->x->f & MP_DESTROYED));
+ assert(!(d->y->f & MP_DESTROYED));
+ assert(!(d->z->f & MP_DESTROYED));
+ }
+ return (d);
+ }
+
+ static ec *ecadd(ec_curve *c, ec *d, const ec *a, const ec *b)
+ {
+ if (a == b)
+ ecdbl(c, d, a);
+ else if (EC_ATINF(a))
+ EC_COPY(d, b);
+ else if (EC_ATINF(b))
+ EC_COPY(d, a);
+ else {
+ field *f = c->f;
+ ecctx *cc = (ecctx *)c;
+ mp *lambda;
+ mp *dx, *dy;
+
+ if (!MP_EQ(a->x, b->x)) {
+ dx = F_ADD(f, MP_NEW, a->x, b->x); /* %$x_0 + x_1$% */
+ dy = F_INV(f, MP_NEW, dx); /* %$(x_0 + x_1)^{-1}$% */
+ dx = F_ADD(f, dx, a->y, b->y); /* %$y_0 + y_1$% */
+ lambda = F_MUL(f, MP_NEW, dy, dx);
+ /* %$\lambda = (y_0 + y_1)/(x_0 + x_1)$% */
+
+ dx = F_SQR(f, dx, lambda); /* %$\lambda^2$% */
+ dx = F_ADD(f, dx, dx, lambda); /* %$\lambda^2 + \lambda$% */
+ dx = F_ADD(f, dx, dx, cc->a); /* %$a + \lambda^2 + \lambda$% */
+ dx = F_ADD(f, dx, dx, a->x); /* %$a + \lambda^2 + \lambda + x_0$% */
+ dx = F_ADD(f, dx, dx, b->x);
+ /* %$x' = a + \lambda^2 + \lambda + x_0 + x_1$% */
+ } else if (!MP_EQ(a->y, b->y) || F_ZEROP(f, a->x)) {
+ EC_SETINF(d);
+ return (d);
+ } else {
+ dx = F_INV(f, MP_NEW, a->x); /* %$x^{-1}$% */
+ dy = F_MUL(f, MP_NEW, dx, a->y); /* %$y/x$% */
+ lambda = F_ADD(f, dy, dy, a->x); /* %$\lambda = x + y/x$% */
+
+ dx = F_SQR(f, dx, lambda); /* %$\lambda^2$% */
+ dx = F_ADD(f, dx, dx, lambda); /* %$\lambda^2 + \lambda$% */
+ dx = F_ADD(f, dx, dx, cc->a); /* %$x' = a + \lambda^2 + \lambda$% */
+ dy = MP_NEW;
+ }
+
+ dy = F_ADD(f, dy, a->x, dx); /* %$ x + x' $% */
+ dy = F_MUL(f, dy, dy, lambda); /* %$ (x + x') \lambda$% */
+ dy = F_ADD(f, dy, dy, a->y); /* %$ (x + x') \lambda + y$% */
+ dy = F_ADD(f, dy, dy, dx); /* %$ y' = (x + x') \lambda + y + x'$% */
+
+ EC_DESTROY(d);
+ d->x = dx;
+ d->y = dy;
+ d->z = 0;
+ MP_DROP(lambda);
+ }
+ return (d);
+ }
+
+ static ec *ecprojadd(ec_curve *c, ec *d, const ec *a, const ec *b)
+ {
+ if (a == b)
+ c->ops->dbl(c, d, a);
+ else if (EC_ATINF(a))
+ EC_COPY(d, b);
+ else if (EC_ATINF(b))
+ EC_COPY(d, a);
+ else {
+ field *f = c->f;
+ ecctx *cc = (ecctx *)c;
+ mp *dx, *dy, *dz, *u, *uu, *v, *t, *s, *ss, *r, *w, *l;
+
+ dz = F_SQR(f, MP_NEW, b->z); /* %$z_1^2$% */
+ u = F_MUL(f, MP_NEW, dz, a->x); /* %$u_0 = x_0 z_1^2$% */
+ t = F_MUL(f, MP_NEW, dz, b->z); /* %$z_1^3$% */
+ s = F_MUL(f, MP_NEW, t, a->y); /* %$s_0 = y_0 z_1^3$% */
+
+ dz = F_SQR(f, dz, a->z); /* %$z_0^2$% */
+ uu = F_MUL(f, MP_NEW, dz, b->x); /* %$u_1 = x_1 z_0^2$% */
+ t = F_MUL(f, t, dz, a->z); /* %$z_0^3$% */
+ ss = F_MUL(f, MP_NEW, t, b->y); /* %$s_1 = y_1 z_0^3$% */
+
+ w = F_ADD(f, u, u, uu); /* %$r = u_0 + u_1$% */
+ r = F_ADD(f, s, s, ss); /* %$w = s_0 + s_1$% */
+ if (F_ZEROP(f, w)) {
+ MP_DROP(w);
+ MP_DROP(uu);
+ MP_DROP(ss);
+ MP_DROP(t);
+ MP_DROP(dz);
+ if (F_ZEROP(f, r)) {
+ MP_DROP(r);
+ return (c->ops->dbl(c, d, a));
+ } else {
+ MP_DROP(r);
+ EC_SETINF(d);
+ return (d);
+ }
+ }
+
+ l = F_MUL(f, t, a->z, w); /* %$l = z_0 w$% */
+
+ dz = F_MUL(f, dz, l, b->z); /* %$z' = l z_1$% */
+
+ ss = F_MUL(f, ss, r, b->x); /* %$r x_1$% */
+ t = F_MUL(f, uu, l, b->y); /* %$l y_1$% */
+ v = F_ADD(f, ss, ss, t); /* %$v = r x_1 + l y_1$% */
+
+ t = F_ADD(f, t, r, dz); /* %$t = r + z'$% */
+
+ uu = F_SQR(f, MP_NEW, dz); /* %$z'^2$% */
+ dx = F_MUL(f, MP_NEW, uu, cc->a); /* %$a z'^2$% */
+ uu = F_MUL(f, uu, t, r); /* %$t r$% */
+ dx = F_ADD(f, dx, dx, uu); /* %$a z'^2 + t r$% */
+ r = F_SQR(f, r, w); /* %$w^2$% */
+ uu = F_MUL(f, uu, r, w); /* %$w^3$% */
+ dx = F_ADD(f, dx, dx, uu); /* %$x' = a z'^2 + t r + w^3$% */
+
+ r = F_SQR(f, r, l); /* %$l^2$% */
+ dy = F_MUL(f, uu, v, r); /* %$v l^2$% */
+ l = F_MUL(f, l, t, dx); /* %$t x'$% */
+ dy = F_ADD(f, dy, dy, l); /* %$y' = t x' + v l^2$% */
+
+ EC_DESTROY(d);
+ d->x = dx;
+ d->y = dy;
+ d->z = dz;
+ MP_DROP(l);
+ MP_DROP(r);
+ MP_DROP(w);
+ MP_DROP(t);
+ MP_DROP(v);
+ }
+ return (d);
+ }
+
+ static int eccheck(ec_curve *c, const ec *p)
+ {
+ ecctx *cc = (ecctx *)c;
+ field *f = c->f;
+ int rc;
+ mp *u, *v;
+
+ v = F_SQR(f, MP_NEW, p->x);
+ u = F_MUL(f, MP_NEW, v, p->x);
+ v = F_MUL(f, v, v, cc->a);
+ u = F_ADD(f, u, u, v);
+ u = F_ADD(f, u, u, cc->b);
+ v = F_MUL(f, v, p->x, p->y);
+ u = F_ADD(f, u, u, v);
+ v = F_SQR(f, v, p->y);
+ u = F_ADD(f, u, u, v);
+ rc = F_ZEROP(f, u);
+ mp_drop(u);
+ mp_drop(v);
+ return (rc);
+ }
+
+ static int ecprojcheck(ec_curve *c, const ec *p)
+ {
+ ec t = EC_INIT;
+ int rc;
+
+ c->ops->fix(c, &t, p);
+ rc = eccheck(c, &t);
+ EC_DESTROY(&t);
+ return (rc);
+ }
+
+ static void ecdestroy(ec_curve *c)
+ {
+ ecctx *cc = (ecctx *)c;
+ MP_DROP(cc->a);
+ MP_DROP(cc->b);
+ if (cc->bb) MP_DROP(cc->bb);
+ DESTROY(cc);
+ }
+
+ /* --- @ec_bin@, @ec_binproj@ --- *
+ *
+ * Arguments: @field *f@ = the underlying field for this elliptic curve
+ * @mp *a, *b@ = the coefficients for this curve
+ *
+ * Returns: A pointer to the curve.
+ *
+ * Use: Creates a curve structure for an elliptic curve defined over
+ * a binary field. The @binproj@ variant uses projective
+ * coordinates, which can be a win.
+ */
+
+ ec_curve *ec_bin(field *f, mp *a, mp *b)
+ {
+ ecctx *cc = CREATE(ecctx);
+ cc->c.ops = &ec_binops;
+ cc->c.f = f;
+ cc->a = F_IN(f, MP_NEW, a);
+ cc->b = F_IN(f, MP_NEW, b);
+ cc->bb = 0;
+ return (&cc->c);
+ }
+
+ ec_curve *ec_binproj(field *f, mp *a, mp *b)
+ {
+ ecctx *cc = CREATE(ecctx);
+ cc->c.ops = &ec_binprojops;
+ cc->c.f = f;
+ cc->a = F_IN(f, MP_NEW, a);
+ cc->b = F_IN(f, MP_NEW, b);
+ cc->bb = F_SQRT(f, MP_NEW, b);
+ cc->bb = F_SQRT(f, cc->bb, cc->bb);
+ return (&cc->c);
+ }
+
+ static const ec_ops ec_binops = {
+ ecdestroy, ec_idin, ec_idout, ec_idfix,
+ 0, ecneg, ecadd, ec_stdsub, ecdbl, eccheck
+ };
+
+ static const ec_ops ec_binprojops = {
+ ecdestroy, ec_projin, ec_projout, ec_projfix,
+ 0, ecprojneg, ecprojadd, ec_stdsub, ecprojdbl, ecprojcheck
+ };
+
+ /*----- Test rig ----------------------------------------------------------*/
+
+ #ifdef TEST_RIG
+
+ #define MP(x) mp_readstring(MP_NEW, #x, 0, 0)
+
+ int main(int argc, char *argv[])
+ {
+ field *f;
+ ec_curve *c;
+ ec g = EC_INIT, d = EC_INIT;
+ mp *p, *a, *b, *r;
+ int i, n = argc == 1 ? 1 : atoi(argv[1]);
+
+ printf("ec-bin: ");
+ fflush(stdout);
+ a = MP(1);
+ b = MP(0x066647ede6c332c7f8c0923bb58213b333b20e9ce4281fe115f7d8f90ad);
+ p = MP(0x20000000000000000000000000000000000000004000000000000000001);
+ r =
+ MP(6901746346790563787434755862277025555839812737345013555379383634485462);
+
+ f = field_binpoly(p);
+ c = ec_binproj(f, a, b);
+
+ g.x = MP(0x0fac9dfcbac8313bb2139f1bb755fef65bc391f8b36f8f8eb7371fd558b);
+ g.y = MP(0x1006a08a41903350678e58528bebf8a0beff867a7ca36716f7e01f81052);
+
+ for (i = 0; i < n; i++) {
+ ec_mul(c, &d, &g, r);
+ if (EC_ATINF(&d)) {
+ fprintf(stderr, "zero too early\n");
+ return (1);
+ }
+ ec_add(c, &d, &d, &g);
+ if (!EC_ATINF(&d)) {
+ fprintf(stderr, "didn't reach zero\n");
+ MP_EPRINTX("d.x", d.x);
+ MP_EPRINTX("d.y", d.y);
+ MP_EPRINTX("d.z", d.y);
+ return (1);
+ }
+ ec_destroy(&d);
+ }
+
+ ec_destroy(&g);
+ ec_destroycurve(c);
+ F_DESTROY(f);
+ MP_DROP(p); MP_DROP(a); MP_DROP(b); MP_DROP(r);
+ assert(!mparena_count(&mparena_global));
+ printf("ok\n");
+ return (0);
+ }
+
+ #endif
+
+ /*----- That's all, folks -------------------------------------------------*/
/* -*-c-*-
*
- * $Id: ec-exp.h,v 1.2 2003/05/15 23:25:59 mdw Exp $
- * $Id: ec-exp.h,v 1.2.4.1 2004/03/20 00:13:31 mdw Exp $
++ * $Id: ec-exp.h,v 1.3 2004/03/21 22:52:06 mdw Exp $
*
* Exponentiation operations for elliptic curves
*
/*----- Revision history --------------------------------------------------*
*
* $Log: ec-exp.h,v $
++ * Revision 1.3 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.2.4.1 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
* Revision 1.2 2003/05/15 23:25:59 mdw
* Make elliptic curve stuff build.
*
/* -*-c-*-
*
- * $Id: ec-prime.c,v 1.3 2003/05/15 23:25:59 mdw Exp $
- * $Id: ec-prime.c,v 1.3.4.3 2004/03/21 22:39:46 mdw Exp $
++ * $Id: ec-prime.c,v 1.4 2004/03/21 22:52:06 mdw Exp $
*
* Elliptic curves over prime fields
*
/*----- Revision history --------------------------------------------------*
*
* $Log: ec-prime.c,v $
++ * Revision 1.4 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.3.4.3 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ * Revision 1.3.4.2 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
+ * Revision 1.3.4.1 2003/06/10 13:43:53 mdw
+ * Simple (non-projective) curves over prime fields now seem to work.
+ *
* Revision 1.3 2003/05/15 23:25:59 mdw
* Make elliptic curve stuff build.
*
/* -*-c-*-
*
- * $Id: ec.c,v 1.4 2003/05/15 23:25:59 mdw Exp $
- * $Id: ec.c,v 1.4.4.2 2004/03/20 00:13:31 mdw Exp $
++ * $Id: ec.c,v 1.5 2004/03/21 22:52:06 mdw Exp $
*
* Elliptic curve definitions
*
/*----- Revision history --------------------------------------------------*
*
* $Log: ec.c,v $
++ * Revision 1.5 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.4.4.2 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
+ * Revision 1.4.4.1 2003/06/10 13:43:53 mdw
+ * Simple (non-projective) curves over prime fields now seem to work.
+ *
* Revision 1.4 2003/05/15 23:25:59 mdw
* Make elliptic curve stuff build.
*
/* -*-c-*-
*
- * $Id: ec.h,v 1.4 2003/05/15 23:25:59 mdw Exp $
- * $Id: ec.h,v 1.4.4.3 2004/03/21 22:39:46 mdw Exp $
++ * $Id: ec.h,v 1.5 2004/03/21 22:52:06 mdw Exp $
*
* Elliptic curve definitions
*
/*----- Revision history --------------------------------------------------*
*
* $Log: ec.h,v $
++ * Revision 1.5 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.4.4.3 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ * Revision 1.4.4.2 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
+ * Revision 1.4.4.1 2003/06/10 13:43:53 mdw
+ * Simple (non-projective) curves over prime fields now seem to work.
+ *
* Revision 1.4 2003/05/15 23:25:59 mdw
* Make elliptic curve stuff build.
*
/* -*-c-*-
*
- * $Id: exp.h,v 1.1 2001/06/16 13:00:59 mdw Exp $
- * $Id: exp.h,v 1.1.4.1 2004/03/20 00:13:31 mdw Exp $
++ * $Id: exp.h,v 1.2 2004/03/21 22:52:06 mdw Exp $
*
* Generalized exponentiation
*
/*----- Revision history --------------------------------------------------*
*
* $Log: exp.h,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.4.1 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
* Revision 1.1 2001/06/16 13:00:59 mdw
* New generic exponentation code. Includes sliding-window simultaneous
* exponentiation.
--- /dev/null
- * $Id: f-binpoly.c,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
+ /* -*-c-*-
+ *
++ * $Id: f-binpoly.c,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Binary fields with polynomial basis representation
+ *
+ * (c) 2004 Straylight/Edgeware
+ */
+
+ /*----- Licensing notice --------------------------------------------------*
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ *
+ * Catacomb 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 Library General Public License for more details.
+ *
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+ /*----- Revision history --------------------------------------------------*
+ *
+ * $Log: f-binpoly.c,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.2.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ */
+
+ /*----- Header files ------------------------------------------------------*/
+
+ #include <mLib/sub.h>
+
+ #include "field.h"
+ #include "gf.h"
+ #include "gfreduce.h"
+
+ /*----- Data structures ---------------------------------------------------*/
+
+ typedef struct fctx {
+ field f;
+ gfreduce r;
+ } fctx;
+
+ /*----- Main code ---------------------------------------------------------*/
+
+ /* --- Field operations --- */
+
+ static void fdestroy(field *ff)
+ {
+ fctx *f = (fctx *)ff;
+ gfreduce_destroy(&f->r);
+ DESTROY(f);
+ }
+
+ static int fzerop(field *ff, mp *x)
+ {
+ return (!MP_LEN(x));
+ }
+
+ static mp *fadd(field *ff, mp *d, mp *x, mp *y)
+ {
+ return (gf_add(d, x, y));
+ }
+
+ static mp *fmul(field *ff, mp *d, mp *x, mp *y)
+ {
+ fctx *f = (fctx *)ff;
+ d = gf_mul(d, x, y);
+ return (gfreduce_do(&f->r, d, d));
+ }
+
+ static mp *fsqr(field *ff, mp *d, mp *x)
+ {
+ fctx *f = (fctx *)ff;
+ d = gf_sqr(d, x);
+ return (gfreduce_do(&f->r, d, d));
+ }
+
+ static mp *finv(field *ff, mp *d, mp *x)
+ {
+ fctx *f = (fctx *)ff;
+ gf_gcd(0, 0, &d, f->r.p, x);
+ return (d);
+ }
+
+ static mp *freduce(field *ff, mp *d, mp *x)
+ {
+ fctx *f = (fctx *)ff;
+ return (gfreduce_do(&f->r, d, x));
+ }
+
+ static mp *fsqrt(field *ff, mp *d, mp *x)
+ {
+ fctx *f = (fctx *)ff;
+ return (gfreduce_sqrt(&f->r, d, x));
+ }
+
+ static mp *fquadsolve(field *ff, mp *d, mp *x)
+ {
+ fctx *f = (fctx *)ff;
+ return (gfreduce_quadsolve(&f->r, d, x));
+ }
+
+ /* --- Field operations table --- */
+
+ static field_ops fops = {
+ fdestroy,
+ freduce, field_id,
+ fzerop, field_id, fadd, fadd, fmul, fsqr, finv, freduce, fsqrt,
+ fquadsolve,
+ 0, 0, 0, 0
+ };
+
+ /* --- @field_binpoly@ --- *
+ *
+ * Arguments: @mp *p@ = the reduction polynomial
+ *
+ * Returns: A pointer to the field.
+ *
+ * Use: Creates a field structure for a binary field mod @p@.
+ */
+
+ field *field_binpoly(mp *p)
+ {
+ fctx *f = CREATE(fctx);
+ f->f.ops = &fops;
+ f->f.zero = MP_ZERO;
+ f->f.one = MP_ONE;
+ gfreduce_create(&f->r, p);
+ return (&f->f);
+ }
+
+ /*----- That's all, folks -------------------------------------------------*/
/* -*-c-*-
*
- * $Id: f-prime.c,v 1.3 2003/05/15 23:25:59 mdw Exp $
- * $Id: f-prime.c,v 1.3.4.3 2004/03/21 22:39:46 mdw Exp $
++ * $Id: f-prime.c,v 1.4 2004/03/21 22:52:06 mdw Exp $
*
* Prime fields with Montgomery arithmetic
*
/*----- Revision history --------------------------------------------------*
*
* $Log: f-prime.c,v $
++ * Revision 1.4 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.3.4.3 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ * Revision 1.3.4.2 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
+ * Revision 1.3.4.1 2003/06/10 13:43:53 mdw
+ * Simple (non-projective) curves over prime fields now seem to work.
+ *
* Revision 1.3 2003/05/15 23:25:59 mdw
* Make elliptic curve stuff build.
*
/* -*-c-*-
*
- * $Id: field.c,v 1.1 2001/05/07 17:30:13 mdw Exp $
- * $Id: field.c,v 1.1.4.1 2003/06/10 13:43:53 mdw Exp $
++ * $Id: field.c,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Abstract field operations
*
- * [Abstract field operations *
* (c) 2001 Straylight/Edgeware
*/
/*----- Revision history --------------------------------------------------*
*
* $Log: field.c,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.4.1 2003/06/10 13:43:53 mdw
+ * Simple (non-projective) curves over prime fields now seem to work.
+ *
* Revision 1.1 2001/05/07 17:30:13 mdw
* Add an internal-representation no-op function.
*
/* -*-c-*-
*
- * $Id: field.h,v 1.3 2002/01/13 13:48:44 mdw Exp $
- * $Id: field.h,v 1.3.4.2 2004/03/21 22:39:46 mdw Exp $
++ * $Id: field.h,v 1.4 2004/03/21 22:52:06 mdw Exp $
*
* Definitions for field arithmetic
*
/*----- Revision history --------------------------------------------------*
*
* $Log: field.h,v $
++ * Revision 1.4 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.3.4.2 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ * Revision 1.3.4.1 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
* Revision 1.3 2002/01/13 13:48:44 mdw
* Further progress.
*
--- /dev/null
- * $Id: gf-arith.c,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
+ /* -*-c-*-
+ *
++ * $Id: gf-arith.c,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Basic arithmetic on binary polynomials
+ *
+ * (c) 2004 Straylight/Edgeware
+ */
+
+ /*----- Licensing notice --------------------------------------------------*
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ *
+ * Catacomb 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 Library General Public License for more details.
+ *
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+ /*----- Revision history --------------------------------------------------*
+ *
+ * $Log: gf-arith.c,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.2.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ */
+
+ /*----- Header files ------------------------------------------------------*/
+
+ #include "gf.h"
+
+ /*----- Macros ------------------------------------------------------------*/
+
+ #define MAX(x, y) ((x) >= (y) ? (x) : (y))
+
+ /*----- Main code ---------------------------------------------------------*/
+
+ /* --- @gf_add@ --- *
+ *
+ * Arguments: @mp *d@ = destination
+ * @mp *a, *b@ = sources
+ *
+ * Returns: Result, @a@ added to @b@.
+ */
+
+ mp *gf_add(mp *d, mp *a, mp *b)
+ {
+ MP_DEST(d, MAX(MP_LEN(a), MP_LEN(b)), (a->f | b->f) & MP_BURN);
+ gfx_add(d->v, d->vl, a->v, a->vl, b->v, b->vl);
+ d->f = (a->f | b->f) & MP_BURN;
+ MP_SHRINK(d);
+ return (d);
+ }
+
+ /* --- @gf_mul@ --- *
+ *
+ * Arguments: @mp *d@ = destination
+ * @mp *a, *b@ = sources
+ *
+ * Returns: Result, @a@ multiplied by @b@.
+ */
+
+ mp *gf_mul(mp *d, mp *a, mp *b)
+ {
+ a = MP_COPY(a);
+ b = MP_COPY(b);
+
+ if (MP_LEN(a) <= MPK_THRESH || MP_LEN(b) <= GFK_THRESH) {
+ MP_DEST(d, MP_LEN(a) + MP_LEN(b), a->f | b->f | MP_UNDEF);
+ gfx_mul(d->v, d->vl, a->v, a->vl, b->v, b->vl);
+ } else {
+ size_t m = MAX(MP_LEN(a), MP_LEN(b));
+ mpw *s;
+ MP_DEST(d, 2 * m, a->f | b->f | MP_UNDEF);
+ s = mpalloc(d->a, 2 * m);
+ gfx_kmul(d->v, d->vl, a->v, a->vl, b->v, b->vl, s, s + 2 * m);
+ mpfree(d->a, s);
+ }
+
+ d->f = (a->f | b->f) & MP_BURN;
+ MP_SHRINK(d);
+ MP_DROP(a);
+ MP_DROP(b);
+ return (d);
+ }
+
+ /* --- @gf_sqr@ --- *
+ *
+ * Arguments: @mp *d@ = destination
+ * @mp *a@ = source
+ *
+ * Returns: Result, @a@ squared.
+ */
+
+ mp *gf_sqr(mp *d, mp *a)
+ {
+ MP_COPY(a);
+ MP_DEST(d, 2 * MP_LEN(a), a->f & MP_BURN);
+ gfx_sqr(d->v, d->vl, a->v, a->vl);
+ d->f = a->f & MP_BURN;
+ MP_SHRINK(d);
+ MP_DROP(a);
+ return (d);
+ }
+
+ /* --- @gf_div@ --- *
+ *
+ * Arguments: @mp **qq, **rr@ = destination, quotient and remainder
+ * @mp *a, *b@ = sources
+ *
+ * Use: Calculates the quotient and remainder when @a@ is divided by
+ * @b@. The destinations @*qq@ and @*rr@ must be distinct.
+ * Either of @qq@ or @rr@ may be null to indicate that the
+ * result is irrelevant. (Discarding both results is silly.)
+ * There is a performance advantage if @a == *rr@.
+ */
+
+ void gf_div(mp **qq, mp **rr, mp *a, mp *b)
+ {
+ mp *r = rr ? *rr : MP_NEW;
+ mp *q = qq ? *qq : MP_NEW;
+
+ /* --- Set the remainder up right --- */
+
+ b = MP_COPY(b);
+ a = MP_COPY(a);
+ if (r)
+ MP_DROP(r);
+ r = a;
+ MP_DEST(r, MP_LEN(b) + 2, a->f | b->f);
+
+ /* --- Fix up the quotient too --- */
+
+ r = MP_COPY(r);
+ MP_DEST(q, MP_LEN(r), r->f | MP_UNDEF);
+ MP_DROP(r);
+
+ /* --- Perform the calculation --- */
+
+ gfx_div(q->v, q->vl, r->v, r->vl, b->v, b->vl);
+
+ /* --- Sort out the sign of the results --- *
+ *
+ * If the signs of the arguments differ, and the remainder is nonzero, I
+ * must add one to the absolute value of the quotient and subtract the
+ * remainder from @b@.
+ */
+
+ q->f = (r->f | b->f) & MP_BURN;
+ r->f = (r->f | b->f) & MP_BURN;
+
+ /* --- Store the return values --- */
+
+ MP_DROP(b);
+
+ if (!qq)
+ MP_DROP(q);
+ else {
+ MP_SHRINK(q);
+ *qq = q;
+ }
+
+ if (!rr)
+ MP_DROP(r);
+ else {
+ MP_SHRINK(r);
+ *rr = r;
+ }
+ }
+
+ /*----- Test rig ----------------------------------------------------------*/
+
+ #ifdef TEST_RIG
+
+ static int verify(const char *op, mp *expect, mp *result, mp *a, mp *b)
+ {
+ if (!MP_EQ(expect, result)) {
+ fprintf(stderr, "\n*** %s failed", op);
+ fputs("\n*** a = ", stderr); mp_writefile(a, stderr, 16);
+ fputs("\n*** b = ", stderr); mp_writefile(b, stderr, 16);
+ fputs("\n*** result = ", stderr); mp_writefile(result, stderr, 16);
+ fputs("\n*** expect = ", stderr); mp_writefile(expect, stderr, 16);
+ fputc('\n', stderr);
+ return (0);
+ }
+ return (1);
+ }
+
+ #define RIG(name, op) \
+ static int t##name(dstr *v) \
+ { \
+ mp *a = *(mp **)v[0].buf; \
+ mp *b = *(mp **)v[1].buf; \
+ mp *r = *(mp **)v[2].buf; \
+ mp *c = op(MP_NEW, a, b); \
+ int ok = verify(#name, r, c, a, b); \
+ mp_drop(a); mp_drop(b); mp_drop(c); mp_drop(r); \
+ assert(mparena_count(MPARENA_GLOBAL) == 0); \
+ return (ok); \
+ }
+
+ RIG(add, gf_add)
+ RIG(mul, gf_mul)
+
+ #undef RIG
+
+ static int tsqr(dstr *v)
+ {
+ mp *a = *(mp **)v[0].buf;
+ mp *r = *(mp **)v[1].buf;
+ mp *c = MP_NEW;
+ int ok = 1;
+ c = gf_sqr(MP_NEW, a);
+ ok &= verify("sqr", r, c, a, MP_ZERO);
+ mp_drop(a); mp_drop(r); mp_drop(c);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
+ return (ok);
+ }
+
+ static int tdiv(dstr *v)
+ {
+ mp *a = *(mp **)v[0].buf;
+ mp *b = *(mp **)v[1].buf;
+ mp *q = *(mp **)v[2].buf;
+ mp *r = *(mp **)v[3].buf;
+ mp *c = MP_NEW, *d = MP_NEW;
+ int ok = 1;
+ gf_div(&c, &d, a, b);
+ ok &= verify("div(quotient)", q, c, a, b);
+ ok &= verify("div(remainder)", r, d, a, b);
+ mp_drop(a); mp_drop(b); mp_drop(c); mp_drop(d); mp_drop(r); mp_drop(q);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
+ return (ok);
+ }
+
+ static test_chunk tests[] = {
+ { "add", tadd, { &type_mp, &type_mp, &type_mp, 0 } },
+ { "mul", tmul, { &type_mp, &type_mp, &type_mp, 0 } },
+ { "sqr", tsqr, { &type_mp, &type_mp, 0 } },
+ { "div", tdiv, { &type_mp, &type_mp, &type_mp, &type_mp, 0 } },
+ { 0, 0, { 0 } },
+ };
+
+ int main(int argc, char *argv[])
+ {
+ sub_init();
+ test_run(argc, argv, tests, SRCDIR "/tests/gf");
+ return (0);
+ }
+
+ #endif
+
+ /*----- That's all, folks -------------------------------------------------*/
--- /dev/null
- * $Id: gf-gcd.c,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
+ /* -*-c-*-
+ *
++ * $Id: gf-gcd.c,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Euclidian algorithm on binary polynomials
+ *
+ * (c) 2004 Straylight/Edgeware
+ */
+
+ /*----- Licensing notice --------------------------------------------------*
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ *
+ * Catacomb 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 Library General Public License for more details.
+ *
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+ /*----- Revision history --------------------------------------------------*
+ *
+ * $Log: gf-gcd.c,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.2.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ */
+
+ /*----- Header files ------------------------------------------------------*/
+
+ #include "gf.h"
+
+ /*----- Main code ---------------------------------------------------------*/
+
+ /* --- @gf_gcd@ --- *
+ *
+ * Arguments: @mp **gcd, **xx, **yy@ = where to write the results
+ * @mp *a, *b@ = sources (must be nonzero)
+ *
+ *
+ * Returns: ---
+ *
+ * Use: Calculates @gcd(a, b)@, and two numbers @x@ and @y@ such that
+ * @ax + by = gcd(a, b)@. This is useful for computing modular
+ * inverses.
+ */
+
+ void gf_gcd(mp **gcd, mp **xx, mp **yy, mp *a, mp *b)
+ {
+ mp *x = MP_ONE, *X = MP_ZERO;
+ mp *y = MP_ZERO, *Y = MP_ONE;
+ mp *u, *v;
+ mp *q = MP_NEW;
+ unsigned f = 0;
+
+ #define f_swap 1u
+ #define f_ext 2u
+
+ /* --- Sort out some initial flags --- */
+
+ if (xx || yy)
+ f |= f_ext;
+
+ /* --- Ensure that @a@ is larger than @b@ --- *
+ *
+ * If they're the same length we don't care which order they're in, so this
+ * unsigned comparison is fine.
+ */
+
+ if (MPX_UCMP(a->v, a->vl, <, b->v, b->vl)) {
+ { mp *t = a; a = b; b = t; }
+ f |= f_swap;
+ }
+
+ /* --- Check for zeroness --- */
+
+ if (MP_EQ(b, MP_ZERO)) {
+
+ /* --- Store %$|a|$% as the GCD --- */
+
+ if (gcd) {
+ if (*gcd) MP_DROP(*gcd);
+ a = MP_COPY(a);
+ *gcd = a;
+ }
+
+ /* --- Store %$1$% and %$0$% in the appropriate bins --- */
+
+ if (f & f_ext) {
+ if (f & f_swap) {
+ mp **t = xx; xx = yy; yy = t;
+ }
+ if (xx) {
+ if (*xx) MP_DROP(*xx);
+ if (MP_EQ(a, MP_ZERO))
+ *xx = MP_ZERO;
+ else
+ *xx = MP_ONE;
+ }
+ if (yy) {
+ if (*yy) MP_DROP(*yy);
+ *yy = MP_ZERO;
+ }
+ }
+ return;
+ }
+
+ /* --- Take a reference to the arguments --- */
+
+ a = MP_COPY(a);
+ b = MP_COPY(b);
+
+ /* --- Make sure @a@ and @b@ are not both even --- */
+
+ MP_SPLIT(a); a->f &= ~MP_NEG;
+ MP_SPLIT(b); b->f &= ~MP_NEG;
+
+ u = MP_COPY(a);
+ v = MP_COPY(b);
+
+ while (MP_LEN(v)) {
+ mp *t;
+ gf_div(&q, &u, u, v);
+ if (f & f_ext) {
+ t = gf_mul(MP_NEW, X, q);
+ t = gf_add(t, x, t);
+ MP_DROP(x); x = X; X = t;
+ t = gf_mul(MP_NEW, Y, q);
+ t = gf_add(t, y, t);
+ MP_DROP(y); y = Y; Y = t;
+ }
+ t = u; u = v; v = t;
+ }
+
+ MP_DROP(q);
+ if (!gcd)
+ MP_DROP(u);
+ else {
+ if (*gcd) MP_DROP(*gcd);
+ u->f &= ~MP_NEG;
+ *gcd = u;
+ }
+
+ /* --- Perform a little normalization --- */
+
+ if (f & f_ext) {
+
+ /* --- If @a@ and @b@ got swapped, swap the coefficients back --- */
+
+ if (f & f_swap) {
+ mp *t = x; x = y; y = t;
+ t = a; a = b; b = t;
+ }
+
+ /* --- Store the results --- */
+
+ if (!xx)
+ MP_DROP(x);
+ else {
+ if (*xx) MP_DROP(*xx);
+ *xx = x;
+ }
+
+ if (!yy)
+ MP_DROP(y);
+ else {
+ if (*yy) MP_DROP(*yy);
+ *yy = y;
+ }
+ }
+
+ MP_DROP(v);
+ MP_DROP(X); MP_DROP(Y);
+ MP_DROP(a); MP_DROP(b);
+ }
+
+ /*----- Test rig ----------------------------------------------------------*/
+
+ #ifdef TEST_RIG
+
+ static int gcd(dstr *v)
+ {
+ int ok = 1;
+ mp *a = *(mp **)v[0].buf;
+ mp *b = *(mp **)v[1].buf;
+ mp *g = *(mp **)v[2].buf;
+ mp *x = *(mp **)v[3].buf;
+ mp *y = *(mp **)v[4].buf;
+
+ mp *gg = MP_NEW, *xx = MP_NEW, *yy = MP_NEW;
+ gf_gcd(&gg, &xx, &yy, a, b);
+ if (!MP_EQ(x, xx)) {
+ fputs("\n*** mp_gcd(x) failed", stderr);
+ fputs("\na = ", stderr); mp_writefile(a, stderr, 16);
+ fputs("\nb = ", stderr); mp_writefile(b, stderr, 16);
+ fputs("\nexpect = ", stderr); mp_writefile(x, stderr, 16);
+ fputs("\nresult = ", stderr); mp_writefile(xx, stderr, 16);
+ fputc('\n', stderr);
+ ok = 0;
+ }
+ if (!MP_EQ(y, yy)) {
+ fputs("\n*** mp_gcd(y) failed", stderr);
+ fputs("\na = ", stderr); mp_writefile(a, stderr, 16);
+ fputs("\nb = ", stderr); mp_writefile(b, stderr, 16);
+ fputs("\nexpect = ", stderr); mp_writefile(y, stderr, 16);
+ fputs("\nresult = ", stderr); mp_writefile(yy, stderr, 16);
+ fputc('\n', stderr);
+ ok = 0;
+ }
+
+ if (!ok) {
+ mp *ax = gf_mul(MP_NEW, a, xx);
+ mp *by = gf_mul(MP_NEW, b, yy);
+ ax = gf_add(ax, ax, by);
+ if (MP_EQ(ax, gg))
+ fputs("\n*** (Alternative result found.)\n", stderr);
+ MP_DROP(ax);
+ MP_DROP(by);
+ }
+
+ if (!MP_EQ(g, gg)) {
+ fputs("\n*** mp_gcd(gcd) failed", stderr);
+ fputs("\na = ", stderr); mp_writefile(a, stderr, 16);
+ fputs("\nb = ", stderr); mp_writefile(b, stderr, 16);
+ fputs("\nexpect = ", stderr); mp_writefile(g, stderr, 16);
+ fputs("\nresult = ", stderr); mp_writefile(gg, stderr, 16);
+ fputc('\n', stderr);
+ ok = 0;
+ }
+ MP_DROP(a); MP_DROP(b); MP_DROP(g); MP_DROP(x); MP_DROP(y);
+ MP_DROP(gg); MP_DROP(xx); MP_DROP(yy);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
+ return (ok);
+ }
+
+ static test_chunk tests[] = {
+ { "gcd", gcd, { &type_mp, &type_mp, &type_mp, &type_mp, &type_mp, 0 } },
+ { 0, 0, { 0 } }
+ };
+
+ int main(int argc, char *argv[])
+ {
+ sub_init();
+ test_run(argc, argv, tests, SRCDIR "/tests/gf");
+ return (0);
+ }
+
+ #endif
+
+ /*----- That's all, folks -------------------------------------------------*/
--- /dev/null
- * $Id: gf.h,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
+ /* -*-c-*-
+ *
++ * $Id: gf.h,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Arithmetic on binary polynomials
+ *
+ * (c) 2004 Straylight/Edgeware
+ */
+
+ /*----- Licensing notice --------------------------------------------------*
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ *
+ * Catacomb 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 Library General Public License for more details.
+ *
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+ /*----- Revision history --------------------------------------------------*
+ *
+ * $Log: gf.h,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.2.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ */
+
+ #ifndef CATACOMB_GF_H
+ #define CATACOMB_GF_H
+
+ #ifdef __cplusplus
+ extern "C" {
+ #endif
+
+ /*----- Header files ------------------------------------------------------*/
+
+ #ifndef CATACOMB_MP_H
+ # include "mp.h"
+ #endif
+
+ #ifndef CATACOMB_GFX_H
+ # include "gfx.h"
+ #endif
+
+ /*----- Functions provided ------------------------------------------------*/
+
+ /* --- @gf_add@ --- *
+ *
+ * Arguments: @mp *d@ = destination
+ * @mp *a, *b@ = sources
+ *
+ * Returns: Result, @a@ added to @b@.
+ */
+
+ extern mp *gf_add(mp */*d*/, mp */*a*/, mp */*b*/);
+ #define gf_sub gf_add
+
+ /* --- @gf_mul@ --- *
+ *
+ * Arguments: @mp *d@ = destination
+ * @mp *a, *b@ = sources
+ *
+ * Returns: Result, @a@ multiplied by @b@.
+ */
+
+ extern mp *gf_mul(mp */*d*/, mp */*a*/, mp */*b*/);
+
+ /* --- @gf_sqr@ --- *
+ *
+ * Arguments: @mp *d@ = destination
+ * @mp *a@ = source
+ *
+ * Returns: Result, @a@ squared.
+ */
+
+ extern mp *gf_sqr(mp */*d*/, mp */*a*/);
+
+ /* --- @gf_div@ --- *
+ *
+ * Arguments: @mp **qq, **rr@ = destination, quotient and remainder
+ * @mp *a, *b@ = sources
+ *
+ * Use: Calculates the quotient and remainder when @a@ is divided by
+ * @b@. The destinations @*qq@ and @*rr@ must be distinct.
+ * Either of @qq@ or @rr@ may be null to indicate that the
+ * result is irrelevant. (Discarding both results is silly.)
+ * There is a performance advantage if @a == *rr@.
+ */
+
+ extern void gf_div(mp **/*qq*/, mp **/*rr*/, mp */*a*/, mp */*b*/);
+
+ /* --- @gf_gcd@ --- *
+ *
+ * Arguments: @mp **gcd, **xx, **yy@ = where to write the results
+ * @mp *a, *b@ = sources (must be nonzero)
+ *
+ *
+ * Returns: ---
+ *
+ * Use: Calculates @gcd(a, b)@, and two numbers @x@ and @y@ such that
+ * @ax + by = gcd(a, b)@. This is useful for computing modular
+ * inverses.
+ */
+
+ extern void gf_gcd(mp **/*gcd*/, mp **/*xx*/, mp **/*yy*/,
+ mp */*a*/, mp */*b*/);
+
+ /*----- That's all, folks -------------------------------------------------*/
+
+ #ifdef __cplusplus
+ }
+ #endif
+
+ #endif
--- /dev/null
- * $Id: gfreduce-exp.h,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
+ /* -*-c-*-
+ *
++ * $Id: gfreduce-exp.h,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Exponentiation operations for binary field reduction
+ *
+ * (c) 2004 Straylight/Edgeware
+ */
+
+ /*----- Licensing notice --------------------------------------------------*
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ *
+ * Catacomb 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 Library General Public License for more details.
+ *
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+ /*----- Revision history --------------------------------------------------*
+ *
+ * $Log: gfreduce-exp.h,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.2.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ */
+
+ #ifndef CATACOMB_GFREDUCE_EXP_H
+ #define CATACOMB_GFREDUCE_EXP_H
+
+ #ifdef __cplusplus
+ extern "C" {
+ #endif
+
+ /*----- Exponentiation definitions ----------------------------------------*/
+
+ #define EXP_TYPE mp *
+
+ #define EXP_COPY(d, x) d = MP_COPY(x)
+ #define EXP_DROP(x) MP_DROP(x)
+
+ #define EXP_MUL(a, x) do { \
+ mp *t = gf_mul(spare, a, x); \
+ spare = a; \
+ a = gfreduce_do(gr, t, t); \
+ } while (0)
+
+ #define EXP_SQR(a) do { \
+ mp *t = gf_sqr(spare, a); \
+ spare = a; \
+ a = gfreduce_do(gr, t, t); \
+ } while (0)
+
+ #define EXP_FIX(x)
+
+ #define EXP_SETMUL(d, x, y) do { \
+ d = gf_mul(MP_NEW, x, y); \
+ d = gfreduce_do(gr, d, d); \
+ } while (0)
+
+ #define EXP_SETSQR(d, x) do { \
+ d = gf_sqr(MP_NEW, x); \
+ d = gfreduce_do(gr, d, d); \
+ } while (0)
+
+ #include "exp.h"
+
+ /*----- That's all, folks -------------------------------------------------*/
+
+ #ifdef __cplusplus
+ }
+ #endif
+
+ #endif
--- /dev/null
- * $Id: gfreduce.c,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
+ /* -*-c-*-
+ *
++ * $Id: gfreduce.c,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Efficient reduction modulo sparse binary polynomials
+ *
+ * (c) 2004 Straylight/Edgeware
+ */
+
+ /*----- Licensing notice --------------------------------------------------*
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ *
+ * Catacomb 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 Library General Public License for more details.
+ *
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+ /*----- Revision history --------------------------------------------------*
+ *
+ * $Log: gfreduce.c,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.2.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ */
+
+ /*----- Header files ------------------------------------------------------*/
+
+ #include <mLib/alloc.h>
+ #include <mLib/darray.h>
+ #include <mLib/macros.h>
+
+ #include "gf.h"
+ #include "gfreduce.h"
+ #include "gfreduce-exp.h"
+ #include "fibrand.h"
+ #include "mprand.h"
+
+ /*----- Data structures ---------------------------------------------------*/
+
+ DA_DECL(instr_v, gfreduce_instr);
+
+ /*----- Main code ---------------------------------------------------------*/
+
+ /* --- What's going on here? --- *
+ *
+ * Let's face it, @gfx_div@ sucks. It works (I hope), but it's not in any
+ * sense fast. Here, we do efficient reduction modulo sparse polynomials.
+ *
+ * Suppose we have a polynomial @X@ we're trying to reduce mod @P@. If we
+ * take the topmost nonzero word of @X@, call it @w@, then we can eliminate
+ * it by subtracting off @w P x^{k}@ for an appropriate value of @k@. The
+ * trick is in observing that if @P@ is sparse we can do this multiplication
+ * and subtraction efficiently, just by XORing appropriate shifts of @w@ into
+ * @X@.
+ *
+ * The first tricky bit is in working out when to stop. I'll use eight-bit
+ * words to demonstrate what I'm talking about.
+ *
+ * xxxxxx xxxxxxxx xxxxxxxx xxxxxxxx xxxxxxxx xxxxxxxx
+ * 001ppppp pppppppp pppppppp pppppppp
+ * |<rp>|
+ * |<------------ bp ------------->|
+ * |<------------ nw --------------->|
+ *
+ * The trick of taking whole words off of @X@ stops working when there are
+ * only @nw@ words left. Then we have to mask off the bottom bits of @w@
+ * before continuing.
+ */
+
+ /* --- @gfreduce_create@ --- *
+ *
+ * Arguments: @gfreduce *r@ = structure to fill in
+ * @mp *x@ = a (hopefully sparse) polynomial
+ *
+ * Returns: ---
+ *
+ * Use: Initializes a context structure for reduction.
+ */
+
+ void gfreduce_create(gfreduce *r, mp *p)
+ {
+ instr_v iv = DA_INIT;
+ unsigned long d, dw;
+ mpscan sc;
+ unsigned long i;
+ gfreduce_instr *ip;
+ unsigned f = 0;
+ size_t w, ww, wi, wl, ll;
+
+ /* --- Sort out the easy stuff --- */
+
+ d = mp_bits(p); assert(d); d--;
+ r->lim = d/MPW_BITS;
+ dw = d%MPW_BITS;
+ if (!dw)
+ r->mask = 0;
+ else {
+ r->mask = MPW(((mpw)-1) << dw);
+ r->lim++;
+ }
+ r->p = mp_copy(p);
+
+ /* --- Stash a new instruction --- */
+
+ #define INSTR(op_, arg_) do { \
+ DA_ENSURE(&iv, 1); \
+ DA(&iv)[DA_LEN(&iv)].op = (op_); \
+ DA(&iv)[DA_LEN(&iv)].arg = (arg_); \
+ DA_EXTEND(&iv, 1); \
+ } while (0)
+
+ #define f_lsr 1u
+
+ w = (d + MPW_BITS - 1)/MPW_BITS;
+ INSTR(GFRI_LOAD, w);
+ wi = DA_LEN(&iv);
+ f = 0;
+ ll = 0;
+ for (i = 0, mp_scan(&sc, p); mp_step(&sc) && i < d; i++) {
+ if (!mp_bit(&sc))
+ continue;
+ ww = (d - i + MPW_BITS - 1)/MPW_BITS;
+ if (ww != w) {
+ wl = DA_LEN(&iv);
+ INSTR(GFRI_STORE, w);
+ if (!ll)
+ ll = DA_LEN(&iv);
+ if (!(f & f_lsr))
+ INSTR(GFRI_LOAD, ww);
+ else {
+ INSTR(GFRI_LOAD, w - 1);
+ for (; wi < wl; wi++) {
+ ip = &DA(&iv)[wi];
+ assert(ip->op == GFRI_LSL);
+ if (ip->arg)
+ INSTR(GFRI_LSR, MPW_BITS - ip->arg);
+ }
+ if (w - 1 != ww) {
+ INSTR(GFRI_STORE, w - 1);
+ INSTR(GFRI_LOAD, ww);
+ }
+ f &= ~f_lsr;
+ }
+ w = ww;
+ wi = DA_LEN(&iv);
+ }
+ INSTR(GFRI_LSL, (i - d)%MPW_BITS);
+ if ((i - d)%MPW_BITS)
+ f |= f_lsr;
+ }
+ wl = DA_LEN(&iv);
+ INSTR(GFRI_STORE, w);
+ if (!ll)
+ ll = DA_LEN(&iv);
+ if (f & f_lsr) {
+ INSTR(GFRI_LOAD, w - 1);
+ for (; wi < wl; wi++) {
+ ip = &DA(&iv)[wi];
+ assert(ip->op == GFRI_LSL);
+ if (ip->arg)
+ INSTR(GFRI_LSR, MPW_BITS - ip->arg);
+ }
+ INSTR(GFRI_STORE, w - 1);
+ }
+
+ #undef INSTR
+
+ r->in = DA_LEN(&iv);
+ r->iv = xmalloc(r->in * sizeof(gfreduce_instr));
+ r->liv = r->iv + ll;
+ memcpy(r->iv, DA(&iv), r->in * sizeof(gfreduce_instr));
+ DA_DESTROY(&iv);
+ }
+
+ /* --- @gfreduce_destroy@ --- *
+ *
+ * Arguments: @gfreduce *r@ = structure to free
+ *
+ * Returns: ---
+ *
+ * Use: Reclaims the resources from a reduction context.
+ */
+
+ void gfreduce_destroy(gfreduce *r)
+ {
+ mp_drop(r->p);
+ xfree(r->iv);
+ }
+
+ /* --- @gfreduce_dump@ --- *
+ *
+ * Arguments: @gfreduce *r@ = structure to dump
+ * @FILE *fp@ = file to dump on
+ *
+ * Returns: ---
+ *
+ * Use: Dumps a reduction context.
+ */
+
+ void gfreduce_dump(gfreduce *r, FILE *fp)
+ {
+ size_t i;
+
+ fprintf(fp, "poly = "); mp_writefile(r->p, fp, 16);
+ fprintf(fp, "\n lim = %lu; mask = %lx\n",
+ (unsigned long)r->lim, (unsigned long)r->mask);
+ for (i = 0; i < r->in; i++) {
+ static const char *opname[] = { "load", "lsl", "lsr", "store" };
+ assert(r->iv[i].op < N(opname));
+ fprintf(fp, " %s %lu\n",
+ opname[r->iv[i].op],
+ (unsigned long)r->iv[i].arg);
+ }
+ }
+
+ /* --- @gfreduce_do@ --- *
+ *
+ * Arguments: @gfreduce *r@ = reduction context
+ * @mp *d@ = destination
+ * @mp *x@ = source
+ *
+ * Returns: Destination, @x@ reduced modulo the reduction poly.
+ */
+
+ static void run(const gfreduce_instr *i, const gfreduce_instr *il,
+ mpw *v, mpw z)
+ {
+ mpw w = 0;
+
+ for (; i < il; i++) {
+ switch (i->op) {
+ case GFRI_LOAD: w = *(v - i->arg); break;
+ case GFRI_LSL: w ^= z << i->arg; break;
+ case GFRI_LSR: w ^= z >> i->arg; break;
+ case GFRI_STORE: *(v - i->arg) = MPW(w); break;
+ default: abort();
+ }
+ }
+ }
+
+ mp *gfreduce_do(gfreduce *r, mp *d, mp *x)
+ {
+ mpw *v, *vl;
+ const gfreduce_instr *il;
+ mpw z;
+
+ /* --- Try to reuse the source's space --- */
+
+ MP_COPY(x);
+ if (d) MP_DROP(d);
+ MP_DEST(x, MP_LEN(x), x->f);
+
+ /* --- Do the reduction --- */
+
+ il = r->iv + r->in;
+ if (MP_LEN(x) >= r->lim) {
+ v = x->v + r->lim;
+ vl = x->vl;
+ while (vl-- > v) {
+ while (*vl) {
+ z = *vl;
+ *vl = 0;
+ run(r->iv, il, vl, z);
+ }
+ }
+ if (r->mask) {
+ while (*vl & r->mask) {
+ z = *vl & r->mask;
+ *vl &= ~r->mask;
+ run(r->iv, il, vl, z);
+ }
+ }
+ }
+
+ /* --- Done --- */
+
+ MP_SHRINK(x);
+ return (x);
+ }
+
+ /* --- @gfreduce_sqrt@ --- *
+ *
+ * Arguments: @gfreduce *r@ = pointer to reduction context
+ * @mp *d@ = destination
+ * @mp *x@ = some polynomial
+ *
+ * Returns: The square root of @x@ modulo @r->p@, or null.
+ */
+
+ mp *gfreduce_sqrt(gfreduce *r, mp *d, mp *x)
+ {
+ mp *y = MP_COPY(x);
+ mp *z, *spare = MP_NEW;
+ unsigned long m = mp_bits(r->p) - 1;
+ unsigned long i;
+
+ for (i = 0; i < m - 1; i++) {
+ mp *t = gf_sqr(spare, y);
+ spare = y;
+ y = gfreduce_do(r, t, t);
+ }
+ z = gf_sqr(spare, y);
+ z = gfreduce_do(r, z, z);
+ if (!MP_EQ(x, z)) {
+ mp_drop(y);
+ y = 0;
+ }
+ mp_drop(z);
+ mp_drop(d);
+ return (y);
+ }
+
+ /* --- @gfreduce_trace@ --- *
+ *
+ * Arguments: @gfreduce *r@ = pointer to reduction context
+ * @mp *x@ = some polynomial
+ *
+ * Returns: The trace of @x@. (%$\Tr(x)=x + x^2 + \cdots + x^{2^{m-1}}$%
+ * if %$x \in \gf{2^m}$%).
+ */
+
+ int gfreduce_trace(gfreduce *r, mp *x)
+ {
+ mp *y = MP_COPY(x);
+ mp *spare = MP_NEW;
+ unsigned long m = mp_bits(r->p) - 1;
+ unsigned long i;
+ int rc;
+
+ for (i = 0; i < m - 1; i++) {
+ mp *t = gf_sqr(spare, y);
+ spare = y;
+ y = gfreduce_do(r, t, t);
+ y = gf_add(y, y, x);
+ }
+ rc = !MP_ISZERO(y);
+ mp_drop(spare);
+ mp_drop(y);
+ return (rc);
+ }
+
+ /* --- @gfreduce_halftrace@ --- *
+ *
+ * Arguments: @gfreduce *r@ = pointer to reduction context
+ * @mp *d@ = destination
+ * @mp *x@ = some polynomial
+ *
+ * Returns: The half-trace of @x@.
+ * (%$\HfTr(x)= x + x^{2^2} + \cdots + x^{2^{m-1}}$%
+ * if %$x \in \gf{2^m}$% with %$m$% odd).
+ */
+
+ mp *gfreduce_halftrace(gfreduce *r, mp *d, mp *x)
+ {
+ mp *y = MP_COPY(x);
+ mp *spare = MP_NEW;
+ unsigned long m = mp_bits(r->p) - 1;
+ unsigned long i;
+
+ mp_drop(d);
+ for (i = 0; i < m - 1; i += 2) {
+ mp *t = gf_sqr(spare, y);
+ spare = y;
+ y = gfreduce_do(r, t, t);
+ t = gf_sqr(spare, y);
+ spare = y;
+ y = gfreduce_do(r, t, t);
+ y = gf_add(y, y, x);
+ }
+ mp_drop(spare);
+ return (y);
+ }
+
+ /* --- @gfreduce_quadsolve@ --- *
+ *
+ * Arguments: @gfreduce *r@ = pointer to reduction context
+ * @mp *d@ = destination
+ * @mp *x@ = some polynomial
+ *
+ * Returns: A polynomial @y@ such that %$y^2 + y = x$%, or null.
+ */
+
+ mp *gfreduce_quadsolve(gfreduce *r, mp *d, mp *x)
+ {
+ unsigned long m = mp_bits(r->p) - 1;
+ mp *t;
+
+ MP_COPY(x);
+ if (m & 1)
+ d = gfreduce_halftrace(r, d, x);
+ else {
+ mp *z, *w, *rho = MP_NEW;
+ mp *spare = MP_NEW;
+ grand *fr = fibrand_create(0);
+ unsigned long i;
+
+ for (;;) {
+ rho = mprand(rho, m, fr, 0);
+ z = MP_ZERO;
+ w = MP_COPY(rho);
+ for (i = 0; i < m - 1; i++) {
+ t = gf_sqr(spare, z); spare = z; z = gfreduce_do(r, t, t);
+ t = gf_sqr(spare, w); spare = w; w = gfreduce_do(r, t, t);
+ t = gf_mul(spare, w, x); t = gfreduce_do(r, t, t); spare = t;
+ z = gf_add(z, z, t);
+ w = gf_add(w, w, rho);
+ }
+ if (!MP_ISZERO(w))
+ break;
+ MP_DROP(z);
+ MP_DROP(w);
+ }
+ if (d) MP_DROP(d);
+ MP_DROP(w);
+ MP_DROP(spare);
+ MP_DROP(rho);
+ fr->ops->destroy(fr);
+ d = z;
+ }
+
+ t = gf_sqr(MP_NEW, d); t = gfreduce_do(r, t, t); t = gf_add(t, t, d);
+ if (!MP_EQ(t, x)) {
+ MP_DROP(d);
+ d = 0;
+ }
+ MP_DROP(t);
+ MP_DROP(x);
+ d->v[0] &= ~(mpw)1;
+ return (d);
+ }
+
+ /* --- @gfreduce_exp@ --- *
+ *
+ * Arguments: @gfreduce *gr@ = pointer to reduction context
+ * @mp *d@ = fake destination
+ * @mp *a@ = base
+ * @mp *e@ = exponent
+ *
+ * Returns: Result, %$a^e \bmod m$%.
+ */
+
+ mp *gfreduce_exp(gfreduce *gr, mp *d, mp *a, mp *e)
+ {
+ mp *x = MP_ONE;
+ mp *spare = (e->f & MP_BURN) ? MP_NEWSEC : MP_NEW;
+
+ MP_SHRINK(e);
+ if (!MP_LEN(e))
+ ;
+ else if (MP_LEN(e) < EXP_THRESH)
+ EXP_SIMPLE(x, a, e);
+ else
+ EXP_WINDOW(x, a, e);
+ mp_drop(d);
+ mp_drop(spare);
+ return (x);
+ }
+
+ /*----- Test rig ----------------------------------------------------------*/
+
+ #ifdef TEST_RIG
+
+ #define MP(x) mp_readstring(MP_NEW, #x, 0, 0)
+
+ static int vreduce(dstr *v)
+ {
+ mp *d = *(mp **)v[0].buf;
+ mp *n = *(mp **)v[1].buf;
+ mp *r = *(mp **)v[2].buf;
+ mp *c;
+ int ok = 1;
+ gfreduce rr;
+
+ gfreduce_create(&rr, d);
+ c = gfreduce_do(&rr, MP_NEW, n);
+ if (!MP_EQ(c, r)) {
+ fprintf(stderr, "\n*** reduction failed\n*** ");
+ gfreduce_dump(&rr, stderr);
+ fprintf(stderr, "\n*** n = "); mp_writefile(n, stderr, 16);
+ fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16);
+ fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16);
+ fprintf(stderr, "\n");
+ ok = 0;
+ }
+ gfreduce_destroy(&rr);
+ mp_drop(n); mp_drop(d); mp_drop(r); mp_drop(c);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
+ return (ok);
+ }
+
+ static int vmodexp(dstr *v)
+ {
+ mp *p = *(mp **)v[0].buf;
+ mp *g = *(mp **)v[1].buf;
+ mp *x = *(mp **)v[2].buf;
+ mp *r = *(mp **)v[3].buf;
+ mp *c;
+ int ok = 1;
+ gfreduce rr;
+
+ gfreduce_create(&rr, p);
+ c = gfreduce_exp(&rr, MP_NEW, g, x);
+ if (!MP_EQ(c, r)) {
+ fprintf(stderr, "\n*** modexp failed\n*** ");
+ fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16);
+ fprintf(stderr, "\n*** g = "); mp_writefile(g, stderr, 16);
+ fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16);
+ fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16);
+ fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16);
+ fprintf(stderr, "\n");
+ ok = 0;
+ }
+ gfreduce_destroy(&rr);
+ mp_drop(p); mp_drop(g); mp_drop(r); mp_drop(x); mp_drop(c);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
+ return (ok);
+ }
+
+ static int vsqrt(dstr *v)
+ {
+ mp *p = *(mp **)v[0].buf;
+ mp *x = *(mp **)v[1].buf;
+ mp *r = *(mp **)v[2].buf;
+ mp *c;
+ int ok = 1;
+ gfreduce rr;
+
+ gfreduce_create(&rr, p);
+ c = gfreduce_sqrt(&rr, MP_NEW, x);
+ if (!MP_EQ(c, r)) {
+ fprintf(stderr, "\n*** sqrt failed\n*** ");
+ fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16);
+ fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16);
+ fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16);
+ fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16);
+ fprintf(stderr, "\n");
+ ok = 0;
+ }
+ gfreduce_destroy(&rr);
+ mp_drop(p); mp_drop(r); mp_drop(x); mp_drop(c);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
+ return (ok);
+ }
+
+ static int vtr(dstr *v)
+ {
+ mp *p = *(mp **)v[0].buf;
+ mp *x = *(mp **)v[1].buf;
+ int r = *(int *)v[2].buf, c;
+ int ok = 1;
+ gfreduce rr;
+
+ gfreduce_create(&rr, p);
+ c = gfreduce_trace(&rr, x);
+ if (c != r) {
+ fprintf(stderr, "\n*** trace failed\n*** ");
+ fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16);
+ fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16);
+ fprintf(stderr, "\n*** c = %d", c);
+ fprintf(stderr, "\n*** r = %d", r);
+ fprintf(stderr, "\n");
+ ok = 0;
+ }
+ gfreduce_destroy(&rr);
+ mp_drop(p); mp_drop(x);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
+ return (ok);
+ }
+
+ static int vhftr(dstr *v)
+ {
+ mp *p = *(mp **)v[0].buf;
+ mp *x = *(mp **)v[1].buf;
+ mp *r = *(mp **)v[2].buf;
+ mp *c;
+ int ok = 1;
+ gfreduce rr;
+
+ gfreduce_create(&rr, p);
+ c = gfreduce_halftrace(&rr, MP_NEW, x);
+ if (!MP_EQ(c, r)) {
+ fprintf(stderr, "\n*** halftrace failed\n*** ");
+ fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16);
+ fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16);
+ fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16);
+ fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16);
+ fprintf(stderr, "\n");
+ ok = 0;
+ }
+ gfreduce_destroy(&rr);
+ mp_drop(p); mp_drop(r); mp_drop(x); mp_drop(c);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
+ return (ok);
+ }
+
+ static int vquad(dstr *v)
+ {
+ mp *p = *(mp **)v[0].buf;
+ mp *x = *(mp **)v[1].buf;
+ mp *r = *(mp **)v[2].buf;
+ mp *c;
+ int ok = 1;
+ gfreduce rr;
+
+ gfreduce_create(&rr, p);
+ c = gfreduce_quadsolve(&rr, MP_NEW, x);
+ if (!MP_EQ(c, r)) {
+ fprintf(stderr, "\n*** quadsolve failed\n*** ");
+ fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16);
+ fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16);
+ fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16);
+ fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16);
+ fprintf(stderr, "\n");
+ ok = 0;
+ }
+ gfreduce_destroy(&rr);
+ mp_drop(p); mp_drop(r); mp_drop(x); mp_drop(c);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
+ return (ok);
+ }
+
+ static test_chunk defs[] = {
+ { "reduce", vreduce, { &type_mp, &type_mp, &type_mp, 0 } },
+ { "modexp", vmodexp, { &type_mp, &type_mp, &type_mp, &type_mp, 0 } },
+ { "sqrt", vsqrt, { &type_mp, &type_mp, &type_mp, 0 } },
+ { "trace", vtr, { &type_mp, &type_mp, &type_int, 0 } },
+ { "halftrace", vhftr, { &type_mp, &type_mp, &type_mp, 0 } },
+ { "quadsolve", vquad, { &type_mp, &type_mp, &type_mp, 0 } },
+ { 0, 0, { 0 } }
+ };
+
+ int main(int argc, char *argv[])
+ {
+ test_run(argc, argv, defs, SRCDIR"/tests/gfreduce");
+ return (0);
+ }
+
+ #endif
+
+ /*----- That's all, folks -------------------------------------------------*/
--- /dev/null
- * $Id: gfreduce.h,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
+ /* -*-c-*-
+ *
++ * $Id: gfreduce.h,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ *
+ * Reduction modulo sparse binary polynomials
+ *
+ * (c) 2004 Straylight/Edgeware
+ */
+
+ /*----- Licensing notice --------------------------------------------------*
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ *
+ * Catacomb 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 Library General Public License for more details.
+ *
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+ /*----- Revision history --------------------------------------------------*
+ *
+ * $Log: gfreduce.h,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.2.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
+ */
+
+ #ifndef CATACOMB_GFREDUCE_H
+ #define CATACOMB_GFREDUCE_H
+
+ #ifdef __cplusplus
+ extern "C" {
+ #endif
+
+ /*----- Header files ------------------------------------------------------*/
+
+ #ifndef CATACOMB_GF_H
+ # include "gf.h"
+ #endif
+
+ /*----- Data structures ---------------------------------------------------*/
+
+ typedef struct gfreduce_instr {
+ unsigned op; /* Instruction opcode */
+ size_t arg; /* Immediate argument */
+ } gfreduce_instr;
+
+ enum {
+ GFRI_LOAD, /* Load @p[arg]@ */
+ GFRI_LSL, /* XOR with @w << arg@ */
+ GFRI_LSR, /* XOR with @w >> arg@ */
+ GFRI_STORE, /* Store @p[arg]@ */
+ GFRI_MAX
+ };
+
+ typedef struct gfreduce {
+ size_t lim; /* Word of degree bit */
+ mpw mask; /* Mask for degree word */
+ mp *p; /* Copy of the polynomial */
+ size_t in; /* Number of instruction words */
+ gfreduce_instr *iv, *liv; /* Vector of instructions */
+ } gfreduce;
+
+ /*----- Functions provided ------------------------------------------------*/
+
+ /* --- @gfreduce_create@ --- *
+ *
+ * Arguments: @gfreduce *r@ = structure to fill in
+ * @mp *x@ = a (hopefully sparse) polynomial
+ *
+ * Returns: ---
+ *
+ * Use: Initializes a context structure for reduction.
+ */
+
+ extern void gfreduce_create(gfreduce */*r*/, mp */*p*/);
+
+ /* --- @gfreduce_destroy@ --- *
+ *
+ * Arguments: @gfreduce *r@ = structure to free
+ *
+ * Returns: ---
+ *
+ * Use: Reclaims the resources from a reduction context.
+ */
+
+ extern void gfreduce_destroy(gfreduce */*r*/);
+
+ /* --- @gfreduce_dump@ --- *
+ *
+ * Arguments: @gfreduce *r@ = structure to dump
+ * @FILE *fp@ = file to dump on
+ *
+ * Returns: ---
+ *
+ * Use: Dumps a reduction context.
+ */
+
+ extern void gfreduce_dump(gfreduce */*r*/, FILE */*fp*/);
+
+ /* --- @gfreduce_do@ --- *
+ *
+ * Arguments: @gfreduce *r@ = reduction context
+ * @mp *d@ = destination
+ * @mp *x@ = source
+ *
+ * Returns: Destination, @x@ reduced modulo the reduction poly.
+ */
+
+ extern mp *gfreduce_do(gfreduce */*r*/, mp */*d*/, mp */*x*/);
+
+ /* --- @gfreduce_sqrt@ --- *
+ *
+ * Arguments: @gfreduce *r@ = pointer to reduction context
+ * @mp *d@ = destination
+ * @mp *x@ = some polynomial
+ *
+ * Returns: The square root of @x@ modulo @r->p@, or null.
+ */
+
+ extern mp *gfreduce_sqrt(gfreduce */*r*/, mp */*d*/, mp */*x*/);
+
+ /* --- @gfreduce_trace@ --- *
+ *
+ * Arguments: @gfreduce *r@ = pointer to reduction context
+ * @mp *x@ = some polynomial
+ *
+ * Returns: The trace of @x@. (%$\Tr(x)=x + x^2 + \cdots + x^{2^{m-1}}$%
+ * if %$x \in \gf{2^m}$%).
+ */
+
+ extern int gfreduce_trace(gfreduce */*r*/, mp */*x*/);
+
+ /* --- @gfreduce_halftrace@ --- *
+ *
+ * Arguments: @gfreduce *r@ = pointer to reduction context
+ * @mp *d@ = destination
+ * @mp *x@ = some polynomial
+ *
+ * Returns: The half-trace of @x@.
+ * (%$\HfTr(x)= x + x^{2^2} + \cdots + x^{2^{m-1}}$%
+ * if %$x \in \gf{2^m}$% with %$m$% odd).
+ */
+
+ extern mp *gfreduce_halftrace(gfreduce */*r*/, mp */*d*/, mp */*x*/);
+
+ /* --- @gfreduce_quadsolve@ --- *
+ *
+ * Arguments: @gfreduce *r@ = pointer to reduction context
+ * @mp *d@ = destination
+ * @mp *x@ = some polynomial
+ *
+ * Returns: A polynomial @y@ such that %$y^2 + y = x$%, or null.
+ */
+
+ extern mp *gfreduce_quadsolve(gfreduce */*r*/, mp */*d*/, mp */*x*/);
+
+ /* --- @gfreduce_exp@ --- *
+ *
+ * Arguments: @gfreduce *gr@ = pointer to reduction context
+ * @mp *d@ = fake destination
+ * @mp *a@ = base
+ * @mp *e@ = exponent
+ *
+ * Returns: Result, %$a^e \bmod m$%.
+ */
+
+ extern mp *gfreduce_exp(gfreduce */*gr*/, mp */*d*/, mp */*a*/, mp */*e*/);
+
+ /*----- That's all, folks -------------------------------------------------*/
+
+ #ifdef __cplusplus
+ }
+ #endif
+
+ #endif
/* -*-c-*-
*
- * $Id: gfx-sqr.c,v 1.1 2000/10/08 15:49:37 mdw Exp $
- * $Id: gfx-sqr.c,v 1.1.4.1 2004/03/21 22:39:46 mdw Exp $
++ * $Id: gfx-sqr.c,v 1.2 2004/03/21 22:52:06 mdw Exp $
*
* Sqaring binary polynomials
*
/*----- Revision history --------------------------------------------------*
*
* $Log: gfx-sqr.c,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.4.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
* Revision 1.1 2000/10/08 15:49:37 mdw
* First glimmerings of binary polynomial arithmetic.
*
/* -*-c-*-
*
- * $Id: gfx.h,v 1.1 2000/10/08 15:49:37 mdw Exp $
- * $Id: gfx.h,v 1.1.4.1 2004/03/21 22:39:46 mdw Exp $
++ * $Id: gfx.h,v 1.2 2004/03/21 22:52:06 mdw Exp $
*
* Low-level arithmetic on binary polynomials
*
/*----- Revision history --------------------------------------------------*
*
* $Log: gfx.h,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.4.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
* Revision 1.1 2000/10/08 15:49:37 mdw
* First glimmerings of binary polynomial arithmetic.
*
/* -*-c-*-
*
- * $Id: mp-gcd.c,v 1.5 2000/10/08 12:02:41 mdw Exp $
- * $Id: mp-gcd.c,v 1.5.4.1 2004/03/21 22:39:46 mdw Exp $
++ * $Id: mp-gcd.c,v 1.6 2004/03/21 22:52:06 mdw Exp $
*
* Extended GCD calculation
*
/*----- Revision history --------------------------------------------------*
*
* $Log: mp-gcd.c,v $
++ * Revision 1.6 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.5.4.1 2004/03/21 22:39:46 mdw
+ * Elliptic curves on binary fields work.
+ *
* Revision 1.5 2000/10/08 12:02:41 mdw
* Use Euclid's algorithm rather than the binary one.
*
/* -*-c-*-
*
- * $Id: mpbarrett-exp.h,v 1.1 2001/06/16 12:58:12 mdw Exp $
- * $Id: mpbarrett-exp.h,v 1.1.4.1 2004/03/20 00:20:05 mdw Exp $
++ * $Id: mpbarrett-exp.h,v 1.2 2004/03/21 22:52:06 mdw Exp $
*
* Exponentiation operations for Barrett reduction
*
/*----- Revision history --------------------------------------------------*
*
* $Log: mpbarrett-exp.h,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.4.1 2004/03/20 00:20:05 mdw
+ * Projective coordinates for prime curves
+ *
* Revision 1.1 2001/06/16 12:58:12 mdw
* Parameters for generic exponentiation.
*
/* -*-c-*-
*
- * $Id: mpmont-exp.h,v 1.1 2001/06/16 12:58:12 mdw Exp $
- * $Id: mpmont-exp.h,v 1.1.4.1 2004/03/20 00:13:31 mdw Exp $
++ * $Id: mpmont-exp.h,v 1.2 2004/03/21 22:52:06 mdw Exp $
*
* Exponentiation operations for Montgomery reduction
*
/*----- Revision history --------------------------------------------------*
*
* $Log: mpmont-exp.h,v $
++ * Revision 1.2 2004/03/21 22:52:06 mdw
++ * Merge and close elliptic curve branch.
++ *
+ * Revision 1.1.4.1 2004/03/20 00:13:31 mdw
+ * Projective coordinates for prime curves
+ *
* Revision 1.1 2001/06/16 12:58:12 mdw
* Parameters for generic exponentiation.
*
--- /dev/null
-# $Id: gf,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
++# $Id: gf,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ #
+ # Test cases for higher-level binary poly arithmetic.
+
+ add {
+ 0 0 0;
+ 1 1 0;
+ 1 2 3;
+ 4 5 1;
+ 0x7fb838a8a0a95046b9d9d9fb4440f7bb
+ 0xc1a7bd3b4e853fc92d4e1588719986aa
+ 0xbe1f8593ee2c6f8f9497cc7335d97111;
+ 0x1e2933215e1c3bba8d2b404d98f43086bfc6198a219b168f214042a5e7df6b21
+ 0x1e2933215e1c3bba8d2b404d98f43086bfc6198a219b168f214042a5e7df6b22 3;
+ }
+
+ mul {
+ 0 0 0;
+ 1 0 0;
+ 0 1 0;
+ 1 1 1;
+ 0x7fb838a8a0a95046b9d9d9fb4440f7bb
+ 0xc1a7bd3b4e853fc92d4e1588719986aa
+ 0x207ccad257b4ed64447158315bfb9aca5cbc5622cfb8fcbb1380eea1bc5c624e;
+ 0xc1a7bd3b4e853fc92d4e1588719986aa
+ 0x283ed59f1226dcefa7ff0ef87ceff5d5
+ 0x1e2933215e1c3bba8d2b404d98f43086bfc6198a219b168f214042a5e7df6b22;
+ 0xbe1f8593ee2c6f8f9497cc7335d97111
+ 0x35a8e33503b3695be00528f8b82db931
+ 0x1e2933215e1c3bba8d2b404d98f43086bfc6198a219b168f214042a5e7df6b21;
+ }
+
+ sqr {
+ 0 0;
+ 1 1;
+ 3 5;
+ 0x7fb838a8a0a95046b9d9d9fb4440f7bb
+ 0x1555454005404440440044411100101445415141514155451010100055154545;
+ 0x01f081e69f45d3254530766ab98d55fa612c7bb27ea31bc2621d894be9c0b196b3
+ 0x0155004001541441551011510504111011050015141444454140511111554414010450154545041554440501455004140401514041104554415000450141144505;
+ }
+
+ div {
+ 0 1 0 0;
+ 0x207ccad257b4ed64447158315bfb9aca5cbc5622cfb8fcbb1380eea1bc5c624e
+ 0x7fb838a8a0a95046b9d9d9fb4440f7bb
+ 0xc1a7bd3b4e853fc92d4e1588719986aa 0;
+ 0x6e0e2a282a5411ae76767ed1103deef069ef4ed3a14ff24b
+ 0x5385621c6661aaa35a24150d2c08332e
+ 0x01c2334cc957151dc7
+ 0x398c4111da6d06cdf3d83704ee403101;
+ }
+
+ gcd {
+ 0xc1a7bd3b4e853fc92d4e1588719986aa
+ 0xbe1f8593ee2c6f8f9497cc7335d97111
+ 3
+ 0x283ed59f1226dcefa7ff0ef87ceff5d5
+ 0x35a8e33503b3695be00528f8b82db931;
+ 0xbe1f8593ee2c6f8f9497cc7335d97111
+ 0xc1a7bd3b4e853fc92d4e1588719986aa
+ 3
+ 0x35a8e33503b3695be00528f8b82db931
+ 0x283ed59f1226dcefa7ff0ef87ceff5d5;
+ 0x800000000000000000000000000000000000000c9
+ 0x3f0eba16286a2d57ea0991168d4994637e8343e36
+ 1
+ 0xa17e704470d80cb5a78f295db0ce543dda16a169
+ 0x3c8c172e24598e90b9542e6b8f6571f54be572b50;
+ }
--- /dev/null
-# $Id: gfreduce,v 1.1.2.1 2004/03/21 22:39:46 mdw Exp $
++# $Id: gfreduce,v 1.2 2004/03/21 22:52:06 mdw Exp $
+ #
+ # Test efficient polynomial reduction
+
+ reduce {
+ 0x10000000
+ 0x4509823098098435
+ 0x8098435;
+ 0x100000000000000050002
+ 0x4509823098098435
+ 0x4509823098098435;
+ 0x100000000000000050002
+ 0x450982309809843545609843098560803495
+ 0x144f98a2f5cbc4773cfd;
+ 0xb2ca471b0867d5fae2e4f27a2d2706da
+ 0xf254423fef93d5d7a76ecf22c656c1352c53257875945d33
+ 0x582f783fc210f72814780e69b0bd29ff;
+ }
+
+ modexp {
+ 0x20000000000000000000000000000000000000000000000000000000000001001
+ 0x02
+ 0x1ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff
+ 1;
+ 0x20000000000000000000000000000000000000000000000000000000000001001
+ 0x435932098459080438094509845
+ 0x1ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff
+ 1;
+ 0x10000000000000000000000000000000000000000003
+ 0x02
+ 0x0fffffffffffffffffffffffffffffffffffffffffff
+ 1;
+ 0x10000000000000000000000000000000000000000003
+ 0x34235950984598345900983409845690805680985
+ 0x0fffffffffffffffffffffffffffffffffffffffffff
+ 1;
+ }
+
+ sqrt {
+ 0x20000000000000000000000000000000000000000000000000000000000001001
+ 0x1f081e69f45d3254530766ab98d55fa612c7bb27ea31bc2621d894be9c0b196b3
+ 0x7fb838a8a0a95046b9d9d9fb4440f7bbc1a7bd3b4e853fc92d4e1588719986aa;
+ 0x10000000000000000000000000000000000000000003
+ 0x4594094509835690805698083560980459903450984
+ 0x820291881a244a02840a2f8ece3f23f88f38bf0b3a;
+ }
+
+ halftrace {
+ 0x20000000000000000000000000000000000000000000000000000000000001001
+ 0x174e65c7d14a8ec286df8c7df17662f13f1d3563f13c8c63f23f5d0bd5d1b45cd
+ 0x8d68905434b020ccb849e17a03a5c441d2a104aaf523699c1cc7a93174d21d9d;
+ }
+
+ quadsolve {
+ 0x20000000000000000000000000000000000000000000000000000000000001001
+ 0x174e65c7d14a8ec286df8c7df17662f13f1d3563f13c8c63f23f5d0bd5d1b45cd
+ 0x8d68905434b020ccb849e17a03a5c441d2a104aaf523699c1cc7a93174d21d9c;
+ 0x10000000000000000000000000000000000000000003
+ 0x3b818b447e90713da04f13c3b07cb5e2681d08e4700
+ 0x27aa17c97dfa80bbdef9f91b243c6e6ddba1a223cac;
+ }
# Test vectors for low-level GF functions
#
- # $Id: gfx,v 1.1 2000/10/08 16:01:48 mdw Exp $
-# $Id: gfx,v 1.1.4.1 2004/03/21 22:39:46 mdw Exp $
++# $Id: gfx,v 1.2 2004/03/21 22:52:06 mdw Exp $
# --- Addition (and subtraction) ---
projects / u / mdw / catacomb / commitdiff