libsec: implement elliptic curve group operations in jacobian coordinate system

1 files changed, 74 insertions, 24 deletions

diff --git a/sys/src/libsec/port/ecc.c b/sys/src/libsec/port/ecc.c
index e4e4df0fd..33038fa24 100644
--- a/sys/src/libsec/port/ecc.c
+++ b/sys/src/libsec/port/ecc.c
@@ -3,12 +3,33 @@
 #include <libsec.h>
 #include <ctype.h>
 
+extern void jacobian_affine(mpint *p,
+	mpint *X, mpint *Y, mpint *Z);
+extern void jacobian_dbl(mpint *p, mpint *a,
+	mpint *X1, mpint *Y1, mpint *Z1,
+	mpint *X3, mpint *Y3, mpint *Z3);
+extern void jacobian_add(mpint *p, mpint *a,
+	mpint *X1, mpint *Y1, mpint *Z1,
+	mpint *X2, mpint *Y2, mpint *Z2,
+	mpint *X3, mpint *Y3, mpint *Z3);
+
 void
-ecassign(ECdomain *, ECpoint *a, ECpoint *b)
+ecassign(ECdomain *dom, ECpoint *a, ECpoint *b)
 {
-	b->inf = a->inf;
+	if((b->inf = a->inf) != 0)
+		return;
 	mpassign(a->x, b->x);
 	mpassign(a->y, b->y);
+	if(b->z != nil){
+		mpassign(a->z != nil ? a->z : mpone, b->z);
+		return;
+	}
+	if(a->z != nil){
+		b->z = mpcopy(a->z);
+		jacobian_affine(dom->p, b->x, b->y, b->z);
+		mpfree(b->z);
+		b->z = nil;
+	}
 }
 
 void
@@ -28,10 +49,26 @@ ecadd(ECdomain *dom, ECpoint *a, ECpoint *b, ECpoint *s)
 		ecassign(dom, a, s);
 		return;
 	}
+
+	if(s->z != nil){
+		if(a == b)
+			jacobian_dbl(dom->p, dom->a,
+				a->x, a->y, a->z != nil ? a->z : mpone,
+				s->x, s->y, s->z);
+		else
+			jacobian_add(dom->p, dom->a,
+				a->x, a->y, a->z != nil ? a->z : mpone,
+				b->x, b->y, b->z != nil ? b->z : mpone,
+				s->x, s->y, s->z);
+		s->inf = mpcmp(s->z, mpzero) == 0;
+		return;
+	}
+
 	if(mpcmp(a->x, b->x) == 0 && (mpcmp(a->y, mpzero) == 0 || mpcmp(a->y, b->y) != 0)){
 		s->inf = 1;
 		return;
 	}
+	s->inf = 0;
 	l = mpnew(0);
 	k = mpnew(0);
 	sx = mpnew(0);
@@ -103,8 +140,10 @@ ecmul(ECdomain *dom, ECpoint *a, mpint *k, ECpoint *s)
 	ns.inf = 1;
 	ns.x = mpnew(0);
 	ns.y = mpnew(0);
+	ns.z = mpnew(0);
 	na.x = mpnew(0);
 	na.y = mpnew(0);
+	na.z = mpnew(0);
 	ecassign(dom, a, &na);
 	l = mpcopy(k);
 	l->sign = 1;
@@ -114,15 +153,17 @@ ecmul(ECdomain *dom, ECpoint *a, mpint *k, ECpoint *s)
 		ecadd(dom, &na, &na, &na);
 		mpright(l, 1, l);
 	}
-	if(k->sign < 0){
+	if(k->sign < 0 && !ns.inf){
 		ns.y->sign = -1;
 		mpmod(ns.y, dom->p, ns.y);
 	}
 	ecassign(dom, &ns, s);
 	mpfree(ns.x);
 	mpfree(ns.y);
+	mpfree(ns.z);
 	mpfree(na.x);
 	mpfree(na.y);
+	mpfree(na.z);
 	mpfree(l);
 }
 
@@ -132,18 +173,17 @@ ecverify(ECdomain *dom, ECpoint *a)
 	mpint *p, *q;
 	int r;
 
+	assert(a->z == nil);	/* need affine coordinates */
 	if(a->inf)
 		return 1;
-	
+
 	p = mpnew(0);
 	q = mpnew(0);
-	mpmul(a->y, a->y, p);
-	mpmod(p, dom->p, p);
-	mpmul(a->x, a->x, q);
-	mpadd(q, dom->a, q);
-	mpmul(a->x, q, q);
-	mpadd(q, dom->b, q);
-	mpmod(q, dom->p, q);
+	mpmodmul(a->y, a->y, dom->p, p);
+	mpmodmul(a->x, a->x, dom->p, q);
+	mpmodadd(q, dom->a, dom->p, q);
+	mpmodmul(q, a->x, dom->p, q);
+	mpmodadd(q, dom->b, dom->p, q);
 	r = mpcmp(p, q);
 	mpfree(p);
 	mpfree(q);
