Coverage for local/lib/python2.7/site-packages/sage/modular/quatalg/brandt.py : 79%

r""" 

Brandt Modules 

Introduction 

============ 

This tutorial outlines the construction of Brandt modules in Sage. The 

importance of this construction is that it provides us with a method 

to compute modular forms on `\Gamma_0(N)` as outlined in Pizer's paper 

[Piz1980]_. In fact there exists a non-canonical Hecke algebra isomorphism 

between the Brandt modules and a certain subspace of 

`S_{2}(\Gamma_0(pM))` which contains all the newforms. 

The Brandt module is the free abelian group on right ideal classes of 

a quaternion order together with a natural Hecke action determined by 

Brandt matrices. 

Quaternion Algebras 

------------------- 

A quaternion algebra over `\QQ` is a central simple algebra of 

dimension 4 over `\QQ`. Such an algebra `A` is said to be 

ramified at a place `v` of `\QQ` if and only if `A_v=A\otimes 

\QQ_v` is a division algebra. Otherwise `A` is said to be split 

at `v`. 

``A = QuaternionAlgebra(p)`` returns the quaternion algebra `A` over 

`\QQ` ramified precisely at the places `p` and `\infty`. 

``A = QuaternionAlgebra(k,a,b)`` returns a quaternion algebra with basis 

`\{1,i,j,j\}` over `\mathbb{K}` such that `i^2=a`, `j^2=b` and `ij=k.` 

An order `R` in a quaternion algebra is a 4-dimensional lattice on `A` 

which is also a subring containing the identity. 

``R = A.maximal_order()`` returns a maximal order `R` in the quaternion 

algebra `A.` 

An Eichler order `\mathcal{O}` in a quaternion algebra is the 

intersection of two maximal orders. The level of `\mathcal{O}` is its 

index in any maximal order containing it. 

``O = A.order_of_level_N`` returns an Eichler order `\mathcal{O}` in `A` 

of level `N` where `p` does not divide `N`. 

A right `\mathcal{O}`-ideal `I` is a lattice on `A` such that 

`I_p=a_p\mathcal{O}` (for some `a_p\in A_p^*`) for all `p<\infty`. Two 

right `\mathcal{O}`-ideals `I` and `J` are said to belong to the same 

class if `I=aJ` for some `a \in A^*`. (Left `\mathcal{O}`-ideals are 

defined in a similar fashion.) 

The right order of `I` is defined to be the set of elements in `A` 

which fix `I` under right multiplication. 

``right_order(R, basis)`` returns the right ideal of `I` in `R` given a 

basis for the right ideal `I` contained in the maximal order `R.` 

``ideal_classes(self)`` returns a tuple of all right ideal classes in self 

which, for the purpose of constructing the Brandt module B(p,M), is 

taken to be an Eichler order of level M. 

The implementation of this method is especially interesting. It 

depends on the construction of a Hecke module defined as a free 

abelian group on right ideal classes of a quaternion algebra with the 

following action 

.. MATH:: 

T_n[I] = \sum_{\phi} [J] 

where `(n,pM)=1` and the sum is over cyclic `\mathcal{O}`-module 

homomorphisms `\phi :I\rightarrow J` of degree `n` up to isomorphism 

of `J`. Equivalently one can sum over the inclusions of the submodules 

`J \rightarrow n^{-1}I`. The rough idea is to start with the trivial 

ideal class containing the order `\mathcal{O}` itself. Using the 

method ``cyclic_submodules(self, I, p)`` one computes `T_p([\mathcal{O}])` 

for some prime integer `p` not dividing the level of the order 

`\mathcal{O}`. Apply this method repeatedly and test for equivalence 

among resulting ideals. A theorem of Serre asserts that one gets a 

complete set of ideal class representatives after a finite number of 

repetitions. 

One can prove that two ideals `I` and `J` are equivalent if and only 

if there exists an element `\alpha \in I \overline{J}` such 

`N(\alpha)=N(I)N(J)`. 

``is_equivalent(I,J)`` returns true if `I` and `J` are equivalent. This 

method first compares the theta series of `I` and `J`. If they are the 

same, it computes the theta series of the lattice `I\overline(J)`. It 

returns true if the `n^{th}` coefficient of this series is nonzero 

where `n=N(J)N(I)`. 

The theta series of a lattice `L` over the quaternion algebra `A` is 

defined as 

.. MATH:: 

\theta_L(q)=\sum_{x \in L} q^{\frac{N(x)}{N(L)}} 

``L.theta_series(T,q)`` returns a power series representing `\theta_L(q)` 

up to a precision of `\mathcal{O}(q^{T+1})`. 

Hecke Structure 

--------------- 

The Hecke structure defined on the Brandt module is given by the 

Brandt matrices which can be computed using the definition of the 

Hecke operators given earlier. 

``hecke_matrix_from_defn(self,n)`` returns the matrix of the n-th Hecke 

operator `B_{0}(n)` acting on self, computed directly from the 

definition. 

However, one can efficiently compute Brandt matrices using theta 

series. In fact, let `\{I_{1},.....,I_{h}\}` be a set of right 

`\mathcal{O}`-ideal class representatives. The (i,j) entry in the 

Brandt matrix `B_{0}(n)` is the product of the `n^{th}` coefficient in 

the theta series of the lattice `I_{i}\overline{I_{j}}` and the first 

coefficient in the theta series of the lattice 

`I_{i}\overline{I_{i}}`. 

``compute_hecke_matrix_brandt(self,n)`` returns the n-th Hecke matrix, 

computed using theta series. 

EXAMPLES:: 

sage: B = BrandtModule(23) 

sage: B.maximal_order() 

Order of Quaternion Algebra (-1, -23) with base ring Rational Field with basis (1/2 + 1/2*j, 1/2*i + 1/2*k, j, k) 

sage: B.right_ideals() 

(Fractional ideal (2 + 2*j, 2*i + 2*k, 4*j, 4*k), Fractional ideal (2 + 2*j, 2*i + 6*k, 8*j, 8*k), Fractional ideal (2 + 10*j + 8*k, 2*i + 8*j + 6*k, 16*j, 16*k)) 

sage: B.hecke_matrix(2) 

[1 2 0] 

[1 1 1] 

[0 3 0] 

sage: B.brandt_series(3) 

[1/4 + q + q^2 + O(q^3) 1/4 + q^2 + O(q^3) 1/4 + O(q^3)] 

[ 1/2 + 2*q^2 + O(q^3) 1/2 + q + q^2 + O(q^3) 1/2 + 3*q^2 + O(q^3)] 

[ 1/6 + O(q^3) 1/6 + q^2 + O(q^3) 1/6 + q + O(q^3)] 

REFERENCES: 

- [Piz1980]_ 

- [Koh2000]_ 

Further Examples 

---------------- 

We decompose a Brandt module over both `\ZZ` and `\QQ`. :: 

sage: B = BrandtModule(43, base_ring=ZZ); B 

Brandt module of dimension 4 of level 43 of weight 2 over Integer Ring 

sage: D = B.decomposition() 

sage: D 

[ 

Subspace of dimension 1 of Brandt module of dimension 4 of level 43 of weight 2 over Integer Ring, 

Subspace of dimension 1 of Brandt module of dimension 4 of level 43 of weight 2 over Integer Ring, 

Subspace of dimension 2 of Brandt module of dimension 4 of level 43 of weight 2 over Integer Ring 

] 

sage: D[0].basis() 

((0, 0, 1, -1),) 

sage: D[1].basis() 

((1, 2, 2, 2),) 

sage: D[2].basis() 

((1, 1, -1, -1), (0, 2, -1, -1)) 

sage: B = BrandtModule(43, base_ring=QQ); B 

Brandt module of dimension 4 of level 43 of weight 2 over Rational Field 

sage: B.decomposition()[2].basis() 

((1, 0, -1/2, -1/2), (0, 1, -1/2, -1/2)) 

AUTHORS: 

- Jon Bober 

- Alia Hamieh 

- Victoria de Quehen 

- William Stein 

- Gonzalo Tornaria 

""" 

