Pari/GP Reference Documentation | Contents - Global index - GP keyboard shortcuts |
msatkinlehner mscosets mscuspidal msdim mseisenstein mseval msfarey msfromcusp msfromell msfromhecke msgetlevel msgetsign msgetweight mshecke msinit msissymbol mslattice msnew msomseval mspadicL mspadicinit mspadicmoments mspadicseries mspathgens mspathlog mspetersson mspolygon msqexpansion mssplit msstar mstooms | |
Let Δ0 := Div0(ℙ1(ℚ)) be the abelian group of divisors of degree 0 on the rational projective line. The standard GL(2,ℚ) action on ℙ1(ℚ) via homographies naturally extends to Δ0. Given * G a finite index subgroup of SL(2,ℤ), * a field F and a finite dimensional representation V/F of GL(2,ℚ), we consider the space of modular symbols M := HomG(Δ0, V). This finite dimensional F-vector space is a G-module, canonically isomorphic to H1c(X(G), V), and allows to compute modular forms for G.
Currently, we only support the groups Γ0(N) (N > 0 an integer)
and the representations Vk = ℚ[X,Y]k-2 (k ≥ 2 an integer) over
ℚ. We represent a space of modular symbols by an ms structure,
created by the function A subspace of M (e.g. the cuspidal or Eisenstein subspaces, the new or old modular symbols, etc.) is given by a structure allowing quick projection and restriction of linear operators; its first component is a matrix whose columns form an F-basis of the subspace.
| |
msatkinlehner(M, Q, {H}) | |
Let M be a full modular symbol space of level N,
as given by
? M = msinit(36,2); \\ M2(Gamma0(36)) ? w = msatkinlehner(M,4); w^2 == 1 %2 = 1 ? #w \\ involution acts on a 13-dimensional space %3 = 13 ? M = msinit(36,2, -1); \\ M2(Gamma0(36))^- ? w = msatkinlehner(M,4); w^2 == 1 %5 = 1 ? #w %6 = 4
The library syntax is
| |
mscosets(gen, inH) | |
Let
? PSL2 = [[0,1;-1,0], [1,1;0,1]]; \\ S and T \\ G = PSL2, H = Gamma0(2) ? [C, M] = mscosets(PSL2, g->g[2,1] % 2 == 0); ? C \\ three cosets %3 = [[1, 0; 0, 1], [0, 1; -1, 0], [0, 1; -1, -1]] ? M %4 = [Vecsmall([2, 1]), Vecsmall([1, 3]), Vecsmall([3, 2])] Looking at M[1] we see that S belongs to the second coset and T to the first (trivial) coset.
The library syntax is
| |
mscuspidal(M, {flag = 0}) | |
M being a full modular symbol space, as given by A subspace is given by a structure allowing quick projection and restriction of linear operators; its first component is a matrix with integer coefficients whose columns form a ℚ-basis of the subspace.
? M = msinit(2,8, 1); \\ M8(Gamma0(2))^+ ? [S,E] = mscuspidal(M, 1); ? E[1] \\ 2-dimensional %3 = [0 -10] [0 -15] [0 -3] [1 0] ? S[1] \\ 1-dimensional %4 = [ 3] [30] [ 6] [-8]
The library syntax is
| |
msdim(M) | |
M being a full modular symbol space or subspace, for instance
as given by
? M = msinit(11,4); msdim(M) %1 = 6 ? M = msinit(11,4,1); msdim(M) %2 = 4 \\ dimension of the '+' part ? [S,E] = mscuspidal(M,1); ? [msdim(S), msdim(E)] %4 = [2, 2]
Note that
The library syntax is
| |
mseisenstein(M) | |
M being a full modular symbol space, as given by
? M = msinit(2,8, 1); \\ M8(Gamma0(2))^+ ? E = mseisenstein(M); ? E[1] \\ 2-dimensional %3 = [0 -10] [0 -15] [0 -3] [1 0] ? E == mscuspidal(M,1)[2] %4 = 1
The library syntax is
| |
mseval(M, s, {p}) | |
Let Δ0 := Div0(ℙ1 (ℚ)).
Let M be a full modular symbol space, as given by
* a
* a
* a
If p is omitted, convert a single symbol s to the second form: a vector
of the s(gi). A
? M = msinit(2,8,1); \\ M8(Gamma0(2))^+ ? g = mspathgens(M)[1] %2 = [[+oo, 0], [0, 1]] ? N = msnew(M)[1]; #N \\ Q-basis of new subspace, dimension 1 %3 = 1 ? s = N[,1] \\ t_COL representation %4 = [-3, 6, -8]~ ? S = mseval(M, s) \\ t_VEC representation %5 = [64*x^6-272*x^4+136*x^2-8, 384*x^5+960*x^4+192*x^3-672*x^2-432*x-72] ? mseval(M,s, g[1]) %6 = 64*x^6 - 272*x^4 + 136*x^2 - 8 ? mseval(M,S, g[1]) %7 = 64*x^6 - 272*x^4 + 136*x^2 - 8 Note that the symbol should have values in V = ℚ[x,y]k-2, we return the de-homogenized values corresponding to y = 1 instead.
The library syntax is
| |
msfarey(F, inH, {&CM}) | |
F being a Farey symbol attached to a group G contained in
PSL2(ℤ) and H a subgroup of G, return a Farey symbol attached
to H. The subgroup H is given by a criterion
*
* or
If present, the argument
\\ Gamma0(N) G0(N) = mspolygon(N); \\ Gamma1(N): direct construction, slow G1(N) = msfarey(mspolygon(1), g -> my(a = g[1,1]%N, c = g[2,1]%N);\ c == 0 && (a == 1 || a == N-1)); \\ Gamma1(N) via Gamma0(N): much faster G1(N) = msfarey(G0(N), g -> my(a=g[1,1]%N); a==1 || a==N-1);
Note that the simpler criterion
\\ Gamma(N) G(N) = msfarey(G1(N), g -> g[1,2]%N==0); G_00(N) = msfarey(G0(N), x -> x[1,2]%N==0); G10(N1,N2) = msfarey(G0(1), x -> x[2,1]%N1==0 && x[1,2]%N2==0); \\ Gamma0(91) has 4 elliptic points of order 3, Gamma1(91) has none D0 = mspolygon(G0(91), 2)[4]; D1 = mspolygon(G1(91), 2)[4]; write("F.tex","\\documentclass{article}\\usepackage{tikz}\\begin{document}",\ D0,"\n",D1,"\\end{document}");
The library syntax is
| |
msfromcusp(M, c) | |
Returns the modular symbol attached to the cusp
c, where M is a modular symbol space of level N, attached to
G = Γ0(N). The cusp c in ℙ1(ℚ)/G is given either as
? M = msinit(2,8); \\ M8(Gamma0(2)) ? E = mseisenstein(M); ? E[1] \\ two-dimensional %3 = [0 -10] [0 -15] [0 -3] [1 0] ? s = msfromcusp(M,oo) %4 = [0, 0, 0, 1]~ ? mseval(M, s) %5 = [1, 0] ? s = msfromcusp(M,1) %6 = [-5/16, -15/32, -3/32, 0]~ ? mseval(M,s) %7 = [-x^6, -6*x^5 - 15*x^4 - 20*x^3 - 15*x^2 - 6*x - 1] In case M was initialized with a nonzero sign, the symbol is given in terms of the fixed basis of the whole symbol space, not the + or - part (to which it need not belong).
? M = msinit(2,8, 1); \\ M8(Gamma0(2))^+ ? E = mseisenstein(M); ? E[1] \\ still two-dimensional, in a smaller space %3 = [ 0 -10] [ 0 3] [-1 0] ? s = msfromcusp(M,oo) \\ in terms of the basis for M8(Gamma0(2)) ! %4 = [0, 0, 0, 1]~ ? mseval(M, s) \\ same symbol as before %5 = [1, 0]
The library syntax is
| |
msfromell(E, {sign = 0}) | |
Let E/ℚ be an elliptic curve of conductor N. For ϵ = ±1, we define the (cuspidal, new) modular symbol xϵ in H1c(X0(N),ℚ)ϵ attached to E. For all primes p not dividing N we have Tp(xϵ) = ap xϵ, where ap = p+1-#E(𝔽p).
Let Ω+ =
This function returns the pair [M, x], where M is
The modular symbols x± are given as a
? E=ellinit([0,-1,1,-10,-20]); \\ X0(11) ? [M,xp]= msfromell(E,1); ? xp %3 = [1/5, -1/2, -1/2]~ ? [M,x]= msfromell(E); ? x \\ x^+, x^- and LE %5 = [[1/5, -1/2, -1/2]~, [0, 1/2, -1/2]~, [1/5, 0; -1, 1; 0, -1]] ? p = 23; (mshecke(M,p) - ellap(E,p))*x[1] %6 = [0, 0, 0]~ \\ true at all primes, including p = 11; same for x[2] ? (mshecke(M,p) - ellap(E,p))*x[3] == 0 %7 = 1 Instead of a single curve E, one may use instead a vector of isogenous curves. The function then returns M and the vector of attached modular symbols.
The library syntax is
| |
msfromhecke(M, v, {H}) | |
Given a msinit M and a vector v of pairs [p, P] (where p is prime and P is a polynomial with integer coefficients), return a basis of all modular symbols such that P(Tp)(s) = 0. If H is present, it must be a Hecke-stable subspace and we restrict to s ∈ H. When Tp has a rational eigenvalue and P(x) = x-ap has degree 1, we also accept the integer ap instead of P.
? E = ellinit([0,-1,1,-10,-20]) \\11a1 ? ellap(E,2) %2 = -2 ? ellap(E,3) %3 = -1 ? M = msinit(11,2); ? S = msfromhecke(M, [[2,-2],[3,-1]]) %5 = [ 1 1] [-5 0] [ 0 -5] ? mshecke(M, 2, S) %6 = [-2 0] [ 0 -2] ? M = msinit(23,4); ? S = msfromhecke(M, [[5, x^4-14*x^3-244*x^2+4832*x-19904]]); ? factor( charpoly(mshecke(M,5,S)) ) %9 = [x^4 - 14*x^3 - 244*x^2 + 4832*x - 19904 2]
The library syntax is
| |
msgetlevel(M) | |
M being a full modular symbol space, as given by
The library syntax is
| |
msgetsign(M) | |
M being a full modular symbol space, as given by
? M = msinit(11,4, 1); ? msgetsign(M) %2 = 1 ? M = msinit(11,4); ? msgetsign(M) %4 = 0
The library syntax is
| |
msgetweight(M) | |
M being a full modular symbol space, as given by
? M = msinit(11,4); ? msgetweight(M) %2 = 4
The library syntax is
| |
mshecke(M, p, {H}) | |
M being a full modular symbol space, as given by
? M = msinit(11,2); \\ M2(Gamma0(11)) ? T2 = mshecke(M,2) %2 = [3 0 0] [1 -2 0] [1 0 -2] ? M = msinit(11,2, 1); \\ M2(Gamma0(11))^+ ? T2 = mshecke(M,2) %4 = [ 3 0] [-1 -2] ? N = msnew(M)[1] \\ Q-basis of new cuspidal subspace %5 = [-2] [-5] ? p = 1009; mshecke(M, p, N) \\ action of T_1009 on N %6 = [-10] ? ellap(ellinit("11a1"), p) %7 = -10
The library syntax is
| |
msinit(G, V, {sign = 0}) | |
Given G a finite index subgroup of SL(2,ℤ)
and a finite dimensional representation V of GL(2,ℚ), creates a
space of modular symbols, the G-module
HomG(Div0(ℙ1(ℚ)), V).
This is canonically isomorphic to H1c(X(G), V), and allows to
compute modular forms for G. If sign is present and nonzero, it
must be ±1 and we consider the subspace defined by Ker (σ -
sign), where σ is induced by
? M = msinit(11,2); msdim(M) \\ Gamma0(11), weight 2 %1 = 3 ? mshecke(M,2) \\ T2 acting on M %2 = [3 1 1] [0 -2 0] [0 0 -2] ? msstar(M) \\ * involution %3 = [1 0 0] [0 0 1] [0 1 0] ? Mp = msinit(11,2, 1); msdim(Mp) \\ + part %4 = 2 ? mshecke(Mp,2) \\ T2 action on M^+ %5 = [3 2] [0 -2] ? msstar(Mp) %6 = [1 0] [0 1]
The library syntax is
| |
msissymbol(M, s) | |
M being a full modular symbol space, as given by
? M = msinit(7,8, 1); \\ M8(Gamma0(7))^+ ? A = msnew(M)[1]; ? s = A[,1]; ? msissymbol(M, s) %4 = 1 ? msissymbol(M, A) %5 = [1, 1, 1] ? S = mseval(M,s); ? msissymbol(M, S) %7 = 1 ? [g,R] = mspathgens(M); g %8 = [[+oo, 0], [0, 1/2], [1/2, 1]] ? #R \\ 3 relations among the generators gi %9 = 3 ? T = S; T[3]++; \\ randomly perturb S(g3) ? msissymbol(M, T) %11 = 0 \\ no longer satisfies the relations
The library syntax is
| |
mslattice(M, {H}) | |
Let Δ0 := Div0(ℙ1(ℚ)) and Vk = ℚ[x,y]k-2.
Let M be a full modular symbol space, as given by
For user convenience, H can be defined by a matrix representing the
ℚ-basis of H (in terms of the canonical ℚ-basis of M fixed by
If omitted, H is the cuspidal part of M as given by
? M=msinit(11,2); ? [S,E] = mscuspidal(M,1); S[1] \\ a primitive Q-basis of S %2 = [ 1 1] [-5 0] [ 0 -5] ? mslattice(M,S) %3 = [-1/5 -1/5] [ 1 0] [ 0 1] ? mslattice(M,E) %4 = [1] [0] [0] ? M=msinit(5,4); ? S=mscuspidal(M); S[1] %6 = [ 7 20] [ 3 3] [-10 -23] [-30 -30] ? mslattice(M,S) %7 = [-1/10 -11/130] [ 0 -1/130] [ 1/10 6/65] [ 0 1/13]
The library syntax is
| |
msnew(M) | |
M being a full modular symbol space, as given by
? M = msinit(11,8, 1); \\ M8(Gamma0(11))^+ ? N = msnew(M); ? #N[1] \\ 6-dimensional %3 = 6
The library syntax is
| |
msomseval(Mp, PHI, path) | |
Return the vectors of moments of the p-adic distribution attached
to the path
? M = msinit(3,6,1); ? Mp= mspadicinit(M,5,10); ? phi = [5,-3,-1]~; ? msissymbol(M,phi) %4 = 1 ? PHI = mstooms(Mp,phi); ? ME = msomseval(Mp,PHI,[oo, 0]);
The library syntax is
| |
mspadicL(mu, {s = 0}, {r = 0}) | |
Returns the value (or r-th derivative)
on a character χs of ℤp* of the p-adic L-function
attached to
Let Φ be the p-adic distribution-valued overconvergent symbol
attached to a modular symbol φ for Γ0(N) (eigenvector for
TN(p) for the eigenvalue ap).
Then Lp(Φ,χs) = Lp(μ,s) is the
p-adic L function defined by
Lp(Φ,χs) = ∫ℤ_{p*} χs(z) dμ(z)
where μ is the distribution on ℤp* defined by the restriction of
Φ([ oo ]-[0]) to ℤp*. The r-th derivative is taken in
direction
*
* s = [s1,s2] with s1 ∈ ℤ ⊂ ℤp and
s2 mod p-1 or
s2 mod 2 for p = 2, encoding the p-adic character χs :=
When ap is a p-adic unit, Lp takes its values in ℚp. When ap is not a unit, it takes its values in the two-dimensional ℚp-vector space Dcris(M(φ)) where M(φ) is the "motive" attached to φ, and we return the two p-adic components with respect to some fixed ℚp-basis.
? M = msinit(3,6,1); phi=[5, -3, -1]~; ? msissymbol(M,phi) %2 = 1 ? Mp = mspadicinit(M, 5, 4); ? mu = mspadicmoments(Mp, phi); \\ no twist \\ End of initializations ? mspadicL(mu,0) \\ Lp(chi^0) %5 = 5 + 2*5^2 + 2*5^3 + 2*5^4 + ... ? mspadicL(mu,1) \\ Lp(chi), zero for parity reasons %6 = [O(5^13)]~ ? mspadicL(mu,2) \\ Lp(chi^2) %7 = 3 + 4*5 + 4*5^2 + 3*5^5 + ... ? mspadicL(mu,[0,2]) \\ Lp(tau^2) %8 = 3 + 5 + 2*5^2 + 2*5^3 + ... ? mspadicL(mu, [1,0]) \\ Lp(<chi>) %9 = 3*5 + 2*5^2 + 5^3 + 2*5^7 + 5^8 + 5^10 + 2*5^11 + O(5^13) ? mspadicL(mu,0,1) \\ Lp'(chi^0) %10 = 2*5 + 4*5^2 + 3*5^3 + ... ? mspadicL(mu, 2, 1) \\ Lp'(chi^2) %11 = 4*5 + 3*5^2 + 5^3 + 5^4 + ...
Now several quadratic twists:
? PHI = mstooms(Mp,phi); ? mu = mspadicmoments(Mp, PHI, 12); \\ twist by 12 ? mspadicL(mu) %14 = 5 + 5^2 + 5^3 + 2*5^4 + ... ? mu = mspadicmoments(Mp, PHI, 8); \\ twist by 8 ? mspadicL(mu) %16 = 2 + 3*5 + 3*5^2 + 2*5^4 + ... ? mu = mspadicmoments(Mp, PHI, -3); \\ twist by -3 < 0 ? mspadicL(mu) %18 = O(5^13) \\ always 0, phi is in the + part and D < 0
One can locate interesting symbols of level N and weight k with
? M=msinit(5,4,1); \\ M4(Gamma0(5))^+ ? L = mssplit(M, msnew(M)); \\ list of irreducible Hecke-subspaces ? phi = L[1]; \\ one Galois orbit of newforms ? #phi[1] \\... this one is rational %4 = 1 ? Mp = mspadicinit(M, 3, 4); ? mu = mspadicmoments(Mp, phi); ? mspadicL(mu) %7 = 1 + 3 + 3^3 + 3^4 + 2*3^5 + 3^6 + O(3^9) ? M = msinit(11,8, 1); \\ M8(Gamma0(11))^+ ? Mp = mspadicinit(M, 3, 4); ? L = mssplit(M, msnew(M)); ? phi = L[1]; #phi[1] \\ ... this one is two-dimensional %11 = 2 ? mu = mspadicmoments(Mp, phi); *** at top-level: mu=mspadicmoments(Mp,ph *** ^ — — — — — — -- *** mspadicmoments: incorrect type in mstooms [dimQ (eigenspace) > 1]
The library syntax is
| |
mspadicinit(M, p, n, {flag}) | |
M being a full modular symbol space, as given by
? E = ellinit("11a1"); ? [M,phi] = msfromell(E,1); ? ellap(E,3) %3 = -1 ? Mp = mspadicinit(M, 3, 10, 0); \\ commit to ordinary symbols ? PHI = mstooms(Mp,phi); If we restrict the range of allowed symbols with flag (for faster initialization), exceptions will occur if vp(ap) violates this bound:
? E = ellinit("15a1"); ? [M,phi] = msfromell(E,1); ? ellap(E,7) %3 = 0 ? Mp = mspadicinit(M,7,5,0); \\ restrict to ordinary symbols ? PHI = mstooms(Mp,phi) *** at top-level: PHI=mstooms(Mp,phi) *** ^ — — — — — *** mstooms: incorrect type in mstooms [vp(ap) > mspadicinit flag] (t_VEC). ? Mp = mspadicinit(M,7,5); \\ no restriction ? PHI = mstooms(Mp,phi); This function uses O(N2(n+k)2p) memory, where N is the level of M.
The library syntax is
| |
mspadicmoments(Mp, PHI, {D = 1}) | |
Given
The returned data μ (p-adic distributions attached to
? M = msinit(3,6, 1); ? phi = [5,-3,-1]~; ? msissymbol(M, phi) %3 = 1 ? p = 5; mshecke(M,p) * phi \\ eigenvector of T5, a5 = 6 %4 = [30, -18, -6]~ ? Mp = mspadicinit(M, p, 10, 0); \\ restrict to ordinary symbols, mod p^10 ? PHI = mstooms(Mp, phi); ? mu = mspadicmoments(Mp, PHI); ? mspadicL(mu) %8 = 5 + 2*5^2 + 2*5^3 + ... ? mu = mspadicmoments(Mp, PHI, 12); \\ twist by 12 ? mspadicL(mu) %10 = 5 + 5^2 + 5^3 + 2*5^4 + ...
The library syntax is
| |
mspadicseries(mu, {i = 0}) | |
Let Φ be the p-adic distribution-valued overconvergent symbol
attached to a modular symbol φ for Γ0(N) (eigenvector for
TN(p) for the eigenvalue ap).
If μ is the distribution on ℤp* defined by the restriction of
Φ([ oo ]-[0]) to ℤp*, let
Lp(μ,τi)(x)
= ∫ℤ_{p*} τi(t) (1+x)logp(t)/logp(u)dμ(t)
Here, τ is the Teichmüller character and u is a specific
multiplicative generator of 1+2pℤp, namely 1+p if p > 2 or 5
if p = 2. To explain
the formula, let G oo := Gal(ℚ(μp oo )/ ℚ),
let χ:G oo → ℤp* be the cyclotomic character (isomorphism)
and γ the element of G oo such that χ(γ) = u;
then
χ(γ)logp(t)/logp(u) = The p-padic precision of individual terms is maximal given the precision of the overconvergent symbol μ.
? [M,phi] = msfromell(ellinit("17a1"),1); ? Mp = mspadicinit(M, 5,7); ? mu = mspadicmoments(Mp, phi,1); \\ overconvergent symbol ? mspadicseries(mu) %4 = (4 + 3*5 + 4*5^2 + 2*5^3 + 2*5^4 + 5^5 + 4*5^6 + 3*5^7 + O(5^9)) \ + (3 + 3*5 + 5^2 + 5^3 + 2*5^4 + 5^6 + O(5^7))*x \ + (2 + 3*5 + 5^2 + 4*5^3 + 2*5^4 + O(5^5))*x^2 \ + (3 + 4*5 + 4*5^2 + O(5^3))*x^3 \ + (3 + O(5))*x^4 + O(x^5) An example with nonzero Teichmüller:
? [M,phi] = msfromell(ellinit("11a1"),1); ? Mp = mspadicinit(M, 3,10); ? mu = mspadicmoments(Mp, phi,1); ? mspadicseries(mu, 2) %4 = (2 + 3 + 3^2 + 2*3^3 + 2*3^5 + 3^6 + 3^7 + 3^10 + 3^11 + O(3^12)) \ + (1 + 3 + 2*3^2 + 3^3 + 3^5 + 2*3^6 + 2*3^8 + O(3^9))*x \ + (1 + 2*3 + 3^4 + 2*3^5 + O(3^6))*x^2 \ + (3 + O(3^2))*x^3 + O(x^4) Supersingular example (not checked)
? E = ellinit("17a1"); ellap(E,3) %1 = 0 ? [M,phi] = msfromell(E,1); ? Mp = mspadicinit(M, 3,7); ? mu = mspadicmoments(Mp, phi,1); ? mspadicseries(mu) %5 = [(2*3^-1 + 1 + 3 + 3^2 + 3^3 + 3^4 + 3^5 + 3^6 + O(3^7)) \ + (2 + 3^3 + O(3^5))*x \ + (1 + 2*3 + O(3^2))*x^2 + O(x^3),\ (3^-1 + 1 + 3 + 3^2 + 3^3 + 3^4 + 3^5 + 3^6 + O(3^7)) \ + (1 + 2*3 + 2*3^2 + 3^3 + 2*3^4 + O(3^5))*x \ + (3^-2 + 3^-1 + O(3^2))*x^2 + O(3^-2)*x^3 + O(x^4)] Example with a twist:
? E = ellinit("11a1"); ? [M,phi] = msfromell(E,1); ? Mp = mspadicinit(M, 3,10); ? mu = mspadicmoments(Mp, phi,5); \\ twist by 5 ? L = mspadicseries(mu) %5 = (2*3^2 + 2*3^4 + 3^5 + 3^6 + 2*3^7 + 2*3^10 + O(3^12)) \ + (2*3^2 + 2*3^6 + 3^7 + 3^8 + O(3^9))*x \ + (3^3 + O(3^6))*x^2 + O(3^2)*x^3 + O(x^4) ? mspadicL(mu) %6 = [2*3^2 + 2*3^4 + 3^5 + 3^6 + 2*3^7 + 2*3^10 + O(3^12)]~ ? ellpadicL(E,3,10,,5) %7 = 2 + 2*3^2 + 3^3 + 2*3^4 + 2*3^5 + 3^6 + 2*3^7 + O(3^10) ? mspadicseries(mu,1) \\ must be 0 %8 = O(3^12) + O(3^9)*x + O(3^6)*x^2 + O(3^2)*x^3 + O(x^4)
The library syntax is
| |
mspathgens(M) | |
Let Δ0 := Div0(ℙ1(ℚ)).
Let M being a full modular symbol space, as given by
An element [v]-[u] in Δ0 is coded by the "path" [u,v],
where
? M = msinit(11,8); \\ M8(Gamma0(11)) ? [g,R] = mspathgens(M); ? g %3 = [[+oo, 0], [0, 1/3], [1/3, 1/2]] \\ 3 paths ? #R \\ a single relation %4 = 1 ? r = R[1]; #r \\ ...involving all 3 generators %5 = 3 ? r[1] %6 = [[1, 1; [1, 1; 0, 1], -1], 1] ? r[2] %7 = [[1, 1; [7, -2; 11, -3], -1], 2] ? r[3] %8 = [[1, 1; [8, -3; 11, -4], -1], 3] The given relation is of the form ∑i (1-γi) gi = 0, with γi ∈ Γ0(11). There will always be a single relation involving all generators (corresponding to a round trip along all cusps), then relations involving a single generator (corresponding to 2 and 3-torsion elements in the group:
? M = msinit(2,8); \\ M8(Gamma0(2)) ? [g,R] = mspathgens(M); ? g %3 = [[+oo, 0], [0, 1]] Note that the output depends only on the group G, not on the representation V.
The library syntax is
| |
mspathlog(M, p) | |
Let Δ0 := Div0(ℙ1(ℚ)).
Let M being a full modular symbol space, as given by Returns (pi) in ℤ[G] such that p = ∑i pi gi.
? M = msinit(2,8); \\ M8(Gamma0(2)) ? [g,R] = mspathgens(M); ? g %3 = [[+oo, 0], [0, 1]] ? p = mspathlog(M, [1/2,2/3]); ? p[1] %5 = [[1, 0; 2, 1] 1] ? p[2] %6 = [[1, 0; 0, 1] 1] [[3, -1; 4, -1] 1] ? mspathlog(M, [1,2;2,3]) == p \\ give path via a 2x2 matrix %7 = 1 Note that the output depends only on the group G, not on the representation V.
The library syntax is
| |
mspetersson(M, {F}, {G = F}) | |
M being a full modular symbol space for Γ = Γ0(N),
as given by Suppose that f1 and f2 are two parabolic forms. Let F1 and F2 be the attached modular symbols Fi(δ) = ∫δ fi(z).(z X + Y)k-2 dz and let Fℝ1, Fℝ2 be the attached real modular symbols Fℝi(δ) = ∫δ Re(fi(z).(z X + Y)k-2 dz) Then we have { Fℝ1, Fℝ2 } = -2 (2i)k-2.
Im( < f1,f2 > Petersson)
and
{ F1, F2 } = (2i)k-2 < f1,f2 > Petersson
In weight 2, the intersection product {F, G} has integer values on the
ℤ-structure on M given by For user convenience, we allow F and G to be matrices and return the attached Gram matrix. If F is omitted: treat it as the full modular space attached to M; if G is omitted, take it equal to F.
? M = msinit(37,2); ? C = mscuspidal(M)[1]; ? mspetersson(M, C) %3 = [ 0 -17 -8 -17] [17 0 -8 -25] [ 8 8 0 -17] [17 25 17 0] ? mspetersson(M, mslattice(M,C)) %4 = [0 -1 0 -1] [1 0 0 -1] [0 0 0 -1] [1 1 1 0] ? E = ellinit("33a1"); ? [M,xpm] = msfromell(E); [xp,xm,L] = xpm; ? mspetersson(M, mslattice(M,L)) %7 = [0 -3] [3 0] ? ellmoddegree(E) %8 = [3, -126] The coefficient 3 in the matrix is the degree of the modular parametrization.
The library syntax is
| |
mspolygon(M, {flag = 0}) | |
M describes a subgroup G of finite index in the modular group
PSL2(ℤ), as given by * Its vertices are an ordered list in ℙ1(ℚ) and contain a representatives of all cusps. * Its edges are hyperbolic arcs joining two consecutive vertices; each edge e is labelled by an integer μ(e) ∈ { oo ,2,3}. * Given a path (a,b) between two elements of ℙ1(ℚ), let (a,b) = (b,a) be the opposite path. There is an involution e → e* on the edges. We have μ(e) = oo if and only if e ! = e*; when μ(e) ! = 3, e is G-equivalent to e*, i.e. there exists γe ∈ G such that e = γe e*; if μ(e) = 3 there exists γe ∈ G of order 3 such that the hyperbolic triangle (e, γe e, γe2 e) is invariant by γe. In all cases, to each edge we have attached γe ∈ G of order μ(e). The polygon is given by a triple [E, A, g] * The list E of its consecutive edges as matrices in M2(ℤ).
* The permutation A attached to the involution: if e = E[i] is the
i-th edge, then * The list g of pairing matrices γe. Remark that γe* = γe-1 if μ(e) ! = 3, i.e., g[i]-1 = g[A[i]] whenever i ! = A[i] (μ(g[i]) = 1) or μ(g[i]) = 2 (g[i]2 = 1). Modulo these trivial relations, the pairing matrices form a system of independant generators of G. Note that γe is elliptic if and only if e* = e. The above data yields a fundamental domain for G acting on Poincaré's half-plane: take the convex hull of the polygon defined by * The edges in E such that e ! = e* or e* = e, where the pairing matrix γe has order 2; * The edges (r,t) and (t,s) where the edge e = (r,s) ∈ E is such that e = e* and γe has order 3 and the triangle (r,t,s) is the image of (0,exp(2iπ/3), oo ) by some element of PSL2(ℚ) formed around the edge. Binary digits of flag mean: 1: return a normalized hyperbolic polygon if set, else a polygon with unimodular edges (matrices of determinant 1). A polygon is normalized in the sense of compact orientable surfaces if the distance d(a,a*) between an edge a and its image by the involution a* is less than 2, with equality if and only if a is linked with another edge b (a, b, a* et b* appear consecutively in E up to cyclic permutation). In particular, the vertices of all edges such that that d(a,a*) ! = 1 (distance is 0 or 2) are all equivalent to 0 modulo G. The external vertices of a a* such that d(a,a*) = 1 are also equivalent to 0; the internal vertices a∩ a* (a single point), together with 0, form a system of representatives of the cusps of G\ℙ1(ℚ). This is useful to compute the homology group H1(G,ℤ) as it gives a symplectic basis for the intersection pairing. In this case, the number of parabolic matrices (trace 2) in the system of generators G is 2(t-1), where t is the number of non equivalent cusps for G. This is currently only implemented for G = Γ0(N).
2: add graphical representations (in LaTeX form) for the hyperbolic polygon
in Poincaré's half-space and the involution a → a* of the Farey symbol.
The corresponding character strings can be included in a LaTeX document
provided the preamble contains
? [V,A,g] = mspolygon(3); ? V %2 = [[-1, 1; -1, 0], [1, 0; 0, 1], [0, 1; -1, 1]] ? A %3 = Vecsmall([2, 1, 3]) ? g %4 = [[-1, -1; 0, -1], [1, -1; 0, 1], [1, -1; 3, -2]] ? [V,A,g, D1,D2] = mspolygon(11,2); \\ D1 and D2 contains pictures ? {write("F.tex", "\\documentclass{article}\\usepackage{tikz}\\begin{document}" D1, "\n", D2, "\\end{document}");} ? [V1,A1] = mspolygon(6,1); \\ normalized ? V1 %8 = [[-1, 1; -1, 0], [1, 0; 0, 1], [0, 1; -1, 3], [1, -2; 3, -5], [-2, 1; -5, 2], [1, -1; 2, -1]] ? A1 %9 = Vecsmall([2, 1, 4, 3, 6, 5]) ? [V0,A0] = mspolygon(6); \\ not normalized V[3]* = V[6], d(V[3],V[6]) = 3 ? A0 %11 = Vecsmall([2, 1, 6, 5, 4, 3]) ? [V,A] = mspolygon(14, 1); ? A %13 = Vecsmall([2, 1, 4, 3, 6, 5, 9, 10, 7, 8]) One can see from this last example that the (normalized) polygon has the form (a1, a1*, a2, a2*, a3, a3*, a4, a5, a4*, a5*), that X0(14) is of genus 1 (in general the genus is the number of blocks of the form aba*b*), has no elliptic points (A has no fixed point) and 4 cusps (number of blocks of the form aa* plus 1). The vertices of edges a4 and a5 all project to 0 in X0(14): the paths a4 and a5 project as loops in X0(14) and give a symplectic basis of the homology H1(X0(14),ℤ).
? [V,A] = mspolygon(15); ? apply(matdet, V) \\ all unimodular %2 = [1, 1, 1, 1, 1, 1, 1, 1, 1, 1] ? [V,A] = mspolygon(15,1); ? apply(matdet, V) \\ normalized polygon but no longer unimodular edges %4 = [1, 1, 1, 1, 2, 2, 47, 11, 47, 11]
The library syntax is
| |
msqexpansion(M, projH, {B = seriesprecision}) | |
M being a full modular symbol space, as given by This function uses a naive O(B2 d3) algorithm, where d = O(kN) is the dimension of Mk(Γ0(N)).
? M = msinit(11,2, 1); \\ M2(Gamma0(11))^+ ? L = mssplit(M, msnew(M)); ? msqexpansion(M,L[1], 20) %3 = [1, -2, -1, 2, 1, 2, -2, 0, -2, -2, 1, -2, 4, 4, -1, -4, -2, 4, 0, 2] ? ellan(ellinit("11a1"), 20) %4 = [1, -2, -1, 2, 1, 2, -2, 0, -2, -2, 1, -2, 4, 4, -1, -4, -2, 4, 0, 2]
The shortcut
? E = ellinit("11a1"); ? [M,S] = msfromell(E); [sp,sm] = S; ? msqexpansion(M,sp,10) \\ in the + eigenspace %3 = [1, -2, -1, 2, 1, 2, -2, 0, -2, -2] ? msqexpansion(M,sm,10) \\ in the - eigenspace %4 = [1, -2, -1, 2, 1, 2, -2, 0, -2, -2] ? ellan(E, 10) %5 = [1, -2, -1, 2, 1, 2, -2, 0, -2, -2]
The library syntax is
| |
mssplit(M, {H}, {dimlim}) | |
Let M denote a full modular symbol space, as given by
? M = msinit(11,8, 1); \\ M8(Gamma0(11))^+ ? L = mssplit(M); \\ split msnew(M) ? #L %3 = 2 ? f = msqexpansion(M,L[1],5); f[1].mod %4 = x^2 + 8*x - 44 ? lift(f) %5 = [1, x, -6*x - 27, -8*x - 84, 20*x - 155] ? g = msqexpansion(M,L[2],5); g[1].mod %6 = x^4 - 558*x^2 + 140*x + 51744 To a Hecke-simple subspace corresponds an orbit of (normalized) newforms, defined over a number field. In the above example, we printed the polynomials defining the said fields, as well as the first 5 Fourier coefficients (at the infinite cusp) of one such form.
The library syntax is
| |
msstar(M, {H}) | |
M being a full modular symbol space, as given by
? M = msinit(11,2); \\ M2(Gamma0(11)) ? w = msstar(M); ? w^2 == 1 %3 = 1
The library syntax is
| |
mstooms(Mp, phi) | |
Given Under the action of TN_{0}(p), φ generates a subspace Wφ of dimension 1 (if p | N) or 2 (if p does not divide N) in the space of modular symbols of level N0.
Let Vp = [p,0;0,1] and Cp = [ap,pk-1;-1,0].
When p does not divide N and ap is divisible by p,
When p does not divide N and ap is not divisible by p,
The resulting overconvergent eigensymbol can then be used in
? M = msinit(3,6, 1); p = 5; ? Tp = mshecke(M, p); factor(charpoly(Tp)) %2 = [x - 3126 2] [ x - 6 1] ? phi = matker(Tp - 6)[,1] \\ generator of p-Eigenspace, ap = 6 %3 = [5, -3, -1]~ ? Mp = mspadicinit(M, p, 10, 0); \\ restrict to ordinary symbols, mod p^10 ? PHI = mstooms(Mp, phi); ? mu = mspadicmoments(Mp, PHI); ? mspadicL(mu) %7 = 5 + 2*5^2 + 2*5^3 + ... A non ordinary symbol.
? M = msinit(4,6,1); p = 3; ? Tp = mshecke(M, p); factor(charpoly(Tp)) %2 = [x - 244 3] [ x + 12 1] ? phi = matker(Tp + 12)[,1] \\ ap = -12 is divisible by p = 3 %3 = [-1/32, -1/4, -1/32, 1]~ ? msissymbol(M,phi) %4 = 1 ? Mp = mspadicinit(M,3,5,0); ? PHI = mstooms(Mp,phi); *** at top-level: PHI=mstooms(Mp,phi) *** ^ — — — — — *** mstooms: incorrect type in mstooms [vp(ap) > mspadicinit flag] (t_VEC). ? Mp = mspadicinit(M,3,5,1); ? PHI = mstooms(Mp,phi);
The library syntax is
| |