@@ -162,10 +202,12 @@ ecpubverify(ECdomain *dom, ECpub *a)
 		return 0;
 	p.x = mpnew(0);
 	p.y = mpnew(0);
+	p.z = mpnew(0);
 	ecmul(dom, a, dom->n, &p);
 	r = p.inf;
 	mpfree(p.x);
 	mpfree(p.y);
+	mpfree(p.z);
 	return r;
 }
 
@@ -272,8 +314,8 @@ mpsqrt(mpint *n, mpint *p, mpint *r)
 	zq = mpnew(0);
 	for(;;){
 		for(;;){
-			mprand(mpsignif(p), genrandom, a);
-			if(mpcmp(a, mpzero) > 0 && mpcmp(a, p) < 0)
+			mpnrand(p, genrandom, a);
+			if(mpcmp(a, mpzero) > 0)
 				break;
 		}
 		mpmul(a, a, t);
@@ -345,11 +387,12 @@ strtoec(ECdomain *dom, char *s, char **rptr, ECpoint *ret)
 		ret->x = mpnew(0);
 		ret->y = mpnew(0);
 	}
+	ret->inf = 0;
 	o = 0;
 	switch(octet(&s)){
 	case 0:
 		ret->inf = 1;
-		return ret;
+		break;
 	case 3:
 		o = 1;
 	case 2:
@@ -370,7 +413,7 @@ strtoec(ECdomain *dom, char *s, char **rptr, ECpoint *ret)
 		mpfree(r);
 		if(!ecverify(dom, ret))
 			goto err;
-		return ret;
+		break;
 	case 4:
 		if(halfpt(dom, s, &s, ret->x) == nil)
 			goto err;
@@ -378,8 +421,12 @@ strtoec(ECdomain *dom, char *s, char **rptr, ECpoint *ret)
 			goto err;
 		if(!ecverify(dom, ret))
 			goto err;
-		return ret;
+		break;
 	}
+	if(ret->z != nil && !ret->inf)
+		mpassign(mpone, ret->z);
+	return ret;
+
 err:
 	if(rptr)
 		*rptr = s;
@@ -403,8 +450,8 @@ ecgen(ECdomain *dom, ECpriv *p)
 		p->d = mpnew(0);
 	}
 	for(;;){
-		mprand(mpsignif(dom->n), genrandom, p->d);
-		if(mpcmp(p->d, mpzero) > 0 && mpcmp(p->d, dom->n) < 0)
+		mpnrand(dom->n, genrandom, p->d);
+		if(mpcmp(p->d, mpzero) > 0)
 			break;
 	}
 	ecmul(dom, &dom->G, p->d, p);
@@ -419,6 +466,7 @@ ecdsasign(ECdomain *dom, ECpriv *priv, uchar *dig, int len, mpint *r, mpint *s)
 
 	tmp.x = mpnew(0);
 	tmp.y = mpnew(0);
+	tmp.z = nil;
 	tmp.d = mpnew(0);
 	E = betomp(dig, len, nil);
 	t = mpnew(0);
@@ -432,8 +480,7 @@ ecdsasign(ECdomain *dom, ECpriv *priv, uchar *dig, int len, mpint *r, mpint *s)
 		mpmul(r, priv->d, s);
 		mpadd(E, s, s);
 		mpinvert(tmp.d, dom->n, t);
-		mpmul(s, t, s);
-		mpmod(s, dom->n, s);
+		mpmodmul(s, t, dom->n, s);
 		if(mpcmp(s, mpzero) != 0)
 			break;
 	}
@@ -461,18 +508,19 @@ ecdsaverify(ECdomain *dom, ECpub *pub, uchar *dig, int len, mpint *r, mpint *s)
 	u2 = mpnew(0);
 	R.x = mpnew(0);
 	R.y = mpnew(0);
+	R.z = mpnew(0);
 	S.x = mpnew(0);
 	S.y = mpnew(0);
+	S.z = mpnew(0);
 	mpinvert(s, dom->n, t);
-	mpmul(E, t, u1);
-	mpmod(u1, dom->n, u1);
-	mpmul(r, t, u2);
-	mpmod(u2, dom->n, u2);
+	mpmodmul(E, t, dom->n, u1);
+	mpmodmul(r, t, dom->n, u2);
 	ecmul(dom, &dom->G, u1, &R);
 	ecmul(dom, pub, u2, &S);
 	ecadd(dom, &R, &S, &R);
 	ret = 0;
 	if(!R.inf){
+		jacobian_affine(dom->p, R.x, R.y, R.z);
 		mpmod(R.x, dom->n, t);
 		ret = mpcmp(r, t) == 0;
 	}
@@ -482,8 +530,10 @@ ecdsaverify(ECdomain *dom, ECpub *pub, uchar *dig, int len, mpint *r, mpint *s)
 	mpfree(u2);
 	mpfree(R.x);
 	mpfree(R.y);
+	mpfree(R.z);
 	mpfree(S.x);
 	mpfree(S.y);
+	mpfree(S.z);
 	return ret;
 }