Fast Γ₀(N) coset methods for SL(3,Z), and extension to SL(4,Z)

lib/CongruenceGroups/congGroups.gi

@@ -101,7 +101,7 @@ InstallMethod( CongruenceSubgroupGamma0, "for integer matrix group and positive
     ##################################################
     ##
     ## Remaining components to be computed in other functions.
-    if DimensionOfMatrixGroup(GG)<4 then ProjectiveSpace(G);fi;
+    if DimensionOfMatrixGroup(GG)<5 then ProjectiveSpace(G);fi;
     AmbientTransversal(G);
     AmbientPosition(G); #A generic method will be used to construct the
     AmbientRepresentation(G); #these two functions except for cases with a
@@ -532,8 +532,10 @@ InstallMethod( IsSubset,
       vertex2word, one, cosetPos, cosetRep, triple2word,
       i,j,u,c,a,b;
 
-        if DimensionOfMatrixGroup(G)>2 then TryNextMethod(); fi;
-        n:=LevelOfCongruenceSubgroup(G);
+    if DimensionOfMatrixGroup(G)>2 and
+       not (DimensionOfMatrixGroup(G) in [3,4] and IsHAPCongruenceSubgroupGamma0(G))
+    then TryNextMethod(); fi;
+    n:=LevelOfCongruenceSubgroup(G);
 
      GG:=AmbientGroupOfCongruenceSubgroup(G);
      ambientGenerators:=GeneratorsOfGroup(GG);

lib/Orru/decsOrru.gd

@@ -4,9 +4,15 @@ DeclareGlobalFunction("HAP_GenericCongruenceSubgroup");
 DeclareOperation("FiniteProjectiveLine",[IsInt]);
 DeclareOperation("FiniteProjectivePlane",[IsInt]);
 DeclareOperation("FiniteProjectiveLine_alt",[IsInt]);
+DeclareOperation("FiniteProjectivePlane_alt",[IsInt]);
+DeclareOperation("FiniteProjectiveSpace3",[IsInt]);
+DeclareOperation("FiniteProjectiveSpace3_alt",[IsInt]);
 DeclareGlobalFunction("HAP_FiniteProjectiveLineIntegers");
 DeclareGlobalFunction("HAP_FiniteProjectiveLineIntegers_alt");
 DeclareGlobalFunction("HAP_FiniteProjectivePlaneIntegers");
+DeclareGlobalFunction("HAP_FiniteProjectivePlaneIntegers_alt");
+DeclareGlobalFunction("HAP_FiniteProjectiveSpace3Integers");
+DeclareGlobalFunction("HAP_FiniteProjectiveSpace3Integers_alt");
 DeclareGlobalFunction("HAP_SL3ZSubgroupTree_fast");
 DeclareOperation("HAPCongruenceSubgroupGamma0",[IsInt,IsInt]);
 DeclareOperation("HAPCongruenceSubgroupGamma0_alt",[IsInt,IsInt]);

lib/Orru/finiteProjectiveSpaces.gi

@@ -27,7 +27,26 @@ function(n)
 return HAP_FiniteProjectiveLineIntegers_alt(n);
 end);
 
+InstallMethod(FiniteProjectivePlane_alt,
+"Finite projective line for the ring Z/nZ",
+[IsInt],
+function(n)
+return HAP_FiniteProjectivePlaneIntegers_alt(n);
+end);
+
+InstallMethod(FiniteProjectiveSpace3_alt,
+"3 dimensional finite projective space for the ring Z/nZ",
+[IsInt],
+function(n)
+return HAP_FiniteProjectiveSpace3Integers_alt(n);
+end);
 
+InstallMethod(FiniteProjectiveSpace3,
+"3 dimensional finite projective space for the ring Z/nZ",
+[IsInt],
+function(n)
+return HAP_FiniteProjectiveSpace3Integers(n);
+end);
 
 
 InstallGlobalFunction(HAP_FiniteProjectiveLineIntegers_alt,
@@ -72,7 +91,7 @@ function(n)
     );
 end);
 