#***************************************************************************** 

# Copyright (C) 2009 William Stein <wstein@gmail.com> 

# 

# This program is free software: you can redistribute it and/or modify 

# it under the terms of the GNU General Public License as published by 

# the Free Software Foundation, either version 2 of the License, or 

# (at your option) any later version. 

# http://www.gnu.org/licenses/ 

#***************************************************************************** 

from __future__ import print_function 

# imports 

from sage.misc.all import prod, verbose 

from sage.rings.all import Integer, ZZ, QQ, PolynomialRing, GF, CommutativeRing 

from sage.algebras.quatalg.quaternion_algebra import QuaternionAlgebra, basis_for_quaternion_lattice 

from sage.algebras.quatalg.quaternion_algebra_cython import rational_matrix_from_rational_quaternions 

from sage.arith.all import gcd, factor, prime_divisors, kronecker, next_prime, lcm 

from sage.modular.hecke.all import (AmbientHeckeModule, HeckeSubmodule, HeckeModuleElement) 

from sage.modular.dirichlet import TrivialCharacter 

from sage.matrix.all import MatrixSpace, matrix 

from sage.misc.mrange import cartesian_product_iterator 

from sage.structure.richcmp import richcmp, richcmp_method 

from sage.misc.cachefunc import cached_method 

