The Tor Browser / file comparison

comparison: security/nss/lib/freebl/ecl/ecp_fpinc.c

security/nss/lib/freebl/ecl/ecp_fpinc.c

branch
TOR_BUG_9701
changeset 15
b8a032363ba2
equal deleted inserted replaced
-1:000000000000 0:e72f209ee978
1 /* This Source Code Form is subject to the terms of the Mozilla Public
2 * License, v. 2.0. If a copy of the MPL was not distributed with this
3 * file, You can obtain one at http://mozilla.org/MPL/2.0/. */
4
5 /* This source file is meant to be included by other source files
6 * (ecp_fp###.c, where ### is one of 160, 192, 224) and should not
7 * constitute an independent compilation unit. It requires the following
8 * preprocessor definitions be made: ECFP_BSIZE - the number of bits in
9 * the field's prime
10 * ECFP_NUMDOUBLES - the number of doubles to store one
11 * multi-precision integer in floating point
12
13 /* Adds a prefix to a given token to give a unique token name. Prefixes
14 * with "ecfp" + ECFP_BSIZE + "_". e.g. if ECFP_BSIZE = 160, then
15 * PREFIX(hello) = ecfp160_hello This optimization allows static function
16 * linking and compiler loop unrolling without code duplication. */
17 #ifndef PREFIX
18 #define PREFIX(b) PREFIX1(ECFP_BSIZE, b)
19 #define PREFIX1(bsize, b) PREFIX2(bsize, b)
20 #define PREFIX2(bsize, b) ecfp ## bsize ## _ ## b
21 #endif
22
23 /* Returns true iff every double in d is 0. (If d == 0 and it is tidied,
24 * this will be true.) */
25 mp_err PREFIX(isZero) (const double *d) {
26 int i;
27
28 for (i = 0; i < ECFP_NUMDOUBLES; i++) {
29 if (d[i] != 0)
30 return MP_NO;
31 }
32 return MP_YES;
33 }
34
35 /* Sets the multi-precision floating point number at t = 0 */
36 void PREFIX(zero) (double *t) {
37 int i;
38
39 for (i = 0; i < ECFP_NUMDOUBLES; i++) {
40 t[i] = 0;
41 }
42 }
43
44 /* Sets the multi-precision floating point number at t = 1 */
45 void PREFIX(one) (double *t) {
46 int i;
47
48 t[0] = 1;
49 for (i = 1; i < ECFP_NUMDOUBLES; i++) {
50 t[i] = 0;
51 }
52 }
53
54 /* Checks if point P(x, y, z) is at infinity. Uses Jacobian coordinates. */
55 mp_err PREFIX(pt_is_inf_jac) (const ecfp_jac_pt * p) {
56 return PREFIX(isZero) (p->z);
57 }
58
59 /* Sets the Jacobian point P to be at infinity. */
60 void PREFIX(set_pt_inf_jac) (ecfp_jac_pt * p) {
61 PREFIX(zero) (p->z);
62 }
63
64 /* Checks if point P(x, y) is at infinity. Uses Affine coordinates. */
65 mp_err PREFIX(pt_is_inf_aff) (const ecfp_aff_pt * p) {
66 if (PREFIX(isZero) (p->x) == MP_YES && PREFIX(isZero) (p->y) == MP_YES)
67 return MP_YES;
68 return MP_NO;
69 }
70
71 /* Sets the affine point P to be at infinity. */
72 void PREFIX(set_pt_inf_aff) (ecfp_aff_pt * p) {
73 PREFIX(zero) (p->x);
74 PREFIX(zero) (p->y);
75 }
76
77 /* Checks if point P(x, y, z, a*z^4) is at infinity. Uses Modified
78 * Jacobian coordinates. */
79 mp_err PREFIX(pt_is_inf_jm) (const ecfp_jm_pt * p) {
80 return PREFIX(isZero) (p->z);
81 }
82
83 /* Sets the Modified Jacobian point P to be at infinity. */
84 void PREFIX(set_pt_inf_jm) (ecfp_jm_pt * p) {
85 PREFIX(zero) (p->z);
86 }
87
88 /* Checks if point P(x, y, z, z^2, z^3) is at infinity. Uses Chudnovsky
89 * Jacobian coordinates */
90 mp_err PREFIX(pt_is_inf_chud) (const ecfp_chud_pt * p) {
91 return PREFIX(isZero) (p->z);
92 }
93
94 /* Sets the Chudnovsky Jacobian point P to be at infinity. */
95 void PREFIX(set_pt_inf_chud) (ecfp_chud_pt * p) {
96 PREFIX(zero) (p->z);
97 }
98
99 /* Copies a multi-precision floating point number, Setting dest = src */
100 void PREFIX(copy) (double *dest, const double *src) {
101 int i;
102
103 for (i = 0; i < ECFP_NUMDOUBLES; i++) {
104 dest[i] = src[i];
105 }
106 }
107
108 /* Sets dest = -src */
109 void PREFIX(negLong) (double *dest, const double *src) {
110 int i;
111
112 for (i = 0; i < 2 * ECFP_NUMDOUBLES; i++) {
113 dest[i] = -src[i];
114 }
115 }
116
117 /* Sets r = -p p = (x, y, z, z2, z3) r = (x, -y, z, z2, z3) Uses
118 * Chudnovsky Jacobian coordinates. */
119 /* TODO reverse order */
120 void PREFIX(pt_neg_chud) (const ecfp_chud_pt * p, ecfp_chud_pt * r) {
121 int i;
122
123 PREFIX(copy) (r->x, p->x);
124 PREFIX(copy) (r->z, p->z);
125 PREFIX(copy) (r->z2, p->z2);
126 PREFIX(copy) (r->z3, p->z3);
127 for (i = 0; i < ECFP_NUMDOUBLES; i++) {
128 r->y[i] = -p->y[i];
129 }
130 }
131
132 /* Computes r = x + y. Does not tidy or reduce. Any combinations of r, x,
133 * y can point to the same data. Componentwise adds first ECFP_NUMDOUBLES
134 * doubles of x and y and stores the result in r. */
135 void PREFIX(addShort) (double *r, const double *x, const double *y) {
136 int i;
137
138 for (i = 0; i < ECFP_NUMDOUBLES; i++) {
139 *r++ = *x++ + *y++;
140 }
141 }
142
143 /* Computes r = x + y. Does not tidy or reduce. Any combinations of r, x,
144 * y can point to the same data. Componentwise adds first
145 * 2*ECFP_NUMDOUBLES doubles of x and y and stores the result in r. */
146 void PREFIX(addLong) (double *r, const double *x, const double *y) {
147 int i;
148
149 for (i = 0; i < 2 * ECFP_NUMDOUBLES; i++) {
150 *r++ = *x++ + *y++;
151 }
152 }
153
154 /* Computes r = x - y. Does not tidy or reduce. Any combinations of r, x,
155 * y can point to the same data. Componentwise subtracts first
156 * ECFP_NUMDOUBLES doubles of x and y and stores the result in r. */
157 void PREFIX(subtractShort) (double *r, const double *x, const double *y) {
158 int i;
159
160 for (i = 0; i < ECFP_NUMDOUBLES; i++) {
161 *r++ = *x++ - *y++;
162 }
163 }
164
165 /* Computes r = x - y. Does not tidy or reduce. Any combinations of r, x,
166 * y can point to the same data. Componentwise subtracts first
167 * 2*ECFP_NUMDOUBLES doubles of x and y and stores the result in r. */
168 void PREFIX(subtractLong) (double *r, const double *x, const double *y) {
169 int i;
170
171 for (i = 0; i < 2 * ECFP_NUMDOUBLES; i++) {
172 *r++ = *x++ - *y++;
173 }
174 }
175
176 /* Computes r = x*y. Both x and y should be tidied and reduced,
177 * r must be different (point to different memory) than x and y.
178 * Does not tidy or reduce. */
179 void PREFIX(multiply)(double *r, const double *x, const double *y) {
180 int i, j;
181
182 for(j=0;j<ECFP_NUMDOUBLES-1;j++) {
183 r[j] = x[0] * y[j];
184 r[j+(ECFP_NUMDOUBLES-1)] = x[ECFP_NUMDOUBLES-1] * y[j];
185 }
186 r[ECFP_NUMDOUBLES-1] = x[0] * y[ECFP_NUMDOUBLES-1];
187 r[ECFP_NUMDOUBLES-1] += x[ECFP_NUMDOUBLES-1] * y[0];
188 r[2*ECFP_NUMDOUBLES-2] = x[ECFP_NUMDOUBLES-1] * y[ECFP_NUMDOUBLES-1];
189 r[2*ECFP_NUMDOUBLES-1] = 0;
190
191 for(i=1;i<ECFP_NUMDOUBLES-1;i++) {
192 for(j=0;j<ECFP_NUMDOUBLES;j++) {
193 r[i+j] += (x[i] * y[j]);
194 }
195 }
196 }
197
198 /* Computes the square of x and stores the result in r. x should be
199 * tidied & reduced, r will be neither tidied nor reduced.
200 * r should point to different memory than x */
201 void PREFIX(square) (double *r, const double *x) {
202 PREFIX(multiply) (r, x, x);
203 }
204
205 /* Perform a point doubling in Jacobian coordinates. Input and output
206 * should be multi-precision floating point integers. */
207 void PREFIX(pt_dbl_jac) (const ecfp_jac_pt * dp, ecfp_jac_pt * dr,
208 const EC_group_fp * group) {
209 double t0[2 * ECFP_NUMDOUBLES], t1[2 * ECFP_NUMDOUBLES],
210 M[2 * ECFP_NUMDOUBLES], S[2 * ECFP_NUMDOUBLES];
211
212 /* Check for point at infinity */
213 if (PREFIX(pt_is_inf_jac) (dp) == MP_YES) {
214 /* Set r = pt at infinity */
215 PREFIX(set_pt_inf_jac) (dr);
216 goto CLEANUP;
217 }
218
219 /* Perform typical point doubling operations */
220
221 /* TODO? is it worthwhile to do optimizations for when pz = 1? */
222
223 if (group->aIsM3) {
224 /* When a = -3, M = 3(px - pz^2)(px + pz^2) */
225 PREFIX(square) (t1, dp->z);
226 group->ecfp_reduce(t1, t1, group); /* 2^23 since the negative
227 * rounding buys another bit */
228 PREFIX(addShort) (t0, dp->x, t1); /* 2*2^23 */
229 PREFIX(subtractShort) (t1, dp->x, t1); /* 2 * 2^23 */
230 PREFIX(multiply) (M, t0, t1); /* 40 * 2^46 */
231 PREFIX(addLong) (t0, M, M); /* 80 * 2^46 */
232 PREFIX(addLong) (M, t0, M); /* 120 * 2^46 < 2^53 */
233 group->ecfp_reduce(M, M, group);
234 } else {
235 /* Generic case */
236 /* M = 3 (px^2) + a*(pz^4) */
237 PREFIX(square) (t0, dp->x);
238 PREFIX(addLong) (M, t0, t0);
239 PREFIX(addLong) (t0, t0, M); /* t0 = 3(px^2) */
240 PREFIX(square) (M, dp->z);
241 group->ecfp_reduce(M, M, group);
242 PREFIX(square) (t1, M);
243 group->ecfp_reduce(t1, t1, group);
244 PREFIX(multiply) (M, t1, group->curvea); /* M = a(pz^4) */
245 PREFIX(addLong) (M, M, t0);
246 group->ecfp_reduce(M, M, group);
247 }
248
249 /* rz = 2 * py * pz */
250 PREFIX(multiply) (t1, dp->y, dp->z);
251 PREFIX(addLong) (t1, t1, t1);
252 group->ecfp_reduce(dr->z, t1, group);
253
254 /* t0 = 2y^2 */
255 PREFIX(square) (t0, dp->y);
256 group->ecfp_reduce(t0, t0, group);
257 PREFIX(addShort) (t0, t0, t0);
258
259 /* S = 4 * px * py^2 = 2 * px * t0 */
260 PREFIX(multiply) (S, dp->x, t0);
261 PREFIX(addLong) (S, S, S);
262 group->ecfp_reduce(S, S, group);
263
264 /* rx = M^2 - 2 * S */
265 PREFIX(square) (t1, M);
266 PREFIX(subtractShort) (t1, t1, S);
267 PREFIX(subtractShort) (t1, t1, S);
268 group->ecfp_reduce(dr->x, t1, group);
269
270 /* ry = M * (S - rx) - 8 * py^4 */
271 PREFIX(square) (t1, t0); /* t1 = 4y^4 */
272 PREFIX(subtractShort) (S, S, dr->x);
273 PREFIX(multiply) (t0, M, S);
274 PREFIX(subtractLong) (t0, t0, t1);
275 PREFIX(subtractLong) (t0, t0, t1);
276 group->ecfp_reduce(dr->y, t0, group);
277
278 CLEANUP:
279 return;
280 }
281
282 /* Perform a point addition using coordinate system Jacobian + Affine ->
283 * Jacobian. Input and output should be multi-precision floating point
284 * integers. */
285 void PREFIX(pt_add_jac_aff) (const ecfp_jac_pt * p, const ecfp_aff_pt * q,
286 ecfp_jac_pt * r, const EC_group_fp * group) {
287 /* Temporary storage */
288 double A[2 * ECFP_NUMDOUBLES], B[2 * ECFP_NUMDOUBLES],
289 C[2 * ECFP_NUMDOUBLES], C2[2 * ECFP_NUMDOUBLES],
290 D[2 * ECFP_NUMDOUBLES], C3[2 * ECFP_NUMDOUBLES];
291
292 /* Check for point at infinity for p or q */
293 if (PREFIX(pt_is_inf_aff) (q) == MP_YES) {
294 PREFIX(copy) (r->x, p->x);
295 PREFIX(copy) (r->y, p->y);
296 PREFIX(copy) (r->z, p->z);
297 goto CLEANUP;
298 } else if (PREFIX(pt_is_inf_jac) (p) == MP_YES) {
299 PREFIX(copy) (r->x, q->x);
300 PREFIX(copy) (r->y, q->y);
301 /* Since the affine point is not infinity, we can set r->z = 1 */
302 PREFIX(one) (r->z);
303 goto CLEANUP;
304 }
305
306 /* Calculates c = qx * pz^2 - px d = (qy * b - py) rx = d^2 - c^3 + 2
307 * (px * c^2) ry = d * (c-rx) - py*c^3 rz = c * pz */
308
309 /* A = pz^2, B = pz^3 */
310 PREFIX(square) (A, p->z);
311 group->ecfp_reduce(A, A, group);
312 PREFIX(multiply) (B, A, p->z);
313 group->ecfp_reduce(B, B, group);
314
315 /* C = qx * A - px */
316 PREFIX(multiply) (C, q->x, A);
317 PREFIX(subtractShort) (C, C, p->x);
318 group->ecfp_reduce(C, C, group);
319
320 /* D = qy * B - py */
321 PREFIX(multiply) (D, q->y, B);
322 PREFIX(subtractShort) (D, D, p->y);
323 group->ecfp_reduce(D, D, group);
324
325 /* C2 = C^2, C3 = C^3 */
326 PREFIX(square) (C2, C);
327 group->ecfp_reduce(C2, C2, group);
328 PREFIX(multiply) (C3, C2, C);
329 group->ecfp_reduce(C3, C3, group);
330
331 /* rz = A = pz * C */
332 PREFIX(multiply) (A, p->z, C);
333 group->ecfp_reduce(r->z, A, group);
334
335 /* C = px * C^2, untidied, unreduced */
336 PREFIX(multiply) (C, p->x, C2);
337
338 /* A = D^2, untidied, unreduced */
339 PREFIX(square) (A, D);
340
341 /* rx = B = A - C3 - C - C = D^2 - (C^3 + 2 * (px * C^2) */
342 PREFIX(subtractShort) (A, A, C3);
343 PREFIX(subtractLong) (A, A, C);
344 PREFIX(subtractLong) (A, A, C);
345 group->ecfp_reduce(r->x, A, group);
346
347 /* B = py * C3, untidied, unreduced */
348 PREFIX(multiply) (B, p->y, C3);
349
350 /* C = px * C^2 - rx */
351 PREFIX(subtractShort) (C, C, r->x);
352 group->ecfp_reduce(C, C, group);
353
354 /* ry = A = D * C - py * C^3 */
355 PREFIX(multiply) (A, D, C);
356 PREFIX(subtractLong) (A, A, B);
357 group->ecfp_reduce(r->y, A, group);
358
359 CLEANUP:
360 return;
361 }
362
363 /* Perform a point addition using Jacobian coordinate system. Input and
364 * output should be multi-precision floating point integers. */
365 void PREFIX(pt_add_jac) (const ecfp_jac_pt * p, const ecfp_jac_pt * q,
366 ecfp_jac_pt * r, const EC_group_fp * group) {
367
368 /* Temporary Storage */
369 double t0[2 * ECFP_NUMDOUBLES], t1[2 * ECFP_NUMDOUBLES],
370 U[2 * ECFP_NUMDOUBLES], R[2 * ECFP_NUMDOUBLES],
371 S[2 * ECFP_NUMDOUBLES], H[2 * ECFP_NUMDOUBLES],
372 H3[2 * ECFP_NUMDOUBLES];
373
374 /* Check for point at infinity for p, if so set r = q */
375 if (PREFIX(pt_is_inf_jac) (p) == MP_YES) {
376 PREFIX(copy) (r->x, q->x);
377 PREFIX(copy) (r->y, q->y);
378 PREFIX(copy) (r->z, q->z);
379 goto CLEANUP;
380 }
381
382 /* Check for point at infinity for p, if so set r = q */
383 if (PREFIX(pt_is_inf_jac) (q) == MP_YES) {
384 PREFIX(copy) (r->x, p->x);
385 PREFIX(copy) (r->y, p->y);
386 PREFIX(copy) (r->z, p->z);
387 goto CLEANUP;
388 }
389
390 /* U = px * qz^2 , S = py * qz^3 */
391 PREFIX(square) (t0, q->z);
392 group->ecfp_reduce(t0, t0, group);
393 PREFIX(multiply) (U, p->x, t0);
394 group->ecfp_reduce(U, U, group);
395 PREFIX(multiply) (t1, t0, q->z);
396 group->ecfp_reduce(t1, t1, group);
397 PREFIX(multiply) (t0, p->y, t1);
398 group->ecfp_reduce(S, t0, group);
399
400 /* H = qx*(pz)^2 - U , R = (qy * pz^3 - S) */
401 PREFIX(square) (t0, p->z);
402 group->ecfp_reduce(t0, t0, group);
403 PREFIX(multiply) (H, q->x, t0);
404 PREFIX(subtractShort) (H, H, U);
405 group->ecfp_reduce(H, H, group);
406 PREFIX(multiply) (t1, t0, p->z); /* t1 = pz^3 */
407 group->ecfp_reduce(t1, t1, group);
408 PREFIX(multiply) (t0, t1, q->y); /* t0 = qy * pz^3 */
409 PREFIX(subtractShort) (t0, t0, S);
410 group->ecfp_reduce(R, t0, group);
411
412 /* U = U*H^2, H3 = H^3 */
413 PREFIX(square) (t0, H);
414 group->ecfp_reduce(t0, t0, group);
415 PREFIX(multiply) (t1, U, t0);
416 group->ecfp_reduce(U, t1, group);
417 PREFIX(multiply) (H3, t0, H);
418 group->ecfp_reduce(H3, H3, group);
419
420 /* rz = pz * qz * H */
421 PREFIX(multiply) (t0, q->z, H);
422 group->ecfp_reduce(t0, t0, group);
423 PREFIX(multiply) (t1, t0, p->z);
424 group->ecfp_reduce(r->z, t1, group);
425
426 /* rx = R^2 - H^3 - 2 * U */
427 PREFIX(square) (t0, R);
428 PREFIX(subtractShort) (t0, t0, H3);
429 PREFIX(subtractShort) (t0, t0, U);
430 PREFIX(subtractShort) (t0, t0, U);
431 group->ecfp_reduce(r->x, t0, group);
432
433 /* ry = R(U - rx) - S*H3 */
434 PREFIX(subtractShort) (t1, U, r->x);
435 PREFIX(multiply) (t0, t1, R);
436 PREFIX(multiply) (t1, S, H3);
437 PREFIX(subtractLong) (t1, t0, t1);
438 group->ecfp_reduce(r->y, t1, group);
439
440 CLEANUP:
441 return;
442 }
443
444 /* Perform a point doubling in Modified Jacobian coordinates. Input and
445 * output should be multi-precision floating point integers. */
446 void PREFIX(pt_dbl_jm) (const ecfp_jm_pt * p, ecfp_jm_pt * r,
447 const EC_group_fp * group) {
448
449 /* Temporary storage */
450 double t0[2 * ECFP_NUMDOUBLES], t1[2 * ECFP_NUMDOUBLES],
451 M[2 * ECFP_NUMDOUBLES], S[2 * ECFP_NUMDOUBLES],
452 U[2 * ECFP_NUMDOUBLES], T[2 * ECFP_NUMDOUBLES];
453
454 /* Check for point at infinity */
455 if (PREFIX(pt_is_inf_jm) (p) == MP_YES) {
456 /* Set r = pt at infinity by setting rz = 0 */
457 PREFIX(set_pt_inf_jm) (r);
458 goto CLEANUP;
459 }
460
461 /* M = 3 (px^2) + a*(pz^4) */
462 PREFIX(square) (t0, p->x);
463 PREFIX(addLong) (M, t0, t0);
464 PREFIX(addLong) (t0, t0, M); /* t0 = 3(px^2) */
465 PREFIX(addShort) (t0, t0, p->az4);
466 group->ecfp_reduce(M, t0, group);
467
468 /* rz = 2 * py * pz */
469 PREFIX(multiply) (t1, p->y, p->z);
470 PREFIX(addLong) (t1, t1, t1);
471 group->ecfp_reduce(r->z, t1, group);
472
473 /* t0 = 2y^2, U = 8y^4 */
474 PREFIX(square) (t0, p->y);
475 group->ecfp_reduce(t0, t0, group);
476 PREFIX(addShort) (t0, t0, t0);
477 PREFIX(square) (U, t0);
478 group->ecfp_reduce(U, U, group);
479 PREFIX(addShort) (U, U, U);
480
481 /* S = 4 * px * py^2 = 2 * px * t0 */
482 PREFIX(multiply) (S, p->x, t0);
483 group->ecfp_reduce(S, S, group);
484 PREFIX(addShort) (S, S, S);
485
486 /* rx = M^2 - 2S */
487 PREFIX(square) (T, M);
488 PREFIX(subtractShort) (T, T, S);
489 PREFIX(subtractShort) (T, T, S);
490 group->ecfp_reduce(r->x, T, group);
491
492 /* ry = M * (S - rx) - U */
493 PREFIX(subtractShort) (S, S, r->x);
494 PREFIX(multiply) (t0, M, S);
495 PREFIX(subtractShort) (t0, t0, U);
496 group->ecfp_reduce(r->y, t0, group);
497
498 /* ra*z^4 = 2*U*(apz4) */
499 PREFIX(multiply) (t1, U, p->az4);
500 PREFIX(addLong) (t1, t1, t1);
501 group->ecfp_reduce(r->az4, t1, group);
502
503 CLEANUP:
504 return;
505 }
506
507 /* Perform a point doubling using coordinates Affine -> Chudnovsky
508 * Jacobian. Input and output should be multi-precision floating point
509 * integers. */
510 void PREFIX(pt_dbl_aff2chud) (const ecfp_aff_pt * p, ecfp_chud_pt * r,
511 const EC_group_fp * group) {
512 double t0[2 * ECFP_NUMDOUBLES], t1[2 * ECFP_NUMDOUBLES],
513 M[2 * ECFP_NUMDOUBLES], twoY2[2 * ECFP_NUMDOUBLES],
514 S[2 * ECFP_NUMDOUBLES];
515
516 /* Check for point at infinity for p, if so set r = O */
517 if (PREFIX(pt_is_inf_aff) (p) == MP_YES) {
518 PREFIX(set_pt_inf_chud) (r);
519 goto CLEANUP;
520 }
521
522 /* M = 3(px)^2 + a */
523 PREFIX(square) (t0, p->x);
524 PREFIX(addLong) (t1, t0, t0);
525 PREFIX(addLong) (t1, t1, t0);
526 PREFIX(addShort) (t1, t1, group->curvea);
527 group->ecfp_reduce(M, t1, group);
528
529 /* twoY2 = 2*(py)^2, S = 4(px)(py)^2 */
530 PREFIX(square) (twoY2, p->y);
531 PREFIX(addLong) (twoY2, twoY2, twoY2);
532 group->ecfp_reduce(twoY2, twoY2, group);
533 PREFIX(multiply) (S, p->x, twoY2);
534 PREFIX(addLong) (S, S, S);
535 group->ecfp_reduce(S, S, group);
536
537 /* rx = M^2 - 2S */
538 PREFIX(square) (t0, M);
539 PREFIX(subtractShort) (t0, t0, S);
540 PREFIX(subtractShort) (t0, t0, S);
541 group->ecfp_reduce(r->x, t0, group);
542
543 /* ry = M(S-rx) - 8y^4 */
544 PREFIX(subtractShort) (t0, S, r->x);
545 PREFIX(multiply) (t1, t0, M);
546 PREFIX(square) (t0, twoY2);
547 PREFIX(subtractLong) (t1, t1, t0);
548 PREFIX(subtractLong) (t1, t1, t0);
549 group->ecfp_reduce(r->y, t1, group);
550
551 /* rz = 2py */
552 PREFIX(addShort) (r->z, p->y, p->y);
553
554 /* rz2 = rz^2 */
555 PREFIX(square) (t0, r->z);
556 group->ecfp_reduce(r->z2, t0, group);
557
558 /* rz3 = rz^3 */
559 PREFIX(multiply) (t0, r->z, r->z2);
560 group->ecfp_reduce(r->z3, t0, group);
561
562 CLEANUP:
563 return;
564 }
565
566 /* Perform a point addition using coordinates: Modified Jacobian +
567 * Chudnovsky Jacobian -> Modified Jacobian. Input and output should be
568 * multi-precision floating point integers. */
569 void PREFIX(pt_add_jm_chud) (ecfp_jm_pt * p, ecfp_chud_pt * q,
570 ecfp_jm_pt * r, const EC_group_fp * group) {
571
572 double t0[2 * ECFP_NUMDOUBLES], t1[2 * ECFP_NUMDOUBLES],
573 U[2 * ECFP_NUMDOUBLES], R[2 * ECFP_NUMDOUBLES],
574 S[2 * ECFP_NUMDOUBLES], H[2 * ECFP_NUMDOUBLES],
575 H3[2 * ECFP_NUMDOUBLES], pz2[2 * ECFP_NUMDOUBLES];
576
577 /* Check for point at infinity for p, if so set r = q need to convert
578 * from Chudnovsky form to Modified Jacobian form */
579 if (PREFIX(pt_is_inf_jm) (p) == MP_YES) {
580 PREFIX(copy) (r->x, q->x);
581 PREFIX(copy) (r->y, q->y);
582 PREFIX(copy) (r->z, q->z);
583 PREFIX(square) (t0, q->z2);
584 group->ecfp_reduce(t0, t0, group);
585 PREFIX(multiply) (t1, t0, group->curvea);
586 group->ecfp_reduce(r->az4, t1, group);
587 goto CLEANUP;
588 }
589 /* Check for point at infinity for q, if so set r = p */
590 if (PREFIX(pt_is_inf_chud) (q) == MP_YES) {
591 PREFIX(copy) (r->x, p->x);
592 PREFIX(copy) (r->y, p->y);
593 PREFIX(copy) (r->z, p->z);
594 PREFIX(copy) (r->az4, p->az4);
595 goto CLEANUP;
596 }
597
598 /* U = px * qz^2 */
599 PREFIX(multiply) (U, p->x, q->z2);
600 group->ecfp_reduce(U, U, group);
601
602 /* H = qx*(pz)^2 - U */
603 PREFIX(square) (t0, p->z);
604 group->ecfp_reduce(pz2, t0, group);
605 PREFIX(multiply) (H, pz2, q->x);
606 group->ecfp_reduce(H, H, group);
607 PREFIX(subtractShort) (H, H, U);
608
609 /* U = U*H^2, H3 = H^3 */
610 PREFIX(square) (t0, H);
611 group->ecfp_reduce(t0, t0, group);
612 PREFIX(multiply) (t1, U, t0);
613 group->ecfp_reduce(U, t1, group);
614 PREFIX(multiply) (H3, t0, H);
615 group->ecfp_reduce(H3, H3, group);
616
617 /* S = py * qz^3 */
618 PREFIX(multiply) (S, p->y, q->z3);
619 group->ecfp_reduce(S, S, group);
620
621 /* R = (qy * z1^3 - s) */
622 PREFIX(multiply) (t0, pz2, p->z);
623 group->ecfp_reduce(t0, t0, group);
624 PREFIX(multiply) (R, t0, q->y);
625 PREFIX(subtractShort) (R, R, S);
626 group->ecfp_reduce(R, R, group);
627
628 /* rz = pz * qz * H */
629 PREFIX(multiply) (t1, q->z, H);
630 group->ecfp_reduce(t1, t1, group);
631 PREFIX(multiply) (t0, p->z, t1);
632 group->ecfp_reduce(r->z, t0, group);
633
634 /* rx = R^2 - H^3 - 2 * U */
635 PREFIX(square) (t0, R);
636 PREFIX(subtractShort) (t0, t0, H3);
637 PREFIX(subtractShort) (t0, t0, U);
638 PREFIX(subtractShort) (t0, t0, U);
639 group->ecfp_reduce(r->x, t0, group);
640
641 /* ry = R(U - rx) - S*H3 */
642 PREFIX(subtractShort) (t1, U, r->x);
643 PREFIX(multiply) (t0, t1, R);
644 PREFIX(multiply) (t1, S, H3);
645 PREFIX(subtractLong) (t1, t0, t1);
646 group->ecfp_reduce(r->y, t1, group);
647
648 if (group->aIsM3) { /* a == -3 */
649 /* a(rz^4) = -3 * ((rz^2)^2) */
650 PREFIX(square) (t0, r->z);
651 group->ecfp_reduce(t0, t0, group);
652 PREFIX(square) (t1, t0);
653 PREFIX(addLong) (t0, t1, t1);
654 PREFIX(addLong) (t0, t0, t1);
655 PREFIX(negLong) (t0, t0);
656 group->ecfp_reduce(r->az4, t0, group);
657 } else { /* Generic case */
658 /* a(rz^4) = a * ((rz^2)^2) */
659 PREFIX(square) (t0, r->z);
660 group->ecfp_reduce(t0, t0, group);
661 PREFIX(square) (t1, t0);
662 group->ecfp_reduce(t1, t1, group);
663 PREFIX(multiply) (t0, group->curvea, t1);
664 group->ecfp_reduce(r->az4, t0, group);
665 }
666 CLEANUP:
667 return;
668 }
669
670 /* Perform a point addition using Chudnovsky Jacobian coordinates. Input
671 * and output should be multi-precision floating point integers. */
672 void PREFIX(pt_add_chud) (const ecfp_chud_pt * p, const ecfp_chud_pt * q,
673 ecfp_chud_pt * r, const EC_group_fp * group) {
674
675 /* Temporary Storage */
676 double t0[2 * ECFP_NUMDOUBLES], t1[2 * ECFP_NUMDOUBLES],
677 U[2 * ECFP_NUMDOUBLES], R[2 * ECFP_NUMDOUBLES],
678 S[2 * ECFP_NUMDOUBLES], H[2 * ECFP_NUMDOUBLES],
679 H3[2 * ECFP_NUMDOUBLES];
680
681 /* Check for point at infinity for p, if so set r = q */
682 if (PREFIX(pt_is_inf_chud) (p) == MP_YES) {
683 PREFIX(copy) (r->x, q->x);
684 PREFIX(copy) (r->y, q->y);
685 PREFIX(copy) (r->z, q->z);
686 PREFIX(copy) (r->z2, q->z2);
687 PREFIX(copy) (r->z3, q->z3);
688 goto CLEANUP;
689 }
690
691 /* Check for point at infinity for p, if so set r = q */
692 if (PREFIX(pt_is_inf_chud) (q) == MP_YES) {
693 PREFIX(copy) (r->x, p->x);
694 PREFIX(copy) (r->y, p->y);
695 PREFIX(copy) (r->z, p->z);
696 PREFIX(copy) (r->z2, p->z2);
697 PREFIX(copy) (r->z3, p->z3);
698 goto CLEANUP;
699 }
700
701 /* U = px * qz^2 */
702 PREFIX(multiply) (U, p->x, q->z2);
703 group->ecfp_reduce(U, U, group);
704
705 /* H = qx*(pz)^2 - U */
706 PREFIX(multiply) (H, q->x, p->z2);
707 PREFIX(subtractShort) (H, H, U);
708 group->ecfp_reduce(H, H, group);
709
710 /* U = U*H^2, H3 = H^3 */
711 PREFIX(square) (t0, H);
712 group->ecfp_reduce(t0, t0, group);
713 PREFIX(multiply) (t1, U, t0);
714 group->ecfp_reduce(U, t1, group);
715 PREFIX(multiply) (H3, t0, H);
716 group->ecfp_reduce(H3, H3, group);
717
718 /* S = py * qz^3 */
719 PREFIX(multiply) (S, p->y, q->z3);
720 group->ecfp_reduce(S, S, group);
721
722 /* rz = pz * qz * H */
723 PREFIX(multiply) (t0, q->z, H);
724 group->ecfp_reduce(t0, t0, group);
725 PREFIX(multiply) (t1, t0, p->z);
726 group->ecfp_reduce(r->z, t1, group);
727
728 /* R = (qy * z1^3 - s) */
729 PREFIX(multiply) (t0, q->y, p->z3);
730 PREFIX(subtractShort) (t0, t0, S);
731 group->ecfp_reduce(R, t0, group);
732
733 /* rx = R^2 - H^3 - 2 * U */
734 PREFIX(square) (t0, R);
735 PREFIX(subtractShort) (t0, t0, H3);
736 PREFIX(subtractShort) (t0, t0, U);
737 PREFIX(subtractShort) (t0, t0, U);
738 group->ecfp_reduce(r->x, t0, group);
739
740 /* ry = R(U - rx) - S*H3 */
741 PREFIX(subtractShort) (t1, U, r->x);
742 PREFIX(multiply) (t0, t1, R);
743 PREFIX(multiply) (t1, S, H3);
744 PREFIX(subtractLong) (t1, t0, t1);
745 group->ecfp_reduce(r->y, t1, group);
746
747 /* rz2 = rz^2 */
748 PREFIX(square) (t0, r->z);
749 group->ecfp_reduce(r->z2, t0, group);
750
751 /* rz3 = rz^3 */
752 PREFIX(multiply) (t0, r->z, r->z2);
753 group->ecfp_reduce(r->z3, t0, group);
754
755 CLEANUP:
756 return;
757 }
758
759 /* Expects out to be an array of size 16 of Chudnovsky Jacobian points.
760 * Fills in Chudnovsky Jacobian form (x, y, z, z^2, z^3), for -15P, -13P,
761 * -11P, -9P, -7P, -5P, -3P, -P, P, 3P, 5P, 7P, 9P, 11P, 13P, 15P */
762 void PREFIX(precompute_chud) (ecfp_chud_pt * out, const ecfp_aff_pt * p,
763 const EC_group_fp * group) {
764
765 ecfp_chud_pt p2;
766
767 /* Set out[8] = P */
768 PREFIX(copy) (out[8].x, p->x);
769 PREFIX(copy) (out[8].y, p->y);
770 PREFIX(one) (out[8].z);
771 PREFIX(one) (out[8].z2);
772 PREFIX(one) (out[8].z3);
773
774 /* Set p2 = 2P */
775 PREFIX(pt_dbl_aff2chud) (p, &p2, group);
776
777 /* Set 3P, 5P, ..., 15P */
778 PREFIX(pt_add_chud) (&out[8], &p2, &out[9], group);
779 PREFIX(pt_add_chud) (&out[9], &p2, &out[10], group);
780 PREFIX(pt_add_chud) (&out[10], &p2, &out[11], group);
781 PREFIX(pt_add_chud) (&out[11], &p2, &out[12], group);
782 PREFIX(pt_add_chud) (&out[12], &p2, &out[13], group);
783 PREFIX(pt_add_chud) (&out[13], &p2, &out[14], group);
784 PREFIX(pt_add_chud) (&out[14], &p2, &out[15], group);
785
786 /* Set -15P, -13P, ..., -P */
787 PREFIX(pt_neg_chud) (&out[8], &out[7]);
788 PREFIX(pt_neg_chud) (&out[9], &out[6]);
789 PREFIX(pt_neg_chud) (&out[10], &out[5]);
790 PREFIX(pt_neg_chud) (&out[11], &out[4]);
791 PREFIX(pt_neg_chud) (&out[12], &out[3]);
792 PREFIX(pt_neg_chud) (&out[13], &out[2]);
793 PREFIX(pt_neg_chud) (&out[14], &out[1]);
794 PREFIX(pt_neg_chud) (&out[15], &out[0]);
795 }
796
797 /* Expects out to be an array of size 16 of Jacobian points. Fills in
798 * Jacobian form (x, y, z), for O, P, 2P, ... 15P */
799 void PREFIX(precompute_jac) (ecfp_jac_pt * precomp, const ecfp_aff_pt * p,
800 const EC_group_fp * group) {
801 int i;
802
803 /* fill precomputation table */
804 /* set precomp[0] */
805 PREFIX(set_pt_inf_jac) (&precomp[0]);
806 /* set precomp[1] */
807 PREFIX(copy) (precomp[1].x, p->x);
808 PREFIX(copy) (precomp[1].y, p->y);
809 if (PREFIX(pt_is_inf_aff) (p) == MP_YES) {
810 PREFIX(zero) (precomp[1].z);
811 } else {
812 PREFIX(one) (precomp[1].z);
813 }
814 /* set precomp[2] */
815 group->pt_dbl_jac(&precomp[1], &precomp[2], group);
816
817 /* set rest of precomp */
818 for (i = 3; i < 16; i++) {
819 group->pt_add_jac_aff(&precomp[i - 1], p, &precomp[i], group);
820 }
821 }

Mercurial Repository: The Tor Browser

Europalab SCM Repository Server

mercurial