##################################################
## Dirichlet Characters
##################################################

> G<a,b,c> := DirichletGroup(8*13, CyclotomicField(12));
> G;
Group of Dirichlet characters of modulus 104 over Cyclotomic 
Field of order 12 and degree 4
> // WARNING: The default ring is Q, not Q(zeta_n).
> DirichletGroup(8*13);   
Group of Dirichlet characters of modulus 104 over Rational Field

> G<a,b,c> := DirichletGroup(8*13, CyclotomicField(12));
> [Order(a), Order(b), Order(c)];
[ 2, 2, 12 ]
> [Conductor(a), Conductor(b), Conductor(c)];
[ 4, 8, 13 ]
> a(3); 
-1
> b(3);
-1
> c(3);
zeta_12^2 - 1

> G<a> := DirichletGroup(5,CyclotomicField(4));
> H<b> := DirichletGroup(7,RationalField());
> Parent(a*b);
Group of Dirichlet characters of modulus 35 over Cyclotomic Field
of order 4 and degree 2

> c := a*b;
> d := Extend(c,70);  // natural extension to character of modulus 70
> Parent(d);
Group of Dirichlet characters of modulus 70 over Cyclotomic Field
of order 4 and degree 2
> Conductor(d);
35
> Modulus(AssociatedPrimitiveCharacter(d));
35

> G := DirichletGroup(35,CyclotomicField(28));  G;
Group of Dirichlet characters of modulus 35 over Cyclotomic Field
of order 28 and degree 12
> e := G!c;
> Parent(e);
Group of Dirichlet characters of modulus 35 over Cyclotomic Field
of order 28 and degree 12
> c(3);
zeta_4
> e(3);
zeta_28^7
> Parent(MinimalBaseRingCharacter(e));
Group of Dirichlet characters of modulus 35 over Cyclotomic Field
of order 4 and degree 2


##################################################
## What Are Modular Symbols?
##################################################

> G<eps> := DirichletGroup(13, CyclotomicField(4));
> eps(-1);
-1
> Order(eps);
4
> M := ModularSymbols(eps, 5);   // eps determines N  !!
> M;
Full modular symbols space of level 13, weight 5, character eps, 
and dimension 8 over Cyclotomic Field of order 4 and degree 2

> [ManinSymbol(x)[1] : x in Basis(M)];
[   <X^3, (0 1)>, <X^3, (1 11)>, <X^3, (1 5)>, <X^3, (1 3)>,
    <X^3, (1 4)>, <X^3, (1 6)>, <X^3, (1 12)>, <X^3, (1 0)>]

> M.1;
X^3*{0, oo}
> M.2;
(1331*X^3 + 363*X^2*Y + 33*X*Y^2 + Y^3)*{-1/11, 0}


##################################################
## Efficiency of Computation of Presentation
##################################################

> time M := ModularSymbols(2004);  M; // on my laptop
Time: 1.440
Full modular symbols space for Gamma_0(2004) of weight 2 and 
dimension 673 over Rational Field
> GetMemoryUsage();
9306624    // about 9.3MB
> time M := ModularSymbols(10000);  
Current total memory usage: 419.4MB, failed memory request: 206.0MB
System error: Out of memory.
> time M := ModularSymbols(10007);  M;   // prime level, so easier
Time: 6.590       
Full modular symbols space for Gamma_0(10007) of weight 2 and 
dimension 1669 over Rational Field
> time t2 := HeckeOperator(M,2);
Time: 3.270

> G<e> := DirichletGroup(389,CyclotomicField(97));
> Order(e);
194
> time M := ModularSymbols(e^2,2);
Time: 0.320
> M;
Full modular symbols space of level 389, weight 2, character e^2,
and dimension 64 over Cyclotomic Field of order 97 and degree 96
> G<e> := DirichletGroup(37,CyclotomicField(36));
> time M := ModularSymbols(e,7);   // weight 7
Time: 3.260
> M;
Full modular symbols space of level 37, weight 7, character e, 
and dimension 38 over Cyclotomic Field of order 36 and degree 12


##################################################
## Efficiency Trick: Work Mod~$p$
##################################################

