ENHANCE: Faster automorphism order computation

by working in suitable hom-stable factors first.

lib/morpheus.gi

@@ -16,36 +16,111 @@
 ##
 MORPHEUSELMS := 50000;
 
+# this method calculates a chief series invariant under `hom` and calculates
+# orders of group elements in factors of this series under action of `hom`.
+# Every time an orbit length is found, `hom` is replaced by the appropriate
+# power. Initially small chief factors are preferred. In the end all
+# generators are used while stepping through the series descendingly, thus
+# ensuring the proper order is found.
 InstallMethod(Order,"for automorphisms",true,[IsGroupHomomorphism],0,
 function(hom)
-local map,phi,o,lo,i,start,img;
+local map,phi,o,lo,i,j,start,img,d,nat,ser,jord,first;
+  d:=Source(hom);
+  if Size(d)<=10000 then
+    ser:=[d,TrivialSubgroup(d)]; # no need to be clever if small
+  else
+    if HasAutomorphismGroup(d) then
+      if IsBound(d!.characteristicSeries) then
+        ser:=d!.characteristicSeries;
+      else
+        ser:=ChiefSeries(d); # could try to be more clever, introduce attribute
+        # `CharacteristicSeries`.
+        ser:=Filtered(ser,x->ForAll(GeneratorsOfGroup(AutomorphismGroup(d)),
+          a->ForAll(GeneratorsOfGroup(x),y->ImageElm(a,y) in x)));
+        d!.characteristicSeries:=ser;
+      fi;
+    else
+      ser:=ChiefSeries(d); # could try to be more clever, introduce attribute
+      # `CharacteristicSeries`.
+      ser:=Filtered(ser,
+            x->ForAll(GeneratorsOfGroup(x),y->ImageElm(hom,y) in x));
+    fi;
+  fi;
+
+  # try to do factors in ascending order in the hope to get short orbits
+  # first
+  jord:=[2..Length(ser)]; # order in which we go through factors
+  if Length(ser)>2 then
+    i:=List(jord,x->Size(ser[x-1])/Size(ser[x]));
+    SortParallel(i,jord);
+  fi;
+
   o:=1;
   phi:=hom;
   map:=MappingGeneratorsImages(phi);
-  i:=1;
-  while i<=Length(map[1]) do
-    lo:=1;
-    start:=map[1][i];
-    img:=map[2][i];
-    while img<>start do
-      img:=ImagesRepresentative(phi,img);
-      lo:=lo+1;
-      # do the bijectivity test only if high local order, then it does not
-      # matter
-      if lo=1000 and not IsBijective(hom) then
-	Error("<hom> must be bijective");
-      fi;
+
+  first:=true;
+  while map[1]<>map[2] do
+    for j in jord do
+      i:=1;
+      while i<=Length(map[1]) do
+        # the first time, do only the generators from prior layer
+        if (not first) 
+           or (map[1][i] in ser[j-1] and not map[1][i] in ser[j]) then
+
+          lo:=1;
+          if j<Length(ser) then
+            nat:=NaturalHomomorphismByNormalSubgroup(d,ser[j]);
+            start:=ImagesRepresentative(nat,map[1][i]);
+            img:=map[2][i];
+            while ImagesRepresentative(nat,img)<>start do
+              img:=ImagesRepresentative(phi,img);
+              lo:=lo+1;
+
+              # do the bijectivity test only if high local order, then it
+              # does not matter. IsBijective is cached, so second test is
+              # cheap.
+              if lo=1000 and not IsBijective(hom) then
+                Error("<hom> must be bijective");
+              fi;
+
+            od;
+
+          else
+            start:=map[1][i];
+            img:=map[2][i];
+            while img<>start do
+              img:=ImagesRepresentative(phi,img);
+              lo:=lo+1;
+
+              # do the bijectivity test only if high local order, then it
+              # does not matter. IsBijective is cached, so second test is
+              # cheap.
+              if lo=1000 and not IsBijective(hom) then
+                Error("<hom> must be bijective");
+              fi;
+
+            od;
+          fi;
+
+          if lo>1 then
+            o:=o*lo;
+            #if i<Length(map[1]) then
+              phi:=phi^lo;
+              map:=MappingGeneratorsImages(phi);
+              i:=0; # restart search, as generator set may have changed.
+            #fi;
+          fi;
+        fi;
+        i:=i+1;
+      od;
     od;
-    if lo>1 then
-      o:=o*lo;
-      if i<Length(map[1]) then
-	phi:=phi^lo;
-	map:=MappingGeneratorsImages(phi);
-        i:=0; # restart search, as generator set may have changed.
-      fi;
-    fi;
-    i:=i+1;
+
+    # if iterating make `jord` standard to we don't skip generators
+    jord:=[2..Length(ser)];
+    first:=false;
   od;
+
   return o;
 end);