-InstallGlobalFunction(HAP_FiniteProjectivePlaneIntegers,
+InstallGlobalFunction(HAP_FiniteProjectivePlaneIntegers_alt,
 function(n)
     local UnitEls, x, y, z, i, c, d, u, UnitsAction, Representatives, 
           RepOf, r, v, w, m, min;
@@ -160,4 +179,148 @@ function(n)
     return Rep;
 end);
 
+InstallGlobalFunction(HAP_FiniteProjectivePlaneIntegers,
+function(n)
+    local divs, Rep, d, q, g, q_g, m, u, U_m, dd, qq, z, uu, zz;
+
+    divs := List(DivisorsInt(n));
+    Remove(divs);
+    divs := Concatenation([n],divs);
+    Rep := [];
+
+    for d in divs do
+      q := n/d;
+      for g in divs do
+        q_g := n/g;
+        m := Gcd(q,q_g);
+        U_m := Filtered([0..m-1], z -> Gcd(z,m) = 1);
+
+        for u in U_m do
+          uu := u;          
+          while Gcd(uu, q_g) <> 1 do
+            uu := uu + m;
+          od;
+
+          dd := Gcd(d,g);
+          qq := n/dd;
+
+          for z in [0..qq-1] do
+            if Gcd(z,dd,qq) = 1 then
+              zz := z;
+              while Gcd(zz, dd) <> 1 do
+                zz := zz + qq;
+              od;
 
+              Add(Rep, [d mod n,g*uu mod n,zz]);
+            fi;
+          od;
+        od;
+      od;
+    od;
+
+    return rec(Reps := Set(List(Rep, Immutable)));;
+end);
+
+InstallGlobalFunction(HAP_FiniteProjectiveSpace3Integers,
+function(n)
+    local divs, Rep, d_x, q_x, d_y, q_y, q_xy, d_z, q_z, Q, y_0, y_00, Y_0, z_0, z_00, Z_0, QQ, t, tt;
+
+    divs := List(DivisorsInt(n));
+    Remove(divs);
+    divs := Concatenation([n],divs);
+    Rep := [];
+
+    for d_x in divs do
+      q_x := n/d_x;
+      for d_y in divs do
+        q_y := n/d_y;
+        q_xy := Gcd(q_x,q_y);
+        Y_0 := Filtered([0..q_xy-1], e -> Gcd(e,q_xy) = 1);
+        Q := Lcm(q_x,q_y);
+
+        for y_0 in Y_0 do
+          y_00 := y_0;          
+          while Gcd(y_00, q_y) <> 1 do
+            y_00 := y_00 + q_xy;
+          od;
+
+          for d_z in divs do
+            q_z := n/d_z;
+            QQ := Lcm(Q,q_z);
+
+            Z_0 := Filtered([0..Gcd(Q,q_z)-1], e -> Gcd(e,Gcd(Q,q_z)) = 1);
+            for z_0 in Z_0 do
+              z_00 := z_0;          
+              while Gcd(z_00, q_z) <> 1 do
+                z_00 := z_00 + Gcd(Q,q_z);
+              od;
+              for t in [0..QQ-1] do
+                if Gcd(t,QQ,n/QQ) = 1 then
+                  tt := t;
+                  while Gcd(tt, n/QQ) <> 1 do
+                    tt := tt + QQ;
+                  od;
+                  Add(Rep, [d_x mod n,d_y*y_00 mod n,d_z*z_00 mod n, tt]);
+                fi;
+              od;
+            od;
+          od;
+        od;
+      od;
+    od;
+
+    return rec(Reps := Set(List(Rep, Immutable)));;
+end);
+
+
+InstallGlobalFunction(HAP_FiniteProjectiveSpace3Integers_alt,
+function(n)
+    local UnitEls, x, y, z, t, i, c, d, u, UnitsAction, Representatives, 
+          RepOf, r, v, w, m, min;
+
+    UnitEls := Units(Integers mod n);
+
+    UnitsAction := function(c, u)
+      local uu;
+      uu :=Int (u);
+      return List(c, x -> (uu * x) mod n);
+    end;
+
+    Representatives := [];
+
+    RepOf := [];
+    for x in [1..n] do
+      RepOf[x] := [];
+      for y in [1..n] do
+        RepOf[x][y] := [];
+        for z in [1..n] do
+          RepOf[x][y][z] := [];
+        od;
+      od;
+    od;
+
+    for x in [0..n-1] do
+      for y in [0..n-1] do
+        for z in [0..n-1] do
+          for t in [0..n-1] do
+            if Gcd(x,y,z,t,n) = 1 then
+              v := [x,y,z,t];
+              if not IsBound(RepOf[x+1][y+1][z+1][t+1]) then
+                m := Orbit(UnitEls,v,UnitsAction);
+                min := Minimum(m);
+                AddSet(Representatives,min);
+                for w in m + 1 do
+                  RepOf[w[1]][w[2]][w[3]][w[4]] := min;
+                od;
+              fi;
+            fi;
+          od;
+        od;
+      od;
+    od;
+
+    return rec(
+        Reps := Set(Representatives),
+        RepOf:= RepOf
+    );
+end);

lib/Orru/initOrru.gi

@@ -4,3 +4,4 @@ ReadPackageHap( "lib/Orru/HAP_SL3ZSubgroupTree_fast.gi");
 ReadPackageHap( "lib/Orru/SL3ZTopHomology.gi");
 ReadPackageHap( "lib/Orru/sl2methodsExtra.gi");
 ReadPackageHap( "lib/Orru/sl3methods.gi");
+ReadPackageHap( "lib/Orru/sl4methods.gi");

lib/Orru/sl2methodsExtra.gi

@@ -1,137 +1,4 @@
-# methods for power of primes
-
 InstallMethod(AmbientPosition,
-    "Returns cosetPos(g) function for the congruence subgroup G",
-    [ IsIntegerMatrixGroup and IsHAPCongruenceSubgroupGamma0 ],
-    function(G)
-        local cosetPos, canonicalRep, n, countRep;
-
-        if DimensionOfMatrixGroup(G) <> 2 then
-            TryNextMethod();
-        fi;
-
-        n := LevelOfCongruenceSubgroup(G);
-
-        if not IsPrimePowerInt(n) then
-            TryNextMethod();
-        fi;
-
-        if IsPrime(n) then
-            TryNextMethod();
-        fi;
-
-        canonicalRep := function(g)
-            local v, vv, U, d, dd, x, y;
-            v := [g[1][1], g[2][1]];
-            vv := List(v, x -> x mod n);
-            U := Filtered([0..n],i->Gcd(i,n)=1);
-            if vv[1] mod n = 0 then
-                return [0,1];
-            elif vv[1] mod n in U then
-                return [1,(Inverse(vv[1]) mod n)*vv[2] mod n];
-            else
-                d := Gcd(vv[1],n);
-                dd := n/d;
-                x := vv[1]/d;
-                y := vv[2]/x mod dd;
-                while not Gcd(d,y) = 1 do
-                    y := y + dd;
-                od;
-                return [d, y];
-            fi;
-        end;
-
-        countRep := function(m)
-            local p, e, i, countp;
-
-            p := Set(Factors(m))[1];
-            e := Length(Factors(m));
-
-            countp := [1,m];
-
-            for i in [2..e] do
-                Add(countp, p^(e-i)*(p-1));
-            od;
-
-            return countp;
-        end;
-
-        cosetPos := function(g)
-            local w, count, e, U;
-            w := canonicalRep(g);
-
-            count := countRep(n);
-
-            if w[1] = 0 then
-                return 1;
-            elif w[1] = 1 then
-                U := [0..n-1];
-                return 1 + Position(U,w[2]);
-            else
-                e := Length(Factors(w[1]));
-                U := Filtered([0..n/w[1]], i-> Gcd(i,n/w[1]) = 1);
-                return Sum(count{[1..1 + e]}) + Position(U, w[2]);
-            fi;
-        end;
-
-        return cosetPos;
-    end);
-
-#InstallMethod(AmbientPosition,
-#    "Returns cosetPos(g) function for the congruence subgroup G",
-#    [ IsIntegerMatrixGroup and IsHAPCongruenceSubgroupGamma0 ],
-#    function(G)
-#        local cosetPos, canonicalRep, n, ProjLine, U;
-#
-#        if DimensionOfMatrixGroup(G) <> 2 then
-#            TryNextMethod();
-#        fi;
-#
-#        n := LevelOfCongruenceSubgroup(G);
-#
-#        if IsPrime(n) then
-#            TryNextMethod();
-#        fi;
-#
-#        if IsPrimePowerInt(n) then
-#            TryNextMethod();
-#        fi;
-#
-#        ProjLine := ProjectiveSpace(G);
-#
-#        U := Filtered([0..n],i -> Gcd(i,n) = 1);
-#
-#        canonicalRep := function(g)
-#            local v, vv, d, dd, x, y;
-#            v := [g[1][1], g[2][1]];
-#            vv := List(v, x -> x mod n);
-#            if vv[1] mod n = 0 then
-#                return [0,1];
-#            elif (vv[1] mod n) in U then
-#                return [1,(Inverse(vv[1]) mod n)*vv[2] mod n];
-#            else
-#                d := Gcd(vv[1],n);
-#                dd := n/d;
-#                x := vv[1]/d;
-#                y := vv[2]/x mod dd;
-#                while not Gcd(d,y) = 1 do
-#                    y := y + dd;
-#                od;
-#                return [d, y];
-#            fi;
-#        end;
-#
-#        cosetPos := function(g)
-#            local w;
-#            w := canonicalRep(g);
-#            return Position(ProjLine,w);
-#        end;
-#
-#        return cosetPos;
-#    end);
-
-## New implementation just to try, with a better formula for the position.
-    InstallMethod(AmbientPosition,
     "Returns cosetPos(g) function for the congruence subgroup G",
     [ IsIntegerMatrixGroup and IsHAPCongruenceSubgroupGamma0 ],
     function(G)
@@ -147,10 +14,6 @@ InstallMethod(AmbientPosition,
             TryNextMethod();
         fi;
 
-        if IsPrimePowerInt(n) then
-            TryNextMethod();
-        fi;
-
         countRep := function(m)
             local d, q, count;