Ticket #11975: trac_11975-part4.patch

File trac_11975-part4.patch, 9.1 KB (added by was, 9 years ago)

part 4, addressing the major issue John Cremona found in equivalence checking.

  • sage/schemes/elliptic_curves/chow_heegner.py

    # HG changeset patch
    # User William Stein <wstein@gmail.com>
    # Date 1327089751 28800
    # Node ID 3c2d896e6b69aeaa6bb1aff3b3429ea806130f68
    # Parent  8df8e78d666e362282d4f0604128662f3c44e9df
    trac 11975 -- Chow-Heegner points; part 4
    
    diff --git a/sage/schemes/elliptic_curves/chow_heegner.py b/sage/schemes/elliptic_curves/chow_heegner.py
    a b  
    453453
    454454        sage: z1 = ComplexField(200)(5,7)
    455455        sage: all(is_gamma0N_equivalent(z1, g.acton(z1), 15, 150) for g in Gamma0(15).gens())
    456         True       
     456        True
     457
     458    These two points have the property that z1 and z2 are equivalent
     459    modulo SL2(Z), and 5*z1 and 5*z2 are *also* equivalent modulo
     460    SL2(Z), but z1 and z2 are not equivalent modulo Gamma0(N).  Note
     461    that both points are equivalent to I, which has extra
     462    automorphisms.  This is interesting because generators for the
     463    modular function field j(z) and j(5*z) agree on z1 and z2, but
     464    the points are inequivalent (this illustrates that this model
     465    for the modular curve X0(5) is singular)::
     466   
     467        sage: from sage.schemes.elliptic_curves.chow_heegner import is_gamma0N_equivalent, is_sl2z_equivalent
     468        sage: z1 = CDF(-2,1)/5                     
     469        sage: z2 = CDF(2,1)/5                       
     470        sage: is_gamma0N_equivalent(z1,z2,5, 30)
     471        False
     472        sage: is_sl2z_equivalent(z1, z2, 30)
     473        True
     474        sage: is_sl2z_equivalent(5*z1, 5*z2, 30)
     475        True
    457476    """
    458     # The algorithm we use below was suggested by Samit Dasgupta at a
    459     # UC Santa Cruz seminar in 2011.
    460     if not is_sl2z_equivalent(z1, z2, prec):
    461         # points are not even sl2z-equivalent, so they can't
    462         # be Gamma_0(N) equivalent
     477    C = ComplexField(prec)
     478    eps = RR(2)**(-prec)
     479    w1, g1 = sl2z_rep_in_fundom(z1, eps)  # g1(z1) = w1 = canonical rep
     480    w2, g2 = sl2z_rep_in_fundom(z2, eps)  # g2(z2) = w2 = canonical rep
     481    a1, g1 = canonicalize_sl2z(w1, g1)
     482    a2, g2 = canonicalize_sl2z(w2, g2)
     483    if abs(a1 - a2) >= eps:
     484        # The points are not even sl2z-equivalent, so they can't be
     485        # Gamma_0(N) equivalent
    463486        return False
    464     # Now we know the two points are SL2Z-equivalent, so j(z) takes
    465     # the same value on z1 and z2.  The modular function j_N(z)=j(N*z)
    466     # on Gamma_0(N) take on the same value on both z1 and z2 if and
    467     # only if N*z1 is equivalent to N*z2.
    468487
    469     # It is a classical theorem that j(z) and j_N(z) generate the
    470     # field of modular functions for Gamma_0(N), i.e., (via standard
    471     # identifications) the field of rational functions on the
    472     # algebraic curve X_0(N)/C.  Thus: if two points are equivalent
    473     # modulo Gamma_0(N), then j(z) and j_N(z) must take on the same
    474     # values on those points; conversely, if j(z) and j_N(z) take the
    475     # same values on two points, then every rational function takes
    476     # the same values on the image of the two points in X_0(N), which
    477     # implies that the points define the same element of X_0(N) (since
    478     # a point on a nonsingular algebraic curve is characterized by the
    479     # set of functions that vanish at it).
     488    # Now we may assume that g1(z1) = g2(z2), because of the adjustments
     489    # made above.  We double check.
     490    # This C2 is purely used for the assert below and nothing else, so this
     491    # "magic constant" not so evil.
     492    C2 = ComplexField(2*prec+10) 
     493    assert abs(g1.acton(C2(z1)) - g2.acton(C2(z2))) < eps
     494   
     495    # So now we have z := g1(z1) = g2(z2), both in the standard
     496    # fundamental domain.
     497    #
     498    # The nontrivial elements of PSL2Z with a fixed point z in the
     499    # standard fundamental domain for the upper half plane are
     500    # Stab(z), where
     501    #
     502    #     * z = i, so Stab(z) generated by S (where S has order 2)
     503    #     * z = rho = exp(2*pi*i/3) so Stab(z) generated by S*T
     504    #     * z = -rhobar = exp(pi*i/3) so Stab(z) generated by T*S
     505    #
     506    # The elements in PSL2Z that send z1 to z2 are the elements
     507    # g2^(-1)*A*g1 for A in Stab(z), so we just check if any are in
     508    # Gamma0(N).
     509   
     510    g2i = g2**(-1)
     511    g = g2i*g1
     512    i = C.gen()
     513    pi = C.pi()
     514    if g[1,0]%N == 0:
     515        return True
     516    S, T = g1.parent().gens()
     517    if a1 == i:
     518        return (g2i*S*g1)[1,0]%N == 0
     519    elif a1 == ((2*pi*i)/3).exp():
     520        return (g2i*S*T*g1)[1,0]%N == 0 or (g2i*S*T*S*T*g1)[1,0]%N == 0
     521    elif a1 == ((pi*i)/3).exp():
     522        return (g2i*T*S*g1)[1,0]%N == 0 or (g2i*T*S*T*S*g1)[1,0]%N == 0
     523    return False
    480524
    481     # Thus: we check whether or not N*z1 is equivalent to N*z2.
    482     return is_sl2z_equivalent(N*z1, N*z2, prec)
     525
    483526
    484527
    485528class X0NPoint(object):
     
    541584        EXAMPLES::
    542585
    543586            sage: from sage.schemes.elliptic_curves.chow_heegner import X0NPoint
    544             sage: z = CDF(1,1); x=X0NPoint(z,11,30); y=X0NPoint((-5*z-1)/(11*z+2),11,30)
    545             sage: x==y
     587            sage: g = SL2Z([-5,-1,11,2])
     588            sage: z = CDF(1,1); x = X0NPoint(z, 11, 30); y = X0NPoint(g.acton(z), 11, 30)
     589            sage: x == y
    546590            True
    547             sage: z = CDF(1,1); x=X0NPoint(z,11,30); y=X0NPoint(z/2,11,30)
    548             sage: x==y
     591            sage: z = CDF(1,1); x = X0NPoint(z, 11, 30); y = X0NPoint(z/2, 11, 30)
     592            sage: x == y
    549593            False
    550594
    551595        Points with different levels::
     
    618662            [-0.26016260 + 0.0081300813*I]
    619663            sage: x.atkin_lehner(3).atkin_lehner(3)
    620664            [1.0000000 + 1.0000000*I]
     665            sage: x.atkin_lehner(3).atkin_lehner(3) == x
     666            True
    621667            sage: x.atkin_lehner(2)
    622668            Traceback (most recent call last):
    623669            ...
     
    635681            [0.11049724 + 0.00055248619*I]
    636682            sage: P.atkin_lehner(10).atkin_lehner(10)
    637683            [0.10554020 + 0.000013888696*I]
     684
     685        This is because applying Atkin-Lehner results in substantial precision loss::
     686       
    638687            sage: P.atkin_lehner(10).atkin_lehner(10) == P
    639             True           
     688            False
     689
     690        Using higher precision for the underlying point fixes the problem::
     691
     692            sage: P = X0NPoint(ComplexField(200)(1,1),90,30)
     693            sage: P.atkin_lehner(10).atkin_lehner(10) == P
     694            True       
    640695        """
    641696        if q is None:
    642697            # Main involution
     
    18851940        low precision (silly example)::
    18861941       
    18871942            sage: phi.points_in_h(RDF(1), equiv_prec=2, min_imag=1e-3)
    1888             [[-0.38 + 0.0020*I], [-0.047 + 0.016*I]]
     1943            [[-0.38 + 0.0020*I], [-0.094 + 0.0015*I], [-0.047 + 0.016*I]]
    18891944
    18901945        Consider points equivalent only if equivalent to high
    18911946        precision, which will result in falsely considering points as
    18921947        inequivalent::
    18931948
    18941949            sage: len(phi.points_in_h(RDF(1), equiv_prec=50, min_imag=1e-3))
    1895             10
    1896             sage: len(phi.points_in_h(RDF(1), equiv_prec=40, min_imag=1e-3))
     1950            9
     1951            sage: len(phi.points_in_h(RDF(1), equiv_prec=40, min_imag=1e-3))  # correct, since moddeg=2
    18971952            2       
    18981953        """
    18991954        C = parent(z)
     
    21842239          the ``numerical_approx`` method.
    21852240       
    21862241        EXAMPLES::
     2242
     2243            sage: P = EllipticCurve('37a').chow_heegner_point(EllipticCurve('37b'))
     2244            sage: P.point_exact(deg1=100)
     2245            (6 : 14 : 1)
     2246
     2247        Another example that illustrates precision issues when testing
     2248        equivalence of points in the upper half plane::
    21872249       
    21882250            sage: P = EllipticCurve('99a').chow_heegner_point(EllipticCurve('99b'))
    2189             sage: P.numerical_approx()
    2190             (1.64062499999... : -1.75195312499... : 1.00000000000000)
    21912251            sage: P.point_exact()
     2252            Traceback (most recent call last):
     2253            ...
     2254            RuntimeError: Found too many points (13 > 12) -- try (good) increasing precision of z (now=53) or (bad) *decreasing* equiv_prec (now=17)
     2255
     2256        As suggested, either increasing the precision (slower, but
     2257        safer) or somewhat decreasing the number equiv_prec of bits
     2258        used for checking equivalence (faster, but less safe) both
     2259        work in this case::
     2260
     2261            sage: P.point_exact(equiv_prec=12)
    21922262            (105/64 : -897/512 : 1)
    2193             sage: P.index()
     2263            sage: P.point_exact(prec=60)
     2264            (105/64 : -897/512 : 1)
     2265            sage: P.index(equiv_prec=12)
    21942266            4
    21952267        """
    21962268        P = self.numerical_approx(**kwds)
  • sage/schemes/elliptic_curves/ell_rational_field.py

    diff --git a/sage/schemes/elliptic_curves/ell_rational_field.py b/sage/schemes/elliptic_curves/ell_rational_field.py
    a b  
    49404940
    49414941            sage: P = EllipticCurve('37a').chow_heegner_point(EllipticCurve('37b')); P
    49424942            Chow-Heegner point on 37a1 associated to 37b1
    4943             sage: P.point_exact()
     4943            sage: P.point_exact(deg1=100)
    49444944            (6 : 14 : 1)
    4945             sage: P.numerical_approx()
     4945            sage: P.numerical_approx(deg1=100)
    49464946            (6.00000000000... : 14.0000000000... : 1.00000000000000)       
    49474947        """
    49484948        import chow_heegner