+++ /dev/null
-/* Copyright (c) 1998,2011,2014 Apple Inc. All Rights Reserved.
- *
- * NOTICE: USE OF THE MATERIALS ACCOMPANYING THIS NOTICE IS SUBJECT
- * TO THE TERMS OF THE SIGNED "FAST ELLIPTIC ENCRYPTION (FEE) REFERENCE
- * SOURCE CODE EVALUATION AGREEMENT" BETWEEN APPLE, INC. AND THE
- * ORIGINAL LICENSEE THAT OBTAINED THESE MATERIALS FROM APPLE,
- * INC. ANY USE OF THESE MATERIALS NOT PERMITTED BY SUCH AGREEMENT WILL
- * EXPOSE YOU TO LIABILITY.
- ***************************************************************************
- *
- * ellipticProj.c - elliptic projective algebra routines.
- *
- * Revision History
- * ----------------
- * 1 Sep 1998 at Apple
- * Created.
- *
- **************************************************************
-
- PROJECTIVE FORMAT
-
- Functions are supplied herein for projective format
- of points. Alternative formats differ in their
- range of applicability, efficiency, and so on.
- Primary advantages of the projective format herein are:
- -- No explicit inversions (until perhaps one such at the end of
- an elliptic multiply operation)
- -- Fairly low operation count (~11 muls for point doubling,
- ~16 muls for point addition)
-
- The elliptic curve is over F_p, with p > 3 prime, and Weierstrass
- parameterization:
-
- y^2 = x^3 + a x + b
-
- The projective-format coordinates are actually stored in
- the form {X, Y, Z}, with true x,y
- coordinates on the curve given by {x,y} = {X/Z^2, Y/Z^3}.
- The function normalizeProj() performs the
- transformation from projective->true.
- (The other direction is trivial, i.e. {x,y} -> {x,y,1} will do.)
- The basic point multiplication function is
-
- ellMulProj()
-
- which obtains the result k * P for given point P and integer
- multiplier k. If true {x,y} are required for a multiple, one
- passes a point P = {X, Y, 1} to ellMulProj(), then afterwards
- calls normalizeProj(),
-
- Projective format is an answer to the typical sluggishness of
- standard elliptic arithmetic, whose explicit inversion in the
- field is, depending of course on the machinery and programmer,
- costly. Projective format is thus especially interesting for
- cryptography.
-
- REFERENCES
-
- perspective," Springer-Verlag, manuscript
-
- Solinas J 1998, IEEE P1363 Annex A (draft standard)
-
-***********************************************************/
-
-#include "ckconfig.h"
-#if CRYPTKIT_ELL_PROJ_ENABLE
-
-#include "ellipticProj.h"
-#include "falloc.h"
-#include "curveParams.h"
-#include "elliptic.h"
-#include "feeDebug.h"
-
-/*
- * convert REC-style smulg to generic imulg
- */
-#define smulg(s, g) imulg((unsigned)s, g)
-
-pointProj newPointProj(unsigned numDigits)
-{
- pointProj pt;
-
- pt = (pointProj) fmalloc(sizeof(pointProjStruct));
- pt->x = newGiant(numDigits);
- pt->y = newGiant(numDigits);
- pt->z = newGiant(numDigits);
- return(pt);
-}
-
-void freePointProj(pointProj pt)
-{
- clearGiant(pt->x); freeGiant(pt->x);
- clearGiant(pt->y); freeGiant(pt->y);
- clearGiant(pt->z); freeGiant(pt->z);
- ffree(pt);
-}
-
-void ptopProj(pointProj pt1, pointProj pt2)
-{
- gtog(pt1->x, pt2->x);
- gtog(pt1->y, pt2->y);
- gtog(pt1->z, pt2->z);
-}
-
-/**************************************************************
- *
- * Elliptic curve operations
- *
- **************************************************************/
-
-/* Begin projective-format functions for
-
- y^2 = x^3 + a x + b.
-
- These are useful in elliptic curve cryptography (ECC).
- A point is kept as a triple {X,Y,Z}, with the true (x,y)
- coordinates given by
-
- {x,y} = {X/Z^2, Y/Z^3}
-
- The function normalizeProj() performs the inverse conversion to get
- the true (x,y) pair.
- */
-
-void ellDoubleProj(pointProj pt, curveParams *cp)
-/* pt := 2 pt on the curve. */
-{
- giant x = pt->x, y = pt->y, z = pt->z;
- giant t1;
- giant t2;
- giant t3;
-
- if(isZero(y) || isZero(z)) {
- int_to_giant(1,x); int_to_giant(1,y); int_to_giant(0,z);
- return;
- }
- t1 = borrowGiant(cp->maxDigits);
- t2 = borrowGiant(cp->maxDigits);
- t3 = borrowGiant(cp->maxDigits);
-
- if((cp->a->sign >= 0) || (cp->a->n[0] != 3)) { /* Path prior to Apr2001. */
- gtog(z,t1); gsquare(t1); feemod(cp, t1);
- gsquare(t1); feemod(cp, t1);
- mulg(cp->a, t1); feemod(cp, t1); /* t1 := a z^4. */
- gtog(x, t2); gsquare(t2); feemod(cp, t2);
- smulg(3, t2); /* t2 := 3x^2. */
- addg(t2, t1); feemod(cp, t1); /* t1 := slope m. */
- } else { /* New optimization for a = -3 (post Apr 2001). */
- gtog(z, t1); gsquare(t1); feemod(cp, t1); /* t1 := z^2. */
- gtog(x, t2); subg(t1, t2); /* t2 := x-z^2. */
- addg(x, t1); smulg(3, t1); /* t1 := 3(x+z^2). */
- mulg(t2, t1); feemod(cp, t1); /* t1 := slope m. */
- }
- mulg(y, z); addg(z,z); feemod(cp, z); /* z := 2 y z. */
- gtog(y, t2); gsquare(t2); feemod(cp, t2); /* t2 := y^2. */
- gtog(t2, t3); gsquare(t3); feemod(cp, t3); /* t3 := y^4. */
- gshiftleft(3, t3); /* t3 := 8 y^4. */
- mulg(x, t2); gshiftleft(2, t2); feemod(cp, t2); /* t2 := 4xy^2. */
- gtog(t1, x); gsquare(x); feemod(cp, x);
- subg(t2, x); subg(t2, x); feemod(cp, x); /* x done. */
- gtog(t1, y); subg(x, t2); mulg(t2, y); subg(t3, y);
- feemod(cp, y);
- returnGiant(t1);
- returnGiant(t2);
- returnGiant(t3);
-}
-
-void ellAddProj(pointProj pt0, pointProj pt1, curveParams *cp)
-/* pt0 := pt0 + pt1 on the curve. */
-{
- giant x0 = pt0->x, y0 = pt0->y, z0 = pt0->z,
- x1 = pt1->x, y1 = pt1->y, z1 = pt1->z;
- giant t1;
- giant t2;
- giant t3;
- giant t4;
- giant t5;
- giant t6;
- giant t7;
-
- if(isZero(z0)) {
- gtog(x1,x0); gtog(y1,y0); gtog(z1,z0);
- return;
- }
- if(isZero(z1)) return;
-
- t1 = borrowGiant(cp->maxDigits);
- t2 = borrowGiant(cp->maxDigits);
- t3 = borrowGiant(cp->maxDigits);
- t4 = borrowGiant(cp->maxDigits);
- t5 = borrowGiant(cp->maxDigits);
- t6 = borrowGiant(cp->maxDigits);
- t7 = borrowGiant(cp->maxDigits);
-
- gtog(x0, t1); gtog(y0,t2); gtog(z0, t3);
- gtog(x1, t4); gtog(y1, t5);
- if(!isone(z1)) {
- gtog(z1, t6);
- gtog(t6, t7); gsquare(t7); feemod(cp, t7);
- mulg(t7, t1); feemod(cp, t1);
- mulg(t6, t7); feemod(cp, t7);
- mulg(t7, t2); feemod(cp, t2);
- }
- gtog(t3, t7); gsquare(t7); feemod(cp, t7);
- mulg(t7, t4); feemod(cp, t4);
- mulg(t3, t7); feemod(cp, t7);
- mulg(t7, t5); feemod(cp, t5);
- negg(t4); addg(t1, t4); feemod(cp, t4);
- negg(t5); addg(t2, t5); feemod(cp, t5);
- if(isZero(t4)) {
- if(isZero(t5)) {
- ellDoubleProj(pt0, cp);
- } else {
- int_to_giant(1, x0); int_to_giant(1, y0);
- int_to_giant(0, z0);
- }
- goto out;
- }
- addg(t1, t1); subg(t4, t1); feemod(cp, t1);
- addg(t2, t2); subg(t5, t2); feemod(cp, t2);
- if(!isone(z1)) {
- mulg(t6, t3); feemod(cp, t3);
- }
- mulg(t4, t3); feemod(cp, t3);
- gtog(t4, t7); gsquare(t7); feemod(cp, t7);
- mulg(t7, t4); feemod(cp, t4);
- mulg(t1, t7); feemod(cp, t7);
- gtog(t5, t1); gsquare(t1); feemod(cp, t1);
- subg(t7, t1); feemod(cp, t1);
- subg(t1, t7); subg(t1, t7); feemod(cp, t7);
- mulg(t7, t5); feemod(cp, t5);
- mulg(t2, t4); feemod(cp, t4);
- gtog(t5, t2); subg(t4,t2); feemod(cp, t2);
- if(t2->n[0] & 1) { /* Test if t2 is odd. */
- addg(cp->basePrime, t2);
- }
- gshiftright(1, t2);
- gtog(t1, x0); gtog(t2, y0); gtog(t3, z0);
-out:
- returnGiant(t1);
- returnGiant(t2);
- returnGiant(t3);
- returnGiant(t4);
- returnGiant(t5);
- returnGiant(t6);
- returnGiant(t7);
-}
-
-
-void ellNegProj(pointProj pt, curveParams *cp)
-/* pt := -pt on the curve. */
-{
- negg(pt->y); feemod(cp, pt->y);
-}
-
-void ellSubProj(pointProj pt0, pointProj pt1, curveParams *cp)
-/* pt0 := pt0 - pt1 on the curve. */
-{
- ellNegProj(pt1, cp);
- ellAddProj(pt0, pt1,cp);
- ellNegProj(pt1, cp);
-}
-
-/*
- * Simple projective multiply.
- *
- * pt := pt * k, result normalized.
- */
-void ellMulProjSimple(pointProj pt0, giant k, curveParams *cp)
-{
- pointProjStruct pt1; // local, giants borrowed
-
- CKASSERT(isone(pt0->z));
- CKASSERT(cp->curveType == FCT_Weierstrass);
-
- /* ellMulProj assumes constant pt0, can't pass as src and dst */
- pt1.x = borrowGiant(cp->maxDigits);
- pt1.y = borrowGiant(cp->maxDigits);
- pt1.z = borrowGiant(cp->maxDigits);
- ellMulProj(pt0, &pt1, k, cp);
- normalizeProj(&pt1, cp);
- CKASSERT(isone(pt1.z));
-
- ptopProj(&pt1, pt0);
- returnGiant(pt1.x);
- returnGiant(pt1.y);
- returnGiant(pt1.z);
-}
-
-void ellMulProj(pointProj pt0, pointProj pt1, giant k, curveParams *cp)
-/* General elliptic multiplication;
- pt1 := k*pt0 on the curve,
- with k an arbitrary integer.
- */
-{
- giant x = pt0->x, y = pt0->y, z = pt0->z,
- xx = pt1->x, yy = pt1->y, zz = pt1->z;
- int ksign, hlen, klen, b, hb, kb;
- giant t0;
-
- CKASSERT(cp->curveType == FCT_Weierstrass);
- if(isZero(k)) {
- int_to_giant(1, xx);
- int_to_giant(1, yy);
- int_to_giant(0, zz);
- return;
- }
- t0 = borrowGiant(cp->maxDigits);
- ksign = k->sign;
- if(ksign < 0) negg(k);
- gtog(x,xx); gtog(y,yy); gtog(z,zz);
- gtog(k, t0); addg(t0, t0); addg(k, t0); /* t0 := 3k. */
- hlen = bitlen(t0);
- klen = bitlen(k);
- for(b = hlen-2; b > 0; b--) {
- ellDoubleProj(pt1,cp);
- hb = bitval(t0, b);
- if(b < klen) kb = bitval(k, b); else kb = 0;
- if((hb != 0) && (kb == 0))
- ellAddProj(pt1, pt0, cp);
- else if((hb == 0) && (kb !=0))
- ellSubProj(pt1, pt0, cp);
- }
- if(ksign < 0) {
- ellNegProj(pt1, cp);
- k->sign = -k->sign;
- }
- returnGiant(t0);
-}
-
-void normalizeProj(pointProj pt, curveParams *cp)
-/* Obtain actual x,y coords via normalization:
- {x,y,z} := {x/z^2, y/z^3, 1}.
- */
-
-{ giant x = pt->x, y = pt->y, z = pt->z;
- giant t1;
-
- CKASSERT(cp->curveType == FCT_Weierstrass);
- if(isZero(z)) {
- int_to_giant(1,x); int_to_giant(1,y);
- return;
- }
- t1 = borrowGiant(cp->maxDigits);
- binvg_cp(cp, z); // was binvaux(p, z);
- gtog(z, t1);
- gsquare(z); feemod(cp, z);
- mulg(z, x); feemod(cp, x);
- mulg(t1, z); mulg(z, y); feemod(cp, y);
- int_to_giant(1, z);
- returnGiant(t1);
-}
-
-static int
-jacobi_symbol(giant a, curveParams *cp)
-/* Standard Jacobi symbol (a/cp->basePrime).
- basePrime must be odd, positive. */
-{
- int t = 1, u;
- giant t5 = borrowGiant(cp->maxDigits);
- giant t6 = borrowGiant(cp->maxDigits);
- giant t7 = borrowGiant(cp->maxDigits);
- int rtn;
-
- gtog(a, t5); feemod(cp, t5);
- gtog(cp->basePrime, t6);
- while(!isZero(t5)) {
- u = (t6->n[0]) & 7;
- while((t5->n[0] & 1) == 0) {
- gshiftright(1, t5);
- if((u==3) || (u==5)) t = -t;
- }
- gtog(t5, t7); gtog(t6, t5); gtog(t7, t6);
- u = (t6->n[0]) & 3;
- if(((t5->n[0] & 3) == 3) && ((u & 3) == 3)) t = -t;
- modg(t6, t5);
- }
- if(isone(t6)) {
- rtn = t;
- }
- else {
- rtn = 0;
- }
- returnGiant(t5);
- returnGiant(t6);
- returnGiant(t7);
-
- return rtn;
-}
-
-static void
-powFp2(giant a, giant b, giant w2, giant n, curveParams *cp)
-/* Perform powering in the field F_p^2:
- a + b w := (a + b w)^n (mod p), where parameter w2 is a quadratic
- nonresidue (formally equal to w^2).
- */
-{
- int j;
- giant t6;
- giant t7;
- giant t8;
- giant t9;
-
- if(isZero(n)) {
- int_to_giant(1,a);
- int_to_giant(0,b);
- return;
- }
- t6 = borrowGiant(cp->maxDigits);
- t7 = borrowGiant(cp->maxDigits);
- t8 = borrowGiant(cp->maxDigits);
- t9 = borrowGiant(cp->maxDigits);
- gtog(a, t8); gtog(b, t9);
- for(j = bitlen(n)-2; j >= 0; j--) {
- gtog(b, t6);
- mulg(a, b); addg(b,b); feemod(cp, b); /* b := 2 a b. */
- gsquare(t6); feemod(cp, t6);
- mulg(w2, t6); feemod(cp, t6);
- gsquare(a); addg(t6, a); feemod(cp, a);
- /* a := a^2 + b^2 w2. */
- if(bitval(n, j)) {
- gtog(b, t6); mulg(t8, b); feemod(cp, b);
- gtog(a, t7); mulg(t9, a); addg(a, b); feemod(cp, b);
- mulg(t9, t6); feemod(cp, t6);
- mulg(w2, t6); feemod(cp, t6);
- mulg(t8, a); addg(t6, a); feemod(cp, a);
- }
- }
- returnGiant(t6);
- returnGiant(t7);
- returnGiant(t8);
- returnGiant(t9);
- return;
-}
-
-static void
-powermodg(
- giant x,
- giant n,
- curveParams *cp
-)
-/* x becomes x^n (mod basePrime). */
-{
- int len, pos;
- giant scratch2 = borrowGiant(cp->maxDigits);
-
- gtog(x, scratch2);
- int_to_giant(1, x);
- len = bitlen(n);
- pos = 0;
- while (1)
- {
- if (bitval(n, pos++))
- {
- mulg(scratch2, x);
- feemod(cp, x);
- }
- if (pos>=len)
- break;
- gsquare(scratch2);
- feemod(cp, scratch2);
- }
- returnGiant(scratch2);
-}
-
-static int sqrtmod(giant x, curveParams *cp)
-/* If Sqrt[x] (mod p) exists, function returns 1, else 0.
- In either case x is modified, but if 1 is returned,
- x:= Sqrt[x] (mod p).
- */
-{
- int rtn;
- giant t0 = borrowGiant(cp->maxDigits);
- giant t1 = borrowGiant(cp->maxDigits);
- giant t2 = borrowGiant(cp->maxDigits);
- giant t3 = borrowGiant(cp->maxDigits);
- giant t4 = borrowGiant(cp->maxDigits);
-
- giant p = cp->basePrime;
-
- feemod(cp, x); /* Justify the argument. */
- gtog(x, t0); /* Store x for eventual validity check on square root. */
- if((p->n[0] & 3) == 3) { /* The case p = 3 (mod 4). */
- gtog(p, t1);
- iaddg(1, t1); gshiftright(2, t1);
- powermodg(x, t1, cp);
- goto resolve;
- }
- /* Next, handle case p = 5 (mod 8). */
- if((p->n[0] & 7) == 5) {
- gtog(p, t1); int_to_giant(1, t2);
- subg(t2, t1); gshiftright(2, t1);
- gtog(x, t2);
- powermodg(t2, t1, cp); /* t2 := x^((p-1)/4) % p. */
- iaddg(1, t1);
- gshiftright(1, t1); /* t1 := (p+3)/8. */
- if(isone(t2)) {
- powermodg(x, t1, cp); /* x^((p+3)/8) is root. */
- goto resolve;
- } else {
- int_to_giant(1, t2); subg(t2, t1);
- /* t1 := (p-5)/8. */
- gshiftleft(2,x);
- powermodg(x, t1, cp);
- mulg(t0, x); addg(x, x); feemod(cp, x);
- /* 2x (4x)^((p-5)/8. */
- goto resolve;
- }
- }
-
- /* Next, handle tougher case: p = 1 (mod 8). */
- int_to_giant(2, t1);
- while(1) { /* Find appropriate nonresidue. */
- gtog(t1, t2);
- gsquare(t2); subg(x, t2); feemod(cp, t2);
- if(jacobi_symbol(t2, cp) == -1) break;
- iaddg(1, t1);
- } /* t2 is now w^2 in F_p^2. */
- int_to_giant(1, t3);
- gtog(p, t4); iaddg(1, t4); gshiftright(1, t4);
- powFp2(t1, t3, t2, t4, cp);
- gtog(t1, x);
-
-resolve:
- gtog(x,t1); gsquare(t1); feemod(cp, t1);
- if(gcompg(t0, t1) == 0) {
- rtn = 1; /* Success. */
- }
- else {
- rtn = 0; /* no square root */
- }
- returnGiant(t0);
- returnGiant(t1);
- returnGiant(t2);
- returnGiant(t3);
- returnGiant(t4);
- return rtn;
-}
-
-
-void findPointProj(pointProj pt, giant seed, curveParams *cp)
-/* Starting with seed, finds a random (projective) point {x,y,1} on curve.
- */
-{
- giant x = pt->x, y = pt->y, z = pt->z;
-
- CKASSERT(cp->curveType == FCT_Weierstrass);
- feemod(cp, seed);
- while(1) {
- gtog(seed, x);
- gsquare(x); feemod(cp, x); // x := seed^2
- addg(cp->a, x); // x := seed^2 + a
- mulg(seed,x); // x := seed^3 + a*seed
- addg(cp->b, x);
- feemod(cp, x); // x := seed^3 + a seed + b.
- /* test cubic form for having root. */
- if(sqrtmod(x, cp)) break;
- iaddg(1, seed);
- }
- gtog(x, y);
- gtog(seed,x);
- int_to_giant(1, z);
-}
-
-#endif /* CRYPTKIT_ELL_PROJ_ENABLE */