1 /* Copyright (C) 2000 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., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */
14 #include "pari.h"
15 #include "paripriv.h"
16
17 /********************************************************************/
18 /** **/
19 /** IGUSA INVARIANTS **/
20 /** (GP2C-generated) **/
21 /** **/
22 /********************************************************************/
23 /*
24 j2(a0,a1,a2,a3,a4,a5,a6) = (-120*a0*a6+20*a1*a5-8*a2*a4+3*a3^2) / 4;
25 */
26 static GEN
igusaj2(GEN a0,GEN a1,GEN a2,GEN a3,GEN a4,GEN a5,GEN a6)27 igusaj2(GEN a0, GEN a1, GEN a2, GEN a3, GEN a4, GEN a5, GEN a6)
28 {
29 pari_sp av = avma;
30 return gerepileupto(av, gmul2n(gadd(gsub(gadd(gmul(gmulsg(-120, a0), a6), gmul(gmulsg(20, a1), a5)), gmul(gmulsg(8, a2), a4)), gmulsg(3, gsqr(a3))), -2));
31 }
32
33 /*
34 j4(a0,a1,a2,a3,a4,a5,a6) = (240*(a0*a3*a4*a5+a1*a2*a3*a6)-400*(a0*a2*a5^2+a1^2*a4*a6)-64*(a0*a4^3+a2^3*a6)+16*(a1*a3*a4^2+a2^2*a3*a5)-672*a0*a3^2*a6+240*a1^2*a5^2-112*a1*a2*a4*a5-8*a1*a3^2*a5+16*a2^2*a4^2-16*a2*a3^2*a4+3*a3^4+2640*a0^2*a6^2-880*a0*a1*a5*a6+1312*a0*a2*a4*a6) / 2^7
35 */
36 static GEN
igusaj4(GEN a0,GEN a1,GEN a2,GEN a3,GEN a4,GEN a5,GEN a6)37 igusaj4(GEN a0, GEN a1, GEN a2, GEN a3, GEN a4, GEN a5, GEN a6)
38 {
39 pari_sp av = avma;
40 return gerepileupto(av,
41 gmul2n(gadd(gsub(gadd(gadd(gsub(gadd(gsub(gsub(gadd(gsub(gadd(gsub(gsub(gmulsg(240,
42 gadd(gmul(gmul(gmul(a0, a3), a4), a5), gmul(gmul(gmul(a1, a2), a3), a6))),
43 gmulsg(400, gadd(gmul(gmul(a0, a2), gsqr(a5)), gmul(gmul(gsqr(a1), a4), a6)))),
44 gmulsg(64, gadd(gmul(a0, gpowgs(a4, 3)), gmul(gpowgs(a2, 3), a6)))), gmulsg(16,
45 gadd(gmul(gmul(a1, a3), gsqr(a4)), gmul(gmul(gsqr(a2), a3), a5)))),
46 gmul(gmul(gmulsg(672, a0), gsqr(a3)), a6)), gmul(gmulsg(240, gsqr(a1)),
47 gsqr(a5))), gmul(gmul(gmul(gmulsg(112, a1), a2), a4), a5)), gmul(gmul(gmulsg(8,
48 a1), gsqr(a3)), a5)), gmul(gmulsg(16, gsqr(a2)), gsqr(a4))),
49 gmul(gmul(gmulsg(16, a2), gsqr(a3)), a4)), gmulsg(3, gpowgs(a3, 4))),
50 gmul(gmulsg(2640, gsqr(a0)), gsqr(a6))), gmul(gmul(gmul(gmulsg(880, a0), a1),
51 a5), a6)), gmul(gmul(gmul(gmulsg(1312, a0), a2), a4), a6)), -7));
52 }
53
54 /*
55 j6(a0,a1,a2,a3,a4,a5,a6) = (1600*(a0^2*a4^2*a5^2+a1^2*a2^2*a6^2)+1600*(a0*a1*a2*a5^3+a1^3*a4*a5*a6)+640*(a0*a1*a3*a4*a5^2+a1^2*a2*a3*a5*a6)-4000*(a0^2*a3*a5^3+a1^3*a3*a6^2)-384*(a0*a1*a4^3*a5+a1*a2^3*a5*a6)-640*(a0*a2^2*a4*a5^2+a1^2*a2*a4^2*a6)+80*(a0*a2*a3^2*a5^2+a1^2*a3^2*a4*a6)+192*(a0*a2*a3*a4^2*a5+a1*a2^2*a3*a4*a6)-48*(a0*a3^3*a4*a5+a1*a2*a3^3*a6)-224*(a1^2*a3*a4^2*a5+a1*a2^2*a3*a5^2)+64*(a1^2*a4^4+a2^4*a5^2)-64*(a1*a2*a3*a4^3+a2^3*a3*a4*a5)+16*(a1*a3^3*a4^2+a2^2*a3^3*a5)-4096*(a0^2*a4^3*a6+a0*a2^3*a6^2)+6400*(a0^2*a2*a5^2*a6+a0*a1^2*a4*a6^2)+10560*(a0^2*a3*a4*a5*a6+a0*a1*a2*a3*a6^2)+2624*(a0*a1*a3*a4^2*a6+a0*a2^2*a3*a5*a6)-4432*a0*a1*a3^2*a5*a6-8*a2*a3^4*a4+a3^6-320*a1^3*a5^3+64*a1^2*a2*a4*a5^2+176*a1^2*a3^2*a5^2+128*a1*a2^2*a4^2*a5+112*a1*a2*a3^2*a4*a5-28*a1*a3^4*a5+16*a2^2*a3^2*a4^2+5120*a0^3*a6^3-2544*a0^2*a3^2*a6^2+312*a0*a3^4*a6-14336*a0^2*a2*a4*a6^2+1024*a0*a2^2*a4^2*a6-2560*a0^2*a1*a5*a6^2-2240*a0*a1^2*a5^2*a6-6528*a0*a1*a2*a4*a5*a6-1568*a0*a2*a3^2*a4*a6) / 2^10
56 */
57 static GEN
igusaj6(GEN a0,GEN a1,GEN a2,GEN a3,GEN a4,GEN a5,GEN a6)58 igusaj6(GEN a0, GEN a1, GEN a2, GEN a3, GEN a4, GEN a5, GEN a6)
59 {
60 pari_sp av = avma;
61 return gerepileupto(av,
62 gmul2n(gsub(gsub(gsub(gsub(gadd(gsub(gadd(gsub(gadd(gadd(gsub(gadd(gadd(gadd(gadd(gsub(gadd(gsub(gsub(gadd(gadd(gadd(gsub(gadd(gsub(gadd(gsub(gsub(gadd(gadd(gsub(gsub(gsub(gadd(gadd(gmulsg(1600,
63 gadd(gmul(gmul(gsqr(a0), gsqr(a4)), gsqr(a5)), gmul(gmul(gsqr(a1), gsqr(a2)),
64 gsqr(a6)))), gmulsg(1600, gadd(gmul(gmul(gmul(a0, a1), a2), gpowgs(a5, 3)),
65 gmul(gmul(gmul(gpowgs(a1, 3), a4), a5), a6)))), gmulsg(640,
66 gadd(gmul(gmul(gmul(gmul(a0, a1), a3), a4), gsqr(a5)),
67 gmul(gmul(gmul(gmul(gsqr(a1), a2), a3), a5), a6)))), gmulsg(4000,
68 gadd(gmul(gmul(gsqr(a0), a3), gpowgs(a5, 3)), gmul(gmul(gpowgs(a1, 3), a3),
69 gsqr(a6))))), gmulsg(384, gadd(gmul(gmul(gmul(a0, a1), gpowgs(a4, 3)), a5),
70 gmul(gmul(gmul(a1, gpowgs(a2, 3)), a5), a6)))), gmulsg(640,
71 gadd(gmul(gmul(gmul(a0, gsqr(a2)), a4), gsqr(a5)), gmul(gmul(gmul(gsqr(a1),
72 a2), gsqr(a4)), a6)))), gmulsg(80, gadd(gmul(gmul(gmul(a0, a2), gsqr(a3)),
73 gsqr(a5)), gmul(gmul(gmul(gsqr(a1), gsqr(a3)), a4), a6)))), gmulsg(192,
74 gadd(gmul(gmul(gmul(gmul(a0, a2), a3), gsqr(a4)), a5), gmul(gmul(gmul(gmul(a1,
75 gsqr(a2)), a3), a4), a6)))), gmulsg(48, gadd(gmul(gmul(gmul(a0, gpowgs(a3, 3)),
76 a4), a5), gmul(gmul(gmul(a1, a2), gpowgs(a3, 3)), a6)))), gmulsg(224,
77 gadd(gmul(gmul(gmul(gsqr(a1), a3), gsqr(a4)), a5), gmul(gmul(gmul(a1,
78 gsqr(a2)), a3), gsqr(a5))))), gmulsg(64, gadd(gmul(gsqr(a1), gpowgs(a4, 4)),
79 gmul(gpowgs(a2, 4), gsqr(a5))))), gmulsg(64, gadd(gmul(gmul(gmul(a1, a2), a3),
80 gpowgs(a4, 3)), gmul(gmul(gmul(gpowgs(a2, 3), a3), a4), a5)))), gmulsg(16,
81 gadd(gmul(gmul(a1, gpowgs(a3, 3)), gsqr(a4)), gmul(gmul(gsqr(a2), gpowgs(a3,
82 3)), a5)))), gmulsg(4096, gadd(gmul(gmul(gsqr(a0), gpowgs(a4, 3)), a6),
83 gmul(gmul(a0, gpowgs(a2, 3)), gsqr(a6))))), gmulsg(6400,
84 gadd(gmul(gmul(gmul(gsqr(a0), a2), gsqr(a5)), a6), gmul(gmul(gmul(a0,
85 gsqr(a1)), a4), gsqr(a6))))), gmulsg(10560, gadd(gmul(gmul(gmul(gmul(gsqr(a0),
86 a3), a4), a5), a6), gmul(gmul(gmul(gmul(a0, a1), a2), a3), gsqr(a6))))),
87 gmulsg(2624, gadd(gmul(gmul(gmul(gmul(a0, a1), a3), gsqr(a4)), a6),
88 gmul(gmul(gmul(gmul(a0, gsqr(a2)), a3), a5), a6)))),
89 gmul(gmul(gmul(gmul(gmulsg(4432, a0), a1), gsqr(a3)), a5), a6)),
90 gmul(gmul(gmulsg(8, a2), gpowgs(a3, 4)), a4)), gpowgs(a3, 6)), gmul(gmulsg(320,
91 gpowgs(a1, 3)), gpowgs(a5, 3))), gmul(gmul(gmul(gmulsg(64, gsqr(a1)), a2), a4),
92 gsqr(a5))), gmul(gmul(gmulsg(176, gsqr(a1)), gsqr(a3)), gsqr(a5))),
93 gmul(gmul(gmul(gmulsg(128, a1), gsqr(a2)), gsqr(a4)), a5)),
94 gmul(gmul(gmul(gmul(gmulsg(112, a1), a2), gsqr(a3)), a4), a5)),
95 gmul(gmul(gmulsg(28, a1), gpowgs(a3, 4)), a5)), gmul(gmul(gmulsg(16, gsqr(a2)),
96 gsqr(a3)), gsqr(a4))), gmul(gmulsg(5120, gpowgs(a0, 3)), gpowgs(a6, 3))),
97 gmul(gmul(gmulsg(2544, gsqr(a0)), gsqr(a3)), gsqr(a6))), gmul(gmul(gmulsg(312,
98 a0), gpowgs(a3, 4)), a6)), gmul(gmul(gmul(gmulsg(14336, gsqr(a0)), a2), a4),
99 gsqr(a6))), gmul(gmul(gmul(gmulsg(1024, a0), gsqr(a2)), gsqr(a4)), a6)),
100 gmul(gmul(gmul(gmulsg(2560, gsqr(a0)), a1), a5), gsqr(a6))),
101 gmul(gmul(gmul(gmulsg(2240, a0), gsqr(a1)), gsqr(a5)), a6)),
102 gmul(gmul(gmul(gmul(gmul(gmulsg(6528, a0), a1), a2), a4), a5), a6)),
103 gmul(gmul(gmul(gmul(gmulsg(1568, a0), a2), gsqr(a3)), a4), a6)), -10));
104 }
105
106 /********************************************************************/
107 /** **/
108 /** A REDUCTION ALGORITHM "A LA TATE" FOR CURVES OF GENUS 2 **/
109 /** **/
110 /********************************************************************/
111 /* Based on genus2reduction-0.3, http://www.math.u-bordeaux.fr/~liu/G2R/
112 * by Qing Liu <liu@math.u-bordeaux.fr>
113 * and Henri Cohen <cohen@math.u-bordeaux.fr>
114
115 * Qing Liu: Modeles minimaux des courbes de genre deux
116 * J. fuer die Reine und Angew. Math., 453 (1994), 137-164.
117 * http://www.math.u-bordeaux.fr/~liu/articles/modregE.ps */
118
119 /* some auxiliary polynomials, gp2c-generated */
120
121 /*
122 apol2(a0,a1,a2) = -5*a1^2+12*a0*a2;
123 */
124 static GEN
apol2(GEN a0,GEN a1,GEN a2)125 apol2(GEN a0, GEN a1, GEN a2)
126 {
127 return gadd(gmulsg(-5, gsqr(a1)), gmul(gmulsg(12, a0), a2));
128 }
129
130 /*
131 apol3(a0,a1,a2,a3) = 5*a1^3+9*a0*(-2*a1*a2+3*a0*a3);
132 */
133 static GEN
apol3(GEN a0,GEN a1,GEN a2,GEN a3)134 apol3(GEN a0, GEN a1, GEN a2, GEN a3)
135 {
136 return gadd(gmulsg(5, gpowgs(a1, 3)), gmul(gmulsg(9, a0), gadd(gmul(gmulsg(-2, a1), a2), gmul(gmulsg(3, a0), a3))));
137 }
138
139 /*
140 apol5(a0,a1,a2,a3,a4,a5) = a1^5+3*a0*(-2*a1^3*a2+9*a0*a1^2*a3-36*a0^2*a1*a4+108*a0^3*a5);
141 */
142 static GEN
apol5(GEN a0,GEN a1,GEN a2,GEN a3,GEN a4,GEN a5)143 apol5(GEN a0, GEN a1, GEN a2, GEN a3, GEN a4, GEN a5)
144 {
145 return gadd(gpowgs(a1, 5), gmul(gmulsg(3, a0), gadd(gsub(gadd(gmul(gmulsg(-2, gpowgs(a1, 3)), a2), gmul(gmul(gmulsg(9, a0), gsqr(a1)), a3)), gmul(gmul(gmulsg(36, gsqr(a0)), a1), a4)), gmul(gmulsg(108, gpowgs(a0, 3)), a5))));
146 }
147
148 /*
149 bpol2(a0,a1,a2,a3,a4) = 2*a2^2-5*a1*a3+10*a0*a4;
150 */
151 static GEN
bpol2(GEN a0,GEN a1,GEN a2,GEN a3,GEN a4)152 bpol2(GEN a0, GEN a1, GEN a2, GEN a3, GEN a4)
153 {
154 return gadd(gsub(gmulsg(2, gsqr(a2)), gmul(gmulsg(5, a1), a3)), gmul(gmulsg(10, a0), a4));
155 }
156
157 static const long VERYBIG = (1L<<20);
158 static long
myval(GEN x,GEN p)159 myval(GEN x, GEN p) { return signe(x)? Z_pval(x,p): VERYBIG; }
160 static long
my3val(GEN x)161 my3val(GEN x) { return signe(x)? Z_lval(x,3): VERYBIG; }
162 /* b in Z[i], return v_3(b) */
163 static long
myval_zi(GEN b)164 myval_zi(GEN b) { return minss(my3val(real_i(b)), my3val(imag_i(b))); }
165 /* b in Z[i, Y]/(Y^2-3), return v_Y(b) */
166 static long
myval_zi2(GEN b)167 myval_zi2(GEN b)
168 {
169 long v0, v1;
170 b = lift_shallow(b);
171 v0 = myval_zi(RgX_coeff(b,0));
172 v1 = myval_zi(RgX_coeff(b,1));
173 return minss(2*v0, 2*v1+1);
174 }
175
176 /* min(a,b,c) */
177 static long
min3(long a,long b,long c)178 min3(long a, long b, long c)
179 {
180 long m = a;
181 if (b < m) m = b;
182 if (c < m) m = c;
183 return m;
184 }
185
186 /* Vector of p-adic factors (over Q_p) to accuracy r of pol. */
187 static GEN
padicfactors(GEN pol,GEN p,long r)188 padicfactors(GEN pol, GEN p, long r) { return gel(factorpadic(pol,p,r),1); }
189
190 /* x(1/t)*t^6, deg x <= 6 */
191 static GEN
RgX_recip6(GEN x)192 RgX_recip6(GEN x)
193 {
194 long lx = lg(x), i, j;
195 GEN y = cgetg(9, t_POL);
196 y[1] = x[1];
197 for (i=8,j=2; j < lx; i--,j++) gel(y,i) = gel(x,j);
198 for ( ; j < 9; i--,j++) gel(y,i) = gen_0;
199 return normalizepol_lg(y, 9);
200 }
201 /* extract coefficients of a polynomial a0 X^6 + ... + a6, of degree <= 6 */
202 static void
RgX_to_06(GEN q,GEN * a0,GEN * a1,GEN * a2,GEN * a3,GEN * a4,GEN * a5,GEN * a6)203 RgX_to_06(GEN q, GEN *a0, GEN *a1, GEN *a2, GEN *a3, GEN *a4, GEN *a5, GEN *a6)
204 {
205 *a0 = gen_0;
206 *a1 = gen_0;
207 *a2 = gen_0;
208 *a3 = gen_0;
209 *a4 = gen_0;
210 *a5 = gen_0;
211 *a6 = gen_0;
212 switch(degpol(q))
213 {
214 case 6: *a0 = gel(q,8); /*fall through*/
215 case 5: *a1 = gel(q,7); /*fall through*/
216 case 4: *a2 = gel(q,6); /*fall through*/
217 case 3: *a3 = gel(q,5); /*fall through*/
218 case 2: *a4 = gel(q,4); /*fall through*/
219 case 1: *a5 = gel(q,3); /*fall through*/
220 case 0: *a6 = gel(q,2); /*fall through*/
221 }
222 }
223 /* extract coefficients a0,...a3 of a polynomial a0 X^6 + ... + a6 */
224 static void
RgX_to_03(GEN q,GEN * a0,GEN * a1,GEN * a2,GEN * a3)225 RgX_to_03(GEN q, GEN *a0, GEN *a1, GEN *a2, GEN *a3)
226 {
227 *a0 = gen_0;
228 *a1 = gen_0;
229 *a2 = gen_0;
230 *a3 = gen_0;
231 switch(degpol(q))
232 {
233 case 6: *a0 = gel(q,8); /*fall through*/
234 case 5: *a1 = gel(q,7); /*fall through*/
235 case 4: *a2 = gel(q,6); /*fall through*/
236 case 3: *a3 = gel(q,5); /*fall through*/
237 }
238 }
239
240 /* deg(H mod p) = 3, return v_p( disc(correspondig p-adic factor) ) */
241 static long
discpart(GEN H,GEN p,long prec)242 discpart(GEN H, GEN p, long prec)
243 {
244 GEN list, prod, dis;
245 long i, j;
246
247 if (degpol(FpX_red(H,p)) != 3)
248 pari_err_BUG("discpart [must not reach]"); /* LCOV_EXCL_LINE */
249 list = padicfactors(H,p,prec);
250 prod = pol_1(varn(H));
251 for(i = 1; i < lg(list); i++)
252 {
253 GEN t = gel(list,i);
254 for(j = 3; j < lg(t); j++) /* include if nonconstant mod p */
255 if (!valp(gel(t,j))) { prod = RgX_mul(prod,t); break; }
256 }
257 if (degpol(prod) != 3) pari_err_BUG("discpart [prod degree]");
258 dis = RgX_disc(prod);
259 return gequal0(dis)? prec+1: valp(dis);
260 }
261
262 /* B = b0 X^6 + ... + b6 a ZX, 0 <= j <= 3.
263 * Let theta_j(H) := min { v_p(b_i) / (i - j), j < i <= 6 } >= 0.
264 * Return 60 theta \in Z */
265 static long
theta_j(GEN B,GEN p,long j)266 theta_j(GEN B, GEN p, long j)
267 {
268 long i, t = 60*myval(RgX_coeff(B,5-j), p);
269 for(i = 2+j; i <= 6; i++)
270 t = minss(t, myval(RgX_coeff(B,6-i), p) * (60 / (i-j)));
271 return t;
272 }
273 /* compute 6 * theta_3 for B in Z[i][X], p = 3 */
274 static long
theta_3_zi(GEN B)275 theta_3_zi(GEN B)
276 {
277 long v2 = myval_zi(RgX_coeff(B,2));
278 long v1 = myval_zi(RgX_coeff(B,1));
279 long v0 = myval_zi(RgX_coeff(B,0));
280 return min3(6*v2, 3*v1, 2*v0);
281 }
282 /* compute 6 * theta_3 for B in (Z[i,Y]/(Y^2-3))[X], p = 3 */
283 static long
theta_3_zi2(GEN B)284 theta_3_zi2(GEN B)
285 {
286 long v2 = myval_zi2(RgX_coeff(B,2));
287 long v1 = myval_zi2(RgX_coeff(B,1));
288 long v0 = myval_zi2(RgX_coeff(B,0));
289 return min3(6*v2, 3*v1, 2*v0);
290 }
291
292 /* Set maxord to the maximal multiplicity of a factor. If there is at least
293 * a triple root (=> maxord >= 3) return it, else return NULL */
294 static GEN
factmz(GEN Q,GEN p,long * maxord)295 factmz(GEN Q, GEN p, long *maxord)
296 {
297 GEN z = FpX_factor_squarefree(Q, p);
298 long m = lg(z)-1; /* maximal multiplicity */
299 *maxord = m;
300 return (m >= 3)? FpX_oneroot(gel(z,m), p): NULL;
301 }
302
303 /* H integral ZX of degree 5 or 6, p > 2. Modify until
304 * y^2 = p^alpha H is minimal over Z_p, alpha = 0,1
305 * Return [H,lambda,60*theta,alpha,quad,beta], where
306 * - quad = 1 if H has a root of order 3 in F_p^2 \ F_p, 0 otherwise
307 * - 0 <= lambda <= 3, index of a coefficient with valuation 0
308 * - theta = theta_j(H(x + r), p, lambda), 60*theta in Z, where r is a root
309 * of H mod p
310 * - beta >= -1 s.t. H = p^n H0(r + p^beta * X) for some n, r in Z, where
311 * H0 is the initial H or polrecip(H) */
312 static GEN
polymini(GEN H,GEN p)313 polymini(GEN H, GEN p)
314 {
315 GEN a0, a1, a2, a3, Hp, rac;
316 long t60, alpha, lambda, maxord, quad = 0, beta = 0;
317
318 alpha = ZX_pvalrem(H, p, &H) & 1;
319 RgX_to_03(H, &a0,&a1,&a2,&a3);
320 if (dvdii(a0,p) && dvdii(a1,p) && dvdii(a2,p) && dvdii(a3,p))
321 {
322 H = RgX_recip6(H);
323 RgX_to_03(H, &a0,&a1,&a2,&a3);
324 }
325 if (!dvdii(a3,p)) lambda = 3;
326 else if (!dvdii(a2,p)) lambda = 2;
327 else if (!dvdii(a1,p)) lambda = 1;
328 else lambda = 0;
329
330 for(;;)
331 { /* lambda <= 3, t60 = 60*theta */
332 long e;
333 t60 = theta_j(H,p,lambda); e = t60 / 60;
334 if (e)
335 {
336 GEN pe = powiu(p,e);
337 /* H <- H(p^e X) / p^(e(6-lambda)) */
338 H = ZX_Z_divexact(ZX_unscale_div(H,pe), powiu(pe,5-lambda));
339 alpha = (alpha + lambda*e)&1;
340 beta += e;
341 t60 -= 60*e;
342 }
343 /* 0 <= t < 60 */
344 Hp = FpX_red(H, p);
345 if (t60) break;
346
347 rac = factmz(Hp,p, &maxord);
348 if (maxord <= 2)
349 {
350 if (degpol(Hp) <= 3) break;
351 goto end;
352 }
353 else
354 { /* maxord >= 3 */
355 if (!rac) { quad = 1; goto end; }
356 if (signe(rac)) H = ZX_translate(H, rac);
357 lambda = 6-maxord;
358 }
359 }
360
361 if (lambda <= 2)
362 {
363 if (myval(RgX_coeff(H,2),p) > 1-alpha &&
364 myval(RgX_coeff(H,1),p) > 2-alpha &&
365 myval(RgX_coeff(H,0),p) > 3-alpha)
366 {
367 H = ZX_unscale(H, p);
368 if (alpha) H = ZX_Z_mul(H, p);
369 return polymini(H, p);
370 }
371 }
372 else if (lambda == 3 && alpha == 1)
373 {
374 if (degpol(Hp) == 3)
375 {
376 if (myval(RgX_coeff(H,6),p) >= 3 &&
377 myval(RgX_coeff(H,5),p) >= 2)
378 { /* too close to root [Kodaira symbol for y^2 = p^alpha*H not
379 implemented when alpha = 1]: go back one step */
380 H = ZX_rescale(H, p); /* H(x/p)p^(deg H) */
381 H = ZX_Z_divexact(H, powiu(p, degpol(H)-3)); /* H(x/p)p^3 */
382 t60 += 60; alpha = 0; beta--;
383 }
384 }
385 else if (degpol(Hp) == 6 && t60)
386 {
387 rac = factmz(RgX_mulXn(Hp, -3), p, &maxord);
388 if (maxord == 3)
389 {
390 GEN T = ZX_unscale(ZX_translate(H,rac),p); /* H(rac + px) */
391 if (ZX_pval(T,p)>= 3)
392 {
393 H = RgX_Rg_div(T, powiu(p,3));
394 alpha = 0; beta++;
395 t60 = theta_j(H,p,3);
396 if (!t60)
397 {
398 Hp = FpX_red(H, p);
399 if (!signe(FpX_disc(Hp,p)))
400 {
401 rac = FpX_oneroot(Hp, p);
402 t60 = theta_j(ZX_translate(H,rac),p,3);
403 }
404 }
405 }
406 }
407 }
408 }
409 end:
410 return mkvec2(H, mkvecsmall5(lambda,t60,alpha,quad,beta));
411 }
412
413 /* a in Q[i], return a^3 mod 3 */
414 static GEN
zi_pow3mod(GEN a)415 zi_pow3mod(GEN a)
416 {
417 GEN x, y;
418 if (typ(a) != t_COMPLEX) return gmodgs(a,3);
419 x = gmodgs(gel(a,1), 3);
420 y = gmodgs(gel(a,2), 3);
421 return mkcomplex(x, negi(y));
422 }
423 static GEN
polymini_zi(GEN pol)424 polymini_zi(GEN pol) /* polynome minimal dans Z[i] */
425 {
426 GEN polh, rac, a0, a1, a2, a3, a4, a5, a6, p = utoipos(3);
427 long alpha, beta = 0, t6;
428
429 alpha = ZX_pval(pol,p) & 1;
430 polh = alpha? RgX_Rg_div(pol, p): pol;
431 rac = mkcomplex(Fp_div(RgX_coeff(polh,3), RgX_coeff(polh,6), p), gen_1);
432 for(;;)
433 {
434 long e;
435 polh = RgX_translate(polh, rac);
436 t6 = theta_3_zi(polh); e = t6 / 6;
437 if (e)
438 {
439 GEN pe = powiu(p,e);
440 polh = RgX_Rg_div(RgX_unscale(polh,pe), powiu(pe,3));
441 alpha = (alpha+e)&1;
442 t6 -= e * 6; beta += e;
443 }
444 RgX_to_06(polh, &a0,&a1,&a2,&a3,&a4,&a5,&a6);
445 if (t6 || !myval_zi(a4) || !myval_zi(a5)) break;
446 rac = zi_pow3mod(gdiv(a6, gneg(a3)));
447 }
448 if (alpha && myval_zi(a0) >= 3 && myval_zi(a1) >= 2 && myval_zi(a2) >= 1)
449 {
450 t6 += 6; beta--; alpha = 0;
451 }
452 if (alpha && beta >= 1) pari_err_BUG("quadratic");
453 return mkvecsmall3(t6, alpha, beta);
454 }
455
456 /* pol is a ZX, minimal polynomial over Z_3[i,Y]/(Y^2-3) */
457 static GEN
polymini_zi2(GEN pol)458 polymini_zi2(GEN pol)
459 {
460 long alpha, beta, t6;
461 GEN a0, a1, a2, a3, a4, a5, a6;
462 GEN polh, rac, y = pol_x(fetch_var()), p = utoipos(3);
463
464 if (ZX_pval(pol,p)) pari_err_BUG("polymini_zi2 [polynomial not minimal]");
465 y = mkpolmod(y, gsubgs(gsqr(y), 3)); /* mod(y,y^2-3) */
466 polh = gdivgs(RgX_unscale(pol, y),27); /* H(y*x) / 27 */
467 if (myval_zi2(RgX_coeff(polh,4)) <= 0 ||
468 myval_zi2(RgX_coeff(polh,2)) <= 0)
469 {
470 (void)delete_var();
471 return mkvecsmall2(0,0);
472 }
473
474 if (myval_zi2(gsub(RgX_coeff(polh,6), RgX_coeff(polh,0))) > 0)
475 rac = gen_I();
476 else
477 rac = gen_1;
478 alpha = 0;
479 beta = 0;
480 for(;;)
481 {
482 long e;
483 polh = RgX_translate(polh, rac);
484 t6 = theta_3_zi2(polh); e = t6 / 6;
485 if (e)
486 {
487 GEN pent = gpowgs(y, e);
488 polh = RgX_Rg_div(RgX_unscale(polh, pent), gpowgs(pent,3));
489 alpha = (alpha+e)&1;
490 t6 -= 6*e; beta += e;
491 }
492 RgX_to_06(polh, &a0,&a1,&a2,&a3,&a4,&a5,&a6);
493 if (t6 || !myval_zi2(a4) || !myval_zi2(a5)) break;
494 a3 = liftpol_shallow(a3); if (typ(a3)==t_POL) a3 = RgX_coeff(a3,0);
495 a6 = liftpol_shallow(a6); if (typ(a6)==t_POL) a6 = RgX_coeff(a6,0);
496 rac = zi_pow3mod(gdiv(a6,gneg(a3)));
497 }
498 if (alpha)
499 {
500 if (myval_zi2(a0) < 3 || myval_zi2(a1) < 2 || myval_zi2(a2) < 1)
501 pari_err_BUG("polymini_zi2 [alpha]");
502 t6 += 6; beta--;
503 }
504 (void)delete_var();
505 if (odd(beta)) pari_err_BUG("quartic [type over Z[i] must be [K-K-(2*m)]]");
506 return mkvecsmall2(t6, beta);
507 }
508
509 struct igusa {
510 GEN j2, i4, j4, j6, j8, j10, i12;
511 GEN a0, A2, A3, A5, B2;
512 };
513 struct igusa_p {
514 long eps, tt, r1, r2, tame;
515 GEN p, stable, val, neron;
516 const char *type;
517 };
518
519 /* initialize Ip */
520 static void
stable_reduction(struct igusa * I,struct igusa_p * Ip,GEN p)521 stable_reduction(struct igusa *I, struct igusa_p *Ip, GEN p)
522 {
523 static const long d[9] = { 0,60,30,30,20,15,12,10 }; /* 120 / deg(X) */
524 GEN j2 = I->j2, i4 = I->i4, j4 = I->j4, j6 = I->j6, j8 = I->j8;
525 GEN val, J, v, Ieps, j10 = I->j10, i12 = I->i12;
526 long s, r1, r2, r3, r4, i, eps;
527
528 Ip->tame = 0;
529 Ip->neron = NULL;
530 Ip->type = NULL;
531 Ip->p = p;
532 Ip->val = val = cgetg(9, t_VECSMALL);
533 val[1] = myval(j2,p);
534 val[2] = myval(j4,p);
535 val[3] = myval(i4,p);
536 val[4] = myval(j6,p);
537 val[5] = myval(j8,p);
538 val[6] = myval(j10,p);
539 val[7] = myval(i12,p);
540 switch(itos_or_0(p))
541 {
542 case 2: eps = 4; val[8] = val[5]; Ieps = j8; break;
543 case 3: eps = 3; val[8] = val[4]; Ieps = j6; break;
544 default: eps = 1; val[8] = val[1]; Ieps = gdivgs(j2,12); break;
545 }
546
547 v = cgetg(8,t_VECSMALL);
548 for(i = 1; i <= 7; i++) v[i] = val[i] * d[i];
549 s = vecsmall_min(v);
550 Ip->eps = eps;
551
552 r1 = 3*eps*val[3];
553 r3 = eps*val[6] + val[8];
554 r2 = eps*val[7];
555 r4 = min3(r1, r2, r3);
556
557 /* s = max(v_p(X) / deg(X)) */
558 J = cgetg(1, t_VEC);
559 if (s == v[6])
560 Ip->tt = 1;
561 else if (s == v[7])
562 {
563 J = mkvec( Fp_to_mod(gmod(gdiv(gpowgs(i4,3),i12), p), p) );
564 Ip->tt = 2;
565 }
566 else if (s == v[3])
567 Ip->tt = (val[2] == val[3] || 2*val[4] == 3*val[3])? 3: 4;
568 else if (r3 == r4)
569 {
570 GEN a,b, P, sj, pj, t = gmul(gpowgs(j10,eps),Ieps);
571 sj = gaddsg(1728, gdiv(gpowgs(i12,eps), t));
572 pj = gdiv(gpowgs(i4,3*eps), t);
573 a = gmod(sj, p);
574 b = gmod(pj, p);
575 P = mkpoln(3, gen_1, Fp_neg(a,p), b, 0); /* X^2 - SX + P: roots j1,j2 */
576 J = FpX_roots(P, p);
577 switch(lg(J)-1)
578 {
579 case 0:
580 P = FpX_to_mod(P, p);
581 a = FpX_to_mod(pol_x(0), p);
582 b = FpX_to_mod(deg1pol_shallow(b, gen_m1,0), p);
583 J = mkvec2(mkpolmod(a,P), mkpolmod(b,P)); break;
584 case 1:
585 a = Fp_to_mod(gel(J,1), p);
586 J = mkvec2(a, a); break;
587 case 2:
588 settyp(J, t_VEC);
589 J = FpV_to_mod(J, p); break;
590 }
591 Ip->tt = 5;
592 }
593 else if (r2 == r4)
594 {
595 J = mkvec( Fp_to_mod(gmod(gdiv(gpowgs(i4,3),i12), p), p) );
596 Ip->tt = 6;
597 }
598 else
599 Ip->tt = 7; /* r1 == r4 */
600 Ip->stable = mkvec2(stoi(Ip->tt), J);
601 }
602
603 struct red {
604 const char *t, *pages;
605 double tnum;
606 GEN g;
607 };
608
609 /* destroy v */
610 static GEN
zv_snf(GEN v)611 zv_snf(GEN v)
612 {
613 long i, l = lg(v);
614 for (i = 1; i < l; i++)
615 {
616 long j, a = v[i];
617 for (j = i+1; j < l; j++)
618 {
619 long b = v[j], d = ugcd(a,b);
620 v[i] = a = a*(b/d);
621 v[j] = d;
622 }
623 }
624 for (i = l-1; i > 0; i--)
625 if (v[i] != 1) { setlg(v, i+1); break; }
626 return zv_to_ZV(v);
627 }
628
629 static GEN
cyclic(long n)630 cyclic(long n)
631 { return (n <= 1)? cgetg(1, t_VECSMALL): mkvecsmall(n); }
632 static GEN
dicyclic(long a,long b)633 dicyclic(long a, long b)
634 {
635 long d;
636 if (!a) a = 1;
637 if (!b) b = 1;
638 if (a < b) lswap(a,b);
639 d = ugcd(a,b);
640 if (d == 1) return cyclic(a*b);
641 return mkvecsmall2(a*b/d, d);
642 }
643 /* Z/2xZ/2, resp Z/4 for n even, resp. odd */
644 static GEN
groupH(long n)645 groupH(long n) { return odd(n)? cyclic(4): dicyclic(2,2); }
646
647 static long
get_red(struct red * S,struct igusa_p * Ip,GEN polh,GEN p,long alpha,long r)648 get_red(struct red *S, struct igusa_p *Ip, GEN polh, GEN p, long alpha, long r)
649 {
650 GEN val = Ip->val;
651 long indice;
652 switch(r)
653 {
654 case 0:
655 indice = FpX_is_squarefree(FpX_red(polh,p), p)
656 ? 0
657 : val[6] - val[7] + val[8]/Ip->eps;
658 S->t = stack_sprintf("I{%ld}", indice);
659 S->tnum = 1;
660 S->pages = "159-177";
661 S->g = cyclic(indice);
662 return indice ? indice: 1;
663 case 6:
664 if (alpha == 0) /* H(px) /p^3 */
665 polh = ZX_Z_divexact(ZX_unscale_div(polh,p), sqri(p));
666 indice = FpX_is_squarefree(FpX_red(polh,p), p)
667 ? 0
668 : val[6] - val[7] + val[8]/Ip->eps;
669 S->t = stack_sprintf("I*{%ld}", indice);
670 S->tnum = 1.5;
671 S->pages = "159-177";
672 S->g = groupH(indice);
673 return indice + 5;
674 case 3:
675 S->t = "III";
676 S->tnum = 3;
677 S->pages = "161-177";
678 S->g = cyclic(2);
679 return 2;
680 case 9:
681 S->t = "III*";
682 S->tnum = 3.5;
683 S->pages = "162-177";
684 S->g = cyclic(2);
685 return 8;
686 case 2:
687 S->t = "II";
688 S->tnum = 2;
689 S->pages = "159-174";
690 S->g = cyclic(1);
691 return 1;
692 case 8:
693 S->t = "IV*";
694 S->tnum = 4.5;
695 S->pages = "160-175";
696 S->g = cyclic(3);
697 return 7;
698 case 4:
699 S->t = "IV";
700 S->tnum = 4;
701 S->pages = "160-174";
702 S->g = cyclic(3);
703 return 3;
704 case 10:
705 S->t = "II*";
706 S->tnum = 2.5;
707 S->pages = "160-174";
708 S->g = cyclic(1);
709 return 9;
710 default: pari_err_BUG("get_red [type]");
711 S->t = "";
712 S->tnum = 0;
713 S->pages = ""; /* gcc -Wall */
714 S->g = NULL;
715 return -1; /*LCOV_EXCL_LINE*/
716 }
717 }
718
719 /* reduce a/b; assume b > 0 */
720 static void
ssQ_red(long a,long b,long * n,long * d)721 ssQ_red(long a, long b, long *n, long *d)
722 {
723 long g = ugcd(labs(a), b);
724 if (g > 1) { a /= g; b /= g; }
725 *n = a; *d = b;
726 }
727 /* denom(a/b); assume b > 0 */
728 static long
ssQ_denom(long a,long b)729 ssQ_denom(long a, long b)
730 {
731 long g = ugcd(labs(a), b);
732 return g == 1? b: b / g;
733 }
734 /* n = lcm(d, denom(a/b)); r = (a/b * n mod n); assume b > 0 and d > 0 */
735 static void
get_nr(long d,long a,long b,long * n,long * r)736 get_nr(long d, long a, long b, long *n, long *r)
737 {
738 long c, A, B;
739 ssQ_red(a, b, &A,&B);
740 c = d / ugcd(d, B);
741 *n = B * c;
742 *r = umodsu(A * c, *n);
743 }
744 /* n = lcm(denom(a/b), denom(c/d)); r = (a/b * n mod n); q = (c/d * n mod n);
745 * assume b > 0 and d > 0 */
746 static void
get_nrq(long a,long b,long c,long d,long * n,long * r,long * q)747 get_nrq(long a, long b, long c, long d, long *n, long *r, long *q)
748 {
749 long g, A, B, C, D;
750 ssQ_red(a, b, &A,&B);
751 ssQ_red(c, d, &C,&D);
752 g = ugcd(B,D);
753 *n = B * (D/g);
754 *r = umodsu(A * (D/g), *n);
755 *q = umodsu(C * (B/g), *n);
756 }
757
758 /* Ip->tt = 1 */
759 static long
tame_1(struct igusa * I,struct igusa_p * Ip)760 tame_1(struct igusa *I, struct igusa_p *Ip)
761 {
762 GEN p = Ip->p, val = Ip->val;
763 long condp = -1, va0, va5, r, n;
764 va0 = myval(I->a0,p);
765 va5 = myval(I->A5,p);
766 if (!gequal0(I->A5) && 20*va0+val[6] > 6*va5)
767 get_nr(ssQ_denom(5*val[6]-6*va5, 40), val[6]-2*va5, 20, &n,&r);
768 else
769 get_nr(ssQ_denom(5*va0-val[6], 10), 10*va0-val[6], 30, &n,&r);
770 switch(n)
771 {
772 case 1:
773 condp = 0;
774 Ip->type = "[I{0-0-0}] page 155";
775 Ip->neron = cyclic(1); break;
776 case 2:
777 switch(r)
778 {
779 case 0:
780 condp = 4;
781 Ip->type = "[I*{0-0-0}] page 155";
782 Ip->neron = mkvecsmall4(2,2,2,2); break;
783 case 1:
784 condp = 2;
785 Ip->type = "[II] page 155";
786 Ip->neron = cyclic(1); break;
787 default: pari_err_BUG("tame_1 [bug1]");
788 }
789 break;
790 case 4:
791 condp = 4;
792 Ip->type = "[VI] page 156";
793 Ip->neron = dicyclic(2,2); break;
794 default: pari_err_BUG("tame_1 [bug8]");
795 }
796 return condp;
797 }
798
799 /* (4.2) */
800 static long
tame_234_init(struct igusa * I,struct igusa_p * Ip,long * n,long * q,long * r)801 tame_234_init(struct igusa *I, struct igusa_p *Ip, long *n, long *q, long *r)
802 {
803 long va0, va5, vb2, v12 = -1, flc = 1;
804 GEN p = Ip->p;
805 switch(Ip->tt)
806 {
807 case 2: v12 = myval(I->i12, Ip->p); break;
808 case 3: v12 = 3*myval(I->i4, Ip->p); break;
809 case 4: v12 = 6*myval(I->j2, Ip->p); break;
810 }
811 va0 = myval(I->a0,p);
812 va5 = myval(I->A5,p);
813 vb2 = myval(I->B2,p);
814 if (9*vb2 >= 6*va0+v12 && 36*va5 >= 120*va0+5*v12)
815 {
816 get_nrq(12*va0-v12,36, 6*va0-v12,12, n, r, q);
817 }
818 else if (120*va0+5*v12 > 36*va5 && 60*vb2 >= 12*va5+5*v12)
819 {
820 ssQ_red(36*va5-25*v12,240, q,n);
821 *r = umodsu(-2* *q, *n);
822 }
823 else /* 6*va0+v12 > 9*vb2 && 12*va5+5*v12 > 60*vb2 */
824 {
825 get_nrq(v12-6*vb2,12, v12-9*vb2,12, n,r,q);
826 flc = 0;
827 }
828 return flc;
829 }
830
831 /* Ip->tt = 2 */
832 static long
tame_2(struct igusa * I,struct igusa_p * Ip)833 tame_2(struct igusa *I, struct igusa_p *Ip)
834 {
835 long condp = -1, d, n, q, r;
836 GEN val = Ip->val;
837 (void)tame_234_init(I, Ip, &n, &q, &r);
838 d = n * (6*val[6]-5*val[7]) / 6;
839 switch(n)
840 {
841 case 1: condp = 1;
842 Ip->type = stack_sprintf("[I{%ld-0-0}] page 170", d);
843 Ip->neron = cyclic(d); break;
844 case 2:
845 switch(r)
846 {
847 case 0: condp = 4;
848 Ip->type = stack_sprintf("[I*{%ld-0-0}] page 171",d/2);
849 Ip->neron = shallowconcat(dicyclic(2,2),groupH(d/2)); break;
850 case 1:
851 switch(q)
852 {
853 case 0: condp = 2;
854 Ip->type = stack_sprintf("[II*{%ld-0}] page 172",d/2);
855 Ip->neron = cyclic(1); break;
856 case 1: condp = 3;
857 Ip->type = stack_sprintf("[II{%ld-0}] page 171",d/2);
858 Ip->neron = cyclic(2*d); break;
859 default: pari_err_BUG("tame2 [bug10]");
860 }
861 break;
862 default: pari_err_BUG("tame2 [bug11]");
863 }
864 break;
865 case 3: condp = 3;
866 Ip->neron = cyclic(d);
867 switch(r)
868 {
869 case 1:
870 Ip->type = stack_sprintf("[II{%ld}-IV] page 175", (d-2)/3);
871 break;
872 case 2:
873 Ip->type = stack_sprintf("[II{%ld}-IV*] page 175", (d-1)/3);
874 break;
875 default: pari_err_BUG("tame2 [bug12]");
876 }
877 break;
878 case 4:
879 switch(r)
880 {
881 case 1:
882 switch(q)
883 {
884 case 1: condp = 3;
885 Ip->type = stack_sprintf("[II{%ld}-III] page 177",(d-2)/4);
886 Ip->neron = cyclic(d/2); break;
887 case 3: condp = 4;
888 Ip->type = stack_sprintf("[II*{%ld}-III*] page 178",(d-2)/4);
889 Ip->neron = cyclic(8); break;
890 default: pari_err_BUG("tame2 [bug13]");
891 }
892 break;
893 case 3:
894 switch(q)
895 {
896 case 1: condp = 4;
897 Ip->type = stack_sprintf("[II*{%ld}-III] page 178",(d-2)/4);
898 Ip->neron = cyclic(8); break;
899 case 3: condp = 3;
900 Ip->type = stack_sprintf("[II{%ld}-III*] page 178",(d-2)/4);
901 Ip->neron = cyclic(d/2); break;
902 default: pari_err_BUG("tame2 [bug14]");
903 }
904 break;
905 default: pari_err_BUG("tame2 [bug15]");
906 }
907 break;
908 case 6:
909 switch(r)
910 {
911 case 2: condp = 4;
912 Ip->type = stack_sprintf("[II*-II*{%ld}] page 176", (d-4)/6);
913 Ip->neron = groupH((d+2)/6); break;
914 case 4: condp = 4;
915 Ip->type = stack_sprintf("[II-II*{%ld}] page 176", (d-2)/6);
916 Ip->neron = groupH((d+4)/6); break;
917 default: pari_err_BUG("tame2 [bug16]");
918 }
919 break;
920 default: pari_err_BUG("tame2 [bug17]");
921 }
922 return condp;
923 }
924
925 /* Ip->tt = 3 */
926 static long
tame_3(struct igusa * I,struct igusa_p * Ip)927 tame_3(struct igusa *I, struct igusa_p *Ip)
928 {
929 long condp = -1, n, q, r, va5, d1, d2;
930 long flc = tame_234_init(I, Ip, &n, &q, &r);
931 GEN val = Ip->val;
932
933 va5 = 2*val[6]-5*val[3];
934 d1 = minss(n * (val[7]-3*val[3]), n * va5 / 4);
935 d2 = n * va5 / 2 - d1;
936 switch(n)
937 {
938 case 1: condp = 2;
939 Ip->type = stack_sprintf("[I{%ld-%ld-0}] page 179", d1,d2);
940 Ip->neron = dicyclic(d1,d2); break;
941 case 2:
942 switch(r)
943 {
944 case 0: condp = 4;
945 Ip->type = stack_sprintf("[I*{%ld-%ld-0}] page 180", d1/2,d2/2);
946 Ip->neron = shallowconcat(groupH(d1/2),groupH(d2/2)); break;
947 case 1: condp = 3;
948 if (flc)
949 {
950 Ip->type = stack_sprintf("[2I{%ld}-0] page 181", d1);
951 Ip->neron = cyclic(d1);
952 }
953 else
954 { /* FIXME: "or" same with d1<->d2 */
955 Ip->type = stack_sprintf("[II{%ld-%ld}] page 182",d1/2,d2/2);
956 Ip->neron = ((d1*d2-4)&7)? cyclic(2*d1): dicyclic(d1,2);
957 }
958 break;
959 default: pari_err_BUG("tame3 [bug20]");
960 }
961 break;
962 case 4: condp = 4;
963 Ip->type = stack_sprintf("[III{%ld}] page 182", d1/2);
964 Ip->neron = groupH(d1/2); break;
965 default: pari_err_BUG("tame3 [bug21]");
966 }
967 return condp;
968 }
969
970 /* Ip->tt = 4 */
971 static long
tame_4(struct igusa * I,struct igusa_p * Ip)972 tame_4(struct igusa *I, struct igusa_p *Ip)
973 {
974 long condp = -1, d1,d2,d3, f1,f2, g, h, n, q, r, vl,vn,vm, e1,e2,e3;
975 GEN val = Ip->val;
976 (void)tame_234_init(I, Ip, &n, &q, &r);
977 vl = val[6]-5*val[1];
978 vn = val[7]-6*val[1];
979 vm = val[2]-2*val[1]; /* all >= 0 */
980 e1 = min3(2*vl, 3*vn, 6*vm);
981 e2 = minss(6*vl - e1, 12*vn - 2*e1); /* >= 0 */
982 e3 = 12*vl - (2*e1+e2); /* >= 0 */
983 d1 = e1*n / 6;
984 d2 = e2*n / 12;
985 d3 = e3*n / 12;
986 g = d1*d2 + d1*d3 + d2*d3;
987 h = ugcd(ugcd(d1,d2),d3);
988 switch(n)
989 {
990 case 1: condp = 2;
991 Ip->type = stack_sprintf("[I{%ld-%ld-%ld}] page 182",d1,d2,d3);
992 Ip->neron = dicyclic(h,g/h); break;
993 case 2:
994 switch(r)
995 {
996 case 0: condp = 4;
997 Ip->type = stack_sprintf("[I*{%ld-%ld-%ld}] page 183",d1/2,d2/2,d3/2);
998 Ip->neron = shallowconcat(groupH(g/4), groupH(2-((h&2)>>1))); break;
999 case 1:
1000 if (d1 == d2 || d1 == d3) f2 = d1;
1001 else if (d2 == d3) f2 = d2;
1002 else {
1003 pari_err_BUG("tame4 [bug23]");
1004 return -1; /*LCOV_EXCL_LINE*/
1005 }
1006 f1 = d1+d2+d3-2*f2;
1007 switch(q)
1008 {
1009 case 0: condp = 3;
1010 Ip->type = stack_sprintf("[II*{%ld-%ld}] page 184", f1/2,f2);
1011 Ip->neron = cyclic(f2); break;
1012 case 1: condp = 3;
1013 Ip->type = stack_sprintf("[II{%ld-%ld}] page 183", f1/2,f2);
1014 Ip->neron = cyclic(2*f1+f2); break;
1015 default: pari_err_BUG("tame4 [bug24]");
1016 }
1017 break;
1018 default: pari_err_BUG("tame4 [bug25]");
1019 }
1020 break;
1021 case 3: condp = 4;
1022 Ip->type = stack_sprintf("[III{%ld}] page 184",d1);
1023 Ip->neron = (d1%3)? cyclic(9): dicyclic(3,3); break;
1024 case 6: condp = 4;
1025 Ip->type = stack_sprintf("[III*{%ld}] page 184",d1/2);
1026 Ip->neron = cyclic(1); break;
1027 default: pari_err_BUG("tame4 [bug26]");
1028 }
1029 return condp;
1030 }
1031
1032 /* p = 3 */
1033 static void
tame_567_init_3(struct igusa_p * Ip,long dk,long * pd,long * pn,long * pdm,long * pr)1034 tame_567_init_3(struct igusa_p *Ip, long dk,
1035 long *pd, long *pn, long *pdm, long *pr)
1036 {
1037 long n = 1 + Ip->r1/6;
1038 *pd = n * dk / 36; /* / (12*Ip->eps) */
1039 *pn = n;
1040 *pr = -1; /* unused */
1041 *pdm = 0;
1042 }
1043
1044 /* (4.3) */
1045 static void
tame_567_init(struct igusa * I,struct igusa_p * Ip,long dk,long * pd,long * pn,long * pdm,long * pr)1046 tame_567_init(struct igusa *I, struct igusa_p *Ip, long dk,
1047 long *pd, long *pn, long *pdm, long *pr)
1048 {
1049 long ndk, ddk;
1050 GEN p = Ip->p, val = Ip->val;
1051
1052 if (equaliu(p,3)) { tame_567_init_3(Ip, dk, pd, pn, pdm, pr); return; }
1053 /* assume p > 3, Ip->eps = 1 */
1054 ssQ_red(dk, 12, &ndk, &ddk);
1055 if (! odd(val[8]))
1056 {
1057 long va0 = myval(I->a0,p), va2 = myval(I->A2,p), va3 = myval(I->A3,p);
1058 long va5 = myval(I->A5,p), vb2 = myval(I->B2,p);
1059 long v1 = 2*va3-4*va0-val[1], v2 = 6*va5-20*va0-5*val[1];
1060 long v3 = 3*vb2-2*va0-2*val[1], v4 = 10*vb2-2*va5-5*val[1];
1061 if (v3 >= 0 && v2 >= 0 && v1 >= 0)
1062 {
1063 if (v1==0 || v2==0) get_nr(ddk, va0+val[1], 6,pn,pr); /* Prop 4.3.1 (a) */
1064 else
1065 { /* Prop 4.3.1 (d) */
1066 long v5 = myval(subii(mulii(I->A2,I->A3),mului(3,I->A5)),p);
1067 if (gequal0(I->A2)) pari_err_BUG("tame567 [bug27]");
1068 get_nr(ddk, 12*va0 + min3(dk, 6*va3-9*va2, 4*v5 - 10*va2), 24, pn,pr);
1069 }
1070 }
1071 else if (v2 < 0 && v4 >= 0)
1072 get_nr(ddk, 2*va5+val[1], 8, pn,pr); /* Prop 4.3.1 (b) */
1073 else /* (v3 < 0 && v4 < 0) */
1074 get_nr(ddk, vb2, 4, pn,pr); /* Prop 4.3.1 (c) */
1075 *pd = (*pn/ddk) * ndk;
1076 }
1077 else
1078 {
1079 *pr = ndk;
1080 *pn = 2*ddk;
1081 *pd = 2*ndk;
1082 }
1083 *pdm = umodsu(*pd, *pn);
1084 }
1085
1086 static long
tame_5(struct igusa * I,struct igusa_p * Ip)1087 tame_5(struct igusa *I, struct igusa_p *Ip)
1088 {
1089 long condp = -1, d, n, dm, r, dk;
1090 GEN val = Ip->val;
1091
1092 dk = Ip->eps*val[6]-5*val[8];
1093 tame_567_init(I, Ip, dk, &d, &n, &dm, &r);
1094 if (! odd(val[8]))
1095 {
1096 switch(n)
1097 {
1098 case 1: condp = 0;
1099 Ip->type = stack_sprintf("[I{0}-I{0}-%ld] page 158", d);
1100 Ip->neron = cyclic(1); break;
1101 case 2:
1102 switch(dm)
1103 {
1104 case 0: condp = 4;
1105 Ip->type = stack_sprintf("[I*{0}-I*{0}-%ld] page 158",(d-2)/2);
1106 Ip->neron = mkvecsmall4(2,2,2,2); break;
1107 case 1: condp = 2;
1108 Ip->type = stack_sprintf("[I{0}-I*{0}-%ld] page 159",(d-1)/2);
1109 Ip->neron = dicyclic(2,2); break;
1110 }
1111 break;
1112 case 3:
1113 switch(dm)
1114 {
1115 case 0: condp = 4;
1116 Ip->type = stack_sprintf("[IV-IV*-%ld] page 165",(d-3)/3);
1117 Ip->neron = dicyclic(3,3); break;
1118 case 1:
1119 switch(r)
1120 {
1121 case 0: case 1: condp = 2;
1122 Ip->type = stack_sprintf("[I{0}-IV-%ld] page 160",(d-1)/3);
1123 Ip->neron = cyclic(3); break;
1124 case 2: condp = 4;
1125 Ip->type = stack_sprintf("[IV*-IV*-%ld] page 166",(d-4)/3);
1126 Ip->neron = dicyclic(3,3); break;
1127 }
1128 break;
1129 case 2:
1130 switch(r)
1131 {
1132 case 0: case 2: condp = 2;
1133 Ip->type = stack_sprintf("[I{0}-IV*-%ld] page 160",(d-2)/3);
1134 Ip->neron = cyclic(3); break;
1135 case 1: condp = 4;
1136 Ip->type = stack_sprintf("[IV-IV-%ld] page 165",(d-2)/3);
1137 Ip->neron = dicyclic(3,3); break;
1138 }
1139 break;
1140 }
1141 break;
1142 case 4:
1143 switch(dm)
1144 {
1145 case 0: condp = 4;
1146 Ip->type = stack_sprintf("[III-III*-%ld] page 169",(d-4)/4);
1147 Ip->neron = dicyclic(2,2); break;
1148 case 1:
1149 switch(r)
1150 {
1151 case 0: case 1: condp = 2;
1152 Ip->type = stack_sprintf("[I{0}-III-%ld] page 161",(d-1)/4);
1153 Ip->neron = cyclic(2); break;
1154 case 2: case 3: condp = 4;
1155 Ip->type = stack_sprintf("[I*{0}-III*-%ld] page 162",(d-5)/4);
1156 Ip->neron = mkvecsmall3(2,2,2); break;
1157 }
1158 break;
1159 case 2: condp = 4;
1160 Ip->neron = dicyclic(2,2);
1161 switch(r)
1162 {
1163 case 1:
1164 Ip->type = stack_sprintf("[III-III-%ld] page 169",(d-2)/4);
1165 break;
1166 case 3:
1167 Ip->type = stack_sprintf("[III*-III*-%ld] page 169",(d-6)/4);
1168 break;
1169 default: pari_err_BUG("tame5 [bug29]");
1170 }
1171 break;
1172 case 3:
1173 switch(r)
1174 {
1175 case 0: case 3: condp = 2;
1176 Ip->type = stack_sprintf("[I{0}-III*-%ld] page 162",(d-3)/4);
1177 Ip->neron = cyclic(2); break;
1178 case 1: case 2: condp = 4;
1179 Ip->type = stack_sprintf("[I*{0}-III-%ld] page 162",(d-3)/4);
1180 Ip->neron = mkvecsmall3(2,2,2); break;
1181 }
1182 break;
1183 }
1184 break;
1185 case 6:
1186 switch(dm)
1187 {
1188 case 0: condp = 4;
1189 Ip->type = stack_sprintf("[II-II*-%ld] page 163",(d-6)/6);
1190 Ip->neron = cyclic(1); break;
1191 case 1:
1192 switch(r)
1193 {
1194 case 0: case 1: condp = 2;
1195 Ip->type = stack_sprintf("[I{0}-II-%ld] page 159",(d-1)/6);
1196 Ip->neron = cyclic(1); break;
1197 case 2: case 5: condp = 4;
1198 Ip->type = stack_sprintf("[II*-IV-%ld] page 164",(d-7)/6);
1199 Ip->neron = cyclic(3); break;
1200 case 3: case 4: condp = 4;
1201 Ip->type = stack_sprintf("[I*{0}-IV*-%ld] page 161",(d-7)/6);
1202 Ip->neron = mkvecsmall2(6,2); break;
1203 }
1204 break;
1205 case 2:
1206 switch(r)
1207 {
1208 case 1: condp = 4;
1209 Ip->type = stack_sprintf("[II-II-%ld] page 163",(d-2)/6);
1210 Ip->neron = cyclic(1); break;
1211 case 3: case 5: condp = 4;
1212 Ip->type = stack_sprintf("[I*{0}-II*-%ld] page 160",(d-8)/6);
1213 Ip->neron = dicyclic(2,2); break;
1214 default: pari_err_BUG("tame5 [bug30]");
1215 }
1216 break;
1217 case 3:
1218 Ip->neron = cyclic(3);
1219 switch(r)
1220 {
1221 case 1: case 2: condp = 4;
1222 Ip->type = stack_sprintf("[II-IV-%ld] page 164",(d-3)/6);
1223 break;
1224 case 4: case 5: condp = 4;
1225 Ip->type = stack_sprintf("[II*-IV*-%ld] page 164",(d-9)/6);
1226 break;
1227 default: pari_err_BUG("tame5 [bug31]");
1228 }
1229 break;
1230 case 4:
1231 switch(r)
1232 {
1233 case 1: case 3: condp = 4;
1234 Ip->type = stack_sprintf("[I*{0}-II-%ld] page 160",(d-4)/6);
1235 Ip->neron = dicyclic(2,2); break;
1236 case 5: condp = 4;
1237 Ip->type = stack_sprintf("[II*-II*-%ld] page 163",(d-10)/6);
1238 Ip->neron = cyclic(1); break;
1239 default: pari_err_BUG("tame5 [bug32]");
1240 }
1241 break;
1242 case 5:
1243 switch(r)
1244 {
1245 case 0: case 5: condp = 2;
1246 Ip->type = stack_sprintf("[I{0}-II*-%ld] page 160",(d-5)/6);
1247 Ip->neron = cyclic(1); break;
1248 case 1: case 4: condp = 4;
1249 Ip->type = stack_sprintf("[II-IV*-%ld] page 164",(d-5)/6);
1250 Ip->neron = cyclic(3); break;
1251 case 2: case 3: condp = 4;
1252 Ip->type = stack_sprintf("[I*{0}-IV-%ld] page 161",(d-5)/6);
1253 Ip->neron = mkvecsmall2(6,2); break;
1254 }
1255 break;
1256 default: pari_err_BUG("tame5 [bug33]");
1257 }
1258 break;
1259 case 12:
1260 condp = 4;
1261 switch(dm)
1262 {
1263 case 1:
1264 switch(r)
1265 {
1266 case 3: case 10:
1267 Ip->type = stack_sprintf("[II*-III-%ld] page 166",(d-13)/12);
1268 Ip->neron = cyclic(2); break;
1269 case 4: case 9:
1270 Ip->type = stack_sprintf("[III*-IV-%ld] page 167",(d-13)/12);
1271 Ip->neron = cyclic(6); break;
1272 default: pari_err_BUG("tame5 [bug34]");
1273 }
1274 break;
1275 case 5:
1276 switch(r)
1277 {
1278 case 2: case 3:
1279 Ip->type = stack_sprintf("[II-III-%ld] page 166",(d-5)/12);
1280 Ip->neron = cyclic(2); break;
1281 case 8: case 9:
1282 Ip->type = stack_sprintf("[III*-IV*-%ld] page 168",(d-17)/12);
1283 Ip->neron = cyclic(6); break;
1284 default: pari_err_BUG("tame5 [bug35]");
1285 }
1286 break;
1287 case 7:
1288 switch(r)
1289 {
1290 case 3: case 4:
1291 Ip->type = stack_sprintf("[III-IV-%ld] page 167",(d-7)/12);
1292 Ip->neron = cyclic(6); break;
1293 case 9: case 10:
1294 Ip->type = stack_sprintf("[II*-III*-%ld] page 167",(d-19)/12);
1295 Ip->neron = cyclic(2); break;
1296 default: pari_err_BUG("tame5 [bug36]");
1297 }
1298 break;
1299 case 11:
1300 switch(r)
1301 {
1302 case 3: case 8:
1303 Ip->type = stack_sprintf("[III-IV*-%ld] page 168",(d-11)/12);
1304 Ip->neron = cyclic(6); break;
1305 case 2: case 9:
1306 Ip->type = stack_sprintf("[II-III*-%ld] page 166",(d-11)/12);
1307 Ip->neron = cyclic(2); break;
1308 default: pari_err_BUG("tame5 [bug37]");
1309 }
1310 break;
1311 default: pari_err_BUG("tame5 [bug38]");
1312 }
1313 break;
1314 default: pari_err_BUG("tame5 [bug39]");
1315 }
1316 }
1317 else
1318 {
1319 r %= (n >> 1);
1320 switch(n)
1321 {
1322 case 2: condp = 2;
1323 Ip->type = stack_sprintf("[2I{0}-%ld] page 159",(d/2));
1324 Ip->neron = cyclic(1); break;
1325 case 4: condp = 4;
1326 Ip->type = stack_sprintf("[2I*{0}-%ld] page 159",(d/2-1)/2);
1327 Ip->neron = dicyclic(2,2); break;
1328 case 6: condp = 4;
1329 Ip->neron = cyclic(3);
1330 switch(r)
1331 {
1332 case 1:
1333 Ip->type = stack_sprintf("[2IV-%ld] page 165",(d/2-1)/3);
1334 break;
1335 case 2:
1336 Ip->type = stack_sprintf("[2IV*-%ld] page 165",(d/2-2)/3);
1337 break;
1338 default: pari_err_BUG("tame5 [bug40]");
1339 }
1340 break;
1341 case 8: condp = 4;
1342 Ip->neron = cyclic(2);
1343 switch(r)
1344 {
1345 case 1:
1346 Ip->type = stack_sprintf("[2III-%ld] page 168",(d/2-1)/4);
1347 break;
1348 case 3:
1349 Ip->type = stack_sprintf("[2III*-%ld] page 168",(d/2-3)/4);
1350 break;
1351 default: pari_err_BUG("tame5 [bug41]");
1352 }
1353 break;
1354 case 12: condp = 4;
1355 Ip->neron = cyclic(1);
1356 switch(r)
1357 {
1358 case 1:
1359 Ip->type = stack_sprintf("[2II-%ld] page 162",(d/2-1)/6);
1360 break;
1361 case 5:
1362 Ip->type = stack_sprintf("[2II*-%ld] page 163",(d/2-5)/6);
1363 break;
1364 default: pari_err_BUG("tame5 [bug42]");
1365 }
1366 break;
1367 default: pari_err_BUG("tame5 [bug43]");
1368 }
1369 }
1370 return condp;
1371 }
1372
1373 static long
tame_6(struct igusa * I,struct igusa_p * Ip)1374 tame_6(struct igusa *I, struct igusa_p *Ip)
1375 {
1376 long condp = -1, d, d1, n, dm, r, dk;
1377 GEN val = Ip->val;
1378
1379 dk = Ip->eps*val[7]-6*val[8];
1380 tame_567_init(I, Ip, dk, &d, &n, &dm, &r);
1381 d1 = n * (Ip->eps*(val[6]-val[7])+val[8]) / Ip->eps;
1382 switch(n)
1383 {
1384 case 1: condp = 1;
1385 Ip->type = stack_sprintf("[I{0}-I{%ld}-%ld] page 170",d1,d);
1386 Ip->neron = cyclic(d1); break;
1387 case 2:
1388 switch(dm)
1389 {
1390 case 0: condp = 4;
1391 Ip->type=stack_sprintf("[I*{0}-I*{%ld}-%ld] page 171", d1/2,(d-2)/2);
1392 Ip->neron = shallowconcat(groupH(d1/2), dicyclic(2,2)); break;
1393 case 1: return -1;
1394 default: pari_err_BUG("tame6 [bug44]");
1395 }
1396 break;
1397 case 3: condp = 3;
1398 Ip->neron = dicyclic(3,d1/3);
1399 switch(dm)
1400 {
1401 case 1:
1402 Ip->type = stack_sprintf("[I{%ld}-IV-%ld] page 173",d1/3,(d-1)/3);
1403 break;
1404 case 2:
1405 Ip->type = stack_sprintf("[I{%ld}-IV*-%ld] page 173",d1/3,(d-2)/3);
1406 break;
1407 default: pari_err_BUG("tame6 [bug45]");
1408 }
1409 break;
1410 case 4:
1411 switch(dm)
1412 {
1413 case 1:
1414 switch(r)
1415 {
1416 case 0: case 1: condp = 3;
1417 Ip->type=stack_sprintf("[I{%ld}-III-%ld] page 176",d1/4,(d-1)/4);
1418 Ip->neron = dicyclic(2,d1/4); break;
1419 case 2: case 3: condp = 4;
1420 Ip->type=stack_sprintf("[I*{%ld}-III*-%ld] page 177",d1/4,(d-5)/4);
1421 Ip->neron = shallowconcat(groupH(d1/4), cyclic(2)); break;
1422 default: pari_err_BUG("tame6 [bug46]");
1423 }
1424 break;
1425 case 3:
1426 switch(r)
1427 {
1428 case 0: case 3: condp = 3;
1429 Ip->type=stack_sprintf("[I{%ld}-III*-%ld] page 176",d1/4,(d-3)/4);
1430 Ip->neron = dicyclic(2,d1/4); break;
1431 case 1: case 2: condp = 4;
1432 Ip->type=stack_sprintf("[I*{%ld}-III-%ld] page 177",d1/4,(d-3)/4);
1433 Ip->neron = shallowconcat(groupH(d1/4), cyclic(2)); break;
1434 default: pari_err_BUG("tame6 [bug47]");
1435 }
1436 break;
1437 default: pari_err_BUG("tame6 [bug48]");
1438 }
1439 break;
1440 case 6:
1441 switch(dm)
1442 {
1443 case 1:
1444 switch(r)
1445 {
1446 case 0: case 1: condp = 3;
1447 Ip->type = stack_sprintf("[I{%ld}-II-%ld] page 172",d1/6,(d-1)/6);
1448 Ip->neron = cyclic(d1/6); break;
1449 case 3: case 4: condp = 4;
1450 Ip->type=stack_sprintf("[I*{%ld}-IV*-%ld] page 174",d1/6,(d-7)/6);
1451 Ip->neron = shallowconcat(groupH(d1/6), cyclic(3)); break;
1452 default: pari_err_BUG("tame6 [bug49]");
1453 }
1454 break;
1455 case 2: condp = 4;
1456 Ip->type = stack_sprintf("[I*{%ld}-II*-%ld] page 174",d1/6,(d-8)/6);
1457 Ip->neron = groupH(d1/6); break;
1458 case 4: condp = 4;
1459 Ip->type = stack_sprintf("[I*{%ld}-II-%ld] page 173",d1/6,(d-4)/6);
1460 Ip->neron = groupH(d1/6); break;
1461 case 5:
1462 switch(r)
1463 {
1464 case 0: case 5: condp = 3;
1465 Ip->type=stack_sprintf("[I{%ld}-II*-%ld] page 172",d1/6,(d-5)/6);
1466 Ip->neron = cyclic(d1/6); break;
1467 case 2: case 3: condp = 4;
1468 Ip->type=stack_sprintf("[I*{%ld}-IV-%ld] page 174",d1/6,(d-5)/6);
1469 Ip->neron = shallowconcat(groupH(d1/6), cyclic(3)); break;
1470 default: pari_err_BUG("tame6 [bug50]");
1471 }
1472 break;
1473 default: pari_err_BUG("tame6 [bug51]");
1474 }
1475 break;
1476 default: pari_err_BUG("tame6 [bug52]");
1477 }
1478 return condp;
1479 }
1480
1481 static long
tame_7(struct igusa * I,struct igusa_p * Ip)1482 tame_7(struct igusa *I, struct igusa_p *Ip)
1483 {
1484 long condp = -1, d, D, d1, d2, n, dm, r, dk;
1485 GEN val = Ip->val;
1486
1487 dk = 3*(Ip->eps*val[3]-2*val[8]);
1488 tame_567_init(I, Ip, dk, &d, &n, &dm, &r);
1489 D = n * (Ip->eps*(val[6]-3*val[3])+val[8]) / Ip->eps;
1490 d1 = minss(n * (val[7]-3*val[3]), D/2);
1491 d2 = D - d1;
1492 /* d1 <= d2 */
1493 switch(n)
1494 {
1495 case 1: condp = 2;
1496 Ip->type = stack_sprintf("[I{%ld}-I{%ld}-%ld] page 179",d1,d2,d);
1497 Ip->neron = dicyclic(d1,d2); break;
1498 case 2:
1499 if (odd(val[8]))
1500 {
1501 condp = 3;
1502 Ip->type = stack_sprintf("[2I{%ld}-%ld] page 181",d1,d/2);
1503 Ip->neron = cyclic(d1);
1504 }
1505 else if (dm == 0)
1506 {
1507 condp = 4;
1508 Ip->type = stack_sprintf("[I*{%ld}-I*{%ld}-%ld] page 180", d1/2,d2/2,(d-2)/2);
1509 Ip->neron = shallowconcat(groupH(d1/2),groupH(d2/2));
1510 }
1511 else
1512 {
1513 GEN H;
1514 if (d1 != d2) return -1;
1515 condp = 3; H = groupH(d1/2);
1516 Ip->type = stack_sprintf("[I{%ld}-I*{%ld}-%ld] page 180", d1/2,d1/2,(d-1)/2);
1517 Ip->neron = shallowconcat(H, H);
1518 }
1519 break;
1520 case 4: condp = 4;
1521 Ip->type = stack_sprintf("[2I*{%ld}-%ld] page 181",d1/2,(d-2)/4);
1522 Ip->neron = groupH(d1/2); break;
1523 default: pari_err_BUG("tame7 [bug55]");
1524 }
1525 return condp;
1526 }
1527
1528 static long labelm3(GEN polh, long t60, long alpha, long Dmin, struct igusa *I, struct igusa_p *Ip);
1529 static long
tame(GEN polh,long t60,long alpha,long Dmin,struct igusa * I,struct igusa_p * Ip)1530 tame(GEN polh, long t60, long alpha, long Dmin, struct igusa *I, struct igusa_p *Ip)
1531 {
1532 long d;
1533 Ip->tame = 1;
1534 switch(Ip->tt)
1535 {
1536 case 1: return tame_1(I,Ip);
1537 case 2: return tame_2(I,Ip);
1538 case 3: return tame_3(I,Ip);
1539 case 4: return tame_4(I,Ip);
1540 case 5: return tame_5(I,Ip);
1541 case 6: d = tame_6(I,Ip); break;
1542 default:d = tame_7(I,Ip); break;
1543 }
1544 if (d < 0) d = labelm3(polh,t60,alpha,Dmin,I,Ip); /* => tt=6 or 7 */
1545 return d;
1546 }
1547
1548 /* maxc = maximum conductor valuation at p */
1549 static long
get_maxc(GEN p)1550 get_maxc(GEN p)
1551 {
1552 switch (itos_or_0(p))
1553 {
1554 case 2: return 20; break;
1555 case 3: return 10; break;
1556 case 5: return 9; break;
1557 default: return 4; break; /* p > 5 */
1558 }
1559 }
1560
1561 /* p = 3 */
1562 static long
quartic(GEN polh,long alpha,long Dmin,struct igusa_p * Ip)1563 quartic(GEN polh, long alpha, long Dmin, struct igusa_p *Ip)
1564 {
1565 GEN val = Ip->val, p = Ip->p;
1566 GEN polf = polymini_zi2(ZX_Z_mul(polh, powiu(p, alpha)));
1567 long condp = -1, d, R, r1, beta;
1568 r1 = polf[1];
1569 beta = polf[2];
1570 R = beta/2;
1571 switch(Ip->tt)
1572 {
1573 case 1: case 5: d = 0;break;
1574 case 3: d = val[6] - 5*val[3]/2;break;
1575 case 7: d = val[6] - 3*val[3] + val[8]/Ip->eps;break;
1576 default: pari_err_BUG("quartic [type choices]");
1577 d = 0; /*LCOV_EXCL_LINE*/
1578 }
1579 switch(r1)
1580 {
1581 case 0:
1582 if (d)
1583 {
1584 condp = 3;
1585 Ip->type = stack_sprintf("[2I{%ld}-%ld] page 181",d,R);
1586 Ip->neron = cyclic(d);
1587 }
1588 else
1589 {
1590 condp = 2;
1591 Ip->neron = cyclic(1);
1592 if (R) Ip->type = stack_sprintf("[2I{0}-%ld] page 159",R);
1593 else Ip->type = "[II] page 155";
1594 }
1595 break;
1596 case 6: condp = 4;
1597 Ip->type = stack_sprintf("[2I*{%ld}-%ld] pages 159, 181",d,R);
1598 Ip->neron = dicyclic(2,2); break;
1599 case 3: condp = 4;
1600 Ip->type = stack_sprintf("[2III-%ld] page 168",R);
1601 Ip->neron = cyclic(2); break;
1602 case 9: condp = 4;
1603 Ip->type = stack_sprintf("[2III*-%ld] page 168",R);
1604 Ip->neron = cyclic(2); break;
1605 case 2: condp = Dmin-12*R-13;
1606 Ip->type = stack_sprintf("[2II-%ld] page 162",R);
1607 Ip->neron = cyclic(1); break;
1608 case 8: condp = Dmin-12*R-19;
1609 Ip->type = stack_sprintf("[2IV*-%ld] page 165",R);
1610 Ip->neron = cyclic(3); break;
1611 case 4: condp = Dmin-12*R-15;
1612 Ip->type = stack_sprintf("[2IV-%ld] page 165",R);
1613 Ip->neron = cyclic(3); break;
1614 case 10: condp = Dmin-12*R-21;
1615 Ip->type = stack_sprintf("[2II*-%ld] page 163",R);
1616 Ip->neron = cyclic(1); break;
1617 default: pari_err_BUG("quartic [type1]");
1618 }
1619 if (condp > get_maxc(p) || condp < 0) pari_err_BUG("quartic [conductor]");
1620 return condp;
1621 }
1622
1623 static long
litredtp(long alpha,long alpha1,long t60,long t60_1,GEN polh,GEN polh1,long Dmin,long R,struct igusa * I,struct igusa_p * Ip)1624 litredtp(long alpha, long alpha1, long t60, long t60_1, GEN polh, GEN polh1,
1625 long Dmin, long R, struct igusa *I, struct igusa_p *Ip)
1626 {
1627 GEN val = Ip->val, p = Ip->p;
1628 long condp = -1, indice, d;
1629
1630 if ((Ip->r1 == 0||Ip->r1 == 6) && (Ip->r2 == 0||Ip->r2 == 6))
1631 { /* (r1,r2) = (0,0), (0,6), (6,0) or (6,6) */
1632 if (Ip->tt == 5)
1633 {
1634 switch(Ip->r1 + Ip->r2)
1635 {
1636 case 0: /* (0,0) */
1637 condp = 0;
1638 Ip->type = stack_sprintf("[I{0}-I{0}-%ld] page 158",R);
1639 Ip->neron = cyclic(1); break;
1640 case 6: /* (0,6) or (6,0) */
1641 condp = 2;
1642 Ip->type = stack_sprintf("[I{0}-I*{0}-%ld] page 159",R);
1643 Ip->neron = dicyclic(2,2); break;
1644 case 12: /* (6,6) */
1645 condp = 4;
1646 Ip->type = stack_sprintf("[I*{0}-I*{0}-%ld] page 158",R);
1647 Ip->neron = mkvecsmall4(2,2,2,2); break;
1648 }
1649 return condp;
1650 }
1651 if (Ip->r1 == Ip->r2) return tame(polh, t60, alpha, Dmin, I, Ip);
1652 if (Ip->tt == 6)
1653 {
1654 d = val[6] - val[7] + val[8]/Ip->eps;
1655 if (Ip->r1 && alpha1 == 0) /* H(px) / p^3 */
1656 polh1 = ZX_Z_divexact(ZX_unscale_div(polh1,p), sqri(p));
1657 if (FpX_is_squarefree(FpX_red(polh1,p),p))
1658 { indice = 0; condp = 3-Ip->r2/6; }
1659 else
1660 { indice = d; condp = 3-Ip->r1/6; }
1661 }
1662 else
1663 { /* Ip->tt == 7 */
1664 long d1;
1665 d = val[6] - 3*val[3] + val[8]/Ip->eps;
1666 if (t60_1 == 60) /* H(px) / p^3 */
1667 polh1 = ZX_Z_divexact(ZX_unscale_div(polh1,p), sqri(p));
1668 d1 = minss(val[7]-3*val[3],d/2);
1669 if (d == 2*d1) indice = d1;
1670 else
1671 {
1672 indice = discpart(polh1,p,d1+1);
1673 if (indice>= d1+1) indice = d-d1; else indice = d1;
1674 }
1675 condp = 3;
1676 }
1677 if (Ip->r1) indice = d - indice; /* (r1,r2) = (6,0) */
1678 Ip->neron = shallowconcat(cyclic(indice),groupH(d-indice));
1679 Ip->type = stack_sprintf("[I{%ld}-I*{%ld}-%ld] page %ld",
1680 indice,d-indice,R, (Ip->tt==6)? 170L: 180L);
1681 return condp;
1682 }
1683 if (Ip->tt == 7) pari_err_BUG("litredtp [switch ri]");
1684 {
1685 struct red __S1, __S2, *S1 = &__S1, *S2 = &__S2;
1686 long f1 = get_red(S1, Ip, polh1, p, alpha1, Ip->r1);
1687 long f2 = get_red(S2, Ip, polh, p, alpha, Ip->r2);
1688 /* reorder to normalize representation */
1689 if (S1->tnum > S2->tnum || (S1->tnum == S2->tnum && f1 > f2))
1690 { struct red *S = S1; S1 = S2; S2 = S; }
1691 Ip->type = stack_sprintf("[%s-%s-%ld] pages %s", S1->t,S2->t, R, S1->pages);
1692 Ip->neron = shallowconcat(S1->g, S2->g);
1693 condp = Dmin - (f1 + f2) + ((R >= 0)? 2-12*R: 4);
1694 }
1695 if (condp > get_maxc(p)) pari_err_BUG("litredtp [conductor]");
1696 return condp;
1697 }
1698
1699 static long
labelm3(GEN h1,long t60_1,long alpha1,long Dmin,struct igusa * I,struct igusa_p * Ip)1700 labelm3(GEN h1, long t60_1, long alpha1, long Dmin, struct igusa *I, struct igusa_p *Ip)
1701 {
1702 GEN h, pm, vs, val = Ip->val, p = Ip->p;
1703 long alpha, t60, lambda, beta, R;
1704
1705 pm = polymini(ZX_Z_mul(RgX_recip6(h1), powiu(p,alpha1)), p);
1706 h = gel(pm,1); vs = gel(pm,2);
1707 lambda= vs[1];
1708 t60 = vs[2];
1709 alpha = vs[3];
1710 beta = vs[5];
1711 if (lambda != 3) pari_err_BUG("labelm3 [lambda != 3]");
1712 R = beta-(alpha1+alpha);
1713 if (odd(R)) pari_err_BUG("labelm3 [R odd]");
1714 R /= 2;
1715 if (R <= -2) pari_err_BUG("labelm3 [R <= -2]");
1716 if (val[8] % (2*Ip->eps)) pari_err_BUG("labelm3 [val(eps2)]");
1717 if (R >= 0 && (alpha+alpha1) >= 1) pari_err_BUG("labelm3 [minimal equation]");
1718 Ip->r1 = t60_1 / 10 + 6*alpha1;
1719 Ip->r2 = t60 / 10 + 6*alpha;
1720 return litredtp(alpha, alpha1, t60, t60_1, h, h1, Dmin, R, I, Ip);
1721 }
1722
1723 /* p = 3 */
1724 static long
quadratic(GEN polh,long alpha,long Dmin,struct igusa * I,struct igusa_p * Ip)1725 quadratic(GEN polh, long alpha, long Dmin, struct igusa *I, struct igusa_p *Ip)
1726 {
1727 long alpha1 = alpha, beta, t6, R;
1728 GEN vs = polymini_zi(ZX_Z_mul(polh, powiu(Ip->p,alpha)));
1729 t6 = vs[1];
1730 alpha = vs[2];
1731 beta = vs[3];
1732 R = beta-alpha;
1733 if (R >= 0 && alpha1)
1734 {
1735 Dmin -= 10;
1736 if (DEBUGLEVEL)
1737 err_printf("(Care: minimal discriminant over Z[i] smaller than over Z)\n");
1738 }
1739 Ip->r2 = Ip->r1 = t6 + 6*alpha;
1740 return litredtp(alpha, alpha, t6*10, t6*10, polh, polh, Dmin, R, I, Ip);
1741 }
1742
1743 static long
genus2localred(struct igusa * I,struct igusa_p * Ip,GEN p,GEN polmini)1744 genus2localred(struct igusa *I, struct igusa_p *Ip, GEN p, GEN polmini)
1745 {
1746 GEN val, vs, polh, list, c1, c2, c3, c4, c5, c6, prod;
1747 long i, vb5, vb6, d, Dmin, alpha, lambda, t60;
1748 long condp = -1, indice, vc6, mm, nb, dism;
1749
1750 stable_reduction(I, Ip, p);
1751 val = Ip->val; Dmin = val[6];
1752 if (Dmin == 0)
1753 {
1754 Ip->tame = 1;
1755 Ip->type = "[I{0-0-0}] page 155";
1756 Ip->neron = cyclic(1); return 0;
1757 }
1758 if (Dmin == 1)
1759 {
1760 Ip->type = "[I{1-0-0}] page 170";
1761 Ip->neron = cyclic(1); return 1;
1762 }
1763 if (Dmin == 2) switch(Ip->tt)
1764 {
1765 case 2:
1766 Ip->type = "[I{2-0-0}] page 170";
1767 Ip->neron = cyclic(2); return 1;
1768 case 3:
1769 Ip->type = "[I{1-1-0}] page 179";
1770 Ip->neron = cyclic(1); return 2;
1771 case 5:
1772 if (cmpis(p,3) <= 0) pari_err_BUG("genus2localred [tt 1]");
1773 Ip->type = "[I{0}-II-0] page 159";
1774 Ip->neron = cyclic(1); return 2;
1775 default: pari_err_BUG("genus2localred [tt 2]");
1776 }
1777 if (absequaliu(p,2)) return -1;
1778 polh = gel(polmini,1); vs = gel(polmini,2);
1779 lambda = vs[1];
1780 t60 = vs[2];
1781 alpha = vs[3];
1782 if (vs[4]) return equaliu(p,3)? quadratic(polh, alpha, Dmin, I, Ip):
1783 tame(polh, t60, alpha, Dmin, I, Ip);
1784 if (!t60 && lambda<= 2)
1785 {
1786 if (Ip->tt >= 5) pari_err_BUG("genus2localred [tt 3]");
1787 return tame(polh, t60, alpha, Dmin, I, Ip);
1788 }
1789 if (Dmin == 3)
1790 {
1791 switch(Ip->tt)
1792 {
1793 case 2: return tame(polh, t60, alpha, Dmin, I, Ip);
1794 case 3: Ip->type = "[I{2-1-0}] page 179"; Ip->neron = cyclic(2); return 2;
1795 case 4: Ip->type = "[I{1-1-1}] page 182"; Ip->neron = cyclic(3); return 2;
1796 case 5:
1797 if (equaliu(p,3) && t60 != 30)
1798 return labelm3(polh,t60,alpha,Dmin,I,Ip);
1799 Ip->type = "[I{0}-III-0] page 161"; Ip->neron = cyclic(2); return 2;
1800 case 6:
1801 if (equaliu(p,3)) pari_err_BUG("genus2localred [conductor]");
1802 Ip->type = "[I{1}-II-0] page 172"; Ip->neron = cyclic(1); return 3;
1803 }
1804 pari_err_BUG("genus2localred [switch tt 4]");
1805 return -1; /* LCOV_EXCL_LINE */
1806 }
1807 switch(lambda)
1808 {
1809 case 0:
1810 switch(t60+alpha)
1811 {
1812 case 10:
1813 condp = Dmin-1;
1814 Ip->type = "[V] page 156";
1815 Ip->neron = cyclic(3); break;
1816 case 11:
1817 condp = Dmin-11;
1818 Ip->type = "[V*] page 156";
1819 Ip->neron = cyclic(3); break;
1820 case 12:
1821 condp = Dmin-2;
1822 Ip->type = "[IX-2] page 157";
1823 Ip->neron = cyclic(5); break;
1824 case 13:
1825 condp = Dmin-12;
1826 Ip->type = "[VIII-4] page 157";
1827 Ip->neron = cyclic(1); break;
1828 case 24:
1829 condp = Dmin-8;
1830 Ip->type = "[IX-4] page 158";
1831 Ip->neron = cyclic(5);
1832 break;
1833 case 15: case 16:
1834 if (Ip->tt>= 5) pari_err_BUG("genus2localred [tt 6]");
1835 return tame(polh, t60, alpha, Dmin, I, Ip);
1836 case 20: case 21:
1837 {
1838 GEN b0, b1, b2, b3, b4, b5, b6, b02, b03, b04, b05;
1839 RgX_to_06(polh, &b0,&b1,&b2,&b3,&b4,&b5,&b6);
1840 vb5 = myval(b5,p);
1841 vb6 = myval(b6,p);
1842 if (vb6 >= 3)
1843 {
1844 if (vb5 < 2) pari_err_BUG("genus2localred [red1]");
1845 if (vb5 >= 3)
1846 {
1847 condp = Dmin-8;
1848 Ip->type = "[II*-IV-(-1)] page 164";
1849 Ip->neron = cyclic(3);
1850 }
1851 else
1852 {
1853 condp = Dmin-7;
1854 Ip->type = "[IV-III*-(-1)] page 167";
1855 Ip->neron = cyclic(6);
1856 }
1857 break;
1858 }
1859 if (dvdii(b0,p)) pari_err_BUG("genus2localred [b0]");
1860 b02 = gsqr(b0);
1861 b03 = gmul(b02, b0);
1862 b04 = gmul(b03, b0);
1863 b05 = gmul(b04, b0);
1864 c1 = gmul2n(b1,-1);
1865 c2 = gmul2n(gsub(gmul(b0,b2), gsqr(c1)),-1);
1866 c3 = gmul2n(gsub(gmul(b02,b3), gmul2n(gmul(c1,c2),1)),-1);
1867 c4 = gsub(gmul(b03,b4), gadd(gmul2n(gmul(c1,c3),1),gsqr(c2)));
1868 c5 = gsub(gmul(b04,b5), gmul2n(gmul(c2,c3),1));
1869 c6 = gsub(gmul(b05,b6), gsqr(c3));
1870 /* b0^5*H(x/b0) = (x^3+c1*x^2+c2*x+c3)^2+c4*x^2+c5*x+c6 */
1871 vc6 = myval(c6,p);
1872 if (vc6 == 2)
1873 {
1874 if (alpha)
1875 {
1876 condp = Dmin-16;
1877 Ip->type = "[IV] page 155";
1878 Ip->neron = cyclic(1);
1879 }
1880 else
1881 {
1882 condp = Dmin-6;
1883 Ip->type = "[III] page 155";
1884 Ip->neron = dicyclic(3,3);
1885 }
1886 }
1887 else
1888 {
1889 if (myval(c3,p) > 1) pari_err_BUG("genus2localred [c3]");
1890 mm = min3(3*myval(c4,p)-4, 3*myval(c5,p)-5, 3*vc6-6);
1891 if (alpha)
1892 {
1893 condp = Dmin-mm-16;
1894 Ip->type = stack_sprintf("[III*{%ld}] page 184", mm);
1895 Ip->neron = cyclic(1);
1896 }
1897 else
1898 {
1899 condp = Dmin-mm-6;
1900 Ip->type = stack_sprintf("[III{%ld}] page 184", mm);
1901 Ip->neron = (mm%3)? cyclic(9): dicyclic(3,3);
1902 }
1903 }
1904 }
1905 break;
1906 case 30:
1907 return equaliu(p,3)? quartic(polh, alpha, Dmin, Ip)
1908 : tame(polh, t60, alpha, Dmin, I, Ip);
1909 default: pari_err_BUG("genus2localred [red2]");
1910 }
1911 break;
1912 case 1:
1913 switch(t60+alpha)
1914 {
1915 case 12:
1916 condp = Dmin;
1917 Ip->type = "[VIII-1] page 156";
1918 Ip->neron = cyclic(1); break;
1919 case 13:
1920 condp = Dmin-10;
1921 Ip->type = "[IX-3] page 157";
1922 Ip->neron = cyclic(5); break;
1923 case 24:
1924 condp = Dmin-4;
1925 Ip->type = "[IX-1] page 157";
1926 Ip->neron = cyclic(5); break;
1927 case 25:
1928 condp = Dmin-14;
1929 Ip->type = "[VIII-3] page 157";
1930 Ip->neron = cyclic(1); break;
1931 case 36:
1932 condp = Dmin-8;
1933 Ip->type = "[VIII-2] page 157";
1934 Ip->neron = cyclic(1); break;
1935 case 15:
1936 condp = Dmin-1;
1937 Ip->type = "[VII] page 156";
1938 Ip->neron = cyclic(2); break;
1939 case 16:
1940 condp = Dmin-11;
1941 Ip->type = "[VII*] page 156";
1942 Ip->neron = cyclic(2); break;
1943 case 20:
1944 if (cmpis(p,3))
1945 {
1946 d = 6*val[6]-5*val[7]-2;
1947 if (d%6) pari_err_BUG("genus2localred [index]");
1948 dism = (d/6);
1949 }
1950 else
1951 {
1952 list = padicfactors(polh,p,Dmin-5);
1953 nb = lg(list);
1954 prod = pol_1(varn(polh));
1955 for(i = 1;i<nb;i++)
1956 {
1957 GEN c = gel(list,i);
1958 if (valp(gel(c,2)) && degpol(c)<= 2) prod = RgX_mul(prod,c);
1959 }
1960 if (degpol(prod) > 2) pari_err_BUG("genus2localred [padicfactors]");
1961 dism = valp(RgX_disc(prod)) - 1;
1962 }
1963 condp = Dmin-dism-3;
1964 Ip->type = stack_sprintf("[II-II*{%ld}] page 176", dism);
1965 Ip->neron = groupH(dism+1); break;
1966 case 21:
1967 vb6 = myval(RgX_coeff(polh,0),p);
1968 if (vb6<2) pari_err_BUG("genus2localred [red3]");
1969 condp = Dmin-14;
1970 Ip->type = "[IV*-II{0}] page 175";
1971 Ip->neron = cyclic(1); break;
1972 case 30:
1973 vb5 = myval(RgX_coeff(polh,1),p);
1974 if (vb5 == 2)
1975 {
1976 if (Ip->tt >= 5) pari_err_BUG("genus2localred [tt 6]");
1977 return tame(polh, t60, alpha, Dmin, I, Ip);
1978 }
1979 condp = Dmin-7;
1980 Ip->type = "[II*-III-(-1)] page 167";
1981 Ip->neron = cyclic(2); break;
1982 }
1983 break;
1984 case 2:
1985 if (ugcd(t60, 60) == 15) /* denom(theta) = 4 */
1986 {
1987 if (Ip->tt>4) pari_err_BUG("genus2localred [tt 5]");
1988 return tame(polh, t60, alpha, Dmin, I, Ip);
1989 }
1990 if (!equaliu(p,3) && ugcd(t60, 60) == 20) /* denom(theta) = 3 */
1991 return tame(polh, t60, alpha, Dmin, I, Ip);
1992 list = padicfactors(polh,p,Dmin-10*alpha);
1993 nb = lg(list); prod = pol_1(varn(polh));
1994 for(i = 1;i<nb;i++)
1995 {
1996 GEN c = gel(list,i);
1997 if (!valp(gel(c,2))) prod = RgX_mul(prod,c);
1998 }
1999 switch(degpol(prod))
2000 {
2001 GEN e0, e1, e2;
2002 case 0:
2003 dism = 0; break;
2004 case 1:
2005 e1 = gel(prod,3);
2006 dism = 2*valp(e1); break;
2007 case 2:
2008 e0 = gel(prod,2);
2009 e1 = gel(prod,3);
2010 e2 = gel(prod,4);
2011 dism = valp(gsub(gsqr(e1),gmul2n(gmul(e0,e2),2))); break;
2012 default:
2013 pari_err_BUG("genus2localred [padicfactors 2]");
2014 dism = 0;
2015 }
2016 switch(t60/5+alpha-4)
2017 {
2018 case 0:
2019 condp = Dmin-dism-1;
2020 Ip->type = stack_sprintf("[IV-II{%ld}] page 175", dism);
2021 Ip->neron = cyclic(3*dism+2); break;
2022 case 1:
2023 condp = Dmin-dism-10;
2024 Ip->type = stack_sprintf("[II*-II*{%ld}] page 176",dism);
2025 Ip->neron = groupH(dism+1); break;
2026 case 2: case 3:
2027 if (myval(RgX_coeff(polh,0),p) == 2)
2028 {
2029 if (Ip->tt>4) pari_err_BUG("genus2localred [tt 5]");
2030 return tame(polh, t60, alpha, Dmin, I, Ip);
2031 }
2032 dism++;
2033 indice = val[6]-(5*val[3]/2)-dism;
2034 condp = Dmin-dism-indice-2;
2035 Ip->type = stack_sprintf("[II{%ld-%ld}] page 182", dism,indice);
2036 Ip->neron = both_odd(dism,indice)? dicyclic(2,2*dism): cyclic(4*dism);
2037 break;
2038 case 4:
2039 condp = Dmin-dism-5;
2040 Ip->type = stack_sprintf("[IV*-II{%ld}] page 175",dism+1);
2041 Ip->neron = cyclic(3*dism+4); break;
2042 }
2043 break;
2044 case 3:
2045 if (!equaliu(p,3) || Ip->tt <= 4)
2046 return tame(polh, t60, alpha, Dmin, I, Ip);
2047 return labelm3(polh,t60,alpha,Dmin,I,Ip); /* p = 3 */
2048 default: pari_err_BUG("genus2localred [switch lambda]");
2049 }
2050 if (condp < 2 || condp > get_maxc(p))
2051 pari_err_BUG("genus2localred [conductor]");
2052 return condp;
2053 }
2054
2055 static long
chk_pol(GEN P)2056 chk_pol(GEN P) {
2057 switch(typ(P))
2058 {
2059 case t_INT: break;
2060 case t_POL: RgX_check_ZX(P,"genus2red"); return varn(P); break;
2061 default: pari_err_TYPE("genus2red", P);
2062 }
2063 return -1;
2064 }
2065
2066 /* P,Q are ZX, study Y^2 + Q(X) Y = P(X) */
2067 GEN
genus2red(GEN PQ,GEN p)2068 genus2red(GEN PQ, GEN p)
2069 {
2070 pari_sp av = avma;
2071 struct igusa I;
2072 GEN P, Q, D;
2073 GEN j22, j42, j2j6, a0,a1,a2,a3,a4,a5,a6, V,polr,facto,factp, vecmini, cond;
2074 long i, l, dd, vP,vQ;
2075
2076 PQ = Q_remove_denom(PQ, &D);
2077 if (typ(PQ) == t_VEC && lg(PQ) == 3)
2078 {
2079 P = gel(PQ,1);
2080 Q = gel(PQ,2);
2081 }
2082 else
2083 {
2084 P = PQ;
2085 Q = gen_0;
2086 }
2087
2088 vP = chk_pol(P);
2089 vQ = chk_pol(Q);
2090 if (vP < 0)
2091 {
2092 if (vQ < 0) pari_err_TYPE("genus2red",mkvec2(P,Q));
2093 P = scalarpol(P,vQ);
2094 }
2095 else if (vQ < 0) Q = scalarpol(Q,vP);
2096 if (p && typ(p) != t_INT) pari_err_TYPE("genus2red", p);
2097 if (D) P = ZX_Z_mul(P,D);
2098
2099 polr = ZX_add(ZX_sqr(Q), gmul2n(P,2)); /* ZX */
2100 switch(degpol(polr))
2101 {
2102 case 5: case 6: break;
2103 default: pari_err_DOMAIN("genus2red","genus","!=", gen_2,mkvec2(P,Q));
2104 }
2105
2106 RgX_to_03(polr, &a0,&a1,&a2,&a3);
2107 I.j10 = !signe(a0)? mulii(sqri(a1), ZX_disc(polr)): ZX_disc(polr);
2108 if (!signe(I.j10))
2109 pari_err_DOMAIN("genus2red","genus","<",gen_2,mkvec2(P,Q));
2110 I.j10 = gmul2n(I.j10, -12); /* t_INT */
2111
2112 if (p == NULL)
2113 {
2114 facto = absZ_factor(I.j10);
2115 factp = gel(facto,1);
2116 }
2117 else
2118 {
2119 factp = mkcol(p);
2120 facto = mkmat2(factp, mkcol(gen_1));
2121 }
2122 l = lg(factp);
2123 vecmini = cgetg(l, t_COL);
2124 for(i = 1; i<l; i++)
2125 {
2126 GEN l = gel(factp,i), pm;
2127 if (i == 1 && absequaliu(l, 2)) { gel(vecmini,1) = gen_0; continue; }
2128 gel(vecmini,i) = pm = polymini(polr, l);
2129 polr = ZX_Q_mul(gel(pm,1), powiu(l, gel(pm,2)[3]));
2130 }
2131 RgX_to_06(polr, &a0,&a1,&a2,&a3,&a4,&a5,&a6);
2132 I.j10 = !signe(a0)? mulii(sqri(a1), ZX_disc(polr)): ZX_disc(polr);
2133 I.j10 = gmul2n(I.j10,-12);
2134
2135 I.a0 = a0;
2136 I.A2 = apol2(a0,a1,a2);
2137 I.A3 = apol3(a0,a1,a2,a3);
2138 I.A5 = apol5(a0,a1,a2,a3,a4,a5);
2139 I.B2 = bpol2(a0,a1,a2,a3,a4);
2140
2141 I.j2 = igusaj2(a0,a1,a2,a3,a4,a5,a6);
2142 I.j4 = igusaj4(a0,a1,a2,a3,a4,a5,a6);
2143 I.i4 = gsub(gsqr(I.j2), gmulsg(24,I.j4));
2144 I.j6 = igusaj6(a0,a1,a2,a3,a4,a5,a6);
2145 j42 = gsqr(I.j4);
2146 j22 = gsqr(I.j2);
2147 j2j6 = gmul(I.j2,I.j6);
2148 I.j8 = gmul2n(gsub(j2j6,j42), -2);
2149 I.i12= gmul2n(gsub(gadd(gmul(j22,j42),gmulsg(36,gmul(j2j6,I.j4))),
2150 gadd(gadd(gmulsg(32,gmul(j42,I.j4)),gmul(j2j6,j22)),gmulsg(108,gsqr(I.j6)))),-2);
2151
2152 for(i = 1; i < l; i++)
2153 gcoeff(facto,i,2) = stoi(Q_pval(I.j10, gel(factp,i)));
2154 dd = ZX_pval(polr,gen_2) & (~1); /* = 2 floor(val/2) */
2155 polr = gmul2n(polr, -dd);
2156
2157 V = cgetg(l, t_VEC);
2158 for (i = 1; i < l; i++)
2159 {
2160 GEN q = gel(factp,i), red, N = NULL;
2161 struct igusa_p Ip;
2162 long f = genus2localred(&I, &Ip, q, gel(vecmini,i));
2163 gcoeff(facto,i,2) = stoi(f);
2164 if (Ip.tame) Ip.type = stack_strcat("(tame) ", Ip.type);
2165 if (f >= 0)
2166 N = zv_snf(Ip.neron);
2167 if (DEBUGLEVEL)
2168 {
2169 if (!p) err_printf("p = %Ps\n", q);
2170 err_printf("(potential) stable reduction: %Ps\n", Ip.stable);
2171 if (f >= 0) {
2172 err_printf("reduction at p: %s, %Ps", Ip.type, N);
2173 err_printf(", f = %ld\n", f);
2174 }
2175 }
2176 red = f >= 0? mkvec2(strtoGENstr(Ip.type), N): cgetg(1, t_VEC);
2177 gel(V, i) = mkvec3(q, Ip.stable, red);
2178 }
2179 if (p) V = gel(V,1);
2180 cond = factorback(facto);
2181 /* remove denominator 2 coming from f = -1 in genuslocalred(, p = 2) */
2182 if (typ(cond) != t_INT) cond = gel(cond,1);
2183 return gerepilecopy(av, mkvec4(cond, facto, polr, V));
2184 }
2185