def weil_representation(self) :
        r"""
        OUTPUT:
        
        - A pair of matrices corresponding to T and S.
        """
        disc_bilinear = lambda a, b: (self._dual_basis * vector(QQ, a.lift())) * self._L * (self._dual_basis * vector(QQ, b.lift()))
        disc_quadratic = lambda a: disc_bilinear(a, a) / ZZ(2)
        
        zeta_order = ZZ(lcm([8, 12, prod(self.invariants())] + map(lambda ed: 2 * ed, self.invariants())))
        K = CyclotomicField(zeta_order); zeta = K.gen()

        R = PolynomialRing(K, 'x'); x = R.gen()
#        sqrt2s = (x**2 - 2).factor()
#        if sqrt2s[0][0][0].complex_embedding().real() > 0 :        
#            sqrt2  = sqrt2s[0][0][0]
#        else : 
#            sqrt2  = sqrt2s[0][1]
        Ldet_rts = (x**2 - prod(self.invariants())).factor()
        if Ldet_rts[0][0][0].complex_embedding().real() > 0 :
            Ldet_rt  = Ldet_rts[0][0][0] 
        else :
            Ldet_rt  = Ldet_rts[0][0][0]
                
        Tmat  = diagonal_matrix( K, [zeta**(zeta_order*disc_quadratic(a)) for a in self] )
        Smat = zeta**(zeta_order / 8 * self._L.nrows()) / Ldet_rt  \
               * matrix( K,  [ [ zeta**ZZ(-zeta_order * disc_bilinear(gamma,delta))
                                 for delta in self ]
                               for gamma in self ])
        
        return (Tmat, Smat)

    def _sage_(self):
        """
        EXAMPLES:
            sage: m = lie('[[1,0,3,3],[12,4,-4,7],[-1,9,8,0],[3,-5,-2,9]]') # optional - lie
            sage: m.sage()  # optional - lie
            [ 1  0  3  3]
            [12  4 -4  7]
            [-1  9  8  0]
            [ 3 -5 -2  9]

        """
        t = self.type()
        if t == "grp":
            raise ValueError, "cannot convert Lie groups to native Sage objects"
        elif t == "mat":
            import sage.matrix.constructor

            return sage.matrix.constructor.matrix(eval(str(self).replace("\n", "").strip()))
        elif t == "pol":
            import sage.misc.misc
            from sage.rings.all import PolynomialRing, QQ

            # Figure out the number of variables
            s = str(self)
            open_bracket = s.find("[")
            close_bracket = s.find("]")
            nvars = len(s[open_bracket:close_bracket].split(","))

            # create the polynomial ring
            R = PolynomialRing(QQ, nvars, "x")
            x = R.gens()
            pol = R(0)

            # Split up the polynomials into terms
            terms = []
            for termgrp in s.split(" - "):
                # The first entry in termgrp has
                # a negative coefficient
                termgrp = "-" + termgrp.strip()
                terms += termgrp.split("+")
            # Make sure we don't accidentally add a negative
            # sign to the first monomial
            if s[0] != "-":
                terms[0] = terms[0][1:]

            # go through all the terms in s
            for term in terms:
                xpos = term.find("X")
                coef = eval(term[:xpos].strip())
                exps = eval(term[xpos + 1 :].strip())
                monomial = sage.misc.misc.prod(map(lambda i: x[i] ** exps[i], range(nvars)))
                pol += coef * monomial

            return pol
        elif t == "tex":
            return repr(self)
        elif t == "vid":
            return None
        else:
            return ExpectElement._sage_(self)

