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