Projective coordinates for prime curves

diff --git a/calc/ecp.cal b/calc/ecp.cal
index 82600f5..43ac1b3 100644 (file)
--- a/calc/ecp.cal
+++ b/calc/ecp.cal
@@ -1,6 +1,6 @@
 /* -*-apcalc-*-
  *
 /* -*-apcalc-*-
  *

- * $Id: ecp.cal,v 1.1.4.1 2003/06/10 13:43:53 mdw Exp $
+ * $Id: ecp.cal,v 1.1.4.2 2004/03/20 00:13:31 mdw Exp $

  *
  * Testbed for elliptic curve arithmetic over prime fields
  *
  *
  * Testbed for elliptic curve arithmetic over prime fields
  *


@@ -30,6 +30,9 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ecp.cal,v $
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ecp.cal,v $

+ * 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.4.1  2003/06/10 13:43:53  mdw
  * Simple (non-projective) curves over prime fields now seem to work.
  *


@@ -42,6 +45,7 @@
 
 obj ecp_curve { a, b, p };
 obj ecp_pt { x, y, e };
 
 obj ecp_curve { a, b, p };
 obj ecp_pt { x, y, e };

+obj ecpp_pt { x, y, z, e };

 
 /*----- Main code ---------------------------------------------------------*/
 
 
 /*----- Main code ---------------------------------------------------------*/
 


@@ -63,6 +67,72 @@ define ecp_pt(x, y, e)
   return (p);
 }
 
   return (p);
 }
 