from copy import copy 

cache = {} 

def BrandtModule(N, M=1, weight=2, base_ring=QQ, use_cache=True): 

""" 

Return the Brandt module of given weight associated to the prime 

power `p^r` and integer `M`, where `p` and `M` are coprime. 

INPUT: 

- `N` -- a product of primes with odd exponents 

- `M` -- an integer coprime to `q` (default: 1) 

- ``weight`` -- an integer that is at least 2 (default: 2) 

- ``base_ring`` -- the base ring (default: ``QQ``) 

- ``use_cache`` -- whether to use the cache (default: ``True``) 

OUTPUT: 

a Brandt module 

EXAMPLES:: 

sage: BrandtModule(17) 

Brandt module of dimension 2 of level 17 of weight 2 over Rational Field 

sage: BrandtModule(17,15) 

Brandt module of dimension 32 of level 17*15 of weight 2 over Rational Field 

sage: BrandtModule(3,7) 

Brandt module of dimension 2 of level 3*7 of weight 2 over Rational Field 

sage: BrandtModule(3,weight=2) 

Brandt module of dimension 1 of level 3 of weight 2 over Rational Field 

sage: BrandtModule(11, base_ring=ZZ) 

Brandt module of dimension 2 of level 11 of weight 2 over Integer Ring 

sage: BrandtModule(11, base_ring=QQbar) 

Brandt module of dimension 2 of level 11 of weight 2 over Algebraic Field 

The ``use_cache`` option determines whether the Brandt module returned 

by this function is cached:: 

sage: BrandtModule(37) is BrandtModule(37) 

True 

sage: BrandtModule(37,use_cache=False) is BrandtModule(37,use_cache=False) 

False 

TESTS: 

Note that `N` and `M` must be coprime:: 

sage: BrandtModule(3,15) 

Traceback (most recent call last): 

... 

ValueError: M must be coprime to N 

Only weight 2 is currently implemented:: 

sage: BrandtModule(3,weight=4) 

Traceback (most recent call last): 

... 

NotImplementedError: weight != 2 not yet implemented 

Brandt modules are cached:: 

sage: B = BrandtModule(3,5,2,ZZ) 

sage: B is BrandtModule(3,5,2,ZZ) 

True 

""" 

N, M, weight = Integer(N), Integer(M), Integer(weight) 

if not N.is_prime(): 

raise NotImplementedError("Brandt modules currently only implemented when N is a prime") 

if M < 1: 

raise ValueError("M must be positive") 

if gcd(M,N) != 1: 

raise ValueError("M must be coprime to N") 

if weight < 2: 

raise ValueError("weight must be at least 2") 

if not isinstance(base_ring, CommutativeRing): 

raise TypeError("base_ring must be a commutative ring") 

key = (N, M, weight, base_ring) 

if use_cache: 

if key in cache: # TODO: re-enable caching! 

return cache[key] 

if weight != 2: 

raise NotImplementedError("weight != 2 not yet implemented") 

B = BrandtModule_class(*key) 

if use_cache: 

cache[key] = B 

return B 

def class_number(p, r, M): 

r""" 

Return the class number of an order of level `N = p^r M` in the 

quaternion algebra over `\QQ` ramified precisely at `p` and infinity. 

This is an implementation of Theorem 1.12 of [Piz1980]_. 

INPUT: 

- `p` -- a prime 

- `r` -- an odd positive integer (default: 1) 

- `M` -- an integer coprime to `q` (default: 1) 

OUTPUT: 

Integer 

EXAMPLES:: 

sage: sage.modular.quatalg.brandt.class_number(389,1,1) 

33 

sage: sage.modular.quatalg.brandt.class_number(389,1,2) # TODO -- right? 

97 

sage: sage.modular.quatalg.brandt.class_number(389,3,1) # TODO -- right? 

4892713 

""" 

N = M * p**r 

D = prime_divisors(M) 

s = 0; t = 0 

if N % 4 != 0: 

s = (1 - kronecker(-4,p))/4 * prod(1 + kronecker(-4,q) for q in D) 

