Line data Source code
1 : /* Copyright (C) 2014 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 : /********************************************************************/
16 : /** **/
17 : /** HYPERELLIPTIC CURVES **/
18 : /** **/
19 : /********************************************************************/
20 : #include "pari.h"
21 : #include "paripriv.h"
22 :
23 : #define DEBUGLEVEL DEBUGLEVEL_hyperell
24 :
25 : /* Implementation of Kedlaya Algorithm for counting point on hyperelliptic
26 : curves by Bill Allombert based on a GP script by Bernadette Perrin-Riou.
27 :
28 : References:
29 : Pierrick Gaudry and Nicolas G\"urel
30 : Counting Points in Medium Characteristic Using Kedlaya's Algorithm
31 : Experiment. Math. Volume 12, Number 4 (2003), 395-402.
32 : http://projecteuclid.org/euclid.em/1087568016
33 :
34 : Harrison, M. An extension of Kedlaya's algorithm for hyperelliptic
35 : curves. Journal of Symbolic Computation, 47 (1) (2012), 89-101.
36 : http://arxiv.org/pdf/1006.4206v3.pdf
37 : */
38 :
39 : /* We use the basis of differentials (x^i*dx/y^k) (i=1 to 2*g-1),
40 : with k either 1 or 3, depending on p and d, see Harrison paper */
41 :
42 : static long
43 1764 : get_basis(long p, long d)
44 : {
45 1764 : if (odd(d))
46 868 : return p < d-1 ? 3 : 1;
47 : else
48 896 : return 2*p <= d-2 ? 3 : 1;
49 : }
50 :
51 : static GEN
52 20265 : FpXXQ_red(GEN S, GEN T, GEN p)
53 : {
54 20265 : pari_sp av = avma;
55 20265 : long i, dS = degpol(S);
56 : GEN A, C;
57 20265 : if (signe(S)==0) return pol_0(varn(T));
58 20265 : A = cgetg(dS+3, t_POL);
59 20265 : C = pol_0(varn(T));
60 1520393 : for(i=dS; i>0; i--)
61 : {
62 1500128 : GEN Si = FpX_add(C, gel(S,i+2), p);
63 1500128 : GEN R, Q = FpX_divrem(Si, T, p, &R);
64 1500128 : gel(A,i+2) = R;
65 1500128 : C = Q;
66 : }
67 20265 : gel(A,2) = FpX_add(C, gel(S,2), p);
68 20265 : A[1] = S[1];
69 20265 : return gerepilecopy(av, FpXX_renormalize(A,dS+3));
70 : }
71 :
72 : static GEN
73 3402 : FpXXQ_sqr(GEN x, GEN T, GEN p)
74 : {
75 3402 : pari_sp av = avma;
76 3402 : long n = degpol(T);
77 3402 : GEN z = FpX_red(ZXX_sqr_Kronecker(x, n), p);
78 3402 : z = Kronecker_to_ZXX(z, n, varn(T));
79 3402 : return gerepileupto(av, FpXXQ_red(z, T, p));
80 : }
81 :
82 : static GEN
83 16863 : FpXXQ_mul(GEN x, GEN y, GEN T, GEN p)
84 : {
85 16863 : pari_sp av = avma;
86 16863 : long n = degpol(T);
87 16863 : GEN z = FpX_red(ZXX_mul_Kronecker(x, y, n), p);
88 16863 : z = Kronecker_to_ZXX(z, n, varn(T));
89 16863 : return gerepileupto(av, FpXXQ_red(z, T, p));
90 : }
91 :
92 : static GEN
93 1309 : ZpXXQ_invsqrt(GEN S, GEN T, ulong p, long e)
94 : {
95 1309 : pari_sp av = avma, av2;
96 : ulong mask;
97 1309 : long v = varn(S), n=1;
98 1309 : GEN a = pol_1(v);
99 1309 : if (e <= 1) return gerepilecopy(av, a);
100 1309 : mask = quadratic_prec_mask(e);
101 1309 : av2 = avma;
102 4676 : for (;mask>1;)
103 : {
104 : GEN q, q2, q22, f, fq, afq;
105 3367 : long n2 = n;
106 3367 : n<<=1; if (mask & 1) n--;
107 3367 : mask >>= 1;
108 3367 : q = powuu(p,n); q2 = powuu(p,n2);
109 3367 : f = RgX_sub(FpXXQ_mul(FpXX_red(S, q), FpXXQ_sqr(a, T, q), T, q), pol_1(v));
110 3367 : fq = ZXX_Z_divexact(f, q2);
111 3367 : q22 = shifti(addiu(q2,1),-1);
112 3367 : afq = FpXX_Fp_mul(FpXXQ_mul(a, fq, T, q2), q22, q2);
113 3367 : a = RgX_sub(a, ZXX_Z_mul(afq, q2));
114 3367 : if (gc_needed(av2,1))
115 : {
116 0 : if(DEBUGMEM>1) pari_warn(warnmem,"ZpXXQ_invsqrt, e = %ld", n);
117 0 : a = gerepileupto(av2, a);
118 : }
119 : }
120 1309 : return gerepileupto(av, a);
121 : }
122 :
123 : static GEN
124 1029007 : to_ZX(GEN a, long v) { return typ(a)==t_INT? scalarpol(a,v): a; }
125 :
126 : static void
127 14 : is_sing(GEN H, ulong p)
128 : {
129 14 : pari_err_DOMAIN("hyperellpadicfrobenius","H","is singular at",utoi(p),H);
130 0 : }
131 :
132 : static void
133 1309 : get_UV(GEN *U, GEN *V, GEN T, ulong p, long e)
134 : {
135 1309 : GEN q = powuu(p,e), d;
136 1309 : GEN dT = FpX_deriv(T, q);
137 1309 : GEN R = polresultantext(T, dT);
138 1309 : long v = varn(T);
139 1309 : if (dvdiu(gel(R,3),p)) is_sing(T, p);
140 1309 : d = Zp_inv(gel(R,3), utoi(p), e);
141 1309 : *U = FpX_Fp_mul(FpX_red(to_ZX(gel(R,1),v),q),d,q);
142 1309 : *V = FpX_Fp_mul(FpX_red(to_ZX(gel(R,2),v),q),d,q);
143 1309 : }
144 :
145 : static GEN
146 133847 : frac_to_Fp(GEN a, GEN b, GEN p)
147 : {
148 133847 : GEN d = gcdii(a, b);
149 133847 : return Fp_div(diviiexact(a, d), diviiexact(b, d), p);
150 : }
151 :
152 : static GEN
153 10094 : ZpXXQ_frob(GEN S, GEN U, GEN V, long k, GEN T, ulong p, long e)
154 : {
155 10094 : pari_sp av = avma, av2;
156 10094 : long i, pr = degpol(S), dT = degpol(T), vT = varn(T);
157 10094 : GEN q = powuu(p,e);
158 10094 : GEN Tp = FpX_deriv(T, q), Tp1 = RgX_shift_shallow(Tp, 1);
159 10094 : GEN M = to_ZX(gel(S,pr+2),vT) , R;
160 10094 : av2 = avma;
161 987868 : for(i = pr-1; i>=k; i--)
162 : {
163 : GEN A, B, H, Bc;
164 : ulong v, r;
165 977774 : H = FpX_divrem(FpX_mul(V,M,q), T, q, &B);
166 977774 : A = FpX_add(FpX_mul(U,M,q), FpX_mul(H, Tp, q),q);
167 977774 : v = u_lvalrem(2*i+1,p,&r);
168 977774 : Bc = ZX_deriv(B);
169 977774 : Bc = FpX_Fp_mul(ZX_divuexact(Bc,upowuu(p,v)),Fp_divu(gen_2, r, q), q);
170 977774 : M = FpX_add(to_ZX(gel(S,i+2),vT), FpX_add(A, Bc, q), q);
171 977774 : if (gc_needed(av2,1))
172 : {
173 0 : if(DEBUGMEM>1) pari_warn(warnmem,"ZpXXQ_frob, step 1, i = %ld", i);
174 0 : M = gerepileupto(av2, M);
175 : }
176 : }
177 10094 : if (degpol(M)<dT-1)
178 5488 : return gerepileupto(av, M);
179 4606 : R = RgX_shift_shallow(M,dT-degpol(M)-2);
180 4606 : av2 = avma;
181 237629 : for(i = degpol(M)-dT+2; i>=1; i--)
182 : {
183 : GEN B, c;
184 233023 : R = RgX_shift_shallow(R, 1);
185 233023 : gel(R,2) = gel(M, i+1);
186 233023 : if (degpol(R) < dT) continue;
187 130935 : B = FpX_add(FpX_mulu(T, 2*i, q), Tp1, q);
188 130935 : c = frac_to_Fp(leading_coeff(R), leading_coeff(B), q);
189 130935 : R = FpX_sub(R, FpX_Fp_mul(B, c, q), q);
190 130935 : if (gc_needed(av2,1))
191 : {
192 0 : if(DEBUGMEM>1) pari_warn(warnmem,"ZpXXQ_frob, step 2, i = %ld", i);
193 0 : R = gerepileupto(av2, R);
194 : }
195 : }
196 4606 : if (degpol(R)==dT-1)
197 : {
198 2912 : GEN c = frac_to_Fp(leading_coeff(R), leading_coeff(Tp), q);
199 2912 : R = FpX_sub(R, FpX_Fp_mul(Tp, c, q), q);
200 2912 : return gerepileupto(av, R);
201 : } else
202 1694 : return gerepilecopy(av, R);
203 : }
204 :
205 : static GEN
206 12026 : revdigits(GEN v)
207 : {
208 12026 : long i, n = lg(v)-1;
209 12026 : GEN w = cgetg(n+2, t_POL);
210 12026 : w[1] = evalsigne(1)|evalvarn(0);
211 168784 : for (i=0; i<n; i++)
212 156758 : gel(w,i+2) = gel(v,n-i);
213 12026 : return FpXX_renormalize(w, n+2);
214 : }
215 :
216 : static GEN
217 10094 : diff_red(GEN s, GEN A, long m, GEN T, GEN p)
218 : {
219 10094 : long v, n, vT = varn(T);
220 : GEN Q, sQ, qS;
221 : pari_timer ti;
222 10094 : if (DEBUGLEVEL>1) timer_start(&ti);
223 10094 : Q = revdigits(FpX_digits(A,T,p));
224 10094 : n = degpol(Q);
225 10094 : if (DEBUGLEVEL>1) timer_printf(&ti,"reddigits");
226 10094 : sQ = FpXXQ_mul(s,Q,T,p);
227 10094 : if (DEBUGLEVEL>1) timer_printf(&ti,"redmul");
228 10094 : qS = RgX_shift_shallow(sQ,m-n);
229 10094 : v = ZX_val(sQ);
230 10094 : if (n > m + v)
231 : {
232 4564 : long i, l = n-m-v;
233 4564 : GEN rS = cgetg(l+1,t_VEC);
234 29190 : for (i = l-1; i >=0 ; i--)
235 24626 : gel(rS,i+1) = to_ZX(gel(sQ, 1+v+l-i), vT);
236 4564 : rS = FpXV_FpX_fromdigits(rS,T,p);
237 4564 : gel(qS,2) = FpX_add(FpX_mul(rS, T, p), gel(qS, 2), p);
238 4564 : if (DEBUGLEVEL>1) timer_printf(&ti,"redadd");
239 : }
240 10094 : return qS;
241 : }
242 :
243 : static GEN
244 10094 : ZC_to_padic(GEN C, GEN q)
245 : {
246 10094 : long i, l = lg(C);
247 10094 : GEN V = cgetg(l,t_COL);
248 102914 : for(i = 1; i < l; i++)
249 92820 : gel(V, i) = gadd(gel(C, i), q);
250 10094 : return V;
251 : }
252 :
253 : static GEN
254 1309 : ZM_to_padic(GEN M, GEN q)
255 : {
256 1309 : long i, l = lg(M);
257 1309 : GEN V = cgetg(l,t_MAT);
258 11403 : for(i = 1; i < l; i++)
259 10094 : gel(V, i) = ZC_to_padic(gel(M, i), q);
260 1309 : return V;
261 : }
262 :
263 : static GEN
264 1743 : ZX_to_padic(GEN P, GEN q)
265 : {
266 1743 : long i, l = lg(P);
267 1743 : GEN Q = cgetg(l, t_POL);
268 1743 : Q[1] = P[1];
269 5978 : for (i=2; i<l ;i++)
270 4235 : gel(Q,i) = gadd(gel(P,i), q);
271 1743 : return normalizepol(Q);
272 : }
273 :
274 : static GEN
275 469 : ZXC_to_padic(GEN x, GEN q)
276 2212 : { pari_APPLY_type(t_COL, ZX_to_padic(gel(x, i), q)) }
277 :
278 : static GEN
279 147 : ZXM_to_padic(GEN x, GEN q)
280 616 : { pari_APPLY_same(ZXC_to_padic(gel(x, i), q)) }
281 :
282 : static GEN
283 1309 : ZlX_hyperellpadicfrobenius(GEN H, ulong p, long n)
284 : {
285 1309 : pari_sp av = avma;
286 : long k, N, i, d;
287 : GEN F, s, Q, pN1, U, V;
288 : pari_timer ti;
289 1309 : if (typ(H) != t_POL) pari_err_TYPE("hyperellpadicfrobenius",H);
290 1309 : if (p == 2) is_sing(H, 2);
291 1309 : d = degpol(H);
292 1309 : if (d <= 0)
293 0 : pari_err_CONSTPOL("hyperellpadicfrobenius");
294 1309 : if (n < 1)
295 0 : pari_err_DOMAIN("hyperellpadicfrobenius","n","<", gen_1, utoi(n));
296 1309 : k = get_basis(p, d);
297 1309 : N = n + ulogint(2*n, p) + 1;
298 1309 : pN1 = powuu(p,N+1);
299 1309 : Q = RgX_to_FpX(H, pN1);
300 1309 : if (dvdiu(leading_coeff(Q),p)) is_sing(H, p);
301 1309 : setvarn(Q,1);
302 1309 : if (DEBUGLEVEL>1) timer_start(&ti);
303 1309 : s = revdigits(FpX_digits(RgX_inflate(Q, p), Q, pN1));
304 1309 : if (DEBUGLEVEL>1) timer_printf(&ti,"s1");
305 1309 : s = ZpXXQ_invsqrt(s, Q, p, N);
306 1309 : if (k==3)
307 35 : s = FpXXQ_mul(s, FpXXQ_sqr(s, Q, pN1), Q, pN1);
308 1309 : if (DEBUGLEVEL>1) timer_printf(&ti,"invsqrt");
309 1309 : get_UV(&U, &V, Q, p, N+1);
310 1309 : F = cgetg(d, t_MAT);
311 11403 : for (i = 1; i < d; i++)
312 : {
313 10094 : pari_sp av2 = avma;
314 : GEN M, D;
315 10094 : D = diff_red(s, monomial(utoipos(p),p*i-1,1),(k*p-1)>>1, Q, pN1);
316 10094 : if (DEBUGLEVEL>1) timer_printf(&ti,"red");
317 10094 : M = ZpXXQ_frob(D, U, V, (k-1)>>1, Q, p, N + 1);
318 10094 : if (DEBUGLEVEL>1) timer_printf(&ti,"frob");
319 10094 : gel(F, i) = gerepilecopy(av2, RgX_to_RgC(M, d-1));
320 : }
321 1309 : return gerepileupto(av, F);
322 : }
323 :
324 : GEN
325 1309 : hyperellpadicfrobenius(GEN H, ulong p, long n)
326 : {
327 1309 : pari_sp av = avma;
328 1309 : GEN M = ZlX_hyperellpadicfrobenius(H, p, n);
329 1309 : GEN q = zeropadic_shallow(utoipos(p),n);
330 1309 : return gerepileupto(av, ZM_to_padic(M, q));
331 : }
332 :
333 : INLINE GEN
334 2247 : FpXXX_renormalize(GEN x, long lx) { return ZXX_renormalize(x,lx); }
335 :
336 : static GEN
337 1806 : ZpXQXXQ_red(GEN F, GEN S, GEN T, GEN q, GEN p, long e)
338 : {
339 1806 : pari_sp av = avma;
340 1806 : long i, dF = degpol(F);
341 : GEN A, C;
342 1806 : if (signe(F)==0) return pol_0(varn(S));
343 1806 : A = cgetg(dF+3, t_POL);
344 1806 : C = pol_0(varn(S));
345 96404 : for(i=dF; i>0; i--)
346 : {
347 94598 : GEN Fi = FpXX_add(C, gel(F,i+2), q);
348 94598 : GEN R, Q = ZpXQX_divrem(Fi, S, T, q, p, e, &R);
349 94598 : gel(A,i+2) = R;
350 94598 : C = Q;
351 : }
352 1806 : gel(A,2) = FpXX_add(C, gel(F,2), q);
353 1806 : A[1] = F[1];
354 1806 : return gerepilecopy(av, FpXXX_renormalize(A,dF+3));
355 : }
356 :
357 : static GEN
358 448 : ZpXQXXQ_sqr(GEN x, GEN S, GEN T, GEN q, GEN p, long e)
359 : {
360 448 : pari_sp av = avma;
361 : GEN z, kx;
362 448 : long n = degpol(S);
363 448 : kx = RgXX_to_Kronecker(x, n);
364 448 : z = Kronecker_to_ZXX(FpXQX_sqr(kx, T, q), n, varn(S));
365 448 : return gerepileupto(av, ZpXQXXQ_red(z, S, T, q, p, e));
366 : }
367 :
368 : static GEN
369 1358 : ZpXQXXQ_mul(GEN x, GEN y, GEN S, GEN T, GEN q, GEN p, long e)
370 : {
371 1358 : pari_sp av = avma;
372 : GEN z, kx, ky;
373 1358 : long n = degpol(S);
374 1358 : kx = RgXX_to_Kronecker(x, n);
375 1358 : ky = RgXX_to_Kronecker(y, n);
376 1358 : z = Kronecker_to_ZXX(FpXQX_mul(ky, kx, T, q), n, varn(S));
377 1358 : return gerepileupto(av, ZpXQXXQ_red(z, S, T, q, p, e));
378 : }
379 :
380 : static GEN
381 441 : FpXXX_red(GEN z, GEN p)
382 : {
383 : GEN res;
384 441 : long i, l = lg(z);
385 441 : res = cgetg(l,t_POL); res[1] = z[1];
386 17388 : for (i=2; i<l; i++)
387 : {
388 16947 : GEN zi = gel(z,i);
389 16947 : if (typ(zi)==t_INT)
390 0 : gel(res,i) = modii(zi,p);
391 : else
392 16947 : gel(res,i) = FpXX_red(zi,p);
393 : }
394 441 : return FpXXX_renormalize(res,lg(res));
395 : }
396 :
397 : static GEN
398 441 : FpXXX_Fp_mul(GEN z, GEN a, GEN p)
399 : {
400 441 : return FpXXX_red(RgX_Rg_mul(z, a), p);
401 : }
402 :
403 : static GEN
404 154 : ZpXQXXQ_invsqrt(GEN F, GEN S, GEN T, ulong p, long e)
405 : {
406 154 : pari_sp av = avma, av2, av3;
407 : ulong mask;
408 154 : long v = varn(F), n=1;
409 : pari_timer ti;
410 154 : GEN a = pol_1(v), pp = utoipos(p);
411 154 : if (DEBUGLEVEL>1) timer_start(&ti);
412 154 : if (e <= 1) return gerepilecopy(av, a);
413 154 : mask = quadratic_prec_mask(e);
414 154 : av2 = avma;
415 595 : for (;mask>1;)
416 : {
417 : GEN q, q2, q22, f, fq, afq;
418 441 : long n2 = n;
419 441 : n<<=1; if (mask & 1) n--;
420 441 : mask >>= 1;
421 441 : q = powuu(p,n); q2 = powuu(p,n2);
422 441 : av3 = avma;
423 441 : f = RgX_sub(ZpXQXXQ_mul(F, ZpXQXXQ_sqr(a, S, T, q, pp, n), S, T, q, pp, n), pol_1(v));
424 441 : fq = gerepileupto(av3, RgX_Rg_divexact(f, q2));
425 441 : q22 = shifti(addiu(q2,1),-1);
426 441 : afq = FpXXX_Fp_mul(ZpXQXXQ_mul(a, fq, S, T, q2, pp, n2), q22, q2);
427 441 : a = RgX_sub(a, RgX_Rg_mul(afq, q2));
428 441 : if (gc_needed(av2,1))
429 : {
430 0 : if(DEBUGMEM>1) pari_warn(warnmem,"ZpXQXXQ_invsqrt, e = %ld", n);
431 0 : a = gerepileupto(av2, a);
432 : }
433 : }
434 154 : return gerepileupto(av, a);
435 : }
436 :
437 : static GEN
438 6573 : frac_to_Fq(GEN a, GEN b, GEN T, GEN q, GEN p, long e)
439 : {
440 6573 : GEN d = gcdii(ZX_content(a), ZX_content(b));
441 6573 : return ZpXQ_div(ZX_Z_divexact(a, d), ZX_Z_divexact(b, d), T, q, p, e);
442 : }
443 :
444 : static GEN
445 469 : ZpXQXXQ_frob(GEN F, GEN U, GEN V, long k, GEN S, GEN T, ulong p, long e)
446 : {
447 469 : pari_sp av = avma, av2;
448 469 : long i, pr = degpol(F), dS = degpol(S), v = varn(T);
449 469 : GEN q = powuu(p,e), pp = utoipos(p);
450 469 : GEN Sp = RgX_deriv(S), Sp1 = RgX_shift_shallow(Sp, 1);
451 469 : GEN M = gel(F,pr+2), R;
452 469 : av2 = avma;
453 52311 : for(i = pr-1; i>=k; i--)
454 : {
455 : GEN A, B, H, Bc;
456 : ulong v, r;
457 51842 : H = ZpXQX_divrem(FpXQX_mul(V, M, T, q), S, T, q, utoipos(p), e, &B);
458 51842 : A = FpXX_add(FpXQX_mul(U, M, T, q), FpXQX_mul(H, Sp, T, q),q);
459 51842 : v = u_lvalrem(2*i+1,p,&r);
460 51842 : Bc = RgX_deriv(B);
461 51842 : Bc = FpXX_Fp_mul(ZXX_Z_divexact(Bc,powuu(p,v)), Fp_divu(gen_2, r, q), q);
462 51842 : M = FpXX_add(gel(F,i+2), FpXX_add(A, Bc, q), q);
463 51842 : if (gc_needed(av2,1))
464 : {
465 0 : if(DEBUGMEM>1) pari_warn(warnmem,"ZpXQXXQ_frob, step 1, i = %ld", i);
466 0 : M = gerepileupto(av2, M);
467 : }
468 : }
469 469 : if (degpol(M)<dS-1)
470 266 : return gerepileupto(av, M);
471 203 : R = RgX_shift_shallow(M,dS-degpol(M)-2);
472 203 : av2 = avma;
473 7175 : for(i = degpol(M)-dS+2; i>=1; i--)
474 : {
475 : GEN B, c;
476 6972 : R = RgX_shift_shallow(R, 1);
477 6972 : gel(R,2) = gel(M, i+1);
478 6972 : if (degpol(R) < dS) continue;
479 6412 : B = FpXX_add(FpXX_mulu(S, 2*i, q), Sp1, q);
480 6412 : c = frac_to_Fq(to_ZX(leading_coeff(R),v), to_ZX(leading_coeff(B),v), T, q, pp, e);
481 6412 : R = FpXX_sub(R, FpXQX_FpXQ_mul(B, c, T, q), q);
482 6412 : if (gc_needed(av2,1))
483 : {
484 0 : if(DEBUGMEM>1) pari_warn(warnmem,"ZpXXQ_frob, step 2, i = %ld", i);
485 0 : R = gerepileupto(av2, R);
486 : }
487 : }
488 203 : if (degpol(R)==dS-1)
489 : {
490 161 : GEN c = frac_to_Fq(to_ZX(leading_coeff(R),v), to_ZX(leading_coeff(Sp),v), T, q, pp, e);
491 161 : R = FpXX_sub(R, FpXQX_FpXQ_mul(Sp, c, T, q), q);
492 161 : return gerepileupto(av, R);
493 : } else
494 42 : return gerepilecopy(av, R);
495 : }
496 :
497 : static GEN
498 469 : Fq_diff_red(GEN s, GEN A, long m, GEN S, GEN T, GEN q, GEN p, long e)
499 : {
500 : long v, n;
501 : GEN Q, sQ, qS;
502 : pari_timer ti;
503 469 : if (DEBUGLEVEL>1) timer_start(&ti);
504 469 : Q = revdigits(ZpXQX_digits(A, S, T, q, p, e));
505 469 : n = degpol(Q);
506 469 : if (DEBUGLEVEL>1) timer_printf(&ti,"reddigits");
507 469 : sQ = ZpXQXXQ_mul(s, Q, S, T, q, p, e);
508 469 : if (DEBUGLEVEL>1) timer_printf(&ti,"redmul");
509 469 : qS = RgX_shift_shallow(sQ,m-n);
510 469 : v = ZX_val(sQ);
511 469 : if (n > m + v)
512 : {
513 189 : long i, l = n-m-v;
514 189 : GEN rS = cgetg(l+1,t_VEC);
515 1547 : for (i = l-1; i >=0 ; i--)
516 1358 : gel(rS,i+1) = gel(sQ, 1+v+l-i);
517 189 : rS = FpXQXV_FpXQX_fromdigits(rS, S, T, q);
518 189 : gel(qS,2) = FpXX_add(FpXQX_mul(rS, S, T, q), gel(qS, 2), q);
519 189 : if (DEBUGLEVEL>1) timer_printf(&ti,"redadd");
520 : }
521 469 : return qS;
522 : }
523 :
524 : static void
525 154 : Fq_get_UV(GEN *U, GEN *V, GEN S, GEN T, ulong p, long e)
526 : {
527 154 : GEN q = powuu(p, e), pp = utoipos(p), d;
528 154 : GEN dS = RgX_deriv(S), R = polresultantext(S, dS), C;
529 154 : long v = varn(S);
530 154 : if (signe(FpX_red(to_ZX(gel(R,3),v), pp))==0) is_sing(S, p);
531 147 : C = FpXQ_red(to_ZX(gel(R, 3),v), T, q);
532 147 : d = ZpXQ_inv(C, T, pp, e);
533 147 : *U = FpXQX_FpXQ_mul(FpXQX_red(to_ZX(gel(R,1),v),T,q),d,T,q);
534 147 : *V = FpXQX_FpXQ_mul(FpXQX_red(to_ZX(gel(R,2),v),T,q),d,T,q);
535 147 : }
536 :
537 : static GEN
538 469 : ZXX_to_FpXC(GEN x, long N, GEN p, long v)
539 : {
540 : long i, l;
541 : GEN z;
542 469 : l = lg(x)-1; x++;
543 469 : if (l > N+1) l = N+1; /* truncate higher degree terms */
544 469 : z = cgetg(N+1,t_COL);
545 2170 : for (i=1; i<l ; i++)
546 : {
547 1701 : GEN xi = gel(x, i);
548 1701 : gel(z,i) = typ(xi)==t_INT? scalarpol(Fp_red(xi, p), v): FpX_red(xi, p);
549 : }
550 511 : for ( ; i<=N ; i++)
551 42 : gel(z,i) = pol_0(v);
552 469 : return z;
553 : }
554 :
555 : GEN
556 154 : ZlXQX_hyperellpadicfrobenius(GEN H, GEN T, ulong p, long n)
557 : {
558 154 : pari_sp av = avma;
559 : long k, N, i, d, N1;
560 : GEN xp, F, s, q, Q, pN1, U, V, pp;
561 : pari_timer ti;
562 154 : if (typ(H) != t_POL) pari_err_TYPE("hyperellpadicfrobenius",H);
563 154 : if (p == 2) is_sing(H, 2);
564 154 : d = degpol(H);
565 154 : if (d <= 0) pari_err_CONSTPOL("hyperellpadicfrobenius");
566 154 : if (n < 1) pari_err_DOMAIN("hyperellpadicfrobenius","n","<", gen_1, utoi(n));
567 154 : k = get_basis(p, d); pp = utoipos(p);
568 154 : N = n + ulogint(2*n, p) + 1;
569 154 : q = powuu(p,n); N1 = N+1;
570 154 : pN1 = powuu(p,N1); T = FpX_get_red(T, pN1);
571 154 : Q = RgX_to_FqX(H, T, pN1);
572 154 : if (signe(FpX_red(to_ZX(leading_coeff(Q),varn(Q)),pp))==0) is_sing(H, p);
573 154 : if (DEBUGLEVEL>1) timer_start(&ti);
574 154 : xp = ZpX_Frobenius(T, pp, N1);
575 154 : s = RgX_inflate(FpXY_FpXQ_evalx(Q, xp, T, pN1), p);
576 154 : s = revdigits(ZpXQX_digits(s, Q, T, pN1, pp, N1));
577 154 : if (DEBUGLEVEL>1) timer_printf(&ti,"s1");
578 154 : s = ZpXQXXQ_invsqrt(s, Q, T, p, N);
579 154 : if (k==3)
580 7 : s = ZpXQXXQ_mul(s, ZpXQXXQ_sqr(s, Q, T, pN1, pp, N1), Q, T, pN1, pp, N1);
581 154 : if (DEBUGLEVEL>1) timer_printf(&ti,"invsqrt");
582 154 : Fq_get_UV(&U, &V, Q, T, p, N+1);
583 147 : if (DEBUGLEVEL>1) timer_printf(&ti,"get_UV");
584 147 : F = cgetg(d, t_MAT);
585 616 : for (i = 1; i < d; i++)
586 : {
587 469 : pari_sp av2 = avma;
588 : GEN M, D;
589 469 : D = Fq_diff_red(s, monomial(pp,p*i-1,1),(k*p-1)>>1, Q, T, pN1, pp, N1);
590 469 : if (DEBUGLEVEL>1) timer_printf(&ti,"red");
591 469 : M = ZpXQXXQ_frob(D, U, V, (k - 1)>>1, Q, T, p, N1);
592 469 : if (DEBUGLEVEL>1) timer_printf(&ti,"frob");
593 469 : gel(F, i) = gerepileupto(av2, ZXX_to_FpXC(M, d-1, q, varn(T)));
594 : }
595 147 : return gerepileupto(av, F);
596 : }
597 :
598 : GEN
599 154 : nfhyperellpadicfrobenius(GEN H, GEN T, ulong p, long n)
600 : {
601 154 : pari_sp av = avma;
602 154 : GEN pp = utoipos(p), q = zeropadic_shallow(pp, n);
603 154 : GEN M = ZlXQX_hyperellpadicfrobenius(lift_shallow(H),T,p,n);
604 147 : GEN MM = ZpXQM_prodFrobenius(M, T, pp, n);
605 147 : GEN m = gmul(ZXM_to_padic(MM, q), gmodulo(gen_1, T));
606 147 : return gerepileupto(av, m);
607 : }
608 :
609 : GEN
610 595 : hyperellpadicfrobenius0(GEN H, GEN Tp, long n)
611 : {
612 : GEN T, p;
613 595 : if (!ff_parse_Tp(Tp, &T,&p,0)) pari_err_TYPE("hyperellpadicfrobenius", Tp);
614 595 : if (lgefint(p) > 3) pari_err_IMPL("large prime in hyperellpadicfrobenius");
615 7 : return T? nfhyperellpadicfrobenius(H, T, itou(p), n)
616 602 : : hyperellpadicfrobenius(H, itou(p), n);
617 : }
618 :
619 : static GEN
620 84 : F2x_genus2charpoly_naive(GEN P, GEN Q)
621 : {
622 84 : long a, b = 1, c = 0;
623 84 : GEN T = mkvecsmall2(P[1], 7);
624 84 : GEN PT = F2x_rem(P, T), QT = F2x_rem(Q, T);
625 84 : long q0 = F2x_eval(Q, 0), q1 = F2x_eval(Q, 1);
626 84 : long dP = F2x_degree(P), dQ = F2x_degree(Q);
627 84 : a= dQ<3 ? 0: dP<=5 ? 1: -1;
628 84 : a += (q0? F2x_eval(P, 0)? -1: 1: 0) + (q1? F2x_eval(P, 1)? -1: 1: 0);
629 84 : b += q0 + q1;
630 84 : if (lgpol(QT))
631 70 : c = (F2xq_trace(F2xq_div(PT, F2xq_sqr(QT, T), T), T)==0 ? 1: -1);
632 84 : return mkvecsmalln(6, 0UL, 4UL, 2*a, (b+2*c+a*a)>>1, a, 1UL);
633 : }
634 :
635 : static GEN
636 252 : Flx_difftable(GEN P, ulong p)
637 : {
638 252 : long i, n = degpol(P);
639 252 : GEN V = cgetg(n+2, t_VEC);
640 252 : gel(V, n+1) = P;
641 1729 : for(i = n; i >= 1; i--)
642 1477 : gel(V, i) = Flx_diff1(gel(V, i+1), p);
643 252 : return V;
644 : }
645 :
646 : static GEN
647 1561 : FlxV_Fl2_eval_pre(GEN V, GEN x, ulong D, ulong p, ulong pi)
648 : {
649 1561 : long i, n = lg(V)-1;
650 1561 : GEN r = cgetg(n+1, t_VEC);
651 11963 : for (i = 1; i <= n; i++)
652 10402 : gel(r, i) = Flx_Fl2_eval_pre(gel(V, i), x, D, p, pi);
653 1561 : return r;
654 : }
655 :
656 : static GEN
657 44590 : Fl2V_next(GEN V, ulong p)
658 : {
659 44590 : long i, n = lg(V)-1;
660 44590 : GEN r = cgetg(n+1, t_VEC);
661 44590 : gel(r, 1) = gel(V, 1);
662 286300 : for (i = 2; i <= n; i++)
663 241710 : gel(r, i) = Flv_add(gel(V, i), gel(V, i-1), p);
664 44590 : return r;
665 : }
666 :
667 : static GEN
668 252 : Flx_genus2charpoly_naive(GEN H, ulong p)
669 : {
670 252 : pari_sp av = avma, av2;
671 252 : ulong pi = get_Fl_red(p);
672 252 : ulong i, j, p2 = p>>1, D = 2, e = ((p&2UL) == 0) ? -1 : 1;
673 252 : long a, b, c = 0, n = degpol(H);
674 252 : GEN t, k = const_vecsmall(p, -1);
675 252 : k[1] = 0;
676 1813 : for (i=1, j=1; i < p; i += 2, j = Fl_add(j, i, p)) k[j+1] = 1;
677 294 : while (k[1+D] >= 0) D++;
678 252 : b = n == 5 ? 0 : 1;
679 252 : a = b ? k[1+Flx_lead(H)]: 0;
680 252 : t = Flx_difftable(H, p);
681 252 : av2 = avma;
682 3626 : for (i=0; i < p; i++)
683 : {
684 3374 : ulong v = Flx_eval(H, i, p);
685 3374 : a += k[1+v];
686 3374 : b += !!v;
687 : }
688 1813 : for (j=1; j <= p2; j++)
689 : {
690 1561 : GEN V = FlxV_Fl2_eval_pre(t, mkvecsmall2(0, j), D, p, pi);
691 1561 : for (i=0;; i++)
692 44590 : {
693 46151 : GEN r2 = gel(V, n+1);
694 92302 : c += uel(r2,2) ?
695 43890 : (uel(r2,1) ? uel(k,1+Fl2_norm_pre(r2, D, p, pi)): e)
696 90041 : : !!uel(r2,1);
697 46151 : if (i == p-1) break;
698 44590 : V = Fl2V_next(V, p);
699 : }
700 1561 : set_avma(av2);
701 : }
702 252 : set_avma(av);
703 252 : return mkvecsmalln(6, 0UL, p*p, a*p, (b+2*c+a*a)>>1, a, 1UL);
704 : }
705 :
706 : static GEN
707 679 : charpoly_funceq(GEN P, GEN q)
708 : {
709 679 : long i, l, g = degpol(P)>>1;
710 679 : GEN R, Q = gpowers0(q, g-1, q); /* Q[i] = q^i, i <= g */
711 679 : R = cgetg_copy(P, &l); R[1] = P[1];
712 3164 : for (i=0; i<g; i++) gel(R, i+2) = mulii(gel(P, 2*g-i+2), gel(Q, g-i));
713 3843 : for (; i<=2*g; i++) gel(R, i+2) = icopy(gel(P, i+2));
714 679 : return R;
715 : }
716 :
717 : static long
718 686 : hyperell_Weil_bound(GEN q, ulong g, GEN p)
719 : {
720 686 : pari_sp av = avma;
721 686 : GEN w = mulii(binomialuu(2*g,g),sqrtint(shifti(powiu(q, g),2)));
722 686 : return gc_long(av, logint(w,p) + 1);
723 : }
724 :
725 : /* return 4P + Q^2 */
726 : static GEN
727 287368 : check_hyperell(GEN PQ)
728 : {
729 : GEN H;
730 287368 : if (is_vec_t(typ(PQ)) && lg(PQ)==3)
731 225651 : H = gadd(gsqr(gel(PQ, 2)), gmul2n(gel(PQ, 1), 2));
732 : else
733 61717 : H = gmul2n(PQ, 2);
734 287368 : return typ(H) == t_POL? H: NULL;
735 : }
736 :
737 : GEN
738 1029 : hyperellcharpoly(GEN PQ)
739 : {
740 1029 : pari_sp av = avma;
741 1029 : GEN M, R, T=NULL, pp=NULL, q;
742 1029 : long d, n, eps = 0;
743 : ulong p;
744 1029 : GEN H = check_hyperell(PQ);
745 1029 : if (!H || !RgX_is_FpXQX(H, &T, &pp) || !pp)
746 0 : pari_err_TYPE("hyperellcharpoly", PQ);
747 1029 : p = itou(pp);
748 1029 : if (!T)
749 : {
750 882 : if (p==2 && is_vec_t(typ(PQ)))
751 : {
752 84 : long dP, dQ, v = varn(H);
753 84 : GEN P = gel(PQ,1), Q = gel(PQ,2);
754 84 : if (typ(P)!=t_POL) P = scalarpol(P, v);
755 84 : if (typ(Q)!=t_POL) Q = scalarpol(Q, v);
756 84 : dP = degpol(P); dQ = degpol(Q);
757 84 : if (dP<=6 && dQ <=3 && (dQ==3 || dP>=5))
758 : {
759 84 : GEN P2 = RgX_to_F2x(P), Q2 = RgX_to_F2x(Q);
760 84 : GEN D = F2x_add(F2x_mul(P2, F2x_sqr(F2x_deriv(Q2))), F2x_sqr(F2x_deriv(P2)));
761 84 : if (F2x_degree(F2x_gcd(D, Q2))) is_sing(PQ, 2);
762 84 : if (dP==6 && dQ<3 && F2x_coeff(P2,5)==F2x_coeff(Q2,2))
763 0 : is_sing(PQ, 2); /* The curve is singular at infinity */
764 84 : R = zx_to_ZX(F2x_genus2charpoly_naive(P2, Q2));
765 84 : return gerepileupto(av, R);
766 : }
767 : }
768 798 : H = RgX_to_FpX(H, pp);
769 798 : d = degpol(H);
770 798 : if (d <= 0) is_sing(H, p);
771 798 : if (p > 2 && ((d == 5 && p < 17500) || (d == 6 && p < 24500)))
772 : {
773 259 : GEN Hp = ZX_to_Flx(H, p);
774 259 : if (!Flx_is_squarefree(Hp, p)) is_sing(H, p);
775 252 : R = zx_to_ZX(Flx_genus2charpoly_naive(Hp, p));
776 252 : return gerepileupto(av, R);
777 : }
778 539 : n = hyperell_Weil_bound(pp, (d-1)>>1, pp);
779 539 : eps = odd(d)? 0: Fp_issquare(leading_coeff(H), pp);
780 539 : M = hyperellpadicfrobenius(H, p, n);
781 539 : R = centerlift(carberkowitz(M, 0));
782 539 : q = pp;
783 : }
784 : else
785 : {
786 : int fixvar;
787 147 : T = typ(T)==t_FFELT? FF_mod(T): RgX_to_FpX(T, pp);
788 147 : q = powuu(p, degpol(T));
789 147 : fixvar = (varncmp(varn(T),varn(H)) <= 0);
790 147 : if (fixvar) setvarn(T, fetch_var());
791 147 : H = RgX_to_FpXQX(H, T, pp);
792 147 : d = degpol(H);
793 147 : if (d <= 0) is_sing(H, p);
794 147 : eps = odd(d)? 0: Fq_issquare(leading_coeff(H), T, pp);
795 147 : n = hyperell_Weil_bound(q, (d-1)>>1, pp);
796 147 : M = nfhyperellpadicfrobenius(H, T, p, n);
797 140 : R = simplify_shallow(centerlift(liftpol_shallow(carberkowitz(M, 0))));
798 140 : if (fixvar) (void)delete_var();
799 : }
800 679 : if (!odd(d))
801 : {
802 301 : GEN b = get_basis(p, d) == 3 ? gen_1 : q;
803 301 : GEN pn = powuu(p, n);
804 301 : R = FpX_div_by_X_x(R, eps? b: negi(b), pn, NULL);
805 301 : R = FpX_center_i(R, pn, shifti(pn,-1));
806 : }
807 679 : return gerepileupto(av, charpoly_funceq(R, q));
808 : }
809 :
810 : int
811 3493 : hyperellisoncurve(GEN W, GEN P)
812 : {
813 3493 : pari_sp av = avma;
814 : long res;
815 : GEN x, y;
816 3493 : if (typ(P)!=t_VEC || lg(P)!=3) pari_err_TYPE("hyperellisoncurve",P);
817 3493 : x = gel(P,1); y = gel(P,2);
818 3493 : if (typ(W)==t_POL)
819 0 : res = gequal(gsqr(y), poleval(W,x));
820 : else
821 : {
822 3493 : if (typ(W)!=t_VEC || lg(W)!=3) pari_err_TYPE("hyperellisoncurve",W);
823 3493 : res = gequal(gmul(y, gadd(y,poleval(gel(W,2), x))), poleval(gel(W,1), x));
824 : }
825 3493 : return gc_int(av, res);
826 : }
827 :
828 : GEN
829 35 : hyperellordinate(GEN W, GEN x)
830 : {
831 35 : pari_sp av = avma;
832 35 : if (typ(W)==t_POL)
833 : {
834 14 : GEN d = poleval(W,x), y;
835 14 : if (gequal0(d)) { return gerepilecopy(av, mkvec(d)); }
836 14 : if (!issquareall(d, &y)) { set_avma(av); return cgetg(1,t_VEC); }
837 7 : return gerepilecopy(av, mkvec2(y, gneg(y)));
838 : }
839 : else
840 : {
841 : GEN b, c, d, rd, y;
842 21 : if (typ(W)!=t_VEC || lg(W)!=3) pari_err_TYPE("hyperellisoncurve",W);
843 21 : b = poleval(gel(W,2), x); c = poleval(gel(W,1), x);
844 21 : d = gadd(gsqr(b), gmul2n(c, 2));
845 21 : if (gequal0(d)) { return gerepilecopy(av, mkvec(gmul2n(gneg(b),-1))); }
846 14 : if (!issquareall(d, &rd)) { set_avma(av); return cgetg(1,t_VEC); }
847 7 : y = gmul2n(gsub(rd, b), -1);
848 7 : return gerepilecopy(av, mkvec2(y, gsub(y,rd)));
849 : }
850 : }
851 :
852 : GEN
853 118068 : hyperelldisc(GEN PQ)
854 : {
855 118068 : pari_sp av = avma;
856 118068 : GEN D, H = check_hyperell(PQ);
857 : long d, g;
858 118068 : if (!H || signe(H)==0) pari_err_TYPE("hyperelldisc",PQ);
859 118068 : d = degpol(H); g = ((d+1)>>1)-1;
860 118068 : D = gmul2n(RgX_disc(H),-4*(g+1));
861 118068 : if (odd(d)) D = gmul(D, gsqr(leading_coeff(H)));
862 118068 : return gerepileupto(av, D);
863 : }
864 :
865 : static long
866 126502 : get_ep(GEN W)
867 : {
868 126502 : GEN P = gel(W,1), Q = gel(W,2);
869 126502 : if (signe(Q)==0) return ZX_lval(P,2);
870 86307 : return minss(ZX_lval(P,2), ZX_lval(Q,2));
871 : }
872 :
873 : static GEN
874 50798 : algo51(GEN W, GEN M)
875 : {
876 50798 : GEN P = gel(W,1), Q = gel(W,2);
877 : for(;;)
878 10640 : {
879 61438 : long vP = ZX_lval(P,2);
880 61438 : long vQ = signe(Q) ? ZX_lval(Q,2): vP+1;
881 : long r;
882 : /* 1 */
883 61438 : if (vQ==0) break;
884 : /* 2 */
885 36120 : if (vP==0)
886 : {
887 : GEN H, H1;
888 : /* a */
889 29602 : RgX_even_odd(FpX_red(P,gen_2),&H, &H1);
890 29602 : if (signe(H1)) break;
891 : /* b */
892 14951 : P = ZX_add(P, ZX_mul(H, ZX_sub(Q, H)));
893 14951 : Q = ZX_sub(Q, ZX_mulu(H, 2));
894 14951 : vP = ZX_lval(P,2);
895 14951 : vQ = signe(Q) ? ZX_lval(Q,2): vP+1;
896 : }
897 : /* 2c */
898 21469 : if (vP==1) break;
899 : /* 2d */
900 10640 : r = minss(2*vQ, vP)>>1;
901 10640 : if (M) gel(M,1) = shifti(gel(M,1), r);
902 10640 : P = ZX_shifti(P, -2*r);
903 10640 : Q = ZX_shifti(Q, -r);
904 : }
905 50798 : return mkvec2(P,Q);
906 : }
907 :
908 : static GEN
909 103343 : algo52(GEN W, GEN c, long *pt_lambda)
910 : {
911 : long lambda;
912 103343 : GEN P = gel(W,1), Q = gel(W,2);
913 : for(;;)
914 116959 : {
915 : GEN H, H1;
916 : /* 1 */
917 220302 : GEN Pc = ZX_affine(P,gen_2,c), Qc = ZX_affine(Q,gen_2,c);
918 220302 : long mP = ZX_lval(Pc,2), mQ = signe(Qc) ? ZX_lval(Qc,2): mP+1;
919 : /* 2 */
920 220302 : if (2*mQ <= mP) { lambda = 2*mQ; break; }
921 : /* 3 */
922 188031 : if (odd(mP)) { lambda = mP; break; }
923 : /* 4 */
924 127718 : RgX_even_odd(FpX_red(ZX_shifti(Pc, -mP),gen_2),&H, &H1);
925 127718 : if (signe(H1)) { lambda = mP; break; }
926 : /* 5 */
927 116959 : P = ZX_add(P, ZX_mul(H, ZX_sub(Q, H)));
928 116959 : Q = ZX_sub(Q, ZX_mulu(H, 2));
929 : }
930 103343 : *pt_lambda = lambda;
931 103343 : return mkvec2(P,Q);
932 : }
933 :
934 : static long
935 146791 : test53(long lambda, long ep, long g)
936 : {
937 146791 : return (lambda <= g+1) || (odd(g) && lambda<g+3 && ep==1);
938 : }
939 :
940 : static long
941 189382 : test55(GEN W, long ep, long g)
942 : {
943 189382 : GEN P = gel(W,1), Q = gel(W,2);
944 189382 : GEN Pe = FpX_red(ep ? ZX_shifti(P,-1): P, gen_2);
945 189382 : GEN Qe = FpX_red(ep ? ZX_shifti(Q,-1): Q, gen_2);
946 189382 : if (ep==0)
947 : {
948 149101 : if (signe(Qe)!=0) return ZX_val(Qe) >= (g + 3)>>1;
949 90723 : else return ZX_val(FpX_deriv(Pe, gen_2)) >= g+1;
950 : }
951 : else
952 40281 : return ZX_val(Qe) >= (g+1)>>1 && ZX_val(Pe) >= g + 1;
953 : }
954 :
955 : static GEN
956 50798 : hyperell_reverse(GEN W, long g)
957 : {
958 50798 : return mkvec2(RgXn_recip_shallow(gel(W,1),2*g+3),
959 50798 : RgXn_recip_shallow(gel(W,2),g+2));
960 : }
961 :
962 : static GEN
963 50784 : algo56(GEN W, long g)
964 : {
965 : long ep;
966 50784 : GEN M = mkvec2(gen_1, matid(2)), Woo;
967 50784 : W = algo51(W, M);
968 50784 : Woo = hyperell_reverse(W, g);
969 50784 : ep = get_ep(Woo);
970 50784 : if (test55(Woo,ep,g))
971 : {
972 : long lambda;
973 11737 : Woo = algo52(Woo, gen_0, &lambda);
974 11737 : if (!test53(lambda,ep,g))
975 : {
976 5969 : long r = lambda>>1;
977 5969 : gel(M,1) = shifti(gel(M,1), r);
978 5969 : gel(M,2) = ZM2_mul(gel(M,2), mkmat22(gen_0, gen_1, gen_2, gen_0));
979 5969 : W = mkvec2(ZX_shifti(ZX_unscale(gel(Woo,1), gen_2), -2*r),
980 5969 : ZX_shifti(ZX_unscale(gel(Woo,2), gen_2), -r));
981 : }
982 : }
983 : for(;;)
984 24892 : {
985 75676 : long j, ep = get_ep(W);
986 189340 : for (j = 0; j < 2; j++)
987 : {
988 : long lambda;
989 138556 : GEN c = utoi(j);
990 138556 : GEN Pc = ZX_affine(gel(W,1), gen_2, c);
991 138556 : GEN Qc = ZX_affine(gel(W,2), gen_2, c);
992 138556 : if (test55(mkvec2(Pc, Qc), ep, g))
993 : {
994 91564 : GEN Wc = algo52(W, c, &lambda);
995 91564 : if (!test53(lambda,ep,g))
996 : {
997 24892 : long r = lambda>>1;
998 24892 : gel(M,1) = shifti(gel(M,1), r);
999 24892 : gel(M,2) = ZM2_mul(gel(M,2), mkmat22(gen_2, c, gen_0, gen_1));
1000 24892 : W = mkvec2(ZX_shifti(ZX_affine(gel(Wc,1), gen_2,c), -2*r),
1001 24892 : ZX_shifti(ZX_affine(gel(Wc,2), gen_2,c), -r));
1002 24892 : break;
1003 : }
1004 : }
1005 : }
1006 75676 : if (j==2) break;
1007 : }
1008 50784 : return mkvec2(W, M);
1009 : }
1010 :
1011 : static GEN
1012 14 : algo56bis(GEN W, long g, long inf, long thr)
1013 : {
1014 14 : pari_sp av = avma;
1015 14 : GEN vl = cgetg(3,t_VEC);
1016 14 : long nl = 1;
1017 14 : W = algo51(W, NULL);
1018 14 : if (inf)
1019 : {
1020 14 : GEN Woo = hyperell_reverse(W, g);
1021 14 : GEN Pc = ZX_unscale(gel(W,1), gen_2);
1022 14 : GEN Qc = ZX_unscale(gel(W,2), gen_2);
1023 14 : long ep = get_ep(Woo);
1024 14 : if (test55(mkvec2(Pc, Qc), ep, g))
1025 : {
1026 : long lambda;
1027 14 : Woo = algo52(Woo, gen_0, &lambda);
1028 14 : if (lambda == thr)
1029 : {
1030 0 : long r = lambda>>1;
1031 0 : gel(vl,nl++) = mkvec2(ZX_shifti(ZX_unscale(gel(Woo,1), gen_2), -2*r),
1032 0 : ZX_shifti(ZX_unscale(gel(Woo,2), gen_2), -r));
1033 : }
1034 : }
1035 : }
1036 : {
1037 14 : long j, ep = get_ep(W);
1038 42 : for (j = 0; j < 2; j++)
1039 : {
1040 : long lambda;
1041 28 : GEN c = utoi(j);
1042 28 : GEN Pc = ZX_affine(gel(W,1), gen_2, c);
1043 28 : GEN Qc = ZX_affine(gel(W,2), gen_2, c);
1044 28 : if (test55(mkvec2(Pc, Qc), ep, g))
1045 : {
1046 28 : GEN Wc = algo52(W, c, &lambda);
1047 28 : if (lambda == thr)
1048 : {
1049 0 : long r = lambda>>1;
1050 0 : gel(vl,nl++) = mkvec2(ZX_shifti(ZX_affine(gel(Wc,1), gen_2,c), -2*r),
1051 0 : ZX_shifti(ZX_affine(gel(Wc,2), gen_2,c), -r));
1052 : }
1053 : }
1054 : }
1055 : }
1056 14 : setlg(vl, nl);
1057 14 : return gerepilecopy(av,vl);
1058 : }
1059 :
1060 : /* return the (degree 2) apolar invariant (the nth transvectant of P and P) */
1061 : static GEN
1062 84 : ZX_apolar(GEN P, long n)
1063 : {
1064 84 : pari_sp av = avma;
1065 84 : long d = degpol(P), i;
1066 84 : GEN s = gen_0, g = cgetg(n+2,t_VEC);
1067 84 : gel(g,1) = gen_1;
1068 588 : for (i = 1; i <= n; i++) gel(g,i+1) = muliu(gel(g,i),i); /* g[i+1] = i! */
1069 658 : for (i = n-d; i <= d; i++)
1070 : {
1071 574 : GEN a = mulii(mulii(gel(g,i+1),gel(g,n-i+1)),
1072 574 : mulii(gel(P,i+2),gel(P,n-i+2)));
1073 574 : s = odd(i)? subii(s, a): addii(s, a);
1074 : }
1075 84 : return gerepileuptoint(av,s);
1076 : }
1077 :
1078 : static GEN
1079 52135 : algo57(GEN F, long g, GEN pr)
1080 : {
1081 : long i, l;
1082 52135 : GEN D, C = content(F);
1083 52135 : GEN e = gel(core2(shifti(C,-vali(C))),2);
1084 52135 : GEN M = mkvec2(e, matid(2));
1085 52135 : long minvd = (2*g+1)>>(odd(g) ? 4:2);
1086 52135 : F = ZX_Z_divexact(F, sqri(e));
1087 52135 : D = absi(hyperelldisc(F));
1088 52135 : if (!pr)
1089 : {
1090 84 : GEN A = gcdii(D, ZX_apolar(F, 2*g+2));
1091 84 : pr = gel(factor(shifti(A, -vali(A))),1);
1092 : }
1093 52135 : l = lg(pr);
1094 311008 : for (i = 1; i < l; i++)
1095 : {
1096 : long ep;
1097 258873 : GEN p = gel(pr, i), ps2 = shifti(p,-1), Fe;
1098 258873 : if (equaliu(p,2) || Z_pval(D,p) < minvd) continue;
1099 196658 : ep = ZX_pvalrem(F,p, &Fe); Fe = FpX_red(Fe, p);
1100 196658 : if (degpol(Fe) < g+1+ep)
1101 : {
1102 6406 : GEN Fi = ZX_unscale(RgXn_recip_shallow(F,2*g+3), p);
1103 6406 : long lambda = ZX_pval(Fi,p);
1104 6406 : if (!test53(lambda,ep,g))
1105 : {
1106 3815 : GEN ppr = powiu(p,lambda>>1);
1107 3815 : F = ZX_Z_divexact(Fi,sqri(ppr));
1108 3815 : gel(M,1) = mulii(gel(M,1), ppr);
1109 3815 : gel(M,2) = ZM2_mul(gel(M,2), mkmat22(gen_0,gen_1,p,gen_0));
1110 : }
1111 : }
1112 : for(;;)
1113 24997 : {
1114 : GEN Fe, R;
1115 221655 : long j, lR, ep = ZX_pvalrem(F,p, &Fe);
1116 221655 : R = FpX_roots_mult(FpX_red(Fe, p), g+2-ep, p); lR = lg(R);
1117 233742 : for (j = 1; j<lR; j++)
1118 : {
1119 37084 : GEN c = Fp_center(gel(R,j), p, ps2);
1120 37084 : GEN Fi = ZX_affine(F,p,c);
1121 37084 : long lambda = ZX_pval(Fi,p);
1122 37084 : if (!test53(lambda,ep,g))
1123 : {
1124 24997 : GEN ppr = powiu(p,lambda>>1);
1125 24997 : F = ZX_Z_divexact(Fi, sqri(ppr));
1126 24997 : gel(M,1) = mulii(gel(M,1), ppr);
1127 24997 : gel(M,2) = ZM2_mul(gel(M,2), mkmat22(p,c,gen_0,gen_1));
1128 24997 : break;
1129 : }
1130 : }
1131 221655 : if (j==lR) break;
1132 : }
1133 : }
1134 52135 : return mkvec2(F, M);
1135 : }
1136 :
1137 : /* if inf=0, ignore point at infinity */
1138 : static GEN
1139 42 : algo57bis(GEN F, long g, GEN p, long inf, long thr)
1140 : {
1141 42 : pari_sp av = avma;
1142 42 : GEN vl = cgetg(3,t_VEC), Fe;
1143 42 : long nl = 1, ep = ZX_pvalrem(F,p, &Fe);
1144 42 : Fe = FpX_red(Fe, p);
1145 : {
1146 42 : GEN R = FpX_roots_mult(Fe, thr-ep, p);
1147 42 : long j, lR = lg(R);
1148 63 : for (j = 1; j<lR; j++)
1149 : {
1150 21 : GEN Fj = ZX_affine(F, p, gel(R,j));
1151 21 : long lambda = ZX_pvalrem(Fj, p, &Fj);
1152 21 : if (lambda == thr) gel(vl,nl++) = odd(lambda)? ZX_Z_mul(Fj, p): Fj;
1153 : }
1154 : }
1155 42 : if (inf==1 && 2*g+2-degpol(Fe) >= thr-ep)
1156 : {
1157 0 : GEN Fj = ZX_unscale(RgXn_recip_shallow(F,2*g+3), p);
1158 0 : long lambda = ZX_pvalrem(Fj, p, &Fj);
1159 0 : if (lambda == thr) gel(vl,nl++) = odd(lambda)? ZX_Z_mul(Fj, p): Fj;
1160 : }
1161 42 : setlg(vl, nl);
1162 42 : return gerepilecopy(av,vl);
1163 : }
1164 :
1165 : static GEN
1166 56 : next_model(GEN G, long g, GEN p, long inf, long thr)
1167 : {
1168 70 : return equaliu(p,2) ? algo56bis(G, g, inf, thr)
1169 70 : : algo57bis(G, g, p, inf, thr);
1170 : }
1171 :
1172 : static GEN
1173 0 : get_extremal(GEN F, GEN G, long g, GEN p)
1174 : {
1175 : while (1)
1176 0 : {
1177 : GEN Wi;
1178 0 : Wi = next_model(G, g, p, 0, g+2);
1179 0 : if (lg(Wi)==1) return F;
1180 0 : F = gel(Wi,1);
1181 0 : Wi = next_model(F, g, p, 0, g+1);
1182 0 : if (lg(Wi)==1) return F;
1183 0 : G = gel(Wi,1);
1184 : }
1185 : }
1186 :
1187 : GEN
1188 63 : hyperellextremalmodels(GEN F, long g, GEN p)
1189 : {
1190 63 : pari_sp av = avma;
1191 : GEN R, W;
1192 : long l;
1193 63 : if (equaliu(p,2))
1194 : {
1195 14 : if (get_ep(F)>0) return mkvec(F);
1196 : } else
1197 49 : if( ZX_pval(F,p) > 0) return mkvec(F);
1198 56 : W = next_model(F, g, p, 1, g+1); l = lg(W);
1199 56 : if (l==1) { set_avma(av); return mkvec(F); }
1200 0 : R = cgetg(3, t_VEC);
1201 0 : gel(R,1) = get_extremal(F, gel(W,1), g, p);
1202 0 : gel(R,2) = l==3 ? get_extremal(F, gel(W,2), g, p) : F;
1203 0 : if (gel(R,2) == gel(R,1)) setlg(R,2);
1204 0 : return gerepilecopy(av, R);
1205 : }
1206 :
1207 : static GEN
1208 301955 : RgX_RgM2_eval(GEN P, GEN A, GEN Bp, long d)
1209 : {
1210 301955 : if (signe(P)==0)
1211 78952 : return P;
1212 : else
1213 : {
1214 223003 : long dP = degpol(P);
1215 223003 : GEN R = RgX_homogenous_evalpow(P, A, Bp);
1216 223003 : if (d > dP)
1217 16948 : R = gmul(R, gel(Bp,1+d-dP));
1218 223003 : return R;
1219 : }
1220 : }
1221 :
1222 : static GEN
1223 52149 : minimalmodel_merge(GEN W2, GEN Modd, long g, long v)
1224 : {
1225 52149 : GEN P = gel(W2,1), Q = gel(W2,2);
1226 52149 : GEN e = gel(Modd,1), M = gel(Modd,2);
1227 52149 : GEN A = deg1pol_shallow(gcoeff(M,1,1), gcoeff(M,1,2), v);
1228 52149 : GEN B = deg1pol_shallow(gcoeff(M,2,1), gcoeff(M,2,2), v);
1229 52149 : GEN Bp = gpowers(B, 2*g+2);
1230 52149 : long f = mod4(e)==1 ? 1: -1;
1231 52149 : GEN m = shifti(f > 0 ? subui(1,e): addui(1,e), -2);
1232 52149 : GEN m24 = subii(shifti(m,1), shifti(sqri(m),2));
1233 52149 : P = RgX_RgM2_eval(P, A, Bp, 2*g+2);
1234 52149 : Q = RgX_RgM2_eval(Q, A, Bp, g+1);
1235 52149 : P = ZX_Z_divexact(ZX_add(P, ZX_Z_mul(ZX_sqr(Q), m24)),sqri(e));
1236 52149 : if (f < 0) Q = ZX_neg(Q);
1237 52149 : return mkvec2(P,Q);
1238 : }
1239 :
1240 : static GEN
1241 104284 : hyperell_redQ(GEN W)
1242 : {
1243 104284 : GEN P = gel(W,1), Q = gel(W,2);
1244 104284 : GEN Pr, Qr = FpX_red(Q, gen_2);
1245 104284 : Pr = ZX_add(P, ZX_shifti(ZX_mul(ZX_sub(Q, Qr),ZX_add(Q, Qr)),-2));
1246 104284 : return mkvec2(Pr, Qr);
1247 : }
1248 :
1249 : static GEN
1250 50735 : minimalmodel_getH(GEN W, GEN Qn, GEN e, GEN M, long g, long v)
1251 : {
1252 50735 : GEN Q = gel(W,2);
1253 50735 : GEN A = deg1pol_shallow(gcoeff(M,1,1), gcoeff(M,1,2), v);
1254 50735 : GEN B = deg1pol_shallow(gcoeff(M,2,1), gcoeff(M,2,2), v);
1255 50735 : GEN Bp = gpowers(B, g+1);
1256 50735 : return ZX_shifti(ZX_sub(ZX_Z_mul(Qn,e),RgX_RgM2_eval(Q, A, Bp, g+1)), -1);
1257 : }
1258 :
1259 : static void
1260 52191 : check_hyperell_Q(const char *fun, GEN *pW, GEN *pF)
1261 : {
1262 52191 : GEN W = *pW, F = check_hyperell(W);
1263 : long v, d;
1264 52191 : if (!F || !signe(F) || !RgX_is_ZX(F)) pari_err_TYPE(fun, W);
1265 52184 : if (!signe(ZX_disc(F))) pari_err_DOMAIN(fun,"disc(W)","==",gen_0,W);
1266 52177 : v = varn(F); d = degpol(F);
1267 52177 : if (typ(W)==t_POL) W = mkvec2(W, pol_0(v));
1268 : else
1269 : {
1270 43603 : GEN P = gel(W, 1), Q = gel(W, 2);
1271 43603 : long g = ((d+1) >> 1) - 1;
1272 43603 : if (typ(P)!=t_POL) P = scalarpol_shallow(P, v);
1273 43603 : if (typ(Q)!=t_POL) Q = scalarpol_shallow(Q, v);
1274 43603 : if (!RgX_is_ZX(P) || !RgX_is_ZX(Q)) pari_err_TYPE(fun,W);
1275 43603 : if (degpol(P) > 2*g+2) pari_err_DOMAIN(fun, "deg(P)", ">", utoi(2*g+2), P);
1276 43603 : if (degpol(Q) > g+1) pari_err_DOMAIN(fun, "deg(Q)", ">", utoi(g+1), Q);
1277 43603 : W = mkvec2(P, Q);
1278 : }
1279 52177 : if (d < 3) pari_err_DOMAIN(fun, "genus", "=", gen_0, gen_0);
1280 52163 : *pW = W; *pF = F;
1281 52163 : }
1282 :
1283 : GEN
1284 52149 : hyperellminimalmodel(GEN W, GEN *pM, GEN pr)
1285 : {
1286 52149 : pari_sp av = avma;
1287 : GEN Wr, F, WM2, F2, W2, M2, Modd, Wf, ef, Mf, Hf;
1288 : long d, g, v;
1289 52149 : check_hyperell_Q("hyperellminimalmodel",&W, &F);
1290 52149 : if (pr && (!is_vec_t(typ(pr)) || !RgV_is_ZV(pr)))
1291 14 : pari_err_TYPE("hyperellminimalmodel",pr);
1292 52135 : d = degpol(F); g = ((d+1)>>1)-1; v = varn(F);
1293 52135 : Wr = hyperell_redQ(W);
1294 52135 : if (!pr || RgV_isin(pr, gen_2))
1295 : {
1296 50784 : WM2 = algo56(Wr,g); W2 = gel(WM2, 1); M2 = gel(WM2, 2);
1297 50784 : F2 = check_hyperell(W2);
1298 : }
1299 : else
1300 : {
1301 1351 : W2 = Wr; F2 = F; M2 = mkvec2(gen_1, matid(2));
1302 : }
1303 52135 : Modd = gel(algo57(F2, g, pr), 2);
1304 52135 : Wf = hyperell_redQ(minimalmodel_merge(W2, Modd, g, v));
1305 52135 : if (!pM) return gerepilecopy(av, Wf);
1306 50721 : ef = mulii(gel(M2,1), gel(Modd,1));
1307 50721 : Mf = ZM2_mul(gel(M2,2), gel(Modd,2));
1308 50721 : Hf = minimalmodel_getH(W, gel(Wf,2), ef, Mf, g, v);
1309 50721 : *pM = mkvec3(ef, Mf, Hf);
1310 50721 : return gc_all(av, 2, &Wf, pM);
1311 : }
1312 :
1313 : GEN
1314 14 : hyperellminimaldisc(GEN W, GEN pr)
1315 : {
1316 14 : pari_sp av = avma;
1317 14 : GEN C = hyperellminimalmodel(W, NULL, pr);
1318 14 : return gerepileuptoint(av, hyperelldisc(C));
1319 : }
1320 :
1321 : static GEN
1322 35 : redqfbsplit(GEN a, GEN b, GEN c, GEN d)
1323 : {
1324 35 : GEN p = subii(d,b), q = shifti(a,1);
1325 35 : GEN U, Q, u, v, w = bezout(p, q, &u, &v);
1326 :
1327 35 : if (!equali1(w)) { p = diviiexact(p, w); q = diviiexact(q, w); }
1328 35 : U = mkmat22(p, negi(v), q, u);
1329 35 : Q = qfb_ZM_apply(mkvec3(a,b,c), U);
1330 35 : b = gel(Q, 2); c = gel(Q,3);
1331 35 : if (signe(b) < 0) gel(U,2) = mkcol2(v, negi(u));
1332 35 : gel(U,2) = ZC_lincomb(gen_1, truedivii(negi(c), d), gel(U,2), gel(U,1));
1333 35 : return U;
1334 : }
1335 :
1336 : static GEN
1337 16386 : polreduce(GEN P, GEN M)
1338 : {
1339 16386 : long v = varn(P), dP = degpol(P), d = odd(dP) ? dP+1: dP;
1340 16386 : GEN A = deg1pol_shallow(gcoeff(M,1,1), gcoeff(M,1,2), v);
1341 16386 : GEN B = deg1pol_shallow(gcoeff(M,2,1), gcoeff(M,2,2), v);
1342 16386 : return RgX_RgM2_eval(P, A, gpowers(B, d), d);
1343 : }
1344 :
1345 : /* assume deg(P) > 2 */
1346 : static GEN
1347 8193 : red_Cremona_Stoll(GEN P, GEN *pM)
1348 : {
1349 : GEN q1, q2, q3, M, R;
1350 8193 : long i, prec = nbits2prec(2*gexpo(P)) + EXTRAPRECWORD, d = degpol(P);
1351 8193 : GEN dP = ZX_deriv(P);
1352 : for (;;)
1353 0 : {
1354 8193 : GEN r = QX_complex_roots(P, prec);
1355 8193 : q1 = gen_0; q2 = gen_0; q3 = gen_0;
1356 41000 : for (i = 1; i <= d; i++)
1357 : {
1358 32807 : GEN ri = gel(r,i);
1359 32807 : GEN s = ginv(gabs(RgX_cxeval(dP,ri,NULL), prec));
1360 32807 : if (d!=4) s = gpow(s, gdivgs(gen_2,d-2), prec);
1361 32807 : q1 = gadd(q1, s);
1362 32807 : q2 = gsub(q2, gmul(real_i(ri), s));
1363 32807 : q3 = gadd(q3, gmul(gnorm(ri), s));
1364 : }
1365 8193 : M = lllgram(mkmat22(q1,q2,q2,q3));
1366 8193 : if (M && lg(M) == 3) break;
1367 0 : prec = precdbl(prec);
1368 : }
1369 8193 : R = polreduce(P, M);
1370 8193 : *pM = M;
1371 8193 : return R;
1372 : }
1373 :
1374 : /* assume deg(P) > 2 */
1375 : GEN
1376 8193 : ZX_hyperellred(GEN P, GEN *pM)
1377 : {
1378 8193 : pari_sp av = avma;
1379 8193 : long d = degpol(P);
1380 : GEN q1, q2, q3, D, vD;
1381 8193 : GEN a = gel(P,d+2), b = gel(P,d+1), c = gel(P, d);
1382 : GEN M, R, M2;
1383 :
1384 8193 : q1 = muliu(sqri(a), d);
1385 8193 : q2 = shifti(mulii(a,b), 1);
1386 8193 : q3 = subii(sqri(b), shifti(mulii(a,c), 1));
1387 8193 : D = gcdii(gcdii(q1, q2), q3);
1388 8193 : if (!equali1(D))
1389 : {
1390 8172 : q1 = diviiexact(q1, D);
1391 8172 : q2 = diviiexact(q2, D);
1392 8172 : q3 = diviiexact(q3, D);
1393 : }
1394 8193 : D = qfb_disc3(q1, q2, q3);
1395 8193 : if (!signe(D))
1396 49 : M = mkmat22(gen_1, truedivii(negi(q2),shifti(q1,1)), gen_0, gen_1);
1397 8144 : else if (issquareall(D,&vD))
1398 35 : M = redqfbsplit(q1, q2, q3, vD);
1399 : else
1400 8109 : M = gel(qfbredsl2(mkqfb(q1,q2,q3,D), NULL), 2);
1401 8193 : R = red_Cremona_Stoll(polreduce(P, M), &M2);
1402 8193 : if (pM) *pM = gmul(M, M2);
1403 8193 : return gc_all(av, pM ? 2: 1, &R, pM);
1404 : }
1405 :
1406 : GEN
1407 42 : hyperellred(GEN W, GEN *pM)
1408 : {
1409 42 : pari_sp av = avma;
1410 : long g, d, v;
1411 : GEN F, M, Wf, Hf;
1412 42 : check_hyperell_Q("hyperellred", &W, &F);
1413 14 : d = degpol(F); g = ((d+1)>>1)-1; v = varn(F);
1414 14 : (void) ZX_hyperellred(F, &M);
1415 14 : Wf = hyperell_redQ(minimalmodel_merge(W, mkvec2(gen_1, M), g, v));
1416 14 : Hf = minimalmodel_getH(W, gel(Wf,2), gen_1, M, g, v);
1417 14 : if (pM) *pM = mkvec3(gen_1, M, Hf);
1418 14 : return gc_all(av, pM ? 2: 1, &Wf, pM);
1419 : }
1420 :
1421 : static void
1422 65296 : check_hyperell_Rg(const char *fun, GEN *pW, GEN *pF)
1423 : {
1424 65296 : GEN W = *pW, F = check_hyperell(W);
1425 : long v;
1426 65296 : if (!F)
1427 7 : pari_err_TYPE(fun, W);
1428 65289 : if (degpol(F) <= 0) pari_err_CONSTPOL(fun);
1429 65282 : v = varn(F);
1430 65282 : if (typ(W)==t_POL) W = mkvec2(W, pol_0(v));
1431 : else
1432 : {
1433 65254 : GEN P = gel(W, 1), Q = gel(W, 2);
1434 65254 : long g = ((degpol(F)+1)>>1)-1;
1435 65254 : if( typ(P)!=t_POL) P = scalarpol(P, v);
1436 65254 : if( typ(Q)!=t_POL) Q = scalarpol(Q, v);
1437 65254 : if (degpol(P) > 2*g+2)
1438 0 : pari_err_DOMAIN(fun, "poldegree(P)", ">", utoi(2*g+2), P);
1439 65254 : if (degpol(Q) > g+1)
1440 0 : pari_err_DOMAIN(fun, "poldegree(Q)", ">", utoi(g+1), Q);
1441 :
1442 65254 : W = mkvec2(P, Q);
1443 : }
1444 65282 : if (pF) *pF = F;
1445 65282 : *pW = W;
1446 65282 : }
1447 :
1448 : static void
1449 65282 : check_hyperell_vc(const char *fun, GEN C, long v, GEN *e, GEN *M, GEN *H)
1450 : {
1451 65282 : if (typ(C) != t_VEC || lg(C) != 4) pari_err_TYPE(fun,C);
1452 65275 : *e = gel(C,1); *M = gel(C,2); *H = gel(C,3);
1453 65275 : if (typ(*M) != t_MAT || lg(*M) != 3 || lgcols(*M) != 3) pari_err_TYPE(fun,C);
1454 65268 : if (typ(*H)!=t_POL || varncmp(varn(*H),v) > 0) *H = scalarpol_shallow(*H,v);
1455 65268 : }
1456 :
1457 : GEN
1458 65296 : hyperellchangecurve(GEN W, GEN C)
1459 : {
1460 65296 : pari_sp av = avma;
1461 : GEN F, P, Q, A, B, Bp, e, M, H;
1462 : long d, g, v;
1463 65296 : check_hyperell_Rg("hyperellchangecurve",&W,&F);
1464 65282 : P = gel(W,1); Q = gel(W,2);
1465 65282 : d = degpol(F); g = ((d+1)>>1)-1; v = varn(F);
1466 65282 : check_hyperell_vc("hyperellchangecurve", C, v, &e, &M, &H);
1467 65268 : if (varncmp(gvar(M),v) <= 0)
1468 0 : pari_err_PRIORITY("hyperellchangecurve",M,"<=",v);
1469 65268 : A = deg1pol_shallow(gcoeff(M,1,1), gcoeff(M,1,2), v);
1470 65268 : B = deg1pol_shallow(gcoeff(M,2,1), gcoeff(M,2,2), v);
1471 65268 : Bp = gpowers(B, 2*g+2);
1472 65268 : P = RgX_RgM2_eval(P, A, Bp, 2*g+2);
1473 65268 : Q = RgX_RgM2_eval(Q, A, Bp, g+1);
1474 65268 : P = RgX_Rg_div(RgX_sub(P, RgX_mul(H,RgX_add(Q,H))), gsqr(e));
1475 65268 : Q = RgX_Rg_div(RgX_add(Q, RgX_mul2n(H,1)), e);
1476 65268 : return gerepilecopy(av, mkvec2(P,Q));
1477 : }
|