Code coverage tests

This page documents the degree to which the PARI/GP source code is tested by our public test suite, distributed with the source distribution in directory src/test/. This is measured by the gcov utility; we then process gcov output using the lcov frond-end.

We test a few variants depending on Configure flags on the pari.math.u-bordeaux.fr machine (x86_64 architecture), and agregate them in the final report:

The target is 90% coverage for all mathematical modules (given that branches depending on DEBUGLEVEL or DEBUGMEM are not covered). This script is run to produce the results below.

LCOV - code coverage report
Current view: top level - modules - genus2red.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.11.0 lcov report (development 22860-5579deb0b) Lines: 1207 1314 91.9 %
Date: 2018-07-18 05:36:42 Functions: 54 54 100.0 %
Legend: Lines: hit not hit

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

Generated by: LCOV version 1.13