def lfsr_connection_polynomial(s):
    """
    INPUT:

    - ``s`` -- a sequence of elements of a finite field of even length

    OUTPUT:

    - ``C(x)`` -- the connection polynomial of the minimal LFSR.

    This implements the algorithm in section 3 of J. L. Massey's article
    [Mas1969]_.

    EXAMPLES::

        sage: F = GF(2)
        sage: F
        Finite Field of size 2
        sage: o = F(0); l = F(1)
        sage: key = [l,o,o,l]; fill = [l,l,o,l]; n = 20
        sage: s = lfsr_sequence(key,fill,n); s
        [1, 1, 0, 1, 0, 1, 1, 0, 0, 1, 0, 0, 0, 1, 1, 1, 1, 0, 1, 0]
        sage: lfsr_connection_polynomial(s)
        x^4 + x + 1
        sage: from sage.matrix.berlekamp_massey import berlekamp_massey
        sage: berlekamp_massey(s)
        x^4 + x^3 + 1

    Notice that ``berlekamp_massey`` returns the reverse of the connection
    polynomial (and is potentially must faster than this implementation).
    """
    # Initialization:
    FF = s[0].base_ring()
    R = PolynomialRing(FF, "x")
    x = R.gen()
    C = R(1); B = R(1); m = 1; b = FF(1); L = 0; N = 0

    while N < len(s):
        if L > 0:
            r = min(L+1,C.degree()+1)
            d = s[N] + sum([(C.list())[i]*s[N-i] for i in range(1,r)])
        if L == 0:
            d = s[N]
        if d == 0:
            m += 1
            N += 1
        if d > 0:
            if 2*L > N:
                C = C - d*b**(-1)*x**m*B
                m += 1
                N += 1
            else:
                T = C
                C = C - d*b**(-1)*x**m*B
                L = N + 1 - L
                m = 1
                b = d
                B = T
                N += 1
    return C

    def __init__(self, R, n, q=None):
        """
        TESTS::
        
            sage: HeckeAlgebraSymmetricGroupT(QQ, 3)
            Hecke algebra of the symmetric group of order 3 on the T basis over Univariate Polynomial Ring in q over Rational Field
        
        ::
        
            sage: HeckeAlgebraSymmetricGroupT(QQ, 3, q=1)
            Hecke algebra of the symmetric group of order 3 with q=1 on the T basis over Rational Field
        """
        self.n = n
        self._basis_keys = permutation.Permutations(n)
        self._name = "Hecke algebra of the symmetric group of order %s"%self.n
        self._one = permutation.Permutation(range(1,n+1))

        if q is None:
            q = PolynomialRing(R, 'q').gen()
            R = q.parent()
        else:
            if q not in R:
                raise ValueError, "q must be in R (= %s)"%R
            self._name += " with q=%s"%q

        self._q = q
        
        CombinatorialAlgebra.__init__(self, R)
        # _repr_ customization: output the basis element indexed by [1,2,3] as [1,2,3]
        self.print_options(prefix="")

def compute_tau0(v0,gamma,wD,return_exact = False):
    r'''
    INPUT:

     - v0: F -> its localization at p
     - gamma: the image of wD (the generator for an order of F) under an optimal embedding

    OUTPUT:

     The element tau_0 such that gamma * [tau_0,1] = wD * [tau_0,1]

    '''
    R = PolynomialRing(QQ,names = 'X')
    X = R.gen()
    F = v0.domain()
    Cp = v0.codomain()
    assert wD.minpoly() == gamma.minpoly()
    a,b,c,d = gamma.list()
    tau0_vec = (c*X**2+(d-a)*X-b).roots(F)
    tau0 = v0(tau0_vec[0][0])
    idx = 0
    if c * tau0 + d != v0(wD):
        tau0 = v0(tau0_vec[1][0])
        idx = 1
    return tau0_vec[idx][0] if return_exact == True else tau0

def is_cm_j_invariant(j):
    """
    Return whether or not this is a CM `j`-invariant.

    INPUT:

    - ``j`` -- an element of a number field `K`

    OUTPUT:

    A pair (bool, (d,f)) which is either (False, None) if `j` is not a
    CM j-invariant or (True, (d,f)) if `j` is the `j`-invariant of the
    imaginary quadratic order of discriminant `D=df^2` where `d` is
    the associated fundamental discriminant and `f` the index.

    .. note::

       The current implementation makes use of the classification of
       all orders of class number up to 100, and hence will raise an
       error if `j` is an algebraic integer of degree greater than
       this.  It would be possible to implement a more general
       version, using the fact that `d` must be supported on the
       primes dividing the discriminant of the minimal polynomial of
       `j`.

    EXAMPLES::

        sage: from sage.schemes.elliptic_curves.cm import is_cm_j_invariant
        sage: is_cm_j_invariant(0)
        (True, (-3, 1))
        sage: is_cm_j_invariant(8000)
        (True, (-8, 1))

        sage: K.<a> = QuadraticField(5)
        sage: is_cm_j_invariant(282880*a + 632000)
        (True, (-20, 1))
        sage: K.<a> = NumberField(x^3 - 2)
        sage: is_cm_j_invariant(31710790944000*a^2 + 39953093016000*a + 50337742902000)
        (True, (-3, 6))

    TESTS::

        sage: from sage.schemes.elliptic_curves.cm import is_cm_j_invariant
        sage: all([is_cm_j_invariant(j) == (True, (d,f)) for d,f,j in cm_j_invariants_and_orders(QQ)])
        True

    """
    from sage.rings.all import NumberFieldElement
    if not isinstance(j, NumberFieldElement) and not j in QQ:
        raise NotImplementedError("is_cm_j_invariant() is only implemented for number field elements")
    if not j.is_integral():
        return False, None
    jpol = PolynomialRing(QQ,'x')([-j,1]) if j in QQ else j.absolute_minpoly()
    h = jpol.degree()
    if h>100:
        raise NotImplementedError("CM data only available for class numbers up to 100")
    for d,f in cm_orders(h):
        if jpol == hilbert_class_polynomial(d*f**2):
            return True, (d,f)
    return False, None

    def _sage_(self):
        """
        EXAMPLES::

            sage: m = lie('[[1,0,3,3],[12,4,-4,7],[-1,9,8,0],[3,-5,-2,9]]') # optional - lie
            sage: m.sage()  # optional - lie
            [ 1  0  3  3]
            [12  4 -4  7]
            [-1  9  8  0]
            [ 3 -5 -2  9]

        """
        t = self.type()
        if t == 'grp':
            raise ValueError("cannot convert Lie groups to native Sage objects")
        elif t == 'mat':
            import sage.matrix.constructor
            return  sage.matrix.constructor.matrix( eval( str(self).replace('\n','').strip())  )
        elif t == 'pol':
            from sage.rings.all import PolynomialRing, QQ

            #Figure out the number of variables
            s = str(self)
            open_bracket = s.find('[')
            close_bracket = s.find(']')
            nvars = len(s[open_bracket:close_bracket].split(','))

            #create the polynomial ring
            R = PolynomialRing(QQ, nvars, 'x')
            x = R.gens()
            pol = R(0)

            #Split up the polynomials into terms
            terms = []
            for termgrp in s.split(' - '):
                #The first entry in termgrp has
                #a negative coefficient
                termgrp = "-"+termgrp.strip()
                terms += termgrp.split('+')
            #Make sure we don't accidentally add a negative
            #sign to the first monomial
            if s[0] != "-":
                terms[0] = terms[0][1:]

            #go through all the terms in s
            for term in terms:
                xpos = term.find('X')
                coef = eval(term[:xpos].strip())
                exps = eval(term[xpos+1:].strip())
                monomial = prod([x[i]**exps[i] for i in range(nvars)])
                pol += coef * monomial

            return pol
        elif t == 'tex':
            return repr(self)
        elif t == 'vid':
            return None
        else:
            return ExpectElement._sage_(self)

