Ticket #9384: trac_9384_update.patch

File trac_9384_update.patch, 4.3 KB (added by ebeyerstedt, 12 years ago)

Replaces previous patch.

  • sage/schemes/elliptic_curves/ell_field.py

    # HG changeset patch
    # User Erin Beyerstedt <ebeyerst@tulane.edu>
    # Date 1278096443 25200
    # Node ID 052083f4e622cfe5647b1f9249a231a0bf92b6df
    # Parent  d48975df81603039025f38745ced56e25bdc3c0c
    [mq]: trac_9384.patch
    
    diff -r d48975df8160 -r 052083f4e622 sage/schemes/elliptic_curves/ell_field.py
    a b  
    547547
    548548        return D
    549549
     550    def descend_to(self, K, f=None):
     551        r"""
     552        Given a subfield `K` and an elliptic curve self defined over a field `L`,
     553        this function determines whether there exists an elliptic curve over `K`
     554        which is isomorphic over `L` to self. If one exists, it finds it.
     555
     556        INPUT:
     557               
     558        - `K` -- a subfield of the base field of self.
     559        - `f` -- an embedding of `K` into the base field of self.
     560
     561        OUTPUT:
     562
     563        Either an elliptic curve defined over `K` which is isomorphic to self
     564        or None if no such curve exists.
     565
     566        .. NOTE::
     567
     568            This only works over number fields and QQ.
     569
     570        EXAMPLES::
     571
     572            sage: E = EllipticCurve([1,2,3,4,5])
     573            sage: E.descend_to(ZZ)
     574            Traceback (most recent call last):
     575            ...
     576            TypeError: Input must be a field.
     577           
     578        ::
     579
     580            sage: F.<b> = QuadraticField(23)
     581            sage: G.<a> = F.extension(x^3+5)
     582            sage: E = EllipticCurve(j=1728*b).change_ring(G)
     583            sage: E.descend_to(F)
     584            Elliptic Curve defined by y^2 = x^3 + (8957952*b-206032896)*x + (-247669456896*b+474699792384) over Number Field in b with defining polynomial x^2 - 23
     585
     586        ::
     587       
     588            sage: L.<a> = NumberField(x^4 - 7)
     589            sage: K.<b> = NumberField(x^2 - 7)
     590            sage: E = EllipticCurve([a^6,0])
     591            sage: E.descend_to(K)
     592            Elliptic Curve defined by y^2 = x^3 + 1296/49*b*x over Number Field in b with defining polynomial x^2 - 7
     593
     594        ::
     595       
     596            sage: K.<a> = QuadraticField(17)
     597            sage: E = EllipticCurve(j = 2*a)
     598            sage: print E.descend_to(QQ)
     599            None
     600        """
     601        if not K.is_field():
     602            raise TypeError, "Input must be a field."
     603        if self.base_field()==K:
     604            return self
     605        j = self.j_invariant()
     606        from sage.rings.all import QQ
     607        if K == QQ:
     608            f = QQ.embeddings(self.base_field())[0]
     609            if j in QQ:
     610                jbase = QQ(j)
     611            else:
     612                return None
     613        elif f == None:
     614            embeddings = K.embeddings(self.base_field())
     615            if len(embeddings) == 0:
     616                raise TypeError, "Input must be a subfield of the base field of the curve."
     617            for g in embeddings:
     618                try:
     619                    jbase = g.preimage(j)
     620                    f = g
     621                    break
     622                except:
     623                    pass
     624            if f == None:
     625                return None
     626        else:
     627            try:
     628                jbase = f.preimage(j)
     629            except:
     630                return None
     631        E = EllipticCurve(j=jbase)
     632        E2 = EllipticCurve(self.base_field(), [f(a) for a in E.a_invariants()])
     633        if jbase==0:
     634            d = self.is_sextic_twist(E2)
     635            if d == 1:
     636                return E
     637            if d == 0:
     638                return None
     639            Etwist = E2.sextic_twist(d)
     640        elif jbase==1728:
     641            d = self.is_quartic_twist(E2)
     642            if d == 1:
     643                return E
     644            if d == 0:
     645                return None
     646            Etwist = E2.quartic_twist(d)
     647        else:
     648            d = self.is_quadratic_twist(E2)
     649            if d == 1:
     650                return E
     651            if d == 0:
     652                return None
     653            Etwist = E2.quadratic_twist(d)
     654        if Etwist.is_isomorphic(self):
     655            try:
     656                Eout = EllipticCurve(K, [f.preimage(a) for a in Etwist.a_invariants()])
     657            except:
     658                return None
     659            else:
     660                return Eout
    550661   
    551662    def isogeny(self, kernel, codomain=None, degree=None, model=None, check=True):
    552663        r"""