> G<a,b,c> := DirichletGroup(2000,GF(5));
> Conductor(c);
5
> Order(c);
4
> time M := ModularSymbols(c,3);
Time: 5.070
> M;
Full modular symbols space of level 2000, weight 3, character c, 
and dimension 1200 over Finite field of size 5

> function f(N,k,p) 
     return Dimension(ModularSymbols(N,k,GF(p))) - 
            Dimension(ModularSymbols(N,k));
  end function;
> [N : N in [2..100] | f(N,2,2) gt 0];
[ 5, 10, 13, 17, 25, 26, 29, 34, 37, 41, 50, 53, 58, 61, 65, 73, 
  74, 82, 85, 89, 97 ]


##################################################
## Hecke Operators on Modular Symbols
##################################################

> HeilbronnMerel(2);
[
[ 1, 0, 0, 2 ],
[ 1, 0, 1, 2 ],
[ 2, 0, 0, 1 ],
[ 2, 1, 0, 1 ]
]
> #HeilbronnMerel(29);
199
> #HeilbronnMerel(10007);   // takes a while
337977
> #HeilbronnCremona(10007);   // in some cases these can be used...
67698

> G<eps> := DirichletGroup(13, CyclotomicField(4)); 
> M := ModularSymbols(eps, 5);   
> T2 := HeckeOperator(M,2);
> Nrows(T2);
8
> T2[1];
(zeta_4 + 16 -3/4 1/4 3/4 -3/4 0 2 -3/2)
> F := CharacteristicPolynomial(T2);
> R<X> := Parent(F);
> Factorization(F);
[   <X - 16*zeta_4 - 1, 1>,
    <X - zeta_4 - 16, 1>,
    <X^3 + (zeta_4 + 1)*X^2 - 23*zeta_4*X - 29*zeta_4 + 29, 2>]


##################################################
## Subspaces of Modular Symbols
##################################################

> G<eps> := DirichletGroup(13, CyclotomicField(4)); 

> M := ModularSymbols(eps, 5);   

> S := CuspidalSubspace(M); S;
Modular symbols space of level 13, weight 5, character eps, and 
dimension 6 over Cyclotomic Field of order 4 and degree 2

> Factorization(CharacteristicPolynomial(HeckeOperator(S,2)));
[  <X^3 + (zeta_4 + 1)*X^2 - 23*zeta_4*X - 29*zeta_4 + 29, 2>  ]

> M := ModularSymbols(33); M;
Full modular symbols space for Gamma_0(33) of weight 2 and 
dimension 9 over Rational Field

> NewSubspace(M);
Modular symbols space for Gamma_0(33) of weight 2 and dimension 3
over Rational Field

> Complement(NewSubspace(M));
Modular symbols space for Gamma_0(33) of weight 2 and dimension 6
over Rational Field


##################################################
## Decomposition of Modular Symbols Spaces
##################################################

> S := CuspidalSubspace(ModularSymbols(33,2));

> NewformDecomposition(S);
[
    Modular symbols space for Gamma_0(33) of weight 2 and 
    dimension 2 over Rational Field,                      // new
    Modular symbols space for Gamma_0(33) of weight 2 and 
    dimension 4 over Rational Field                       // old
]

> time M := ModularSymbols(700,2,+1);
Time: 0.270
> time S := CuspidalSubspace(M);
Time: 0.190
> time D := NewformDecomposition(S);
Time: 10.990
> #D;
34
> D;
[
    Modular symbols space for Gamma_0(700) of weight 2 and 
    dimension 1 over Rational Field,
    Modular symbols space for Gamma_0(700) of weight 2 and 
    dimension 1 over Rational Field,
    ...
    Modular symbols space for Gamma_0(700) of weight 2 and 
    dimension 6 over Rational Field
]

% > HeckeBound(M);
% 8