def eisenstein_basis(N, k, verbose=False):
    r"""
    Find spanning list of 'easy' generators for the subspace of
    `M_k(\Gamma_0(N))` generated by level 1 Eisenstein series and
    their images of even integer weights up to `k`.

    INPUT:
        - N -- positive integer
        - k -- positive integer
        - ``verbose`` -- bool (default: False)

    OUTPUT:
        - list of monomials in images of level 1 Eisenstein series
        - prec of q-expansions needed to determine element of
          `M_k(\Gamma_0(N))`.

    EXAMPLES::

        sage: from psage.modform.rational.special import eisenstein_basis
        sage: eisenstein_basis(5,4)
        ([E4(q^5)^1, E4(q^1)^1, E2^*(q^5)^2], 3)
        sage: eisenstein_basis(11,2,verbose=True)  # warning below because of verbose
        Warning -- not enough series.
        ([E2^*(q^11)^1], 2)
        sage: eisenstein_basis(11,2,verbose=False)
        ([E2^*(q^11)^1], 2)
    """
    assert N > 1
    if k % 2 != 0:
        return []
    # Make list E of Eisenstein series, to enough precision to
    # determine them, until we span space.
    M = ModularForms(N, k)
    prec = M.echelon_basis()[-1].valuation() + 1
    
    gens = eisenstein_gens(N, k, prec)
    R = PolynomialRing(ZZ, len(gens), ['E%sq%s'%(g[1],g[0]) for g in gens])
    z = [(R.gen(i), g[1]) for i, g in enumerate(gens)]
    m = monomials(z, k)
    
    A = QQ**prec
    V = A.zero_subspace()
    E = []
    for i, z in enumerate(m):
        d = z.degrees()
        f = prod(g[2]**d[i] for i, g in enumerate(gens) if d[i])
        v = A(f.padded_list(prec))
        if v not in V:
            V = V + A.span([v])
            w = [(gens[i][0],gens[i][1],d[i]) for i in range(len(d)) if d[i]]
            E.append(EisensteinMonomial(w))
            if V.dimension() == M.dimension():
                 return E, prec

    if verbose: print "Warning -- not enough series."
    return E, prec

    def local_coordinates_at_infinity(self, prec = 20, name = 't'):
        """
        For the genus `g` hyperelliptic curve `y^2 = f(x)`, return
        `(x(t), y(t))` such that `(y(t))^2 = f(x(t))`, where `t = x^g/y` is
        the local parameter at infinity

        INPUT:

        - ``prec`` -- desired precision of the local coordinates
        - ``name`` -- generator of the power series ring (default: ``t``)

        OUTPUT:

        `(x(t),y(t))` such that `y(t)^2 = f(x(t))` and `t = x^g/y`
        is the local parameter at infinity

        EXAMPLES::

            sage: R.<x> = QQ['x']
            sage: H = HyperellipticCurve(x^5-5*x^2+1)
            sage: x,y = H.local_coordinates_at_infinity(10)
            sage: x
            t^-2 + 5*t^4 - t^8 - 50*t^10 + O(t^12)
            sage: y
            t^-5 + 10*t - 2*t^5 - 75*t^7 + 50*t^11 + O(t^12)

        ::

            sage: R.<x> = QQ['x']
            sage: H = HyperellipticCurve(x^3-x+1)
            sage: x,y = H.local_coordinates_at_infinity(10)
            sage: x
            t^-2 + t^2 - t^4 - t^6 + 3*t^8 + O(t^12)
            sage: y
            t^-3 + t - t^3 - t^5 + 3*t^7 - 10*t^11 + O(t^12)

        AUTHOR:
            - Jennifer Balakrishnan (2007-12)
        """
        g = self.genus()
        pol = self.hyperelliptic_polynomials()[0]
        K = LaurentSeriesRing(self.base_ring(), name, default_prec=prec+2)
        t = K.gen()
        L = PolynomialRing(self.base_ring(),'x')
        x = L.gen()
        i = 0
        w = (x**g/t)**2-pol
        wprime = w.derivative(x)
        x = t**-2
        for i in range((RR(log(prec+2)/log(2))).ceil()):
            x = x-w(x)/wprime(x)
        y = x**g/t
        return x+O(t**(prec+2)) , y+O(t**(prec+2))

    def tiny_integrals_on_basis(self, P, Q):
        r"""
        Evaluate the integrals `\{\int_P^Q x^i dx/2y \}_{i=0}^{2g-1}`
        by formally integrating a power series in a local parameter `t`.
        `P` and `Q` MUST be in the same residue disc for this result to make sense.

        INPUT:

        - P a point on self
        - Q a point on self (in the same residue disc as P)

        OUTPUT:

        The integrals `\{\int_P^Q x^i dx/2y \}_{i=0}^{2g-1}`

        EXAMPLES::

            sage: K = pAdicField(17, 5)
            sage: E = EllipticCurve(K, [-31/3, -2501/108]) # 11a
            sage: P = E(K(14/3), K(11/2))
            sage: TP = E.teichmuller(P);
            sage: E.tiny_integrals_on_basis(P, TP)
            (17 + 14*17^2 + 17^3 + 8*17^4 + O(17^5), 16*17 + 5*17^2 + 8*17^3 + 14*17^4 + O(17^5))

        ::

            sage: K = pAdicField(11, 5)
            sage: x = polygen(K)
            sage: C = HyperellipticCurve(x^5 + 33/16*x^4 + 3/4*x^3 + 3/8*x^2 - 1/4*x + 1/16)
            sage: P = C.lift_x(11^(-2))
            sage: Q = C.lift_x(3*11^(-2))
            sage: C.tiny_integrals_on_basis(P,Q)
            (3*11^3 + 7*11^4 + 4*11^5 + 7*11^6 + 5*11^7 + O(11^8), 3*11 + 10*11^2 + 8*11^3 + 9*11^4 + 7*11^5 + O(11^6), 4*11^-1 + 2 + 6*11 + 6*11^2 + 7*11^3 + O(11^4), 11^-3 + 6*11^-2 + 2*11^-1 + 2 + O(11^2))


        Note that this fails if the points are not in the same residue disc::

            sage: S = C(0,1/4)
            sage: C.tiny_integrals_on_basis(P,S)
            Traceback (most recent call last):
            ...
            ValueError: (11^-2 + O(11^3) : 11^-5 + 8*11^-2 + O(11^0) : 1 + O(11^5)) and (0 : 3 + 8*11 + 2*11^2 + 8*11^3 + 2*11^4 + O(11^5) : 1 + O(11^5)) are not in the same residue disc

        """
        if P == Q:
            V = VectorSpace(self.base_ring(), 2*self.genus())
            return V(0)
        R = PolynomialRing(self.base_ring(), ['x', 'y'])
        x, y = R.gens()
        return self.tiny_integrals([x**i for i in range(2*self.genus())], P, Q)

    def local_coordinates_at_infinity(self, prec=20, name="t"):
        """
        For the genus g hyperelliptic curve y^2 = f(x), returns (x(t), y(t)) such that
        (y(t))^2 = f(x(t)), where t = x^g/y is the local parameter at infinity

        INPUT:
            - prec: desired precision of the local coordinates
            - name: gen of the power series ring (default: 't')

        OUTPUT:
        (x(t),y(t)) such that y(t)^2 = f(x(t)) and t = x^g/y
        is the local parameter at infinity


        EXAMPLES:
            sage: R.<x> = QQ['x']
            sage: H = HyperellipticCurve(x^5-5*x^2+1)
            sage: x,y = H.local_coordinates_at_infinity(10)
            sage: x
            t^-2 + 5*t^4 - t^8 - 50*t^10 + O(t^12)
            sage: y 
            t^-5 + 10*t - 2*t^5 - 75*t^7 + 50*t^11 + O(t^12)

            sage: R.<x> = QQ['x']
            sage: H = HyperellipticCurve(x^3-x+1)
            sage: x,y = H.local_coordinates_at_infinity(10)
            sage: x
            t^-2 + t^2 - t^4 - t^6 + 3*t^8 + O(t^12)
            sage: y
            t^-3 + t - t^3 - t^5 + 3*t^7 - 10*t^11 + O(t^12)


        AUTHOR:
            - Jennifer Balakrishnan (2007-12)
        """
        g = self.genus()
        pol = self.hyperelliptic_polynomials()[0]
        K = LaurentSeriesRing(self.base_ring(), name)
        t = K.gen()
        K.set_default_prec(prec + 2)
        L = PolynomialRing(self.base_ring(), "x")
        x = L.gen()
        i = 0
        w = (x ** g / t) ** 2 - pol
        wprime = w.derivative(x)
        x = t ** -2
        for i in range((RR(log(prec + 2) / log(2))).ceil()):
            x = x - w(x) / wprime(x)
        y = x ** g / t
        return x + O(t ** (prec + 2)), y + O(t ** (prec + 2))

    def __init__(self, params, asym=False):

        (self.n, self.q, sigma, self.sigma_prime, self.k) = params

        S, x = PolynomialRing(ZZ, 'x').objgen()
        self.R = S.quotient_ring(S.ideal(x**self.n + 1))

        Sq = PolynomialRing(Zmod(self.q), 'x')
        self.Rq = Sq.quotient_ring(Sq.ideal(x**self.n + 1))

        # draw z_is uniformly from Rq and compute its inverse in Rq
        if asym:
            z = [self.Rq.random_element() for i in range(self.k)]
            self.zinv = [z_i**(-1) for z_i in z]
        else: # or do symmetric version
            z = self.Rq.random_element()
            zinv = z**(-1)
            z, self.zinv = zip(*[(z,zinv) for i in range(self.k)])

        # set up some discrete Gaussians
        DGSL_sigma = DGSL(ZZ**self.n, sigma)
        self.D_sigma = lambda: self.Rq(list(DGSL_sigma()))

        # discrete Gaussian in ZZ^n with stddev sigma_prime, yields random level-0 encodings
        DGSL_sigmap_ZZ = DGSL(ZZ**self.n, self.sigma_prime)
        self.D_sigmap_ZZ = lambda: self.Rq(list(DGSL_sigmap_ZZ()))

        # draw g repeatedly from a Gaussian distribution of Z^n (with param sigma)
        # until g^(-1) in QQ[x]/<x^n + 1> is small (< n^2)
        Sk = PolynomialRing(QQ, 'x')
        K = Sk.quotient_ring(Sk.ideal(x**self.n + 1)) 
        while True:
            l = self.D_sigma()
            ginv_K = K(mod_near_poly(l, self.q))**(-1)
            ginv_size = vector(ginv_K).norm()

            if ginv_size < self.n**2:
                g = self.Rq(l)
                self.ginv = g**(-1)
                break

        # discrete Gaussian in I = <g>, yields random encodings of 0
        short_g = vector(ZZ, mod_near_poly(g,self.q))
        DGSL_sigmap_I = DGSL(short_g, self.sigma_prime)
        self.D_sigmap_I = lambda: self.Rq(list(DGSL_sigmap_I()))

        # compute zero-testing parameter p_zt
        # randomly draw h (in Rq) from a discrete Gaussian with param q^(1/2)
        self.h = self.Rq(list(DGSL(ZZ**self.n, round(sqrt(self.q)))()))

        # create p_zt
        self.p_zt = self.ginv * self.h * prod(z)

