1 /* Copyright (C) 2012 The PARI group.
2
3 This file is part of the PARI/GP package.
4
5 PARI/GP is free software; you can redistribute it and/or modify it under the
6 terms of the GNU General Public License as published by the Free Software
7 Foundation; either version 2 of the License, or (at your option) any later
8 version. It is distributed in the hope that it will be useful, but WITHOUT
9 ANY WARRANTY WHATSOEVER.
10
11 Check the License for details. You should have received a copy of it, along
12 with the package; see the file 'COPYING'. If not, write to the Free Software
13 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
14
15 #include "pari.h"
16 #include "paripriv.h"
17
18 /* Not so fast arithmetic with points over elliptic curves over F_2^n */
19
20 /***********************************************************************/
21 /** **/
22 /** F2xqE **/
23 /** **/
24 /***********************************************************************/
25
26 /* Theses functions deal with point over elliptic curves over F_2^n defined
27 * by an equation of the form:
28 ** y^2+x*y=x^3+a_2*x^2+a_6 if the curve is ordinary.
29 ** y^2+a_3*y=x^3+a_4*x+a_6 if the curve is supersingular.
30 * Most of the time a6 is omitted since it can be recovered from any point
31 * on the curve.
32 * For supersingular curves, the parameter a2 is replaced by [a3,a4,a3^-1].
33 */
34
35 GEN
RgE_to_F2xqE(GEN x,GEN T)36 RgE_to_F2xqE(GEN x, GEN T)
37 {
38 if (ell_is_inf(x)) return x;
39 retmkvec2(Rg_to_F2xq(gel(x,1),T),Rg_to_F2xq(gel(x,2),T));
40 }
41
42 GEN
F2xqE_changepoint(GEN x,GEN ch,GEN T)43 F2xqE_changepoint(GEN x, GEN ch, GEN T)
44 {
45 pari_sp av = avma;
46 GEN p1,z,u,r,s,t,v,v2,v3;
47 if (ell_is_inf(x)) return x;
48 u = gel(ch,1); r = gel(ch,2);
49 s = gel(ch,3); t = gel(ch,4);
50 v = F2xq_inv(u, T); v2 = F2xq_sqr(v, T); v3 = F2xq_mul(v,v2, T);
51 p1 = F2x_add(gel(x,1),r);
52 z = cgetg(3,t_VEC);
53 gel(z,1) = F2xq_mul(v2, p1, T);
54 gel(z,2) = F2xq_mul(v3, F2x_add(gel(x,2), F2x_add(F2xq_mul(s, p1, T),t)), T);
55 return gerepileupto(av, z);
56 }
57
58 GEN
F2xqE_changepointinv(GEN x,GEN ch,GEN T)59 F2xqE_changepointinv(GEN x, GEN ch, GEN T)
60 {
61 GEN u, r, s, t, X, Y, u2, u3, u2X, z;
62 if (ell_is_inf(x)) return x;
63 X = gel(x,1); Y = gel(x,2);
64 u = gel(ch,1); r = gel(ch,2);
65 s = gel(ch,3); t = gel(ch,4);
66 u2 = F2xq_sqr(u, T); u3 = F2xq_mul(u,u2, T);
67 u2X = F2xq_mul(u2,X, T);
68 z = cgetg(3, t_VEC);
69 gel(z,1) = F2x_add(u2X,r);
70 gel(z,2) = F2x_add(F2xq_mul(u3,Y, T), F2x_add(F2xq_mul(s,u2X, T), t));
71 return z;
72 }
73
74 static GEN
nonzerotrace_F2xq(GEN T)75 nonzerotrace_F2xq(GEN T)
76 {
77 pari_sp av = avma;
78 long n = F2x_degree(T), vs = T[1];
79 GEN a;
80 if (odd(n))
81 return pol1_F2x(vs);
82 do
83 {
84 set_avma(av);
85 a = random_F2x(n, vs);
86 } while (F2xq_trace(a, T)==0);
87 return a;
88 }
89
90 void
F2xq_elltwist(GEN a,GEN a6,GEN T,GEN * pt_a,GEN * pt_a6)91 F2xq_elltwist(GEN a, GEN a6, GEN T, GEN *pt_a, GEN *pt_a6)
92 {
93 pari_sp av = avma;
94 GEN n = nonzerotrace_F2xq(T);
95 if (typ(a)==t_VECSMALL)
96 {
97 *pt_a = gerepileuptoleaf(av, F2x_add(n, a));
98 *pt_a6 = vecsmall_copy(a6);
99 } else
100 {
101 GEN a6t = F2x_add(a6, F2xq_mul(n, F2xq_sqr(gel(a,1), T), T));
102 *pt_a6 = gerepileuptoleaf(av, a6t);
103 *pt_a = vecsmall_copy(a);
104 }
105 }
106
107 static GEN
F2xqE_dbl_slope(GEN P,GEN a,GEN T,GEN * slope)108 F2xqE_dbl_slope(GEN P, GEN a, GEN T, GEN *slope)
109 {
110 GEN x, y, Q;
111 if (ell_is_inf(P)) return ellinf();
112 x = gel(P,1); y = gel(P,2);
113 if (typ(a)==t_VECSMALL)
114 {
115 GEN a2 = a;
116 if (!lgpol(gel(P,1))) return ellinf();
117 *slope = F2x_add(x, F2xq_div(y, x, T));
118 Q = cgetg(3,t_VEC);
119 gel(Q, 1) = F2x_add(F2xq_sqr(*slope, T), F2x_add(*slope, a2));
120 gel(Q, 2) = F2x_add(F2xq_mul(*slope, F2x_add(x, gel(Q, 1)), T), F2x_add(y, gel(Q, 1)));
121 }
122 else
123 {
124 GEN a3 = gel(a,1), a4 = gel(a,2), a3i = gel(a,3);
125 *slope = F2xq_mul(F2x_add(a4, F2xq_sqr(x, T)), a3i, T);
126 Q = cgetg(3,t_VEC);
127 gel(Q, 1) = F2xq_sqr(*slope, T);
128 gel(Q, 2) = F2x_add(F2xq_mul(*slope, F2x_add(x, gel(Q, 1)), T), F2x_add(y, a3));
129 }
130 return Q;
131 }
132
133 GEN
F2xqE_dbl(GEN P,GEN a,GEN T)134 F2xqE_dbl(GEN P, GEN a, GEN T)
135 {
136 pari_sp av = avma;
137 GEN slope;
138 return gerepileupto(av, F2xqE_dbl_slope(P, a, T,&slope));
139 }
140
141 static GEN
F2xqE_add_slope(GEN P,GEN Q,GEN a,GEN T,GEN * slope)142 F2xqE_add_slope(GEN P, GEN Q, GEN a, GEN T, GEN *slope)
143 {
144 GEN Px, Py, Qx, Qy, R;
145 if (ell_is_inf(P)) return Q;
146 if (ell_is_inf(Q)) return P;
147 Px = gel(P,1); Py = gel(P,2);
148 Qx = gel(Q,1); Qy = gel(Q,2);
149 if (F2x_equal(Px, Qx))
150 {
151 if (F2x_equal(Py, Qy))
152 return F2xqE_dbl_slope(P, a, T, slope);
153 else
154 return ellinf();
155 }
156 *slope = F2xq_div(F2x_add(Py, Qy), F2x_add(Px, Qx), T);
157 R = cgetg(3,t_VEC);
158 if (typ(a)==t_VECSMALL)
159 {
160 GEN a2 = a;
161 gel(R, 1) = F2x_add(F2x_add(F2x_add(F2x_add(F2xq_sqr(*slope, T), *slope), Px), Qx), a2);
162 gel(R, 2) = F2x_add(F2xq_mul(*slope, F2x_add(Px, gel(R, 1)), T), F2x_add(Py, gel(R, 1)));
163 }
164 else
165 {
166 GEN a3 = gel(a,1);
167 gel(R, 1) = F2x_add(F2x_add(F2xq_sqr(*slope, T), Px), Qx);
168 gel(R, 2) = F2x_add(F2xq_mul(*slope, F2x_add(Px, gel(R, 1)), T), F2x_add(Py, a3));
169 }
170 return R;
171 }
172
173 GEN
F2xqE_add(GEN P,GEN Q,GEN a,GEN T)174 F2xqE_add(GEN P, GEN Q, GEN a, GEN T)
175 {
176 pari_sp av = avma;
177 GEN slope;
178 return gerepileupto(av, F2xqE_add_slope(P, Q, a, T, &slope));
179 }
180
181 static GEN
F2xqE_neg_i(GEN P,GEN a)182 F2xqE_neg_i(GEN P, GEN a)
183 {
184 GEN LHS;
185 if (ell_is_inf(P)) return P;
186 LHS = typ(a)==t_VECSMALL ? gel(P,1): gel(a,1);
187 return mkvec2(gel(P,1), F2x_add(LHS, gel(P,2)));
188 }
189
190 GEN
F2xqE_neg(GEN P,GEN a,GEN T)191 F2xqE_neg(GEN P, GEN a, GEN T)
192 {
193 GEN LHS;
194 (void) T;
195 if (ell_is_inf(P)) return ellinf();
196 LHS = typ(a)==t_VECSMALL ? gel(P,1): gel(a,1);
197 return mkvec2(gcopy(gel(P,1)), F2x_add(LHS, gel(P,2)));
198 }
199
200 GEN
F2xqE_sub(GEN P,GEN Q,GEN a2,GEN T)201 F2xqE_sub(GEN P, GEN Q, GEN a2, GEN T)
202 {
203 pari_sp av = avma;
204 GEN slope;
205 return gerepileupto(av, F2xqE_add_slope(P, F2xqE_neg_i(Q, a2), a2, T, &slope));
206 }
207
208 struct _F2xqE
209 {
210 GEN a2, a6;
211 GEN T;
212 };
213
214 static GEN
_F2xqE_dbl(void * E,GEN P)215 _F2xqE_dbl(void *E, GEN P)
216 {
217 struct _F2xqE *ell = (struct _F2xqE *) E;
218 return F2xqE_dbl(P, ell->a2, ell->T);
219 }
220
221 static GEN
_F2xqE_add(void * E,GEN P,GEN Q)222 _F2xqE_add(void *E, GEN P, GEN Q)
223 {
224 struct _F2xqE *ell=(struct _F2xqE *) E;
225 return F2xqE_add(P, Q, ell->a2, ell->T);
226 }
227
228 static GEN
_F2xqE_mul(void * E,GEN P,GEN n)229 _F2xqE_mul(void *E, GEN P, GEN n)
230 {
231 pari_sp av = avma;
232 struct _F2xqE *e=(struct _F2xqE *) E;
233 long s = signe(n);
234 if (!s || ell_is_inf(P)) return ellinf();
235 if (s<0) P = F2xqE_neg(P, e->a2, e->T);
236 if (is_pm1(n)) return s>0? gcopy(P): P;
237 return gerepilecopy(av, gen_pow_i(P, n, e, &_F2xqE_dbl, &_F2xqE_add));
238 }
239
240 GEN
F2xqE_mul(GEN P,GEN n,GEN a2,GEN T)241 F2xqE_mul(GEN P, GEN n, GEN a2, GEN T)
242 {
243 struct _F2xqE E;
244 E.a2 = a2; E.T = T;
245 return _F2xqE_mul(&E, P, n);
246 }
247
248 /* Finds a random nonsingular point on E */
249 GEN
random_F2xqE(GEN a,GEN a6,GEN T)250 random_F2xqE(GEN a, GEN a6, GEN T)
251 {
252 pari_sp ltop = avma;
253 GEN x, y, rhs, u;
254 do
255 {
256 set_avma(ltop);
257 x = random_F2x(F2x_degree(T),T[1]);
258 if (typ(a) == t_VECSMALL)
259 {
260 GEN a2 = a, x2;
261 if (!lgpol(x))
262 { set_avma(ltop); retmkvec2(pol0_Flx(T[1]), F2xq_sqrt(a6,T)); }
263 u = x; x2 = F2xq_sqr(x, T);
264 rhs = F2x_add(F2xq_mul(x2,F2x_add(x,a2),T),a6);
265 rhs = F2xq_div(rhs,x2,T);
266 }
267 else
268 {
269 GEN a3 = gel(a,1), a4 = gel(a,2), a3i = gel(a,3), u2i;
270 u = a3; u2i = F2xq_sqr(a3i,T);
271 rhs = F2x_add(F2xq_mul(x,F2x_add(F2xq_sqr(x,T),a4),T),a6);
272 rhs = F2xq_mul(rhs,u2i,T);
273 }
274 } while (F2xq_trace(rhs,T));
275 y = F2xq_mul(F2xq_Artin_Schreier(rhs, T), u, T);
276 return gerepilecopy(ltop, mkvec2(x, y));
277 }
278
279 static GEN
_F2xqE_rand(void * E)280 _F2xqE_rand(void *E)
281 {
282 struct _F2xqE *ell=(struct _F2xqE *) E;
283 return random_F2xqE(ell->a2, ell->a6, ell->T);
284 }
285
286 static const struct bb_group F2xqE_group={_F2xqE_add,_F2xqE_mul,_F2xqE_rand,hash_GEN,zvV_equal,ell_is_inf, NULL};
287
288 const struct bb_group *
get_F2xqE_group(void ** pt_E,GEN a2,GEN a6,GEN T)289 get_F2xqE_group(void ** pt_E, GEN a2, GEN a6, GEN T)
290 {
291 struct _F2xqE *e = (struct _F2xqE *) stack_malloc(sizeof(struct _F2xqE));
292 e->a2 = a2; e->a6 = a6; e->T = T;
293 *pt_E = (void *) e;
294 return &F2xqE_group;
295 }
296
297 GEN
F2xqE_order(GEN z,GEN o,GEN a2,GEN T)298 F2xqE_order(GEN z, GEN o, GEN a2, GEN T)
299 {
300 pari_sp av = avma;
301 struct _F2xqE e;
302 e.a2=a2; e.T=T;
303 return gerepileuptoint(av, gen_order(z, o, (void*)&e, &F2xqE_group));
304 }
305
306 GEN
F2xqE_log(GEN a,GEN b,GEN o,GEN a2,GEN T)307 F2xqE_log(GEN a, GEN b, GEN o, GEN a2, GEN T)
308 {
309 pari_sp av = avma;
310 struct _F2xqE e;
311 e.a2=a2; e.T=T;
312 return gerepileuptoint(av, gen_PH_log(a, b, o, (void*)&e, &F2xqE_group));
313 }
314
315 /***********************************************************************/
316 /** **/
317 /** Pairings **/
318 /** **/
319 /***********************************************************************/
320
321 /* Derived from APIP from and by Jerome Milan, 2012 */
322
323 static long
is_2_torsion(GEN Q,GEN a)324 is_2_torsion(GEN Q, GEN a)
325 {
326 return (typ(a)==t_VEC || lgpol(gel(Q, 1))) ? 0: 1;
327 }
328
329 static GEN
F2xqE_vert(GEN P,GEN Q,GEN a,GEN T)330 F2xqE_vert(GEN P, GEN Q, GEN a, GEN T)
331 {
332 long vT = T[1];
333 if (ell_is_inf(P))
334 return pol1_F2x(T[1]);
335 if (!F2x_equal(gel(Q, 1), gel(P, 1)))
336 return F2x_add(gel(Q, 1), gel(P, 1));
337 if (is_2_torsion(Q, a))
338 return F2xq_inv(gel(Q,2), T);
339 return pol1_F2x(vT);
340 }
341
342 static GEN
F2xqE_Miller_line(GEN R,GEN Q,GEN slope,GEN a,GEN T)343 F2xqE_Miller_line(GEN R, GEN Q, GEN slope, GEN a, GEN T)
344 {
345 long vT = T[1];
346 GEN x = gel(Q, 1), y = gel(Q, 2);
347 GEN tmp1 = F2x_add(x, gel(R, 1));
348 GEN tmp2 = F2x_add(F2xq_mul(tmp1, slope, T), gel(R, 2));
349 GEN s1, s2, ix;
350 if (!F2x_equal(y, tmp2))
351 return F2x_add(y, tmp2);
352 if (is_2_torsion(Q, a)) return pol1_F2x(vT);
353 if (typ(a)==t_VEC)
354 {
355 GEN a4 = gel(a,2), a3i = gel(a,3);
356 s1 = F2xq_mul(F2x_add(a4, F2xq_sqr(x, T)), a3i, T);
357 if (!F2x_equal(s1, slope))
358 return F2x_add(s1, slope);
359 s2 = F2xq_mul(F2x_add(x, F2xq_sqr(s1, T)), a3i, T);
360 if (lgpol(s2)) return s2;
361 return zv_copy(a3i);
362 } else
363 {
364 GEN a2 = a ;
365 ix = F2xq_inv(x, T);
366 s1 = F2x_add(x, F2xq_mul(y, ix, T));
367 if (!F2x_equal(s1, slope))
368 return F2x_add(s1, slope);
369 s2 =F2x_add(a2, F2x_add(F2xq_sqr(s1,T), s1));
370 if (!F2x_equal(s2, x))
371 return F2x_add(pol1_F2x(vT), F2xq_mul(s2, ix, T));
372 return ix;
373 }
374 }
375
376 /* Computes the equation of the line tangent to R and returns its
377 evaluation at the point Q. Also doubles the point R.
378 */
379
380 static GEN
F2xqE_tangent_update(GEN R,GEN Q,GEN a2,GEN T,GEN * pt_R)381 F2xqE_tangent_update(GEN R, GEN Q, GEN a2, GEN T, GEN *pt_R)
382 {
383 if (ell_is_inf(R))
384 {
385 *pt_R = ellinf();
386 return pol1_F2x(T[1]);
387 }
388 else if (is_2_torsion(R, a2))
389 {
390 *pt_R = ellinf();
391 return F2xqE_vert(R, Q, a2, T);
392 } else {
393 GEN slope;
394 *pt_R = F2xqE_dbl_slope(R, a2, T, &slope);
395 return F2xqE_Miller_line(R, Q, slope, a2, T);
396 }
397 }
398
399 /* Computes the equation of the line through R and P, and returns its
400 evaluation at the point Q. Also adds P to the point R.
401 */
402
403 static GEN
F2xqE_chord_update(GEN R,GEN P,GEN Q,GEN a2,GEN T,GEN * pt_R)404 F2xqE_chord_update(GEN R, GEN P, GEN Q, GEN a2, GEN T, GEN *pt_R)
405 {
406 if (ell_is_inf(R))
407 {
408 *pt_R = gcopy(P);
409 return F2xqE_vert(P, Q, a2, T);
410 }
411 else if (ell_is_inf(P))
412 {
413 *pt_R = gcopy(R);
414 return F2xqE_vert(R, Q, a2, T);
415 }
416 else if (F2x_equal(gel(P, 1), gel(R, 1)))
417 {
418 if (F2x_equal(gel(P, 2), gel(R, 2)))
419 return F2xqE_tangent_update(R, Q, a2, T, pt_R);
420 else
421 {
422 *pt_R = ellinf();
423 return F2xqE_vert(R, Q, a2, T);
424 }
425 } else {
426 GEN slope;
427 *pt_R = F2xqE_add_slope(P, R, a2, T, &slope);
428 return F2xqE_Miller_line(R, Q, slope, a2, T);
429 }
430 }
431
432 struct _F2xqE_miller { GEN T, a2, P; };
433
434 static GEN
F2xqE_Miller_dbl(void * E,GEN d)435 F2xqE_Miller_dbl(void* E, GEN d)
436 {
437 struct _F2xqE_miller *m = (struct _F2xqE_miller *)E;
438 GEN T = m->T, a2 = m->a2, P = m->P;
439 GEN v, line, point = gel(d,3);
440 GEN N = F2xq_sqr(gel(d,1), T);
441 GEN D = F2xq_sqr(gel(d,2), T);
442 line = F2xqE_tangent_update(point, P, a2, T, &point);
443 N = F2xq_mul(N, line, T);
444 v = F2xqE_vert(point, P, a2, T);
445 D = F2xq_mul(D, v, T); return mkvec3(N, D, point);
446 }
447 static GEN
F2xqE_Miller_add(void * E,GEN va,GEN vb)448 F2xqE_Miller_add(void* E, GEN va, GEN vb)
449 {
450 struct _F2xqE_miller *m = (struct _F2xqE_miller *)E;
451 GEN T = m->T, a2 = m->a2, P = m->P;
452 GEN v, line, point;
453 GEN na = gel(va,1), da = gel(va,2), pa = gel(va,3);
454 GEN nb = gel(vb,1), db = gel(vb,2), pb = gel(vb,3);
455 GEN N = F2xq_mul(na, nb, T);
456 GEN D = F2xq_mul(da, db, T);
457 line = F2xqE_chord_update(pa, pb, P, a2, T, &point);
458 N = F2xq_mul(N, line, T);
459 v = F2xqE_vert(point, P, a2, T);
460 D = F2xq_mul(D, v, T); return mkvec3(N, D, point);
461 }
462 /* Returns the Miller function f_{m, Q} evaluated at the point P using
463 * the standard Miller algorithm. */
464 static GEN
F2xqE_Miller(GEN Q,GEN P,GEN m,GEN a2,GEN T)465 F2xqE_Miller(GEN Q, GEN P, GEN m, GEN a2, GEN T)
466 {
467 pari_sp av = avma;
468 struct _F2xqE_miller d;
469 GEN v, N, D, g1;
470
471 d.a2 = a2; d.T = T; d.P = P;
472 g1 = pol1_F2x(T[1]);
473 v = gen_pow_i(mkvec3(g1,g1,Q), m, (void*)&d, F2xqE_Miller_dbl, F2xqE_Miller_add);
474 N = gel(v,1); D = gel(v,2);
475 return gerepileupto(av, F2xq_div(N, D, T));
476 }
477
478 GEN
F2xqE_weilpairing(GEN P,GEN Q,GEN m,GEN a2,GEN T)479 F2xqE_weilpairing(GEN P, GEN Q, GEN m, GEN a2, GEN T)
480 {
481 pari_sp av = avma;
482 GEN N, D;
483 if (ell_is_inf(P) || ell_is_inf(Q)
484 || (F2x_equal(gel(P,1),gel(Q,1)) && F2x_equal(gel(P,2),gel(Q,2))))
485 return pol1_F2x(T[1]);
486 N = F2xqE_Miller(P, Q, m, a2, T);
487 D = F2xqE_Miller(Q, P, m, a2, T);
488 return gerepileupto(av, F2xq_div(N, D, T));
489 }
490
491 GEN
F2xqE_tatepairing(GEN P,GEN Q,GEN m,GEN a2,GEN T)492 F2xqE_tatepairing(GEN P, GEN Q, GEN m, GEN a2, GEN T)
493 {
494 if (ell_is_inf(P) || ell_is_inf(Q))
495 return pol1_F2x(T[1]);
496 return F2xqE_Miller(P, Q, m, a2, T);
497 }
498
499 /***********************************************************************/
500 /** **/
501 /** Point counting **/
502 /** **/
503 /***********************************************************************/
504
505 static GEN
Z2x_rshift(GEN y,long x)506 Z2x_rshift(GEN y, long x)
507 {
508 GEN z;
509 long i, l;
510 if (!x) return pol0_Flx(y[1]);
511 z = cgetg_copy(y, &l); z[1] = y[1];
512 for(i=2; i<l; i++) z[i] = y[i]>>x;
513 return Flx_renormalize(z, l);
514 }
515
516 /* Solve the linear equation approximation in the Newton algorithm */
517
518 static GEN
gen_Z2x_Dixon(GEN F,GEN V,long N,void * E,GEN lin (void * E,GEN F,GEN d,long N),GEN invl (void * E,GEN d))519 gen_Z2x_Dixon(GEN F, GEN V, long N, void *E, GEN lin(void *E, GEN F, GEN d, long N), GEN invl(void *E, GEN d))
520 {
521 pari_sp av = avma;
522 long N2, M;
523 GEN VN2, V2, VM, bil;
524 ulong q = 1UL<<N;
525 if (N == 1) return invl(E, V);
526 V = Flx_red(V, q);
527 N2 = (N + 1)>>1; M = N - N2;
528 F = FlxT_red(F, q);
529 VN2 = gen_Z2x_Dixon(F, V, N2, E, lin, invl);
530 bil = lin(E, F, VN2, N);
531 V2 = Z2x_rshift(Flx_sub(V, bil, q), N2);
532 VM = gen_Z2x_Dixon(F, V2, M, E, lin, invl);
533 return gerepileupto(av, Flx_add(VN2, Flx_Fl_mul(VM, 1UL<<N2, q), q));
534 }
535
536 /* Solve F(X) = V mod 2^N
537 F(Xn) = V [mod 2^n]
538 Vm = (V-F(Xn))/(2^n)
539 F(Xm) = Vm
540 X = Xn + 2^n*Xm
541 */
542
543 static GEN
gen_Z2X_Dixon(GEN F,GEN V,long N,void * E,GEN lin (void * E,GEN F,GEN d,long N),GEN lins (void * E,GEN F,GEN d,long N),GEN invls (void * E,GEN d))544 gen_Z2X_Dixon(GEN F, GEN V, long N, void *E,
545 GEN lin(void *E, GEN F, GEN d, long N),
546 GEN lins(void *E, GEN F, GEN d, long N),
547 GEN invls(void *E, GEN d))
548 {
549 pari_sp av = avma;
550 long n, m;
551 GEN Xn, Xm, FXn, Vm;
552 if (N<BITS_IN_LONG)
553 {
554 ulong q = 1UL<<N;
555 return Flx_to_ZX(gen_Z2x_Dixon(ZXT_to_FlxT(F,q), ZX_to_Flx(V,q),N,E,lins,invls));
556 }
557 V = ZX_remi2n(V, N);
558 n = (N + 1)>>1; m = N - n;
559 F = ZXT_remi2n(F, N);
560 Xn = gen_Z2X_Dixon(F, V, n, E, lin, lins, invls);
561 FXn = lin(E, F, Xn, N);
562 Vm = ZX_shifti(ZX_sub(V, FXn), -n);
563 Xm = gen_Z2X_Dixon(F, Vm, m, E, lin, lins, invls);
564 return gerepileupto(av, ZX_remi2n(ZX_add(Xn, ZX_shifti(Xm, n)), N));
565 }
566
567 /* H -> H mod 2*/
568
_can_invls(void * E,GEN V)569 static GEN _can_invls(void *E, GEN V) {(void) E; return V; }
570
571 /* H -> H-(f0*H0-f1*H1) */
572
_can_lin(void * E,GEN F,GEN V,long N)573 static GEN _can_lin(void *E, GEN F, GEN V, long N)
574 {
575 pari_sp av=avma;
576 GEN d0, d1, z;
577 (void) E;
578 RgX_even_odd(V, &d0, &d1);
579 z = ZX_sub(V, ZX_sub(ZX_mul(gel(F,1), d0), ZX_mul(gel(F,2), d1)));
580 return gerepileupto(av, ZX_remi2n(z, N));
581 }
582
_can_lins(void * E,GEN F,GEN V,long N)583 static GEN _can_lins(void *E, GEN F, GEN V, long N)
584 {
585 GEN D=Flx_splitting(V, 2), z;
586 ulong q = 1UL<<N;
587 (void) E;
588 z = Flx_sub(Flx_mul(gel(F,1), gel(D,1), q), Flx_mul(gel(F,2), gel(D,2), q), q);
589 return Flx_sub(V, z, q);
590 }
591
592 /* P -> P-(P0^2-X*P1^2) */
593
594 static GEN
_can_iter(void * E,GEN f2,GEN q)595 _can_iter(void *E, GEN f2, GEN q)
596 {
597 GEN f0, f1, z;
598 (void) E;
599 RgX_even_odd(f2, &f0, &f1);
600 z = ZX_add(ZX_sub(f2, FpX_sqr(f0, q)), RgX_shift_shallow(FpX_sqr(f1, q), 1));
601 return mkvec3(z,f0,f1);
602 }
603
604 /* H -> H-(2*P0*H0-2*X*P1*H1) */
605
606 static GEN
_can_invd(void * E,GEN V,GEN v,GEN q,long M)607 _can_invd(void *E, GEN V, GEN v, GEN q, long M)
608 {
609 GEN F;
610 (void)E; (void)q;
611 F = mkvec2(ZX_shifti(gel(v,2),1), ZX_shifti(RgX_shift_shallow(gel(v,3),1),1));
612 return gen_Z2X_Dixon(F, V, M, NULL, _can_lin, _can_lins, _can_invls);
613 }
614
615 /* Lift P to Q such that Q(x^2)=Q(x)*Q(-x) mod 2^n
616 if Q = Q0(X^2)+X*Q1(X^2), solve Q(X^2) = Q0^2-X*Q1^2
617 */
618 GEN
F2x_Teichmuller(GEN P,long n)619 F2x_Teichmuller(GEN P, long n)
620 { return gen_ZpX_Newton(F2x_to_ZX(P),gen_2, n, NULL, _can_iter, _can_invd); }
621
622 static GEN
Z2XQ_frob(GEN x,GEN T,GEN q)623 Z2XQ_frob(GEN x, GEN T, GEN q)
624 {
625 return FpX_rem(RgX_inflate(x, 2), T, q);
626 }
627
628 static GEN
Z2xq_frob(GEN x,GEN T,ulong q)629 Z2xq_frob(GEN x, GEN T, ulong q)
630 {
631 return Flx_rem(Flx_inflate(x, 2), T, q);
632 }
633
634 struct _frob_lift
635 {
636 GEN T, sqx;
637 };
638
639 /* H -> S^-1(H) mod 2 */
640
_frob_invls(void * E,GEN V)641 static GEN _frob_invls(void *E, GEN V)
642 {
643 struct _frob_lift *F = (struct _frob_lift*) E;
644 GEN sqx = F->sqx;
645 return F2x_to_Flx(F2xq_sqrt_fast(Flx_to_F2x(V), gel(sqx,1), gel(sqx,2)));
646 }
647
648 /* H -> f1*S(H) + f2*H */
649
_frob_lin(void * E,GEN F,GEN x2,long N)650 static GEN _frob_lin(void *E, GEN F, GEN x2, long N)
651 {
652 GEN T = gel(F,3);
653 GEN q = int2n(N);
654 GEN y2 = Z2XQ_frob(x2, T, q);
655 GEN lin = ZX_add(ZX_mul(gel(F,1), y2), ZX_mul(gel(F,2), x2));
656 (void) E;
657 return FpX_rem(ZX_remi2n(lin, N), T, q);
658 }
659
_frob_lins(void * E,GEN F,GEN x2,long N)660 static GEN _frob_lins(void *E, GEN F, GEN x2, long N)
661 {
662 GEN T = gel(F,3);
663 ulong q = 1UL<<N;
664 GEN y2 = Z2xq_frob(x2, T, q);
665 GEN lin = Flx_add(Flx_mul(gel(F,1), y2,q), Flx_mul(gel(F,2), x2,q),q);
666 (void) E;
667 return Flx_rem(lin, T, q);
668 }
669
670 /* X -> P(X,S(X)) */
671
672 static GEN
_lift_iter(void * E,GEN x2,GEN q)673 _lift_iter(void *E, GEN x2, GEN q)
674 {
675 struct _frob_lift *F = (struct _frob_lift*) E;
676 long N = expi(q);
677 GEN TN = ZXT_remi2n(F->T, N);
678 GEN y2 = Z2XQ_frob(x2, TN, q);
679 GEN x2y2 = FpX_rem(ZX_remi2n(ZX_mul(x2, y2), N), TN, q);
680 GEN s = ZX_add(ZX_add(x2, ZX_shifti(y2, 1)), ZX_shifti(x2y2, 3));
681 GEN V = ZX_add(ZX_add(ZX_sqr(s), y2), ZX_shifti(x2y2, 2));
682 return mkvec4(FpX_rem(ZX_remi2n(V, N), TN, q),x2,y2,s);
683 }
684
685 /* H -> Dx*H+Dy*S(H) */
686
687 static GEN
_lift_invd(void * E,GEN V,GEN v,GEN qM,long M)688 _lift_invd(void *E, GEN V, GEN v, GEN qM, long M)
689 {
690 struct _frob_lift *F = (struct _frob_lift*) E;
691 GEN TM = ZXT_remi2n(F->T, M);
692 GEN x2 = gel(v,2), y2 = gel(v,3), s = gel(v,4), r;
693 GEN Dx = ZX_add(ZX_mul(ZX_Z_add(ZX_shifti(y2, 4), gen_2), s),
694 ZX_shifti(y2, 2));
695 GEN Dy = ZX_add(ZX_Z_add(ZX_mul(ZX_Z_add(ZX_shifti(x2, 4), utoi(4)), s),
696 gen_1), ZX_shifti(x2, 2));
697 Dx = FpX_rem(ZX_remi2n(Dx, M), TM, qM);
698 Dy = FpX_rem(ZX_remi2n(Dy, M), TM, qM);
699 r = mkvec3(Dy, Dx, TM);
700 return gen_Z2X_Dixon(r, V, M, E, _frob_lin, _frob_lins, _frob_invls);
701 }
702
703 /*
704 Let P(X,Y)=(X+2*Y+8*X*Y)^2+Y+4*X*Y
705 Solve P(x,S(x))=0 [mod 2^n,T]
706 assuming x = x0 [mod 2,T]
707
708 we set s = X+2*Y+8*X*Y, P = s^2+Y+4*X*Y
709 Dx = dP/dx = (16*s+4)*x+(4*s+1)
710 Dy = dP/dy = (16*y+2)*s+4*y
711 */
712
713 static GEN
solve_AGM_eqn(GEN x0,long n,GEN T,GEN sqx)714 solve_AGM_eqn(GEN x0, long n, GEN T, GEN sqx)
715 {
716 struct _frob_lift F;
717 F.T=T; F.sqx=sqx;
718 return gen_ZpX_Newton(x0, gen_2, n, &F, _lift_iter, _lift_invd);
719 }
720
721 static GEN
Z2XQ_invnorm_pcyc(GEN a,GEN T,long e)722 Z2XQ_invnorm_pcyc(GEN a, GEN T, long e)
723 {
724 GEN q = int2n(e);
725 GEN z = ZpXQ_norm_pcyc(a, T, q, gen_2);
726 return Fp_inv(z, q);
727 }
728
729 /* Assume a = 1 [4] */
730 static GEN
Z2XQ_invnorm(GEN a,GEN T,long e)731 Z2XQ_invnorm(GEN a, GEN T, long e)
732 {
733 pari_timer ti;
734 GEN pe = int2n(e), s;
735 if (degpol(a)==0)
736 return Fp_inv(Fp_powu(gel(a,2), get_FpX_degree(T), pe), pe);
737 if (DEBUGLEVEL>=3) timer_start(&ti);
738 s = ZpXQ_log(a, T, gen_2, e);
739 if (DEBUGLEVEL>=3) timer_printf(&ti,"Z2XQ_log");
740 s = Fp_neg(FpXQ_trace(s, T, pe), pe);
741 if (DEBUGLEVEL>=3) timer_printf(&ti,"FpXQ_trace");
742 s = modii(gel(Qp_exp(cvtop(s, gen_2, e-2)),4),pe);
743 if (DEBUGLEVEL>=3) timer_printf(&ti,"Qp_exp");
744 return s;
745 }
746
747 /* Assume a2==0, so 4|E(F_p): if t^4 = a6 then (t,t^2) is of order 4
748 8|E(F_p) <=> trace(a6)==0
749 */
750
751 static GEN
F2xq_elltrace_Harley(GEN a6,GEN T2)752 F2xq_elltrace_Harley(GEN a6, GEN T2)
753 {
754 pari_sp ltop = avma;
755 pari_timer ti;
756 GEN T, sqx;
757 GEN x, x2, t;
758 long n = F2x_degree(T2), N = ((n + 1)>>1) + 2;
759 long ispcyc;
760 if (n==1) return gen_m1;
761 if (n==2) return F2x_degree(a6) ? gen_1 : stoi(-3);
762 if (n==3) return F2x_degree(a6) ? (F2xq_trace(a6,T2) ? stoi(-3): gen_1)
763 : stoi(5);
764 timer_start(&ti);
765 sqx = mkvec2(F2xq_sqrt(polx_F2x(T2[1]),T2), T2);
766 if (DEBUGLEVEL>1) timer_printf(&ti,"Sqrtx");
767 ispcyc = zx_is_pcyc(F2x_to_Flx(T2));
768 T = ispcyc? F2x_to_ZX(T2): F2x_Teichmuller(T2, N-2);
769 if (DEBUGLEVEL>1) timer_printf(&ti,"Teich");
770 T = FpX_get_red(T, int2n(N));
771 if (DEBUGLEVEL>1) timer_printf(&ti,"Barrett");
772 x = solve_AGM_eqn(F2x_to_ZX(a6), N-2, T, sqx);
773 if (DEBUGLEVEL>1) timer_printf(&ti,"Lift");
774 x2 = ZX_Z_add_shallow(ZX_shifti(x,2), gen_1);
775 t = (ispcyc? Z2XQ_invnorm_pcyc: Z2XQ_invnorm)(x2, T, N);
776 if (DEBUGLEVEL>1) timer_printf(&ti,"Norm");
777 if (cmpii(sqri(t), int2n(n + 2)) > 0)
778 t = subii(t, int2n(N));
779 return gerepileuptoint(ltop, t);
780 }
781
782 static GEN
get_v(GEN u,GEN a,GEN T)783 get_v(GEN u, GEN a, GEN T)
784 {
785 GEN a4 = gel(a,2), a3i = gel(a,3);
786 GEN ui2 = F2xq_mul(u, a3i, T), ui4 = F2xq_sqr(ui2, T);
787 return F2xq_mul(a4, ui4, T);
788 }
789
790 static GEN
get_r(GEN w,GEN u,GEN T)791 get_r(GEN w, GEN u, GEN T)
792 {
793 return F2xq_sqr(F2xq_mul(F2xq_Artin_Schreier(w, T), u, T), T);
794 }
795
796 static GEN
F2xq_elltracej(GEN a,GEN a6,GEN T,GEN q,long n)797 F2xq_elltracej(GEN a, GEN a6, GEN T, GEN q, long n)
798 {
799 GEN a3 = gel(a,1), a4 = gel(a,2), a3i = gel(a,3);
800 GEN r, pi, rhs;
801 long t, s, m = n >> 1;
802 if (odd(n))
803 {
804 GEN u, v, w;
805 u = F2xq_pow(a3,diviuexact(subiu(shifti(q,1), 1), 3),T);
806 v = get_v(u, a, T);
807 if (F2xq_trace(v, T)==0) return gen_0;
808 w = F2xq_Artin_Schreier(F2x_1_add(v), T);
809 if (F2xq_trace(w, T)) w = F2x_1_add(w);
810 r = get_r(w, u, T);
811 pi = int2n(m+1);
812 s = (((m+3)&3L) <= 1) ? -1: 1;
813 }
814 else
815 {
816 if (F2x_degree(F2xq_pow(a3,diviuexact(subiu(q, 1), 3),T))==0)
817 {
818 GEN u, v, w;
819 u = F2xq_sqrtn(a3, utoi(3), T, NULL);
820 v = get_v(u, a, T);
821 if (F2xq_trace(v, T)==1) return gen_0;
822 w = F2xq_Artin_Schreier(v, T);
823 if (F2xq_trace(w, T)==1) return gen_0;
824 r = get_r(w, u, T);
825 pi = int2n(m+1);
826 s = odd(m+1) ? -1: 1;
827 }
828 else
829 {
830 long sv = T[1];
831 GEN P = mkpoln(5, pol1_F2x(sv), pol0_F2x(sv), pol0_F2x(sv), a3, a4);
832 r = F2xq_sqr(gel(F2xqX_roots(P,T), 1), T);
833 pi = int2n(m);
834 s = odd(m) ? -1: 1;
835 }
836 }
837 rhs = F2x_add(F2xq_mul(F2x_add(F2xq_sqr(r, T), a4), r, T), a6);
838 t = F2xq_trace(F2xq_mul(rhs, F2xq_sqr(a3i, T), T), T);
839 return (t==0)^(s==1)? pi: negi(pi);
840 }
841
842 GEN
F2xq_ellcard(GEN a,GEN a6,GEN T)843 F2xq_ellcard(GEN a, GEN a6, GEN T)
844 {
845 pari_sp av = avma;
846 long n = F2x_degree(T);
847 GEN c;
848 if (typ(a)==t_VECSMALL)
849 {
850 GEN t = F2xq_elltrace_Harley(a6, T);
851 c = addii(int2u(n), F2xq_trace(a,T) ? addui(1,t): subui(1,t));
852 } else if (n==1)
853 {
854 long a4i = lgpol(gel(a,2)), a6i = lgpol(a6);
855 return utoi(a4i? (a6i? 1: 5): 3);
856 }
857 else
858 {
859 GEN q = int2u(n);
860 c = subii(addiu(q,1), F2xq_elltracej(a, a6, T, q, n));
861 }
862 return gerepileuptoint(av, c);
863 }
864
865 /***********************************************************************/
866 /** **/
867 /** Group structure **/
868 /** **/
869 /***********************************************************************/
870
871 static GEN
_F2xqE_pairorder(void * E,GEN P,GEN Q,GEN m,GEN F)872 _F2xqE_pairorder(void *E, GEN P, GEN Q, GEN m, GEN F)
873 {
874 struct _F2xqE *e = (struct _F2xqE *) E;
875 return F2xq_order(F2xqE_weilpairing(P,Q,m,e->a2,e->T), F, e->T);
876 }
877
878 GEN
F2xq_ellgroup(GEN a2,GEN a6,GEN N,GEN T,GEN * pt_m)879 F2xq_ellgroup(GEN a2, GEN a6, GEN N, GEN T, GEN *pt_m)
880 {
881 struct _F2xqE e;
882 GEN q = int2u(F2x_degree(T));
883 e.a2=a2; e.a6=a6; e.T=T;
884 return gen_ellgroup(N, subiu(q,1), pt_m, (void*)&e, &F2xqE_group,
885 _F2xqE_pairorder);
886 }
887
888 GEN
F2xq_ellgens(GEN a2,GEN a6,GEN ch,GEN D,GEN m,GEN T)889 F2xq_ellgens(GEN a2, GEN a6, GEN ch, GEN D, GEN m, GEN T)
890 {
891 GEN P;
892 pari_sp av = avma;
893 struct _F2xqE e;
894 e.a2=a2; e.a6=a6; e.T=T;
895 switch(lg(D)-1)
896 {
897 case 0:
898 return cgetg(1,t_VEC);
899 case 1:
900 P = gen_gener(gel(D,1), (void*)&e, &F2xqE_group);
901 P = mkvec(F2xqE_changepoint(P, ch, T));
902 break;
903 default:
904 P = gen_ellgens(gel(D,1), gel(D,2), m, (void*)&e, &F2xqE_group,
905 _F2xqE_pairorder);
906 gel(P,1) = F2xqE_changepoint(gel(P,1), ch, T);
907 gel(P,2) = F2xqE_changepoint(gel(P,2), ch, T);
908 break;
909 }
910 return gerepilecopy(av, P);
911 }
912