if N % 9 != 0: 

t = (1 - kronecker(-3,p))/3 * prod(1 + kronecker(-3,q) for q in D) 

h = (N/Integer(12))*(1 - 1/p)*prod(1+1/q for q in D) + s + t 

return Integer(h) 

def maximal_order(A): 

""" 

Return a maximal order in the quaternion algebra ramified 

at `p` and infinity. 

This is an implementation of Proposition 5.2 of [Piz1980]_. 

INPUT: 

- `A` -- quaternion algebra ramified precisely at `p` and infinity 

OUTPUT: 

a maximal order in `A` 

EXAMPLES:: 

sage: A = BrandtModule(17).quaternion_algebra() 

sage: sage.modular.quatalg.brandt.maximal_order(A) 

Order of Quaternion Algebra (-17, -3) with base ring Rational Field with basis (1/2 + 1/2*j, 1/2*i + 1/2*k, -1/3*j - 1/3*k, k) 

sage: A = QuaternionAlgebra(17,names='i,j,k') 

sage: A.maximal_order() 

Order of Quaternion Algebra (-3, -17) with base ring Rational Field with basis (1/2 + 1/2*i, 1/2*j - 1/2*k, -1/3*i + 1/3*k, -k) 

""" 

return A.maximal_order() 

def basis_for_left_ideal(R, gens): 

""" 

Return a basis for the left ideal of `R` with given generators. 

INPUT: 

- `R` -- quaternion order 

- ``gens`` -- list of elements of `R` 

OUTPUT: 

list of four elements of `R` 

EXAMPLES:: 

sage: B = BrandtModule(17); A = B.quaternion_algebra(); i,j,k = A.gens() 

sage: sage.modular.quatalg.brandt.basis_for_left_ideal(B.maximal_order(), [i+j,i-j,2*k,A(3)]) 

[1/2 + 1/6*j + 2/3*k, 1/2*i + 1/2*k, 1/3*j + 1/3*k, k] 

sage: sage.modular.quatalg.brandt.basis_for_left_ideal(B.maximal_order(), [3*(i+j),3*(i-j),6*k,A(3)]) 

[3/2 + 1/2*j + 2*k, 3/2*i + 3/2*k, j + k, 3*k] 

""" 

return basis_for_quaternion_lattice([b*g for b in R.basis() for g in gens]) 

def right_order(R, basis): 

""" 

Given a basis for a left ideal `I`, return the right order in the 

quaternion order `R` of elements `x` such that `I x` is contained in `I`. 

INPUT: 

- `R` -- order in quaternion algebra 

- ``basis`` -- basis for an ideal `I` 

OUTPUT: 

order in quaternion algebra 

EXAMPLES: 

We do a consistency check with the ideal equal to a maximal order:: 

sage: B = BrandtModule(17); basis = sage.modular.quatalg.brandt.basis_for_left_ideal(B.maximal_order(), B.maximal_order().basis()) 

sage: sage.modular.quatalg.brandt.right_order(B.maximal_order(), basis) 

Order of Quaternion Algebra (-17, -3) with base ring Rational Field with basis (1/2 + 1/6*j + 2/3*k, 1/2*i + 1/2*k, 1/3*j + 1/3*k, k) 

sage: basis 

[1/2 + 1/6*j + 2/3*k, 1/2*i + 1/2*k, 1/3*j + 1/3*k, k] 

sage: B = BrandtModule(17); A = B.quaternion_algebra(); i,j,k = A.gens() 

sage: basis = sage.modular.quatalg.brandt.basis_for_left_ideal(B.maximal_order(), [i*j-j]) 

sage: sage.modular.quatalg.brandt.right_order(B.maximal_order(), basis) 

Order of Quaternion Algebra (-17, -3) with base ring Rational Field with basis (1/2 + 1/2*i + 1/2*j + 17/2*k, i, j + 8*k, 9*k) 

""" 

# Compute matrix of multiplication by each element of the basis. 

B = R.basis() 

Z = R.quaternion_algebra() 

M = MatrixSpace(QQ, 4) 

# I = matrix with rows the given basis for I 

I = M([list(f) for f in basis]) 

# psi = matrix of right multiplication on each basis element 

psi = [M([list(f*x) for x in Z.basis()]) for f in basis] 

# invert them 

psi_inv = [x**(-1) for x in psi] 

# apply the four inverses to I 

W = [I*x for x in psi_inv] 

# The right order is the intersection of the row span of the W with the row span of B. 

