Ticket #5307: trac_5307.patch

File trac_5307.patch, 8.8 KB (added by John Cremona, 13 years ago)
  • sage/schemes/elliptic_curves/ell_local_data.py

    # HG changeset patch
    # User John Cremona <john.cremona@gmail.com>
    # Date 1244824557 -3600
    # Node ID c99ab9fd4f6aa44fdbdd67faf36e6b39b1f4a1b2
    # Parent  6c7354ace3b6e0bf4c33ddb730531bd23a9fb92c
    [mq]: tate
    
    diff -r 6c7354ace3b6 -r c99ab9fd4f6a sage/schemes/elliptic_curves/ell_local_data.py
    a b  
    490490        t = verbose("Running Tate's algorithm with P = %s"%P, level=1)
    491491        F = OK.residue_field(P)
    492492        p = F.characteristic()
    493         if P.is_principal():
     493
     494        # In case P is not principal we mostly use a uniformiser which
     495        # is globally integral (with positive valuation at some other
     496        # primes); for this to work, it is essential that we can
     497        # reduce (mod P) elements of K which are not integral (but are
     498        # P-integral).  However, if the model is non-minimal and we
     499        # end up dividing a_i by pi^i then at that point we use a
     500        # uniformiser pi which has non-positive valuation at all other
     501        # primes, so that we can divide by it without losing
     502        # integrality at other primes.
     503
     504        principal_flag = P.is_principal()
     505        if principal_flag:
    494506            pi = P.gens_reduced()[0]
    495507            verbose("P is principal, generator pi = %s"%pi, t, 1)
    496508        else:
    497             pi = K.uniformizer(P, 'negative')
     509            pi = K.uniformizer(P, 'positive')
    498510            verbose("P is not principal, uniformizer pi = %s"%pi, t, 1)
     511        pi2 = pi*pi; pi3 = pi*pi2; pi4 = pi*pi3
    499512        prime = pi if K is QQ else P
    500513
    501514        pval = lambda x: x.valuation(prime)
     
    630643                break #return
    631644            if pval(b6) < 3:
    632645                ## Type IV
    633                 if _pquadroots(1, a3 / pi, -a6/(pi*pi)):
    634                     cp = 3
    635                 else:
    636                     cp = 1
     646                cp = 1
     647                a3t = preduce(a3/pi)
     648                a6t = preduce(a6/pi2)
     649                if _pquadroots(1, a3t, -a6t): cp = 3
    637650                KS = KodairaSymbol("IV")
    638651                fp = val_disc - 2
    639652                break #return
    640653
    641             # If our curve is none of these types, we change types so that p | a1, a2 and p^2 | a3, a4 and p^3 | a6
     654            # If our curve is none of these types, we change coords so that
     655            # p | a1, a2;  p^2 | a3, a4;  p^3 | a6
    642656            if p == 2:
    643                 s = proot(a2, 2)
    644                 t = pi*proot(a6/(pi*pi), 2)
     657                s = proot(a2, 2)        # so s^2=a2 (mod pi)
     658                t = pi*proot(a6/pi2, 2) # so t^2=a6 (mod pi^3)
    645659            elif p == 3:
    646                 s = a1
    647                 t = a3
     660                s = a1       # so a1'=2s+a1=3a1=0 (mod pi)
     661                t = a3       # so a3'=2t+a3=3a3=0 (mod pi^2)
    648662            else:
    649                 s = -a1*halfmodp
    650                 t = -a3*halfmodp
     663                s = -a1*halfmodp   # so a1'=2s+a1=0 (mod pi)
     664                t = -a3*halfmodp   # so a3'=2t+a3=0 (mod pi^2)
    651665            C = C.rst_transform(0, s, t)
    652666            (a1, a2, a3, a4, a6) = C.a_invariants()
    653667            (b2, b4, b6, b8) = C.b_invariants()
     
    665679            if min(pval(a1), pval(a2), pval(a3), pval(a4), pval(a6)) < 0:
    666680                raise RuntimeError, "Non-integral model after second transform!"
    667681
    668             # Analyze roots of the cubic T^3 + bT^2 + cT + d = 0, where b = a2/p, c = a4/p^2, d = a6/p^3
    669             b = a2/pi
    670             c = a4/(pi*pi)
    671             d = a6/(pi**3)
     682            # Analyze roots of the cubic T^3 + bT^2 + cT + d = 0 mod P, where
     683            # b = a2/p, c = a4/p^2, d = a6/p^3
     684            b = preduce(a2/pi)
     685            c = preduce(a4/pi2)
     686            d = preduce(a6/pi3)
    672687            bb = b*b
    673688            cc = c*c
    674689            bc = b*c
     
    703718                C = C.rst_transform(r, 0, 0)
    704719                (a1, a2, a3, a4, a6) = C.a_invariants()
    705720                (b2, b4, b6, b8) = C.b_invariants()
    706                 ix = 3; iy = 3; mx = pi*pi; my = pi*pi
     721                # The rest of this branch is just to compute cp, fp, KS.
     722                # We use pi to keep transforms integral.
     723                ix = 3; iy = 3; mx = pi2; my = mx
    707724                while True:
    708                     a2t = a2 / pi
    709                     a3t = a3 / my
    710                     a4t = a4 / (pi*mx)
    711                     a6t = a6 / (mx*my)
     725                    a2t = preduce(a2 / pi)
     726                    a3t = preduce(a3 / my)
     727                    a4t = preduce(a4 / (pi*mx))
     728                    a6t = preduce(a6 / (mx*my))
    712729                    if pdiv(a3t*a3t + 4*a6t):
    713730                        if p == 2:
    714731                            t = my*proot(a6t, 2)
     
    717734                        C = C.rst_transform(0, 0, t)
    718735                        (a1, a2, a3, a4, a6) = C.a_invariants()
    719736                        (b2, b4, b6, b8) = C.b_invariants()
    720                         my = my*pi
     737                        my *= pi
    721738                        iy += 1
    722                         a2t = a2/pi
    723                         a3t = a3/my
    724                         a4t = a4/(pi*mx)
    725                         a6t = a6/(mx*my)
     739                        a2t = preduce(a2 / pi)
     740                        a3t = preduce(a3/my)
     741                        a4t = preduce(a4/(pi*mx))
     742                        a6t = preduce(a6/(mx*my))
    726743                        if pdiv(a4t*a4t - 4*a6t*a2t):
    727744                            if p == 2:
    728745                                r = mx*proot(a6t*pinv(a2t), 2)
     
    731748                            C = C.rst_transform(r, 0, 0)
    732749                            (a1, a2, a3, a4, a6) = C.a_invariants()
    733750                            (b2, b4, b6, b8) = C.b_invariants()
    734                             mx = mx*pi
     751                            mx *= pi
    735752                            ix += 1 # and stay in loop
    736753                        else:
    737754                            if _pquadroots(a2t, a4t, a6t):
     
    767784                    raise RuntimeError, "Non-integral model after third transform!"
    768785                if pval(a2) < 2 or pval(a4) < 3 or pval(a6) < 4:
    769786                    raise RuntimeError, "Cubic after transform does not have a triple root at 0"
    770                 a3t = a3/(pi*pi)
    771                 a6t = a6/(pi**4)
     787                a3t = preduce(a3/pi2)
     788                a6t = preduce(a6/pi4)
    772789                # We test for Type IV*
    773790                if not pdiv(a3t*a3t + 4*a6t):
    774791                    cp = 3 if _pquadroots(1, a3t, -a6t) else 1
     
    776793                    fp = val_disc - 6
    777794                    break #return
    778795                # Now change coordinates so that p^3|a3, p^5|a6
    779                 t =        -pi*pi*proot(a6t, 2) if p==2 \
    780                       else  pi*pi*preduce(-a3t*halfmodp)
     796                if p==2:
     797                    t = -pi2*proot(a6t, 2)
     798                else: 
     799                    t = pi2*preduce(-a3t*halfmodp)
    781800                C = C.rst_transform(0, 0, t)
    782801                (a1, a2, a3, a4, a6) = C.a_invariants()
    783802                (b2, b4, b6, b8) = C.b_invariants()
     
    794813                    fp = val_disc - 8
    795814                    cp = 1
    796815                    break #return
    797                 a1 /= pi
    798                 a2 /= pi**2
    799                 a3 /= pi**3
    800                 a4 /= pi**4
    801                 a6 /= pi**6
     816                if not principal_flag:
     817                    pi = K.uniformizer(P, 'negative')
     818                pie = pi
     819                a1 /= pie
     820                pie *= pi # now pie=pi^2
     821                a2 /= pie
     822                pie *= pi # now pie=pi^3
     823                a3 /= pie
     824                pie *= pi # now pie=pi^4
     825                a4 /= pie
     826                pie *= pi # now pie=pi^5
     827                pie *= pi # now pie=pi^6
     828                a6 /= pie
    802829                verbose("Non-minimal equation, dividing out...\nNew model is %s"%([a1, a2, a3, a4, a6]), t, 1)
    803830        C = C._tidy_model()
    804831        return (C, p, val_disc, fp, KS, cp, split)
  • sage/schemes/elliptic_curves/ell_number_field.py

    diff -r 6c7354ace3b6 -r c99ab9fd4f6a sage/schemes/elliptic_curves/ell_number_field.py
    a b  
    848848            sage: E=EllipticCurve([0,0,0, 21796814856932765568243810*a - 134364590724198567128296995, 121774567239345229314269094644186997594*a - 750668847495706904791115375024037711300])
    849849            sage: E.conductor()
    850850            Fractional ideal (1)
     851
     852        An example which used to fail (see trac #5307)::
     853
     854            sage: K.<w>=NumberField(x^2+x+6)
     855            sage: E=EllipticCurve([w,-1,0,-w-6,0])
     856            sage: E.conductor()
     857            Fractional ideal (86304, w + 5898)
    851858        """
    852859        try:
    853860            return self._conductor