attachment:root_number.sage of days22/dokchitser

   1 def root_number(self, P=None):
   2     r"""
   3     Returns the global or local root number of this elliptic curve.
   4 
   5     The computations are based on:
   6     "Galois theory, elliptic curves, and root numbers" -- D. Rohrlich
   7     (P|3) "The local root number of elliptic curves with wild ramification" -- S. Kobayashi
   8     (P|2) "Root numbers of elliptic curves in residue characteristic 2" -- T. and V. Dokchitser
   9 
  10     INPUT::
  11 
  12          - `P` -- optional, default (None); if given, return the local
  13                   root number at the prime ideal `P`
  14 
  15     EXAMPLES::
  16     
  17         sage: K.<a>=NumberField(x^4+2)
  18         sage: E = EllipticCurve(K, [1, a, 0, 1+a, 0])
  19         sage: E.root_number()
  20         1
  21         sage: E.root_number(K.ideal(a+1))
  22         1
  23     """
  24     if P is None:
  25         # contribution from places at infinity
  26         s = self.base_field().signature()
  27         rn = (-1)^(s[0] + s[1])
  28         # contribution from places where E has bad reduction
  29         for d in self.local_data():
  30             rn = _root_number_local(self, d.prime()) * rn
  31         return rn
  32     else:
  33         return _root_number_local(self, P)
  34 
  35 
  36 def _root_number_local(E, P):
  37     r"""
  38     Returns the local root number of this elliptic curve at `P`.
  39 
  40     The computation for primes not dividing 2 or 3 is based on:
  41     "Galois theory, elliptic curves, and root numbers" -- D. Rohrlich
  42 
  43     INPUT::
  44 
  45          - `P` -- return the local root number at the prime ideal `P`
  46     """
  47     K = E.base_field()
  48     d = E.local_data(P)
  49     if d.has_good_reduction() or d.has_nonsplit_multiplicative_reduction():
  50         return 1
  51     elif d.has_split_multiplicative_reduction():
  52         return -1
  53     else: # additive reduction
  54         j = E.j_invariant()
  55         jv = j.valuation(P) if not K is QQ else valuation(j, P.gen())
  56         if jv < 0:
  57             # potential multiplicative reduction
  58             if P.base_ring() is ZZ:
  59                 return hilbert_symbol(-1, -E.c6(), P.gen())
  60             else:
  61                 return generalized_hilbert_symbol(-1, -E.c6(), P)
  62                 # return hilbert_symbol_magma(-1, -E.c6(), P)
  63         else:
  64             # potential good reduction
  65             p = P.gen() if K is QQ else P.smallest_integer()
  66             if p==2: # prime | 2
  67                 return _root_number_local_2(E, P)
  68             if p==3: # prime | 3
  69                 return _root_number_local_3(E, P)
  70 
  71             if K is QQ:
  72                 f = 1
  73                 e = 12 / gcd(12, valuation(E.discriminant(), p))
  74             else:
  75                 f = P.residue_class_degree()
  76                 e = 12 / gcd(12, E.discriminant().valuation(P))
  77 
  78             if f%2==0 or e==1:
  79                 eps = 1
  80             elif e==2 or e==6:
  81                 eps = legendre_symbol(-1, p)
  82             elif e==3:
  83                 eps = legendre_symbol(-3, p)
  84             elif e==4:
  85                 eps = legendre_symbol(-2, p)
  86             return eps
  87 
  88 
  89 def _root_number_local_3(E, P):
  90     r"""
  91     Returns the local root number of this elliptic curve at `P`
  92     where `P` is a prime ideal dividing 3 at which this elliptic
  93     curve has (bad) potential good reduction.
  94 
  95     The computation for primes dividing 3 is based on:
  96     "The local root number of elliptic curves with wild ramification" -- S. Kobayashi
  97 
  98     INPUT::
  99 
 100          - `P` -- return the local root number at the prime ideal `P`
 101     """
 102     def _val(a, P):
 103         if P.base_ring() is ZZ:
 104             return valuation(a, P.gen())
 105         return a.valuation(P)
 106     def quadr_residue_symbol(a):
 107         return Kr(a).is_square() and 1 or -1
 108     def quadr_symbol(a):
 109         av =  _val(a, P)
 110         if av%2==1:
 111             return -1
 112         if P.base_ring() is ZZ:
 113             pi = P.gen()
 114         else:
 115             pi = P.number_field().uniformizer(P)
 116         b = a / pi^av
 117         return quadr_residue_symbol(b)
 118     # assume potential good reduction
 119     Es = E.local_data(P).minimal_model(tidy=False)
 120     Es = Es.short_weierstrass_model(complete_cube=False)
 121 
 122     disc = Es.discriminant()
 123     vD = _val(disc, P)
 124     Kr = P.residue_field()
 125 
 126     ks = str(Es.kodaira_symbol(P))
 127     if ks=='I0' or ks=='I0*':
 128         assert vD%2==0, "error -- valuation of discriminant should be even"
 129         return quadr_residue_symbol(-1)^(vD/2)
 130     elif ks=='III' or ks=='III*':
 131         return quadr_residue_symbol(-2)
 132     else: 
 133         assert ks=='II' or ks=='II*' or ks=='IV' or ks=='IV*', "error -- unexpected kodeira symbol"
 134         d = quadr_symbol(disc)
 135         c = Es.a6()
 136         cv = _val(c, P)
 137         assert cv%3!=0, "error -- 3 should not divide the valuation of the constant term anymore"
 138         if P.base_ring() is ZZ:
 139             hs = hilbert_symbol(disc, c, P.gen())
 140         else:
 141             hs = generalized_hilbert_symbol(disc, c, P)
 142         qr1 = quadr_residue_symbol(cv)^vD
 143         qr2 = quadr_residue_symbol(-1)^(vD*(vD-1)/2)
 144         return d*hs*qr1*qr2
 145 
 146 
 147 def _root_number_local_2(E, P):
 148     r"""
 149     Returns the local root number of this elliptic curve at `P`
 150     where `P` is a prime ideal dividing 2 at which this elliptic
 151     curve has (bad) potential good reduction.
 152 
 153     The computation for primes dividing 2 is based on:
 154     "Root numbers of elliptic curves in residue characteristic 2" -- T. and V. Dokchitser
 155 
 156     INPUT::
 157 
 158          - `P` -- return the local root number at the prime ideal `P`
 159     """
 160     K = E.base_field()
 161 
 162     if K is QQ:
 163         rn = -1
 164     else:
 165         rn = (-1)^(P.residue_class_degree() * P.ramification_index())
 166 
 167     t = polygen(K)
 168     assert E.a1() == 0
 169     assert E.a3() == 0
 170     f = t^3+E.a2()*t^2+E.a4()*t+E.a6()
 171 
 172     d = 1
 173     roots = []
 174     for g,_ in f.factor():
 175         if g.degree() == 1:
 176             roots += g.roots(multiplicities=False)
 177         if g.degree() == 2:
 178             d = 2
 179             K2.<a2> = K.extension(g)
 180             # roots += g.change_ring(K2).roots(multiplicities=False)
 181         if g.degree() == 3:
 182             d = 3
 183             K3.<a3> = K.extension(g)
 184             for h,_ in g.change_ring(K3).factor():
 185                 if h.degree() == 1:
 186                     # roots += h.roots(multiplicities=False)
 187                     pass
 188                 if h.degree() == 2:
 189                     d = 6
 190                     K6.<a6> = K3.extension(h)
 191                     # roots += h.change_ring(K6).roots(multiplicities=False)
 192 
 193     # print "d = ", d
 194 
 195     def _val(a, P):
 196         if P.base_ring() is ZZ:
 197             return valuation(a, P.gen())
 198         return a.valuation(P)
 199 
 200     def _change_ring(E, K, toK):
 201         return EllipticCurve(K, [0, toK(E.a2()), 0, toK(E.a4()), toK(E.a6())])
 202 
 203     def _translate_w_2_isogeny(E, r):
 204         E1 = E.change_weierstrass_model(1,r,0,0)
 205         E2 = EllipticCurve(E.base_field(), [0,-2*E1.a2(),0,E1.a2()^2-4*E1.a4(),0])
 206         # print E.is_isogenous(E2)
 207         return E2
 208 
 209     def _compute_H(E, P, ord2=True, ord3=False):
 210         # print E
 211         t = E.tamagawa_number(P)
 212         # print t
 213         if ord2:
 214             cv = valuation(t,2) + _u_of_E(E, P)
 215             H2 = 1 if cv%2==0 else -1
 216         if ord3:
 217             cv = valuation(t,3)
 218             H3 = 1 if cv%2==0 else -1
 219         if ord2 and ord3: return (H2, H3)
 220         if ord2: return H2
 221         return H3
 222 
 223     def _u_of_E(E, P):
 224         disc = E.discriminant()
 225         Emin = E.local_data(P).minimal_model(tidy=False)
 226         discmin = Emin.discriminant()
 227         assert _val(discmin, P) == E.local_data(P)._val_disc
 228         # print "quot = ", disc/discmin
 229         # return (disc/discmin).norm()^(1/12)
 230         # print "quot = ", (disc/discmin).valuation(P) / 12
 231         # print "rd = ", P.residue_class_degree()
 232         # if not E.base_field() is QQ:
 233             # print P.ramification_index()
 234         if E.base_field() is QQ:
 235             return Integer(_val(disc/discmin, P) / 12)
 236         else:
 237             # print P, disc/discmin, P.residue_class_degree()
 238             return Integer(P.residue_class_degree() * _val(disc/discmin, P) / 12)
 239             # return Integer(_val(disc/discmin, P) / 12)
 240 
 241     if d==1:
 242         # print roots
 243         H = _compute_H(E, P)
 244         E1 = _translate_w_2_isogeny(E, roots[0])
 245         H1 = _compute_H(E1, P)
 246         E2 = _translate_w_2_isogeny(E, roots[1])
 247         H2 = _compute_H(E2, P)
 248         E3 = _translate_w_2_isogeny(E, roots[2])
 249         H3 = _compute_H(E3, P)
 250         return rn * H * H1 * H2 * H3
 251     elif d==2:
 252         H = _compute_H(E, P)
 253         # 2-isogenous curve over base field
 254         # print roots[0]
 255         E1 = _translate_w_2_isogeny(E, roots[0]) # root[0] is the root over the base field
 256         # print E1
 257         H1 = _compute_H(E1, P)
 258 
 259         # a2 in K2 is root not in base field
 260         #TODO: over QQ not necessary
 261         K2b.<b2> = K2.absolute_field()
 262         from2, to2 = K2b.structure()
 263         # Unfairly picking one prime factor
 264         # Pfactors = K2b.ideal([ to2(pgen) for pgen in P.gens() ]).factor()
 265         # print Pfactors
 266         # P2b = K2b.ideal([ to2(pgen) for pgen in P.gens() ]).factor()[0][0]
 267         for P2b,_ in K2b.ideal([ to2(pgen) for pgen in P.gens() ]).factor():
 268             # print P2b
 269             E2b = _change_ring(E1, K2b, to2)
 270             H2 = _compute_H(E2b, P2b)
 271             # print "H2 = ", H2
 272 
 273             E2b = _change_ring(E, K2b, to2)
 274             E3b = _translate_w_2_isogeny(E2b, to2(a2))
 275             assert E3b.a6() == 0
 276             H3 = _compute_H(E3b, P2b)
 277             rn = rn * H2 * H3
 278 
 279         # print H, H1, H2, H3
 280 
 281         return rn * H * H1
 282     elif d==3:
 283         # a3 in K3
 284         K3b.<b3> = K3.absolute_field()
 285         from3, to3 = K3b.structure()
 286         # Unfairly picking one prime factor
 287         # print K3b.ideal([ to3(pgen) for pgen in P.gens() ]).factor()
 288         for P3b,_ in K3b.ideal([ to3(pgen) for pgen in P.gens() ]).factor():
 289             # P3b = K3b.ideal([ to3(pgen) for pgen in P.gens() ]).factor()[0][0]
 290             E3b = _change_ring(E, K3b, to3)
 291             H = _compute_H(E3b, P3b)
 292             # print "H = ", H
 293             E3bt = _translate_w_2_isogeny(E3b, to3(a3))
 294             H2 = _compute_H(E3bt, P3b)
 295             # print H, H2
 296             rn = rn * H * H2
 297         # print "H2 = ", H2
 298         return rn
 299     elif d==6:
 300         # a3 in K3
 301         K3b.<b3> = K3.absolute_field()
 302         from3, to3 = K3b.structure()
 303         # Unfairly picking one prime factor
 304         P3b = K3b.ideal([ to3(pgen) for pgen in P.gens() ]).factor()[0][0]
 305         E3b = _change_ring(E, K3b, to3)
 306         H = _compute_H(E3b, P3b)
 307 
 308         E3bt = _translate_w_2_isogeny(E3b, to3(a3))
 309         H2 = _compute_H(E3bt, P3b)
 310 
 311         # a6 in K6
 312         K6b.<b6> = K6.absolute_field()
 313         from6, to6 = K6b.structure()
 314         # Unfairly picking one prime factor
 315         P6b = K6b.ideal([ to6(pgen) for pgen in P.gens() ]).factor()[0][0]
 316         E6b = _change_ring(E, K6b, to6)
 317         H3, H33 = _compute_H(E6b, P6b, ord2=True, ord3=True)
 318 
 319         E6bt = _translate_w_2_isogeny(E6b, to6(a3))
 320         H4 = _compute_H(E6bt, P6b)
 321 
 322         x = polygen(K)
 323         Kd.<d> = K.extension(x^2 - E.discriminant())
 324         KDb.<bd> = Kd.absolute_field()
 325         fromD, toD = KDb.structure()
 326         # Unfairly picking one prime factor
 327         PDb = KDb.ideal([ toD(pgen) for pgen in P.gens() ]).factor()[0][0]
 328         EDb = _change_ring(E, KDb, toD)
 329         HD = _compute_H(EDb, PDb, ord2=False, ord3=True)
 330 
 331         return rn * H * H2 * H3 * H33 * H4 * HD
 332 
 333     raise NotImplementedError("additive reduction over a prime dividing 2")