X = M([list(b) for b in B]).row_module(ZZ) 

for A in W: 

X = X.intersection(A.row_module(ZZ)) 

C = [Z(list(b)) for b in X.basis()] 

return Z.quaternion_order(C) 

def quaternion_order_with_given_level(A, level): 

""" 

Return an order in the quaternion algebra A with given level. 

(Implemented only when the base field is the rational numbers.) 

INPUT: 

- ``level`` -- The level of the order to be returned. Currently this 

is only implemented when the level is divisible by at 

most one power of a prime that ramifies in this 

quaternion algebra. 

EXAMPLES:: 

sage: from sage.modular.quatalg.brandt import quaternion_order_with_given_level, maximal_order 

sage: A.<i,j,k> = QuaternionAlgebra(5) 

sage: level = 2 * 5 * 17 

sage: O = quaternion_order_with_given_level(A, level) 

sage: M = maximal_order(A) 

sage: L = O.free_module() 

sage: N = M.free_module() 

sage: L.index_in(N) == level/5 #check that the order has the right index in the maximal order 

True 

""" 

if A.base_ring() is not QQ: 

raise NotImplementedError("base field must be rational numbers") 

from sage.modular.quatalg.brandt import maximal_order 

if len(A.ramified_primes()) > 1: 

raise NotImplementedError("Currently this algorithm only works when the quaternion algebra is only ramified at one finite prime.") 

# (The algorithm we use is similar to that in Magma (by David Kohel).) 

# in the following magma code, M denotes is the level 

level = abs(level) 

N = A.discriminant() 

N1 = gcd(level, N) 

M1 = level/N1 

O = maximal_order(A) 

if 0 and N1 != 1: # we don't know why magma does the following, so we don't do it. 

for p in A.ramified_primes(): 

if level % p**2 == 0: 

raise NotImplementedError("Currently sage can only compute orders whose level is divisible by at most one power of any prime that ramifies in the quaternion algebra") 

P = basis_for_left_ideal(O, [N1] + [x*y - y*x for x, y in cartesian_product_iterator([A.basis(), A.basis()]) ]) 

O = A.quaternion_order(P) 

fact = factor(M1) 

B = O.basis() 

for (p, r) in fact: 

a = int(-p) // 2 

for v in GF(p)**4: 

x = sum([int(v[i]+a)*B[i] for i in range(4)]) 

D = x.reduced_trace()**2 - 4 * x.reduced_norm() 

#x = O.random_element((-p/2).floor(), (p/2).ceil()) 

if kronecker(D, p) == 1: break 

X = PolynomialRing(GF(p), 'x').gen() 

a = ZZ((X**2 - ZZ(x.reduced_trace()) * X + ZZ(x.reduced_norm())).roots()[0][0]) 

I = basis_for_left_ideal(O, [p**r, (x-a)**r] ) 

O = right_order(O, I) # right_order returns the RightOrder of I inside O, so we don't need to do another intersection 

return O 

class BrandtSubmodule(HeckeSubmodule): 

def _repr_(self): 

""" 

Return string representation of this Brandt submodule. 

EXAMPLES:: 

sage: BrandtModule(11)[0]._repr_() 

'Subspace of dimension 1 of Brandt module of dimension 2 of level 11 of weight 2 over Rational Field' 

""" 

return "Subspace of dimension %s of %s" % (self.dimension(), self.ambient_module()) 

class BrandtModuleElement(HeckeModuleElement): 

def __init__(self, parent, x): 

""" 

EXAMPLES:: 

sage: B = BrandtModule(37) 

sage: x = B([1,2,3]); x 

(1, 2, 3) 

sage: parent(x) 

Brandt module of dimension 3 of level 37 of weight 2 over Rational Field 

""" 

if isinstance(x, HeckeModuleElement): 

x = x.element() 

HeckeModuleElement.__init__(self, parent, parent.free_module()(x)) 

def _richcmp_(self, other, op): 

""" 

EXAMPLES:: 

sage: B = BrandtModule(13,5) 

sage: B.0 

(1, 0, 0, 0, 0, 0) 

sage: B.0 == B.1 

False 

sage: B.0 == 0 

False 

sage: B(0) == 0 

True 

sage: B.0 + 2*B.1 == 2*B.1 + B.0 

True 

sage: loads(dumps(B.0)) == B.0 

True 

""" 

return richcmp(self.element(), other.element(), op) 

def monodromy_pairing(self, x): 