def _check_muqt(mu, q, t, pi=None):
    """
    EXAMPLES::

        sage: from sage.combinat.sf.ns_macdonald import _check_muqt
        sage: P, q, t, n, R, x = _check_muqt([0,0,1],None,None)
        sage: P
        Fraction Field of Multivariate Polynomial Ring in q, t over Rational Field
        sage: q
        q
        sage: t
        t
        sage: n
        Nonattacking fillings of [0, 0, 1]
        sage: R
        Multivariate Polynomial Ring in x0, x1, x2 over Fraction Field of Multivariate Polynomial Ring in q, t over Rational Field
        sage: x
        (x0, x1, x2)

    ::

        sage: q,t = var('q,t')
        sage: P, q, t, n, R, x = _check_muqt([0,0,1],q,None)
        Traceback (most recent call last):
        ...
        ValueError: you must specify either both q and t or neither of them

    ::

        sage: P, q, t, n, R, x = _check_muqt([0,0,1],q,2)
        Traceback (most recent call last):
        ...
        ValueError: the parents of q and t must be the same
    """
    if q is None and t is None:
        P = PolynomialRing(QQ,'q,t').fraction_field()
        q,t = P.gens()
    elif q is not None and t is not None:
        if q.parent() != t.parent():
            raise ValueError("the parents of q and t must be the same")
        P = q.parent()
    else:
        raise ValueError("you must specify either both q and t or neither of them")
    n = NonattackingFillings(mu, pi)
    R = PolynomialRing(P, len(n._shape), 'x')
    x = R.gens()
    return P, q, t, n, R, x

    def __init__(self, group, base_ring, red_hom):
        r"""
        Abstract (Hecke) forms ring.

        INPUT:

        - ``group``       - The Hecke triangle group (default: ``HeckeTriangleGroup(3)``)
        - ``base_ring``   - The base_ring (default: ``ZZ``).
        - ``red_hom``     - If True then results of binary operations are considered
                            homogeneous whenever it makes sense (default: False).
                            This is mainly used by the (Hecke) forms.

        OUTPUT:

        The corresponding abstract (Hecke) forms ring.

        EXAMPLES::

            sage: from graded_ring import ModularFormsRing
            sage: MR = ModularFormsRing(group=5, base_ring=ZZ, red_hom=True)
            sage: MR
            ModularFormsRing(n=5) over Integer Ring
            sage: MR.group()
            Hecke triangle group for n = 5
            sage: MR.base_ring()
            Integer Ring
            sage: MR.has_reduce_hom()
            True
            sage: MR.is_homogeneous()
            False
        """

        from graded_ring import canonical_parameters
        (group, base_ring, red_hom) = canonical_parameters(group, base_ring, red_hom)

        if (group == infinity):
            raise NotImplementedError

        #if (not group.is_arithmetic() and base_ring.characteristic()>0):
        #    raise NotImplementedError
        #if (base_ring.characteristic().divides(2*group.n()*(group.n()-2))):
        #    raise NotImplementedError
        if (base_ring.characteristic() > 0):
            raise NotImplementedError
        self._group               = group
        self._red_hom             = red_hom
        self._base_ring           = base_ring
        self._coeff_ring          = FractionField(PolynomialRing(base_ring,'d'))
        self._pol_ring            = PolynomialRing(base_ring,'x,y,z,d')
        self._rat_field           = FractionField(self._pol_ring)

        # default values
        self._weight              = None
        self._ep                  = None
        self._analytic_type       = self.AT(["quasi", "mero"])

        self.default_prec(10)
        self.disp_prec(5)
        self.default_num_prec(53)

    def curve_over_ram_extn(self, deg):
        r"""
        Returns self over $\Q_p(p^(1/deg))$

        INPUT:

        - deg: the degree of the ramified extension

        OUTPUT:

        self over $\Q_p(p^(1/deg))$

        EXAMPLES::

            sage: R.<x> = QQ['x']
            sage: H = HyperellipticCurve(x^5-23*x^3+18*x^2+40*x)
            sage: K = Qp(11,5)
            sage: HK = H.change_ring(K)
            sage: HL = HK.curve_over_ram_extn(2)
            sage: HL
            Hyperelliptic Curve over Eisenstein Extension of 11-adic Field with capped relative precision 5 in a defined by (1 + O(11^5))*x^2 + (O(11^6))*x + (10*11 + 10*11^2 + 10*11^3 + 10*11^4 + 10*11^5 + O(11^6)) defined by (1 + O(a^10))*y^2 = (1 + O(a^10))*x^5 + (10 + 8*a^2 + 10*a^4 + 10*a^6 + 10*a^8 + O(a^10))*x^3 + (7 + a^2 + O(a^10))*x^2 + (7 + 3*a^2 + O(a^10))*x

        AUTHOR:

        - Jennifer Balakrishnan

        """
        from sage.schemes.hyperelliptic_curves.constructor import HyperellipticCurve

        K = self.base_ring()
        p = K.prime()
        A = PolynomialRing(QQ, "x")
        x = A.gen()
        J = K.extension(x ** deg - p, names="a")
        pol = self.hyperelliptic_polynomials()[0]
        H = HyperellipticCurve(A(pol))
        HJ = H.change_ring(J)
        self._curve_over_ram_extn = HJ
        self._curve_over_ram_extn._curve_over_Qp = self
        return HJ

