X-Git-Url: https://git.saurik.com/apple/security.git/blobdiff_plain/80e2389990082500d76eb566d4946be3e786c3ef..d8f41ccd20de16f8ebe2ccc84d47bf1cb2b26bbb:/Security/libsecurity_cryptkit/lib/CurveParamDocs/ellproj.c diff --git a/Security/libsecurity_cryptkit/lib/CurveParamDocs/ellproj.c b/Security/libsecurity_cryptkit/lib/CurveParamDocs/ellproj.c new file mode 100644 index 00000000..b98b2573 --- /dev/null +++ b/Security/libsecurity_cryptkit/lib/CurveParamDocs/ellproj.c @@ -0,0 +1,448 @@ +/************************************************************** + * + * ellproj.c + * + Fast algorithms for fundamental elliptic curve arithmetic, + projective format. Such algorithms apply in domains such as: + -- factoring + -- primality studies (e.g. rigorous primality proofs) + -- elliptic curve cryptography (ECC) + + 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 normalize_proj() 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 + + ell_mul_proj() + + 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 ell_mul_proj(), then afterwards + calls normalize_proj(), + + 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 + + Crandall R and Pomerance C 1998, "Prime numbers: a computational + perspective," Springer-Verlag, manuscript + + Solinas J 1998, IEEE P1363 Annex A (draft standard) + + LEGAL AND PATENT NOTICE + + This and related PSI library source code is intended solely for + educational and research applications, and should not be used + for commercial purposes without explicit permission from PSI + (not to mention proper clearance of legal restrictions beyond + the purview of PSI). + The code herein will not perform cryptography per se, + although the number-theoretical algorithms herein -- all of which + are in the public domain -- can be used in principle to effect + what is known as elliptic curve cryptography (ECC). Note that + there are strict constraints on how cryptography may be used, + especially in regard to exportability. + Therefore one should avoid any casual distribution of actual + cryptographic software. Along similar lines, because of various + patents, proprietary to Apple Computer, Inc., and perhaps to other + organizations, one should not tacitly assume that an ECC scheme is + unconstrained. For example,the commercial use of certain fields + F_p^k (i.e., fixation of certain primes p) is covered in Apple + patents. + + * Updates: + * 3 Apr 98 REC Creation + * + * c. 1998 Perfectly Scientific, Inc. + * All Rights Reserved. + * + * + *************************************************************/ + +/* include files */ + +#include +#include +#include +#include +#ifdef _WIN32 + +#include + +#endif + +#include +#include "giants.h" +#include "ellproj.h" +#include "tools.h" + +/* global variables */ + +static giant t0 = NULL, t1 = NULL, t2 = NULL, t3 = NULL, t4 = NULL, + t5 = NULL, t6 = NULL, t7 = NULL; + +/************************************************************** + * + * Maintenance functions + * + **************************************************************/ + +point_proj +new_point_proj(int shorts) +{ + point_proj pt; + + if(t0 == NULL) init_ell_proj(shorts); + pt = (point_proj) malloc(sizeof(point_struct_proj)); + pt->x = newgiant(shorts); + pt->y = newgiant(shorts); + pt->z = newgiant(shorts); + return(pt); +} + +void +free_point_proj(point_proj pt) +{ + free(pt->x); free(pt->y); free(pt->z); + free(pt); +} + +void +ptop_proj(point_proj pt1, point_proj pt2) +{ + gtog(pt1->x, pt2->x); + gtog(pt1->y, pt2->y); + gtog(pt1->z, pt2->z); +} + +void +init_ell_proj(int shorts) +/* Called by new_point_proj(), to set up giant registers. */ +{ + t0 = newgiant(shorts); + t1 = newgiant(shorts); + t2 = newgiant(shorts); + t3 = newgiant(shorts); + t4 = newgiant(shorts); + t5 = newgiant(shorts); + t6 = newgiant(shorts); + t7 = newgiant(shorts); +} + + +/************************************************************** + * + * 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 normalize_proj() performs the inverse conversion to get + the true (x,y) pair. + */ + +void +ell_double_proj(point_proj pt, giant a, giant p) +/* pt := 2 pt on the curve. */ +{ + giant x = pt->x, y = pt->y, z = pt->z; + + if(isZero(y) || isZero(z)) { + itog(1,x); itog(1,y); itog(0,z); + return; + } + gtog(z,t1); squareg(t1); modg(p, t1); + squareg(t1); modg(p, t1); + mulg(a, t1); modg(p, t1); /* t1 := a z^4. */ + gtog(x, t2); squareg(t2); smulg(3, t2); modg(p, t2); /* t2 := 3x^2. */ + addg(t2, t1); modg(p, t1); /* t1 := slope m. */ + mulg(y, z); addg(z,z); modg(p, z); /* z := 2 y z. */ + gtog(y, t2); squareg(t2); modg(p, t2); /* t2 := y^2. */ + gtog(t2, t3); squareg(t3); modg(p, t3); /* t3 := y^4. */ + gshiftleft(3, t3); /* t3 := 8 y^4. */ + mulg(x, t2); gshiftleft(2, t2); modg(p, t2); /* t2 := 4xy^2. */ + gtog(t1, x); squareg(x); modg(p, x); + subg(t2, x); subg(t2, x); modg(p, x); /* x done. */ + gtog(t1, y); subg(x, t2); mulg(t2, y); subg(t3, y); + modg(p, y); +} +/* +elldouble[pt_] := Block[{x,y,z,m,y2,s}, + x = pt[[1]]; y = pt[[2]]; z = pt[[3]]; + If[(y==0) || (z==0), Return[{1,1,0}]]; + m = Mod[3 x^2 + a Mod[Mod[z^2,p]^2,p],p]; + z = Mod[2 y z, p]; + y2 = Mod[y^2,p]; + s = Mod[4 x y2,p]; + x = Mod[m^2 - 2s,p]; + y = Mod[m(s - x) - 8 y2^2,p]; + Return[{x,y,z}]; +]; +*/ + +void +ell_add_proj(point_proj pt0, point_proj pt1, giant a, giant p) +/* 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; + + if(isZero(z0)) { + gtog(x1,x0); gtog(y1,y0); gtog(z1,z0); + return; + } + if(isZero(z1)) return; + gtog(x0, t1); gtog(y0,t2); gtog(z0, t3); + gtog(x1, t4); gtog(y1, t5); + if(!isone(z1)) { + gtog(z1, t6); + gtog(t6, t7); squareg(t7); modg(p, t7); + mulg(t7, t1); modg(p, t1); + mulg(t6, t7); modg(p, t7); + mulg(t7, t2); modg(p, t2); + } + gtog(t3, t7); squareg(t7); modg(p, t7); + mulg(t7, t4); modg(p, t4); + mulg(t3, t7); modg(p, t7); + mulg(t7, t5); modg(p, t5); + negg(t4); addg(t1, t4); modg(p, t4); + negg(t5); addg(t2, t5); modg(p, t5); + if(isZero(t4)) { + if(isZero(t5)) { + ell_double_proj(pt0, a, p); + } else { + itog(1, x0); itog(1, y0); itog(0, z0); + } + return; + } + addg(t1, t1); subg(t4, t1); modg(p, t1); + addg(t2, t2); subg(t5, t2); modg(p, t2); + if(!isone(z1)) { + mulg(t6, t3); modg(p, t3); + } + mulg(t4, t3); modg(p, t3); + gtog(t4, t7); squareg(t7); modg(p, t7); + mulg(t7, t4); modg(p, t4); + mulg(t1, t7); modg(p, t7); + gtog(t5, t1); squareg(t1); modg(p, t1); + subg(t7, t1); modg(p, t1); + subg(t1, t7); subg(t1, t7); modg(p, t7); + mulg(t7, t5); modg(p, t5); + mulg(t2, t4); modg(p, t4); + gtog(t5, t2); subg(t4,t2); modg(p, t2); + if(t2->n[0] & 1) { /* Test if t2 is odd. */ + addg(p, t2); + } + gshiftright(1, t2); + gtog(t1, x0); gtog(t2, y0); gtog(t3, z0); +} + +/* +elladd[pt0_, pt1_] := Block[ + {x0,y0,z0,x1,y1,z1, + t1,t2,t3,t4,t5,t6,t7}, + x0 = pt0[[1]]; y0 = pt0[[2]]; z0 = pt0[[3]]; + x1 = pt1[[1]]; y1 = pt1[[2]]; z1 = pt1[[3]]; + If[z0 == 0, Return[pt1]]; + If[z1 == 0, Return[pt0]]; + + t1 = x0; + t2 = y0; + t3 = z0; + t4 = x1; + t5 = y1; + If[(z1 != 1), + t6 = z1; + t7 = Mod[t6^2, p]; + t1 = Mod[t1 t7, p]; + t7 = Mod[t6 t7, p]; + t2 = Mod[t2 t7, p]; + ]; + t7 = Mod[t3^2, p]; + t4 = Mod[t4 t7, p]; + t7 = Mod[t3 t7, p]; + t5 = Mod[t5 t7, p]; + t4 = Mod[t1-t4, p]; + t5 = Mod[t2 - t5, p]; + If[t4 == 0, If[t5 == 0, + Return[elldouble[pt0]], + Return[{1,1,0}] + ] + ]; + t1 = Mod[2t1 - t4,p]; + t2 = Mod[2t2 - t5, p]; + If[z1 != 1, t3 = Mod[t3 t6, p]]; + t3 = Mod[t3 t4, p]; + t7 = Mod[t4^2, p]; + t4 = Mod[t4 t7, p]; + t7 = Mod[t1 t7, p]; + t1 = Mod[t5^2, p]; + t1 = Mod[t1-t7, p]; + t7 = Mod[t7 - 2t1, p]; + t5 = Mod[t5 t7, p]; + t4 = Mod[t2 t4, p]; + t2 = Mod[t5-t4, p]; + If[EvenQ[t2], t2 = t2/2, t2 = (p+t2)/2]; + Return[{t1, t2, t3}]; +]; +*/ + +void +ell_neg_proj(point_proj pt, giant p) +/* pt := -pt on the curve. */ +{ + negg(pt->y); modg(p, pt->y); +} + +void +ell_sub_proj(point_proj pt0, point_proj pt1, giant a, giant p) +/* pt0 := pt0 - pt1 on the curve. */ +{ + ell_neg_proj(pt1, p); + ell_add_proj(pt0, pt1,a,p); + ell_neg_proj(pt1,p); +} + +void +ell_mul_proj(point_proj pt0, point_proj pt1, giant k, giant a, giant p) +/* 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; + + if(isZero(k)) { + itog(1, xx); + itog(1, yy); + itog(0, zz); + return; + } + 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--) { + ell_double_proj(pt1,a,p); + hb = bitval(t0, b); + if(b < klen) kb = bitval(k, b); else kb = 0; + if((hb != 0) && (kb == 0)) + ell_add_proj(pt1, pt0, a, p); + else if((hb == 0) && (kb !=0)) + ell_sub_proj(pt1, pt0, a, p); + } + if(ksign < 0) { + ell_neg_proj(pt1, p); + k->sign = -k->sign; + } +} + +/* +elliptic[pt_, k_] := Block[{pt2, hh, kk, hb, kb, lenh, lenk}, + If[k==0, Return[{1,1,0}]]; + hh = Reverse[bitList[3k]]; + kk = Reverse[bitList[k]]; + pt2 = pt; + lenh = Length[hh]; + lenk = Length[kk]; + Do[ + pt2 = elldouble[pt2]; + hb = hh[[b]]; + If[b <= lenk, kb = kk[[b]], kb = 0]; + If[{hb,kb} == {1,0}, + pt2 = elladd[pt2, pt], + If[{hb, kb} == {0,1}, + pt2 = ellsub[pt2, pt]] + ] + ,{b, lenh-1, 2,-1} + ]; + Return[pt2]; +]; +*/ + +void +normalize_proj(point_proj pt, giant p) +/* 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; + + if(isZero(z)) { + itog(1,x); itog(1,y); + return; + } + binvaux(p, z); gtog(z, t1); + squareg(z); modg(p, z); + mulg(z, x); modg(p, x); + mulg(t1, z); mulg(z, y); modg(p, y); + itog(1, z); +} + +/* +normalize[pt_] := Block[{z,z2,z3}, + If[pt[[3]] == 0, Return[pt]]; + z = ellinv[pt[[3]]]; + z2 = Mod[z^2,p]; + z3 = Mod[z z2,p]; + Return[{Mod[pt[[1]] z2, p], Mod[pt[[2]] z3, p], 1}]; + ]; +*/ + + +void +find_point_proj(point_proj pt, giant seed, giant a, giant b, giant p) +/* Starting with seed, finds a random (projective) point {x,y,1} on curve. + */ +{ giant x = pt->x, y = pt->y, z = pt->z; + + modg(p, seed); + while(1) { + gtog(seed, x); + squareg(x); modg(p, x); + addg(a, x); + mulg(seed,x); addg(b, x); + modg(p, x); /* x := seed^3 + a seed + b. */ + if(sqrtmod(p, x)) break; /* Test if cubic form has root. */ + iaddg(1, seed); + } + gtog(x, y); + gtog(seed,x); + itog(1, z); +}