Spracherkennung für: .gi vermutete Sprache: Unknown {[0] [0] [0]} [Methode: Schwerpunktbildung, einfache Gewichte, sechs Dimensionen]
#############################################################################
##
## genss.gi genss package
## Max Neunhoeffer
## Felix Noeske
##
## Copyright 2006 Lehrstuhl D für Mathematik, RWTH Aachen
##
## Implementation stuff for generic Schreier-Sims
##
#############################################################################
#############################################################################
# The following global record contains default values for options for the
# main function "StabilizerChain":
#############################################################################
# Initial hash size for orbits:
GENSS.InitialHashSize := NextPrimeInt(1000);
# Number of points to process before reporting:
GENSS.Report := 30000;
# Number of random elements to consider for the determination of short orbits:
GENSS.ShortOrbitsNrRandoms := 10;
# Absolute limit for orbit length in search for short orbits:
GENSS.ShortOrbitsOrbLimit := 80000;
# First limit for the parallel enumeration of orbits looking for short orbits:
GENSS.ShortOrbitsInitialLimit := 400;
# Absolute limit for single orbit length:
GENSS.OrbitLengthLimit := 10000000;
# Number of points in the previous orbit to consider for the next base point:
GENSS.NumberPrevOrbitPoints := 10;
# Number of (evenly distributed) random generators for the stabilizer:
GENSS.RandomStabGens := 3;
# Product replacement parameters for the stabilizer element generation:
# Now actually used for the generation of all random elements on top
# level. Random elements further down are created on top and then
# sifted.
GENSS.StabGenScramble := 30;
GENSS.StabGenScrambleFactor := 6;
GENSS.StabGenAddSlots := 3;
GENSS.StabGenMaxDepth := 400;
# Number of random elements used for verification,
# note that this is changed by the "random" and "ErrorBound" options!
GENSS.VerifyElements := 10; # this amounts to an error bound of 1/1024
GENSS.DeterministicVerification := false;
# Number of random elements used to do immediate verification:
GENSS.ImmediateVerificationElements := 3;
# Are we working projectively?
GENSS.Projective := false;
# To find a very short orbit we try two basis vectors and a random vector:
GENSS.VeryShortOrbLimit := 500;
# Never consider more than this number of candidates for short orbits:
GENSS.LimitShortOrbCandidates := 50;
# Do not throw Errors but return fail:
GENSS.FailInsteadOfError := false;
# Number of Schreier generators to create in TC verification:
GENSS.NumberSchreierGens := 20;
# Maximal number of Schreier generators to create in TC verification:
GENSS.MaxNumberSchreierGens := 100;
# By default do 3 rounds of birthday paradox method:
GENSS.TryBirthdayParadox := 3;
# By default do 1 short orbit tries:
GENSS.TryShortOrbit := 1;
# Limit for orbit length during orbit estimation using birthday paradox:
GENSS.OrbitLimitBirthdayParadox := 1000000;
# We immediately take an orbit if its estimate is lower than this limit:
GENSS.OrbitLimitImmediatelyTake := 10000;
# Limit for number of random elements during orbit estimation:
GENSS.NrRandElsBirthdayParadox := 6000;
# Set this to false to allow orbits of length 1:
GENSS.Reduced := true;
# Product replacement parameters for Stab:
GENSS.StabScramble := 10;
GENSS.StabScrambleFactor := 1;
GENSS.StabAddSlots := 1;
GENSS.StabMaxDepth := 400;
# Initial limit for orbit length for Stab computation:
GENSS.StabInitialLimit := 1000;
# Patience to find random elements mapping the orbit into itself:
GENSS.StabInitialPatience := 10;
# Maximal amount of memory used for the orbit:
GENSS.StabOrbitLimit := 1000000;
# If the probability of a wrong stabiliser is smaller than this, do no
# longer try to create stabiliser elements:
GENSS.StabAssumeCompleteLimit := 1/(10^7);
# The following switches off the log used in all orbits for stabilizer
# chains if set to false. Do not do this unless you know what you are
# doing since since could prevent Schreier trees from being shallow.
GENSS.OrbitsWithLog := true;
#############################################################################
# A few helper functions needed elsewhere:
#############################################################################
InstallGlobalFunction( GENSS_CopyDefaultOptions,
function( defopt, opt )
local n;
for n in RecNames(defopt) do
if not(IsBound(opt.(n))) then
opt.(n) := defopt.(n);
fi;
od;
end );
InstallGlobalFunction( GENSS_MapBaseImage,
function( bi, el, S )
# bi must be a base image belonging to S, el a group element
local i,l,res;
l := Length(bi);
res := 0*[1..l];
i := 1;
while true do
res[i] := S!.orb!.op(bi[i],el);
i := i + 1;
if i = l+1 then break; fi;
S := S!.stab;
od;
return res;
end );
InstallGlobalFunction( GENSS_FindVectorsWithShortOrbit,
# This implements Murray/O'Brien-like heuristics.
# It produces new random elements using GENSS_RandomElementFromAbove
# and stores them in opt.FindBasePointCandidatesData.randpool for
# later usage.
function(g,opt,parentS)
# Needs opt components "ShortOrbitsNrRandoms"
local l, f, data, x, c, onlydegs, v, vv, w, i, nw, inters, sortfun,
wb, j, ww;
Info(InfoGenSS,4,"Trying Murray/O'Brien heuristics...");
l := ShallowCopy(GeneratorsOfGroup(g));
f := DefaultFieldOfMatrixGroup(g);
data := opt.FindBasePointCandidatesData;
for i in [1..opt.ShortOrbitsNrRandoms] do
if parentS = false then
x := GENSS_RandomElementFromAbove(opt,0);
else
x := GENSS_RandomElementFromAbove(parentS,parentS!.layer);
fi;
Add(l,x);
Add(data.randpool,x);
od;
ForgetMemory(l);
c := List(l,x->Set(Factors(CharacteristicPolynomial(x,1):
onlydegs := [1..3])));
v := [];
for i in [1..Length(l)] do
for j in [1..Length(c[i])] do
vv := EmptyPlist(Length(c[i]));
Add(vv,[NullspaceMat(Value(c[i][j],l[i])),
Degree(c[i][j]),
WeightVecFFE(CoefficientsOfLaurentPolynomial(c[i][j])[1]),
1]);
od;
Add(v,vv);
od;
Info(InfoGenSS,4,"Have eigenspaces.");
# Now collect a list of all those spaces together with all
# possible intersects:
w := [];
i := 1;
while i <= Length(l) and Length(w) < GENSS.LimitShortOrbCandidates do
nw := [];
for j in [1..Length(v[i])] do
for ww in w do
inters := SumIntersectionMat(ww[1],v[i][j][1])[2];
if Length(inters) > 0 then
Add(nw,[inters,Minimum(ww[2],v[i][j][2]),
Minimum(ww[3],v[i][j][3]),ww[4]+v[i][j][4]]);
fi;
od;
Add(nw,v[i][j]);
od;
Append(w,nw);
i := i + 1;
od;
sortfun := function(a,b)
if a[2] < b[2] then return true;
elif a[2] > b[2] then return false;
elif a[3] < b[3] then return true;
elif a[3] > b[3] then return false;
elif a[4] < b[4] then return true;
elif a[4] > b[4] then return false;
elif Length(a[1]) < Length(b[1]) then return true;
else return false;
fi;
end;
Sort(w,sortfun);
wb := List(w,ww->ww[1][1]);
Info(InfoGenSS,3,"Have ",Length(wb)," vectors for possibly short orbits.");
for ww in [1..Length(wb)] do
if not(IsMutable(wb[ww])) then
wb[ww] := ShallowCopy(wb[ww]);
fi;
ORB_NormalizeVector(wb[ww]);
od;
return wb;
end );
InstallGlobalFunction( GENSS_FindShortOrbit,
function( g, opt, parentS )
# Needs opt components:
# "ShortOrbitsNrRandoms" (because it uses GENSS_FindVectorsWithShortOrbit)
# "ShortOrbitsOrbLimit"
# "ShortOrbitsInitialartLimit"
local ThrowAwayOrbit,found,gens,hashlen,i,j,limit,newnrorbs,nrorbs,o,wb;
wb := GENSS_FindVectorsWithShortOrbit(g,opt,parentS);
if Length(wb) = 0 then return fail; fi;
# Now we have a list of vectors with (hopefully) short orbits.
# We start enumerating all those orbits, but first only 50 elements:
nrorbs := Minimum(Length(wb),64); # take only the 64 first
gens := GeneratorsOfGroup(g);
o := [];
hashlen := NextPrimeInt(QuoInt(opt.ShortOrbitsOrbLimit,2));
for i in [1..nrorbs] do
Add(o,Orb(gens,ShallowCopy(wb[i]),OnLines,
rec( treehashsize := hashlen )));
od;
limit := opt.ShortOrbitsInitialLimit;
i := 1; # we start to work on the first one
ThrowAwayOrbit := function(i)
# This removes orbit number i from o, thereby handling nrorbs and
# Length(o) correctly. If you want to use o[i] further, please
# make a copy (of the reference) before calling this function.
if Length(o) > nrorbs then
o[i] := o[nrorbs+1];
o{[nrorbs+1..Length(o)-1]} := o{[nrorbs+2..Length(o)]};
Unbind(o[Length(o)]);
else
o{[i..nrorbs-1]} := o{[i+1..nrorbs]};
Unbind(o[nrorbs]);
nrorbs := nrorbs-1;
fi;
end;
repeat
Enumerate(o[i],limit);
found := IsClosed(o[i]);
if Length(o[i]) = 1 then
Info(InfoGenSS,3,"Orbit Number ",i," has length 1.");
found := false;
# Now throw away this orbit:
ThrowAwayOrbit(i);
# we intentionally do not increase i here!
elif not(found) then
i := i + 1;
fi;
if i > nrorbs then
Info(InfoGenSS,3,"Done ",nrorbs,
" orbit(s) to limit ",limit,".");
limit := limit * 2;
if limit > opt.ShortOrbitsOrbLimit then
Info(InfoGenSS,3,"Limit reached, giving up.");
return fail;
fi;
i := 1;
if nrorbs < i then
Info(InfoGenSS,3,"No orbits left, giving up.");
return fail;
fi;
if nrorbs > 1 then
newnrorbs := QuoInt((nrorbs+1),2);
for j in [newnrorbs+1..nrorbs] do
Unbind(o[j]);
od;
nrorbs := newnrorbs;
fi;
fi;
until found;
Info(InfoGenSS,2,"Found short orbit of length ",Length(o[i])," (#",i,").");
return o[i];
end );
InstallGlobalFunction( GENSS_IsOneProjective,
function(el)
local s;
s := el[1][1];
if IsZero(s) then return false; fi;
if not(IsOne(s)) then
s := s^-1;
el := s * el;
fi;
return IsOne( el );
end );
#############################################################################
# Now to the heart of the method, the Schreier-Sims:
#############################################################################
InstallMethod( FindBasePointCandidates,
"for a scalar matrix group over a FF",
[ IsGroup and IsMatrixGroup and IsFinite, IsRecord, IsInt, IsObject ],
50, # highest weight,
function( grp, opt, i, parentS )
local F, q, d, gens, v;
Info( InfoGenSS, 3, "Finding nice base points (scalar)..." );
F := DefaultFieldOfMatrixGroup(grp);
q := Size(F);
d := DimensionOfMatrixGroup(grp);
gens := GeneratorsOfGroup(grp);
if IsObjWithMemory(gens[1]) then
gens := StripMemory(gens);
fi;
if ForAny(gens,x->not(GENSS_IsOneProjective(x))) then
opt.FindBasePointCandidatesData := rec( randpool := [], vecs := [] );
# This is needed for communication between different methods
TryNextMethod();
fi;
v := ZeroMutable(gens[1][1]);
v[1] := PrimitiveRoot(F); # to indicate the field!
return rec( points := [v], ops := [OnRight], used := 0 );
end );
InstallMethod( FindBasePointCandidates,
"for a matrix group over a FF, very short orbit",
[ IsGroup and IsMatrixGroup and IsFinite, IsRecord, IsInt, IsObject ],
40, # highest weight,
function( grp, opt, i, parentS )
local F, q, d, gens, op, v, vv, k, kk, o, cand, j;
Info( InfoGenSS, 3, "Finding nice base points (very short)..." );
F := DefaultFieldOfMatrixGroup(grp);
q := Size(F);
d := DimensionOfMatrixGroup(grp);
gens := GeneratorsOfGroup(grp);
if IsObjWithMemory(gens[1]) then
gens := StripMemory(gens);
fi;
# Try two standard basis vectors and a random vector to find a very
# short orbit:
if q = 2 then
op := OnPoints;
else
op := OnLines;
fi;
v := [];
# Find first standard basis vector that is moved:
vv := ZeroMutable( gens[1][1] );
k := 1;
while k <= d do
vv[k] := One(F);
if ForAny(gens,x->op(vv,x) <> vv) then
Add(v,vv);
break;
fi;
vv[k] := Zero(F);
k := k + 1;
od;
# Find last standard basis vector that is moved:
vv := ZeroMutable( gens[1][1] );
kk := d;
while kk >= k do
vv[kk] := One(F);
if ForAny(gens,x->op(vv,x) <> vv) then
Add(v,vv);
break;
fi;
vv[kk] := Zero(F);
kk := kk - 1;
od;
# Pick a random vector:
vv := ZeroMutable( gens[1][1] );
if IsPlistRep(vv) then
for j in [1..Length(vv)] do
vv[j] := Random(F);
od;
else
Randomize(vv);
fi;
ORB_NormalizeVector(vv);
Add(v,vv);
# Now investigate these up to a certain limit:
for j in [1..Length(v)] do
o := Orb(gens,v[j],op,
rec(treehashsize := QuoInt(opt.VeryShortOrbLimit,2)+1));
Enumerate(o,opt.VeryShortOrbLimit);
if Length(o) > 1 and Length(o) < opt.VeryShortOrbLimit then
Info( InfoGenSS, 3, "Found orbit of length ",Length(o) );
cand := rec( points := [v[j]], ops := [op], used := 0 );
# Note that if we work non-projectively, then the same
# point will be taken next with OnRight automatically!
return cand;
fi;
od;
Info( InfoGenSS, 3, "Found no very short orbit up to limit ",
opt.VeryShortOrbLimit );
Append(opt.FindBasePointCandidatesData.vecs,v); # hand on vectors
TryNextMethod();
end );
# Method with rank 30 taking some commutators and invariant spaces
InstallMethod( FindBasePointCandidates,
"for a matrix group over a FF, using birthday paradox method",
[ IsGroup and IsMatrixGroup and IsFinite, IsRecord, IsInt, IsObject ], 20,
function( grp, opt, mode, parentS )
local F, q, d, randels, immorblimit, orblimit, data, op, v, l, c, e, ht,
val, x, w, cand, minest, minpos, round, i, j, gens;
F := DefaultFieldOfMatrixGroup(grp);
q := Size(F);
d := DimensionOfMatrixGroup(grp);
randels := opt.NrRandElsBirthdayParadox;
immorblimit := opt.OrbitLimitImmediatelyTake;
orblimit := opt.OrbitLimitBirthdayParadox;
Info( InfoGenSS, 3, "Finding base points (birthday paradox, limit=",
orblimit,", randels=",randels,")..." );
data := opt.FindBasePointCandidatesData; # this we get from earlier methods
if q = 2 then
op := OnPoints;
else
op := OnLines;
fi;
gens := GeneratorsOfGroup(grp);
if IsObjWithMemory(gens[1]) then
gens := StripMemory(gens);
fi;
for round in [1..opt.TryBirthdayParadox] do
v := Set(GENSS_FindVectorsWithShortOrbit(grp,opt,parentS));
if round = 1 then
Append(v,data.vecs); # take previously tried ones as well
fi;
v := Filtered(v,vv->ForAny(gens,x-> vv <> op(vv,x)));
l := Length(v);
if l = 0 then
# Find a vector on which grp acts non-trivially.
# At least one basis vector is guaranteed to be moved,
# unless grp is diagonal and op is OnLines, in which case
# they might all be fixed. However, in that case, the
# vector [1,1,...,1] is moved.
v := OneMutable(gens[1]); # List of basis vectors
Add(v, Sum(v)); # add vector [1,1,...,1]
v := Filtered(v,vv->ForAny(gens,x-> vv <> op(vv,x)));
l := Length(v);
fi;
c := 0*[1..l]; # the number of coincidences
e := ListWithIdenticalEntries(l,infinity); # the current estimates
ht := HTCreate(v[1]*PrimitiveRoot(F),
rec(hashlen := NextPrimeInt(l * randels * 4)));
for i in [1..l] do
val := HTValue(ht,v[i]);
if val = fail then
HTAdd(ht,v[i],[i]);
else
AddSet(val,i);
fi;
od;
for i in [1..randels] do
if parentS = false then
x := GENSS_RandomElementFromAbove(opt,0);
else
x := GENSS_RandomElementFromAbove(parentS,parentS!.layer);
fi;
Add(data.randpool,x);
for j in [1..l] do
if IsObjWithMemory(x) then
w := op(v[j],x!.el);
else
w := op(v[j],x);
fi;
val := HTValue(ht,w);
if val <> fail then # we know this point!
if j in val then # a coincidence!
c[j] := c[j] + 1;
e[j] := QuoInt(i^2,2*c[j]);
if (c[j] >= 3 and e[j] <= immorblimit) or
(c[j] >= 15 and e[j] <= orblimit) then
Info( InfoGenSS, 2, "Found orbit with estimated ",
"length ",e[j]," (coinc=",c[j],")" );
cand := rec(points := [v[j]], ops := [op],
used := 0);
for i in [1..l] do
if i <> j and c[i] >= 10 and
e[i] <= orblimit then
Add(cand.points,v[i]);
Add(cand.ops,op);
fi;
od;
if Length(cand.points) > 1 then
Info( InfoGenSS, 2, "Adding ",
Length(cand.points)-1, " more vectors",
" to candidates.");
fi;
return cand;
fi;
else
AddSet(val,j);
fi;
else
HTAdd(ht,w,[j]);
fi;
od;
od;
minest := Minimum(e);
minpos := Position(e,minest);
Info( InfoGenSS,2,"Birthday #", round, ": no small enough estimate. ",
"MinEst=",minest," Coinc=",c[j] );
randels := randels * 2;
orblimit := orblimit * 4;
od;
TryNextMethod();
end );
InstallMethod( FindBasePointCandidates,
"for a matrix group over a FF, original try short orbit method",
[ IsGroup and IsMatrixGroup and IsFinite, IsRecord, IsInt, IsObject ], 10,
function( grp, opt, i, parentS )
local F, d, data, cand, res;
F := DefaultFieldOfMatrixGroup(grp);
d := DimensionOfMatrixGroup(grp);
data := opt.FindBasePointCandidatesData;
Info( InfoGenSS, 3, "Finding nice base points (TryShortOrbit)..." );
# Next possibly "TryShortOrbit":
cand := rec( points := [], used := 0 );
if IsBound(opt.TryShortOrbit) and opt.TryShortOrbit > 0 then
repeat
opt.TryShortOrbit := opt.TryShortOrbit - 1;
Info(InfoGenSS,1,"Looking for short orbit (",opt.TryShortOrbit,
")...");
res := GENSS_FindShortOrbit(grp,opt,parentS);
until res <> fail or opt.TryShortOrbit = 0;
if res <> fail then
if Size(F) > 2 then
cand.points := [res[1]];
cand.ops := [OnLines];
# Note that if we work non-projectively, then the same
# point will be taken next with OnRight automatically!
else
cand.points := [res[1]];
cand.ops := [OnRight];
fi;
return cand;
fi;
fi;
TryNextMethod();
end );
InstallMethod( FindBasePointCandidates,
"for a matrix group over a FF, traditional Murray/O'Brien",
[ IsGroup and IsMatrixGroup and IsFinite, IsRecord, IsInt, IsObject ],
function( grp, opt, i, parentS )
local F, d, bv, cand, w, v;
F := DefaultFieldOfMatrixGroup(grp);
d := DimensionOfMatrixGroup(grp);
Info( InfoGenSS, 3, "Finding nice base points (Murray/O'Brien)..." );
# Standard Murray/O'Brien heuristics:
if i = 0 and
((opt!.Projective = false and Size(F)^d > 300000) or
(opt!.Projective = true and Size(F)^(d-1) > 300000)) then
bv := GENSS_FindVectorsWithShortOrbit(grp,opt,parentS);
if Length(bv) < 3 then
bv := MutableCopyMat(One(grp));
fi;
bv := bv{[1..3]}; # just take 3 of them
else
bv := One(grp);
fi;
cand := rec( points := [], ops := [], used := 0 );
for v in bv do
w := ORB_NormalizeVector(ShallowCopy(v));
if Size(F) = 2 then
Add(cand.points,w);
Add(cand.ops,OnRight);
else
Add(cand.points,w);
Add(cand.ops,OnLines);
# Note that if we work non-projectively, then the same
# point will be taken next with OnRight automatically!
fi;
od;
return cand;
end );
InstallMethod( FindBasePointCandidates, "for a permutation group",
[ IsGroup and IsPermGroup, IsRecord, IsInt, IsObject ],
function( grp, opt, i, parentS )
local ops,points;
if i = 0 then
points := [1..Minimum(20,LargestMovedPoint(grp))];
else
points := [1..LargestMovedPoint(grp)];
fi;
ops := List([1..Length(points)],x->OnPoints);
return rec( points := points, ops := ops, used := 0 );
end );
GENSS_HACK := OnPoints;
InstallGlobalFunction( GENSS_OpFunctionMaker, function(op,index)
local name,s,f;
name := NAME_FUNC(op);
if name = "unknown" or name = "GENSS_HACK" then return fail; fi;
s := Concatenation( "GENSS_HACK := function(x,el) return ",
name, "(x,el[",String(index),"]); end;" );
f := InputTextString(s);
Read(f);
return GENSS_HACK;
end);
InstallMethod( FindBasePointCandidates, "for a direct product",
[ IsGroup, IsRecord, IsInt, IsObject ],
function( grp, opt, i, parentS )
local gens,l,cand,j,factgens,fac,cand2,k,op,S2,op2;
gens := GeneratorsOfGroup(grp);
if not(ForAll(gens,IsDirectProductElement)) or Length(gens) = 0 then
TryNextMethod();
fi;
l := Length(gens[1]);
cand := rec( points := [], ops := [] );
for j in [1..l] do
factgens := List(gens,x->x[j]);
fac := Group(factgens);
S2 := StabilizerChain(fac);
cand2 := BaseStabilizerChain(S2);
for k in [1..Length(cand2.ops)] do
op := cand2.ops[k];
op2 := GENSS_OpFunctionMaker(op,j);
if op2 <> fail then
Add(cand.points,cand2.points[k]);
Add(cand.ops,op2);
fi;
od;
od;
cand.used := 0;
return cand;
end );
InstallGlobalFunction( GENSS_NextBasePoint,
function( gens, cand, opt, S )
local NotFixedUnderAllGens,i,notfixed;
if IsBound(opt.FindBasePointCandidatesData) then
Unbind(opt.FindBasePointCandidatesData);
fi;
NotFixedUnderAllGens := function( gens, x, op )
if IsObjWithMemory(gens[1]) then
return ForAny( gens, g->op(x,g!.el) <> x );
else
return ForAny( gens, g->op(x,g) <> x );
fi;
end;
# S can be false or a stabilizer chain record
if S <> false and not(opt.StrictlyUseCandidates) then
# try points in previous orbit
for i in [2..Minimum(Length(S!.orb),opt.NumberPrevOrbitPoints)] do
if NotFixedUnderAllGens(gens,S!.orb[i],S!.orb!.op) then
Info(InfoGenSS,3,"Taking another point in same orbit.");
return rec( point := S!.orb[i], op := S!.orb!.op,
cand := cand );
fi;
od;
# Maybe we can take the last base point now acting non-projectively?
if opt.Projective = false and IsIdenticalObj(S!.orb!.op,OnLines) and
NotFixedUnderAllGens(gens,S!.orb[1],OnRight) then
Info(InfoGenSS,3,"Taking same point in non-projective action.");
return rec( point := S!.orb[1], op := OnRight, cand := cand );
fi;
fi;
# Now use the candidates we already have:
repeat
if cand.used >= Length(cand.points) then
if IsBound(cand.points2) and Length(cand.points2) > 0 then
cand := rec( points := cand.points2, ops := cand.ops2,
used := 0 );
else
cand := FindBasePointCandidates(Group(gens),opt,1,S);
opt.StrictlyUseCandidates := false;
opt.Reduced := true;
fi;
fi;
cand.used := cand.used + 1;
notfixed := NotFixedUnderAllGens(gens,
cand.points[cand.used],cand.ops[cand.used]);
if notfixed = false and opt.Reduced = true then # everything is fixed!
if IsBound(cand.points2) then
# Maybe our current stabilizer is too small, we still
# might want to use these points again:
Add(cand.points2,cand.points[cand.used]);
Add(cand.ops2,cand.ops[cand.used]);
fi;
fi;
until opt.Reduced = false or notfixed;
return rec( point := cand.points[cand.used], op := cand.ops[cand.used],
cand := cand );
end );
InstallGlobalFunction( GENSS_CreateStabChainRecord,
function( parentS, gens, size, nextpoint, nextop, cand, opt )
# parentS can be false or the parent in the stabiliser chain
local base, layer, stronggens, layergens, nr, orb, S, hashsize;
Info( InfoGenSS, 4, "Creating new stab chain record..." );
if parentS = false then
base := [];
layer := 1;
stronggens := ShallowCopy(gens);
layergens := [1..Length(gens)]; # indices in stronggens
else
base := parentS!.base;
layer := parentS!.layer + 1;
stronggens := parentS!.stronggens;
nr := Length(stronggens);
Append(stronggens,gens);
layergens := [nr+1..Length(stronggens)];
fi;
gens := ShallowCopy(gens);
# Note that we do ShallowCopy such that the original list and the
# one in the orbit record are different from each other.
if IsInt(size) then
hashsize := NextPrimeInt(Minimum(size,opt.InitialHashSize));
else
hashsize := opt.InitialHashSize;
fi;
orb := Orb( gens, nextpoint, nextop,
rec( treehashsize := hashsize, schreier := true,
log := opt.OrbitsWithLog,
report := opt.Report ) );
S := rec( stab := false, orb := orb, cand := cand, base := base,
opt := opt, layer := layer, parentS := parentS,
stronggens := stronggens, layergens := layergens,
size := size, randpool := [], IsOne := opt.IsOne );
if parentS = false then
S!.topS := S;
else
S!.topS := parentS!.topS;
fi;
if IsBound(opt.FindBasePointCandidatesData) then
S.randpool := opt.FindBasePointCandidatesData.randpool;
fi;
Add(base,nextpoint);
S.orb!.stabilizerchain := S;
Objectify( StabChainByOrbType, S );
return S;
end );
InstallGlobalFunction( GENSS_RandomElementFromAbove,
function( S, i )
# This function provides a random element for the i-th stabiliser
# in the chain "from above". These elements eventually come from
# the one product replacer object in the opt record for the whole
# group. They are sifted through i layers of the stabiliser chain
# infrastructure until they stabilise the first i points (i=0 simply
# means a random element in the whole group). If we find out underways
# that an orbit is too small (that is, a stabiliser was not complete),
# we fix the stabiliser chain as we go. Random elements that have
# been generated in a layer j < i previously and are not yet used
# as strong generators can be taken from a pool that was kept in
# the stabiliser chain. S must either be layer i of the stabiliser
# or (for i=0) it can be equal to the options record opt chain.
local SS, x, topS, j, o, p, po;
if i = 0 then
return Next(S.pr); # in this case we got the options record
fi;
if S!.layer <> i then
Error("i must be equal to the layer");
fi;
SS := S;
while true do # will be left by break
if Length(SS!.randpool) > 0 then
x := Remove(SS!.randpool,Length(SS!.randpool));
topS := SS!.topS;
break;
fi;
if SS!.parentS = false then # we have reached the top
if IsBound(SS!.opt.RandomElmFunc) then
x := SS!.opt.RandomElmFunc();
else
x := Next(SS!.opt.pr);
fi;
topS := SS; # remember top
break;
fi;
SS := SS!.parentS;
od;
# We now have a random element x in layer SS!.layer, that is,
# it is contained in the SS!.layer-1-th stabiliser, sift it down:
j := SS!.layer-1;
while j < i do
o := SS!.orb;
if IsObjWithMemory(x) then
p := o!.op(o[1],x!.el);
else
p := o!.op(o[1],x);
fi;
po := Position(o,p);
if po = fail then # not in current stabilizer
Info(InfoGenSS,3,"Random element from top found error in layer ",
SS!.layer);
AddGeneratorToStabilizerChain(topS,x);
# We add it at the top to add it in every orbit above except
# the first one!
return GENSS_RandomElementFromAbove(S,i);
fi;
# Now sift through Schreier tree:
while po > 1 do
x := x * SS!.orb!.gensi[o!.schreiergen[po]];
po := o!.schreierpos[po];
od;
SS := SS!.stab;
j := j + 1;
od;
# After this, we have successfully reached i-th stabilizer
return x;
end );
InstallGlobalFunction( GENSS_StabilizerChainInner,
function( gens, size, layer, cand, opt, parentS )
# Computes a stabilizer chain for the group generated by gens
# with known size size (can be false if not known). This will be
# layer layer in the final stabilizer chain. cand is a (shared)
# record for base point candidates and opt the (shared) option
# record. This is called in StabilizerChain and calls itself.
# It also can be called if a new layer is needed.
local base,gen,S,i,merk,merk2,next,pr,r,stabgens,x;
Info(InfoGenSS,4,"Entering GENSS_StabilizerChainInner layer=",layer);
next := GENSS_NextBasePoint(gens,cand,opt,parentS);
cand := next.cand; # This could have changed
S := GENSS_CreateStabChainRecord(parentS,gens,size,
next.point,next.op,next.cand,opt);
base := S!.base;
Info( InfoGenSS, 3, "Entering orbit enumeration layer ",layer,"..." );
repeat
Enumerate(S!.orb,opt.OrbitLengthLimit);
if not(IsClosed(S!.orb)) then
if opt.FailInsteadOfError then
return "Orbit too long, increase opt.OrbitLengthLimit";
else
Error("Orbit too long, increase opt.OrbitLengthLimit");
fi;
fi;
until IsClosed(S!.orb);
Info(InfoGenSS, 2, "Layer ", layer, ": Orbit length is ", Length(S!.orb));
if layer > 1 then
parentS!.stab := S; # such that from now on random element
# generation works!
else
if (Length(S!.orb) > 50 or S!.orb!.depth > 5) and
S!.opt.OrbitsWithLog then
Info(InfoGenSS, 3, "Trying to make Schreier tree shallower (depth=",
S!.orb!.depth,")...");
merk := Length(S!.orb!.gens);
merk2 := Length(S!.stronggens);
MakeSchreierTreeShallow(S!.orb);
Append(S!.stronggens,S!.orb!.gens{[merk+1..Length(S!.orb!.gens)]});
Append(S!.layergens,[merk2+1..Length(S!.stronggens)]);
Info(InfoGenSS, 3, "Depth is now ",S!.orb!.depth);
fi;
fi;
S!.orb!.gensi := List(S!.orb!.gens,x->x^-1);
# Are we done?
if size <> false and Length(S!.orb) = size then
S!.proof := true;
Info(InfoGenSS,4,"Leaving GENSS_StabilizerChainInner layer=",layer);
return S;
fi;
# Now create a few random stabilizer elements:
stabgens := EmptyPlist(opt.RandomStabGens);
for i in [1..opt.RandomStabGens] do
x := GENSS_RandomElementFromAbove(S,layer);
if not(S!.IsOne(x)) then
Add(stabgens,x);
fi;
od;
Info(InfoGenSS,3,"Created ",opt.RandomStabGens,
" random stab els, ",
Length(stabgens)," non-trivial.");
if Length(stabgens) > 0 then # there is a non-trivial stabiliser
Info(InfoGenSS,3,"Found ",Length(stabgens)," non-trivial ones.");
if size <> false then
S!.stab := GENSS_StabilizerChainInner(stabgens,size/Length(S!.orb),
layer+1,cand,opt,S);
else
S!.stab := GENSS_StabilizerChainInner(stabgens,false,
layer+1,cand,opt,S);
fi;
if IsString(S!.stab) then return S!.stab; fi;
if opt.ImmediateVerificationElements > 0 then
Info(InfoGenSS,2,"Doing immediate verification in layer ",
S!.layer," (",opt.ImmediateVerificationElements,
" elements)...");
i := 0;
#Error("foobar");
while i < opt.ImmediateVerificationElements do
i := i + 1;
x := GENSS_RandomElementFromAbove(S,layer);
if AddGeneratorToStabilizerChain(S!.topS,x) then
Info( InfoGenSS, 2, "Immediate verification found error ",
"(layer ",S!.layer,")..." );
i := 0;
fi;
od;
fi;
S!.proof := S!.stab!.proof; # hand up information
else
# We are not sure that the next stabiliser is trivial, but we believe!
Info(InfoGenSS,3,"Found no non-trivial ones.");
S!.proof := false;
fi;
Info(InfoGenSS,4,"Leaving GENSS_StabilizerChainInner layer=",layer);
return S;
end );
InstallGlobalFunction( GENSS_DeriveCandidatesFromStabChain,
function( S )
local cand;
cand := rec( points := [], ops := [], used := 0 );
while S <> false do
Add( cand.points, S!.orb[1] );
Add( cand.ops, S!.orb!.op );
S := S!.stab;
od;
return cand;
end );
InstallGlobalFunction( GENSS_TrivialOp,
function( p, x )
return p;
end );
InstallMethod( StabilizerChain, "for a group object", [ IsGroup ],
function( grp )
return StabilizerChain( grp, rec() );
end );
InstallMethod(VerifyStabilizerChainMC,
"for a stabilizer chain and a positive integer",
[ IsStabilizerChainByOrb, IsInt ],
function( S, nrels )
local i,x;
Info(InfoGenSS,2,"Doing randomized verification...");
i := 0;
while i < nrels do
i := i + 1;
if IsBound(S!.opt.RandomElmFunc) then
x := S!.opt.RandomElmFunc();
else
x := Next(S!.opt.pr);
fi;
if AddGeneratorToStabilizerChain(S,x) then
Info( InfoGenSS, 2, "Verification found error ... ",
"new size ", Size(S) );
i := 0;
fi;
od;
end );
InstallMethod( StabilizerChain, "for a group object and a record",
[ IsGroup, IsRecord ],
function( grp, opt )
# Computes a stabilizer chain for the group grp
# Possible options:
# random: compatibility to StabChain, works as there
# ErrorBound: rational giving the error bound, probability of an
# error will be proven to be lower than this
# if given, this sets VerifyElements and
# DeterministicVerification
# VerifyElements: number of random elements for verification
# DeterministicVerification: flag, whether to do or not
# (not yet implemented)
# Base: specify a known base, this can either be the
# StabilizerChain of a known supergroup or a list of
# points, if "BaseOps" is not given, it is either
# set to OnPoints if Projective is not set or false,
# or is set to OnLines if Projective is set to true,
# if Base is set
# BaseOps: a list of operation functions for the base points
# Projective: if set to true indicates that the matrices given
# as generators are to be thought as projective
# elements, that is, the whole StabilizerChain will
# be one for a projective group!
# StrictlyUseCandidates: if set to true exactly the points in
# opt.cand will be used as base points and no other
# points from previous orbits are used first
# this is set when Base was specified
# Reduced: if set to true no orbits of length 1 will occur
# (except for the trivial group), is true by default
# Cand: initializer for cand, which are base points
# candidates
# must be a record with components "points", "ops",
# "used", the latter indicates the largest index
# that has already been used. "used" is automatically
# set to 0 if not set.
# TryShortOrbit: Number of tries for the short orbit finding alg.
# StabGenScramble,
# StabGenScrambleFactor,
# StabGenAddSlots,
# StabGenMaxDepth: parameters for product replacer for generating
# random elements
# VeryShortOrbLimit: when looking for short orbits try a random
# vector and enumerate its orbit until this limit
#
# ... to be continued
#
local S,cand,i,pr,prob,x,gens,SS;
# First a few preparations, then we delegate to GENSS_StabilizerChainInner:
# Add some default options:
if (HasSize(grp) and not(IsBound(opt.Projective) and opt.Projective))
or IsBound(opt.Size) then
if not(IsBound(opt.ImmediateVerificationElements)) then
opt.ImmediateVerificationElements := 0;
fi;
fi;
if IsBound(opt.Projective) and opt.Projective then
opt.IsOne := GENSS_IsOneProjective;
elif not(IsBound(opt.IsOne)) then
opt.IsOne := IsOne;
fi;
# Now opt.IsOne is set to a function to check whether or not a group
# element is equal to the identity.
GENSS_CopyDefaultOptions(GENSS,opt);
# Check for the identity group:
gens := GeneratorsOfGroup(grp);
if Length(gens) = 0 or
ForAll(gens,opt.IsOne) then
# Set up a trivial stabilizer chain record:
S := GENSS_CreateStabChainRecord(false,gens,1,1,GENSS_TrivialOp,
rec( points := [], ops := [],
used := 0),opt);
Enumerate(S!.orb);
S!.orb!.gensi := List(S!.orb!.gens,x->x^-1);
S!.proof := true;
S!.trivialgroup := true;
if (not(IsBound(opt.Projective)) or
opt.Projective = false) and
IsIdenticalObj(opt.IsOne,IsOne) then
SetStoredStabilizerChain(grp,S);
fi;
return S;
fi;
# Setup a random element generator for the whole group and store
# it in the opt record:
opt.pr := ProductReplacer( gens,
rec( scramble := opt.StabGenScramble,
scramblefactor := opt.StabGenScrambleFactor,
addslots := opt.StabGenAddSlots,
maxdepth := opt.StabGenMaxDepth ));
# Old style error probability for compatibility:
if IsBound(opt.random) then
if opt.random = 0 then
opt.VerifyElements := 0;
elif opt.random = 1000 then
opt.DeterministicVerification := true;
Info(InfoGenSS,1,"Warning: Deterministic verification not yet ",
"implemented!");
else
prob := 1/2;
opt.VerifyElements := 1;
while prob * (1000-opt.random) >= 1 do
prob := prob / 2;
opt.VerifyElements := opt.VerifyElements + 1;
od;
fi;
fi;
# The new style error bounds:
if IsBound(opt.ErrorBound) then
prob := 1/2;
opt.VerifyElements := 1;
while prob >= opt.ErrorBound do
prob := prob / 2;
opt.VerifyElements := opt.VerifyElements + 1;
od;
fi;
# Find base point candidates:
if IsBound(opt.Base) then
if IsStabilizerChain(opt.Base) then
cand := GENSS_DeriveCandidatesFromStabChain(opt.Base);
else # directly take the base points:
cand := rec( points := opt.Base, used := 0, points2 := [],
ops2 := [] );
if not(IsBound(opt.BaseOps)) then
# Let's guess the ops:
if IsBound(opt.Projective) and opt.Projective = true then
cand.ops := ListWithIdenticalEntries(Length(opt.Base),
OnLines);
else
cand.ops := ListWithIdenticalEntries(Length(opt.Base),
OnPoints);
fi;
else
cand.ops := opt.BaseOps;
fi;
fi;
opt.StrictlyUseCandidates := true;
elif IsBound(opt.Cand) then
cand := opt.Cand;
if not(IsBound(cand.used)) then
cand.used := 0;
fi;
cand.points2 := [];
cand.ops2 := [];
else
# Otherwise try different things later using generic methods:
cand := rec( points := [], ops := [], used := 0 );
fi;
if not(IsBound(opt.StrictlyUseCandidates)) then
opt.StrictlyUseCandidates := false;
fi;
if HasSize(grp) and not(opt.Projective) then
S := GENSS_StabilizerChainInner(GeneratorsOfGroup(grp), Size(grp), 1,
cand, opt, false);
elif IsBound(opt.Size) then
S := GENSS_StabilizerChainInner(GeneratorsOfGroup(grp), opt.Size, 1,
cand, opt, false);
else
S := GENSS_StabilizerChainInner(GeneratorsOfGroup(grp), false, 1,
cand, opt, false);
fi;
if IsString(S) then return S; fi;
# Do we already have a proof?
if S!.proof then
if not(IsBound(opt.Projective)) or opt.Projective = false then
SetStoredStabilizerChain(grp,S);
fi;
return S;
fi;
if opt.VerifyElements = 0 then return S; fi;
Info(InfoGenSS,2,"Current size found: ",Size(S));
# Now a possible verification phase:
if S!.size <> false then # we knew the size in advance
Info(InfoGenSS,2,"Doing verification via known size...");
while Size(S) < S!.size do
Info(InfoGenSS,2,"Known size not reached, throwing in a random ",
"element...");
if IsBound(opt.RandomElmFunc) then
x := opt.RandomElmFunc();
else
x := Next(opt.pr);
fi;
if AddGeneratorToStabilizerChain(S,x) then
Info( InfoGenSS, 2, "Increased size to ",Size(S) );
fi;
od;
S!.proof := true;
if (not(IsBound(opt.Projective)) or opt.Projective = false) and
IsIdenticalObj(opt.IsOne,IsOne) then
SetStoredStabilizerChain(grp,S);
fi;
else
# Do some verification here:
VerifyStabilizerChainMC(S,opt.VerifyElements);
fi;
# Now clean up the random element pools to save memory:
SS := S;
while SS <> false do
if IsBound(SS!.randpool) and Length(SS!.randpool) > 0 then
SS!.randpool := [];
fi;
SS := SS!.stab;
od;
return S;
end );
InstallMethod( AddGeneratorToStabilizerChain,
"for a stabilizer chain and a new generator",
[ IsStabilizerChain and IsStabilizerChainByOrb, IsObject ],
function( S, el )
# Increases the set represented by S by the generator el.
local SS, r, n, pr, i, newstrongnr;
if IsBound(S!.trivialgroup) and S!.trivialgroup then
if S!.IsOne(el) then
return false;
fi;
SS := StabilizerChain(Group(el),S!.opt);
if IsString(SS) then return SS; fi;
for n in NamesOfComponents(SS) do
S!.(n) := SS!.(n);
od;
Unbind(S!.trivialgroup);
return true;
fi;
r := SiftGroupElement( S, el );
# if this is one then el is already contained in the stabilizer chain
if r.isone then # already in the group!
return false;
fi;
# Now there remain two cases:
# (1) the sift stopped somewhere and we have to add a generator there
# (2) the sift ran all through the chain and the element still was not
# the identity, then we have to prolong the chain
if r.S <> false then # case (1)
SS := r.S;
Info( InfoGenSS, 2, "Adding new generator to stabilizer chain ",
"in layer ", SS!.layer, "..." );
Add(SS!.stronggens,r.rem);
Add(SS!.layergens,Length(SS!.stronggens));
AddGeneratorsToOrbit(SS!.orb,[r.rem]);
Add(SS!.orb!.gensi,r.rem^-1);
newstrongnr := Length(SS!.stronggens);
Info( InfoGenSS, 4, "Entering orbit enumeration layer ",SS!.layer,
"..." );
repeat
Enumerate(SS!.orb,S!.opt.OrbitLengthLimit);
if not(IsClosed(SS!.orb)) then
if S!.opt.FailInsteadOfError then
return "Orbit too long, increase S!.opt.OrbitLengthLimit";
else
Error("Orbit too long, increase S!.opt.OrbitLengthLimit!");
fi;
fi;
until IsClosed(SS!.orb);
Info( InfoGenSS, 4, "Done orbit enumeration layer ",SS!.layer,"..." );
SS!.proof := false;
else # case (2)
# Note that we do not create a pr instance here for one
# generator, this will be done later on as needed...
SS := r.preS;
newstrongnr := Length(SS!.stronggens)+1; # r.rem will end up there !
SS!.stab := GENSS_StabilizerChainInner([r.rem],false,
SS!.layer+1,SS!.cand, SS!.opt, SS );
if IsString(SS!.stab) then return SS!.stab; fi;
SS := SS!.stab;
fi;
# Now we have added a new generator (or a new layer) at layer SS,
# the new gen came from layer S (we were called here, after all),
# thus we have to check, whether all the orbits between S (inclusively)
# and SS (exclusively) are also closed under the new generator r.rem,
# we add it to all these orbits, thereby also making the Schreier trees
# shallower:
while S!.layer < SS!.layer do
Info(InfoGenSS,2,"Adding new generator to orbit in layer ",S!.layer);
Add(S!.layergens,newstrongnr);
AddGeneratorsToOrbit(S!.orb,[r.rem]);
Add(S!.orb!.gensi,r.rem^-1);
S := S!.stab;
od;
# Finally, we have to add it to the product replacer!
AddGeneratorToProductReplacer(S!.opt!.pr,r.rem);
return true;
end );
InstallMethod( SiftGroupElement, "for a stabilizer chain and a group element",
[ IsStabilizerChain and IsStabilizerChainByOrb, IsObject ],
function( S, x )
local o,p,po,preS,r,isone;
isone := S!.IsOne;
preS := false;
while S <> false do
o := S!.orb;
if IsObjWithMemory(x) then
p := o!.op(o[1],x!.el);
else
p := o!.op(o[1],x);
fi;
po := Position(o,p);
if po = fail then # not in current stabilizer
return rec( isone := false, rem := x, S := S, preS := preS );
fi;
# Now sift through Schreier tree:
while po > 1 do
x := x * S!.orb!.gensi[o!.schreiergen[po]];
po := o!.schreierpos[po];
od;
preS := S;
S := S!.stab;
od;
r := rec( rem := x, S := false, preS := preS, isone := isone(x) );
return r;
end );
InstallMethod( SiftGroupElementSLP,
"for a stabilizer chain and a group element",
[ IsStabilizerChain and IsStabilizerChainByOrb, IsObject ],
function( S, x )
local preS, nrstrong, slp, o, p, po, r, isone;
preS := false;
isone := S!.IsOne;
nrstrong := Length(S!.stronggens);
slp := []; # will be reversed in the end
while S <> false do
o := S!.orb;
p := o!.op(o[1],x);
po := Position(o,p);
if po = fail then # not in current stabilizer
return rec( isone := false, rem := x, S := S, preS := preS,
slp := fail );
fi;
# Now sift through Schreier tree:
while po > 1 do
x := x * S!.orb!.gensi[o!.schreiergen[po]];
Add(slp,1);
Add(slp,S!.layergens[o!.schreiergen[po]]);
po := o!.schreierpos[po];
od;
preS := S;
S := S!.stab;
od;
r := rec( rem := x, S := false, preS := preS, isone := isone(x) );
if r.isone then
if Length(slp) = 0 then # element was the identity!
r.slp := StraightLineProgramNC([[1,0]],nrstrong);
else
r.slp := StraightLineProgramNC(
[slp{[Length(slp),Length(slp)-1..1]}],nrstrong);
fi;
else
r.slp := fail;
fi;
return r;
end );
InstallMethod( StrongGenerators, "for a stabilizer chain",
[ IsStabilizerChain and IsStabilizerChainByOrb ],
function( S )
return S!.stronggens;
end );
InstallMethod( NrStrongGenerators, "for a stabilizer chain",
[ IsStabilizerChain and IsStabilizerChainByOrb ],
function( S )
return Length(S!.stronggens);
end );
InstallMethod( ForgetMemory, "for a stabilizer chain",
[ IsStabilizerChain and IsStabilizerChainByOrb ],
function( S )
ForgetMemory(S!.stronggens);
while S <> false do
ForgetMemory(S!.orb);
ForgetMemory(S!.orb!.gensi);
S := S!.stab;
od;
end );
InstallMethod( AddNormalizingGenToLayer,
"for a stabilizer chain, a group element, and a prime number",
[ IsStabilizerChain and IsStabilizerChainByOrb, IsObject, IsPosInt ],
function( S, x, p )
# This is a very specialised function which is not officially documented.
# It is used for the matrix PCGS code where a stabiliser chain is
# built up using a known base in a generator by generator fashion
# where the generators are the generators of the PCGS. Thus a lot
# is known about them.
# This function assumes that x is a new generator normalising the
# group generated by the previous generators in the layer described
# by S. So in particular, S is not necessarily the top of the
# stabiliser chain but the layer where actually something happens.
# Namely, since it is known that the new generator generates
# a group which is p times as big as the previous one (p a prime)
# it is known that the orbit in the given layer will become
# p times as big and we can enumerate it relatively cheaply.
local lmp,o,oldnrgens;
o := S!.orb;
oldnrgens := Length(o!.gens);
Add(o!.gens,x);
Add(o!.gensi,x^-1);
if IsPermOnIntOrbitRep(o) then
lmp := LargestMovedPoint(o!.gens);
if lmp > Length(o!.tab) then
Append(o!.tab,ListWithIdenticalEntries(lmp-Length(o!.tab),0));
fi;
fi;
ResetFilterObj(o,IsClosed);
o!.pos := 1;
o!.genstoapply := [oldnrgens+1..Length(o!.gens)];
repeat
Enumerate(o,S!.opt.OrbitLengthLimit);
if not(IsClosed(o)) then
if S!.opt.FailInsteadOfError then
return "Orbit too long, increase S!.opt.OrbitLengthLimit";
else
Error("Orbit too long, increase S!.opt.OrbitLengthLimit!");
fi;
fi;
until IsClosed(S!.orb);
o!.genstoapply := [1..Length(o!.gens)];
# Now fix up the stabilizer chain record:
if not x in S!.stronggens then
Add(S!.stronggens,x);
Add(S!.layergens,Length(S!.stronggens));
fi;
Unbind(S!.size); # whatever we knew before, we know no longer
end );
InstallMethod( GENSS_CreateSchreierGenerator,
"for a stabilizer chain, a position in the first orbit and a gen number",
[ IsStabilizerChain and IsStabilizerChainByOrb, IsPosInt, IsPosInt ],
function( S, i, j )
local o,p,w1,w2;
o := S!.orb;
p := Position(o,o!.op(o[i],o!.gens[j]));
if o!.schreierpos[p] = i and o!.schreiergen[p] = j then
return fail;
fi;
w1 := TraceSchreierTreeForward(o,i);
w2 := TraceSchreierTreeBack(o,p);
return [w1,j,w2];
end );
InstallGlobalFunction( GENSS_Prod,
function( gens, word )
if Length(word) = 0 then
return gens[1]^0;
else
return Product(gens{word});
fi;
end );
InstallGlobalFunction( VerifyStabilizerChainTC,
function( S ) Error("Currently not functional!"); return fail; end );
## function( S )
## local Grels,Hrels,Prels,MakeSchreierGens,ct,f,gens,gensi,i,j,k,l,li,max,
## newpres,nrgens,nrschr,o,ords,pres,sb,sgs,slp,subgens,v,w,x,
## cosetnrlimitfactor;
##
## allstrong := function(S)
## st := Set(S!.layergens);
## while S!.stab <> false do
## S := S!.stab;
## UniteSet(st,S!.layergens);
## od;
## return st;
## end;
## if S!.stab <> false then
## pres := VerifyStabilizerChainTC(S!.stab);
## if IsList(pres) then return pres; fi;
## strongbelow := allstrong(S!.stab);
## else
## pres := StraightLineProgram([[]],0);
## strongbelow := [];
## fi;
## strong := allstrong(S);
## Info(InfoGenSS,1,"Verifying stabilizer chain in layer ",S!.layer);
##
## # First create a few Schreier generators:
## sgs := [];
## i := 1;
## j := 1;
## o := S!.orb;
## nrgens := Length(o!.gens);
## MakeSchreierGens := function(n)
## local sg;
## Info(InfoGenSS,3,"Creating ",n," Schreier generators...");
## while Length(sgs) < n and
## i <= Length(o) do
## sg := GENSS_CreateSchreierGenerator(S,i,j);
## j := j + 1;
## if j > nrgens then
## j := 1;
## i := i + 1;
## fi;
## if sg <> fail then
## Add(sgs,sg);
## fi;
## od;
## end;
##
## nrschr := S!.opt.NumberSchreierGens;
## MakeSchreierGens(nrschr);
## f := FreeGroup(Length(strong));
## gens := GeneratorsOfGroup(f);
## gensi := List(gens,x->x^-1);
## subgens := gens{List(strongbelow,i->Position(strong,i))};
## Hrels := ResultOfStraightLineProgram(pres,subgens);
## if S!.opt.Projective then
## ords := List([1..nrgens],i->ProjectiveOrder(o!.gens[i]));
## else
## ords := List([1..nrgens],i->Order(o!.gens[i]));
## fi;
## Prels := List([1..nrgens],i->gens[i+sb]^ords[i]);
## Grels := [];
## cosetnrlimitfactor := 4;
## while true do # will be left by return eventually
## for k in [Length(Grels)+1..Length(sgs)] do
## Grels[k] := GENSS_Prod(gens,sgs[k][1]+sb) * gens[sgs[k][2]+sb] *
## GENSS_Prod(gensi,sgs[k][3]+sb);
## x := GENSS_Prod(o!.gens,sgs[k][1]) * o!.gens[sgs[k][2]] *
## GENSS_Prod(o!.gensi,sgs[k][3]);
## if S!.stab <> false then
## slp := SiftGroupElementSLP(S!.stab,x);
## if not(slp.isone) then
## return [fail,S!.layer];
## fi;
## Grels[k] := Grels[k] / ResultOfStraightLineProgram(slp.slp,
## subgens);
## sgs[k][4] := slp.slp;
## else
## if not(S!.IsOne(x)) then
## return [fail,S!.layer];
## fi;
## sgs[k][4] := false;
## fi;
## od;
## Info(InfoGenSS,2,"Doing coset enumeration with limit ",
## cosetnrlimitfactor*Length(o));
## ct := CosetTableFromGensAndRels(gens,Concatenation(Prels,Hrels,Grels),
## subgens:max := cosetnrlimitfactor*Length(o),silent);
## if ct = fail then # did not close!
## cosetnrlimitfactor := QuoInt(cosetnrlimitfactor*3,2);
## Info(InfoGenSS,2,"Coset enumeration did not finish!");
## #Error(1);
## if nrschr > Length(sgs) # or
## # nrschr > S!.opt.MaxNumberSchreierGens
## then # we are done!
## # Something is wrong!
## return [fail, S!.layer];
## fi;
## else
## Info(InfoGenSS,2,"Coset enumeration found ",Length(ct[1]),
## " cosets.");
## #Error(2);
## if Length(ct[1]) = Length(o) then
## # Verification is OK, now build a presentation:
## l := GeneratorsWithMemory(
## ListWithIdenticalEntries(S!.nrstrong,()));
## li := List(l,x->x^-1);
## newpres := ResultOfStraightLineProgram(pres,
## l{[1..S!.strongbelow]});
## for k in [1..nrgens] do
## Add(newpres,l[k+sb]^ords[k]);
## od;
## for k in [1..Length(sgs)] do
## if sgs[k][4] <> false then
## Add(newpres,
## GENSS_Prod(l,sgs[k][1]+sb)*l[sgs[k][2]+sb]*
## GENSS_Prod(li,sgs[k][3]+sb)*
## ResultOfStraightLineProgram(sgs[k][4],l)^-1);
## else
## Add(newpres,GENSS_Prod(l,sgs[k][1]+sb)*l[sgs[k][2]+sb]*
## GENSS_Prod(li,sgs[k][3]+sb));
## fi;
## od;
## Info(InfoGenSS,2,"Found presentation for layer ",S!.layer,
## " using ",Length(newpres)," relators.");
## return SLPOfElms(newpres);
## fi;
## fi;
## nrschr := QuoInt(nrschr*4,2);
## MakeSchreierGens(nrschr);
## od;
## end);
GENSS_CosetTableFromGensAndRelsInit :=
function(fgens,grels,fsgens,limit)
local next, prev, # next and previous coset on lists
firstFree, lastFree, # first and last free coset
firstDef, lastDef, # first and last defined coset
table, # columns in the table for gens
rels, # representatives of the relators
relsGen, # relators sorted by start generator
subgroup, # rows for the subgroup gens
i, gen, inv, # loop variables for generator
g, # loop variable for generator col
rel, # loop variables for relation
p, p1, p2, # generator position numbers
app, # arguments list for 'MakeConsequences'
j, # integer variable
length, length2, # length of relator (times 2)
cols,
nums,
l,
nrdef, # number of defined cosets
nrmax, # maximal value of the above
nrdel, # number of deleted cosets
nrinf; # number for next information message
# give some information
Info( InfoFpGroup, 3, " defined deleted alive maximal");
nrdef := 1;
nrmax := 1;
nrdel := 0;
nrinf := 1000;
# define one coset (1)
firstDef := 1; lastDef := 1;
firstFree := 2; lastFree := limit;
# make the lists that link together all the cosets
next := [ 2 .. limit + 1 ]; next[1] := 0; next[limit] := 0;
prev := [ 0 .. limit - 1 ]; prev[2] := 0;
# compute the representatives for the relators
rels := RelatorRepresentatives( grels );
# make the columns for the generators
table := [];
for gen in fgens do
g := ListWithIdenticalEntries( limit, 0 );
Add( table, g );
if not ( gen^2 in rels or gen^-2 in rels ) then
g := ListWithIdenticalEntries( limit, 0 );
fi;
Add( table, g );
od;
# make the rows for the relators and distribute over relsGen
relsGen := RelsSortedByStartGen( fgens, rels, table, true );
# make the rows for the subgroup generators
subgroup := [];
for rel in fsgens do
#T this code should use ExtRepOfObj -- its faster
# cope with SLP elms
if IsStraightLineProgElm(rel) then
rel:=EvalStraightLineProgElm(rel);
fi;
length := Length( rel );
length2 := 2 * length;
nums := [ ]; nums[length2] := 0;
cols := [ ]; cols[length2] := 0;
# compute the lists.
i := 0; j := 0;
while i < length do
i := i + 1; j := j + 2;
gen := Subword( rel, i, i );
p := Position( fgens, gen );
if p = fail then
p := Position( fgens, gen^-1 );
p1 := 2 * p;
p2 := 2 * p - 1;
else
p1 := 2 * p - 1;
p2 := 2 * p;
fi;
nums[j] := p1; cols[j] := table[p1];
nums[j-1] := p2; cols[j-1] := table[p2];
od;
Add( subgroup, [ nums, cols ] );
od;
# make the structure that is passed to 'MakeConsequences'
app := [ table, next, prev, relsGen, subgroup,
firstFree, lastFree, firstDef, lastDef, 0, 0 ];
# we do not want minimal gaps to be marked in the coset table
app[12] := 0;
return rec( nrdef := nrdef, nrmax := nrmax, nrdel := nrdel, nrinf := nrinf,
app := app, closed := false, table := table, limit := limit,
fgens := fgens );
end;
GENSS_DoToddCoxeter := function( r )
local app,fgens,firstDef,firstFree,gen,i,inv,lastDef,lastFree,next,nrdef,
nrdel,nrinf,nrmax,prev,relsGen,subgroup,table;
fgens := r.fgens;
nrdef := r.nrdef;
nrmax := r.nrmax;
nrdel := r.nrdel;
nrinf := r.nrinf;
app := r.app;
table := app[1];
next := app[2];
prev := app[3];
relsGen := app[4];
subgroup := app[5];
firstFree := app[6];
lastFree := app[7];
firstDef := app[8];
lastDef := app[9];
# run over all the cosets:
while firstDef <> 0 do
# run through all the rows and look for undefined entries
for i in [ 1 .. Length( table ) ] do
gen := table[i];
if gen[firstDef] <= 0 then
inv := table[i + 2*(i mod 2) - 1];
# if necessary expand the table
if firstFree = 0 then
app[6] := firstFree;
app[7] := lastFree;
app[8] := firstDef;
app[9] := lastDef;
r.nrdef := nrdef;
r.nrmax := nrmax;
r.nrdel := nrdel;
r.nrinf := nrinf;
return fail;
fi;
# update the debugging information
nrdef := nrdef + 1;
if nrmax <= firstFree then
nrmax := firstFree;
fi;
# define a new coset
gen[firstDef] := firstFree;
inv[firstFree] := firstDef;
next[lastDef] := firstFree;
prev[firstFree] := lastDef;
lastDef := firstFree;
firstFree := next[firstFree];
next[lastDef] := 0;
# set up the deduction queue and run over it until it's empty
app[6] := firstFree;
app[7] := lastFree;
app[8] := firstDef;
app[9] := lastDef;
app[10] := i;
app[11] := firstDef;
nrdel := nrdel + MakeConsequences( app );
firstFree := app[6];
lastFree := app[7];
firstDef := app[8];
lastDef := app[9];
# give some information
if nrinf <= nrdef+nrdel then
Info( InfoFpGroup, 3, "\t", nrdef, "\t", nrinf-nrdef,
"\t", 2*nrdef-nrinf, "\t", nrmax );
nrinf := ( Int(nrdef+nrdel)/1000 + 1 ) * 1000;
fi;
fi;
od;
firstDef := next[firstDef];
od;
Info( InfoFpGroup, 2, "\t", nrdef, "\t", nrdel, "\t", nrdef-nrdel, "\t",
nrmax );
# separate pairs of identical table columns.
for i in [ 1 .. Length( fgens ) ] do
if IsIdenticalObj( table[2*i-1], table[2*i] ) then
table[2*i] := StructuralCopy( table[2*i-1] );
fi;
od;
# standardize the table
StandardizeTable( table );
# return the success result
r.closed := true;
return true;
end;
GENSS_ToddCoxeterExpandLimit := function(r,newlimit)
local app,g,l,limit,next,prev,table;
if newlimit <= r.limit then return r; fi;
app := r.app;
table := app[1];
next := app[2];
prev := app[3];
limit := r.limit;
next[newlimit] := 0;
prev[newlimit] := newlimit-1;
for g in table do g[newlimit] := 0; od;
for l in [ limit+2 .. newlimit-1 ] do
next[l] := l+1;
prev[l] := l-1;
for g in table do g[l] := 0; od;
od;
next[limit+1] := limit+2;
prev[limit+1] := 0;
for g in table do g[limit+1] := 0; od;
app[6] := limit+1; # this is firstFree
r.limit := newlimit;
app[7] := newlimit; # this is lastFree
end;
GENSS_TCAddRelators := function(r,newrels)
local i,relsGen;
newrels := RelatorRepresentatives(newrels);
newrels := RelsSortedByStartGen( r.fgens, newrels, r.table, true );
relsGen := r.app[4];
for i in [1..Length(relsGen)] do
Append(relsGen[i],newrels[i]);
od;
end;
## VerifyStabilizerChainTC5 :=
## function( S )
## local FindTwoWords,Grels,Hrels,Prels,allgens,cosetlimit,done,el,f,
## gens,gensi,hom,k,l,li,newpres,newrel,nrcosets,nrgens,o,opgens,
## ords,pres,r,rels,sb,stronggens,strongi,subgens,tc,words;
##
## if S!.stab <> false then
## pres := VerifyStabilizerChainTC5(S!.stab);
## if IsList(pres) then return pres; fi;
## else
## pres := StraightLineProgram([[]],0);
## fi;
## Info(InfoGenSS,1,"Verifying stabilizer chain in layer ",S!.layer);
##
## o := S!.orb;
## nrgens := Length(o!.gens);
## sb := S!.strongbelow;
##
## f := FreeGroup(sb+nrgens);
## gens := GeneratorsOfGroup(f);
## gensi := List(gens,x->x^-1);
## allgens := 0*[1..2*Length(gens)];
## allgens{[1,3..2*Length(gens)-1]} := gens;
## allgens{[2,4..2*Length(gens)]} := gensi;
## subgens := gens{[1..S!.strongbelow]};
## Hrels := ResultOfStraightLineProgram(pres,subgens);
## if S!.opt.Projective then
## ords := List([1..nrgens],i->ProjectiveOrder(o!.gens[i]));
## else
## ords := List([1..nrgens],i->Order(o!.gens[i]));
## fi;
## Prels := List([1..nrgens],i->gens[i+sb]^ords[i]);
## Grels := [];
## stronggens := StrongGenerators(S);
## strongi := List(stronggens,x->x^-1);
## opgens := 0*[1..2*Length(stronggens)];
## opgens{[1,3..2*Length(stronggens)-1]} := stronggens;
## opgens{[2,4..2*Length(stronggens)]} := strongi;
##
## # Now start up a coset enumeration:
## cosetlimit := Maximum(QuoInt(7 * Length(o),6),Length(o)+5);
## Info(InfoGenSS,2,"Starting coset enumeration with limit ",
## cosetlimit," and ",Length(Hrels),
## "+",Length(Prels),"+",Length(Grels)," relations...");
## rels := Concatenation(Hrels,Prels);
## tc := GENSS_CosetTableFromGensAndRelsInit(gens,rels,subgens,cosetlimit);
## done := GENSS_DoToddCoxeter(tc);
##
## FindTwoWords := function(o,opgens,table)
## local TraceWord,cosets,cosetsrevtab,i,j,new,nrcosets,pt,pts,
## ptsrev,schgen,schpt,w1,w2,x,y;
## nrcosets := Length(table[1]);
## cosets := [1];
## cosetsrevtab := 0*[1..nrcosets];
## cosetsrevtab[1] := 1;
## pts := [o[1]];
## ptsrev := 0*[1..Length(o)];
## ptsrev[1] := 1;
## schpt := [fail]; # the Schreier tree
## schgen := [fail];
## i := 1;
## TraceWord := function(pos)
## local w;
## w := [];
## while pos > 1 do
## Add(w,schgen[pos]);
## pos := schpt[pos];
## od;
## return Reversed(w);
## end;
## while i <= Length(cosets) do
## for j in [1..Length(opgens)] do
## x := table[j][cosets[i]];
## if x <> 0 then # image is defined:
## if cosetsrevtab[x] = 0 then # not visited
## Add(cosets,x);
## new := Length(cosets);
## cosetsrevtab[x] := new;
## schpt[new] := i;
## schgen[new] := j;
## pt := o!.op(pts[i],opgens[j]);
## y := Position(o,pt);
## if ptsrev[y] = 0 then
## Add(pts,pt);
## ptsrev[y] := new;
## else
## # We have reached a new coset by a word that
## # maps the starting point of the orbit to the
## # same point as the one of another coset!
## w1 := TraceWord(ptsrev[y]);
## w2 := TraceWord(new);
## return [w1,w2];
## fi;
## fi;
## fi;
## od;
## i := i + 1;
## od;
## Error("Bad, this should never have been reached!");
## return fail;
## end;
##
## while true do # will be left by return eventually
## if done = true then
## nrcosets := Length(tc.table[1]);
## Info(InfoGenSS,2,"Coset enumeration found ",nrcosets," cosets.");
## if nrcosets = Length(o) then
## # Verification is OK, now build a presentation:
## l := GeneratorsWithMemory(
## ListWithIdenticalEntries(S!.nrstrong,()));
## newpres := ResultOfStraightLineProgram(pres,
## l{[1..S!.strongbelow]});
## for k in [1..nrgens] do
## Add(newpres,l[k+sb]^ords[k]);
## od;
## hom := GroupHomomorphismByImagesNC(f,Group(l),gens,l);
## for k in Grels do
## Add(newpres,ImageElm(hom,k));
## od;
## Info(InfoGenSS,2,"Found presentation for layer ",S!.layer,
## " using ",Length(newpres)," relators.");
## return SLPOfElms(newpres);
## elif nrcosets < Length(o) then
## Error("This cannot possibly have happened!");
## return [fail,S!.layer];
## else # nrcosets > Length(o)
## Info(InfoGenSS,2,"Too many cosets, we must have forgotten ",
## "another relation!");
## Error("This cannot possibly have happened2!");
## return [fail,S!.layer];
## fi;
## fi;
## # Now we have to find another relation, we do a breadth-first
## # search through the already defined cosets to find two cosets
## # that are still different but ought to be equal because the
## # corresponding orbit points are equal:
## words := FindTwoWords(o,opgens,tc.table);
## el := EvaluateWord(opgens,words[1])/EvaluateWord(opgens,words[2]);
## r := SiftGroupElementSLP(S,el);
## if not(r.isone) then
## # Error, we found a new stabilizer element!
## return [fail,S!.layer,el];
## fi;
## newrel := [(EvaluateWord(allgens,words[1])
## /EvaluateWord(allgens,words[2]))
## / ResultOfStraightLineProgram(r.slp,gens)];
## GENSS_TCAddRelators(tc,newrel);
## Add(Grels,newrel[1]);
## if Length(words[1]) > 0 then
## tc.app[10] := words[1][1];
## else
## # More difficult:
## words := ExtRepOfObj(newrel[1]);
## if words[2] > 0 then
## tc.app[10] := 2*words[1]-1;
## else
## tc.app[10] := 2*words[1];
## fi;
## fi;
## #Print("<\c");
## #for i in [1..tc.limit] do
## # for j in [1..Length(allgens)] do
## # if tc.table[j][i] <> 0 then
## # tc.app[11] := i;
## # tc.app[10] := j;
## # tc.nrdel := tc.nrdel + MakeConsequences( tc.app );
## # fi;
## # od;
## #od;
## tc.app[11] := 1;
## tc.nrdel := tc.nrdel + MakeConsequences( tc.app );
## #Print("-\c");
## done := GENSS_DoToddCoxeter(tc);
## #Print(">\c");
## if Length(Grels) mod 100 = 0 then
## #Print("\n");
## Info(InfoGenSS,2,"Currently using ",Length(Hrels),"+",
## Length(Prels),"+",Length(Grels)," relations.");
## fi;
## od;
## end;
##
## VerifyStabilizerChainMax :=
## function( S )
## local bl,count,d,gens,i,invtab,j,k,o,p,r,res,s,w1,w2,x,xi,y,isone;
##
## isone := S!.IsOne;
## if S!.stab <> false then
## res := VerifyStabilizerChainMax(S!.stab);
## if res <> true then return res; fi;
## fi;
## Info(InfoGenSS,1,"Verifying stabilizer chain in layer ",S!.layer,
## " (orbit of length ",Length(S!.orb),")");
##
## count := 0;
## gens := ShallowCopy(Set(S!.orb!.gens));
## invtab := [];
## k := Length(gens);
## for i in [1..k] do
## x := gens[i];
## xi := x^-1;
## p := Position(gens,xi);
## if p = fail then
## Add(gens,xi);
## invtab[i] := Length(gens);
## invtab[Length(gens)] := i;
## else
## invtab[i] := p;
## fi;
## od;
## o := Orb(gens,S!.orb[1],S!.orb!.op,rec(schreier := true, log := true));
## Enumerate(o);
##
## if Length(o) <> Length(S!.orb) then
## Error("something is fishy, orbits do not have the same length");
## return fail;
## fi;
##
## # Now we have a nice breadth-first orbit such that all inverses of
## # generators are again generators.
## # First check whether the generators fix the first point:
## for x in gens do
## if o!.op(o[1],x) = o[1] then
## count := count + 1;
## if S!.stab = false then
## r := rec( isone := isone(x) );
## else
## r := SiftGroupElement(S!.stab,x);
## fi;
## if not(r.isone) then
## return [fail,S!.layer,x,"generator"];
## fi;
## fi;
## od;
##
## bl := BlistList([1..Length(o)],[]);
## for d in [1..o!.depth] do
## Info(InfoGenSS,2,"Testing in depth ",d," (of ",o!.depth,") checking ",
## o!.depthmarks[d+1]-o!.depthmarks[d]," points (of ",Length(o),")");
## for i in [o!.depthmarks[d]..o!.depthmarks[d+1]-1] do
## # These indices contain points in depth d
## if bl[o!.schreierpos[i]] then
## # this subtree does not need to be done!
## bl[i] := true; # mark subtree as done
## else
## x := o[i];
## for j in [1..Length(gens)] do
## if j <> invtab[o!.schreiergen[i]] then
## y := o!.op(x,gens[j]);
## p := Position(o,y);
## # if p > i then this is a new point, do nothing
## if p <= i then
## count := count + 1;
## w1 := TraceSchreierTreeForward(o,i);
## w2 := TraceSchreierTreeForward(o,p);
## s := (EvaluateWord(gens,w1)*gens[j])
## /EvaluateWord(gens,w2);
## #Print(Length(gens),w1,j,w2,"\n");
## if S!.stab = false then
## r := rec( isone := isone(s) );
## else
## r := SiftGroupElement(S!.stab,s);
## fi;
## if not(r.isone) then
## return [fail,S!.layer,s,"schreier gen"];
## fi;
## fi;
## if p < o!.depthmarks[d] then
## # this goes up in the tree, this means, if this
## # Schreier generator is OK, we do not have to
## # look below!
## bl[i] := true;
## break;
## fi;
## fi;
## od;
## fi;
## od;
## od;
## Info(InfoGenSS,2,"Have sifted ",count," elements.");
## return true;
## end;
#############################################################################
# The following operations apply to stabilizer chains:
#############################################################################
InstallMethod( Size, "for a stabilizer chain",
[ IsStabilizerChain and IsStabilizerChainByOrb ],
function( S )
local size;
size := 1;
while S <> false do
size := size * Length(S!.orb);
S := S!.stab;
od;
return size;
end );
InstallMethod( \in, "for a group element and a stabilizer chain",
[ IsObject, IsStabilizerChain and IsStabilizerChainByOrb ],
function( el, S )
local r;
r := SiftGroupElement(S,el);
if r.isone then
return true;
else
return false;
fi;
end );
GENSS_VIEWDEPTH := 0; # to please the interpreter
InstallMethod( ViewObj, "for a stabilizer chain",
[ IsStabilizerChain and IsStabilizerChainByOrb ],
function( S )
local i;
if not(IsBound(GENSS_VIEWDEPTH)) then
GENSS_VIEWDEPTH := 0;
else
GENSS_VIEWDEPTH := GENSS_VIEWDEPTH + 1;
fi;
for i in [1..GENSS_VIEWDEPTH] do Print(" "); od;
Print("<stabchain size=",Size(S)," orblen=",Length(S!.orb),
" layer=",S!.layer," SchreierDepth=",S!.orb!.depth,">");
if S!.stab <> false then
Print("\n");
ViewObj(S!.stab);
fi;
GENSS_VIEWDEPTH := GENSS_VIEWDEPTH - 1;
if GENSS_VIEWDEPTH < 0 then
Unbind(GENSS_VIEWDEPTH);
fi;
end );
Unbind(GENSS_VIEWDEPTH);
InstallMethod( BaseStabilizerChain,
"generic method for stabilizer chains",
[IsStabilizerChain],
function( S )
local SS,pts,ops;
pts := [];
ops := [];
SS := S;
while SS <> false do
Add(pts,SS!.orb[1]);
Add(ops,SS!.orb!.op);
SS := SS!.stab;
od;
return rec( points := pts, ops := ops );
end );
InstallMethod( SiftBaseImage,
"generic method for genss stabilizer chains",
[IsStabilizerChain, IsList],
function(S,bi)
local l, SS, i, o, pos;
l := Length(bi);
SS := S;
i := 1;
while i <= Length(bi) do
o := SS!.orb;
pos := Position(o,bi[i]);
if pos = fail then
return false;
fi;
while pos > 1 do
if o!.memorygens then
bi{[i..l]} := GENSS_MapBaseImage(bi{[i..l]},
o!.gensi[o!.schreiergen[pos]]!.el,SS);
else
bi{[i..l]} := GENSS_MapBaseImage(bi{[i..l]},
o!.gensi[o!.schreiergen[pos]],SS);
fi;
pos := o!.schreierpos[pos];
od;
if o[1] <> bi[i] then
Error("this should not have happened, tell the authors");
fi;
i := i + 1;
SS := SS!.stab;
od;
return true;
end );
InstallMethod( IsProved, "for a stabilizer chain",
[ IsStabilizerChain and IsStabilizerChainByOrb ],
function( S )
return S!.proof;
end );
InstallMethod( StabChainOp, "for a permutation group and a stabilizer chain",
[ IsPermGroup, IsStabilizerChain and IsStabilizerChainByOrb ],
function( g, S )
local base;
if HasSize(g) and Size(g) <> Size(S) then
Error("known size of group does not match stabiliser chain");
return fail;
else
SetSize(g,Size(S));
fi;
base := BaseStabilizerChain(S);
if not(ForAll(base.points,IsInt) and ForAll(base.ops,x->x=OnPoints)) then
return StabChainOp(g);
fi;
return StabChainOp(g,rec( base := base.points, knownBase := base.points,
size := Size(S) ));
end );
InstallGlobalFunction( GENSS_FindGensStabilizer,
function( S, pt, opt )
# Assumes that S is a correct stabilizer chain for the group generated by
# its gens (in S!.orb!.gens) and that pt lies in the first orbit S!.orb.
# Finds an SLP to express generators of the stabilizer of pt in
# the group in terms of the gens.
# Assumes that the group and the first stabilizer in the chain are
# non-trivial.
local gens,g,pos,mgens,word,cosetrep,invgens,pr,stabgens,size,stabsize,
randel,newpt,ptpos,ptcosetrep,el,SS,stab,i;
gens := S!.orb!.gens;
g := GroupWithGenerators(gens);
SetSize(g,Size(S));
pos := Position(S!.orb,pt);
if pos = fail then
Error("point pt must lie in first orbit");
return fail;
fi;
if not(IsBound(opt.scramble)) then opt.scramble := 30; fi;
if not(IsBound(opt.scramblefactor)) then opt.scramblefactors := 0; fi;
if not(IsBound(opt.maxdepth)) then opt.maxdepth := 200; fi;
if not(IsBound(opt.addslots)) then
opt.addslots := Maximum(0,5-Length(gens));
fi;
# Give the gens a (new) memory:
if IsObjWithMemory(gens[1]) then
mgens := GeneratorsWithMemory(List(gens,x->x!.el));
else
mgens := GeneratorsWithMemory(gens);
fi;
# FIXME: use GENSS_Prod
if pos = 1 then
cosetrep := mgens[1]^0;
else
word := TraceSchreierTreeForward(S!.orb,pos);
cosetrep := Product(mgens{word});
fi;
invgens := List(mgens,x->x^-1);
# Set up the product replacer:
pr := ProductReplacer( mgens, opt );
stabgens := []; # here we collect the stabilizer generators
size := 1;
stabsize := Size(S!.stab); # this has to exist!
repeat
randel := Next(pr);
newpt := S!.orb!.op(pt,randel!.el);
ptpos := Position(S!.orb,newpt);
word := TraceSchreierTreeBack(S!.orb,ptpos);
if Length(word) > 0 then
ptcosetrep := Product(invgens{word});
el := randel * ptcosetrep;
else
el := randel;
fi;
if pos <> 1 then
el := el * cosetrep;
fi;
# now el is an element in the stabilizer of pt
if S!.IsOne(el) then
Info(InfoGenSS,3,"Found trivial stabilizer element.");
else
Add(stabgens,el);
stab := GroupWithGenerators(stabgens);
SS := StabilizerChain(stab,rec( Base := S, IsOne := S!.IsOne ));
Info(InfoGenSS,1,"Have group size ",Size(SS)," (of ",
stabsize,")");
if Size(SS) = size then
Remove(stabgens);
else
size := Size(SS);
fi;
fi;
until size = stabsize;
# Try leaving out the first generators found:
i := 1;
repeat
Info(InfoGenSS,1,"Need ",Length(stabgens)+1-i," generators.");
i := i + 1;
if i < Length(stabgens) then
stab := Group(stabgens{[i..Length(stabgens)]});
SS := StabilizerChain(stab,rec(Base := S, IsOne := S!.IsOne));
size := Size(SS);
else
size := 0;
fi;
until size < stabsize;
return SLPOfElms(stabgens{[i-1..Length(stabgens)]});
end );
InstallGlobalFunction( GENSS_FindShortGensStabilizerOld,
function( S, pt )
# Assumes that S is a correct stabilizer chain for the group generated by
# its gens (in S!.orb!.gens) and that pt lies in the first orbit S!.orb.
# Finds an SLP to express generators of the stabilizer of pt in
# the group in terms of the gens.
# Assumes that the group and the first stabilizer in the chain are
# non-trivial.
local gens,g,pos,mgens,word,cosetrep,invgens,pr,stabgens,size,stabsize,
randel,newpt,ptpos,ptcosetrep,el,SS,stab,i,j;
gens := S!.orb!.gens;
g := GroupWithGenerators(gens);
SetSize(g,Size(S));
pos := Position(S!.orb,pt);
if pos = fail then
Error("point pt must lie in first orbit");
return fail;
fi;
# Give the gens a (new) memory:
if IsObjWithMemory(gens[1]) then
mgens := GeneratorsWithMemory(List(gens,x->x!.el));
else
mgens := GeneratorsWithMemory(gens);
fi;
if pos = 1 then
cosetrep := mgens[1]^0;
else
word := TraceSchreierTreeForward(S!.orb,pos);
cosetrep := Product(mgens{word});
fi;
invgens := List(mgens,x->x^-1);
# Set up the search:
pr := ShallowCopy(mgens);
i := 1; # counts points
j := 1; # counts gens
stabgens := []; # here we collect the stabilizer generators
size := 1;
stabsize := Size(S!.stab); # this has to exist!
repeat
Add(pr,pr[i]*mgens[j]);
j := j + 1;
if j > Length(mgens) then
j := 1;
i := i + 1;
fi;
randel := pr[Length(pr)];
newpt := S!.orb!.op(pt,randel!.el);
ptpos := Position(S!.orb,newpt);
word := TraceSchreierTreeBack(S!.orb,ptpos);
if Length(word) > 0 then
ptcosetrep := Product(invgens{word});
el := randel * ptcosetrep;
else
el := randel;
fi;
if pos <> 1 then
el := el * cosetrep;
fi;
# now el is an element in the stabilizer of pt
if S!.IsOne(el) then
Info(InfoGenSS,2,"Found trivial stabilizer element.");
else
Add(stabgens,el);
stab := GroupWithGenerators(stabgens);
SS := StabilizerChain(stab,rec( Base := S, IsOne := S!.IsOne ));
Info(InfoGenSS,1,"Have group size ",Size(SS)," (of ",
stabsize,")");
if Size(SS) = size then
Remove(stabgens);
else
size := Size(SS);
fi;
fi;
until size = stabsize;
# Try leaving out the first generators found:
i := 1;
repeat
Info(InfoGenSS,1,"Need ",Length(stabgens)+1-i," generators.");
i := i + 1;
if i < Length(stabgens) then
stab := Group(stabgens{[i..Length(stabgens)]});
SS := StabilizerChain(stab,rec(Base := S, IsOne := S!.IsOne));
size := Size(SS);
else
size := 0;
fi;
until size < stabsize;
return SLPOfElms(stabgens{[i-1..Length(stabgens)]});
end );
InstallGlobalFunction( GENSS_FindShortGensStabilizer,
function( gens, pt, op, grpsize, stabsize, S )
# Assumes that S is a correct stabilizer chain a supergroup of the
# group g generated by gens. pt and op is an action for g.
# grpsize is the Size of g and stabsize is the size of the stabilizer.
# Finds an SLP to express generators of the stabilizer of pt in
# g in terms of the gens.
# Assumes that g and the orbit of pt under g are non-trivial.
local SS,el,g,i,invgens,j,mgens,newpt,orb,pr,ptcosetrep,ptpos,randel,size,
stab,stabgens,word;
g := GroupWithGenerators(gens);
Info(InfoGenSS,2,"Enumerating orbit of length ",grpsize/stabsize);
orb := Orb(gens,pt,op,rec(schreier := true,stabsizebound := stabsize));
Enumerate(orb);
# Give the gens a (new) memory:
if IsObjWithMemory(gens[1]) then
mgens := GeneratorsWithMemory(List(gens,x->x!.el));
else
mgens := GeneratorsWithMemory(gens);
fi;
invgens := List(mgens,x->x^-1);
# Set up the search:
pr := ShallowCopy(mgens);
i := 1; # counts points
j := 1; # counts gens
stabgens := []; # here we collect the stabilizer generators
size := 1;
repeat
Add(pr,pr[i]*mgens[j]);
j := j + 1;
if j > Length(mgens) then
j := 1;
i := i + 1;
fi;
randel := pr[Length(pr)];
newpt := op(pt,randel!.el);
ptpos := Position(orb,newpt);
word := TraceSchreierTreeBack(orb,ptpos);
if Length(word) > 0 then
ptcosetrep := Product(invgens{word});
el := randel * ptcosetrep;
else
el := randel;
fi;
# now el is an element in the stabilizer of pt
if S!.IsOne(el) then
Info(InfoGenSS,3,"Found trivial stabilizer element.");
else
Add(stabgens,el);
stab := GroupWithGenerators(stabgens);
SS := StabilizerChain(stab,rec( Base := S, IsOne := S!.IsOne ));
Info(InfoGenSS,1,"Have group size ",Size(SS)," (of ",
stabsize,")");
if Size(SS) = size then
Remove(stabgens);
else
size := Size(SS);
fi;
fi;
until size >= stabsize;
if size > stabsize then
Error("Something is wrong, stabsize limit was too small");
return fail;
fi;
# Try leaving out the first generators found:
i := 1;
Info(InfoGenSS,1,"Need ",Length(stabgens)," generators.");
while i < Length(stabgens) do
stab := Group(stabgens{Concatenation([1..i-1],[i+1..Length(stabgens)])});
SS := StabilizerChain(stab,rec(Base := S, IsOne := S!.IsOne));
if Size(SS) < stabsize then
i := i + 1; # keep generator
else
Remove(stabgens,i);
Info(InfoGenSS,1,"Need ",Length(stabgens)," generators.");
fi;
od;
return SLPOfElms(stabgens);
end );
InstallGlobalFunction( SLPChainStabilizerChain,
function( S, gens )
# S a stabilizer chain assumed to be correct for the group g generated
# by generators gens.
# Returns a list of straight line programs expressing successively
# the stabilizers in the chain, each in terms of the generators of the
# previous, the first in terms of gens.
local l,ll,slp;
l := [];
ll := [Size(S)];
while S!.stab <> false do
Info(InfoGenSS,1,"Working on group with size ",Size(S),"...");
slp := GENSS_FindShortGensStabilizer(
gens, S!.orb[1], S!.orb!.op, Size(S), Size(S!.stab), S );
Add(l,slp);
Add(ll,Size(S!.stab));
gens := ResultOfStraightLineProgram(slp,gens);
Info(InfoGenSS,1,"Found SLP with ",
Length(LinesOfStraightLineProgram(slp))," lines to compute ",
Length(gens)," generators.");
S := S!.stab;
od;
return rec(slps := l,sizes := ll);
end );
InstallGlobalFunction( GENSS_MakeIterRecord,
function( S )
local SS,Slist,it;
Slist := [S];
SS := S!.stab;
while SS <> false do
Add(Slist,SS);
SS := SS!.stab;
od;
it := rec( NextIterator := GENSS_GroupNextIterator,
IsDoneIterator := GENSS_GroupIsDoneIterator,
ShallowCopy := GENSS_GroupShallowCopy,
pos := ListWithIdenticalEntries( Length( S!.base ), 1 ),
Slist := Slist,
cache := ListWithIdenticalEntries( Length(S!.base),
One(S!.orb!.gens[1]) ),
one := One(S!.orb!.gens[1]),
lens := List(Slist,x->Length(x!.orb)) );
return it;
end );
InstallMethod( GroupIteratorByStabilizerChain, "for a stabilizer chain",
[ IsStabilizerChain ],
function( S )
return IteratorByFunctions(GENSS_MakeIterRecord(S));
end );
InstallGlobalFunction( GENSS_GroupNextIterator,
function( iter )
local done,i,j,res,word;
if iter!.pos[1] > iter!.lens[1] then return fail; fi;
res := iter!.cache[Length(iter!.cache)];
# Now step forwards:
i := Length(iter!.Slist);
repeat
iter!.pos[i] := iter!.pos[i] + 1;
if iter!.pos[i] > iter!.lens[i] then
iter!.pos[i] := 1;
i := i - 1;
done := false;
else
done := true;
fi;
until done or i = 0;
if i = 0 then
iter!.pos[1] := iter!.lens[1]+1;
return res;
fi; # we are done
word := TraceSchreierTreeForward(iter!.Slist[i]!.orb,iter!.pos[i]);
iter!.cache[i] := Product(iter!.Slist[i]!.orb!.gens{word});
if i > 1 then iter!.cache[i] := iter!.cache[i] * iter!.cache[i-1]; fi;
for j in [i+1..Length(iter!.Slist)] do
iter!.cache[j] := iter!.cache[j-1];
od;
return res;
end );
InstallGlobalFunction( GENSS_GroupIsDoneIterator,
function( iter )
return iter!.pos[1] > iter!.lens[1];
end );
InstallGlobalFunction( GENSS_GroupShallowCopy,
function( iter )
return GENSS_MakeIterRecord( iter!.Slist[1] );
end );
#############################################################################
# We can store a stabilizer chain in the attribute
# StoredStabilizerChain, then the following methods for groups apply:
#############################################################################
InstallMethod( SetStabilizerChain, "for a group and a stabilizer chain",
[IsGroup, IsStabilizerChain],
function(g,S)
if IsIdenticalObj(S!.IsOne,IsOne) and
HasSize(g) and Size(g) <> Size(S) then
Error("you try to set a stabilizer chain for the wrong group");
return;
fi;
SetStoredStabilizerChain(g,S);
end );
InstallMethod( Size, "for a group with a stored stabilizer chain",
[IsGroup and HasStoredStabilizerChain],
function(g)
return Size(StoredStabilizerChain(g));
end );
InstallMethod( Size, "for a perm. group with a stored stabilizer chain",
[IsPermGroup and HasStoredStabilizerChain],
function(g)
return Size(StoredStabilizerChain(g));
end );
InstallMethod( Size, "for a group with a stored stabilizer chain",
[IsGroup and HasStoredStabilizerChain and IsHandledByNiceMonomorphism],
function(g)
return Size(StoredStabilizerChain(g));
end );
InstallMethod( \in, "for a group elm and a group with stored stabilizer chain",
[IsPerm, IsPermGroup and HasStoredStabilizerChain], 1,
function(x, g)
local S,r;
S := StoredStabilizerChain(g);
r := SiftGroupElement(S,x);
return r.isone;
end );
InstallMethod( \in, "for a group elm and a group with stored stabilizer chain",
[IsObject, IsMatrixGroup and HasStoredStabilizerChain and
IsHandledByNiceMonomorphism],
function(x, g)
local S,r;
S := StoredStabilizerChain(g);
r := SiftGroupElement(S,x);
return r.isone;
end );
InstallMethod( \in, "for a group elm and a group with stored stabilizer chain",
[IsObject, IsGroup and HasStoredStabilizerChain],
function(x, g)
local S,r;
S := StoredStabilizerChain(g);
r := SiftGroupElement(S,x);
return r.isone;
end );
InstallOtherMethodWithRandomSource( Random, "for a random source and a stabilizer chain",
[ IsRandomSource, IsStabilizerChainByOrb ],
function(rs, S)
local x,pos,w;
x := One(S!.orb!.gens[1]);
while S <> false do
pos := Random(rs, 1,Length(S!.orb));
w := TraceSchreierTreeForward(S!.orb,pos);
x := GENSS_Prod(S!.orb!.gens,w) * x;
S := S!.stab;
od;
return x;
end );
InstallMethodWithRandomSource( Random, "for a random source and a group with a stored stabilizer chain",
[ IsRandomSource, IsGroup and HasStoredStabilizerChain ],
function(rs, g) return Random(rs, StoredStabilizerChain(g)); end );
InstallMethodWithRandomSource( Random, "for a random source and a permgroup with a stored stabilizer chain",
[ IsRandomSource, IsPermGroup and HasStoredStabilizerChain ], 10,
function(rs, g) return Random(rs, StoredStabilizerChain(g)); end );
InstallMethodWithRandomSource( Random, "for a random source and a matrixgroup with a stored stabilizer chain",
[ IsRandomSource, IsMatrixGroup and IsHandledByNiceMonomorphism and
HasStoredStabilizerChain ],
function(rs, g) return Random(rs, StoredStabilizerChain(g)); end );
InstallMethod( SizeMC, "for a group and an error bound",
[IsGroup, IsRat],
function( G, err )
local S;
S := StabilizerChain(G,rec( ErrorBound := err ));
return Size(S);
end );
InstallMethod( SizeMC, "for a permutation group and an error bound, genss",
[IsGroup and IsPermGroup, IsRat], 1,
function( G, err )
local S;
S := StabilizerChain(G,rec( ErrorBound := err ));
return Size(S);
end );
InstallGlobalFunction( GENSS_ImageElm,
function( data, x )
local r;
r := SiftGroupElementSLP(data!.Sg,x);
if not(r.isone) then return fail; fi;
return ResultOfStraightLineProgram(r.slp,data!.stronggims);
end );
InstallGlobalFunction( GENSS_PreImagesRepresentative,
function( data, x )
local r;
r := SiftGroupElementSLP(data!.Si,x);
if not(r.isone) then return fail; fi;
return ResultOfStraightLineProgram(r.slp,data!.strongipre);
end );
InstallGlobalFunction( GroupHomomorphismByImagesNCStabilizerChain,
function( g, h, images, opt1, opt2 )
local Sg,Si,data,gm,im,slpstrongg,slpstrongi,strongg,stronggims,
strongi,strongipre;
gm := GroupWithMemory(g);
Sg := StabilizerChain(gm,opt1);
strongg := StrongGenerators(Sg);
slpstrongg := SLPOfElms(strongg);
ForgetMemory(Sg);
stronggims := ResultOfStraightLineProgram(slpstrongg,images);
im := GroupWithMemory(images);
Si := StabilizerChain(im,opt2);
strongi := StrongGenerators(Si);
slpstrongi := SLPOfElms(strongi);
ForgetMemory(Si);
strongipre := ResultOfStraightLineProgram(slpstrongi,GeneratorsOfGroup(g));
data := rec( Sg := Sg, stronggims := stronggims,
Si := Si, strongipre := strongipre,
slpstrongg := slpstrongg, slpstrongi := slpstrongi );
return GroupHomByFuncWithData( g, h, GENSS_ImageElm, false,
GENSS_PreImagesRepresentative, data );
end );
InstallMethod( ORB_StabilizerChainKnownSize,
"GENSS method for arbitrary groups",
[IsGroup,IsPosInt],
function(g,size)
if HasStoredStabilizerChain(g) and
size = Size(StoredStabilizerChain(g)) then
return StoredStabilizerChain(g);
fi;
return StabilizerChain(g,rec( Size := size ));
end );
InstallMethod( ORB_BaseStabilizerChain,
"GENSS method for arbitary groups",
[IsStabilizerChain],
BaseStabilizerChain );
InstallMethod( ORB_StabilizerChainKnownBase,
"GENSS method for arbitrary groups",
[IsGroup,IsObject],
function(g,base)
if HasStoredStabilizerChain(g) then
return StoredStabilizerChain(g);
fi;
if base = fail then
return StabilizerChain( g );
else
return StabilizerChain( g, rec( Cand := ShallowCopy(base) ) );
fi;
end );
InstallMethod( ORB_SizeStabilizerChain,
"GENSS method for arbitrary groups",
[IsStabilizerChain], Size );
InstallMethod( ORB_IsWordInStabilizerChain,
"GENSS method for arbitrary groups",
[IsList,IsList,IsList,IsStabilizerChain],
function(word, gens, gensi, S)
local x, b, ops, i;
if Size(S) = 1 then
x := ORB_ApplyWord(gens[1]^0,word,gens,gensi,OnRight);
return S!.IsOne(x);
fi;
b := BaseStabilizerChain(S);
ops := b.ops;
b := b.points;
for i in [1..Length(word)] do
if word[i] > 0 then
b := List([1..Length(b)],j->ops[j](b[j],gens[word[i]]));
else
b := List([1..Length(b)],j->ops[j](b[j],gensi[-word[i]]));
fi;
od;
return SiftBaseImage(S,b);
end );
InstallMethod( ORB_IsElementInStabilizerChain,
"GENSS method for arbitrary groups",
[IsObject, IsStabilizerChain],
function(el, S)
return el in S;
end );
#############################################################################
# The following methods are about methods to compute stabilisers:
#############################################################################
InstallMethod( Stab, "add empty options record",
[IsGroup, IsObject, IsFunction],
function( g, x, op )
return Stab(g,x,op,rec());
end );
InstallMethod( Stab, "add empty options record",
[IsList, IsObject, IsFunction],
function( l, x, op )
return Stab(l,x,op,rec());
end );
InstallMethod( Stab, "create group from list of generators",
[IsList, IsObject, IsFunction, IsRecord],
function( l, x, op, opt )
return Stab(GroupWithGenerators(l),x,op,opt);
end );
# Now the real one:
InstallMethod( Stab, "by Orb orbit enumeration",
[IsGroup, IsObject, IsFunction, IsRecord],
function( g, x, op, opt )
local S,count,el,errorprob,found,gens,i,j,limit,memperpt,nrrand,o,
pat,pos,pr,res,stab,stabchain,stabel,stabgens,stabsizeest,w1,w2,y;
GENSS_CopyDefaultOptions(GENSS,opt);
if HasStoredStabilizerChain(g) then
S := StoredStabilizerChain(g);
if IsIdenticalObj(S!.orb!.op,op) and x in S!.orb then
if S!.stab = false then
y := rec( stab := TrivialSubgroup(g), size := 1, proof := true );
y.stabilizerchain := StabilizerChain(y.stab);
return y;
fi;
pos := Position(S!.orb,x);
w1 := TraceSchreierTreeForward(S!.orb,pos);
y := GENSS_Prod(S!.orb!.gens,w1);
stab := Group(List(S!.stab!.orb!.gens,a->y^-1*a*y));
return rec( stab := stab, size := Size(S!.stab),
stabilizerchain := StabilizerChain(stab,rec( Base := S)),
proof := true );
fi;
fi;
# Now we have to do it ourselves:
if not(IsBound(opt.ErrorBound)) then # we are working deterministically
if not(HasSize(g)) then
# We are basically stuffed, unless we want to check all
# Schreier generators of an orbit!
Error("I do not want to check all Schreier generators, ",
"need size of group or ErrorBound option!");
return fail;
fi;
else
if opt.ErrorBound < opt.StabAssumeCompleteLimit then
opt.StabAssumeCompleteLimit := opt.ErrorBound;
fi;
fi;
# Now we either know the group size or work Monte Carlo:
errorprob := 1;
limit := opt.StabInitialLimit;
pat := opt.StabInitialPatience;
o := Orb(g,x,op,rec( report := opt.Report,
treehashsize := opt.InitialHashSize,
schreier := true ) );
stabgens := [];
stabsizeest := 1;
pr := ProductReplacer(GeneratorsOfGroup(g),
rec( scramble := opt.StabScramble,
scramblefactor := opt.StabScrambleFactor,
addslots := opt.StabAddSlots,
maxdepth := opt.StabMaxDepth ));
gens := GeneratorsOfGroup(g);
# Now go through a cycle of orbit enumeration and stabilizer generation:
repeat
if not(IsClosed(o)) then
if HasSize(g) and limit > QuoInt(Size(g),stabsizeest*2)+1 then
limit := QuoInt(Size(g),stabsizeest*2)+2;
fi;
Info(InfoGenSS,2,"Enumerating orbit with limit ",limit);
Enumerate(o,limit);
if IsClosed(o) then
Info(InfoGenSS,2,"Orbit closed, size is ",Length(o));
fi;
fi;
if errorprob < opt.StabAssumeCompleteLimit then
if HasSize(g) and
2 * Length(o) * stabsizeest > Size(g) then
# Done!
#stab := Subgroup(g,stabgens);
if Length(stabgens) > 0 then
stab := Group(stabgens);
SetParent(stab,g);
SetSize(stab,stabsizeest);
SetStabilizerChain(stab,stabchain);
else
stab := TrivialSubgroup(g);
stabchain := StabilizerChain(stab);
fi;
return rec( stab := stab, size := stabsizeest,
stabilizerchain := stabchain,
proof := true );
fi;
limit := 2*limit;
if not IsClosed(o) then continue; fi;
fi;
count := 0;
found := false;
repeat
el := Next(pr);
i := 1;
while i <= Length(o) and not(found) do
y := op(o[i],el);
j := Position(o,y);
if j <> fail then
found := true;
else
i := i + 1;
fi;
od;
count := count + 1;
until found or count >= pat;
if not(IsClosed(o)) then
if not(found) then
limit := QuoInt(3*limit,2);
else
if count < 10 then
limit := limit + 100;
else
limit := QuoInt(3*limit,2);
fi;
fi;
fi;
pat := QuoInt(pat*5,4);
if found then # Have a potential stabiliser element:
Info(InfoGenSS,3,"Found a potential stabilising element...");
w1 := TraceSchreierTreeForward(o,i);
w2 := TraceSchreierTreeForward(o,j);
stabel := EvaluateWord(gens,w1)*el/EvaluateWord(gens,w2);
if IsOne(stabel) then
errorprob := errorprob / 2;
Info(InfoGenSS,3,"... which was the identity.");
else
Add(stabgens,stabel);
if Length(stabgens) < 2 then
Info(InfoGenSS,3,"Waiting for second stabiliser element.");
else # Length(stabgens) >= 2 then
if Length(stabgens) = 2 then
Info(InfoGenSS,2,"Estimating stabilizer...");
stabchain := StabilizerChain(Group(stabgens));
errorprob := 1;
else
res := AddGeneratorToStabilizerChain(stabchain,stabel);
if not(res) then
Remove(stabgens,Length(stabgens));
errorprob := errorprob / 2;
Info(InfoGenSS,2,"Error probablity now < ",
errorprob);
else
Info(InfoGenSS,2,
"Added generator to stabilizer...");
VerifyStabilizerChainMC(stabchain,10);
errorprob := 1;
fi;
fi;
if Size(stabchain) > stabsizeest then
stabsizeest := Size(stabchain);
Info(InfoGenSS,2,"New stabilizer estimate: ",
stabsizeest);
else
Info(InfoGenSS,2,"Stabiliser estimate unchanged.");
fi;
if HasSize(g) and
2 * Length(o) * stabsizeest > Size(g) then
# Done!
if Length(stabgens) > 0 then
stab := Group(stabgens);
SetParent(stab,g);
SetSize(stab,stabsizeest);
SetStabilizerChain(stab,stabchain);
else
stab := TrivialSubgroup(g);
stabchain := StabilizerChain(stab);
fi;
return rec( stab := stab, size := stabsizeest,
stabilizerchain := stabchain,
proof := true );
fi;
if IsBound(opt.ErrorBound) and
errorprob < opt.ErrorBound then
#stab := Subgroup(g,stabgens);
if Length(stabgens) > 0 then
stab := Group(stabgens);
SetParent(stab,g);
else
stab := TrivialSubgroup(g);
stabchain := StabilizerChain(stab);
fi;
return rec( stab := stab, size := stabsizeest,
stabilizerchain := stabchain,
proof := false );
fi;
fi;
fi;
else # no element found
Info(InfoGenSS,3,"No stabiliser element found!");
fi;
until Length(o) > opt.StabOrbitLimit;
return fail;
end );
InstallMethod( StabMC, "add empty options record",
[IsGroup, IsObject, IsFunction],
function( g, x, op )
return StabMC(g,x,op,rec());
end );
InstallMethod( StabMC, "add empty options record",
[IsList, IsObject, IsFunction],
function( l, x, op )
return StabMC(l,x,op,rec());
end );
InstallMethod( StabMC, "create group from list of generators",
[IsList, IsObject, IsFunction, IsRecord],
function( l, x, op, opt )
return StabMC(GroupWithGenerators(l),x,op,opt);
end );
InstallMethod( StabMC, "call Stab with errorbound",
[IsGroup, IsObject, IsFunction, IsRecord],
function( g, x, op, opt )
if not(IsBound(opt.ErrorBound)) then
opt.ErrorBound := 1/1025;
fi;
if not(IsBound(opt.DoEstimate)) then
opt.DoEstimate := 1000; # try for 1 second
fi;
return Stab(g,x,op,opt);
end );
InstallGlobalFunction( BacktrackSearchStabilizerChainElement,
function( S, # a stabilizer chain
P, # a property function G -> {true,false}
g, # a group element
pruner )
# Let G be a group and U being the subgroup defined by S.
# Does a backtrack search in U using S for an element x of U
# for which P(xg) is true. g not necessarily in U!
# pruner is either false or a function taking arguments as seen below.
local cache,gen,i,ii,res,t,w,x,dotree;
cache := WeakPointerObj([[S!.orb!.gens[1]^0,[]]]);
for i in [1..Length(S!.orb)] do
ii := S!.orb!.orbind[i];
if ii > 1 then
t := ElmWPObj(cache,S!.orb!.schreierpos[ii]);
if t <> fail then
gen := S!.orb!.schreiergen[ii];
w := EmptyPlist(Length(t[2])+1);
Append(w,t[2]);
Add(w,gen);
t := t[1] * S!.orb!.gens[gen];
else
w := TraceSchreierTreeForward(S!.orb,ii);
t := Product(S!.orb!.gens{w});
fi;
SetElmWPObj(cache,ii,[t,w]);
x := t*g;
else
w := [];
t := S!.orb!.gens[1]^0;
x := g;
fi;
if S!.stab = false then
if P(x) then
return rec( elm := x, vec := [i], word := S!.layergens{w} );
fi;
else
if pruner = false or pruner(S,ii,x,t,w) then
res := BacktrackSearchStabilizerChainElement(S!.stab,P,x,
pruner);
if res <> fail then
Add(res.vec,ii,1);
res.word := Concatenation(res.word,S!.layergens{w});
return res;
fi;
fi;
fi;
od;
return fail;
end );
InstallGlobalFunction( ComputeSuborbitsForStabilizerChain,
function( S )
local SS,gens;
SS := S;
while SS!.stab <> false do
gens := StrongGenerators(S){SS!.stab!.layergens};
SS!.suborbs := FindSuborbits(S!.orb,gens);
SS := SS!.stab;
od;
end );
InstallGlobalFunction( BacktrackSearchStabilizerChainSubgroup,
function(S,P,pruner)
# Let G be a group and U being the subgroup defined by S.
# Does a backtrack search in U using S for the subgroup H of U
# with H := {h \in U | P(h)=true}. P must be such that H is
# a subgroup of U.
# pruner is either false or a function taking arguments as seen below.
local SS,SSS,cache,dotree,gen,i,ii,newgens,o,res,t,w,T;
cache := WeakPointerObj([[S!.orb!.gens[1]^0,[]]]);
if S!.stab = false then # lowest level
newgens := [];
o := false;
for i in [2..Length(S!.orb)] do
ii := S!.orb!.orbind[i];
if o = false or not(S!.orb[ii] in o) then
if ii > 1 then
t := ElmWPObj(cache,S!.orb!.schreierpos[ii]);
if t <> fail then
gen := S!.orb!.schreiergen[ii];
w := EmptyPlist(Length(t[2])+1);
Append(w,t[2]);
Add(w,gen);
t := t[1] * S!.orb!.gens[gen];
else
w := TraceSchreierTreeForward(S!.orb,ii);
t := Product(S!.orb!.gens{w});
fi;
SetElmWPObj(cache,ii,[t,w]);
else
w := [];
t := S!.orb!.gens[1]^0;
fi;
if P(t) then
Add(newgens,t);
if Length(newgens) = 1 then
o := Orb(ShallowCopy(newgens),S!.orb[1],S!.orb!.op,
rec( treehashsize := 200, schreier := true,
log := true ));
else
AddGeneratorsToOrbit(o,[t]);
fi;
Enumerate(o);
fi;
fi;
od;
if o = false then # No solution found
newgens := [S!.orb!.gens[1]^0];
o := Orb(newgens,S!.orb[1],S!.orb!.op,
rec(schreier := true, log := true));
Enumerate(o);
fi;
SSS := rec( stab := false, size := Length(o), base := [o[1]],
opt := ShallowCopy(S!.opt), layer := S!.layer,
orb := o, stronggens := newgens,
layergens := [1..Length(newgens)], randpool := [],
IsOne := S!.opt.IsOne, proof := true,
parentS := false);
Objectify( StabChainByOrbType, SSS );
SSS!.opt.pr := ProductReplacer(SSS!.stronggens);
SSS!.topS := SSS; # this will be changed further up!
o!.stabilizerchain := SSS;
o!.gensi := List(o!.gens,x->x^-1);
Info(InfoGenSS,2,"Layer ",SSS!.layer," completed, size ",Size(SSS));
return SSS;
else # S!.stab <> false
# First go down in tree, essentially trying coset rep 1 first:
SS := BacktrackSearchStabilizerChainSubgroup(S!.stab,P,pruner);
newgens := [];
o := false;
for i in [2..Length(S!.orb)] do
ii := S!.orb!.orbind[i];
if o = false or not(S!.orb[ii] in o) then
if ii > 1 then
t := ElmWPObj(cache,S!.orb!.schreierpos[ii]);
if t <> fail then
gen := S!.orb!.schreiergen[ii];
w := EmptyPlist(Length(t[2])+1);
Append(w,t[2]);
Add(w,gen);
t := t[1] * S!.orb!.gens[gen];
else
w := TraceSchreierTreeForward(S!.orb,ii);
t := Product(S!.orb!.gens{w});
fi;
SetElmWPObj(cache,ii,[t,w]);
else
w := [];
t := S!.orb!.gens[1]^0;
fi;
if pruner = false or pruner(S,ii,t,t,w) then
res := BacktrackSearchStabilizerChainElement(S!.stab,P,t,
pruner);
if res <> fail then
Add(newgens,res.elm);
if Length(newgens) = 1 then
o := Orb(Concatenation(StrongGenerators(SS),newgens),
S!.orb[1],S!.orb!.op,
rec( treehashsize :=
Maximum(100,QuoInt(Length(S!.orb),3)),
schreier := true, log := true ));
else
AddGeneratorsToOrbit(o,[res.elm]);
fi;
Enumerate(o);
fi;
fi;
fi;
od;
if o = false then # No solution found
o := Orb([S!.orb!.gens[1]^0],S!.orb[1],S!.orb!.op,
rec(schreier := true, log := true));
Enumerate(o);
fi;
SSS := rec( stab := SS, size := Length(o)*Size(SS),
base := SS!.base, opt := ShallowCopy(S!.opt),
layer := S!.layer,
orb := o, stronggens := SS!.stronggens,
layergens := [1..Length(SS!.stronggens)+Length(newgens)],
randpool := [], IsOne := SS!.opt.IsOne, proof := true,
parentS := false );
Add(SSS.base,o[1],1);
Append(SSS.stronggens,newgens);
SSS.opt.pr := ProductReplacer(ShallowCopy(SSS.stronggens));
Objectify( StabChainByOrbType, SSS );
SS!.parentS := SSS;
T := SSS;
while T <> false do
T!.topS := SSS;
T := T!.stab;
od;
o!.stabilizerchain := SSS;
o!.gensi := List(o!.gens,x->x^-1);
Info(InfoGenSS,2,"Layer ",SSS!.layer," completed, size ",Size(SSS));
return SSS;
fi;
end );
InstallMethod( FindShortGeneratorsOfSubgroup,
"without option rec or func, permutation group method",
[ IsPermGroup, IsPermGroup ],
function(G,U)
local S;
S := StabilizerChain(U);
return FindShortGeneratorsOfSubgroup(G,U,
rec( membershiptest :=
function(a,b)
return SiftGroupElement(S,a).isone;
end,
sizetester :=
function(a)
return Size(StabilizerChain(a));
end ));
end );
InstallMethod( FindShortGeneratorsOfSubgroup,
"without option rec or func, matrix group method",
[ IsMatrixGroup, IsMatrixGroup ],
function(G,U)
local S;
S := StabilizerChain(U);
return FindShortGeneratorsOfSubgroup(G,U,
rec( membershiptest :=
function(a,b)
return SiftGroupElement(S,a).isone;
end,
sizetester :=
function(a)
return Size(StabilizerChain(a));
end ));
end );