def ubs(f):
    r"""
    Given a sextic form `f`, return a dictionary of the invariants of Mestre, p 317 [M]_.

    `f` may be homogeneous in two variables or inhomogeneous in one.

    EXAMPLES::

        sage: from sage.schemes.hyperelliptic_curves.invariants import ubs
        sage: x = QQ['x'].0
        sage: ubs(x^6 + 1)
        {'A': 2, 'C': -2/9, 'B': 2/3, 'D': 0, 'f': x^6 + h^6, 'i': 2*x^2*h^2, 'Delta': -2/3*x^2*h^2, 'y1': 0, 'y3': 0, 'y2': 0}

        sage: R.<u, v> = QQ[]
        sage: ubs(u^6 + v^6)
        {'A': 2, 'C': -2/9, 'B': 2/3, 'D': 0, 'f': u^6 + v^6, 'i': 2*u^2*v^2, 'Delta': -2/3*u^2*v^2, 'y1': 0, 'y3': 0, 'y2': 0}

        sage: R.<t> = GF(31)[]
        sage: ubs(t^6 + 2*t^5 + t^2 + 3*t + 1)
        {'A': 0, 'C': -15, 'B': -12, 'D': -15, 'f': t^6 + 2*t^5*h + t^2*h^4 + 3*t*h^5 + h^6, 'i': -4*t^4 + 10*t^3*h + 2*t^2*h^2 - 9*t*h^3 - 7*h^4, 'Delta': -10*t^4 + 12*t^3*h + 7*t^2*h^2 - 5*t*h^3 + 2*h^4, 'y1': 4*t^2 - 10*t*h - 13*h^2, 'y3': 4*t^2 - 4*t*h - 9*h^2, 'y2': 6*t^2 - 4*t*h + 2*h^2}
    """
    ub = Ueberschiebung
    if f.parent().ngens() == 1:
        f = PolynomialRing(f.parent().base_ring(), 1, f.parent().variable_name())(f)
        x1, x2 = f.homogenize().parent().gens()
        f = sum([ f[i]*x1**i*x2**(6-i) for i in range(7) ])
    U = {}
    U['f'] = f
    U['i'] = ub(f, f, 4)
    U['Delta'] = ub(U['i'], U['i'], 2)
    U['y1'] = ub(f, U['i'], 4)
    U['y2'] = ub(U['i'], U['y1'], 2)
    U['y3'] = ub(U['i'], U['y2'], 2)
    U['A'] = ub(f, f, 6)
    U['B'] = ub(U['i'], U['i'], 4)
    U['C'] = ub(U['i'], U['Delta'], 4)
    U['D'] = ub(U['y3'], U['y1'], 2)
    return U

 def generator_relations(self, K) :
     """
     An ideal `I` in a polynomial ring `R` over `K`, such that the
     associated ring is `R / I` surjects onto the ring of modular forms
     with coefficients in `K`.
     
     INPUT:
         - `K` -- A ring.
         
     OUTPUT:
         An ideal in a polynomial ring.
     
     TESTS::
         sage: from psage.modform.fourier_expansion_framework.modularforms.modularform_testtype import *
         sage: t = ModularFormTestType_vectorvalued()
         sage: t.generator_relations(QQ)
         Ideal (g1^2 - g2, g1^3 - g3, g1^4 - g4, g1^5 - g5) of Multivariate Polynomial Ring in g1, g2, g3, g4, g5, v1, v2, v3 over Rational Field
     """
     if K.has_coerce_map_from(ZZ) :
         R = PolynomialRing(K, self.non_vector_valued()._generator_names(K) + self._generator_names(K))
         return R.ideal().parent()(self.non_vector_valued().generator_relations(K))
         
     raise NotImplementedError