+define ecpp_pt(p)
+{
+  local obj ecpp_pt pp;
+  if (istype(p, 1))
+    return (0);
+  pp.x = p.x;
+  pp.y = p.y;
+  pp.z = 1;
+  pp.e = p.e;
+  return (pp);
+}
+
+define ecpp_fix(pp)
+{
+  local obj ecp_pt p;
+  local e, zi, z2, z3;
+  if (istype(pp, 1) || pp.z == 0)
+    return (0);
+  e = pp.e;
+  zi = minv(pp.z, e.p);
+  z2 = zi * zi;
+  z3 = zi * z2;
+  p.x = pp.x * z2 % e.p;
+  p.y = pp.y * z3 % e.p;
+  p.e = e;
+  return (p);
+}
+
+define ecpp_dbl(a)
+{
+  local m, s, t, y2;
+  local e;
+  local obj ecpp_pt d;
+  if (istype(a, 1) || a.y == 0)
+    return (0);
+  e = a.e;
+  if (e.a % e.p == e.p - 3) {
+    m = a.z^3 % e.p;
+    m = 3 * (a.x + t4) * (a.x - t4) % e.p;
+  } else {
+    m = (3 * a.x^2 - e.a * a.z^4) % e.p;
+  }
+  d.z = 2 * a.y * a.z % e.p;
+  y2 = a.y^2 % e.p;
+  s = 4 * a.x * a.y % e.p;
+  d.x = (m^2 - 2 * s) % e.p;
+  d.y = (m * (s - d.x) - y * y2^2) % e.p;
+  d.e = e;
+  return (d);
+}
+
+define ecpp_add(a, b)
+{
+  if (a == 0)
+    d = b;
+  else if (b == 0)
+    d = a;
+  else if (!istype(a, b))
+    quit "bad type arguments to ecp_pt_add";
+  else if (a.e != b.e)
+    quit "points from different curves in ecp_pt_add";
+  else {
+    e = a.e;
+    
+}
+

 define ecp_pt_print(a)
 {
   print "(" : a.x : ", " : a.y : ")" :;
 define ecp_pt_print(a)
 {
   print "(" : a.x : ", " : a.y : ")" :;


@@ -103,6 +173,8 @@ define ecp_pt_dbl(a)
 {
   local e, alpha;
   local obj ecp_pt d;
 {
   local e, alpha;
   local obj ecp_pt d;

+  if (istype(a, 1))
+    return (0);

   e = a.e;
   alpha = (3 * a.x^2 + e.a) * minv(2 * a.y, e.p) % e.p;
   d.x = (alpha^2 - 2 * a.x) % e.p;
   e = a.e;
   alpha = (3 * a.x^2 + e.a) * minv(2 * a.y, e.p) % e.p;
   d.x = (alpha^2 - 2 * a.x) % e.p;

diff --git a/ec-exp.h b/ec-exp.h
index 9ba7bc5..bb4e08a 100644 (file)
--- a/ec-exp.h
+++ b/ec-exp.h
@@ -1,6 +1,6 @@
 /* -*-c-*-
  *
 /* -*-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 $

  *
  * Exponentiation operations for elliptic curves
  *
  *
  * Exponentiation operations for elliptic curves
  *


@@ -30,6 +30,9 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ec-exp.h,v $
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ec-exp.h,v $

+ * 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.
  *
  * Revision 1.2  2003/05/15 23:25:59  mdw
  * Make elliptic curve stuff build.
  *


@@ -58,6 +61,7 @@
 
 #define EXP_MUL(a, x) EC_ADD(c, &(a), &(a), &(x))
 #define EXP_SQR(a) EC_DBL(c, &(a), &(a));
 
 #define EXP_MUL(a, x) EC_ADD(c, &(a), &(a), &(x))
 #define EXP_SQR(a) EC_DBL(c, &(a), &(a));

+#define EXP_FIX(x) EC_FIX(c, &(x), &(x));

 
 #define EXP_SETMUL(d, x, y) do {                                       \
   EC_CREATE(&(d));                                                     \
 
 #define EXP_SETMUL(d, x, y) do {                                       \
   EC_CREATE(&(d));                                                     \

diff --git a/ec-prime.c b/ec-prime.c
index b2bfd52..14e4c16 100644 (file)
--- a/ec-prime.c
+++ b/ec-prime.c
@@ -1,6 +1,6 @@
 /* -*-c-*-
  *
 /* -*-c-*-
  *

- * $Id: ec-prime.c,v 1.3.4.1 2003/06/10 13:43:53 mdw Exp $
+ * $Id: ec-prime.c,v 1.3.4.2 2004/03/20 00:13:31 mdw Exp $

  *
  * Elliptic curves over prime fields
  *
  *
  * Elliptic curves over prime fields
  *


@@ -30,6 +30,9 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ec-prime.c,v $
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ec-prime.c,v $

+ * 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.4.1  2003/06/10 13:43:53  mdw
  * Simple (non-projective) curves over prime fields now seem to work.
  *


@@ -59,7 +62,7 @@ typedef struct ecctx {
 
 /*----- Simple prime curves -----------------------------------------------*/
 
 
 /*----- Simple prime curves -----------------------------------------------*/
 

-static const ec_ops ec_primeops;
+static const ec_ops ec_primeops, ec_primeprojops, ec_primeprojxops;

 
 static ec *ecneg(ec_curve *c, ec *d, const ec *p)
 {
 
 static ec *ecneg(ec_curve *c, ec *d, const ec *p)
 {


@@ -68,11 +71,33 @@ static ec *ecneg(ec_curve *c, ec *d, const ec *p)
   return (d);
 }
 
   return (d);
 }
 

+static ec *ecfind(ec_curve *c, ec *d, mp *x)
+{
+  mp *p, *q;
+  ecctx *cc = (ecctx *)c;
+  field *f = c->f;
+
+  q = F_SQR(f, MP_NEW, x);
+  p = F_MUL(f, MP_NEW, x, q);
+  q = F_MUL(f, q, x, cc->a);
+  p = F_ADD(f, p, p, q);
+  p = F_ADD(f, p, p, cc->b);
+  MP_DROP(q);
+  p = F_SQRT(f, p, p);
+  if (!p)
+    return (0);
+  EC_DESTROY(d);
+  d->x = MP_COPY(x);
+  d->y = p;
+  d->z = MP_COPY(f->one);
+  return (d);
+}
+

 static ec *ecdbl(ec_curve *c, ec *d, const ec *a)
 {
   if (EC_ATINF(a))
     EC_SETINF(d);
 static ec *ecdbl(ec_curve *c, ec *d, const ec *a)
 {
   if (EC_ATINF(a))
     EC_SETINF(d);

-  else if (!MP_LEN(a->y))
+  else if (F_ZEROP(c->f, a->y))

     EC_COPY(d, a);
   else {
     field *f = c->f;
     EC_COPY(d, a);
   else {
     field *f = c->f;


@@ -80,19 +105,19 @@ static ec *ecdbl(ec_curve *c, ec *d, const ec *a)
     mp *lambda;
     mp *dy, *dx;
 
     mp *lambda;
     mp *dy, *dx;
 

-    dx = F_SQR(f, MP_NEW, a->x);
-    dy = F_DBL(f, MP_NEW, a->y);
-    dx = F_TPL(f, dx, dx);
-    dx = F_ADD(f, dx, dx, cc->a);
-    dy = F_INV(f, dy, dy);
-    lambda = F_MUL(f, MP_NEW, dx, dy);
+    dx = F_SQR(f, MP_NEW, a->x);       /* %$x^2$% */
+    dy = F_DBL(f, MP_NEW, a->y);       /* %$2 y$% */
+    dx = F_TPL(f, dx, dx);             /* %$3 x^2$% */
+    dx = F_ADD(f, dx, dx, cc->a);      /* %$3 x^2 + A$% */
+    dy = F_INV(f, dy, dy);             /* %$(2 y)^{-1}$% */
+    lambda = F_MUL(f, MP_NEW, dx, dy); /* %$\lambda = (3 x^2 + A)/(2 y)$% */

 
 

-    dx = F_SQR(f, dx, lambda);
-    dy = F_DBL(f, dy, a->x);
-    dx = F_SUB(f, dx, dx, dy);
-    dy = F_SUB(f, dy, a->x, dx);
-    dy = F_MUL(f, dy, lambda, dy);
-    dy = F_SUB(f, dy, dy, a->y);
+    dx = F_SQR(f, dx, lambda);         /* %$\lambda^2$% */
+    dy = F_DBL(f, dy, a->x);           /* %$2 x$% */
+    dx = F_SUB(f, dx, dx, dy);         /* %$x' = \lambda^2 - 2 x */
+    dy = F_SUB(f, dy, a->x, dx);       /* %$x - x'$% */
+    dy = F_MUL(f, dy, lambda, dy);     /* %$\lambda (x - x')$% */
+    dy = F_SUB(f, dy, dy, a->y);       /* %$y' = \lambda (x - x') - y$% */

 
     EC_DESTROY(d);
     d->x = dx;
 
     EC_DESTROY(d);
     d->x = dx;


@@ -103,6 +128,94 @@ static ec *ecdbl(ec_curve *c, ec *d, const ec *a)
   return (d);
 }
 
   return (d);
 }
 

+static ec *ecprojdbl(ec_curve *c, ec *d, const ec *a)
+{
+  if (EC_ATINF(a))
+    EC_SETINF(d);
+  else if (F_ZEROP(c->f, a->y))
+    EC_COPY(d, a);
+  else {
+    field *f = c->f;
+    ecctx *cc = (ecctx *)c;
+    mp *p, *q, *m, *s, *dx, *dy, *dz;
+
+    p = F_SQR(f, MP_NEW, a->z);                /* %$z^2$% */
+    q = F_SQR(f, MP_NEW, p);           /* %$z^4$% */
+    p = F_MUL(f, p, q, cc->a);         /* %$A z^4$% */
+    m = F_SQR(f, MP_NEW, a->x);                /* %$x^2$% */
+    m = F_TPL(f, m, m);                        /* %$3 x^2$% */
+    m = F_ADD(f, m, m, p);             /* %$m = 3 x^2 + A z^4$% */
+
+    q = F_DBL(f, q, a->y);             /* %$2 y$% */
+    dz = F_MUL(f, MP_NEW, q, a->z);    /* %$z' = 2 y z$% */
+
+    p = F_SQR(f, p, q);                        /* %$4 y^2$% */
+    s = F_MUL(f, MP_NEW, p, a->x);     /* %$s = 4 x y^2$% */
+    q = F_SQR(f, q, p);                        /* %$16 y^4$% */
+    q = F_HLV(f, q, q);                        /* %$t = 8 y^4$% */
+
+    p = F_DBL(f, p, s);                        /* %$2 s$% */
+    dx = F_SQR(f, MP_NEW, m);          /* %$m^2$% */
+    dx = F_SUB(f, dx, dx, p);          /* %$x' = m^2 - 2 s$% */
+
+    s = F_SUB(f, s, s, dx);            /* %$s - x'$% */
+    dy = F_MUL(f, p, m, s);            /* %$m (s - x')$% */
+    dy = F_SUB(f, dy, dy, q);          /* %$y' = m (s - x') - t$% */
+
+    EC_DESTROY(d);
+    d->x = dx;
+    d->y = dy;
+    d->z = dz;
+    MP_DROP(m);
+    MP_DROP(q);
+    MP_DROP(s);
+  }
+  return (d);
+}
+
+static ec *ecprojxdbl(ec_curve *c, ec *d, const ec *a)
+{
+  if (EC_ATINF(a))
+    EC_SETINF(d);
+  else if (F_ZEROP(c->f, a->y))
+    EC_COPY(d, a);
+  else {
+    field *f = c->f;
+    mp *p, *q, *m, *s, *dx, *dy, *dz;
+
+    m = F_SQR(f, MP_NEW, a->z);                /* %$z^2$% */
+    p = F_SUB(f, MP_NEW, a->x, m);     /* %$x - z^2$% */
+    q = F_ADD(f, MP_NEW, a->x, m);     /* %$x + z^2$% */
+    m = F_MUL(f, m, p, q);             /* %$x^2 - z^4$% */
+    m = F_TPL(f, m, m);                        /* %$m = 3 x^2 - 3 z^4$% */
+
+    q = F_DBL(f, q, a->y);             /* %$2 y$% */
+    dz = F_MUL(f, MP_NEW, q, a->z);    /* %$z' = 2 y z$% */
+
+    p = F_SQR(f, p, q);                        /* %$4 y^2$% */
+    s = F_MUL(f, MP_NEW, p, a->x);     /* %$s = 4 x y^2$% */
+    q = F_SQR(f, q, p);                        /* %$16 y^4$% */
+    q = F_HLV(f, q, q);                        /* %$t = 8 y^4$% */
+
+    p = F_DBL(f, p, s);                        /* %$2 s$% */
+    dx = F_SQR(f, MP_NEW, m);          /* %$m^2$% */
+    dx = F_SUB(f, dx, dx, p);          /* %$x' = m^2 - 2 s$% */
+
+    s = F_SUB(f, s, s, dx);            /* %$s - x'$% */
+    dy = F_MUL(f, p, m, s);            /* %$m (s - x')$% */
+    dy = F_SUB(f, dy, dy, q);          /* %$y' = m (s - x') - t$% */
+
+    EC_DESTROY(d);
+    d->x = dx;
+    d->y = dy;
+    d->z = dz;
+    MP_DROP(m);
+    MP_DROP(q);
+    MP_DROP(s);
+  }
+  return (d);
+}
+

 static ec *ecadd(ec_curve *c, ec *d, const ec *a, const ec *b)
 {
   if (a == b)
 static ec *ecadd(ec_curve *c, ec *d, const ec *a, const ec *b)
 {
   if (a == b)


@@ -117,29 +230,32 @@ static ec *ecadd(ec_curve *c, ec *d, const ec *a, const ec *b)
     mp *dy, *dx;
 
     if (!MP_EQ(a->x, b->x)) {
     mp *dy, *dx;
 
     if (!MP_EQ(a->x, b->x)) {

-      dy = F_SUB(f, MP_NEW, a->y, b->y);
-      dx = F_SUB(f, MP_NEW, a->x, b->x);
-      dx = F_INV(f, dx, dx);
+      dy = F_SUB(f, MP_NEW, a->y, b->y); /* %$y_0 - y_1$% */
+      dx = F_SUB(f, MP_NEW, a->x, b->x); /* %$x_0 - x_1$% */
+      dx = F_INV(f, dx, dx);           /* %$(x_0 - x_1)^{-1}$% */

       lambda = F_MUL(f, MP_NEW, dy, dx);
       lambda = F_MUL(f, MP_NEW, dy, dx);

-    } else if (!MP_LEN(a->y) || !MP_EQ(a->y, b->y)) {
+                                  /* %$\lambda = (y_0 - y1)/(x_0 - x_1)$% */
+    } else if (F_ZEROP(c->f, a->y) || !MP_EQ(a->y, b->y)) {

       EC_SETINF(d);
       return (d);
     } else {
       ecctx *cc = (ecctx *)c;
       EC_SETINF(d);
       return (d);
     } else {
       ecctx *cc = (ecctx *)c;

-      dx = F_SQR(f, MP_NEW, a->x);
-      dx = F_TPL(f, dx, dx);
-      dx = F_ADD(f, dx, dx, cc->a);
-      dy = F_DBL(f, MP_NEW, a->y);
-      dy = F_INV(f, dy, dy);
+      dx = F_SQR(f, MP_NEW, a->x);     /* %$x_0^2$% */
+      dx = F_TPL(f, dx, dx);           /* %$3 x_0^2$% */
+      dx = F_ADD(f, dx, dx, cc->a);    /* %$3 x_0^2 + A$% */
+      dy = F_DBL(f, MP_NEW, a->y);     /* %$2 y_0$% */
+      dy = F_INV(f, dy, dy);           /* %$(2 y_0)^{-1}$% */

       lambda = F_MUL(f, MP_NEW, dx, dy);
       lambda = F_MUL(f, MP_NEW, dx, dy);

+                                   /* %$\lambda = (3 x_0^2 + A)/(2 y_0)$% */

     }
 
     }
 

-    dx = F_SQR(f, dx, lambda);
-    dx = F_SUB(f, dx, dx, a->x);
-    dx = F_SUB(f, dx, dx, b->x);
-    dy = F_SUB(f, dy, b->x, dx);
-    dy = F_MUL(f, dy, lambda, dy);
+    dx = F_SQR(f, dx, lambda);         /* %$\lambda^2$% */
+    dx = F_SUB(f, dx, dx, a->x);       /* %$\lambda^2 - x_0$% */
+    dx = F_SUB(f, dx, dx, b->x);       /* %$x' = \lambda^2 - x_0 - x_1$% */
+    dy = F_SUB(f, dy, b->x, dx);       /* %$x_1 - x'$% */
+    dy = F_MUL(f, dy, lambda, dy);     /* %$\lambda (x_1 - x')$% */

     dy = F_SUB(f, dy, dy, b->y);
     dy = F_SUB(f, dy, dy, b->y);

+                                      /* %$y' = \lambda (x_1 - x') - y_1$% */

 
     EC_DESTROY(d);
     d->x = dx;
 
     EC_DESTROY(d);
     d->x = dx;


@@ -150,6 +266,101 @@ static ec *ecadd(ec_curve *c, ec *d, const ec *a, const ec *b)
   return (d);
 }
 
   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;
+    mp *p, *q, *r, *w, *u, *s, *dx, *dy, *dz;
+
+    q = F_SQR(f, MP_NEW, a->z);                /* %$z_0^2$% */
+    u = F_MUL(f, MP_NEW, q, b->x);     /* %$u = x_1 z_0^2$% */
+    p = F_MUL(f, MP_NEW, q, b->y);     /* %$y_1 z_0^2$% */
+    s = F_MUL(f, q, p, a->z);          /* %$s = y_1 z_0^3$% */
+
+    w = F_SUB(f, p, a->x, u);          /* %$w = x_0 - u$% */
+    r = F_SUB(f, MP_NEW, a->y, s);     /* %$r = y_0 - s$% */
+    if (F_ZEROP(f, w)) {
+      if (F_ZEROP(f, r)) {
+       MP_DROP(w);
+       MP_DROP(r);
+       MP_DROP(u);
+       MP_DROP(s);
+       return (c->ops->dbl(c, d, a));
+      } else {
+       MP_DROP(w);
+       MP_DROP(r);
+       MP_DROP(u);
+       MP_DROP(s);
+       EC_SETINF(d);
+       return (d);
+      }
+    }
+    u = F_ADD(f, u, u, a->x);          /* %$t = x_0 + u$% */
+    s = F_ADD(f, s, s, a->y);          /* %$m = y_0 + r$% */
+
+    dz = F_MUL(f, MP_NEW, a->z, w);    /* %$z' = z_0 w$% */
+
+    p = F_SQR(f, MP_NEW, w);           /* %$w^2$% */
+    q = F_MUL(f, MP_NEW, p, u);                /* %$t w^2$% */
+    u = F_MUL(f, u, p, w);             /* %$w^3$% */
+    p = F_MUL(f, p, u, s);             /* %$m w^3$% */
+    
+    dx = F_SQR(f, u, r);               /* %$r^2$% */
+    dx = F_SUB(f, dx, dx, q);          /* %$x' = r^2 - t w^2$% */
+
+    s = F_DBL(f, s, dx);               /* %$2 x'$% */
+    q = F_SUB(f, q, q, s);             /* %$v = t w^2 - 2 x'$% */
+    dy = F_MUL(f, s, q, r);            /* %$v r$% */
+    dy = F_SUB(f, dy, dy, p);          /* %$v r - m w^3$% */
+    dy = F_HLV(f, dy, dy);             /* %$y' = (v r - m w^3)/2$% */
+
+    EC_DESTROY(d);
+    d->x = dx;
+    d->y = dy;
+    d->z = dz;
+    MP_DROP(p);
+    MP_DROP(q);
+    MP_DROP(r);
+    MP_DROP(w);
+  }
+  return (d);
+}
+
+static int eccheck(ec_curve *c, const ec *p)
+{
+  ecctx *cc = (ecctx *)c;
+  field *f = c->f;
+  int rc;
+  mp *l = F_SQR(f, MP_NEW, p->y);
+  mp *x = F_SQR(f, MP_NEW, p->x);
+  mp *r = F_MUL(f, MP_NEW, x, p->x);
+  x = F_MUL(f, x, cc->a, p->x);
+  r = F_ADD(f, r, r, x);
+  r = F_ADD(f, r, r, cc->b);
+  rc = MP_EQ(l, r) ? 0 : -1;
+  mp_drop(l);
+  mp_drop(x);
+  mp_drop(r);
+  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;
 static void ecdestroy(ec_curve *c)
 {
   ecctx *cc = (ecctx *)c;


@@ -180,8 +391,37 @@ extern ec_curve *ec_prime(field *f, mp *a, mp *b)
   return (&cc->c);
 }
 
   return (&cc->c);
 }
 

+extern ec_curve *ec_primeproj(field *f, mp *a, mp *b)
+{
+  ecctx *cc = CREATE(ecctx);
+  mp *ax;
+
+  ax = mp_add(MP_NEW, a, MP_THREE);
+  ax = F_IN(f, ax, ax);
+  if (F_ZEROP(f, ax))
+    cc->c.ops = &ec_primeprojxops;
+  else
+    cc->c.ops = &ec_primeprojops;
+  MP_DROP(ax);
+  cc->c.f = f;
+  cc->a = F_IN(f, MP_NEW, a);
+  cc->b = F_IN(f, MP_NEW, b);
+  return (&cc->c);
+}
+

 static const ec_ops ec_primeops = {
 static const ec_ops ec_primeops = {

-  ecdestroy, ec_idin, ec_idout, 0, ecneg, ecadd, ec_stdsub, ecdbl
+  ecdestroy, ec_idin, ec_idout, ec_idfix,
+  0, ecneg, ecadd, ec_stdsub, ecdbl, eccheck
+};
+
+static const ec_ops ec_primeprojops = {
+  ecdestroy, ec_projin, ec_projout, ec_projfix,
+  0, ecneg, ecprojadd, ec_stdsub, ecprojdbl, ecprojcheck
+};
+
+static const ec_ops ec_primeprojxops = {
+  ecdestroy, ec_projin, ec_projout, ec_projfix,
+  0, ecneg, ecprojadd, ec_stdsub, ecprojxdbl, ecprojcheck

 };
 
 /*----- Test rig ----------------------------------------------------------*/
 };
 
 /*----- Test rig ----------------------------------------------------------*/

diff --git a/ec.c b/ec.c
index 4111928..a2b229f 100644 (file)
--- a/ec.c
+++ b/ec.c
@@ -1,6 +1,6 @@
 /* -*-c-*-
  *
 /* -*-c-*-
  *

- * $Id: ec.c,v 1.4.4.1 2003/06/10 13:43:53 mdw Exp $
+ * $Id: ec.c,v 1.4.4.2 2004/03/20 00:13:31 mdw Exp $

  *
  * Elliptic curve definitions
  *
  *
  * Elliptic curve definitions
  *


@@ -30,6 +30,9 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ec.c,v $
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ec.c,v $

+ * 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.4.1  2003/06/10 13:43:53  mdw
  * Simple (non-projective) curves over prime fields now seem to work.
  *


@@ -113,7 +116,7 @@ ec *ec_copy(ec *d, const ec *p) { EC_COPY(d, p); return (d); }
 
 /*----- Standard curve operations -----------------------------------------*/
 
 
 /*----- Standard curve operations -----------------------------------------*/
 

-/* --- @ec_idin@, @ec_idout@ --- *
+/* --- @ec_idin@, @ec_idout@, @ec_idfix@ --- *

  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
  *             @ec *d@ = pointer to the destination
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
  *             @ec *d@ = pointer to the destination


@@ -152,6 +155,12 @@ ec *ec_idout(ec_curve *c, ec *d, const ec *p)
   return (d);
 }
 
   return (d);
 }
 

+ec *ec_idfix(ec_curve *c, ec *d, const ec *p)
+{
+  EC_COPY(d, p);
+  return (d);
+}
+

 /* --- @ec_projin@, @ec_projout@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
 /* --- @ec_projin@, @ec_projout@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve


@@ -182,12 +191,15 @@ ec *ec_projout(ec_curve *c, ec *d, const ec *p)
   if (EC_ATINF(p))
     EC_SETINF(d);
   else {
   if (EC_ATINF(p))
     EC_SETINF(d);
   else {

-    mp *x, *y, *z;
+    mp *x, *y, *z, *zz;

     field *f = c->f;
     z = F_INV(f, MP_NEW, p->z);
     field *f = c->f;
     z = F_INV(f, MP_NEW, p->z);

-    x = F_MUL(f, d->x, p->x, z);
+    zz = F_SQR(f, MP_NEW, z);
+    z = F_MUL(f, z, zz, z);
+    x = F_MUL(f, d->x, p->x, zz);

     y = F_MUL(f, d->y, p->y, z);
     mp_drop(z);
     y = F_MUL(f, d->y, p->y, z);
     mp_drop(z);

+    mp_drop(zz);

     mp_drop(d->z);
     d->x = F_OUT(f, x, x);
     d->y = F_OUT(f, y, y);
     mp_drop(d->z);
     d->x = F_OUT(f, x, x);
     d->y = F_OUT(f, y, y);


@@ -196,6 +208,28 @@ ec *ec_projout(ec_curve *c, ec *d, const ec *p)
   return (d);
 }
 
   return (d);
 }
 

+ec *ec_projfix(ec_curve *c, ec *d, const ec *p)
+{
+  if (EC_ATINF(p))
+    EC_SETINF(d);
+  else if (d->z == c->f->one)
+    EC_COPY(d, p);
+  else {
+    mp *z, *zz;
+    field *f = c->f;
+    z = F_INV(f, MP_NEW, p->z);
+    zz = F_SQR(f, MP_NEW, z);
+    z = F_MUL(f, z, zz, z);
+    d->x = F_MUL(f, d->x, p->x, zz);
+    d->y = F_MUL(f, d->y, p->y, z);
+    mp_drop(z);
+    mp_drop(zz);
+    mp_drop(d->z);
+    d->z = MP_COPY(f->one);
+  }
+  return (d);  
+}
+

 /* --- @ec_stdsub@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
 /* --- @ec_stdsub@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve


@@ -213,6 +247,7 @@ ec *ec_stdsub(ec_curve *c, ec *d, const ec *p, const ec *q)
 {
   ec t = EC_INIT;
   EC_NEG(c, &t, q);
 {
   ec t = EC_INIT;
   EC_NEG(c, &t, q);

+  EC_FIX(c, &t, &t);

   EC_ADD(c, d, p, &t);
   EC_DESTROY(&t);
   return (d);
   EC_ADD(c, d, p, &t);
   EC_DESTROY(&t);
   return (d);


@@ -249,7 +284,7 @@ ec *ec_find(ec_curve *c, ec *d, mp *x)
   x = F_IN(c->f, MP_NEW, x);
   if ((d = EC_FIND(c, d, x)) != 0)
     EC_OUT(c, d, d);
   x = F_IN(c->f, MP_NEW, x);
   if ((d = EC_FIND(c, d, x)) != 0)
     EC_OUT(c, d, d);

-  mp_drop(x);
+  MP_DROP(x);

   return (d);
 }
 
   return (d);
 }
 


@@ -335,6 +370,29 @@ ec *ec_dbl(ec_curve *c, ec *d, const ec *p)
   return (EC_OUT(c, d, d));
 }
 
   return (EC_OUT(c, d, d));
 }
 

+/* --- @ec_check@ --- *
+ *
+ * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
+ *             @const ec *p@ = pointer to the point
+ *
+ * Returns:    Zero if OK, nonzero if this is an invalid point.
+ *
+ * Use:                Checks that a point is actually on an elliptic curve.
+ */
+
+int ec_check(ec_curve *c, const ec *p)
+{
+  ec t = EC_INIT;
+  int rc;
+
+  if (EC_ATINF(p))
+    return (0);
+  EC_IN(c, &t, p);
+  rc = EC_CHECK(c, &t);
+  EC_DESTROY(&t);
+  return (rc);
+}
+

 /* --- @ec_imul@, @ec_mul@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
 /* --- @ec_imul@, @ec_mul@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve


@@ -360,10 +418,14 @@ ec *ec_imul(ec_curve *c, ec *d, const ec *p, mp *n)
   EC_SETINF(d);
   if (MP_LEN(n) == 0)
     ;
   EC_SETINF(d);
   if (MP_LEN(n) == 0)
     ;

-  else if (MP_LEN(n) < EXP_THRESH)
-    EXP_SIMPLE(*d, t, n);
-  else
-    EXP_WINDOW(*d, t, n);
+  else {
+    if (n->f & MP_NEG)
+      EC_NEG(c, &t, &t);
+    if (MP_LEN(n) < EXP_THRESH)
+      EXP_SIMPLE(*d, t, n);
+    else
+      EXP_WINDOW(*d, t, n);
+  }

   EC_DESTROY(&t);
   return (d);
 }
   EC_DESTROY(&t);
   return (d);
 }

diff --git a/ec.h b/ec.h
index 72ee7a9..476d898 100644 (file)
--- a/ec.h
+++ b/ec.h
@@ -1,6 +1,6 @@
 /* -*-c-*-
  *
 /* -*-c-*-
  *

- * $Id: ec.h,v 1.4.4.1 2003/06/10 13:43:53 mdw Exp $
+ * $Id: ec.h,v 1.4.4.2 2004/03/20 00:13:31 mdw Exp $

  *
  * Elliptic curve definitions
  *
  *
  * Elliptic curve definitions
  *


@@ -30,6 +30,9 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ec.h,v $
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: ec.h,v $

+ * 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.4.1  2003/06/10 13:43:53  mdw
  * Simple (non-projective) curves over prime fields now seem to work.
  *


@@ -83,27 +86,37 @@ typedef struct ec_mulfactor {
   mp *exp;                             /* The exponent */
 } ec_mulfactor;
 
   mp *exp;                             /* The exponent */
 } ec_mulfactor;
 

-/* --- Elliptic curve operations --- */
+/* --- Elliptic curve operations --- *
+ *
+ * All operations (apart from @destroy@ and @in@) are guaranteed to be
+ * performed on internal representations of points.  Moreover, the second
+ * argument to @add@ and @mul@ is guaranteed to be the output of @in@ or
+ * @fix@.
+ */

 
 typedef struct ec_ops {
   void (*destroy)(ec_curve */*c*/);
   ec *(*in)(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
   ec *(*out)(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 
 typedef struct ec_ops {
   void (*destroy)(ec_curve */*c*/);
   ec *(*in)(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
   ec *(*out)(ec_curve */*c*/, ec */*d*/, const ec */*p*/);

+  ec *(*fix)(ec_curve */*c*/, ec */*d*/, const ec */*p*/);

   ec *(*find)(ec_curve */*c*/, ec */*d*/, mp */*x*/);
   ec *(*neg)(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
   ec *(*add)(ec_curve */*c*/, ec */*d*/, const ec */*p*/, const ec */*q*/);
   ec *(*sub)(ec_curve */*c*/, ec */*d*/, const ec */*p*/, const ec */*q*/);
   ec *(*dbl)(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
   ec *(*find)(ec_curve */*c*/, ec */*d*/, mp */*x*/);
   ec *(*neg)(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
   ec *(*add)(ec_curve */*c*/, ec */*d*/, const ec */*p*/, const ec */*q*/);
   ec *(*sub)(ec_curve */*c*/, ec */*d*/, const ec */*p*/, const ec */*q*/);
   ec *(*dbl)(ec_curve */*c*/, ec */*d*/, const ec */*p*/);

+  int (*check)(ec_curve */*c*/, const ec */*p*/);

 } ec_ops;
 
 #define EC_IN(c, d, p)         (c)->ops->in((c), (d), (p))
 #define EC_OUT(c, d, p)                (c)->ops->out((c), (d), (p))
 } ec_ops;
 
 #define EC_IN(c, d, p)         (c)->ops->in((c), (d), (p))
 #define EC_OUT(c, d, p)                (c)->ops->out((c), (d), (p))

+#define EC_FIX(c, d, p)                (c)->ops->fix((c), (d), (p))

 
 #define EC_FIND(c, d, x)       (c)->ops->find((c), (d), (x))
 #define EC_NEG(c, d, x)                (c)->ops->neg((c), (d), (x))
 #define EC_ADD(c, d, p, q)     (c)->ops->add((c), (d), (p), (q))
 #define EC_SUB(c, d, p, q)     (c)->ops->sub((c), (d), (p), (q))
 #define EC_DBL(c, d, p)                (c)->ops->dbl((c), (d), (p))
 
 #define EC_FIND(c, d, x)       (c)->ops->find((c), (d), (x))
 #define EC_NEG(c, d, x)                (c)->ops->neg((c), (d), (x))
 #define EC_ADD(c, d, p, q)     (c)->ops->add((c), (d), (p), (q))
 #define EC_SUB(c, d, p, q)     (c)->ops->sub((c), (d), (p), (q))
 #define EC_DBL(c, d, p)                (c)->ops->dbl((c), (d), (p))

+#define EC_CHECK(c, p)         (c)->ops->check((c), (p))

 
 /*----- Simple memory management things -----------------------------------*/
 
 
 /*----- Simple memory management things -----------------------------------*/
 


@@ -210,32 +223,6 @@ extern ec *ec_copy(ec */*d*/, const ec */*p*/);
 
 /*----- Interesting arithmetic --------------------------------------------*/
 
 
 /*----- Interesting arithmetic --------------------------------------------*/
 

-/* --- @ec_in@ --- *
- *
- * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
- *             @ec *d@ = pointer to the destination point
- *             @const ec *p@ = pointer to the source point
- *
- * Returns:    The destination point.
- *
- * Use:                Converts a point to internal representation.
- */
-
-extern ec *ec_in(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
-
-/* --- @ec_out@ --- *
- *
- * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
- *             @ec *d@ = pointer to the destination point
- *             @const ec *p@ = pointer to the source point
- *
- * Returns:    The destination point.
- *
- * Use:                Converts a point to external representation.
- */
-
-extern ec *ec_out(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
-

 /* --- @ec_find@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
 /* --- @ec_find@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve


@@ -306,6 +293,18 @@ extern ec *ec_sub(ec_curve */*c*/, ec */*d*/,
 
 extern ec *ec_dbl(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 
 
 extern ec *ec_dbl(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 

+/* --- @ec_check@ --- *
+ *
+ * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
+ *             @const ec *p@ = pointer to the point
+ *
+ * Returns:    Zero if OK, nonzero if this is an invalid point.
+ *
+ * Use:                Checks that a point is actually on an elliptic curve.
+ */
+
+extern int ec_check(ec_curve */*c*/, const ec */*p*/);
+

 /* --- @ec_mul@, @ec_imul@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
 /* --- @ec_mul@, @ec_imul@ --- *
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve


@@ -343,7 +342,7 @@ extern ec *ec_immul(ec_curve */*c*/, ec */*d*/,
 
 /*----- Standard curve operations -----------------------------------------*/
 
 
 /*----- Standard curve operations -----------------------------------------*/
 

-/* --- @ec_idin@, @ec_idout@ --- *
+/* --- @ec_idin@, @ec_idout@, @ec_idfix@ --- *

  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
  *             @ec *d@ = pointer to the destination
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
  *             @ec *d@ = pointer to the destination


@@ -358,8 +357,9 @@ extern ec *ec_immul(ec_curve */*c*/, ec */*d*/,
 
 extern ec *ec_idin(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 extern ec *ec_idout(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 
 extern ec *ec_idin(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 extern ec *ec_idout(ec_curve */*c*/, ec */*d*/, const ec */*p*/);

+extern ec *ec_idfix(ec_curve */*c*/, ec */*d*/, const ec */*p*/);

 
 

-/* --- @ec_projin@, @ec_projout@ --- *
+/* --- @ec_projin@, @ec_projout@, @ec_projfix@ --- *

  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
  *             @ec *d@ = pointer to the destination
  *
  * Arguments:  @ec_curve *c@ = pointer to an elliptic curve
  *             @ec *d@ = pointer to the destination


@@ -373,6 +373,7 @@ extern ec *ec_idout(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 
 extern ec *ec_projin(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 extern ec *ec_projout(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 
 extern ec *ec_projin(ec_curve */*c*/, ec */*d*/, const ec */*p*/);
 extern ec *ec_projout(ec_curve */*c*/, ec */*d*/, const ec */*p*/);

+extern ec *ec_projfix(ec_curve */*c*/, ec */*d*/, const ec */*p*/);

 
 /* --- @ec_stdsub@ --- *
  *
 
 /* --- @ec_stdsub@ --- *
  *

diff --git a/exp.h b/exp.h
index 6cfdfd8..fc9e3a9 100644 (file)
--- a/exp.h
+++ b/exp.h
@@ -1,6 +1,6 @@
 /* -*-c-*-
  *
 /* -*-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 $

  *
  * Generalized exponentiation
  *
  *
  * Generalized exponentiation
  *


@@ -30,6 +30,9 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: exp.h,v $
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: exp.h,v $

+ * 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.
  * Revision 1.1  2001/06/16 13:00:59  mdw
  * New generic exponentation code.  Includes sliding-window simultaneous
  * exponentiation.


@@ -99,6 +102,10 @@ typedef struct exp_simul {
  * @EXP_MUL(a, x)@             Multiplies @a@ by @x@ (writing the result
  *                             back to @a@).
  *
  * @EXP_MUL(a, x)@             Multiplies @a@ by @x@ (writing the result
  *                             back to @a@).
  *

+ * @EXP_FIX(x)@                        Makes @x@ be a canonical representation of
+ *                             its value.  All multiplications have the
+ *                             right argument canonical.
+ *

  * @EXP_SQR(a)@                        Multiplies @a@ by itself.
  *
  * @EXP_SETMUL(d, x, y)@       Sets @d@ to be the product of @x@ and @y@.
  * @EXP_SQR(a)@                        Multiplies @a@ by itself.
  *
  * @EXP_SETMUL(d, x, y)@       Sets @d@ to be the product of @x@ and @y@.


@@ -140,6 +147,7 @@ typedef struct exp_simul {
                                                                        \
   /* --- Do the main body of the work --- */                           \
                                                                        \
                                                                        \
   /* --- Do the main body of the work --- */                           \
                                                                        \

+  EXP_FIX(g);                                                          \

   for (;;) {                                                           \
     EXP_MUL(a, g);                                                     \
     sq = 0;                                                            \
   for (;;) {                                                           \
     EXP_MUL(a, g);                                                     \
     sq = 0;                                                            \


@@ -184,11 +192,15 @@ exp_simple_exit:;                                                 \
                                                                        \
   /* --- Do the precomputation --- */                                  \
                                                                        \
                                                                        \
   /* --- Do the precomputation --- */                                  \
                                                                        \

+  EXP_FIX(g);                                                          \

   EXP_SETSQR(g2, g);                                                   \
   EXP_SETSQR(g2, g);                                                   \

+  EXP_FIX(g2);                                                         \

   v = xmalloc(EXP_TABSZ * sizeof(EXP_TYPE));                           \
   EXP_COPY(v[0], g);                                                   \
   v = xmalloc(EXP_TABSZ * sizeof(EXP_TYPE));                           \
   EXP_COPY(v[0], g);                                                   \

-  for (i = 1; i < EXP_TABSZ; i++)                                      \
+  for (i = 1; i < EXP_TABSZ; i++) {                                    \

     EXP_SETMUL(v[i], v[i - 1], g2);                                    \
     EXP_SETMUL(v[i], v[i - 1], g2);                                    \

+    EXP_FIX(v[i]);                                                     \
+  }                                                                    \

   EXP_DROP(g2);                                                                \
                                                                        \
   /* --- Skip top-end zero bits --- *                                  \
   EXP_DROP(g2);                                                                \
                                                                        \
   /* --- Skip top-end zero bits --- *                                  \


@@ -286,17 +298,21 @@ exp_window_exit:;                                                 \
   j = 1;                                                               \
   for (i = 0; i < n; i++) {                                            \
     EXP_COPY(v[j], f[n - 1 - i].base);                                 \
   j = 1;                                                               \
   for (i = 0; i < n; i++) {                                            \
     EXP_COPY(v[j], f[n - 1 - i].base);                                 \

+    EXP_FIX(v[j]);                                                     \

     j <<= 1;                                                           \
   }                                                                    \
   k = n * EXP_WINSZ;                                                   \
   jj = 1;                                                              \
   for (; i < k; i++) {                                                 \
     EXP_SETSQR(v[j], v[jj]);                                           \
     j <<= 1;                                                           \
   }                                                                    \
   k = n * EXP_WINSZ;                                                   \
   jj = 1;                                                              \
   for (; i < k; i++) {                                                 \
     EXP_SETSQR(v[j], v[jj]);                                           \

+    EXP_FIX(v[j]);                                                     \

     j <<= 1; jj <<= 1;                                                 \
   }                                                                    \
   for (i = 1; i < vn; i <<= 1) {                                       \
     j <<= 1; jj <<= 1;                                                 \
   }                                                                    \
   for (i = 1; i < vn; i <<= 1) {                                       \

-    for (j = 1; j < i; j++)                                            \
+    for (j = 1; j < i; j++) {                                          \

       EXP_SETMUL(v[j + i], v[j], v[i]);                                        \
       EXP_SETMUL(v[j + i], v[j], v[i]);                                        \

+      EXP_FIX(v[j + i]);                                               \
+    }                                                                  \

   }                                                                    \
                                                                        \
   /* --- Set up the bitscanners --- *                                  \
   }                                                                    \
                                                                        \
   /* --- Set up the bitscanners --- *                                  \


@@ -381,7 +397,7 @@ exp_window_exit:;                                                   \
                                                                        \
 exp_simul_done:                                                                \
   while (sq--) EXP_SQR(a);                                             \
                                                                        \
 exp_simul_done:                                                                \
   while (sq--) EXP_SQR(a);                                             \

-  for (i = 1; i < vn; i++)                                     \
+  for (i = 1; i < vn; i++)                                             \

     EXP_DROP(v[i]);                                                    \
   xfree(v);                                                            \
 } while (0)
     EXP_DROP(v[i]);                                                    \
   xfree(v);                                                            \
 } while (0)

diff --git a/f-prime.c b/f-prime.c
index 0d76da0..7c1dae5 100644 (file)
--- a/f-prime.c
+++ b/f-prime.c
@@ -1,6 +1,6 @@
 /* -*-c-*-
  *
 /* -*-c-*-
  *

- * $Id: f-prime.c,v 1.3.4.1 2003/06/10 13:43:53 mdw Exp $
+ * $Id: f-prime.c,v 1.3.4.2 2004/03/20 00:13:31 mdw Exp $

  *
  * Prime fields with Montgomery arithmetic
  *
  *
  * Prime fields with Montgomery arithmetic
  *


@@ -30,6 +30,9 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: f-prime.c,v $
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: f-prime.c,v $

+ * 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.4.1  2003/06/10 13:43:53  mdw
  * Simple (non-projective) curves over prime fields now seem to work.
  *


@@ -82,6 +85,11 @@ static mp *fout(field *ff, mp *d, mp *x)
   return (mpmont_reduce(&f->mm, d, x));
 }
 
   return (mpmont_reduce(&f->mm, d, x));
 }
 

+static int fzerop(field *ff, mp *x)
+{
+  return (!MP_LEN(x));
+}
+

 static mp *fneg(field *ff, mp *d, mp *x)
 {
   fctx *f = (fctx *)ff;
 static mp *fneg(field *ff, mp *d, mp *x)
 {
   fctx *f = (fctx *)ff;


@@ -157,11 +165,37 @@ static mp *ftpl(field *ff, mp *d, mp *x)
   return (d);
 }
 
   return (d);
 }
 

+static mp *fqdl(field *ff, mp *d, mp *x)
+{
+  fctx *f = (fctx *)ff;
+  d = mp_lsl(d, x, 2);
+  while (MP_CMP(d, >, f->mm.m))
+    d = mp_sub(d, d, f->mm.m);
+  return (d);
+}
+
+static mp *fhlv(field *ff, mp *d, mp *x)
+{
+  fctx *f = (fctx *)ff;
+  if (!MP_LEN(x)) {
+    MP_COPY(x);
+    MP_DROP(d);
+    return (x);
+  }
+  if (x->v[0] & 1) {
+    d = mp_add(d, x, f->mm.m);
+    x = d;
+  }
+  return (mp_lsr(d, x, 1));
+}
+

 static mp *fsqrt(field *ff, mp *d, mp *x)
 {
   fctx *f = (fctx *)ff;
   d = mpmont_reduce(&f->mm, d, x);
   d = mp_modsqrt(d, d, f->mm.m);
 static mp *fsqrt(field *ff, mp *d, mp *x)
 {
   fctx *f = (fctx *)ff;
   d = mpmont_reduce(&f->mm, d, x);
   d = mp_modsqrt(d, d, f->mm.m);

+  if (!d)
+    return (d);

   return (mpmont_mul(&f->mm, d, d, f->mm.r2));
 }
 
   return (mpmont_mul(&f->mm, d, d, f->mm.r2));
 }
 


@@ -170,8 +204,8 @@ static mp *fsqrt(field *ff, mp *d, mp *x)
 static field_ops fops = {
   fdestroy,
   fin, fout,
 static field_ops fops = {
   fdestroy,
   fin, fout,

-  fneg, fadd, fsub, fmul, fsqr, finv, freduce,
-  fdbl, ftpl, fsqrt
+  fzerop, fneg, fadd, fsub, fmul, fsqr, finv, freduce,
+  fdbl, ftpl, fqdl, fhlv, fsqrt

 };
 
 /* --- @field_prime@ --- *
 };
 
 /* --- @field_prime@ --- *

diff --git a/field.h b/field.h
index 5a70649..8096700 100644 (file)
--- a/field.h
+++ b/field.h
@@ -1,6 +1,6 @@
 /* -*-c-*-
  *
 /* -*-c-*-
  *

- * $Id: field.h,v 1.3 2002/01/13 13:48:44 mdw Exp $
+ * $Id: field.h,v 1.3.4.1 2004/03/20 00:13:31 mdw Exp $

  *
  * Definitions for field arithmetic
  *
  *
  * Definitions for field arithmetic
  *


@@ -30,6 +30,9 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: field.h,v $
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: field.h,v $

+ * 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.
  *
  * Revision 1.3  2002/01/13 13:48:44  mdw
  * Further progress.
  *


@@ -70,6 +73,7 @@ typedef struct field_ops {
   mp *(*in)(field */*f*/, mp */*d*/, mp */*x*/);
   mp *(*out)(field */*f*/, mp */*d*/, mp */*x*/);
 
   mp *(*in)(field */*f*/, mp */*d*/, mp */*x*/);
   mp *(*out)(field */*f*/, mp */*d*/, mp */*x*/);
 

+  int (*zerop)(field */*f*/, mp */*x*/);

   mp *(*neg)(field */*f*/, mp */*d*/, mp */*x*/);
   mp *(*add)(field */*f*/, mp */*d*/, mp */*x*/, mp */*y*/);
   mp *(*sub)(field */*f*/, mp */*d*/, mp */*x*/, mp */*y*/);
   mp *(*neg)(field */*f*/, mp */*d*/, mp */*x*/);
   mp *(*add)(field */*f*/, mp */*d*/, mp */*x*/, mp */*y*/);
   mp *(*sub)(field */*f*/, mp */*d*/, mp */*x*/, mp */*y*/);


@@ -82,6 +86,8 @@ typedef struct field_ops {
 
   mp *(*dbl)(field */*f*/, mp */*d*/, mp */*x*/);
   mp *(*tpl)(field */*f*/, mp */*d*/, mp */*x*/);
 
   mp *(*dbl)(field */*f*/, mp */*d*/, mp */*x*/);
   mp *(*tpl)(field */*f*/, mp */*d*/, mp */*x*/);

+  mp *(*qdl)(field */*f*/, mp */*d*/, mp */*x*/);
+  mp *(*hlv)(field */*f*/, mp */*d*/, mp */*x*/);

   mp *(*sqrt)(field */*f*/, mp */*d*/, mp */*x*/);
 
 } field_ops;
   mp *(*sqrt)(field */*f*/, mp */*d*/, mp */*x*/);
 
 } field_ops;


@@ -90,6 +96,8 @@ typedef struct field_ops {
 
 #define F_IN(f, d, x)          (f)->ops->in((f), (d), (x))
 #define F_OUT(f, d, x)         (f)->ops->out((f), (d), (x))
 
 #define F_IN(f, d, x)          (f)->ops->in((f), (d), (x))
 #define F_OUT(f, d, x)         (f)->ops->out((f), (d), (x))

+
+#define F_ZEROP(f, x)          (f)->ops->zerop((f), (x))

 #define F_NEG(f, d, x)         (f)->ops->neg((f), (d), (x))
 #define F_ADD(f, d, x, y)      (f)->ops->add((f), (d), (x), (y))
 #define F_SUB(f, d, x, y)      (f)->ops->sub((f), (d), (x), (y))
 #define F_NEG(f, d, x)         (f)->ops->neg((f), (d), (x))
 #define F_ADD(f, d, x, y)      (f)->ops->add((f), (d), (x), (y))
 #define F_SUB(f, d, x, y)      (f)->ops->sub((f), (d), (x), (y))


@@ -100,6 +108,8 @@ typedef struct field_ops {
 
 #define F_DBL(f, d, x)         (f)->ops->dbl((f), (d), (x))
 #define F_TPL(f, d, x)         (f)->ops->tpl((f), (d), (x))
 
 #define F_DBL(f, d, x)         (f)->ops->dbl((f), (d), (x))
 #define F_TPL(f, d, x)         (f)->ops->tpl((f), (d), (x))

+#define F_QDL(f, d, x)         (f)->ops->qdl((f), (d), (x))
+#define F_HLV(f, d, x)         (f)->ops->hlv((f), (d), (x))

 #define F_SQRT(f, d, x)                (f)->ops->sqrt((f), (d), (x))
 
 /*----- Helpful field operations ------------------------------------------*/
 #define F_SQRT(f, d, x)                (f)->ops->sqrt((f), (d), (x))
 
 /*----- Helpful field operations ------------------------------------------*/

diff --git a/mpmont-exp.h b/mpmont-exp.h
index e732e6f..0f82b9c 100644 (file)
--- a/mpmont-exp.h
+++ b/mpmont-exp.h
@@ -1,6 +1,6 @@
 /* -*-c-*-
  *
 /* -*-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 $

  *
  * Exponentiation operations for Montgomery reduction
  *
  *
  * Exponentiation operations for Montgomery reduction
  *


@@ -30,6 +30,9 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: mpmont-exp.h,v $
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: mpmont-exp.h,v $

+ * 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.
  *
  * Revision 1.1  2001/06/16 12:58:12  mdw
  * Parameters for generic exponentiation.
  *


@@ -61,6 +64,8 @@
   a = mpmont_reduce(mm, t, t);                                         \
 } while (0)
 
   a = mpmont_reduce(mm, t, t);                                         \
 } while (0)
 

+#define EXP_FIX(x)
+

 #define EXP_SETMUL(d, x, y) d = mpmont_mul(mm, MP_NEW, x, y)
 
 #define EXP_SETSQR(d, x) do {                                          \
 #define EXP_SETMUL(d, x, y) d = mpmont_mul(mm, MP_NEW, x, y)
 
 #define EXP_SETSQR(d, x) do {                                          \