1 #include <stdarg.h>
2 #include <stdio.h>
3 #include <stdint.h> // for intptr_t
4 #include <stdlib.h> //for rand, pbc_malloc, pbc_free
5 #include <string.h> //for strcmp
6 #include <gmp.h>
7 #include "pbc_utils.h"
8 #include "pbc_field.h"
9 #include "pbc_fp.h"
10 #include "pbc_fieldquadratic.h"
11 #include "pbc_param.h"
12 #include "pbc_pairing.h"
13 #include "pbc_curve.h"
14 #include "pbc_random.h"
15 #include "pbc_memory.h"
16 #include "ecc/param.h"
17 #include "pbc_a_param.h"
18 #include "pbc_a1_param.h"
19
20 typedef struct {
21 int exp2;
22 int exp1;
23 int sign1;
24 int sign0;
25 mpz_t r; // r = 2^exp2 + sign1 * 2^exp1 + sign0 * 1
26 mpz_t q; // we work in E(F_q) (and E(F_q^2))
27 mpz_t h; // r * h = q + 1
28 } *a_param_ptr;
29
30 typedef struct {
31 field_t Fq, Fq2, Eq;
32 int exp2, exp1;
33 int sign1;
34 } *a_pairing_data_ptr;
35
a_out_str(FILE * stream,void * data)36 static void a_out_str(FILE *stream, void *data) {
37 a_param_ptr p = data;
38 param_out_type(stream, "a");
39 param_out_mpz(stream, "q", p->q);
40 param_out_mpz(stream, "h", p->h);
41 param_out_mpz(stream, "r", p->r);
42 param_out_int(stream, "exp2", p->exp2);
43 param_out_int(stream, "exp1", p->exp1);
44 param_out_int(stream, "sign1", p->sign1);
45 param_out_int(stream, "sign0", p->sign0);
46 }
47
a_clear(void * data)48 static void a_clear(void *data) {
49 a_param_ptr sp = data;
50 mpz_clear(sp->r);
51 mpz_clear(sp->q);
52 mpz_clear(sp->h);
53 pbc_free(data);
54 }
55
phi_identity(element_ptr out,element_ptr in,pairing_ptr pairing)56 static void phi_identity(element_ptr out, element_ptr in, pairing_ptr pairing) {
57 UNUSED_VAR(pairing);
58 element_set(out, in);
59 }
60
compute_abc_tangent(element_ptr a,element_ptr b,element_ptr c,element_ptr Vx,element_ptr Vy,element_ptr e0)61 static void compute_abc_tangent(element_ptr a, element_ptr b, element_ptr c,
62 element_ptr Vx, element_ptr Vy, element_ptr e0) {
63 //a = -slope_tangent(V.x, V.y);
64 //b = 1;
65 //c = -(V.y + aV.x);
66 //but we multiply by -2*V.y to avoid division so:
67 //a = -(3 Vx^2 + cc->a)
68 //b = 2 * Vy
69 //c = -(2 Vy^2 + a Vx);
70 element_square(a, Vx);
71 //element_mul_si(a, a, 3);
72 element_add(e0, a, a);
73 element_add(a, e0, a);
74 element_set1(b);
75 element_add(a, a, b);
76 element_neg(a, a);
77
78 element_double(b, Vy);
79
80 element_mul(e0, b, Vy);
81 element_mul(c, a, Vx);
82 element_add(c, c, e0);
83 element_neg(c, c);
84 }
85
compute_abc_tangent_proj(element_ptr a,element_ptr b,element_ptr c,element_ptr Vx,element_ptr Vy,element_ptr z,element_ptr z2,element_ptr e0)86 static void compute_abc_tangent_proj(element_ptr a, element_ptr b, element_ptr c,
87 element_ptr Vx, element_ptr Vy,
88 element_ptr z, element_ptr z2, element_ptr e0) {
89 //a = -(3x^2 + cca z^4)
90 //for this case cca = 1
91 //b = 2 y z^3
92 //c = -(2 y^2 + x a)
93 //a = z^2 a
94 element_square(a, z2);
95 element_square(b, Vx);
96 ////element_mul_si(b, b, 3);
97 element_double(e0, b);
98 element_add(b, e0, b);
99 element_add(a, a, b);
100 element_neg(a, a);
101
102 ////element_mul_si(e0, Vy, 2);
103 element_double(e0, Vy);
104 element_mul(b, e0, z2);
105 element_mul(b, b, z);
106
107 element_mul(c, Vx, a);
108 element_mul(a, a, z2);
109 element_mul(e0, e0, Vy);
110 element_add(c, c, e0);
111 element_neg(c, c);
112 }
113
compute_abc_line(element_ptr a,element_ptr b,element_ptr c,element_ptr Vx,element_ptr Vy,element_ptr V1x,element_ptr V1y,element_ptr e0)114 static void compute_abc_line(element_ptr a, element_ptr b, element_ptr c,
115 element_ptr Vx, element_ptr Vy,
116 element_ptr V1x, element_ptr V1y,
117 element_ptr e0) {
118 //a = -(B.y - A.y) / (B.x - A.x);
119 //b = 1;
120 //c = -(A.y + a * A.x);
121 //but we'll multiply by B.x - A.x to avoid division, so
122 //a = -(By - Ay)
123 //b = Bx - Ax
124 //c = -(Ay b + a Ax);
125 element_sub(a, Vy, V1y);
126 element_sub(b, V1x, Vx);
127 element_mul(c, Vx, V1y);
128 element_mul(e0, Vy, V1x);
129 element_sub(c, c, e0);
130 }
131
132 struct pp_coeff_s {
133 element_t a;
134 element_t b;
135 element_t c;
136 };
137 typedef struct pp_coeff_s pp_coeff_t[1];
138 typedef struct pp_coeff_s *pp_coeff_ptr;
139
pp_coeff_set(pp_coeff_ptr p,element_t a,element_t b,element_t c)140 static void pp_coeff_set(pp_coeff_ptr p, element_t a, element_t b, element_t c) {
141 element_init(p->a, a->field);
142 element_init(p->b, b->field);
143 element_init(p->c, c->field);
144 element_set(p->a, a);
145 element_set(p->b, b);
146 element_set(p->c, c);
147 }
148
a_pairing_pp_init(pairing_pp_t p,element_ptr in1,pairing_t pairing)149 static void a_pairing_pp_init(pairing_pp_t p, element_ptr in1, pairing_t pairing) {
150 int i, n;
151 a_pairing_data_ptr ainfo = pairing->data;
152 p->data = pbc_malloc(sizeof(pp_coeff_t) * (ainfo->exp2 + 1));
153 pp_coeff_t *coeff = (pp_coeff_t *) p->data;
154 element_t V, V1;
155 element_t a, b, c;
156 element_t e0;
157 element_ptr Vx, Vy;
158 element_ptr V1x, V1y;
159
160 #define do_tangent() \
161 compute_abc_tangent(a, b, c, Vx, Vy, e0); \
162 pp_coeff_set(coeff[i], a, b, c);
163
164 #define do_line() \
165 compute_abc_line(a, b, c, Vx, Vy, V1x, V1y, e0); \
166 pp_coeff_set(coeff[i], a, b, c);
167
168 element_init(V, ainfo->Eq);
169 element_init(V1, ainfo->Eq);
170 element_set(V, in1);
171 Vx = curve_x_coord(V);
172 Vy = curve_y_coord(V);
173 V1x = curve_x_coord(V1);
174 V1y = curve_y_coord(V1);
175 element_init(e0, ainfo->Fq);
176 element_init(a, ainfo->Fq);
177 element_init(b, ainfo->Fq);
178 element_init(c, ainfo->Fq);
179
180 n = ainfo->exp1;
181 for (i=0; i<n; i++) {
182 do_tangent();
183 element_double(V, V);
184 }
185
186 if (ainfo->sign1 < 0) {
187 element_neg(V1, V);
188 } else {
189 element_set(V1, V);
190 }
191 n = ainfo->exp2;
192 for (; i<n; i++) {
193 do_tangent();
194 element_double(V, V);
195 }
196
197 do_line();
198
199 element_clear(e0);
200 element_clear(a);
201 element_clear(b);
202 element_clear(c);
203 element_clear(V);
204 element_clear(V1);
205 #undef do_tangent
206 #undef do_line
207 }
208
a_pairing_pp_clear(pairing_pp_t p)209 static void a_pairing_pp_clear(pairing_pp_t p) {
210 a_pairing_data_ptr ainfo = p->pairing->data;
211 pp_coeff_t *coeff = (pp_coeff_t *) p->data;
212 int i, n = ainfo->exp2 + 1;
213 for (i=0; i<n; i++) {
214 pp_coeff_ptr pp = coeff[i];
215 element_clear(pp->a);
216 element_clear(pp->b);
217 element_clear(pp->c);
218 }
219 pbc_free(p->data);
220 }
221
222 // Requires cofactor to be odd.
223 // Overwrites in and temp, out != in.
224 // Luckily this touchy routine is only used internally.
225 // TODO: rewrite to allow (out == in)? would simplify a_finalpow()
lucas_odd(element_ptr out,element_ptr in,element_ptr temp,mpz_t cofactor)226 static void lucas_odd(element_ptr out, element_ptr in, element_ptr temp, mpz_t cofactor) {
227 element_ptr in0 = element_x(in);
228 element_ptr in1 = element_y(in);
229 element_ptr v0 = element_x(out);
230 element_ptr v1 = element_y(out);
231 element_ptr t0 = element_x(temp);
232 element_ptr t1 = element_y(temp);
233 int j;
234
235 element_set_si(t0, 2);
236 element_double(t1, in0);
237
238 element_set(v0, t0);
239 element_set(v1, t1);
240
241 j = mpz_sizeinbase(cofactor, 2) - 1;
242 for (;;) {
243 if (!j) {
244 element_mul(v1, v0, v1);
245 element_sub(v1, v1, t1);
246 element_square(v0, v0);
247 element_sub(v0, v0, t0);
248 break;
249 }
250 if (mpz_tstbit(cofactor, j)) {
251 element_mul(v0, v0, v1);
252 element_sub(v0, v0, t1);
253 element_square(v1, v1);
254 element_sub(v1, v1, t0);
255 } else {
256 element_mul(v1, v0, v1);
257 element_sub(v1, v1, t1);
258 element_square(v0, v0);
259 element_sub(v0, v0, t0);
260 }
261 j--;
262 }
263
264 //assume cofactor = (q + 1) / r is even
265 //(r should be odd and q + 1 is always even)
266 //thus v0 = V_k, v1 = V_{k+1}
267 //and V_{k-1} = P v0 - v1
268
269 //so U_k = (P V_k - 2 V_{k-1}) / (P^2 - 4)
270 // = (2 v1 - P v0) / (P^2 - 4)
271
272 element_mul(in0, v0, t1);
273 element_double(v1, v1);
274 element_sub(v1, v1, in0);
275
276 element_square(t1, t1);
277 element_sub(t1, t1, t0);
278 element_sub(t1, t1, t0);
279 element_div(v1, v1, t1);
280
281 element_halve(v0, v0);
282 element_mul(v1, v1, in1);
283 }
284
a_tateexp(element_ptr out,element_ptr in,element_ptr temp,mpz_t cofactor)285 static inline void a_tateexp(element_ptr out, element_ptr in, element_ptr temp, mpz_t cofactor) {
286 element_ptr in1 = element_y(in);
287 //simpler but slower:
288 //element_pow_mpz(out, f, tateexp);
289
290 //1. Exponentiate by q-1
291 //which is equivalent to the following
292
293 element_invert(temp, in);
294 element_neg(in1, in1);
295 element_mul(in, in, temp);
296
297 //2. Exponentiate by (q+1)/r
298
299 //Instead of:
300 // element_pow_mpz(out, in, cofactor);
301 //we use Lucas sequences (see "Compressed Pairings", Scott and Barreto)
302 lucas_odd(out, in, temp, cofactor);
303 }
304
305 //computes a Qx + b Qy + c for type A pairing
a_miller_evalfn(element_ptr out,element_ptr a,element_ptr b,element_ptr c,element_ptr Qx,element_ptr Qy)306 static inline void a_miller_evalfn(element_ptr out,
307 element_ptr a, element_ptr b, element_ptr c,
308 element_ptr Qx, element_ptr Qy) {
309 //we'll map Q via (x,y) --> (-x, iy)
310 //hence Re(a Qx + b Qy + c) = -a Q'x + c and
311 //Im(a Qx + b Qy + c) = b Q'y
312 element_mul(element_y(out), a, Qx);
313 element_sub(element_x(out), c, element_y(out));
314 element_mul(element_y(out), b, Qy);
315 }
316
a_pairing_pp_apply(element_ptr out,element_ptr in2,pairing_pp_t p)317 static void a_pairing_pp_apply(element_ptr out, element_ptr in2, pairing_pp_t p) {
318 //TODO: use proj coords here too to shave off a little time
319 element_ptr Qx = curve_x_coord(in2);
320 element_ptr Qy = curve_y_coord(in2);
321 element_t f, f0;
322 int i, n;
323 a_pairing_data_ptr ainfo = p->pairing->data;
324 pp_coeff_t *coeff = p->data;
325 element_init(f, ainfo->Fq2);
326 element_init(f0, ainfo->Fq2);
327
328 element_set1(f);
329 n = ainfo->exp1;
330 for (i=0; i<n; i++) {
331 pp_coeff_ptr pp = coeff[i];
332 element_square(f, f);
333 a_miller_evalfn(f0, pp->a, pp->b, pp->c, Qx, Qy);
334 element_mul(f, f, f0);
335 }
336 if (ainfo->sign1 < 0) {
337 element_invert(out, f);
338 } else {
339 element_set(out, f);
340 }
341 n = ainfo->exp2;
342 for (; i<n; i++) {
343 element_square(f, f);
344 pp_coeff_ptr pp = coeff[i];
345 a_miller_evalfn(f0, pp->a, pp->b, pp->c, Qx, Qy);
346 element_mul(f, f, f0);
347 }
348
349 element_mul(f, f, out);
350 {
351 pp_coeff_ptr pp = coeff[i];
352 a_miller_evalfn(f0, pp->a, pp->b, pp->c, Qx, Qy);
353 element_mul(f, f, f0);
354 }
355
356 a_tateexp(out, f, f0, p->pairing->phikonr);
357
358 element_clear(f);
359 element_clear(f0);
360 }
361
362 // in1, in2 are from E(F_q), out from F_q^2.
363 // Pairing via elliptic nets (see Stange).
a_pairing_ellnet(element_ptr out,element_ptr in1,element_ptr in2,pairing_t pairing)364 static void a_pairing_ellnet(element_ptr out, element_ptr in1, element_ptr in2,
365 pairing_t pairing) {
366 element_ptr x = curve_x_coord(in1);
367 element_ptr y = curve_y_coord(in1);
368
369 element_ptr x2 = curve_x_coord(in2);
370 element_ptr y2 = curve_y_coord(in2);
371
372 //we map (x2,y2) to (-x2, i y2) before pairing
373 //notation: cmi means c_{k-i}, ci means c_{k+i}
374 element_t cm3, cm2, cm1, c0, c1, c2, c3, c4;
375 element_t dm1, d0, d1;
376 element_t A, B, C;
377
378 element_init_same_as(cm3, x);
379 element_init_same_as(cm2, x);
380 element_init_same_as(cm1, x);
381 element_init_same_as(c0, x);
382 element_init_same_as(c1, x);
383 element_init_same_as(c2, x);
384 element_init_same_as(c3, x);
385 element_init_same_as(c4, x);
386 element_init_same_as(C, x);
387
388 element_init_same_as(dm1, out);
389 element_init_same_as(d0, out);
390 element_init_same_as(d1, out);
391 element_init_same_as(A, x);
392 element_init_same_as(B, out);
393
394 // c1 = 2y
395 // c0 = 1
396 // cm2 = -1
397 // cm3 = -2y
398 element_double(c1, y);
399 element_set1(c0);
400 element_neg(cm3, c1);
401 element_neg(cm2, c0);
402
403 // a = 1, b = 0 for Y^2 = X^3 + X
404 //hence c3 = c_{k+3} = c_4 = 4y(x^6 + 5(x^4 - x^2) - 1)
405 //use cm1, C, c2 as temp variables for now
406 element_square(cm1, x);
407 element_square(C, cm1);
408 element_sub(c2, C, cm1);
409 element_double(c3, c2);
410 element_double(c3, c3);
411 element_add(c3, c3, c2);
412 element_mul(c2, C, cm1);
413 element_add(c3, c3, c2);
414 element_add(c3, c3, cm2);
415 element_mul(c3, c3, c1);
416 element_double(c3, c3);
417
418 // c2 = c_3 = 3x^4 + 6x^2 - 1
419 element_double(cm1, cm1);
420 element_add(cm1, cm1, C);
421 element_double(C, cm1);
422 element_add(C, C, cm1);
423 element_add(c2, C, cm2);
424
425 // c4 = c_5 = c_2^3 c_4 - c_3^3 = c1^3 c3 - c2^3
426 element_square(C, c1);
427 element_mul(c4, C, c1);
428 element_mul(c4, c4, c3);
429 element_square(C, c2);
430 element_mul(C, C, c2);
431 element_sub(c4, c4, C);
432
433 //compute A, B, d1 (which is d_2 since k = 1)
434 //(recall phi takes x2 to -x2, y2 to i y2)
435 element_add(A, x, x2);
436 element_double(C, x);
437 element_sub(C, C, x2);
438 element_square(cm1, A);
439 element_mul(cm1, C, cm1);
440 element_set(element_x(d1), y);
441 element_set(element_y(d1), y2);
442 element_square(d1, d1);
443 element_sub(element_x(d1), element_x(d1), cm1);
444 element_neg(B, d1);
445 element_invert(B, B);
446 element_invert(A, A);
447 element_mul(element_x(d1), y, A);
448 element_neg(element_x(d1), element_x(d1));
449 element_mul(element_y(d1), y2, A);
450 element_square(d1, d1);
451 element_sub(element_x(d1), C, element_x(d1));
452 element_neg(element_y(d1), element_y(d1));
453
454 // cm1 = 0
455 // C = (2y)^-1
456 element_set0(cm1);
457 element_invert(C, c1);
458
459 element_set1(dm1);
460 element_set1(d0);
461
462 element_t sm2, sm1;
463 element_t s0, s1, s2, s3;
464 element_t tm2, tm1;
465 element_t t0, t1, t2, t3;
466 element_t e0, e1;
467 element_t u, v;
468
469 element_init_same_as(sm2, x);
470 element_init_same_as(sm1, x);
471 element_init_same_as(s0, x);
472 element_init_same_as(s1, x);
473 element_init_same_as(s2, x);
474 element_init_same_as(s3, x);
475
476 element_init_same_as(tm2, x);
477 element_init_same_as(tm1, x);
478 element_init_same_as(t0, x);
479 element_init_same_as(t1, x);
480 element_init_same_as(t2, x);
481 element_init_same_as(t3, x);
482
483 element_init_same_as(e0, x);
484 element_init_same_as(e1, x);
485
486 element_init_same_as(u, d0);
487 element_init_same_as(v, d0);
488
489 int m = mpz_sizeinbase(pairing->r, 2) - 2;
490 for (;;) {
491 element_square(sm2, cm2);
492 element_square(sm1, cm1);
493 element_square(s0, c0);
494 element_square(s1, c1);
495 element_square(s2, c2);
496 element_square(s3, c3);
497
498 element_mul(tm2, cm3, cm1);
499 element_mul(tm1, cm2, c0);
500 element_mul(t0, cm1, c1);
501 element_mul(t1, c0, c2);
502 element_mul(t2, c1, c3);
503 element_mul(t3, c2, c4);
504
505 element_square(u, d0);
506 element_mul(v, dm1, d1);
507
508 if (mpz_tstbit(pairing->r, m)) {
509 //double-and-add
510 element_mul(e0, t0, sm2);
511 element_mul(e1, tm2, s0);
512 element_sub(cm3, e0, e1);
513 element_mul(cm3, cm3, C);
514
515 element_mul(e0, t0, sm1);
516 element_mul(e1, tm1, s0);
517 element_sub(cm2, e0, e1);
518
519 element_mul(e0, t1, sm1);
520 element_mul(e1, tm1, s1);
521 element_sub(cm1, e0, e1);
522 element_mul(cm1, cm1, C);
523
524 element_mul(e0, t1, s0);
525 element_mul(e1, t0, s1);
526 element_sub(c0, e0, e1);
527
528 element_mul(e0, t2, s0);
529 element_mul(e1, t0, s2);
530 element_sub(c1, e0, e1);
531 element_mul(c1, c1, C);
532
533 element_mul(e0, t2, s1);
534 element_mul(e1, t1, s2);
535 element_sub(c2, e0, e1);
536
537 element_mul(e0, t3, s1);
538 element_mul(e1, t1, s3);
539 element_sub(c3, e0, e1);
540 element_mul(c3, c3, C);
541
542 element_mul(e0, t3, s2);
543 element_mul(e1, t2, s3);
544 element_sub(c4, e0, e1);
545
546 element_mul(element_x(out), element_x(u), t0);
547 element_mul(element_y(out), element_y(u), t0);
548 element_mul(element_x(dm1), element_x(v), s0);
549 element_mul(element_y(dm1), element_y(v), s0);
550 element_sub(dm1, dm1, out);
551
552 element_mul(element_x(out), element_x(u), t1);
553 element_mul(element_y(out), element_y(u), t1);
554 element_mul(element_x(d0), element_x(v), s1);
555 element_mul(element_y(d0), element_y(v), s1);
556 element_sub(d0, d0, out);
557 element_mul(element_x(d0), element_x(d0), A);
558 element_mul(element_y(d0), element_y(d0), A);
559
560 element_mul(element_x(out), element_x(u), t2);
561 element_mul(element_y(out), element_y(u), t2);
562 element_mul(element_x(d1), element_x(v), s2);
563 element_mul(element_y(d1), element_y(v), s2);
564 element_sub(d1, d1, out);
565 element_mul(d1, d1, B);
566 } else {
567 //double
568 element_mul(e0, tm1, sm2);
569 element_mul(e1, tm2, sm1);
570 element_sub(cm3, e0, e1);
571
572 element_mul(e0, t0, sm2);
573 element_mul(e1, tm2, s0);
574 element_sub(cm2, e0, e1);
575 element_mul(cm2, cm2, C);
576
577 element_mul(e0, t0, sm1);
578 element_mul(e1, tm1, s0);
579 element_sub(cm1, e0, e1);
580
581 element_mul(e0, t1, sm1);
582 element_mul(e1, tm1, s1);
583 element_sub(c0, e0, e1);
584 element_mul(c0, c0, C);
585
586 element_mul(e0, t1, s0);
587 element_mul(e1, t0, s1);
588 element_sub(c1, e0, e1);
589
590 element_mul(e0, t2, s0);
591 element_mul(e1, t0, s2);
592 element_sub(c2, e0, e1);
593 element_mul(c2, c2, C);
594
595 element_mul(e0, t2, s1);
596 element_mul(e1, t1, s2);
597 element_sub(c3, e0, e1);
598
599 element_mul(e0, t3, s1);
600 element_mul(e1, t1, s3);
601 element_sub(c4, e0, e1);
602 element_mul(c4, c4, C);
603
604 element_mul(element_x(out), element_x(u), tm1);
605 element_mul(element_y(out), element_y(u), tm1);
606 element_mul(element_x(dm1), element_x(v), sm1);
607 element_mul(element_y(dm1), element_y(v), sm1);
608 element_sub(dm1, dm1, out);
609
610 element_mul(element_x(out), element_x(u), t0);
611 element_mul(element_y(out), element_y(u), t0);
612 element_mul(element_x(d0), element_x(v), s0);
613 element_mul(element_y(d0), element_y(v), s0);
614 element_sub(d0, d0, out);
615
616 element_mul(element_x(out), element_x(u), t1);
617 element_mul(element_y(out), element_y(u), t1);
618 element_mul(element_x(d1), element_x(v), s1);
619 element_mul(element_y(d1), element_y(v), s1);
620 element_sub(d1, d1, out);
621 element_mul(element_x(d1), element_x(d1), A);
622 element_mul(element_y(d1), element_y(d1), A);
623 }
624 if (!m) break;
625 m--;
626 }
627 // since c_k lies base field
628 // it gets killed by the final powering
629 //element_invert(c1, c1);
630 //element_mul(element_x(d1), element_x(d1), c1);
631 //element_mul(element_y(d1), element_y(d1), c1);
632
633 a_tateexp(out, d1, d0, pairing->phikonr);
634
635 element_clear(dm1);
636 element_clear(d0);
637 element_clear(d1);
638
639 element_clear(cm3);
640 element_clear(cm2);
641 element_clear(cm1);
642 element_clear(c0);
643 element_clear(c1);
644 element_clear(c2);
645 element_clear(c3);
646 element_clear(c4);
647
648 element_clear(sm2);
649 element_clear(sm1);
650 element_clear(s0);
651 element_clear(s1);
652 element_clear(s2);
653 element_clear(s3);
654
655 element_clear(tm2);
656 element_clear(tm1);
657 element_clear(t0);
658 element_clear(t1);
659 element_clear(t2);
660 element_clear(t3);
661
662 element_clear(e0);
663 element_clear(e1);
664 element_clear(A);
665 element_clear(B);
666 element_clear(C);
667 element_clear(u);
668 element_clear(v);
669 }
670
671 struct ellnet_pp_st_s {
672 element_t sm1, s0, s1, s2;
673 element_t tm1, t0, t1, t2;
674 };
675 typedef struct ellnet_pp_st_s ellnet_pp_st_t[1];
676 typedef struct ellnet_pp_st_s *ellnet_pp_st_ptr;
677
678 struct ellnet_pp_s {
679 element_t x;
680 element_t y;
681 ellnet_pp_st_t *seq;
682 };
683 typedef struct ellnet_pp_s ellnet_pp_t[1];
684 typedef struct ellnet_pp_s *ellnet_pp_ptr;
685
a_pairing_ellnet_pp_init(pairing_pp_t p,element_ptr in1,pairing_t pairing)686 static void a_pairing_ellnet_pp_init(pairing_pp_t p, element_ptr in1, pairing_t pairing) {
687 element_ptr x = curve_x_coord(in1);
688 element_ptr y = curve_y_coord(in1);
689 int i, rbits = mpz_sizeinbase(pairing->r, 2);
690 ellnet_pp_ptr pp = p->data = pbc_malloc(sizeof(ellnet_pp_t));
691 pp->seq = pbc_malloc(sizeof(ellnet_pp_st_t) * rbits);
692 element_init_same_as(pp->x, x);
693 element_init_same_as(pp->y, y);
694 element_set(pp->x, x);
695 element_set(pp->y, y);
696 for (i=0; i<rbits; i++) {
697 ellnet_pp_st_ptr seq = pp->seq[i];
698 element_init_same_as(seq->sm1, x);
699 element_init_same_as(seq->s0, x);
700 element_init_same_as(seq->s1, x);
701 element_init_same_as(seq->s2, x);
702 element_init_same_as(seq->tm1, x);
703 element_init_same_as(seq->t0, x);
704 element_init_same_as(seq->t1, x);
705 element_init_same_as(seq->t2, x);
706 }
707
708 //we map (x2,y2) to (-x2, i y2) before pairing
709 //notation: cmi means c_{k-i}, ci means c_{k+i}
710 element_t cm3, cm2, cm1, c0, c1, c2, c3, c4;
711 element_t C;
712
713 element_init_same_as(cm3, x);
714 element_init_same_as(cm2, x);
715 element_init_same_as(cm1, x);
716 element_init_same_as(c0, x);
717 element_init_same_as(c1, x);
718 element_init_same_as(c2, x);
719 element_init_same_as(c3, x);
720 element_init_same_as(c4, x);
721 element_init_same_as(C, x);
722
723 // c1 = 2y
724 // c0 = 1
725 // cm2 = -1
726 // cm3 = -2y
727 element_double(c1, y);
728 element_set1(c0);
729 element_neg(cm3, c1);
730 element_neg(cm2, c0);
731
732 // a = 1, b = 0 for Y^2 = X^3 + X
733 //hence c3 = c_{k+3} = c_4 = 4y(x^6 + 5(x^4 - x^2) - 1)
734 //use cm1, C, c2 as temp variables for now
735 element_square(cm1, x);
736 element_square(C, cm1);
737 element_sub(c2, C, cm1);
738 element_double(c3, c2);
739 element_double(c3, c3);
740 element_add(c3, c3, c2);
741 element_mul(c2, C, cm1);
742 element_add(c3, c3, c2);
743 element_add(c3, c3, cm2);
744 element_mul(c3, c3, c1);
745 element_double(c3, c3);
746
747 // c2 = c_3 = 3x^4 + 6x^2 - 1
748 element_double(cm1, cm1);
749 element_add(cm1, cm1, C);
750 element_double(C, cm1);
751 element_add(C, C, cm1);
752 element_add(c2, C, cm2);
753
754 // c4 = c_5 = c_2^3 c_4 - c_3^3 = c1^3 c3 - c2^3
755 element_square(C, c1);
756 element_mul(c4, C, c1);
757 element_mul(c4, c4, c3);
758 element_square(C, c2);
759 element_mul(C, C, c2);
760 element_sub(c4, c4, C);
761
762 // cm1 = 0
763 // C = (2y)^-1
764 element_set0(cm1);
765 element_invert(C, c1);
766
767 int k = 0;
768 element_t sm2, s3;
769 element_t tm2, t3;
770 element_ptr sm1, s0, s1, s2;
771 element_ptr tm1, t0, t1, t2;
772 element_t e0, e1;
773
774 element_init_same_as(sm2, x);
775 element_init_same_as(s3, x);
776
777 element_init_same_as(tm2, x);
778 element_init_same_as(t3, x);
779
780 element_init_same_as(e0, x);
781 element_init_same_as(e1, x);
782
783 int m = rbits - 2;
784 for (;;) {
785 ellnet_pp_st_ptr seq = pp->seq[k];
786 sm1 = seq->sm1;
787 s0 = seq->s0;
788 s1 = seq->s1;
789 s2 = seq->s2;
790 tm1 = seq->tm1;
791 t0 = seq->t0;
792 t1 = seq->t1;
793 t2 = seq->t2;
794
795 element_square(sm2, cm2);
796 element_square(sm1, cm1);
797 element_square(s0, c0);
798 element_square(s1, c1);
799 element_square(s2, c2);
800 element_square(s3, c3);
801
802 element_mul(tm2, cm3, cm1);
803 element_mul(tm1, cm2, c0);
804 element_mul(t0, cm1, c1);
805 element_mul(t1, c0, c2);
806 element_mul(t2, c1, c3);
807 element_mul(t3, c2, c4);
808
809 if (!m) break;
810 k++;
811
812 if (mpz_tstbit(pairing->r, m)) {
813 //double-and-add
814 element_mul(e0, t0, sm2);
815 element_mul(e1, tm2, s0);
816 element_sub(cm3, e0, e1);
817 element_mul(cm3, cm3, C);
818
819 element_mul(e0, t0, sm1);
820 element_mul(e1, tm1, s0);
821 element_sub(cm2, e0, e1);
822
823 element_mul(e0, t1, sm1);
824 element_mul(e1, tm1, s1);
825 element_sub(cm1, e0, e1);
826 element_mul(cm1, cm1, C);
827
828 element_mul(e0, t1, s0);
829 element_mul(e1, t0, s1);
830 element_sub(c0, e0, e1);
831
832 element_mul(e0, t2, s0);
833 element_mul(e1, t0, s2);
834 element_sub(c1, e0, e1);
835 element_mul(c1, c1, C);
836
837 element_mul(e0, t2, s1);
838 element_mul(e1, t1, s2);
839 element_sub(c2, e0, e1);
840
841 element_mul(e0, t3, s1);
842 element_mul(e1, t1, s3);
843 element_sub(c3, e0, e1);
844 element_mul(c3, c3, C);
845
846 element_mul(e0, t3, s2);
847 element_mul(e1, t2, s3);
848 element_sub(c4, e0, e1);
849
850 } else {
851 //double
852 element_mul(e0, tm1, sm2);
853 element_mul(e1, tm2, sm1);
854 element_sub(cm3, e0, e1);
855
856 element_mul(e0, t0, sm2);
857 element_mul(e1, tm2, s0);
858 element_sub(cm2, e0, e1);
859 element_mul(cm2, cm2, C);
860
861 element_mul(e0, t0, sm1);
862 element_mul(e1, tm1, s0);
863 element_sub(cm1, e0, e1);
864
865 element_mul(e0, t1, sm1);
866 element_mul(e1, tm1, s1);
867 element_sub(c0, e0, e1);
868 element_mul(c0, c0, C);
869
870 element_mul(e0, t1, s0);
871 element_mul(e1, t0, s1);
872 element_sub(c1, e0, e1);
873
874 element_mul(e0, t2, s0);
875 element_mul(e1, t0, s2);
876 element_sub(c2, e0, e1);
877 element_mul(c2, c2, C);
878
879 element_mul(e0, t2, s1);
880 element_mul(e1, t1, s2);
881 element_sub(c3, e0, e1);
882
883 element_mul(e0, t3, s1);
884 element_mul(e1, t1, s3);
885 element_sub(c4, e0, e1);
886 element_mul(c4, c4, C);
887 }
888 m--;
889 }
890
891 element_clear(cm3);
892 element_clear(cm2);
893 element_clear(cm1);
894 element_clear(c0);
895 element_clear(c1);
896 element_clear(c2);
897 element_clear(c3);
898 element_clear(c4);
899
900 element_clear(sm2);
901 element_clear(s3);
902
903 element_clear(tm2);
904 element_clear(t3);
905
906 element_clear(e0);
907 element_clear(e1);
908 element_clear(C);
909 }
910
a_pairing_ellnet_pp_clear(pairing_pp_t p)911 static void a_pairing_ellnet_pp_clear(pairing_pp_t p) {
912 ellnet_pp_ptr pp = p->data;
913 int i, rbits = mpz_sizeinbase(p->pairing->r, 2);
914 for (i=0; i<rbits; i++) {
915 ellnet_pp_st_ptr seq = pp->seq[i];
916 element_clear(seq->sm1);
917 element_clear(seq->s0);
918 element_clear(seq->s1);
919 element_clear(seq->s2);
920 element_clear(seq->tm1);
921 element_clear(seq->t0);
922 element_clear(seq->t1);
923 element_clear(seq->t2);
924 }
925 element_clear(pp->x);
926 element_clear(pp->y);
927 pbc_free(pp->seq);
928 pbc_free(p->data);
929 }
930
a_pairing_ellnet_pp_apply(element_ptr out,element_ptr in2,pairing_pp_t p)931 static void a_pairing_ellnet_pp_apply(element_ptr out, element_ptr in2, pairing_pp_t p) {
932 element_ptr x2 = curve_x_coord(in2);
933 element_ptr y2 = curve_y_coord(in2);
934 ellnet_pp_ptr pp = p->data;
935 int rbits = mpz_sizeinbase(p->pairing->r, 2);
936 int k = 0;
937 int m = rbits - 2;
938 element_t A, B;
939 element_t e0, e1;
940 element_t dm1, d0, d1;
941 element_t u, v;
942
943 element_init_same_as(A, x2);
944 element_init_same_as(B, out);
945 element_init_same_as(e0, x2);
946 element_init_same_as(e1, x2);
947 element_init_same_as(dm1, out);
948 element_init_same_as(d0, out);
949 element_init_same_as(d1, out);
950 element_init_same_as(u, out);
951 element_init_same_as(v, out);
952
953 element_add(A, pp->x, x2);
954 element_double(e0, pp->x);
955 element_sub(e0, e0, x2);
956 element_square(e1, A);
957 element_mul(e1, e0, e1);
958 element_set(element_x(d1), pp->y);
959 element_set(element_y(d1), y2);
960 element_square(d1, d1);
961 element_sub(element_x(d1), element_x(d1), e1);
962 element_neg(B, d1);
963 element_invert(B, B);
964 element_invert(A, A);
965 element_mul(element_x(d1), pp->y, A);
966 element_neg(element_x(d1), element_x(d1));
967 element_mul(element_y(d1), y2, A);
968 element_square(d1, d1);
969 element_sub(element_x(d1), e0, element_x(d1));
970 element_neg(element_y(d1), element_y(d1));
971
972 element_set1(dm1);
973 element_set1(d0);
974 for (;;) {
975 element_ptr sm1, s0, s1, s2;
976 element_ptr tm1, t0, t1, t2;
977 ellnet_pp_st_ptr seq = pp->seq[k];
978 sm1 = seq->sm1;
979 s0 = seq->s0;
980 s1 = seq->s1;
981 s2 = seq->s2;
982 tm1 = seq->tm1;
983 t0 = seq->t0;
984 t1 = seq->t1;
985 t2 = seq->t2;
986 k++;
987
988 element_square(u, d0);
989 element_mul(v, dm1, d1);
990
991 if (mpz_tstbit(p->pairing->r, m)) {
992 //double-and-add
993 element_mul(element_x(out), element_x(u), t0);
994 element_mul(element_y(out), element_y(u), t0);
995 element_mul(element_x(dm1), element_x(v), s0);
996 element_mul(element_y(dm1), element_y(v), s0);
997 element_sub(dm1, dm1, out);
998
999 element_mul(element_x(out), element_x(u), t1);
1000 element_mul(element_y(out), element_y(u), t1);
1001 element_mul(element_x(d0), element_x(v), s1);
1002 element_mul(element_y(d0), element_y(v), s1);
1003 element_sub(d0, d0, out);
1004 element_mul(element_x(d0), element_x(d0), A);
1005 element_mul(element_y(d0), element_y(d0), A);
1006
1007 element_mul(element_x(out), element_x(u), t2);
1008 element_mul(element_y(out), element_y(u), t2);
1009 element_mul(element_x(d1), element_x(v), s2);
1010 element_mul(element_y(d1), element_y(v), s2);
1011 element_sub(d1, d1, out);
1012 element_mul(d1, d1, B);
1013 } else {
1014 //double
1015 element_mul(element_x(out), element_x(u), tm1);
1016 element_mul(element_y(out), element_y(u), tm1);
1017 element_mul(element_x(dm1), element_x(v), sm1);
1018 element_mul(element_y(dm1), element_y(v), sm1);
1019 element_sub(dm1, dm1, out);
1020
1021 element_mul(element_x(out), element_x(u), t0);
1022 element_mul(element_y(out), element_y(u), t0);
1023 element_mul(element_x(d0), element_x(v), s0);
1024 element_mul(element_y(d0), element_y(v), s0);
1025 element_sub(d0, d0, out);
1026
1027 element_mul(element_x(out), element_x(u), t1);
1028 element_mul(element_y(out), element_y(u), t1);
1029 element_mul(element_x(d1), element_x(v), s1);
1030 element_mul(element_y(d1), element_y(v), s1);
1031 element_sub(d1, d1, out);
1032 element_mul(element_x(d1), element_x(d1), A);
1033 element_mul(element_y(d1), element_y(d1), A);
1034 }
1035 if (!m) break;
1036 m--;
1037 }
1038 a_tateexp(out, d1, d0, p->pairing->phikonr);
1039
1040 element_clear(A);
1041 element_clear(B);
1042 element_clear(e0);
1043 element_clear(e1);
1044 element_clear(dm1);
1045 element_clear(d0);
1046 element_clear(d1);
1047 element_clear(u);
1048 element_clear(v);
1049 }
1050
1051 //in1, in2 are from E(F_q), out from F_q^2
a_pairing_proj(element_ptr out,element_ptr in1,element_ptr in2,pairing_t pairing)1052 static void a_pairing_proj(element_ptr out, element_ptr in1, element_ptr in2,
1053 pairing_t pairing) {
1054 a_pairing_data_ptr p = pairing->data;
1055 element_t V, V1;
1056 element_t z, z2;
1057 element_t f, f0, f1;
1058 element_t a, b, c;
1059 element_t e0;
1060 const element_ptr e1 = a, e2 = b, e3 = c;
1061 int i, n;
1062 element_ptr Vx, Vy;
1063 element_ptr V1x, V1y;
1064 element_ptr Qx = curve_x_coord(in2);
1065 element_ptr Qy = curve_y_coord(in2);
1066
1067 //could save a couple of inversions by avoiding
1068 //this function and rewriting do_line() to handle projective coords
1069 //convert V from weighted projective (Jacobian) to affine
1070 //i.e. (X, Y, Z) --> (X/Z^2, Y/Z^3)
1071 //also sets z to 1
1072 #define point_to_affine() \
1073 element_invert(z, z); \
1074 element_square(e0, z); \
1075 element_mul(Vx, Vx, e0); \
1076 element_mul(e0, e0, z); \
1077 element_mul(Vy, Vy, e0); \
1078 element_set1(z); \
1079 element_set1(z2);
1080
1081 #define proj_double() { \
1082 /* e0 = 3x^2 + (cc->a) z^4 */ \
1083 /* for this case a = 1 */ \
1084 element_square(e0, Vx); \
1085 /*element_mul_si(e0, e0, 3);*/ \
1086 element_double(e1, e0); \
1087 element_add(e0, e1, e0); \
1088 element_square(e1, z2); \
1089 element_add(e0, e0, e1); \
1090 \
1091 /* z_out = 2 y z */ \
1092 element_mul(z, Vy, z); \
1093 /*element_mul_si(z, z, 2);*/ \
1094 element_double(z, z); \
1095 element_square(z2, z); \
1096 \
1097 /* e1 = 4 x y^2 */ \
1098 element_square(e2, Vy); \
1099 element_mul(e1, Vx, e2); \
1100 /*element_mul_si(e1, e1, 4);*/ \
1101 element_double(e1, e1); \
1102 element_double(e1, e1); \
1103 \
1104 /* x_out = e0^2 - 2 e1 */ \
1105 element_double(e3, e1); \
1106 element_square(Vx, e0); \
1107 element_sub(Vx, Vx, e3); \
1108 \
1109 /* e2 = 8y^4 */ \
1110 element_square(e2, e2); \
1111 /*element_mul_si(e2, e2, 8);*/ \
1112 element_double(e2, e2); \
1113 element_double(e2, e2); \
1114 element_double(e2, e2); \
1115 \
1116 /*y_out = e0(e1 - x_out) - e2*/\
1117 element_sub(e1, e1, Vx); \
1118 element_mul(e0, e0, e1); \
1119 element_sub(Vy, e0, e2); \
1120 }
1121
1122 #define do_tangent() \
1123 compute_abc_tangent_proj(a, b, c, Vx, Vy, z, z2, e0); \
1124 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
1125 element_mul(f, f, f0);
1126
1127 #define do_line() \
1128 compute_abc_line(a, b, c, Vx, Vy, V1x, V1y, e0); \
1129 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
1130 element_mul(f, f, f0);
1131
1132 element_init(V, p->Eq);
1133 element_init(V1, p->Eq);
1134 element_set(V, in1);
1135
1136 Vx = curve_x_coord(V);
1137 Vy = curve_y_coord(V);
1138 V1x = curve_x_coord(V1);
1139 V1y = curve_y_coord(V1);
1140
1141 element_init(f, p->Fq2);
1142 element_init(f0, p->Fq2);
1143 element_init(f1, p->Fq2);
1144 element_set1(f);
1145 element_init(a, p->Fq);
1146 element_init(b, p->Fq);
1147 element_init(c, p->Fq);
1148 element_init(e0, p->Fq);
1149 element_init(z, p->Fq);
1150 element_init(z2, p->Fq);
1151 element_set1(z);
1152 element_set1(z2);
1153 n = p->exp1;
1154 for (i=0; i<n; i++) {
1155 //f = f^2 g_V,V(Q)
1156 //where g_V,V = tangent at V
1157 element_square(f, f);
1158 do_tangent();
1159 proj_double();
1160 }
1161 point_to_affine();
1162 if (p->sign1 < 0) {
1163 element_neg(V1, V);
1164 element_invert(f1, f);
1165 } else {
1166 element_set(V1, V);
1167 element_set(f1, f);
1168 }
1169 n = p->exp2;
1170 for (; i<n; i++) {
1171 element_square(f, f);
1172 do_tangent();
1173 proj_double();
1174 }
1175
1176 element_mul(f, f, f1);
1177 point_to_affine();
1178 do_line();
1179
1180 a_tateexp(out, f, f0, pairing->phikonr);
1181
1182 element_clear(f);
1183 element_clear(f0);
1184 element_clear(f1);
1185 element_clear(z);
1186 element_clear(z2);
1187 element_clear(V);
1188 element_clear(V1);
1189 element_clear(a);
1190 element_clear(b);
1191 element_clear(c);
1192 element_clear(e0);
1193 #undef point_to_affine
1194 #undef proj_double
1195 #undef do_tangent
1196 #undef do_line
1197 }
1198
1199 //in1, in2 are from E(F_q), out from F_q^2
a_pairing_affine(element_ptr out,element_ptr in1,element_ptr in2,pairing_t pairing)1200 static void a_pairing_affine(element_ptr out, element_ptr in1, element_ptr in2,
1201 pairing_t pairing) {
1202 a_pairing_data_ptr p = pairing->data;
1203 element_t V, V1;
1204 element_t f, f0, f1;
1205 element_t a, b, c;
1206 element_t e0;
1207 int i, n;
1208 element_ptr Qx = curve_x_coord(in2);
1209 element_ptr Qy = curve_y_coord(in2);
1210 element_ptr Vx, Vy;
1211 element_ptr V1x, V1y;
1212
1213 #define do_tangent() \
1214 compute_abc_tangent(a, b, c, Vx, Vy, e0); \
1215 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
1216 element_mul(f, f, f0);
1217
1218 #define do_line() \
1219 compute_abc_line(a, b, c, Vx, Vy, V1x, V1y, e0); \
1220 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
1221 element_mul(f, f, f0);
1222
1223 element_init(V, p->Eq);
1224 element_init(V1, p->Eq);
1225 Vx = curve_x_coord(V);
1226 Vy = curve_y_coord(V);
1227
1228 V1x = curve_x_coord(V1);
1229 V1y = curve_y_coord(V1);
1230
1231 element_set(V, in1);
1232 element_init(f, p->Fq2);
1233 element_init(f0, p->Fq2);
1234 element_init(f1, p->Fq2);
1235 element_set1(f);
1236 element_init(a, p->Fq);
1237 element_init(b, p->Fq);
1238 element_init(c, p->Fq);
1239 element_init(e0, p->Fq);
1240 n = p->exp1;
1241 for (i=0; i<n; i++) {
1242 //f = f^2 g_V,V(Q)
1243 //where g_V,V = tangent at V
1244 element_square(f, f);
1245 do_tangent();
1246 element_double(V, V);
1247 }
1248 if (p->sign1 < 0) {
1249 element_neg(V1, V);
1250 element_invert(f1, f);
1251 } else {
1252 element_set(V1, V);
1253 element_set(f1, f);
1254 }
1255 n = p->exp2;
1256 for (; i<n; i++) {
1257 element_square(f, f);
1258 do_tangent();
1259 element_double(V, V);
1260 }
1261
1262 element_mul(f, f, f1);
1263 do_line();
1264
1265 a_tateexp(out, f, f0, pairing->phikonr);
1266
1267 element_clear(f);
1268 element_clear(f0);
1269 element_clear(f1);
1270 element_clear(V);
1271 element_clear(V1);
1272 element_clear(a);
1273 element_clear(b);
1274 element_clear(c);
1275 element_clear(e0);
1276 #undef do_tangent
1277 #undef do_line
1278 }
1279
1280 // On Computing Products of Pairing
1281 //in1, in2 are from E(F_q), out from F_q^2
a_pairings_affine(element_ptr out,element_t in1[],element_t in2[],int n_prod,pairing_t pairing)1282 void a_pairings_affine(element_ptr out, element_t in1[], element_t in2[],
1283 int n_prod, pairing_t pairing) {
1284 a_pairing_data_ptr p = pairing->data;
1285 element_t* V = pbc_malloc(sizeof(element_t)*n_prod);
1286 element_t* V1 = pbc_malloc(sizeof(element_t)*n_prod);
1287 element_t f, f0, f1;
1288 element_t a, b, c;
1289 element_t e0;
1290 int i, j, n;
1291 element_ptr Qx, Qy;
1292 element_ptr Vx, Vy;
1293 element_ptr V1x, V1y;
1294
1295 #define do_tangents() \
1296 for(j=0; j<n_prod; j++){ \
1297 Vx = curve_x_coord(V[j]); \
1298 Vy = curve_y_coord(V[j]); \
1299 Qx = curve_x_coord(in2[j]); \
1300 Qy = curve_y_coord(in2[j]); \
1301 \
1302 compute_abc_tangent(a, b, c, Vx, Vy, e0); \
1303 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
1304 element_mul(f, f, f0); \
1305 }
1306
1307 #define do_lines() \
1308 for(j=0;j<n_prod;j++){ \
1309 Vx = curve_x_coord(V[j]); \
1310 Vy = curve_y_coord(V[j]); \
1311 V1x = curve_x_coord(V1[j]); \
1312 V1y = curve_y_coord(V1[j]); \
1313 Qx = curve_x_coord(in2[j]); \
1314 Qy = curve_y_coord(in2[j]); \
1315 \
1316 compute_abc_line(a, b, c, Vx, Vy, V1x, V1y, e0); \
1317 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
1318 element_mul(f, f, f0); \
1319 }
1320
1321 for(i=0; i<n_prod; i++){
1322 element_init(V[i],p->Eq);
1323 element_init(V1[i],p->Eq);
1324 element_set(V[i],in1[i]);
1325 }
1326
1327
1328 element_init(f, p->Fq2);
1329 element_init(f0, p->Fq2);
1330 element_init(f1, p->Fq2);
1331 element_set1(f);
1332 element_init(a, p->Fq);
1333 element_init(b, p->Fq);
1334 element_init(c, p->Fq);
1335 element_init(e0, p->Fq);
1336 n = p->exp1;
1337 for (i=0; i<n; i++) {
1338 //f = f^2 g_V,V(Q)
1339 //where g_V,V = tangent at V
1340 element_square(f, f);
1341 do_tangents();
1342 element_multi_double(V, V, n_prod); //V_i = V_i + V_i for all i at one time.
1343 }
1344 if (p->sign1 < 0) {
1345 for(j=0; j<n_prod; j++){
1346 element_neg(V1[j], V[j]);
1347 }
1348 element_invert(f1, f);
1349 } else {
1350 for(j=0; j<n_prod; j++){
1351 element_set(V1[j], V[j]);
1352 }
1353 element_set(f1, f);
1354 }
1355 n = p->exp2;
1356 for (; i<n; i++) {
1357 element_square(f, f);
1358 do_tangents();
1359 element_multi_double(V, V, n_prod);
1360 }
1361
1362 element_mul(f, f, f1);
1363 do_lines();
1364
1365 a_tateexp(out, f, f0, pairing->phikonr);
1366
1367 element_clear(f);
1368 element_clear(f0);
1369 element_clear(f1);
1370 for(j=0;j<n_prod;j++){
1371 element_clear(V[j]);
1372 element_clear(V1[j]);
1373 }
1374 pbc_free(V);
1375 pbc_free(V1);
1376 element_clear(a);
1377 element_clear(b);
1378 element_clear(c);
1379 element_clear(e0);
1380 #undef do_tangents
1381 #undef do_lines
1382 }
1383
a_pairing_clear(pairing_t pairing)1384 static void a_pairing_clear(pairing_t pairing) {
1385 field_clear(pairing->GT);
1386
1387 a_pairing_data_ptr p = pairing->data;
1388 field_clear(p->Eq);
1389 field_clear(p->Fq);
1390 field_clear(p->Fq2);
1391 pbc_free(p);
1392
1393 mpz_clear(pairing->r);
1394 mpz_clear(pairing->phikonr);
1395 field_clear(pairing->Zr);
1396 }
1397
a_pairing_option_set(pairing_t pairing,char * key,char * value)1398 static void a_pairing_option_set(pairing_t pairing, char *key, char *value) {
1399 if (!strcmp(key, "method")) {
1400 if (!strcmp(value, "miller")) {
1401 pairing->map = a_pairing_proj;
1402 pairing->pp_init = a_pairing_pp_init;
1403 pairing->pp_clear = a_pairing_pp_clear;
1404 pairing->pp_apply = a_pairing_pp_apply;
1405 } else if (!strcmp(value, "miller-affine")) {
1406 pairing->map = a_pairing_affine;
1407 pairing->pp_init = a_pairing_pp_init;
1408 pairing->pp_clear = a_pairing_pp_clear;
1409 pairing->pp_apply = a_pairing_pp_apply;
1410 } else if (!strcmp(value, "shipsey-stange")) {
1411 pairing->map = a_pairing_ellnet;
1412 pairing->pp_init = a_pairing_ellnet_pp_init;
1413 pairing->pp_clear = a_pairing_ellnet_pp_clear;
1414 pairing->pp_apply = a_pairing_ellnet_pp_apply;
1415 }
1416 }
1417 }
1418
a_finalpow(element_t e)1419 static void a_finalpow(element_t e) {
1420 pairing_ptr pairing = e->field->pairing;
1421 element_t t0, t1;
1422 element_init_same_as(t0, e->data);
1423 element_init_same_as(t1, e->data);
1424 a_tateexp(t0, e->data, t1, pairing->phikonr);
1425 element_set(e->data, t0);
1426 element_clear(t0);
1427 element_clear(t1);
1428 }
1429
a_init_pairing(pairing_ptr pairing,void * data)1430 static void a_init_pairing(pairing_ptr pairing, void *data) {
1431 a_param_ptr param = data;
1432 element_t a, b;
1433 a_pairing_data_ptr p;
1434
1435 p = pairing->data = pbc_malloc(sizeof(*p));
1436 p->exp2 = param->exp2;
1437 p->exp1 = param->exp1;
1438 p->sign1 = param->sign1;
1439 mpz_init(pairing->r);
1440 mpz_set(pairing->r, param->r);
1441 field_init_fp(pairing->Zr, pairing->r);
1442 pairing->map = a_pairing_proj;
1443 pairing->prod_pairings = a_pairings_affine;
1444
1445 field_init_fp(p->Fq, param->q);
1446 element_init(a, p->Fq);
1447 element_init(b, p->Fq);
1448 element_set1(a);
1449 element_set0(b);
1450 field_init_curve_ab(p->Eq, a, b, pairing->r, param->h);
1451 element_clear(a);
1452 element_clear(b);
1453
1454 field_init_fi(p->Fq2, p->Fq);
1455
1456 //k=2, hence phi_k(q) = q + 1, phikonr = (q+1)/r
1457 mpz_init(pairing->phikonr);
1458 mpz_set(pairing->phikonr, param->h);
1459
1460 pairing->G1 = p->Eq;
1461 pairing->G2 = pairing->G1;
1462 pairing->phi = phi_identity;
1463 pairing_GT_init(pairing, p->Fq2);
1464 pairing->finalpow = a_finalpow;
1465
1466 pairing->clear_func = a_pairing_clear;
1467 pairing->option_set = a_pairing_option_set;
1468 pairing->pp_init = a_pairing_pp_init;
1469 pairing->pp_clear = a_pairing_pp_clear;
1470 pairing->pp_apply = a_pairing_pp_apply;
1471 }
1472
a_param_init(pbc_param_ptr par)1473 static void a_param_init(pbc_param_ptr par) {
1474 static pbc_param_interface_t interface = {{
1475 a_clear,
1476 a_init_pairing,
1477 a_out_str,
1478 }};
1479 par->api = interface;
1480 a_param_ptr p = par->data = pbc_malloc(sizeof(*p));
1481 mpz_init(p->r);
1482 mpz_init(p->q);
1483 mpz_init(p->h);
1484 }
1485
1486 // Public interface for type A pairings:
1487
pbc_param_init_a(pbc_param_ptr par,struct symtab_s * tab)1488 int pbc_param_init_a(pbc_param_ptr par, struct symtab_s *tab) {
1489 a_param_init(par);
1490 a_param_ptr p = par->data;
1491
1492 int err = 0;
1493 err += lookup_mpz(p->q, tab, "q");
1494 err += lookup_mpz(p->r, tab, "r");
1495 err += lookup_mpz(p->h, tab, "h");
1496 err += lookup_int(&p->exp2, tab, "exp2");
1497 err += lookup_int(&p->exp1, tab, "exp1");
1498 err += lookup_int(&p->sign1, tab, "sign1");
1499 err += lookup_int(&p->sign0, tab, "sign0");
1500 return err;
1501 }
1502
pbc_param_init_a_gen(pbc_param_ptr par,int rbits,int qbits)1503 void pbc_param_init_a_gen(pbc_param_ptr par, int rbits, int qbits) {
1504 a_param_init(par);
1505 a_param_ptr sp = par->data;
1506 int found = 0;
1507
1508 mpz_ptr q = sp->q;
1509 mpz_ptr r = sp->r;
1510 mpz_ptr h = sp->h;
1511
1512 do {
1513 int i;
1514 mpz_set_ui(r, 0);
1515
1516 if (rand() % 2) {
1517 sp->exp2 = rbits - 1;
1518 sp->sign1 = 1;
1519 } else {
1520 sp->exp2 = rbits;
1521 sp->sign1 = -1;
1522 }
1523 mpz_setbit(r, sp->exp2);
1524
1525 //use q as a temp variable
1526 mpz_set_ui(q, 0);
1527 sp->exp1 = (rand() % (sp->exp2 - 1)) + 1;
1528 mpz_setbit(q, sp->exp1);
1529 if (sp->sign1 > 0) {
1530 mpz_add(r, r, q);
1531 } else {
1532 mpz_sub(r, r, q);
1533 }
1534
1535 if (rand() % 2) {
1536 sp->sign0 = 1;
1537 mpz_add_ui(r, r, 1);
1538 } else {
1539 sp->sign0 = -1;
1540 mpz_sub_ui(r, r, 1);
1541 }
1542 if (!mpz_probab_prime_p(r, 10)) continue;
1543 for (i=0; i<10; i++) {
1544 int bit;
1545 //use q as a temp variable
1546 mpz_set_ui(q, 0);
1547 bit = qbits - rbits - 4 + 1;
1548 if (bit < 3) bit = 3;
1549 mpz_setbit(q, bit);
1550 pbc_mpz_random(h, q);
1551 mpz_mul_ui(h, h, 12);
1552 //finally q takes the value it should
1553 mpz_mul(q, h, r);
1554 mpz_sub_ui(q, q, 1);
1555 if (mpz_probab_prime_p(q, 10)) {
1556 found = 1;
1557 break;
1558 }
1559 }
1560 } while (!found);
1561 }
1562
1563 // Type A1 pairings:
1564
1565 struct a1_param_s {
1566 mpz_t p;
1567 mpz_t n;
1568 int l;
1569 };
1570 typedef struct a1_param_s a1_param_t[1];
1571 typedef struct a1_param_s *a1_param_ptr;
1572
1573 struct a1_pairing_data_s {
1574 field_t Fp, Fp2, Ep;
1575 };
1576 typedef struct a1_pairing_data_s a1_pairing_data_t[1];
1577 typedef struct a1_pairing_data_s *a1_pairing_data_ptr;
1578
a1_clear(void * data)1579 static void a1_clear(void *data) {
1580 a1_param_ptr param = data;
1581 mpz_clear(param->p);
1582 mpz_clear(param->n);
1583 pbc_free(data);
1584 }
1585
a1_out_str(FILE * stream,void * data)1586 static void a1_out_str(FILE *stream, void *data) {
1587 a1_param_ptr p = data;
1588 param_out_type(stream, "a1");
1589 param_out_mpz(stream, "p", p->p);
1590 param_out_mpz(stream, "n", p->n);
1591 param_out_int(stream, "l", p->l);
1592 }
1593
1594 struct pp2_coeff_s {
1595 element_t cx2;
1596 element_t cy2;
1597 element_t cxy;
1598 element_t cx;
1599 element_t cy;
1600 element_t c;
1601 };
1602 typedef struct pp2_coeff_s pp2_coeff_t[1];
1603 typedef struct pp2_coeff_s *pp2_coeff_ptr;
1604
pp2_coeff_set(pp2_coeff_ptr p,element_t cx2,element_t cy2,element_t cxy,element_t cx,element_t cy,element_t c)1605 static void pp2_coeff_set(pp2_coeff_ptr p,
1606 element_t cx2, element_t cy2, element_t cxy,
1607 element_t cx, element_t cy, element_t c) {
1608 element_init(p->cx2, cx2->field);
1609 element_init(p->cy2, cy2->field);
1610 element_init(p->cxy, cxy->field);
1611 element_init(p->cx, cx->field);
1612 element_init(p->cy, cy->field);
1613 element_init(p->c, c->field);
1614 element_set(p->cx2, cx2);
1615 element_set(p->cy2, cy2);
1616 element_set(p->cxy, cxy);
1617 element_set(p->cx, cx);
1618 element_set(p->cy, cy);
1619 element_set(p->c, c);
1620 }
1621
a1_pairing_pp_clear(pairing_pp_t p)1622 static void a1_pairing_pp_clear(pairing_pp_t p) {
1623 void **pp = p->data;
1624 while (*pp) {
1625 pbc_free(*pp);
1626 pp++;
1627 }
1628 pbc_free(p->data);
1629 }
1630
a1_pairing_pp_init(pairing_pp_t p,element_ptr in1,pairing_t pairing)1631 static void a1_pairing_pp_init(pairing_pp_t p, element_ptr in1, pairing_t pairing) {
1632 int m;
1633 element_ptr Px = curve_x_coord(in1);
1634 element_ptr Py = curve_y_coord(in1);
1635 a1_pairing_data_ptr a1info = pairing->data;
1636 p->data = pbc_malloc(sizeof(void *) * mpz_sizeinbase(pairing->r, 2));
1637 void **pp = p->data;
1638 element_t V;
1639 element_t a, b, c;
1640 element_t a2, b2, c2;
1641 element_t e0, e1, e2;
1642 element_ptr Vx, Vy;
1643
1644 #define do_tangent() compute_abc_tangent(a, b, c, Vx, Vy, e0);
1645
1646 #define do_line() compute_abc_line(a2, b2, c2, Vx, Vy, Px, Py, e0);
1647
1648 element_init(V, a1info->Ep);
1649 element_set(V, in1);
1650 Vx = curve_x_coord(V);
1651 Vy = curve_y_coord(V);
1652
1653 element_init(a, a1info->Fp);
1654 element_init(b, a1info->Fp);
1655 element_init(c, a1info->Fp);
1656 element_init(e0, a1info->Fp);
1657 element_init(e1, a1info->Fp);
1658 element_init(e2, a1info->Fp);
1659 element_init(a2, a1info->Fp);
1660 element_init(b2, a1info->Fp);
1661 element_init(c2, a1info->Fp);
1662
1663 m = mpz_sizeinbase(pairing->r, 2) - 2;
1664
1665 for(;;) {
1666 do_tangent();
1667 if (!m) break;
1668 element_double(V, V);
1669
1670 if (mpz_tstbit(pairing->r, m)) {
1671 do_line();
1672 element_add(V, V, in1);
1673 //preprocess two at once
1674 //e0 = coeff of x
1675 element_mul(e0, a, c2);
1676 element_mul(e1, a2, c);
1677 element_add(e0, e0, e1);
1678
1679 //e1 = coeff of y
1680 element_mul(e1, b2, c);
1681 element_mul(e2, b, c2);
1682 element_add(e1, e1, e2);
1683
1684 //c = constant term
1685 element_mul(c, c, c2);
1686
1687 //c2 = coeff of xy
1688 element_mul(c2, a, b2);
1689 element_mul(e2, a2, b);
1690 element_add(c2, c2, e2);
1691
1692 //a = coeff of x^2
1693 element_mul(a, a, a2);
1694
1695 //b = coeff of y^2
1696 element_mul(b, b, b2);
1697
1698 *pp = pbc_malloc(sizeof(pp2_coeff_t));
1699 pp2_coeff_set(*pp, a, b, c2, e0, e1, c);
1700 } else {
1701 *pp = pbc_malloc(sizeof(pp_coeff_t));
1702 pp_coeff_set(*pp, a, b, c);
1703 }
1704 pp++;
1705 m--;
1706 }
1707 *pp = pbc_malloc(sizeof(pp_coeff_t));
1708 pp_coeff_set(*pp, a, b, c);
1709 pp++;
1710 *pp = NULL;
1711
1712 element_clear(a2);
1713 element_clear(b2);
1714 element_clear(c2);
1715 element_clear(e2);
1716 element_clear(e1);
1717 element_clear(e0);
1718 element_clear(a);
1719 element_clear(b);
1720 element_clear(c);
1721 element_clear(V);
1722 #undef do_tangent
1723 #undef do_line
1724 }
1725
a1_pairing_pp_apply(element_ptr out,element_ptr in2,pairing_pp_t p)1726 static void a1_pairing_pp_apply(element_ptr out, element_ptr in2, pairing_pp_t p) {
1727 void **pp = p->data;
1728 a1_pairing_data_ptr a1info = p->pairing->data;
1729 element_t f, f0;
1730 element_t e0, e1;
1731 int m;
1732 element_ptr Qx = curve_x_coord(in2);
1733 element_ptr Qy = curve_y_coord(in2);
1734 element_t Qx2, Qy2, Qxy;
1735
1736 #define do_tangent() \
1737 pp_coeff_ptr ppp = *pp; \
1738 a_miller_evalfn(f0, ppp->a, ppp->b, ppp->c, Qx, Qy);
1739
1740 #define do_line() { \
1741 pp2_coeff_ptr ppp = *pp; \
1742 /*we'll map Q via (x,y) --> (-x, iy) */ \
1743 /*hence Qx^2 = x^2, Qy^2 = -y^2, Qx Qy = -ixy */\
1744 /*where x = Q'x, y = Q'y */ \
1745 \
1746 /* Re = cx2 x^2 - cy2 y^2 - cx x + c */ \
1747 /* Im = -cxy xy + cy y */ \
1748 element_mul(e0, ppp->cx2, Qx2); \
1749 element_mul(e1, ppp->cy2, Qy2); \
1750 element_sub(e0, e0, e1); \
1751 element_mul(e1, ppp->cx, Qx); \
1752 element_sub(e0, e0, e1); \
1753 element_add(element_x(f0), e0, ppp->c); \
1754 \
1755 element_mul(e0, ppp->cy, Qy); \
1756 element_mul(e1, ppp->cxy, Qxy); \
1757 element_sub(element_y(f0), e0, e1); \
1758 }
1759
1760 element_init(f, out->field);
1761 element_init(f0, out->field);
1762
1763 element_set1(f);
1764
1765 element_init(e0, a1info->Fp);
1766 element_init(e1, a1info->Fp);
1767 element_init(Qx2, a1info->Fp);
1768 element_init(Qy2, a1info->Fp);
1769 element_init(Qxy, a1info->Fp);
1770
1771 element_square(Qx2, Qx);
1772 element_square(Qy2, Qy);
1773 element_mul(Qxy, Qx, Qy);
1774
1775 m = mpz_sizeinbase(p->pairing->r, 2) - 2;
1776
1777 while (m > 0) {
1778 if (mpz_tstbit(p->pairing->r, m)) {
1779 do_line();
1780 } else {
1781 do_tangent();
1782 }
1783 element_mul(f, f, f0);
1784 pp++;
1785 m--;
1786 element_square(f, f);
1787 }
1788 do_tangent();
1789 element_mul(f, f, f0);
1790
1791 //Tate exponentiation
1792 //simpler but slower:
1793 //element_pow_mpz(out, f, p->tateexp);
1794 //use this trick instead:
1795 element_invert(f0, f);
1796 element_neg(element_y(f), element_y(f));
1797 element_mul(f, f, f0);
1798 element_pow_mpz(out, f, p->pairing->phikonr);
1799
1800 /* We could use this instead but p->h is small so this does not help much
1801 a_tateexp(out, f, f0, p->h);
1802 */
1803
1804 element_clear(Qx2);
1805 element_clear(Qy2);
1806 element_clear(Qxy);
1807 element_clear(f);
1808 element_clear(f0);
1809 element_clear(e1);
1810 element_clear(e0);
1811 #undef do_tangent
1812 #undef do_line
1813 }
1814
1815 // e0 is a temp var.
1816 // Mixed coordinates.
compute_abc_line_proj(element_ptr a,element_ptr b,element_ptr c,element_ptr Vx,element_ptr Vy,element_ptr z,element_ptr z2,element_ptr V1x,element_ptr V1y,element_ptr e0)1817 static void compute_abc_line_proj(element_ptr a, element_ptr b, element_ptr c,
1818 element_ptr Vx, element_ptr Vy, element_ptr z, element_ptr z2,
1819 element_ptr V1x, element_ptr V1y, element_ptr e0) {
1820 //temporally used to store Z1^3
1821 element_mul(c,z,z2);
1822 //a = Y1-Y2*Z1^3
1823 element_mul(e0,V1y,c);
1824 element_sub(a,Vy,e0);
1825 //b = -(X1*Z1-X2*Z1^3)
1826 element_mul(b,c,V1x);
1827 element_mul(e0,Vx,z);
1828 element_sub(b,b,e0);
1829 //c = -(Y2*b+X2*a)
1830 element_mul(c,b,V1y);
1831 element_mul(e0,a,V1x);
1832 element_add(c,c,e0);
1833 element_neg(c,c);
1834 }
1835
1836 // in1, in2 are from E(F_q), out from F_q^2
a1_pairing_proj(element_ptr out,element_ptr in1,element_ptr in2,pairing_t pairing)1837 static void a1_pairing_proj(element_ptr out, element_ptr in1, element_ptr in2,
1838 pairing_t pairing) {
1839 a1_pairing_data_ptr p = pairing->data;
1840 element_t V;
1841 element_t z, z2;
1842 element_t f, f0;
1843 element_t a, b, c;
1844 element_t e0;
1845 const element_ptr e1 = a, e2 = b, e3 = c; // used in point_to_affine() etc.
1846 int m;
1847 element_ptr Px = curve_x_coord(in1);
1848 element_ptr Py = curve_y_coord(in1);
1849 element_ptr Qx = curve_x_coord(in2);
1850 element_ptr Qy = curve_y_coord(in2);
1851 element_ptr Vx;
1852 element_ptr Vy;
1853
1854 #define point_to_affine() \
1855 element_invert(z, z); \
1856 element_square(e0, z); \
1857 element_mul(Vx, Vx, e0); \
1858 element_mul(e0, e0, z); \
1859 element_mul(Vy, Vy, e0); \
1860 element_set1(z); \
1861 element_set1(z2);
1862
1863 //TODO: do I need to check if V=-in1?
1864 //Where V=(Vx,Vy,z) and in1=(Px,Py,1), a mixed coordinates.
1865 #define proj_add() { \
1866 /* H=X2*Z1^2-X1 */ \
1867 element_mul(e0,Px,z2); \
1868 element_sub(e0,e0,Vx); \
1869 /* H^2 */ \
1870 element_square(e1,e0); \
1871 /* r=Y2*Z1^3-Y1 */ \
1872 element_mul(e2,z,z2); \
1873 element_mul(e2,e2,Py); \
1874 element_sub(e2,e2,Vy); \
1875 \
1876 /* X3=r^2-H^3-2X1*H^2 */ \
1877 element_set(z2,Vx); /* use z2 to store X1 and update Vx=X3 */ \
1878 element_square(Vx,e2); \
1879 element_mul(e3,e0,e1); /* e3=H^3 */ \
1880 element_sub(Vx,Vx,e3); \
1881 element_double(e3,z2); \
1882 element_mul(e3,e3,e1); /* 2X1*H^2 */ \
1883 element_sub(Vx,Vx,e3); \
1884 /* Y3=r(X1*H^2-X3)-Y1*H^3 */ \
1885 element_mul(e3,z2,e1); \
1886 element_sub(e3,e3,Vx); \
1887 element_mul(e3,e3,e2); \
1888 element_mul(e2,e0,e1); /* e2 no longer used. */ \
1889 element_mul(e2,e2,Vy); \
1890 element_sub(Vy,e3,e2); \
1891 /* Z3=Z1*H */ \
1892 element_mul(z,z,e0); \
1893 element_square(z2,z); \
1894 }
1895
1896 #define proj_double() { \
1897 /* e0 = 3x^2 + (cc->a) z^4 */ \
1898 /* for this case a = 1 */ \
1899 element_square(e0, Vx); \
1900 /* element_mul_si(e0, e0, 3); */ \
1901 element_double(e1, e0); \
1902 element_add(e0, e1, e0); \
1903 element_square(e1, z2); \
1904 element_add(e0, e0, e1); \
1905 \
1906 /* z_out = 2 y z */ \
1907 element_mul(z, Vy, z); \
1908 /* element_mul_si(z, z, 2); */ \
1909 element_double(z, z); \
1910 element_square(z2, z); \
1911 \
1912 /* e1 = 4 x y^2 */ \
1913 element_square(e2, Vy); \
1914 element_mul(e1, Vx, e2); \
1915 /* element_mul_si(e1, e1, 4); */ \
1916 element_double(e1, e1); \
1917 element_double(e1, e1); \
1918 \
1919 /* x_out = e0^2 - 2 e1 */ \
1920 element_double(e3, e1); \
1921 element_square(Vx, e0); \
1922 element_sub(Vx, Vx, e3); \
1923 \
1924 /* e2 = 8y^4 */ \
1925 element_square(e2, e2); \
1926 /* element_mul_si(e2, e2, 8); */ \
1927 element_double(e2, e2); \
1928 element_double(e2, e2); \
1929 element_double(e2, e2); \
1930 \
1931 /* y_out = e0(e1 - x_out) - e2 */ \
1932 element_sub(e1, e1, Vx); \
1933 element_mul(e0, e0, e1); \
1934 element_sub(Vy, e0, e2); \
1935 }
1936
1937 #define do_tangent() { \
1938 compute_abc_tangent_proj(a, b, c, Vx, Vy, z, z2, e0); \
1939 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
1940 element_mul(f, f, f0); \
1941 }
1942
1943 #define do_line() { \
1944 compute_abc_line_proj(a, b, c, Vx, Vy, z, z2, Px, Py, e0); \
1945 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
1946 element_mul(f, f, f0); \
1947 }
1948
1949 element_init(V, p->Ep);
1950 element_set(V, in1);
1951 Vx = curve_x_coord(V);
1952 Vy = curve_y_coord(V);
1953
1954 element_init(f, p->Fp2);
1955 element_init(f0, p->Fp2);
1956 element_set1(f);
1957 element_init(a, p->Fp);
1958 element_init(b, p->Fp);
1959 element_init(c, p->Fp);
1960 element_init(e0, p->Fp);
1961 element_init(z, p->Fp);
1962 element_init(z2, p->Fp);
1963 element_set1(z);
1964 element_set1(z2);
1965
1966 m = mpz_sizeinbase(pairing->r, 2) - 2;
1967 //TODO: sliding NAF
1968 for(;;) {
1969 do_tangent();
1970 if (!m) break;
1971
1972 proj_double(); //V=2V
1973 if (mpz_tstbit(pairing->r, m)) {
1974 // point_to_affine();
1975 do_line();
1976 proj_add(); //V=V+in1
1977 }
1978
1979 m--;
1980 element_square(f, f);
1981 }
1982
1983 // Tate exponentiation.
1984 // Simpler but slower:
1985 // element_pow_mpz(out, f, p->tateexp);
1986 // Use this trick instead:
1987 element_invert(f0, f);
1988 element_neg(element_y(f), element_y(f));
1989 element_mul(f, f, f0);
1990 element_pow_mpz(out, f, pairing->phikonr);
1991
1992 /* We could use this instead but p->h is small so this does not help much
1993 a_tateexp(out, f, f0, p->h);
1994 */
1995
1996 element_clear(f);
1997 element_clear(f0);
1998 element_clear(z);
1999 element_clear(z2);
2000 element_clear(V);
2001 element_clear(a);
2002 element_clear(b);
2003 element_clear(c);
2004 element_clear(e0);
2005 #undef point_to_affine
2006 #undef proj_add
2007 #undef proj_double
2008 #undef do_tangent
2009 #undef do_line
2010 }
2011
2012 //in1, in2 are from E(F_q), out from F_q^2
a1_pairing(element_ptr out,element_ptr in1,element_ptr in2,pairing_t pairing)2013 static void a1_pairing(element_ptr out, element_ptr in1, element_ptr in2,
2014 pairing_t pairing) {
2015 a1_pairing_data_ptr p = pairing->data;
2016 element_t V;
2017 element_t f, f0;
2018 element_t a, b, c;
2019 element_t e0;
2020 int m;
2021 element_ptr Px = curve_x_coord(in1);
2022 element_ptr Py = curve_y_coord(in1);
2023 element_ptr Qx = curve_x_coord(in2);
2024 element_ptr Qy = curve_y_coord(in2);
2025 element_ptr Vx;
2026 element_ptr Vy;
2027
2028 #define do_tangent() { \
2029 compute_abc_tangent(a, b, c, Vx, Vy, e0); \
2030 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
2031 element_mul(f, f, f0); \
2032 }
2033
2034 #define do_line() { \
2035 compute_abc_line(a, b, c, Vx, Vy, Px, Py, e0); \
2036 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
2037 element_mul(f, f, f0); \
2038 }
2039
2040 element_init(V, p->Ep);
2041 element_set(V, in1);
2042 Vx = curve_x_coord(V);
2043 Vy = curve_y_coord(V);
2044
2045 element_init(f, p->Fp2);
2046 element_init(f0, p->Fp2);
2047 element_set1(f);
2048 element_init(a, p->Fp);
2049 element_init(b, p->Fp);
2050 element_init(c, p->Fp);
2051 element_init(e0, p->Fp);
2052
2053 m = mpz_sizeinbase(pairing->r, 2) - 2;
2054
2055 //TODO: sliding NAF
2056 for(;;) {
2057 do_tangent();
2058 if (!m) break;
2059
2060 element_double(V, V);
2061 if (mpz_tstbit(pairing->r, m)) {
2062 do_line();
2063 element_add(V, V, in1);
2064 }
2065
2066 m--;
2067 element_square(f, f);
2068 }
2069
2070 // Tate exponentiation.
2071 // Simpler but slower:
2072 // element_pow_mpz(out, f, p->tateexp);
2073 // Use this trick instead:
2074 element_invert(f0, f);
2075 element_neg(element_y(f), element_y(f));
2076 element_mul(f, f, f0);
2077 element_pow_mpz(out, f, pairing->phikonr);
2078
2079 /* We could use this instead but p->h is small so this does not help much
2080 a_tateexp(out, f, f0, p->h);
2081 */
2082
2083 element_clear(f);
2084 element_clear(f0);
2085 element_clear(V);
2086 element_clear(a);
2087 element_clear(b);
2088 element_clear(c);
2089 element_clear(e0);
2090 #undef do_tangent
2091 #undef do_line
2092 }
2093
2094 //in1, in2 are from E(F_q), out from F_q^2
a1_pairings_affine(element_ptr out,element_t in1[],element_t in2[],int n_prod,pairing_t pairing)2095 void a1_pairings_affine(element_ptr out, element_t in1[], element_t in2[],
2096 int n_prod, pairing_t pairing) {
2097 a1_pairing_data_ptr p = pairing->data;
2098 element_t* V = pbc_malloc(sizeof(element_t)*n_prod);
2099 element_t f, f0;
2100 element_t a, b, c;
2101 element_t e0;
2102 int m, i;
2103 element_ptr Px, Py;
2104 element_ptr Qx, Qy;
2105 element_ptr Vx, Vy;
2106
2107 #define do_tangents() { \
2108 for(i=0; i<n_prod; i++){ \
2109 Vx = curve_x_coord(V[i]); \
2110 Vy = curve_y_coord(V[i]); \
2111 Qx = curve_x_coord(in2[i]); \
2112 Qy = curve_y_coord(in2[i]); \
2113 compute_abc_tangent(a, b, c, Vx, Vy, e0); \
2114 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
2115 element_mul(f, f, f0); \
2116 } \
2117 }
2118
2119 #define do_lines() { \
2120 for(i=0; i<n_prod; i++){ \
2121 Vx = curve_x_coord(V[i]); \
2122 Vy = curve_y_coord(V[i]); \
2123 Px = curve_x_coord(in1[i]); \
2124 Py = curve_y_coord(in1[i]); \
2125 Qx = curve_x_coord(in2[i]); \
2126 Qy = curve_y_coord(in2[i]); \
2127 compute_abc_line(a, b, c, Vx, Vy, Px, Py, e0); \
2128 a_miller_evalfn(f0, a, b, c, Qx, Qy); \
2129 element_mul(f, f, f0); \
2130 } \
2131 }
2132
2133 for(i=0; i<n_prod; i++){
2134 element_init(V[i], p->Ep);
2135 element_set(V[i], in1[i]);
2136 }
2137 element_init(f, p->Fp2);
2138 element_init(f0, p->Fp2);
2139 element_set1(f);
2140 element_init(a, p->Fp);
2141 element_init(b, p->Fp);
2142 element_init(c, p->Fp);
2143 element_init(e0, p->Fp);
2144
2145 m = mpz_sizeinbase(pairing->r, 2) - 2;
2146
2147 //TODO: sliding NAF
2148 for(;;) {
2149 do_tangents();
2150 if (!m) break;
2151 element_multi_double(V, V, n_prod);
2152 if (mpz_tstbit(pairing->r, m)) {
2153 do_lines();
2154 element_multi_add(V, V, in1, n_prod);
2155 }
2156
2157 m--;
2158 element_square(f, f);
2159 }
2160
2161 // Tate exponentiation.
2162 // Simpler but slower:
2163 // element_pow_mpz(out, f, p->tateexp);
2164 // Use this trick instead:
2165 element_invert(f0, f);
2166 element_neg(element_y(f), element_y(f));
2167 element_mul(f, f, f0);
2168 element_pow_mpz(out, f, pairing->phikonr);
2169
2170 /* We could use this instead but p->h is small so this does not help much
2171 a_tateexp(out, f, f0, p->h);
2172 */
2173
2174 element_clear(f);
2175 element_clear(f0);
2176 for(i=0; i<n_prod; i++){
2177 element_clear(V[i]);
2178 }
2179 pbc_free(V);
2180 element_clear(a);
2181 element_clear(b);
2182 element_clear(c);
2183 element_clear(e0);
2184 #undef do_tangents
2185 #undef do_lines
2186 }
2187
a1_pairing_clear(pairing_t pairing)2188 static void a1_pairing_clear(pairing_t pairing) {
2189 field_clear(pairing->GT);
2190
2191 a1_pairing_data_ptr p = pairing->data;
2192 field_clear(p->Ep);
2193 field_clear(p->Fp2);
2194 field_clear(p->Fp);
2195 pbc_free(p);
2196
2197 mpz_clear(pairing->phikonr);
2198 mpz_clear(pairing->r);
2199 field_clear(pairing->Zr);
2200 }
2201
a1_pairing_option_set(pairing_t pairing,char * key,char * value)2202 static void a1_pairing_option_set(pairing_t pairing, char *key, char *value) {
2203 if (!strcmp(key, "method")) {
2204 if (!strcmp(value, "miller")) {
2205 pairing->map = a1_pairing_proj;
2206 pairing->pp_init = a1_pairing_pp_init;
2207 pairing->pp_clear = a1_pairing_pp_clear;
2208 pairing->pp_apply = a1_pairing_pp_apply;
2209 } else if (!strcmp(value, "miller-affine")){
2210 pairing->map = a1_pairing;
2211 pairing->pp_init = a1_pairing_pp_init;
2212 pairing->pp_clear = a1_pairing_pp_clear;
2213 pairing->pp_apply = a1_pairing_pp_apply;
2214 } else if (!strcmp(value, "shipsey-stange")) {
2215 pairing->map = a_pairing_ellnet;
2216 pairing->pp_init = a_pairing_ellnet_pp_init;
2217 pairing->pp_clear = a_pairing_ellnet_pp_clear;
2218 pairing->pp_apply = a_pairing_ellnet_pp_apply;
2219 }
2220 }
2221 }
2222
a1_init_pairing(pairing_t pairing,void * data)2223 static void a1_init_pairing(pairing_t pairing, void *data) {
2224 a1_param_ptr param = data;
2225 element_t a, b;
2226 mpz_init(pairing->r);
2227 mpz_set(pairing->r, param->n);
2228 field_init_fp(pairing->Zr, pairing->r);
2229
2230 a1_pairing_data_ptr p;
2231
2232 p = pairing->data = pbc_malloc(sizeof(a1_pairing_data_t));
2233
2234 //k=2, hence phi_k(q) = q + 1, phikonr = (q+1)/r
2235 mpz_init(pairing->phikonr);
2236 mpz_set_ui(pairing->phikonr, param->l);
2237
2238 field_init_fp(p->Fp, param->p);
2239 element_init(a, p->Fp);
2240 element_init(b, p->Fp);
2241 element_set1(a);
2242 element_set0(b);
2243 field_init_curve_ab(p->Ep, a, b, pairing->r, pairing->phikonr);
2244
2245 // Turns out to be faster.
2246 field_curve_use_random_solvefory(p->Ep);
2247
2248 element_clear(a);
2249 element_clear(b);
2250 field_init_fi(p->Fp2, p->Fp);
2251
2252 pairing->finalpow = a_finalpow;
2253 pairing->G1 = pbc_malloc(sizeof(field_t));
2254 pairing->G2 = pairing->G1 = p->Ep;
2255 pairing_GT_init(pairing, p->Fp2);
2256
2257 pairing->map = a1_pairing_proj; //default uses projective coordinates.
2258 pairing->phi = phi_identity;
2259 pairing->prod_pairings = a1_pairings_affine;
2260
2261 pairing->clear_func = a1_pairing_clear;
2262
2263 pairing->pp_init = a1_pairing_pp_init;
2264 pairing->pp_clear = a1_pairing_pp_clear;
2265 pairing->pp_apply = a1_pairing_pp_apply;
2266 pairing->option_set = a1_pairing_option_set;
2267 }
2268
a1_init(pbc_param_t p)2269 static void a1_init(pbc_param_t p) {
2270 static pbc_param_interface_t interface = {{
2271 a1_clear,
2272 a1_init_pairing,
2273 a1_out_str,
2274 }};
2275 p->api = interface;
2276 a1_param_ptr param = p->data = pbc_malloc(sizeof(*param));
2277 mpz_init(param->p);
2278 mpz_init(param->n);
2279 }
2280
2281 // Public interface:
2282
pbc_param_init_a1(pbc_param_ptr par,struct symtab_s * tab)2283 int pbc_param_init_a1(pbc_param_ptr par, struct symtab_s *tab) {
2284 a1_init(par);
2285 a1_param_ptr p = par->data;
2286
2287 int err = 0;
2288 err += lookup_mpz(p->p, tab, "p");
2289 err += lookup_mpz(p->n, tab, "n");
2290 err += lookup_int(&p->l, tab, "l");
2291 return err;
2292 }
2293
pbc_param_init_a1_gen(pbc_param_ptr par,mpz_t order)2294 void pbc_param_init_a1_gen(pbc_param_ptr par, mpz_t order) {
2295 a1_init(par);
2296 a1_param_ptr param = par->data;
2297 // If order is even, ideally check all even l, not just multiples of 4
2298 // but I don't see a good reason for having an even order.
2299 unsigned int l = 4;
2300 mpz_t n;
2301 mpz_ptr p = param->p;
2302 mpz_init(n);
2303 mpz_mul_ui(n, order, 4);
2304 mpz_sub_ui(p, n, 1);
2305 for (;;) {
2306 if (mpz_probab_prime_p(p, 20)) {
2307 break;
2308 }
2309 mpz_add(p, p, n);
2310 l += 4;
2311 }
2312 param->l = l;
2313 mpz_set(param->n, order);
2314 mpz_clear(n);
2315 }
2316