/* MAGMA code used for computations in Section 6 of
"Automorphic forms and rational homology 3-spheres"
by Frank Calegari and Nathan M. Dunfield
see that section for some of the basic setup.
If you just want to _use_ the code as is, skip ahead to
"SuperTestFile"
*/
/* Basic setup for twist knot orbifolds.
n = number of twists, k = orbifold amount
n = 0: unknot, n = 1: trefoil, n = -1: figure eight
Here is the fundamental group:
*/
TwistOrbifoldGroup := function(n, k)
F := Group;
w := b * (b * a)^-1 * a;
return F/sub;
end function;
// Helper function for moving a multivariable polynomial from
// one ring to another.
FlattenLiftPolyonomial := function(P, F, z)
R := PolynomialRing(F);
if BaseRing(P) eq Rationals() then
return elt;
else
return elt, z): x in Coefficients(P)]>;
end if;
end function;
// We need to tell MAGMA how to make a
// permutation representation after it has constructed
// the matrices in PSL(2, F).
// The points in P^1(F)
PointsInP1 := function(F)
x := PrimitiveElement(F);
return [ [1, 0] ] cat [ [0,1] ] cat [ [x^k, 1] : k in [0..(#F - 2)]];
end function;
IndexOfPointInP1 := function(x)
if x[2] eq 0 then
return 1;
end if;
if x[1] eq 0 then
return 2;
end if;
return Log(x[1]/x[2]) + 3;
end function;
// Action of a matrix on P1
ProjectiveAction := function(M, F)
V := VectorSpace(F, 2);
return [ IndexOfPointInP1( (V ! x) * M) : x in PointsInP1(F) ];
end function;
// Computes the first betti number mod a suitable prime
HomologyOfProjCover := function(f)
return #AbelianQuotientInvariants( sub, 31991);
end function;
// alternate version, with different prime.
HomologyOfProjCover2 := function(f)
return #AbelianQuotientInvariants( sub, 32003);
end function;
// Code used to format the output when sending to a file
CarRetToSpace := function(x)
if x eq "\n" then return " "; else return x; end if;
end function;
RemoveCarReturns := function(str)
return &cat [ CarRetToSpace(s) : s in Eltseq(str) ];
end function;
FormatLine := function(p, ranks, t)
return RemoveCarReturns(Sprint([* p, ranks, t*])) cat "," ;
end function;
/*
Now we move on to the equations defining the character variety,
i.e. the which traces give rise to representations. Below is the
recursive function; the one you want to call follows.
*/
HSPolyRec := function(n, R)
x := R.1;
z := R.2;
if n eq 0 then return One(R); end if;
if n eq 1 then return z - 1; end if;
if n gt 0 then
return (2 + 2*x - 2*z -x*z + z^2)*$$(n-1, R) - $$(n-2, R);
else
return (2 + 2*x - 2*z -x*z + z^2)*$$(n+1, R) - $$(n+2, R);
end if;
end function;
/*
Polynomial from Hoste-Shanahan. The original variables are
x = tr(a^2), z = tr(ab).
I want to use a different normalization than HS.
The two final two variables
are "a" which is the trace of "a" (not the trace of a^2).
ans "c" which is the trace of "ab" (denoted "z" in HS).
My conventions are implements by the substitution
at the last stage.
*/
HSPoly := function(n)
R:= PolynomialRing(Integers(),2);
return R, Evaluate(HSPolyRec(n, R), a, a^2 - 2);
end function;
// Return tr(a) for the matrix with eigenvalues
// e^(pi I /k) , e^(-pi I /k)
FieldWithTraceOfRootOfOne := function(k)
if (k mod 2) eq 0 then m := 2*k; else m := k; end if;
L := CyclotomicField(m);
p := MinimalPolynomial(z + 1/z);
K := NumberField(p);
return K.1, K;
end function;
// Now, we will use the HS polynomial to write
// down the def. polynomial for the character variety.
TraceVarietyEquation := function(n,k)
x, K := FieldWithTraceOfRootOfOne(k);
R := PolynomialRing(K);
S := PolynomialRing(K,2);
J, r := HSPoly(n);
return x, K, R, R ! UnivariatePolynomial(Evaluate(S ! r, X, x));
end function;
/*
Unfortunately, the character variety is not always irreducible over Q.
This is important because we want to look at just the congruence
quotients associated to the holonomy representation, not from some
other random rep.
For the range of n,k constrained to here, the non-geometric factor is
easy to figure out because it corresponds to representations with all
real traces --- thus reps into PSL(2,R) or SU(2). For these examples,
it is is always a linear factor over Q(tr(a)) and so it's easy to remove.
*/
TraceHolonomyEquation := function(n, k)
assert (Abs(n) le 7) and (k le 10);
z, K, R, h := TraceVarietyEquation(n, k);
factors := [x[1] : x in Factorization(h)];
d, i := Maximum([Degree(x) : x in factors]);
return K, R, factors[i];
end function;
/*
Now on to finding the actual reps. The first function is used to mod
out by Galois automorphisms of the target PSL(2, F) when F is not a
prime field.
*/
// assumes the elements of x are Sbar/S Galois invariant
ModOutByFrob := function(S, x)
if #x eq 0 then return []; end if;
d := Degree(S);
p := Characteristic(S);
G := SymmetricGroup(#x);
H := sub;
return [ x[y[2]] : y in OrbitRepresentatives(H)];
end function;
// The possible traces for rho(a) in F
PossATraces := function(k, F)
if (k mod 2) eq 0 then m := 2*k; else m := k; end if;
L := CyclotomicField(m);
R := PolynomialRing(F);
f := elt< R | MinimalPolynomial(z + 1/z)>;
roots := ModOutByFrob(PrimeField(F), [r[1] : r in Roots(f)]);
return [roots[j] : j in [1..#roots] | not -roots[j] in [roots[i]: i in [1..(j-1)]]];
end function;
// The possible traces of rho(ab) in F
TParamsForATrace := function(F, x, h)
return ModOutByFrob(sub, [z[1] : z in Roots(FlattenLiftPolyonomial(h, F, x))]);
end function;
/*
Here's all the points on the (geometric component) of character
variety T(n,k) over the field F. So we don't repeat ourselves, we
require that tr(rho(a)) and tr(rho(ab)) generate F.
*/
TraceCoorOfFiniteReps := function(n,k,F)
K, R, h := TraceHolonomyEquation(n,k);
ans := &cat [[[a, c] : c in TParamsForATrace(F, a, h)] : a in PossATraces(k,F)];
return [x : x in ans | #CommonOverfield(sub, sub) eq #F];
end function;
/*
The following function builds the matrices of the rep from the given
trace data. Over the algebraic closure of F, this realization is
essentially unique. (Unless the representation is irreducible, that
is, which corresponds to some polynomial in a and t vanishes. I don't
bother to check for this here since later we check if the rep surjects
onto PSL(2,F).)
*/
RepresentationFromTraces := function(a, t, F)
R := PolynomialRing(F);
for b1 in F do
roots := Roots(1 - a*b1 + b1^2 + b3^2 - b3*(-(a*b1) + t));
if #roots gt 0 then
r3 := roots[1][1];
b2 := a*b1 + r3 - t;
return Matrix([[a, One(F)],[-1, 0]]), Matrix([[b1, b2], [r3, a - b1]]);
end if;
end for;
assert false;
end function;
/*
Given the same info on traces, constructs the homomorphism from
TwistOrbifoldGroup to Sym(P^1(F)).
*/
TraceCreateRepresentation := function(G, F, a, t)
S := SymmetricGroup(#F + 1);
A, B := RepresentationFromTraces(a,t,F);
imagegens := [S ! ProjectiveAction(A, F), S ! ProjectiveAction(B,F)];
assert IsSatisfied(Relations(G), imagegens);
return hom< G -> S | imagegens>;
end function;
/*
Create all the _surjective_ reps TwistOrbifoldGroup -> PSL(2, F)
the test for surjectivity can take a while. Probably you could use
some theorem to essentially eliminate this check. I didn't bother
since computing the homology of the cover is the bottleneck, anyhow.
*/
TraceCreateRepresentations := function(n,k,q)
G := TwistOrbifoldGroup(n, k);
F := GF(q);
N := Order(PSL(2, F));
ans := [ TraceCreateRepresentation(G, F, x[1], x[2]) : x in TraceCoorOfFiniteReps(n,k,F)];
return [f : f in ans | Order(Image(f)) eq N ];
end function;
/*
Ok finally, here's the main function.
It computes the homology of the Gamma_0 congruence covers of T(n,k)
for primes with norm between minP and maxP (inclusive). The output is
concatenated to the file "filename" in the format:
[* norm of prime, [ list of bettis of congruence covers ], time spent *],
*/
SuperTestFile := procedure(n, k, minP, maxP, filename)
for p in [p : p in [minP..maxP] | IsPrimePower(p)] do
t := Cputime();
ranks := [HomologyOfProjCover(f) : f in TraceCreateRepresentations(n, k, p)];
t := Cputime(t);
if #ranks gt 0 then
Write(filename, FormatLine(p, ranks, t)) ;
end if;
end for;
end procedure;
/* Code for looking at other congruence covers of the example itself */
M0Group := function()
return Group;
end function;
M0PrimeGroup := function()
G := M0Group();
S := CyclicGroup(2);
imagegens := [One(S), One(S), S.1, S.1];
h := hom S| imagegens>;
return Simplify(Kernel(h));
end function;
M1Group := function()
return Group < a,b,c,d | a*b^-1*c^-1*b*a^-1*d*c*d^-1, a*b*a^-1*d*c*d*a^-2*b*c, a*d*c*d*a^-1*b*d^-2*c^-1*b^-1, c*d^2*c*d^2*c*d^2, c^3, a*c*b*c*a*b*d^-2 >;
end function;
M0MatrixRep := function()
R := PolynomialRing(Integers());
K := ext;
K := AbsoluteField(K);
pi := 1 - r;
pibar := 1 + r;
one := Matrix([[1,0], [0,1]]);
I := Matrix([[i,0], [0, -i]]);
J := Matrix([[0,1],[-3,0]]);
S := (1/2)*(I + J);
T := (1/2)*(one + I*J);
u := I;
v := -2 * I + pi * S;
x := r*one + pi*I + pibar * S - r*T;
y := r*one + S - r*T;
return K, [u,v,x,y];
end function;
M0ConguenceQuotient := function(G, mats, P)
F, f := ResidueClassField(P);
half := F ! 1/2;
mapmat := function(M)
return half*Matrix( [[ f(2*M[1,1]), f(2*M[1,2])], [f(2*M[2,1]), f(2*M[2,2])]]);
end function;
images := [mapmat(M) : M in mats];
S := SymmetricGroup(#F + 1);
imagegens := [S ! ProjectiveAction(A,F) : A in images];
assert IsSatisfied(Relations(G), imagegens);
return hom< G -> S | imagegens>;
end function;
/*
Prints beta_1 of Gamma_0(M_0, P) for one of the primes P in
Z[sqrt(-2)] sitting over each rational prime in [minP, maxP].
Only does for primes which split completely in the larger ring
Z[sqrt(-2), sqrt(-1) for technical reasons.
*/
FrankTest := procedure(minP, maxP)
G := M0Group();
K, mats := M0MatrixRep();
O := MaximalOrder(K);
for p in [p : p in [minP..maxP] | IsPrime(p)] do
facs := Decomposition(O, p);
if #facs eq 4 then
print p, HomologyOfProjCover(M0ConguenceQuotient(G, mats, facs[1][1]));
end if;
end for;
end procedure;
/* same as above, but does the primes P omitted there */
FrankTest2 := procedure(minP, maxP)
L := QuadraticField(-2);
G := M0Group();
K, mats := M0MatrixRep();
O := MaximalOrder(K);
for p in [p : p in [minP..maxP] | IsPrime(p)] do
facs := Decomposition(O, p);
if (#Decomposition(L, p) eq 2) and (#facs eq 2) then
f := M0ConguenceQuotient(G, mats, facs[1][1]);
Q := Image(f);
orbs := [o : o in Orbits(Q) | #o eq p + 1];
assert #orbs eq 1;
g := Action(Q, orbs[1]);
print p, HomologyOfProjCover(f*g);
end if;
end for;
end procedure;
/* Ok, now lets look at the Bianchi group for PGL(Z[sqrt(-2)]) */
BianchiGroup := function()
return Group;
end function;
BianchiMatrixRep := function()
K := QuadraticField(-2);
return K, [
Matrix( [[0, 1], [-1, 1]]),
Matrix( [[0,1],[1,0]]),
Matrix( [[0,1],[-1,0]]),
Matrix( [[r, 1], [1, 0]])];
end function;
BianchiConguenceQuotient := function(G, mats, P)
F, f := ResidueClassField(P);
mapmat := function(M)
return Matrix( [[ f(M[1,1]), f(M[1,2])], [f(M[2,1]), f(M[2,2])]]);
end function;
images := [mapmat(M) : M in mats];
S := SymmetricGroup(#F + 1);
imagegens := [S ! ProjectiveAction(A,F) : A in images];
assert IsSatisfied(Relations(G), imagegens);
return hom< G -> S | imagegens>;
end function;
FrankBianchiTest := procedure(minP, maxP)
G := BianchiGroup();
K, mats := BianchiMatrixRep();
O := MaximalOrder(K);
for p in [p : p in [minP..maxP] | IsPrime(p)] do
facs := Decomposition(O, p);
if #facs eq 2 then
print p,
HomologyOfProjCover(BianchiConguenceQuotient(G, mats, facs[1][1]));
end if;
end for;
end procedure;
FrankBianchi3Test := procedure(minP, maxP)
G := BianchiGroup();
K, mats := BianchiMatrixRep();
O := MaximalOrder(K);
f := BianchiConguenceQuotient(G, mats, Decomposition(O, 3)[1][1]);
g:= BianchiConguenceQuotient(G, mats, Decomposition(O, 3)[2][1]);
H1:= sub;
H2:= sub;
H := H1 meet H2;
for p in [p : p in [minP..maxP] | IsPrime(p)] do
facs := Decomposition(O, p);
if #facs eq 2 then
h := BianchiConguenceQuotient(G, mats, facs[1][1]);
N := sub;
print p, #AQInvariants(N, 31991), #AQInvariants(H1 meet N, 31991), #AQInvariants(H2 meet N, 31991), #AQInvariants(H meet N, 31991);
end if;
end for;
end procedure;