""" 

Return the monodromy pairing of ``self`` and ``x``. 

EXAMPLES:: 

sage: B = BrandtModule(5,13) 

sage: B.monodromy_weights() 

(1, 3, 1, 1, 1, 3) 

sage: (B.0 + B.1).monodromy_pairing(B.0 + B.1) 

4 

TESTS: 

One check for :trac:`12866`:: 

sage: Br = BrandtModule(2,7) 

sage: g1, g2 = Br.basis() 

sage: g = g1 - g2 

sage: g.monodromy_pairing(g) 

6 

""" 

B = self.parent() 

w = B.monodromy_weights() 

x = B(x).element() 

v = self.element() 

return sum(x[i] * v[i] * w[i] for i in range(len(v))) 

def __mul__(self, right): 

""" 

Return the monodromy pairing of ``self`` and ``right``. 

EXAMPLES:: 

sage: B = BrandtModule(7,10) 

sage: B.monodromy_weights() 

(1, 1, 1, 2, 1, 1, 2, 1, 1, 1) 

sage: B.0 * B.0 

1 

sage: B.3 * B.3 

2 

sage: (B.0+B.3) * (B.0 + B.1 + 2*B.3) 

5 

""" 

return self.monodromy_pairing(right) 

def _add_(self, right): 

""" 

Return the sum of ``self`` and ``right``. 

EXAMPLES:: 

sage: B = BrandtModule(11) 

sage: B.0 + B.1 # indirect doctest 

(1, 1) 

""" 

return BrandtModuleElement(self.parent(), self.element() + right.element()) 

def _sub_(self, right): 

""" 

Return the difference of ``self`` and ``right``. 

EXAMPLES:: 

sage: B = BrandtModule(11) 

sage: B.0 - B.1 # indirect doctest 

(1, -1) 

""" 

return BrandtModuleElement(self.parent(), self.element() - right.element()) 

def _neg_(self): 

""" 

Return the opposite of ``self``. 

EXAMPLES:: 

sage: B = BrandtModule(11) 

sage: -B.0 # indirect doctest 

(-1, 0) 

""" 

return BrandtModuleElement(self.parent(), -self.element()) 

@richcmp_method 

class BrandtModule_class(AmbientHeckeModule): 

""" 

A Brandt module. 

EXAMPLES:: 

sage: BrandtModule(3, 10) 

Brandt module of dimension 4 of level 3*10 of weight 2 over Rational Field 

""" 

def __init__(self, N, M, weight, base_ring): 

""" 

INPUT: 

- N -- ramification number (coprime to M) 

- M -- auxiliary level 

- weight -- integer 2 

- base_ring -- the base ring 

EXAMPLES:: 

sage: BrandtModule(3, 5, weight=2, base_ring=ZZ) 

Brandt module of dimension 2 of level 3*5 of weight 2 over Integer Ring 

""" 

assert weight == 2 

self.__N = N 

self.__M = M 

if not N.is_prime(): 

raise NotImplementedError("right now N must be prime") 

rank = class_number(N, 1, M) 

self.__key = (N, M, weight, base_ring) 

AmbientHeckeModule.__init__(self, base_ring, rank, N * M, weight=2) 

self._populate_coercion_lists_(coerce_list=[self.free_module()]) 

Element = BrandtModuleElement 

def _submodule_class(self): 

""" 

Return the Python class of submodules of this ambient Brandt module. 

EXAMPLES:: 

sage: BrandtModule(37)._submodule_class() 

<class 'sage.modular.quatalg.brandt.BrandtSubmodule'> 

""" 

return BrandtSubmodule 

@cached_method 

def free_module(self): 

""" 

Return the underlying free module of the Brandt module. 

EXAMPLES:: 

sage: B = BrandtModule(10007,389) 

sage: B.free_module() 

Vector space of dimension 325196 over Rational Field 

""" 

return self.base_ring() ** self.dimension() 

def N(self): 

""" 

Return ramification level `N`. 

EXAMPLES:: 

sage: BrandtModule(7,5,2,ZZ).N() 

7 

""" 

return self.__N 

def M(self): 

""" 

Return the auxiliary level (prime to `p` part) of the quaternion 

order used to compute this Brandt module. 

EXAMPLES:: 

sage: BrandtModule(7,5,2,ZZ).M() 

5 

""" 

return self.__M 

def character(self): 

r""" 

The character of this space. 

Always trivial. 

EXAMPLES:: 

sage: BrandtModule(11,5).character() 

Dirichlet character modulo 55 of conductor 1 mapping 12 |--> 1, 46 |--> 1 

""" 