$ ls magma/package/Geometry/ModSym/*.m
analytic.m       cusps.m       inner_twists.m          operators.m
arith.m          decomp.m      intersection_pairing.m  period.m
boundary.m       derivative.m  linalg.m                qexpansion.m
calc.m           dims.m        maps.m                  representation.m
charpolyhecke.m  dirichlet.m   mestre.m                subspace.m
compgrp.m        eisenstein.m  modsym.m                tests.m
core.m           elliptic.m    multichar.m             verbose.m


##################################################
## Computing Modular Forms Using MAGMA
##################################################

> G<eps> := DirichletGroup(13,CyclotomicField(6));
> M := ModularSymbols(eps,2, +1);
> S := CuspidalSubspace(M);
> S;
Modular symbols space of level 13, weight 2, character eps, and 
dimension 1 over Cyclotomic Field of order 6 and degree 2
> qExpansionBasis(S,4);
[   q + (-zeta_6 - 1)*q^2 + (2*zeta_6 - 2)*q^3 + O(q^4) ]

> M := ModularSymbols(33,2);
> S := CuspidalSubspace(M);

> qExpansionBasis(S,10);
[   q - q^5 - 2*q^6 + 2*q^7 - 2*q^8 - q^9 + O(q^10),
    q^2 - q^4 - q^5 - q^6 + 2*q^7 - q^8 + q^9 + O(q^10),
    q^3 - 2*q^6 - q^9 + O(q^10)   ]

> qExpansionBasis(OldSubspace(S),10);
[   q - 2*q^2 + 2*q^4 + q^5 - 2*q^7 - 3*q^9 + O(q^10),
    q^3 - 2*q^6 - q^9 + O(q^10) ]

> qExpansionBasis(NewSubspace(S),10);
[   q + q^2 - q^3 - q^4 - 2*q^5 - q^6 + 4*q^7 - 3*q^8 + q^9 + 
        O(q^10) ]

> M := ModularSymbols(389,2, 1);
> S := CuspidalSubspace(M);
> D := Decomposition(S,2);
> V := D[3]; V;
Modular symbols space for Gamma_0(389) of weight 2 and dimension 
3 over Rational Field
> qExpansionBasis(V,10);
[   q - q^5 - 2*q^6 - q^7 + 2*q^8 - q^9 + O(q^10),
    q^2 - q^3 + O(q^10),
    q^4 - q^5 - q^6 + q^9 + O(q^10) ]
> qEigenform(V,6);    // eigenform in span of above q-expansions
q + a*q^2 - a*q^3 + (a^2 - 2)*q^4 + (-a^2 + 1)*q^5 + O(q^6)
> BaseRing(Modulus(Parent($1)));
Univariate Quotient Polynomial Algebra in a over Rational Field
with modulus a^3 - 4*a - 2

> M := ModularForms(33,2); M;
Space of modular forms on Gamma_0(33) of weight 2 and 
dimension 6 over Integer Ring.
> Basis(M);
[   1 + O(q^8),
    q - q^5 + 2*q^7 + O(q^8),
    q^2 + 2*q^7 + O(q^8),
    q^3 + O(q^8),
    q^4 + q^5 + O(q^8),
    q^6 + O(q^8)   ]
> SetPrecision(M,15);
> Basis(M);
[   1 + 12*q^11 + O(q^15),
    q - q^5 + 2*q^7 - 2*q^8 + q^9 - 2*q^10 - q^11 + 4*q^12 + 4*q^14 + O(q^15),
    q^2 + 2*q^7 + q^8 + q^9 + 2*q^10 + q^11 + q^12 + 2*q^14 + O(q^15),
    q^3 + q^9 - 2*q^11 + 4*q^12 + O(q^15),
    ... ]

> M := ModularForms(1,12);
> Newforms(M);
[* [*
q - 24*q^2 + 252*q^3 - 1472*q^4 + 4830*q^5 - 6048*q^6 - 16744*q^7
+ O(q^8)
*], [*
691/65520 + q + 2049*q^2 + 177148*q^3 + 4196353*q^4 + 
48828126*q^5 + 362976252*q^6 + 1977326744*q^7 + O(q^8)
*] *]

> M := ModularForms(23,2);
> S := CuspidalSubspace(M);
> S;
Space of modular forms on Gamma_0(23) of weight 2 and dimension 2
over Integer Ring.
> Newforms(S);
[* [*
q + a*q^2 + (-2*a - 1)*q^3 + (-a - 1)*q^4 + 2*a*q^5 + (a - 2)*q^6
  + (2*a + 2)*q^7 + O(q^8),
q + b*q^2 + (-2*b - 1)*q^3 + (-b - 1)*q^4 + 2*b*q^5 + (b - 2)*q^6
  + (2*b + 2)*q^7 + O(q^8)
*] *]
> Parent($1[1][1]);
Space of modular forms on Gamma_0(23) of weight 2 and dimension 2
over Number Field with defining polynomial x^2 + x - 1 over the 
Rational Field.

> f := Newforms(S)[1][1];  g := Newforms(S)[1][2];
> f + g;
>> f + g;
     ^
Runtime error in '+': Arguments 1 and 2 have incompatible 
coefficient rings.

> ComplexEmbeddings(f);
[* [* q - 1.618033988749894848204586834365638117720*q^2 +  ...
      q + 0.618033988749894848204586834365638117720*q^2 -  ... *] *]
> $1[1][1] + $1[1][2];
2*q - q^2 - q^4 - 2.0...*q^5 - 
> pAdicEmbeddings(f,2);
[* [* O(2^20) + (1 + O(2^20))*q + ((1 + O(2^20))*a + O(2^20))*q^2 + ...
      O(2^20) + (1 + O(2^20))*q + ((1 + O(2^20))*b + O(2^20))*q^2 + ... *] *]
> pAdicEmbeddings(f,11);
[* [*  O(11^20) + (1 + O(11^20))*q + (273946294811098331671 + ...
*], [* O(11^20) + (1 + O(11^20))*q - (273946294811098331672 + ... *] *]
> $1[1][1] + $1[2][1];
O(11^20) + (2 + O(11^20))*q - (1 + O(11^20))*q^2 + ...

> Reductions(f,2);
[* [*  q + $.1*q^2 + q^3 + $.1^2*q^4 + $.1*q^6 + O(q^8),
       q + $.1^2*q^2 + q^3 + $.1*q^4 + $.1^2*q^6 + O(q^8) *] *]
> Reductions(f,11);
[* [* q + 7*q^2 + 7*q^3 + 3*q^4 + 3*q^5 + 5*q^6 + 5*q^7 + O(q^8) *], 
   [* q + 3*q^2 + 4*q^3 + 7*q^4 + 6*q^5 + q^6 + 8*q^7 + O(q^8) *] *]
> f11 := Reductions(f,11)[1][1];
> Type(f11);
ModFrmElt
> f11;
q + 7*q^2 + 7*q^3 + 3*q^4 + 3*q^5 + 5*q^6 + 5*q^7 + O(q^8)
> PowerSeries(f11,15);
q + 7*q^2 + 7*q^3 + 3*q^4 + 3*q^5 + 5*q^6 + 5*q^7 + 7*q^8 + 2*q^9
    + 10*q^10 + 4*q^11 + 10*q^12 + 3*q^13 + 2*q^14 + O(q^15)

> G<eps> := DirichletGroup(13, CyclotomicField(6));
> M := ModularForms(eps);
> BaseRing(M);
Integer Ring
> S := CuspidalSubspace(M);
> S;
Space of modular forms on Gamma_1(13) with character all 
conjugates of [eps], weight 2, and dimension 2 over Integer Ring.
> Basis(S);
[
    q - 4*q^3 - q^4 + 3*q^5 + 6*q^6 + O(q^8),
    q^2 - 2*q^3 - q^4 + 2*q^5 + 2*q^6 + O(q^8)
]

%> Newforms(S);
% [* [*
% q + (-a - 1)*q^2 + (2*a - 2)*q^3 + a*q^4 + (-2*a + 1)*q^5 + (-2*a
% + 4)*q^6 + O(q^8),
% q + (-b - 1)*q^2 + (2*b - 2)*q^3 + b*q^4 + (-2*b + 1)*q^5 + (-2*b
% + 4)*q^6 + O(q^8)
% *] *]

$ ls magma/package/Geometry/ModFrm/*.m
abelian_varieties.m  eisenstein.m      modular_symbols.m  relations.m
arithmetic.m         elliptic_curve.m  newforms.m         subspaces.m
bases.m              hecke_algebras.m  operators.m        tests.m
categories.m         input_output.m    p-adic.m           verbose.m
congruences.m        l_series.m        predicates.m       weight1table.m
creation.m           level1.m          q-expansions.m
degeneracy_maps.m    misc.m            qexp_mappings.m


##################################################
## Modular Abelian Varieties and Modular Symbols
##################################################

> M := ModularSymbols(389);    
> S := CuspidalSubspace(M);
> D := NewformDecomposition(S);
> [Dimension(A)/2 : A in D];          // dimensions of abvars A_f
[ 1, 2, 3, 6, 20 ]
> [ModularDegree(D[i]) : i in [1..#D]];   
[ 40, 144, 992, 17856, 20480 ]
> [LRatio(D[i],1) : i in [1..#D]];    // BSD Ratios L(A_f,1)/Omega
[ 0, 0, 0, 0, 51200/97 ]
> Factorization(51200);
[ <2, 11>, <5, 2> ]
> LSeriesLeadingCoefficient(D[1],1,100);
0.75931650029224679065762600319 2
> E := EllipticCurve(A);  AnalyticRank(E);   // Watkin's new code
2 0.7593000000
> LSeriesLeadingCoefficient(D[2],1,100);
1.487184621319346836916654326667 1

> TamagawaNumber(D[1],389);          // c_{389} = 1 for elliptic curve
1
> TamagawaNumber(D[5],389);          // c_{389} = 97 for 20-dim quotient
97
> TorsionBound(D[5],13);             // multiple of order of torsion 
97 
> #RationalCuspidalSubgroup(D[5]);   // divisor of order of torsion
97
> Invariants(IntersectionGroup(D[1],D[2]));
[ 2, 2 ]
> Invariants(IntersectionGroup(D[1],D[5]));
[ 20, 20 ]


##################################################
## The Modular Abelian Varieties Package
##################################################

> J := JZero(389); J;
Modular abelian variety JZero(389) of dimension 32 and level 389 
over Q
> D := Decomposition(J);
> [Dimension(A) : A in D];
[ 1, 2, 3, 6, 20 ]
> [ModularDegree(A) : A in D];   
[ 40, 144, 992, 17856, 20480 ]
> [LRatio(A,1) : A in D];
[ 0, 0, 0, 0, 51200/97 ]
> L := LSeries(D[1]); L;
L(389A,s): L-series of Modular abelian variety 389A of dimension 
1, level 389 and conductor 389 over Q
> LeadingCoefficient(L,1,200);    // slow, since doesn't use Watkins (but *general*)
0.75931650029224679065762600319 2
> TamagawaNumber(D[1],389);
1 1 true

> TamagawaNumber(D[5],389);  
97 97 true   
> TorsionLowerBound(D[5]);
97
> TorsionMultiple(D[5]);
97
> G := RationalCuspidalSubgroup(D[5]); G;
Finitely generated subgroup ... with invariants [ 97 ]
> B := D[5]/G; B;     // quotients by anything are defined.
Modular abelian variety of dimension 20 and level 389 over Q
> H := D[1] meet D[5];  H; // takes a while
Finitely generated subgroup ... with invariants [ 20, 20 ]

> J := JZero(22);
> [Matrix(phi) : phi in Basis(End(J))];
...

> function f(N) 
      J := JZero(N); 
      T := HeckeAlgebra(J); 
      return Index(Saturation(T),T); 
  end function;
> for N in [1..120] do print N, f(N); end for;
...

$ ls magma/package/Geometry/ModAbVar/*.m
arithabvar.m   elt.m       homspace.m      misc.m       periods.m
compgrp.m      endo_alg.m  inner_twists.m  modabvar.m   rings.m
complements.m  fields.m    linalg.m        morphisms.m  subgrp.m
decomp.m       heegner.m   lser.m          new_old.m    test.m
ellcrv.m       homology.m  map.m           operators.m  torsion.m



##########################################################
# The following are all the examples not in functions.   #
##########################################################

def examples():
    """
    
    """


if __name__ ==  '__main__':
    import doctest, sys
    doctest.testmod(sys.modules[__name__])