1 #include "sympow.h"
2 #define DEBUG (FALSE || GLOBAL_DEBUG)
3 
4 #define tlc tame_local_conductor
5 #define wlc wild_local_conductor
6 
tame_local_conductor(int type,int sympow)7 int tame_local_conductor(int type,int sympow)
8 {if (sympow<=0) return 0;
9  switch(type)
10  {case 0: errorit("case 0 in tame_local_conductor"); /* good */
11   case 1: return sympow; /* multiplicative */
12   case 29: case 30: case 31: /* potentially multiplicative */
13    if (sympow&1) return sympow+1; else return sympow;
14   case 2: case 16: case 17: /* C2 */
15    if (sympow&1) return sympow+1; else return 0;
16   case 3: case 14: case 18:  /* C3 */
17   switch(sympow%3)
18   {case 0: return (sympow*2)/3; case 1: return sympow+1-(sympow-1)/3;
19    case 2: return sympow+1-(sympow+1)/3;}
20   case 4: case 22: case 23: /* C4 */
21   switch(sympow&3)
22   {case 1: case 3: return sympow+1;
23    case 0: return sympow/2; case 2: return 1+sympow/2;}
24   case 6: case 15: case 19: case 20: case 21: /* C6 */
25   switch(sympow%6)
26   {case 1: case 3: case 5: return sympow+1;
27    case 0: return (sympow*2)/3; case 4: return sympow+1-(sympow-1)/3;
28    case 2: return sympow+1-(sympow+1)/3;}
29   case 8: case 9: case 10: /* Q8 */
30   switch(sympow&3)
31   {case 1: case 3: return (sympow+1);
32    case 0: return (3*sympow)/4; case 2: return (3*sympow+6)/4;}
33   case 12: case 13: /* C3x|C4 */
34   switch(sympow%12)
35   {case 1: case 3: case 5: case 7: case 9: case 11: return sympow+1;
36    case 0: return (5*sympow)/6; case 2: return (5*sympow+8)/6;
37    case 4: return (5*sympow+4)/6; case 6: return (5*sympow)/6+1;
38    case 8: return (5*sympow+2)/6; case 10: return (5*sympow+10)/6;}
39   case 24: case 25: case 26: case 27: /* SL2F3 */
40   switch(sympow%12)
41   {case 1: case 3: case 5: case 7: case 9: case 11: return sympow+1;
42    case 0: return (11*sympow)/12; case 2: return (11*sympow+14)/12;
43    case 4: return (11*sympow+16)/12; case 6: return (11*sympow+6)/12;
44    case 8: return (11*sympow+8)/12; case 10: return (11*sympow+22)/12;}
45   default: printf("%i ",type); errorit("Bad type in tame_local_conductor");}
46  return -1;}
47 
wild_local_conductor(int type,int m)48 int wild_local_conductor(int type,int m)
49 {if (m<=0) return 0;
50  switch(type)
51  {case 0: case 1: case 2: case 3: case 4: case 6: case 31: return 0; /* p>5 */
52   case 8: return tlc(8,m)+tlc(2,m)/2;
53   case 9: return tlc(8,m)+tlc(2,m);
54   case 10: return tlc(8,m)+tlc(4,m)+tlc(2,m);
55   case 12: return tlc(3,m)/2;
56   case 13: return 3*tlc(3,m)/2;
57   case 14: case 15: case 18: case 19: return tlc(3,m);
58   case 16: case 20: return tlc(2,m);
59   case 17: case 21: return 2*tlc(2,m);
60   case 22: case 23: return 2*tlc(4,m)+tlc(2,m);
61   case 24: return (2*tlc(8,m)+tlc(2,m))/6;
62   case 25: return (tlc(8,m)+2*tlc(2,m))/3;
63   case 26: return (tlc(8,m)+5*tlc(2,m))/3;
64   case 27: return (10*tlc(8,m)+5*tlc(2,m))/6;
65   case 29: return 2*tlc(2,m);
66   case 30: return tlc(2,m);
67   default: printf("%i ",type); errorit("Bad type in wild_local_conductor");}
68  return -1;}
69 
badprimetype_output(int x,llint p)70 void badprimetype_output(int x,llint p)
71 {switch(x)
72  {case 0: printf("C1 ADDITIVE REDUCTION\n"); break;
73   case 1: printf("C1 MULTIPLICATIVE REDUCTION\n"); break;
74   case 2: case 3: case 4: case 6:
75   {if ((p%x)==1) printf("C%i abelian\n",x);
76    else printf("C%i nonabelian\n",x); break;}
77   case 8: printf("Q8 v2(N)=5\n"); break;
78   case 9: printf("Q8 v2(N)=6\n"); break;
79   case 10: printf("Q8 v2(N)=8\n"); break;
80   case 12: printf("C3x|C4 v3(N)=3\n"); break;
81   case 13: printf("C3x|C4 v3(N)=5\n"); break;
82   case 14: printf("C3 p=3 abelian\n"); break;
83   case 15: printf("C6 p=3 abelian\n"); break;
84   case 16: printf("C2 p=2 v2(N)=4\n"); break;
85   case 17: printf("C2 p=2 v2(N)=6\n"); break;
86   case 18: printf("C3 p=3 nonabelian\n"); break;
87   case 19: printf("C6 p=3 nonabelian\n"); break;
88   case 20: printf("C6 p=2 v2(N)=4\n"); break;
89   case 21: printf("C6 p=2 v2(N)=6\n"); break;
90   case 22: printf("C4 p=2 abelian\n"); break;
91   case 23: printf("C4 p=2 nonabelian\n"); break;
92   case 24: printf("SL2F3 v2(N)=3\n"); break;
93   case 25: printf("SL2F3 v2(N)=4\n"); break;
94   case 26: printf("SL2F3 v2(N)=6\n"); break;
95   case 27: printf("SL2F3 v2(N)=7\n"); break;
96   case 29: printf("C2 p=2 MULTIPLICATIVE v2(N)=6\n"); break;
97   case 30: printf("C2 p=2 MULTIPLICATIVE v2(N)=4\n"); break;
98   case 31: printf("C2 MULTIPLICATIVE\n"); break;
99   default: errorit("Unknown badprime type");}}
100 
bpt_at(double p)101 static int bpt_at(double p)/* p>3.0 */
102 {int v4,v6,vd; if (p<=3.0) errorit("prime too small in bpt_at");
103  if (!QD_is_divisible(Edisc,p)) return 0;
104  if (!QD_is_divisible(Ec4,p)) return 1;
105  v4=QD_valuation(Ec4,p); v6=QD_valuation(Ec6,p); vd=QD_valuation(Edisc,p);
106  if ((v4==2) && (v6==3) && (vd>6)) return 31;
107  if (vd==6) return 2; if ((vd==2) || (vd==10)) return 6;
108  if ((vd==3) || (vd==9)) return 4; if ((vd==4) || (vd==8)) return 3;
109  errorit("bad valuation in bpt_ap"); return -1;}
110 
ecfp3(int c4t,int c6t,int Dt)111 static int ecfp3(int c4t,int c6t,int Dt)
112 {int c4r,c6r; if (DEBUG) printf("ecfp3 %i %i %i\n",c4t,c6t,Dt);
113  if ((Dt%3)!=0) return 0; if ((c6t%27)!=0) return 1;
114  c4r=(c4t/9)%3; c6r=(c6t/27)%9;
115  if (c4r==0)
116  {if (c6r==0) {if ((c4t%81)==0) return 5; return 4;}
117   switch (c6r)
118   {case 1: case 2: case 7: case 8: return 3;
119    case 3: case 6: return 5; case 4: case 5: return 2;}}
120  if (c4r==1)
121  {switch(c6r)
122   {case 0: return 2; case 1: case 3: case 6: case 8: return 3;
123    case 2: case 4: case 5: case 7: return 4;}}
124  if (c4r==2)
125  {switch(c6r)
126   {case 0: return 2; case 2: case 7: return 2;
127   case 1: case 3: case 4: case 5: case 6: case 8: return 3;}}
128  return 0;
129 }
130 
bpt_at3()131 static int bpt_at3()
132 {int v4,v6,vd,t,u,c4,c6,D,c4t,c6t,Dt;
133  if (DEBUG) printf("bpt_at3\n");
134  if (!QD_is_divisible(Edisc,3.0)) return 0;
135  if (!QD_is_divisible(Ec4,3.0)) {fp3=1; return 1;}
136  v4=QD_valuation(Ec4,3.0); v6=QD_valuation(Ec6,3.0);
137  vd=QD_valuation(Edisc,3.0); c4=QD_modi(Ec4,387420489);
138  c6=QD_modi(Ec6,387420489); D=QD_modi(Edisc,387420489);
139  if (DEBUG) printf("%i %i %i\n",c4,c6,D);
140  if ((v4==2) && (v6==3) && (v6!=5) && (vd>6)) {fp3=2; return 31;}
141  if ((v4>=2) && ((v6>=3) && (v6!=5)) && (vd==6)) {fp3=2; return 2;}
142  if ((v4<2) || (v6<3) || (v6==5) || (vd<6)) {c4t=c4; c6t=c6; Dt=D;}
143  else {c4t=c4/9; c6t=(c6-387420489)/(-27); Dt=D/729;}
144  fp3=ecfp3(c4t,c6t,Dt); ASSERT(fp3>=2);
145  if (fp3==2) return 4; if (fp3==3) return 12; if (fp3==5) return 13;
146  if (vd&3) t=15; else t=14;
147  if (((c6t%243)==0) && ((c4t%81)==54)) t+=4;
148  if ((c4t%27)==9) {u=c6t%243; if ((u==54) || (u==189)) t+=4;}
149  return t;
150 }
151 
ecfp2(int c4,int c6)152 static int ecfp2(int c4,int c6)
153 {int c4r,c6r; c4r=c4&15; c6r=c6&15;
154  if ((c4r&3)==1)
155  {switch(c6r&3)
156   {case 0: return 6; case 1: return 4; case 3: return 3;
157   case 2:
158   {if (c4r==1) switch(c6r)
159    {case 2: return 4; case 6: case 10: return 5; case 14: return 3;}
160    if (c4r==5) switch(c6r)
161    {case 2: return 3; case 6: return 2; case 10: case 14: return 4;}
162    if (c4r==9) switch(c6r)
163    {case 2: case 14: return 5; case 6: return 3; case 10: return 4;}
164    if (c4r==13) switch(c6r)
165    {case 2: case 6: return 4; case 10: return 3; case 14: return 2;}}}}
166  if ((c4r&3)==2)
167  {switch(c6r&3)
168   {case 1: return 3; case 2: return 6; case 3: return 4;
169    case 0: if (c6r&7) return 7; else return 8;}}
170  if ((c4r&3)==3)
171  {switch(c6r&3)
172   {case 0: return 5; case 1: return 2; case 2: return 7; case 3: return 4;}}
173  if (((c4r&3)==0) && ((c6r&3)!=0))
174  {switch(c6r&3) {case 1: return 2; case 3: return 4;}}
175  c4r=(c4/4)&1; c6r=(c6/4)&3;
176  if (c6r==3) return 4; if (c4r) return 3; return 2;
177 }
178 
bpt_at2()179 static int bpt_at2()
180 {int v4,v6,vd,e,c4,c6,D,c4t,c6t,Dt,t1=FALSE;
181  if (!QD_is_divisible(Edisc,2.0)) {fp2=0; return 0;}
182  if (!QD_is_divisible(Ec4,2.0)) {fp2=1; return 1;}
183  v6=QD_valuation(Ec6,2.0); if (v6==3) {fp2=0; return 0;}
184  vd=QD_valuation(Edisc,2.0); e=minimum(v6/3,vd/6); v4=QD_valuation(Ec4,2.0);
185  c4=QD_modi(Ec4,1<<27); c6=QD_modi(Ec6,1<<27); D=QD_modi(Edisc,1<<27);
186  v4-=2*e; v6-=3*e; if ((v4==2) || (v4==3) || (v6==4)) e--;
187  c4t=c4>>(2*e); c6t=c6>>(3*e); Dt=D>>(6*e);
188  if ((c6t&31)!=0)
189  {if ((Dt&1)==0) {if (e&1) {fp2=6; return 29;} fp2=4; return 30;}
190   else {if (e&1) {fp2=6; return 17;} fp2=4; return 16;}}
191  fp2=ecfp2(c4t/16,c6t/32); if (fp2==4) {fp2=ecfp2(c4t/16,-c6t/32); t1=TRUE;}
192  ASSERT(fp2>=2); ASSERT(fp2!=4);
193  if (fp2==2)
194  {if ((e&1)==0) {if (t1) {fp2=4; return 20;} return 3;} fp2=6; return 21;}
195  if (fp2==3)
196  {if ((e&1)==0) {if (t1) {fp2=4; return 25;} return 24;} fp2=6; return 26;}
197  if (fp2==5) {if ((e&1)==0) return 8; fp2=6; return 9;}
198  if (fp2==6) {if ((e&1)==1) {fp2=5; return 8;} return 9;}
199  if (fp2==7) return 27;
200  if (fp2==8)
201  {if ((c6t&511)==0) return 10; if ((c4t&127)==32) return 23; return 22;}
202  return -1;
203 }
204 
CMUL(llint p)205 static void CMUL(llint p) /* multiply COND0 by p and check if overflow */
206 {if (((COND0*p)/p)!=COND0) errorit("Curve conductor is too large"); COND0*=p;}
207 
do_badprime(llint p)208 int do_badprime(llint p) /* E is a minimal model */
209 {int f=0,i;
210  if (p>3) {f=bpt_at((double) p); if (f==1) CMUL(p); else CMUL(p*p);}
211  if (p==3) {f=bpt_at3(); for (i=0;i<fp3;i++) CMUL(3);}
212  if (p==2) {f=bpt_at2(); for (i=0;i<fp2;i++) CMUL(2);} return f;}
213 
lkup(int m)214 static int lkup(int m)
215 {if ((m>=10) && (m<=20)) return 0; m=m%10; if (m>=4) return 0; return m;}
216 
compute_conductor(int m)217 void compute_conductor(int m)
218 {int i,T,W; llint P; QD A; QDpoly v; char STR[4][4]={"th","st","nd","rd"};
219  QD_copy(wmax,QD_one,COND[m]);
220  for (i=0;badprimes[i]!=0;i++)
221  {T=tlc(badprimetype[i],m); P=badprimes[i]; W=wlc(badprimetype[i],m);
222   if (VERBOSE) printf("sp %i: Conductor at %lli is %i+%i, root number is %i\n",
223 		      m,P,T,W,local_rootno(i,badprimes[i],m));
224   QD_copy(wmax,QD_zero,A); A[0]=(double) P; QD_powi(wmax,A,T+W,A);
225   QD_mul(wmax,COND[m],A,COND[m]);
226   if (VERBOSE)
227   {printf("sp %i: Euler factor at %lli is ",m,P);
228    euler_factor(P,m,&v); QDpoly_intout(v); delQDpoly(&v);}}
229  if (VERBOSE)
230  {printf("%i%s sym power conductor is ",m,STR[lkup(m)]);
231   QD_intout_noplus(COND[m]);
232   printf(", global root number is %i\n",global_rootno(m));}}
233 
prim_part(int CM)234 static int prim_part(int CM)
235 {if ((CM==-27) || (CM== -12)) return 3; if (CM==-28) return 7;
236  if ((CM==-16) || (CM==-4) || (CM==-8)) return 2; return -CM;}
237 
f(int CM)238 static int f(int CM)
239 {if ((CM==-16) || (CM==-4)) return 2; if (CM==-8) return 3; return 1;}
240 
compute_conductor_hecke(int m)241 void compute_conductor_hecke(int m)
242 {int i,T,W,u; llint P=0; QD A; QDpoly v; char STR[4][4]={"th","st","nd","rd"};
243  if (!CM_CASE) errorit("not CM case in compute_conductor_hecke\n");
244  QD_copy(wmax,QD_one,COND[m]);
245  for (i=0;badprimes[i]!=0;i++)
246  {T=tlc(badprimetype[i],m)-tlc(badprimetype[i],m-2);
247   W=wlc(badprimetype[i],m)-wlc(badprimetype[i],m-2); P=badprimes[i];
248   if ((m&3)==0) {u=prim_part(CM_CASE); if (P==u) {T+=1; W+=f(CM_CASE)-1;}}
249   if ((m&3)==2) {u=prim_part(CM_CASE); if (P==u) {T-=1; W-=f(CM_CASE)-1;}}
250   if (VERBOSE) printf("sp %i: Conductor at %lli is %i+%i\n",m,P,T,W);
251   QD_copy(wmax,QD_zero,A); A[0]=(double) P; QD_powi(wmax,A,T+W,A);
252   QD_mul(wmax,COND[m],A,COND[m]);
253   if (VERBOSE)
254   {printf("sp %i: Euler factor at %lli is ",m,P);
255    euler_factor(P,m,&v); QDpoly_intout(v); delQDpoly(&v);}}
256  if (VERBOSE)
257  {printf("%i%s sym power conductor is ",m,STR[lkup(m)]);
258   QD_intout_noplus(COND[m]); printf("\n");}}
259 
get_tame_conductor(llint p,int m)260 int get_tame_conductor(llint p,int m)
261 {int i; if ((COND0%p)!=0) return 0;
262  if (HECKE) errorit("Tame conductor for Hecke?");
263  for (i=0;badprimes[i]!=p;i++); return(tlc(badprimetype[i],m));}
264 
get_wild_conductor(llint p,int m)265 int get_wild_conductor(llint p,int m)
266 {int i; if ((COND0%p)!=0) return 0;
267  if (HECKE) errorit("Wild conductor for Hecke?");
268  for (i=0;badprimes[i]!=p;i++); return(wlc(badprimetype[i],m));}
269 
get_conductor(llint p,int m)270 int get_conductor(llint p,int m)
271 {int T,W,i,u; if ((COND0%p)!=0) return 0;
272  for (i=0;badprimes[i]!=p;i++);
273  if (!HECKE) return(tlc(badprimetype[i],m)+wlc(badprimetype[i],m));
274  T=tlc(badprimetype[i],m)-tlc(badprimetype[i],m-2);
275  W=wlc(badprimetype[i],m)-wlc(badprimetype[i],m-2);
276  if ((m&3)==0) {u=prim_part(CM_CASE); if (p==u) T+=f(CM_CASE);}
277  if ((m&3)==2) {u=prim_part(CM_CASE); if (p==u) T-=f(CM_CASE);} return T+W;}
278