def _load_euler_table(self, n, force=False, verbose=False):
    r"""
    Load the euler table for self over the degree n extension of
    $\mathbb{F}_q$ to disk. If self is the versal j-curve, the table
    is pulled from
    
    SAGE_ROOT/data/jcurve_euler_tables .
    
    Otherwise, the table is pulled from the `user` table

    SAGE_ROOT/data/local_euler_tables .

    This should eventually be implemented using MongoDB.

    It currently doesn't check if the key exist. If the key doesn't
    exist, a RuntimeError is raised by sage.database.db. This
    RuntimeError should be sufficient, so key checking may not be
    necessary.

    INPUT:

        - n -- the degree of the extension of F_q
        - force -- boolean that overwrites self's euler table with
                   one from database

    EXAMPLES::

        sage: import psage
        sage: K.<t> = psage.FunctionField(GF(11))
        sage: E = psage.ellff_EllipticCurve(K,[0,0,0,-27*t/(t-1728),54*t/(t-1728)])
        sage: E._euler_table(1)
        Traceback (most recent call last):
            ...
        RuntimeError: table is empty
        sage: E._load_euler_table(1)
        sage: E._euler_table(1)
        [0, 0, 4, -6, 3, 5, 1, -2, 4, -2, 3, 1]
            
    """
        
    import os
    SAGE_ROOT = os.environ['SAGE_ROOT']
        
    K = self.K
    R = self.R
    t = K.gens()[0]
    p = self.p
    d = self.d
    q = self.q
    R2 = PolynomialRing(GF(q), 's')
    s = R2.gens()[0]
    a1n = R2(0)
    a1d = R2(1)
    a2n = R2(0)
    a2d = R2(1)
    a3n = R2(0)
    a3d = R2(1)
    a4n = R2(self.a4.numerator().coeffs())
    a4d = R2(self.a4.denominator().coeffs())
    a6n = R2(self.a6.numerator().coeffs())
    a6d = R2(self.a6.denominator().coeffs())

    ainvs = [0, 0, 0, self.a4, self.a6]
    ainvs_pairs = ((a1n, a1d), (a2n, a2d), (a3n, a3d), (a4n, a4d), (a6n, a6d))

    # recognize if self is j-curve and use special repository
    if ainvs == [0,0,0,-27*t*(t-1728)**3,54*t*(t-1728)**5]:
        if verbose:
            print 'j-curve recognized; saving euler table to database'
        if not os.path.exists(SAGE_ROOT + '/data/jcurve_euler_tables/jcurve_euler_tables'):
            print 'Database does not exist; cannot load from it'
        else:
            euler_db = jCurveEulerTables()
            # check that keys exist?
            self._set_euler_table(n, euler_db[q][n], force)
            
    # work with user's repository of euler tables
    else:
        if not os.path.exists(SAGE_ROOT + '/data/local_euler_tables/local_euler_tables'):
            print 'Database does not exist; cannot load from it'
        else:                
            local_euler_db = LocalEulerTables()
            # check that keys exist?
            self._set_euler_table(n, local_euler_db[ainvs_pairs][q][n], force)

    def parametrization(self, point=None, morphism=True):
        r"""
        Return a parametrization `f` of ``self`` together with the
        inverse of `f`.

        If ``point`` is specified, then that point is used
        for the parametrization. Otherwise, use ``self.rational_point()``
        to find a point.

        If ``morphism`` is True, then `f` is returned in the form
        of a Scheme morphism. Otherwise, it is a tuple of polynomials
        that gives the parametrization.

        ALGORITHM:

        Uses the PARI/GP function ``qfparam``.

        EXAMPLES ::

            sage: c = Conic([1,1,-1])
            sage: c.parametrization()
            (Scheme morphism:
              From: Projective Space of dimension 1 over Rational Field
              To:   Projective Conic Curve over Rational Field defined by x^2 + y^2 - z^2
              Defn: Defined on coordinates by sending (x : y) to
                    (2*x*y : x^2 - y^2 : x^2 + y^2),
             Scheme morphism:
              From: Projective Conic Curve over Rational Field defined by x^2 + y^2 - z^2
              To:   Projective Space of dimension 1 over Rational Field
              Defn: Defined on coordinates by sending (x : y : z) to
                    (1/2*x : -1/2*y + 1/2*z))

        An example with ``morphism = False`` ::

            sage: R.<x,y,z> = QQ[]
            sage: C = Curve(7*x^2 + 2*y*z + z^2)
            sage: (p, i) = C.parametrization(morphism = False); (p, i)
            ([-2*x*y, x^2 + 7*y^2, -2*x^2], [-1/2*x, 1/7*y + 1/14*z])
            sage: C.defining_polynomial()(p)
            0
            sage: i[0](p) / i[1](p)
            x/y

        A ``ValueError`` is raised if ``self`` has no rational point ::

            sage: C = Conic(x^2 + 2*y^2 + z^2)
            sage: C.parametrization()
            Traceback (most recent call last):
            ...
            ValueError: Conic Projective Conic Curve over Rational Field defined by x^2 + 2*y^2 + z^2 has no rational points over Rational Field!

        A ``ValueError`` is raised if ``self`` is not smooth ::

            sage: C = Conic(x^2 + y^2)
            sage: C.parametrization()
            Traceback (most recent call last):
            ...
            ValueError: The conic self (=Projective Conic Curve over Rational Field defined by x^2 + y^2) is not smooth, hence does not have a parametrization.
        """
        if (not self._parametrization is None) and not point:
            par = self._parametrization
        else:
            if not self.is_smooth():
                raise ValueError("The conic self (=%s) is not smooth, hence does not have a parametrization." % self)
            if point is None:
                point = self.rational_point()
            point = Sequence(point)
            Q = PolynomialRing(QQ, 'x,y')
            [x, y] = Q.gens()
            gens = self.ambient_space().gens()
            M = self.symmetric_matrix()
            M *= lcm([ t.denominator() for t in M.list() ])
            par1 = qfparam(M, point)
            B = Matrix([[par1[i][j] for j in range(3)] for i in range(3)])
            # self is in the image of B and does not lie on a line,
            # hence B is invertible
            A = B.inverse()
            par2 = [sum([A[i,j]*gens[j] for j in range(3)]) for i in [1,0]]
            par = ([Q(pol(x/y)*y**2) for pol in par1], par2)
            if self._parametrization is None:
                self._parametrization = par
        if not morphism:
            return par
        P1 = ProjectiveSpace(self.base_ring(), 1, 'x,y')
        return P1.hom(par[0],self), self.Hom(P1)(par[1], check = False)

class FormsRing_abstract(Parent):
    r"""
    Abstract (Hecke) forms ring.

    This should never be called directly. Instead one should
    instantiate one of the derived classes of this class.
    """

    from graded_ring_element import FormsRingElement
    Element = FormsRingElement

    from analytic_type import AnalyticType
    AT = AnalyticType()

    def __init__(self, group, base_ring, red_hom):
        r"""
        Abstract (Hecke) forms ring.

        INPUT:

        - ``group``       - The Hecke triangle group (default: ``HeckeTriangleGroup(3)``)
        - ``base_ring``   - The base_ring (default: ``ZZ``).
        - ``red_hom``     - If True then results of binary operations are considered
                            homogeneous whenever it makes sense (default: False).
                            This is mainly used by the (Hecke) forms.

        OUT