return TrivialCharacter(self.__N * self.__M) 

def _repr_(self): 

""" 

Return string representation of this Brandt module. 

EXAMPLES:: 

sage: BrandtModule(7,5,2,ZZ)._repr_() 

'Brandt module of dimension 4 of level 7*5 of weight 2 over Integer Ring' 

""" 

aux = '' if self.__M == 1 else '*%s'%self.__M 

return "Brandt module of dimension %s of level %s%s of weight %s over %s"%( 

self.rank(), self.__N, aux, self.weight(), self.base_ring()) 

def __richcmp__(self, other, op): 

r""" 

Compare ``self`` to ``other``. 

EXAMPLES:: 

sage: BrandtModule(37, 5, 2, ZZ) == BrandtModule(37, 5, 2, QQ) 

False 

sage: BrandtModule(37, 5, 2, ZZ) == BrandtModule(37, 5, 2, ZZ) 

True 

sage: BrandtModule(37, 5, 2, ZZ) == loads(dumps(BrandtModule(37, 5, 2, ZZ))) 

True 

""" 

if not isinstance(other, BrandtModule_class): 

return NotImplemented 

return richcmp((self.__M, self.__N, self.weight(), self.base_ring()), 

(other.__M, other.__N, other.weight(), other.base_ring()), 

op) 

@cached_method 

def quaternion_algebra(self): 

r""" 

Return the quaternion algebra `A` over `\QQ` ramified precisely at 

`p` and infinity used to compute this Brandt module. 

EXAMPLES:: 

sage: BrandtModule(997).quaternion_algebra() 

Quaternion Algebra (-2, -997) with base ring Rational Field 

sage: BrandtModule(2).quaternion_algebra() 

Quaternion Algebra (-1, -1) with base ring Rational Field 

sage: BrandtModule(3).quaternion_algebra() 

Quaternion Algebra (-1, -3) with base ring Rational Field 

sage: BrandtModule(5).quaternion_algebra() 

Quaternion Algebra (-2, -5) with base ring Rational Field 

sage: BrandtModule(17).quaternion_algebra() 

Quaternion Algebra (-17, -3) with base ring Rational Field 

""" 

p = self.N() 

assert p.is_prime(), "we have only implemented the prime case" 

if p == 2: 

QA = -1 

QB = -1 

elif p % 4 == 3: 

QA = -1 

QB = -p 

elif p % 8 == 5: 

QA = -2 

QB = -p 

elif p % 8 == 1: 

q = 3 

while q % 4 != 3 or kronecker(p, q) != -1: 

q = next_prime(q) 

QA = -p 

QB = -q 

return QuaternionAlgebra(QQ, QA, QB) 

@cached_method 

def maximal_order(self): 

""" 

Return a maximal order in the quaternion algebra associated to this Brandt module. 

EXAMPLES:: 

sage: BrandtModule(17).maximal_order() 

Order of Quaternion Algebra (-17, -3) with base ring Rational Field with basis (1/2 + 1/2*j, 1/2*i + 1/2*k, -1/3*j - 1/3*k, k) 

sage: BrandtModule(17).maximal_order() is BrandtModule(17).maximal_order() 

True 

""" 

return maximal_order(self.quaternion_algebra()) 

@cached_method 

def order_of_level_N(self): 

""" 

Return Eichler order of level `N = p^{2 r + 1} M` in the 

quaternion algebra. 

EXAMPLES:: 

sage: BrandtModule(7).order_of_level_N() 

Order of Quaternion Algebra (-1, -7) with base ring Rational Field with basis (1/2 + 1/2*j, 1/2*i + 1/2*k, j, k) 

sage: BrandtModule(7,13).order_of_level_N() 

Order of Quaternion Algebra (-1, -7) with base ring Rational Field with basis (1/2 + 1/2*j + 12*k, 1/2*i + 9/2*k, j + 11*k, 13*k) 

sage: BrandtModule(7,3*17).order_of_level_N() 

Order of Quaternion Algebra (-1, -7) with base ring Rational Field with basis (1/2 + 1/2*j + 35*k, 1/2*i + 65/2*k, j + 19*k, 51*k) 

""" 

return quaternion_order_with_given_level(self.quaternion_algebra(), self.level()) 

def cyclic_submodules(self, I, p): 

