Coverage for local/lib/python2.7/site-packages/sage/modular/etaproducts.py: 76%

return 1/ZZ(24)/gcd(cusp.width(), self.level()//cusp.width()) * sum( [ell*self.r(ell)/cusp.width() * (gcd(cusp.width(), self.level()//ell))**2 for ell in self._keys] )

def divisor(self):

r"""

Return the divisor of self, as a formal sum of CuspFamily objects.

EXAMPLES::

sage: e = EtaProduct(12, {1:-336, 2:576, 3:696, 4:-216, 6:-576, 12:-144})

sage: e.divisor() # FormalSum seems to print things in a random order?

-131*(Inf) - 50*(c_{2}) + 11*(0) + 50*(c_{6}) + 169*(c_{4}) - 49*(c_{3})

sage: e = EtaProduct(2^8, {8:1,32:-1})

sage: e.divisor() # random

-(c_{2}) - (Inf) - (c_{8,2}) - (c_{8,3}) - (c_{8,4}) - (c_{4,2}) - (c_{8,1}) - (c_{4,1}) + (c_{32,4}) + (c_{32,3}) + (c_{64,1}) + (0) + (c_{32,2}) + (c_{64,2}) + (c_{128}) + (c_{32,1})

"""

return FormalSum([ (self.order_at_cusp(c), c) for c in AllCusps(self.level())])

def degree(self):

r"""

Return the degree of self as a map

`X_0(N) \to \mathbb{P}^1`, which is equal to the sum of

all the positive coefficients in the divisor of self.

EXAMPLES::

sage: e = EtaProduct(12, {1:-336, 2:576, 3:696, 4:-216, 6:-576, 12:-144})

sage: e.degree()

230

"""

return sum( [self.order_at_cusp(c) for c in AllCusps(self.level()) if self.order_at_cusp(c) > 0])

# def plot(self):

# r""" Returns an error as it's not clear what plotting an eta product means. """

# raise NotImplementedError

def r(self, d):

r"""

Return the exponent `r_d` of `\eta(q^d)` in self.

EXAMPLES::

sage: e = EtaProduct(12, {2:24, 3:-24})

sage: e.r(3)

-24

sage: e.r(4)

"""

return self._rdict.get(d, 0)

# def __call__(self, cusp):

# r""" Calculate the value of self at the given cusp. """

# assert self.level() == cusp.level()

# if self.order_at_cusp(cusp) < 0:

# return Infinity

# if self.order_at_cusp(cusp) > 0:

# return 0

# else:

# s = ZZ(1)

# for ell in divisors(self.level()):

# s *= 1/ZZ(cusp.width())*gcd(cusp.width(), self.level() // ell)**(self.r(ell) / ZZ(2))

# return s

def num_cusps_of_width(N, d):

r"""

Return the number of cusps on `X_0(N)` of width d.

INPUT:

- ``N`` - (integer): the level

- ``d`` - (integer): an integer dividing N, the cusp

width

EXAMPLES::

sage: [num_cusps_of_width(18,d) for d in divisors(18)]

[1, 1, 2, 2, 1, 1]

sage: num_cusps_of_width(4,8)

Traceback (most recent call last):

...

ValueError: N and d must be positive integers with d|N

"""

N = ZZ(N)

d = ZZ(d)

if N <= 0 or d <= 0 or (N % d) != 0:

raise ValueError("N and d must be positive integers with d|N")

return euler_phi(gcd(d, N//d))

def AllCusps(N):

r"""

Return a list of CuspFamily objects corresponding to the cusps of

`X_0(N)`.

INPUT:

- ``N`` - (integer): the level

EXAMPLES::

sage: AllCusps(18)

[(Inf), (c_{2}), (c_{3,1}), (c_{3,2}), (c_{6,1}), (c_{6,2}), (c_{9}), (0)]

sage: AllCusps(0)

Traceback (most recent call last):

...

ValueError: N must be positive

"""

from six.moves import range

N = ZZ(N)

if N <= 0:

raise ValueError("N must be positive")

c = []

for d in divisors(N):

n = num_cusps_of_width(N, d)

if n == 1:

c.append(CuspFamily(N, d))

elif n > 1:

for i in range(n):

c.append(CuspFamily(N, d, label=str(i + 1)))

return c

class CuspFamily(SageObject):

r"""

A family of elliptic curves parametrising a region of

`X_0(N)`.

"""

def __init__(self, N, width, label = None):

r"""

Create the cusp of width d on X_0(N) corresponding to the family

of Tate curves

`(\CC_p/q^d, \langle \zeta q\rangle)`. Here

`\zeta` is a primitive root of unity of order `r`

with `\mathrm{lcm}(r,d) = N`. The cusp doesn't store zeta,

so we store an arbitrary label instead.

EXAMPLES::

sage: CuspFamily(8, 4)

(c_{4})

sage: CuspFamily(16, 4, '1')

(c_{4,1})

"""

N = ZZ(N)

if N <= 0:

raise ValueError("N must be positive")

self._N = N

self._width = width

if (N % width):

raise ValueError("Bad width")

if num_cusps_of_width(N, width) > 1 and label is None:

raise ValueError("There are %s > 1 cusps of width %s on X_0(%s): specify a label" % (num_cusps_of_width(N,width), width, N))

if num_cusps_of_width(N, width) == 1 and label is not None:

raise ValueError("There is only one cusp of width %s on X_0(%s): no need to specify a label" % (width, N))

self.label = label

def width(self):

r"""

The width of this cusp.

EXAMPLES::

sage: e = CuspFamily(10, 1)

sage: e.width()

"""

return self._width

def level(self):

r"""

The level of this cusp.

EXAMPLES::

sage: e = CuspFamily(10, 1)

sage: e.level()

"""

return self._N

def sage_cusp(self):

"""

Return the corresponding element of

`\mathbb{P}^1(\QQ)`.

EXAMPLES::

sage: CuspFamily(10, 1).sage_cusp() # not implemented

Infinity

"""

raise NotImplementedError

def _repr_(self):

r"""

Return a string representation of self.

EXAMPLES::

sage: CuspFamily(16, 4, "1")._repr_()

'(c_{4,1})'

"""

if self.width() == 1:

return "(Inf)"

elif self.width() == self.level():

return "(0)"

else:

return "(c_{%s%s})" % (self.width(), ((self.label and (","+self.label)) or ""))

def qexp_eta(ps_ring, prec):

r"""

Return the q-expansion of `\eta(q) / q^{1/24}`, where

`\eta(q)` is Dedekind's function

.. MATH::

\eta(q) = q^{1/24}\prod_{n=1}^\infty (1-q^n),

as an element of ps_ring, to precision prec.

INPUT:

- ``ps_ring`` - (PowerSeriesRing): a power series ring

- ``prec`` - (integer): the number of terms to compute.

OUTPUT: An element of ps_ring which is the q-expansion of

`\eta(q)/q^{1/24}` truncated to prec terms.

ALGORITHM: We use the Euler identity

.. MATH::

\eta(q) = q^{1/24}( 1 + \sum_{n \ge 1} (-1)^n (q^{n(3n+1)/2} + q^{n(3n-1)/2})

to compute the expansion.

EXAMPLES::

sage: qexp_eta(ZZ[['q']], 100)

1 - q - q^2 + q^5 + q^7 - q^12 - q^15 + q^22 + q^26 - q^35 - q^40 + q^51 + q^57 - q^70 - q^77 + q^92 + O(q^100)

"""

prec = Integer(prec)

assert prec>0, "prec must be a positive integer"

v = [Integer(0)] * prec

pm = Integer(1)

v[0] = pm

try:

n = 1

while True:

pm = -pm

v[n*(3*n-1)//2] = pm

v[n*(3*n+1)//2] = pm

n += 1

except IndexError:

pass

return ps_ring(v, prec=prec)

def eta_poly_relations(eta_elements, degree, labels=['x1','x2'], verbose=False):

r"""

Find polynomial relations between eta products.

INPUT:

- ``eta_elements`` - (list): a list of EtaGroupElement objects.

Not implemented unless this list has precisely two elements. degree

- ``degree`` - (integer): the maximal degree of polynomial to look for.

- ``labels`` - (list of strings): labels to use for the polynomial returned.

- ``verbose``` - (boolean, default False): if True, prints information as

it goes.

OUTPUT: a list of polynomials which is a Groebner basis for the

part of the ideal of relations between eta_elements which is

generated by elements up to the given degree; or None, if no

relations were found.

ALGORITHM: An expression of the form

`\sum_{0 \le i,j \le d} a_{ij} x^i y^j` is zero if and

only if it vanishes at the cusp infinity to degree at least

`v = d(deg(x) + deg(y))`. For all terms up to `q^v`

in the `q`-expansion of this expression to be zero is a

system of `v + k` linear equations in `d^2`

coefficients, where `k` is the number of nonzero negative

coefficients that can appear.

Solving these equations and calculating a basis for the solution

space gives us a set of polynomial relations, but this is generally

far from a minimal generating set for the ideal, so we calculate a

Groebner basis.

As a test, we calculate five extra terms of `q`-expansion

and check that this doesn't change the answer.

EXAMPLES::

sage: t = EtaProduct(26, {2:2,13:2,26:-2,1:-2})

sage: u = EtaProduct(26, {2:4,13:2,26:-4,1:-2})

sage: eta_poly_relations([t, u], 3)

sage: eta_poly_relations([t, u], 4)

[x1^3*x2 - 13*x1^3 - 4*x1^2*x2 - 4*x1*x2 - x2^2 + x2]

Use verbose=True to see the details of the computation::

sage: eta_poly_relations([t, u], 3, verbose=True)

Trying to find a relation of degree 3

Lowest order of a term at infinity = -12

Highest possible degree of a term = 15

Trying all coefficients from q^-12 to q^15 inclusive

No polynomial relation of order 3 valid for 28 terms

Check:

Trying all coefficients from q^-12 to q^20 inclusive

No polynomial relation of order 3 valid for 33 terms

sage: eta_poly_relations([t, u], 4, verbose=True)

Trying to find a relation of degree 4

Lowest order of a term at infinity = -16

Highest possible degree of a term = 20

Trying all coefficients from q^-16 to q^20 inclusive

Check:

Trying all coefficients from q^-16 to q^25 inclusive

[x1^3*x2 - 13*x1^3 - 4*x1^2*x2 - 4*x1*x2 - x2^2 + x2]

"""

if len(eta_elements) > 2:

raise NotImplementedError("Don't know how to find relations between more than two elements")

eta1, eta2 = eta_elements

if verbose:

print("Trying to find a relation of degree %s" % degree)

inf = CuspFamily(eta1.level(), 1)

loterm = -(min([0, eta1.order_at_cusp(inf)]) + min([0,eta2.order_at_cusp(inf)]))*degree

if verbose:

print("Lowest order of a term at infinity = %s" % -loterm)

maxdeg = sum([eta1.degree(), eta2.degree()])*degree

if verbose:

print("Highest possible degree of a term = %s" % maxdeg)

m = loterm + maxdeg + 1

oldgrob = _eta_relations_helper(eta1, eta2, degree, m, labels, verbose)

if verbose:

print("Check:")

newgrob = _eta_relations_helper(eta1, eta2, degree, m+5, labels, verbose)

if oldgrob != newgrob:

if verbose:

raise ArithmeticError("Answers different!")

else:

raise ArithmeticError("Check: answers different!")

return newgrob

def _eta_relations_helper(eta1, eta2, degree, qexp_terms, labels, verbose):

r"""

Helper function used by eta_poly_relations. Finds a basis for the

space of linear relations between the first qexp_terms of the

`q`-expansions of the monomials

`\eta_1^i * \eta_2^j` for `0 \le i,j < degree`,

and calculates a Groebner basis for the ideal generated by these

relations.

Liable to return meaningless results if qexp_terms isn't at least

`1 + d*(m_1,m_2)` where

.. MATH::

m_i = min(0, {\text degree of the pole of $\eta_i$ at $\infty$})

as then 1 will be in the ideal.

EXAMPLES::

sage: from sage.modular.etaproducts import _eta_relations_helper

sage: r,s = EtaGroup(4).basis()

sage: _eta_relations_helper(r,s,4,100,['a','b'],False)

[a*b - a + 16]

sage: _eta_relations_helper(EtaProduct(26, {2:2,13:2,26:-2,1:-2}),EtaProduct(26, {2:4,13:2,26:-4,1:-2}),3,12,['a','b'],False) # not enough terms, will return rubbish

[1]

"""

from six.moves import range

indices = [(i,j) for j in range(degree) for i in range(degree)]

inf = CuspFamily(eta1.level(), 1)

pole_at_infinity = -(min([0, eta1.order_at_cusp(inf)]) + min([0,eta2.order_at_cusp(inf)]))*degree

if verbose:

print("Trying all coefficients from q^%s to q^%s inclusive" % (-pole_at_infinity, -pole_at_infinity + qexp_terms - 1))

rows = []

for j in range(qexp_terms):

rows.append([])

for i in indices:

func = (eta1**i[0]*eta2**i[1]).qexp(qexp_terms)

for j in range(qexp_terms):

rows[j].append(func[j - pole_at_infinity])

M = matrix(rows)

V = M.right_kernel()

if V.dimension() == 0:

if verbose:

print("No polynomial relation of order %s valid for %s terms" % (degree, qexp_terms))

return None

if V.dimension() >= 1:

R = PolynomialRing(QQ, 2, labels)

x,y = R.gens()

relations = []

for c in V.basis():

relations.append(sum( [ c[v] * x**indices[v][0] * y**indices[v][1] for v in range(len(indices))]))

id = R.ideal(relations)

return id.groebner_basis()

Coverage for local/lib/python2.7/site-packages/sage/modular/etaproducts.py : 76%

288 statements 219 run 69 missing 0 excluded