Line data Source code
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 : #define DEBUGLEVEL DEBUGLEVEL_ellcard
19 :
20 : /* Not so fast arithmetic with points over elliptic curves over F_2^n */
21 :
22 : /***********************************************************************/
23 : /** **/
24 : /** F2xqE **/
25 : /** **/
26 : /***********************************************************************/
27 : /* These functions deal with points over elliptic curves over F_2^n defined
28 : * by an equation of the form:
29 : * y^2+x*y = x^3+a_2*x^2+a_6 if the curve is ordinary.
30 : * y^2+a_3*y = x^3+a_4*x+a_6 if the curve is supersingular.
31 : * Most of the time a6 is omitted since it can be recovered from any point
32 : * on the curve.
33 : * For supersingular curves, the parameter a2 is replaced by [a3,a4,a3^-1]. */
34 :
35 : GEN
36 10360 : RgE_to_F2xqE(GEN x, GEN T)
37 : {
38 10360 : if (ell_is_inf(x)) return x;
39 9877 : retmkvec2(Rg_to_F2xq(gel(x,1),T),Rg_to_F2xq(gel(x,2),T));
40 : }
41 :
42 : GEN
43 171794 : F2xqE_changepoint(GEN x, GEN ch, GEN T)
44 : {
45 171794 : pari_sp av = avma;
46 : GEN p1,z,u,r,s,t,v,v2,v3;
47 171794 : if (ell_is_inf(x)) return x;
48 171794 : u = gel(ch,1); r = gel(ch,2);
49 171794 : s = gel(ch,3); t = gel(ch,4);
50 171794 : v = F2xq_inv(u, T); v2 = F2xq_sqr(v, T); v3 = F2xq_mul(v,v2, T);
51 171794 : p1 = F2x_add(gel(x,1),r);
52 171794 : z = cgetg(3,t_VEC);
53 171794 : gel(z,1) = F2xq_mul(v2, p1, T);
54 171794 : gel(z,2) = F2xq_mul(v3, F2x_add(gel(x,2), F2x_add(F2xq_mul(s, p1, T),t)), T);
55 171794 : return gc_upto(av, z);
56 : }
57 :
58 : GEN
59 10360 : F2xqE_changepointinv(GEN x, GEN ch, GEN T)
60 : {
61 : GEN u, r, s, t, X, Y, u2, u3, u2X, z;
62 10360 : if (ell_is_inf(x)) return x;
63 9877 : X = gel(x,1); Y = gel(x,2);
64 9877 : u = gel(ch,1); r = gel(ch,2);
65 9877 : s = gel(ch,3); t = gel(ch,4);
66 9877 : u2 = F2xq_sqr(u, T); u3 = F2xq_mul(u,u2, T);
67 9877 : u2X = F2xq_mul(u2,X, T);
68 9877 : z = cgetg(3, t_VEC);
69 9877 : gel(z,1) = F2x_add(u2X,r);
70 9877 : gel(z,2) = F2x_add(F2xq_mul(u3,Y, T), F2x_add(F2xq_mul(s,u2X, T), t));
71 9877 : return z;
72 : }
73 :
74 : static GEN
75 3514 : nonzerotrace_F2xq(GEN T)
76 : {
77 3514 : pari_sp av = avma;
78 3514 : long n = F2x_degree(T), vs = T[1];
79 : GEN a;
80 3514 : if (odd(n))
81 1162 : return pol1_F2x(vs);
82 : do
83 : {
84 4668 : set_avma(av);
85 4668 : a = random_F2x(n, vs);
86 4668 : } while (F2xq_trace(a, T)==0);
87 2352 : return a;
88 : }
89 :
90 : void
91 3514 : F2xq_elltwist(GEN a, GEN a6, GEN T, GEN *pt_a, GEN *pt_a6)
92 : {
93 3514 : pari_sp av = avma;
94 3514 : GEN n = nonzerotrace_F2xq(T);
95 3514 : if (typ(a)==t_VECSMALL)
96 : {
97 3514 : *pt_a = gc_leaf(av, F2x_add(n, a));
98 3514 : *pt_a6 = vecsmall_copy(a6);
99 : } else
100 : {
101 0 : GEN a6t = F2x_add(a6, F2xq_mul(n, F2xq_sqr(gel(a,1), T), T));
102 0 : *pt_a6 = gc_leaf(av, a6t);
103 0 : *pt_a = vecsmall_copy(a);
104 : }
105 3514 : }
106 :
107 : static GEN
108 120660 : F2xqE_dbl_slope(GEN P, GEN a, GEN T, GEN *slope)
109 : {
110 : GEN x, y, Q;
111 120660 : if (ell_is_inf(P)) return ellinf();
112 119723 : x = gel(P,1); y = gel(P,2);
113 119723 : if (typ(a)==t_VECSMALL)
114 : {
115 99661 : GEN a2 = a;
116 99661 : if (!lgpol(gel(P,1))) { *slope = NULL; return ellinf(); }
117 95747 : *slope = F2x_add(x, F2xq_div(y, x, T));
118 95747 : Q = cgetg(3,t_VEC);
119 95747 : gel(Q, 1) = F2x_add(F2xq_sqr(*slope, T), F2x_add(*slope, a2));
120 95747 : 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 20062 : GEN a3 = gel(a,1), a4 = gel(a,2), a3i = gel(a,3);
125 20062 : *slope = F2xq_mul(F2x_add(a4, F2xq_sqr(x, T)), a3i, T);
126 20062 : Q = cgetg(3,t_VEC);
127 20062 : gel(Q, 1) = F2xq_sqr(*slope, T);
128 20062 : gel(Q, 2) = F2x_add(F2xq_mul(*slope, F2x_add(x, gel(Q, 1)), T), F2x_add(y, a3));
129 : }
130 115809 : return Q;
131 : }
132 :
133 : GEN
134 103686 : F2xqE_dbl(GEN P, GEN a, GEN T)
135 : {
136 103686 : pari_sp av = avma;
137 : GEN slope;
138 103686 : return gc_upto(av, F2xqE_dbl_slope(P, a, T,&slope));
139 : }
140 :
141 : static GEN
142 37656 : F2xqE_add_slope(GEN P, GEN Q, GEN a, GEN T, GEN *slope)
143 : {
144 : GEN Px, Py, Qx, Qy, R;
145 37656 : if (ell_is_inf(P)) { *slope = NULL; return Q; }
146 37160 : if (ell_is_inf(Q)) { *slope = NULL; return P; }
147 37160 : Px = gel(P,1); Py = gel(P,2);
148 37160 : Qx = gel(Q,1); Qy = gel(Q,2);
149 37160 : if (F2x_equal(Px, Qx))
150 : {
151 2108 : if (F2x_equal(Py, Qy))
152 77 : return F2xqE_dbl_slope(P, a, T, slope);
153 : else
154 2031 : { *slope = NULL; return ellinf(); }
155 : }
156 35052 : *slope = F2xq_div(F2x_add(Py, Qy), F2x_add(Px, Qx), T);
157 35052 : R = cgetg(3,t_VEC);
158 35052 : if (typ(a)==t_VECSMALL)
159 : {
160 33687 : GEN a2 = a;
161 33687 : gel(R, 1) = F2x_add(F2x_add(F2x_add(F2x_add(F2xq_sqr(*slope, T), *slope), Px), Qx), a2);
162 33687 : 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 1365 : GEN a3 = gel(a,1);
167 1365 : gel(R, 1) = F2x_add(F2x_add(F2xq_sqr(*slope, T), Px), Qx);
168 1365 : gel(R, 2) = F2x_add(F2xq_mul(*slope, F2x_add(Px, gel(R, 1)), T), F2x_add(Py, a3));
169 : }
170 35052 : return R;
171 : }
172 :
173 : GEN
174 37621 : F2xqE_add(GEN P, GEN Q, GEN a, GEN T)
175 : {
176 37621 : pari_sp av = avma;
177 : GEN slope;
178 37621 : return gc_upto(av, F2xqE_add_slope(P, Q, a, T, &slope));
179 : }
180 :
181 : static GEN
182 0 : F2xqE_neg_i(GEN P, GEN a)
183 : {
184 : GEN LHS;
185 0 : if (ell_is_inf(P)) return P;
186 0 : LHS = typ(a)==t_VECSMALL ? gel(P,1): gel(a,1);
187 0 : return mkvec2(gel(P,1), F2x_add(LHS, gel(P,2)));
188 : }
189 :
190 : GEN
191 84 : F2xqE_neg(GEN P, GEN a, GEN T)
192 : {
193 : GEN LHS;
194 : (void) T;
195 84 : if (ell_is_inf(P)) return ellinf();
196 84 : LHS = typ(a)==t_VECSMALL ? gel(P,1): gel(a,1);
197 84 : return mkvec2(gcopy(gel(P,1)), F2x_add(LHS, gel(P,2)));
198 : }
199 :
200 : GEN
201 0 : F2xqE_sub(GEN P, GEN Q, GEN a2, GEN T)
202 : {
203 0 : pari_sp av = avma;
204 : GEN slope;
205 0 : return gc_upto(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
215 103686 : _F2xqE_dbl(void *E, GEN P)
216 : {
217 103686 : struct _F2xqE *ell = (struct _F2xqE *) E;
218 103686 : return F2xqE_dbl(P, ell->a2, ell->T);
219 : }
220 :
221 : static GEN
222 37621 : _F2xqE_add(void *E, GEN P, GEN Q)
223 : {
224 37621 : struct _F2xqE *ell=(struct _F2xqE *) E;
225 37621 : return F2xqE_add(P, Q, ell->a2, ell->T);
226 : }
227 :
228 : static GEN
229 29520 : _F2xqE_mul(void *E, GEN P, GEN n)
230 : {
231 29520 : pari_sp av = avma;
232 29520 : struct _F2xqE *e=(struct _F2xqE *) E;
233 29520 : long s = signe(n);
234 29520 : if (!s || ell_is_inf(P)) return ellinf();
235 29499 : if (s<0) P = F2xqE_neg(P, e->a2, e->T);
236 29499 : if (is_pm1(n)) return s>0? gcopy(P): P;
237 28097 : return gc_GEN(av, gen_pow_i(P, n, e, &_F2xqE_dbl, &_F2xqE_add));
238 : }
239 :
240 : GEN
241 0 : F2xqE_mul(GEN P, GEN n, GEN a2, GEN T)
242 : {
243 : struct _F2xqE E;
244 0 : E.a2 = a2; E.T = T;
245 0 : return _F2xqE_mul(&E, P, n);
246 : }
247 :
248 : /* Finds a random nonsingular point on E */
249 : GEN
250 181796 : random_F2xqE(GEN a, GEN a6, GEN T)
251 : {
252 181796 : pari_sp ltop = avma;
253 : GEN x, y, rhs, u;
254 : do
255 : {
256 375121 : set_avma(ltop);
257 375121 : x = random_F2x(F2x_degree(T),T[1]);
258 375121 : if (typ(a) == t_VECSMALL)
259 : {
260 61977 : GEN a2 = a, x2;
261 61977 : if (!lgpol(x))
262 856 : { set_avma(ltop); retmkvec2(pol0_Flx(T[1]), F2xq_sqrt(a6,T)); }
263 61121 : u = x; x2 = F2xq_sqr(x, T);
264 61121 : rhs = F2x_add(F2xq_mul(x2,F2x_add(x,a2),T),a6);
265 61121 : rhs = F2xq_div(rhs,x2,T);
266 : }
267 : else
268 : {
269 313144 : GEN a3 = gel(a,1), a4 = gel(a,2), a3i = gel(a,3), u2i;
270 313144 : u = a3; u2i = F2xq_sqr(a3i,T);
271 313144 : rhs = F2x_add(F2xq_mul(x,F2x_add(F2xq_sqr(x,T),a4),T),a6);
272 313144 : rhs = F2xq_mul(rhs,u2i,T);
273 : }
274 374265 : } while (F2xq_trace(rhs,T));
275 180940 : y = F2xq_mul(F2xq_Artin_Schreier(rhs, T), u, T);
276 180940 : return gc_GEN(ltop, mkvec2(x, y));
277 : }
278 :
279 : static GEN
280 13824 : _F2xqE_rand(void *E)
281 : {
282 13824 : struct _F2xqE *ell=(struct _F2xqE *) E;
283 13824 : 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 *
289 0 : get_F2xqE_group(void ** pt_E, GEN a2, GEN a6, GEN T)
290 : {
291 0 : struct _F2xqE *e = (struct _F2xqE *) stack_malloc(sizeof(struct _F2xqE));
292 0 : e->a2 = a2; e->a6 = a6; e->T = T;
293 0 : *pt_E = (void *) e;
294 0 : return &F2xqE_group;
295 : }
296 :
297 : GEN
298 280 : F2xqE_order(GEN z, GEN o, GEN a2, GEN T)
299 : {
300 280 : pari_sp av = avma;
301 : struct _F2xqE e;
302 280 : e.a2=a2; e.T=T;
303 280 : return gc_INT(av, gen_order(z, o, (void*)&e, &F2xqE_group));
304 : }
305 :
306 : GEN
307 42 : F2xqE_log(GEN a, GEN b, GEN o, GEN a2, GEN T)
308 : {
309 42 : pari_sp av = avma;
310 : struct _F2xqE e;
311 42 : e.a2=a2; e.T=T;
312 42 : return gc_INT(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
324 19373 : is_2_torsion(GEN Q, GEN a)
325 : {
326 19373 : return (typ(a)==t_VEC || lgpol(gel(Q, 1))) ? 0: 1;
327 : }
328 :
329 : static GEN
330 35102 : F2xqE_vert(GEN P, GEN Q, GEN a, GEN T)
331 : {
332 35102 : long vT = T[1];
333 35102 : if (ell_is_inf(P))
334 9085 : return pol1_F2x(T[1]);
335 26017 : if (!F2x_equal(gel(Q, 1), gel(P, 1)))
336 24255 : return F2x_add(gel(Q, 1), gel(P, 1));
337 1762 : if (is_2_torsion(Q, a))
338 0 : return F2xq_inv(gel(Q,2), T);
339 1762 : return pol1_F2x(vT);
340 : }
341 :
342 : static GEN
343 16932 : F2xqE_Miller_line(GEN R, GEN Q, GEN slope, GEN a, GEN T)
344 : {
345 16932 : long vT = T[1];
346 16932 : GEN x = gel(Q, 1), y = gel(Q, 2);
347 16932 : GEN tmp1 = F2x_add(x, gel(R, 1));
348 16932 : GEN tmp2 = F2x_add(F2xq_mul(tmp1, slope, T), gel(R, 2));
349 : GEN s1, s2, ix;
350 16932 : if (!F2x_equal(y, tmp2))
351 16218 : return F2x_add(y, tmp2);
352 714 : if (is_2_torsion(Q, a)) return pol1_F2x(vT);
353 714 : if (typ(a)==t_VEC)
354 : {
355 700 : GEN a4 = gel(a,2), a3i = gel(a,3);
356 700 : s1 = F2xq_mul(F2x_add(a4, F2xq_sqr(x, T)), a3i, T);
357 700 : if (!F2x_equal(s1, slope))
358 350 : return F2x_add(s1, slope);
359 350 : s2 = F2xq_mul(F2x_add(x, F2xq_sqr(s1, T)), a3i, T);
360 350 : if (lgpol(s2)) return s2;
361 0 : return zv_copy(a3i);
362 : } else
363 : {
364 14 : GEN a2 = a ;
365 14 : ix = F2xq_inv(x, T);
366 14 : s1 = F2x_add(x, F2xq_mul(y, ix, T));
367 14 : if (!F2x_equal(s1, slope))
368 7 : return F2x_add(s1, slope);
369 7 : s2 =F2x_add(a2, F2x_add(F2xq_sqr(s1,T), s1));
370 7 : if (!F2x_equal(s2, x))
371 7 : return F2x_add(pol1_F2x(vT), F2xq_mul(s2, ix, T));
372 0 : 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
381 16897 : F2xqE_tangent_update(GEN R, GEN Q, GEN a2, GEN T, GEN *pt_R)
382 : {
383 16897 : if (ell_is_inf(R))
384 : {
385 0 : *pt_R = ellinf();
386 0 : return pol1_F2x(T[1]);
387 : }
388 16897 : else if (is_2_torsion(R, a2))
389 : {
390 0 : *pt_R = ellinf();
391 0 : return F2xqE_vert(R, Q, a2, T);
392 : } else {
393 : GEN slope;
394 16897 : *pt_R = F2xqE_dbl_slope(R, a2, T, &slope);
395 16897 : 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
404 9120 : F2xqE_chord_update(GEN R, GEN P, GEN Q, GEN a2, GEN T, GEN *pt_R)
405 : {
406 9120 : if (ell_is_inf(R))
407 : {
408 0 : *pt_R = gcopy(P);
409 0 : return F2xqE_vert(P, Q, a2, T);
410 : }
411 9120 : else if (ell_is_inf(P))
412 : {
413 0 : *pt_R = gcopy(R);
414 0 : return F2xqE_vert(R, Q, a2, T);
415 : }
416 9120 : else if (F2x_equal(gel(P, 1), gel(R, 1)))
417 : {
418 9085 : if (F2x_equal(gel(P, 2), gel(R, 2)))
419 0 : return F2xqE_tangent_update(R, Q, a2, T, pt_R);
420 : else
421 : {
422 9085 : *pt_R = ellinf();
423 9085 : return F2xqE_vert(R, Q, a2, T);
424 : }
425 : } else {
426 : GEN slope;
427 35 : *pt_R = F2xqE_add_slope(P, R, a2, T, &slope);
428 35 : return F2xqE_Miller_line(R, Q, slope, a2, T);
429 : }
430 : }
431 :
432 : struct _F2xqE_miller { GEN T, a2, P; };
433 :
434 : static GEN
435 16897 : F2xqE_Miller_dbl(void* E, GEN d)
436 : {
437 16897 : struct _F2xqE_miller *m = (struct _F2xqE_miller *)E;
438 16897 : GEN T = m->T, a2 = m->a2, P = m->P;
439 16897 : GEN v, line, point = gel(d,3);
440 16897 : GEN N = F2xq_sqr(gel(d,1), T);
441 16897 : GEN D = F2xq_sqr(gel(d,2), T);
442 16897 : line = F2xqE_tangent_update(point, P, a2, T, &point);
443 16897 : N = F2xq_mul(N, line, T);
444 16897 : v = F2xqE_vert(point, P, a2, T);
445 16897 : D = F2xq_mul(D, v, T); return mkvec3(N, D, point);
446 : }
447 : static GEN
448 9120 : F2xqE_Miller_add(void* E, GEN va, GEN vb)
449 : {
450 9120 : struct _F2xqE_miller *m = (struct _F2xqE_miller *)E;
451 9120 : GEN T = m->T, a2 = m->a2, P = m->P;
452 : GEN v, line, point;
453 9120 : GEN na = gel(va,1), da = gel(va,2), pa = gel(va,3);
454 9120 : GEN nb = gel(vb,1), db = gel(vb,2), pb = gel(vb,3);
455 9120 : GEN N = F2xq_mul(na, nb, T);
456 9120 : GEN D = F2xq_mul(da, db, T);
457 9120 : line = F2xqE_chord_update(pa, pb, P, a2, T, &point);
458 9120 : N = F2xq_mul(N, line, T);
459 9120 : v = F2xqE_vert(point, P, a2, T);
460 9120 : 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
465 9120 : F2xqE_Miller(GEN Q, GEN P, GEN m, GEN a2, GEN T)
466 : {
467 9120 : pari_sp av = avma;
468 : struct _F2xqE_miller d;
469 : GEN v, N, D, g1;
470 :
471 9120 : d.a2 = a2; d.T = T; d.P = P;
472 9120 : g1 = pol1_F2x(T[1]);
473 9120 : v = gen_pow_i(mkvec3(g1,g1,Q), m, (void*)&d, F2xqE_Miller_dbl, F2xqE_Miller_add);
474 9120 : N = gel(v,1); D = gel(v,2);
475 9120 : return gc_upto(av, F2xq_div(N, D, T));
476 : }
477 :
478 : GEN
479 5371 : F2xqE_weilpairing(GEN P, GEN Q, GEN m, GEN a2, GEN T)
480 : {
481 5371 : pari_sp av = avma;
482 : GEN N, D;
483 5371 : if (ell_is_inf(P) || ell_is_inf(Q)
484 4821 : || (F2x_equal(gel(P,1),gel(Q,1)) && F2x_equal(gel(P,2),gel(Q,2))))
485 825 : return pol1_F2x(T[1]);
486 4546 : N = F2xqE_Miller(P, Q, m, a2, T);
487 4546 : D = F2xqE_Miller(Q, P, m, a2, T);
488 4546 : return gc_upto(av, F2xq_div(N, D, T));
489 : }
490 :
491 : GEN
492 28 : F2xqE_tatepairing(GEN P, GEN Q, GEN m, GEN a2, GEN T)
493 : {
494 28 : if (ell_is_inf(P) || ell_is_inf(Q))
495 0 : return pol1_F2x(T[1]);
496 28 : return F2xqE_Miller(P, Q, m, a2, T);
497 : }
498 :
499 : /***********************************************************************/
500 : /** **/
501 : /** Point counting **/
502 : /** **/
503 : /***********************************************************************/
504 :
505 : static GEN
506 69228 : Z2x_rshift(GEN y, long x)
507 : {
508 : GEN z;
509 : long i, l;
510 69228 : if (!x) return pol0_Flx(y[1]);
511 69228 : z = cgetg_copy(y, &l); z[1] = y[1];
512 7719563 : for(i=2; i<l; i++) z[i] = y[i]>>x;
513 69228 : return Flx_renormalize(z, l);
514 : }
515 :
516 : /* Solve the linear equation approximation in the Newton algorithm */
517 :
518 : static GEN
519 171575 : 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 171575 : pari_sp av = avma;
522 : long N2, M;
523 : GEN VN2, V2, VM, bil;
524 171575 : ulong q = 1UL<<N;
525 171575 : if (N == 1) return invl(E, V);
526 69228 : V = Flx_red(V, q);
527 69228 : N2 = (N + 1)>>1; M = N - N2;
528 69228 : F = FlxT_red(F, q);
529 69228 : VN2 = gen_Z2x_Dixon(F, V, N2, E, lin, invl);
530 69228 : bil = lin(E, F, VN2, N);
531 69228 : V2 = Z2x_rshift(Flx_sub(V, bil, q), N2);
532 69228 : VM = gen_Z2x_Dixon(F, V2, M, E, lin, invl);
533 69228 : return gc_upto(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
544 33139 : 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 33139 : pari_sp av = avma;
550 : long n, m;
551 : GEN Xn, Xm, FXn, Vm;
552 33139 : if (N<BITS_IN_LONG)
553 : {
554 33119 : ulong q = 1UL<<N;
555 33119 : return Flx_to_ZX(gen_Z2x_Dixon(ZXT_to_FlxT(F,q), ZX_to_Flx(V,q),N,E,lins,invls));
556 : }
557 20 : V = ZX_remi2n(V, N);
558 20 : n = (N + 1)>>1; m = N - n;
559 20 : F = ZXT_remi2n(F, N);
560 20 : Xn = gen_Z2X_Dixon(F, V, n, E, lin, lins, invls);
561 20 : FXn = lin(E, F, Xn, N);
562 20 : Vm = ZX_shifti(ZX_sub(V, FXn), -n);
563 20 : Xm = gen_Z2X_Dixon(F, Vm, m, E, lin, lins, invls);
564 20 : return gc_upto(av, ZX_remi2n(ZX_add(Xn, ZX_shifti(Xm, n)), N));
565 : }
566 :
567 : /* H -> H mod 2*/
568 :
569 : static GEN
570 50918 : _can_invls(void *E, GEN V) {(void) E; return V; }
571 :
572 : /* H -> H-(f0*H0-f1*H1) */
573 :
574 : static GEN
575 10 : _can_lin(void *E, GEN F, GEN V, long N)
576 : {
577 10 : pari_sp av=avma;
578 : GEN d0, d1, z;
579 : (void) E;
580 10 : RgX_even_odd(V, &d0, &d1);
581 10 : z = ZX_sub(V, ZX_sub(ZX_mul(gel(F,1), d0), ZX_mul(gel(F,2), d1)));
582 10 : return gc_upto(av, ZX_remi2n(z, N));
583 : }
584 :
585 : static GEN
586 34481 : _can_lins(void *E, GEN F, GEN V, long N)
587 : {
588 34481 : GEN D=Flx_splitting(V, 2), z;
589 34481 : ulong q = 1UL<<N;
590 : (void) E;
591 34481 : z = Flx_sub(Flx_mul(gel(F,1), gel(D,1), q), Flx_mul(gel(F,2), gel(D,2), q), q);
592 34481 : return Flx_sub(V, z, q);
593 : }
594 :
595 : /* P -> P-(P0^2-X*P1^2) */
596 :
597 : static GEN
598 16427 : _can_iter(void *E, GEN f2, GEN q)
599 : {
600 : GEN f0, f1, z;
601 : (void) E;
602 16427 : RgX_even_odd(f2, &f0, &f1);
603 16427 : z = ZX_add(ZX_sub(f2, FpX_sqr(f0, q)), RgX_shift_shallow(FpX_sqr(f1, q), 1));
604 16427 : return mkvec3(z,f0,f1);
605 : }
606 :
607 : /* H -> H-(2*P0*H0-2*X*P1*H1) */
608 :
609 : static GEN
610 16427 : _can_invd(void *E, GEN V, GEN v, GEN q, long M)
611 : {
612 : GEN F;
613 : (void)E; (void)q;
614 16427 : F = mkvec2(ZX_shifti(gel(v,2),1), ZX_shifti(RgX_shift_shallow(gel(v,3),1),1));
615 16427 : return gen_Z2X_Dixon(F, V, M, NULL, _can_lin, _can_lins, _can_invls);
616 : }
617 :
618 : /* Lift P to Q such that Q(x^2)=Q(x)*Q(-x) mod 2^n
619 : if Q = Q0(X^2)+X*Q1(X^2), solve Q(X^2) = Q0^2-X*Q1^2
620 : */
621 : GEN
622 7191 : F2x_Teichmuller(GEN P, long n)
623 7191 : { return gen_ZpX_Newton(F2x_to_ZX(P),gen_2, n, NULL, _can_iter, _can_invd); }
624 :
625 : static GEN
626 16682 : Z2XQ_frob(GEN x, GEN T, GEN q)
627 : {
628 16682 : return FpX_rem(RgX_inflate(x, 2), T, q);
629 : }
630 :
631 : static GEN
632 34747 : Z2xq_frob(GEN x, GEN T, ulong q)
633 : {
634 34747 : return Flx_rem(Flx_inflate(x, 2), T, q);
635 : }
636 :
637 : struct _frob_lift
638 : {
639 : GEN T, sqx;
640 : };
641 :
642 : /* H -> S^-1(H) mod 2 */
643 :
644 : static GEN
645 51429 : _frob_invls(void *E, GEN V)
646 : {
647 51429 : struct _frob_lift *F = (struct _frob_lift*) E;
648 51429 : GEN sqx = F->sqx;
649 51429 : return F2x_to_Flx(F2xq_sqrt_fast(Flx_to_F2x(V), gel(sqx,1), gel(sqx,2)));
650 : }
651 :
652 : /* H -> f1*S(H) + f2*H */
653 :
654 : static GEN
655 10 : _frob_lin(void *E, GEN F, GEN x2, long N)
656 : {
657 10 : GEN T = gel(F,3);
658 10 : GEN q = int2n(N);
659 10 : GEN y2 = Z2XQ_frob(x2, T, q);
660 10 : GEN lin = ZX_add(ZX_mul(gel(F,1), y2), ZX_mul(gel(F,2), x2));
661 : (void) E;
662 10 : return FpX_rem(ZX_remi2n(lin, N), T, q);
663 : }
664 :
665 : static GEN
666 34747 : _frob_lins(void *E, GEN F, GEN x2, long N)
667 : {
668 34747 : GEN T = gel(F,3);
669 34747 : ulong q = 1UL<<N;
670 34747 : GEN y2 = Z2xq_frob(x2, T, q);
671 34747 : GEN lin = Flx_add(Flx_mul(gel(F,1), y2,q), Flx_mul(gel(F,2), x2,q),q);
672 : (void) E;
673 34747 : return Flx_rem(lin, T, q);
674 : }
675 :
676 : /* X -> P(X,S(X)) */
677 :
678 : static GEN
679 16672 : _lift_iter(void *E, GEN x2, GEN q)
680 : {
681 16672 : struct _frob_lift *F = (struct _frob_lift*) E;
682 16672 : long N = expi(q);
683 16672 : GEN TN = ZXT_remi2n(F->T, N);
684 16672 : GEN y2 = Z2XQ_frob(x2, TN, q);
685 16672 : GEN x2y2 = FpX_rem(ZX_remi2n(ZX_mul(x2, y2), N), TN, q);
686 16672 : GEN s = ZX_add(ZX_add(x2, ZX_shifti(y2, 1)), ZX_shifti(x2y2, 3));
687 16672 : GEN V = ZX_add(ZX_add(ZX_sqr(s), y2), ZX_shifti(x2y2, 2));
688 16672 : return mkvec4(FpX_rem(ZX_remi2n(V, N), TN, q),x2,y2,s);
689 : }
690 :
691 : /* H -> Dx*H+Dy*S(H) */
692 : static GEN
693 16672 : _lift_invd(void *E, GEN V, GEN v, GEN qM, long M)
694 : {
695 16672 : struct _frob_lift *F = (struct _frob_lift*) E;
696 16672 : GEN TM = ZXT_remi2n(F->T, M);
697 16672 : GEN x2 = gel(v,2), y2 = gel(v,3), s = gel(v,4), r;
698 16672 : GEN Dx = ZX_add(ZX_mul(ZX_Z_add(ZX_shifti(y2, 4), gen_2), s),
699 : ZX_shifti(y2, 2));
700 16672 : GEN Dy = ZX_add(ZX_Z_add(ZX_mul(ZX_Z_add(ZX_shifti(x2, 4), utoi(4)), s),
701 : gen_1), ZX_shifti(x2, 2));
702 16672 : Dx = FpX_rem(ZX_remi2n(Dx, M), TM, qM);
703 16672 : Dy = FpX_rem(ZX_remi2n(Dy, M), TM, qM);
704 16672 : r = mkvec3(Dy, Dx, TM);
705 16672 : return gen_Z2X_Dixon(r, V, M, E, _frob_lin, _frob_lins, _frob_invls);
706 : }
707 :
708 : /* Let P(X,Y)=(X+2*Y+8*X*Y)^2+Y+4*X*Y
709 : * Solve P(x,S(x))=0 [mod 2^n,T]
710 : * assuming x = x0 [mod 2,T]
711 : *
712 : * we set s = X+2*Y+8*X*Y, P = s^2+Y+4*X*Y
713 : * Dx = dP/dx = (16*s+4)*x+(4*s+1)
714 : * Dy = dP/dy = (16*y+2)*s+4*y */
715 : static GEN
716 7422 : solve_AGM_eqn(GEN x0, long n, GEN T, GEN sqx)
717 : {
718 : struct _frob_lift F;
719 7422 : F.T=T; F.sqx=sqx;
720 7422 : return gen_ZpX_Newton(x0, gen_2, n, &F, _lift_iter, _lift_invd);
721 : }
722 :
723 : static GEN
724 238 : Z2XQ_invnorm_pcyc(GEN a, GEN T, long e)
725 : {
726 238 : GEN q = int2n(e);
727 238 : GEN z = ZpXQ_norm_pcyc(a, T, q, gen_2);
728 238 : return Zp_inv(z, gen_2, e);
729 : }
730 :
731 : /* Assume a = 1 [4] */
732 : static GEN
733 7184 : Z2XQ_invnorm(GEN a, GEN T, long e)
734 : {
735 : pari_timer ti;
736 7184 : GEN pe = int2n(e), s;
737 7184 : if (degpol(a)==0)
738 56 : return Zp_inv(Fp_powu(gel(a,2), get_FpX_degree(T), pe), gen_2, e);
739 7128 : if (DEBUGLEVEL>=3) timer_start(&ti);
740 7128 : s = ZpXQ_log(a, T, gen_2, e);
741 7128 : if (DEBUGLEVEL>=3) timer_printf(&ti,"Z2XQ_log");
742 7128 : s = Fp_neg(FpXQ_trace(s, T, pe), pe);
743 7128 : if (DEBUGLEVEL>=3) timer_printf(&ti,"FpXQ_trace");
744 7128 : s = Qp_exp(cvtop(s, gen_2, e-2));
745 7128 : s = modii(padic_u(s), pe);
746 7128 : if (DEBUGLEVEL>=3) timer_printf(&ti,"Qp_exp");
747 7128 : return s;
748 : }
749 :
750 : static int
751 7422 : F2x_is_pcyc(GEN T)
752 : {
753 : long i, n;
754 : ulong k;
755 7422 : if (!uisprime(F2x_degree(T) + 1)) return 0;
756 1127 : n = lg(T)-1;
757 1136 : for (i = 2; i < n; i++)
758 16 : if (uel(T,i) != ~0UL) return 0;
759 1120 : k = uel(T,i); return !(k & (k+1)); /* test if k+1 is a power of 2 */
760 : }
761 :
762 : /* Assume a2==0, so 4|E(F_p): if t^4 = a6 then (t,t^2) is of order 4
763 : 8|E(F_p) <=> trace(a6)==0 */
764 : static GEN
765 7765 : F2xq_elltrace_Harley(GEN a6, GEN T2)
766 : {
767 7765 : pari_sp ltop = avma;
768 : pari_timer ti;
769 : GEN T, sqx;
770 : GEN x, x2, t;
771 7765 : long n = F2x_degree(T2), N = ((n + 1)>>1) + 2;
772 : long ispcyc;
773 7765 : if (n==1) return gen_m1;
774 7576 : if (n==2) return F2x_degree(a6) ? gen_1 : stoi(-3);
775 7625 : if (n==3) return F2x_degree(a6) ? (F2xq_trace(a6,T2) ? stoi(-3): gen_1)
776 203 : : stoi(5);
777 7422 : timer_start(&ti);
778 7422 : sqx = mkvec2(F2xq_sqrt(polx_F2x(T2[1]),T2), T2);
779 7422 : if (DEBUGLEVEL>1) timer_printf(&ti,"Sqrtx");
780 7422 : ispcyc = F2x_is_pcyc(T2);
781 7422 : T = ispcyc? F2x_to_ZX(T2): F2x_Teichmuller(T2, N-2);
782 7422 : if (DEBUGLEVEL>1) timer_printf(&ti,"Teich");
783 7422 : T = FpX_get_red(T, int2n(N));
784 7422 : if (DEBUGLEVEL>1) timer_printf(&ti,"Barrett");
785 7422 : x = solve_AGM_eqn(F2x_to_ZX(a6), N-2, T, sqx);
786 7422 : if (DEBUGLEVEL>1) timer_printf(&ti,"Lift");
787 7422 : x2 = ZX_Z_add_shallow(ZX_shifti(x,2), gen_1);
788 7422 : t = (ispcyc? Z2XQ_invnorm_pcyc: Z2XQ_invnorm)(x2, T, N);
789 7422 : if (DEBUGLEVEL>1) timer_printf(&ti,"Norm");
790 7422 : if (cmpii(sqri(t), int2n(n + 2)) > 0)
791 3501 : t = subii(t, int2n(N));
792 7422 : return gc_INT(ltop, t);
793 : }
794 :
795 : static GEN
796 23855 : get_v(GEN u, GEN a, GEN T)
797 : {
798 23855 : GEN a4 = gel(a,2), a3i = gel(a,3);
799 23855 : GEN ui2 = F2xq_mul(u, a3i, T), ui4 = F2xq_sqr(ui2, T);
800 23855 : return F2xq_mul(a4, ui4, T);
801 : }
802 :
803 : static GEN
804 10465 : get_r(GEN w, GEN u, GEN T)
805 : {
806 10465 : return F2xq_sqr(F2xq_mul(F2xq_Artin_Schreier(w, T), u, T), T);
807 : }
808 :
809 : static GEN
810 37734 : F2xq_elltracej(GEN a, GEN a6, GEN T, GEN q, long n)
811 : {
812 37734 : GEN a3 = gel(a,1), a4 = gel(a,2), a3i = gel(a,3);
813 : GEN r, pi, rhs;
814 37734 : long t, s, m = n >> 1;
815 37734 : if (odd(n))
816 : {
817 : GEN u, v, w;
818 16877 : u = F2xq_pow(a3,diviuexact(subiu(shifti(q,1), 1), 3),T);
819 16877 : v = get_v(u, a, T);
820 16877 : if (F2xq_trace(v, T)==0) return gen_0;
821 8692 : w = F2xq_Artin_Schreier(F2x_1_add(v), T);
822 8692 : if (F2xq_trace(w, T)) w = F2x_1_add(w);
823 8692 : r = get_r(w, u, T);
824 8692 : pi = int2n(m+1);
825 8692 : s = (((m+3)&3L) <= 1) ? -1: 1;
826 : }
827 : else
828 : {
829 20857 : if (F2x_degree(F2xq_pow(a3,diviuexact(subiu(q, 1), 3),T))==0)
830 : {
831 : GEN u, v, w;
832 6978 : u = F2xq_sqrtn(a3, utoi(3), T, NULL);
833 6978 : v = get_v(u, a, T);
834 6978 : if (F2xq_trace(v, T)==1) return gen_0;
835 3556 : w = F2xq_Artin_Schreier(v, T);
836 3556 : if (F2xq_trace(w, T)==1) return gen_0;
837 1773 : r = get_r(w, u, T);
838 1773 : pi = int2n(m+1);
839 1773 : s = odd(m+1) ? -1: 1;
840 : }
841 : else
842 : {
843 13879 : long sv = T[1];
844 13879 : GEN P = mkpoln(5, pol1_F2x(sv), pol0_F2x(sv), pol0_F2x(sv), a3, a4);
845 13879 : r = F2xq_sqr(gel(F2xqX_roots(P,T), 1), T);
846 13879 : pi = int2n(m);
847 13879 : s = odd(m) ? -1: 1;
848 : }
849 : }
850 24344 : rhs = F2x_add(F2xq_mul(F2x_add(F2xq_sqr(r, T), a4), r, T), a6);
851 24344 : t = F2xq_trace(F2xq_mul(rhs, F2xq_sqr(a3i, T), T), T);
852 24344 : return (t==0)^(s==1)? pi: negi(pi);
853 : }
854 :
855 : GEN
856 45699 : F2xq_ellcard(GEN a, GEN a6, GEN T)
857 : {
858 45699 : pari_sp av = avma;
859 45699 : long n = F2x_degree(T);
860 : GEN c;
861 45699 : if (typ(a)==t_VECSMALL)
862 : {
863 7765 : GEN t = F2xq_elltrace_Harley(a6, T);
864 7765 : c = addii(int2u(n), F2xq_trace(a,T) ? addui(1,t): subui(1,t));
865 37934 : } else if (n==1)
866 : {
867 200 : long a4i = lgpol(gel(a,2)), a6i = lgpol(a6);
868 200 : return utoi(a4i? (a6i? 1: 5): 3);
869 : }
870 : else
871 : {
872 37734 : GEN q = int2u(n);
873 37734 : c = subii(addiu(q,1), F2xq_elltracej(a, a6, T, q, n));
874 : }
875 45499 : return gc_INT(av, c);
876 : }
877 :
878 : /***********************************************************************/
879 : /** **/
880 : /** Group structure **/
881 : /** **/
882 : /***********************************************************************/
883 :
884 : static GEN
885 401 : _F2xqE_pairorder(void *E, GEN P, GEN Q, GEN m, GEN F)
886 : {
887 401 : struct _F2xqE *e = (struct _F2xqE *) E;
888 401 : return F2xq_order(F2xqE_weilpairing(P,Q,m,e->a2,e->T), F, e->T);
889 : }
890 :
891 : GEN
892 3808 : F2xq_ellgroup(GEN a2, GEN a6, GEN N, GEN T, GEN *pt_m)
893 : {
894 : struct _F2xqE e;
895 3808 : GEN q = int2u(F2x_degree(T));
896 3808 : e.a2=a2; e.a6=a6; e.T=T;
897 3808 : return gen_ellgroup(N, subiu(q,1), pt_m, (void*)&e, &F2xqE_group,
898 : _F2xqE_pairorder);
899 : }
900 :
901 : GEN
902 3738 : F2xq_ellgens(GEN a2, GEN a6, GEN ch, GEN D, GEN m, GEN T)
903 : {
904 : GEN P;
905 3738 : pari_sp av = avma;
906 : struct _F2xqE e;
907 3738 : e.a2=a2; e.a6=a6; e.T=T;
908 3738 : switch(lg(D)-1)
909 : {
910 28 : case 0:
911 28 : return cgetg(1,t_VEC);
912 3598 : case 1:
913 3598 : P = gen_gener(gel(D,1), (void*)&e, &F2xqE_group);
914 3598 : P = mkvec(F2xqE_changepoint(P, ch, T));
915 3598 : break;
916 112 : default:
917 112 : P = gen_ellgens(gel(D,1), gel(D,2), m, (void*)&e, &F2xqE_group,
918 : _F2xqE_pairorder);
919 112 : gel(P,1) = F2xqE_changepoint(gel(P,1), ch, T);
920 112 : gel(P,2) = F2xqE_changepoint(gel(P,2), ch, T);
921 112 : break;
922 : }
923 3710 : return gc_GEN(av, P);
924 : }
925 :
926 : GEN
927 46168 : ell_to_a4a6_F2xq(GEN E, GEN T)
928 : {
929 46168 : long v = T[1];
930 46168 : GEN a1 = Rg_to_F2xq(ell_get_a1(E),T);
931 46168 : GEN a2 = Rg_to_F2xq(ell_get_a2(E),T);
932 46168 : GEN a3 = Rg_to_F2xq(ell_get_a3(E),T);
933 46168 : GEN a4 = Rg_to_F2xq(ell_get_a4(E),T);
934 46168 : GEN a6 = Rg_to_F2xq(ell_get_a6(E),T);
935 46168 : if (lgpol(a1))
936 : {
937 7961 : GEN a1i = F2xq_inv(a1,T);
938 7961 : GEN a1i2 = F2xq_sqr(a1i,T);
939 7961 : GEN a1i3 = F2xq_mul(a1i,a1i2,T);
940 7961 : GEN a1i6 = F2xq_sqr(a1i3,T);
941 7961 : GEN d = F2xq_mul(a3,a1i,T);
942 7961 : GEN dd = F2xq_mul(d,a1i2,T);
943 7961 : GEN e = F2xq_mul(F2x_add(a4,F2xq_sqr(d,T)),a1i,T);
944 7961 : GEN ee = F2xq_mul(e,a1i3,T);
945 7961 : GEN da2 = F2x_add(a2,d);
946 7961 : GEN d2 = F2xq_mul(da2,a1i2,T);
947 7961 : GEN d6 = F2xq_mul(F2x_add(F2x_add(F2xq_mul(F2x_add(F2xq_mul(da2,d,T),a4),d,T),a6),F2xq_sqr(e,T)),a1i6,T);
948 7961 : retmkvec3(d2, d6, mkvec4(a1i,dd,pol0_F2x(v),ee));
949 : }
950 : else
951 : { /* allow a1 = a3 = 0: singular curve */
952 38207 : GEN d4 = F2x_add(F2xq_sqr(a2,T),a4);
953 38207 : GEN d6 = F2x_add(F2xq_mul(a2,a4,T),a6);
954 38207 : GEN a3i = lgpol(a3)? F2xq_inv(a3,T): a3;
955 38207 : retmkvec3(mkvec3(a3,d4,a3i), d6, mkvec4(pol1_F2x(v),a2,pol0_F2x(T[1]),pol0_F2x(T[1])));
956 : }
957 : }
|