""" 

Return a list of rescaled versions of the fractional right 

ideals `J` such that `J` contains `I` and the quotient has 

group structure the product of two cyclic groups of order `p`. 

We emphasize again that `J` is rescaled to be integral. 

INPUT: 

- `I` -- ideal I in R = self.order_of_level_N() 

- `p` -- prime `p` coprime to self.level() 

OUTPUT: 

list of the `p+1` fractional right R-ideals that contain I 

such that J/I is GF(p) x GF(p). 

EXAMPLES:: 

sage: B = BrandtModule(11) 

sage: I = B.order_of_level_N().unit_ideal() 

sage: B.cyclic_submodules(I, 2) 

[Fractional ideal (1/2 + 3/2*j + k, 1/2*i + j + 1/2*k, 2*j, 2*k), 

Fractional ideal (1/2 + 1/2*i + 1/2*j + 1/2*k, i + k, j + k, 2*k), 

Fractional ideal (1/2 + 1/2*j + k, 1/2*i + j + 3/2*k, 2*j, 2*k)] 

sage: B.cyclic_submodules(I, 3) 

[Fractional ideal (1/2 + 1/2*j, 1/2*i + 5/2*k, 3*j, 3*k), 

Fractional ideal (1/2 + 3/2*j + 2*k, 1/2*i + 2*j + 3/2*k, 3*j, 3*k), 

Fractional ideal (1/2 + 3/2*j + k, 1/2*i + j + 3/2*k, 3*j, 3*k), 

Fractional ideal (1/2 + 5/2*j, 1/2*i + 1/2*k, 3*j, 3*k)] 

sage: B.cyclic_submodules(I, 11) 

Traceback (most recent call last): 

... 

ValueError: p must be coprime to the level 

""" 

if not Integer(p).is_prime(): 

raise ValueError("p must be a prime") 

if self.level() % p == 0: 

raise ValueError("p must be coprime to the level") 

R = self.order_of_level_N() 

A = R.quaternion_algebra() 

B = R.basis() 

V = GF(p)**4 

# step 1: Compute alpha, beta, and the matrix of their action on I/pI. 

# NOTE: Move this code to orders once we have it all working... 

try: 

alpha, beta = self.__cyclic_submodules[p] 

compute = False 

except AttributeError: 

self.__cyclic_submodules = {} 

compute = True 

except KeyError: 

compute = True 

if compute: 

d = R.free_module().basis_matrix().determinant() 

S = None 

for v in V: 

if not v: continue 

alpha = sum(Integer(v[i])*B[i] for i in range(4)) 

# If the quadratic polynomial over GF(p) given by 

# X^2 - alpha.reduced_trace() * X + alpha.reduced_norm() 

# is not irreducible, we try again with a new element. 

if p == 2: 

# special case p == 2, since there is a unique quadratic irreducible poly. 

if alpha.reduced_trace()%2 == 0 or alpha.reduced_norm()%2 == 0: 

continue 

else: 

# check if the discriminant is a square -- if so, poly is reducible 

b = alpha.reduced_trace(); c = alpha.reduced_norm() 

if kronecker(b*b - 4*c, p) != -1: 

continue 

for w in V: 

if not w: continue 

beta = sum(Integer(w[i])*B[i] for i in range(4)) 

v = [A(1), alpha, beta, alpha*beta] 

M = rational_matrix_from_rational_quaternions(v) 

e = M.determinant() 

if e and (d/e).valuation(p) == 0: 

S = A.quaternion_order(v) 

break 

if S is not None: 

break 

self.__cyclic_submodules[p] = (alpha, beta) 

# right multiplication by X changes something to be written 

# in terms of the basis for I. 

Y = I.basis_matrix() 

X = Y**(-1) 

# Compute the matrix of right multiplication by alpha acting on 

# our fixed choice of basis for this ideal. 

M_alpha = (matrix([(i*alpha).coefficient_tuple() for i in I.basis()]) * X).change_ring(GF(p)).change_ring(GF(p)) 

M_beta = (matrix([(i*beta).coefficient_tuple() for i in I.basis()]) * X).change_ring(GF(p)).change_ring(GF(p)) 

# step 2: Find j such that if f=I[j], then mod 2 we have span(I[0],alpha*I[i]) 

# has trivial intersection with span(I[j],alpha*I[j]). 

# 

# In terms of our matrices alpha, beta, we can now think of I/p*I 

# as being the GF(p)^4 that M_alpha and M_beta naturally act on, 

# and I[0], I[1], I[2], I[3] correspond to the standard basis. 

# 

# We try each of the standard basis vectors. 

W0 = V.span([V.gen(0), V.gen(0)*M_alpha]) 

assert W0.dimension() == 2 

j = None