|
#############################################################################
##
## This file is part of GAP, a system for computational discrete algebra.
## This file's authors include Alexander Hulpke.
##
## Copyright of GAP belongs to its developers, whose names are too numerous
## to list here. Please refer to the COPYRIGHT file for details.
##
## SPDX-License-Identifier: GPL-2.0-or-later
##
## This file contains the generic operations for cosets.
##
#############################################################################
##
#R IsRightCosetDefaultRep
##
DeclareRepresentation( "IsRightCosetDefaultRep",
IsComponentObjectRep and IsAttributeStoringRep and IsRightCoset, [] );
#############################################################################
##
#M Enumerator
##
BindGlobal( "NumberElement_RightCoset", function( enum, elm )
return Position( enum!.groupEnumerator, elm / enum!.representative, 0 );
end );
BindGlobal( "ElementNumber_RightCoset", function( enum, pos )
return enum!.groupEnumerator[ pos ] * enum!.representative;
end );
InstallMethod( Enumerator,
"for a right coset",
[ IsRightCoset ],
function( C )
local enum;
enum:= EnumeratorByFunctions( C, rec(
NumberElement := NumberElement_RightCoset,
ElementNumber := ElementNumber_RightCoset,
groupEnumerator := Enumerator( ActingDomain( C ) ),
representative := Representative( C ) ) );
SetLength( enum, Size( ActingDomain( C ) ) );
return enum;
end );
#############################################################################
##
#R IsDoubleCosetDefaultRep
##
DeclareRepresentation( "IsDoubleCosetDefaultRep",
IsComponentObjectRep and IsAttributeStoringRep and IsDoubleCoset, [] );
InstallMethod(ComputedAscendingChains,"init",true,[IsGroup],0,G->[]);
#############################################################################
##
#F AscendingChain(<G>,<U>) . . . . . . . chain of subgroups G=G_1>...>G_n=U
##
InstallGlobalFunction( AscendingChain, function(G,U)
local c,i;
if not IsSubgroup(G,U) then
Error("not subgroup");
fi;
c:=ComputedAscendingChains(U);
i:=PositionProperty(c,i->i[1]=G);
if i=fail then
i:=AscendingChainOp(G,U);
Add(c,[G,i]);
return i;
else
return c[i][2];
fi;
end );
# Find element in G to conjugate B into A
# call with G,A,B;
InstallGlobalFunction(DoConjugateInto,function(g,a,b,onlyone)
local cla,clb,i,j,k,bd,r,rep,b2,dc,
gens,conjugate;
Info(InfoCoset,2,"call DoConjugateInto ",Size(g)," ",Size(a)," ",Size(b));
conjugate:=function(act,asub,genl,nr)
local i,dc,j,z,r,r2,found;
found:=[];
Info(InfoCoset,2,"conjugate ",Size(act)," ",Size(asub)," ",nr);
z:=Centralizer(act,genl[nr]);
if Index(act,z)<Maximum(List(cla[nr],Size)) then
Info(InfoCoset,2,"!orbsize ",Index(act,z));
# asub orbits on the act-class of genl[nr]
dc:=DoubleCosetRepsAndSizes(act,z,asub);
for j in dc do
z:=genl[nr]^j[1];
if z in a then
r:=j[1];
if nr=Length(genl) then
Add(found,r);
if onlyone then return found; fi;
else
r2:=conjugate(Centralizer(act,z),Centralizer(asub,z),
List(genl,x->x^r),nr+1);
if Length(r2)>0 then
Append(found,r*r2);
if onlyone then return found; fi;
fi;
fi;
fi;
od;
else
for i in cla[nr] do
Info(InfoCoset,2,"!classize ",Size(i)," ",
Index(act,Centralizer(act,genl[nr]))," ",
QuoInt(Size(a),Size(Centralizer(i))*Size(asub)));
# split up a-classes to asub-classes
dc:=DoubleCosetRepsAndSizes(a,Centralizer(i),asub);
Info(InfoCoset,2,Length(dc)," double cosets");
for j in dc do
z:=Representative(i)^j[1];
r:=RepresentativeAction(act,genl[nr],z);
if r<>fail then
if nr=Length(genl) then
Add(found,r);
if onlyone then return found; fi;
else
r2:=conjugate(Centralizer(act,z),Centralizer(asub,z),
List(genl,x->x^r),nr+1);
if Length(r2)>0 then
Append(found,r*r2);
if onlyone then return found; fi;
fi;
fi;
fi;
od;
od;
fi;
return found;
end;
if onlyone and IsSubset(a,b) then return One(g);fi;
# match points of perm group
if IsPermGroup(g) and IsSubset(g,a) and IsSubset(g,b) then
# how can we map orbits into orbits?
cla:=List(Orbits(a,MovedPoints(g)),Set);
clb:=List(Orbits(b,MovedPoints(g)),Set);
# no improvement if all orbits of a are fixed
if ForAny(cla,x->ForAny(GeneratorsOfGroup(g),y->OnSets(x,y)<>x)) then
r:=AllSubsetSummations(List(cla,Length),List(clb,Length),10^5);
if r=fail then
Info(InfoCoset,1,"Too many subset combinations");
else
Info(InfoCoset,1,"Testing ",Length(r)," combinations");
dc:=[];
for i in r do
k:=List(i,x->Union(clb{x}));
k:=RepresentativeAction(g,k,cla,OnTuplesSets);
if k<>fail then
Add(dc,[i,k]);
fi;
od;
if Length(dc)>0 then g:=Stabilizer(g,cla,OnTuplesSets);fi;
rep:=[];
for i in dc do
r:=DoConjugateInto(g,a,b^i[2],onlyone);
if onlyone then
if r<>fail then return i[2]*r;fi;
else
if r<>fail then Append(rep,List(r,x->i[2]*x));fi;
fi;
od;
if onlyone then return fail; #otherwise would have found and stopped
else return rep;fi;
fi;
else
# orbits are fixed. Make sure b is so
if ForAny(clb,x->not ForAny(cla,y->IsSubset(y,x))) then
if onlyone then return fail;else return [];fi;
fi;
fi;
fi;
# don't try the `MorGen...` search for more than two generators if
# generator number seems OK
if Length(SmallGeneratingSet(b))=AbelianRank(b) and
Length(SmallGeneratingSet(b))>2 then
gens:=SmallGeneratingSet(b);
elif IsPermGroup(b) and Size(b)<RootInt(NrMovedPoints(b)^3,2) then
r:=SmallerDegreePermutationRepresentation(b:cheap);
k:=Image(r,b);
gens:=MorFindGeneratingSystem(k,MorMaxFusClasses(MorRatClasses(k)));
gens:=List(gens,x->PreImagesRepresentative(r,x));
else
gens:=MorFindGeneratingSystem(b,MorMaxFusClasses(MorRatClasses(b)));
fi;
clb:=ConjugacyClasses(a);
cla:=[];
r:=[];
for i in gens do
b2:=Centralizer(g,i);
bd:=Size(Centralizer(b,i));
k:=Order(i);
rep:=[];
for j in [1..Length(clb)] do
if Order(Representative(clb[j]))=k
and (Size(a)/Size(clb[j])) mod bd=0 then
if not IsBound(r[j]) then
r[j]:=Size(Centralizer(g,Representative(clb[j])));
fi;
if r[j]=Size(b2) then
Add(rep,clb[j]);
fi;
fi;
od;
if Length(rep)=0 then
return []; # cannot have any
fi;
Add(cla,rep);
od;
r:=List(cla,x->-Maximum(List(x,Size)));
r:=Sortex(r);
gens:=Permuted(gens,r);
cla:=Permuted(cla,r);
r:=conjugate(g,a,gens,1);
if onlyone then
# get one
if Length(r)=0 then
return fail;
else
return r[1];
fi;
fi;
Info(InfoCoset,2,"Found ",Length(r)," reps");
# remove duplicate groups
rep:=[];
b2:=[];
for i in r do
bd:=b^i;
if ForAll(b2,x->RepresentativeAction(a,x,bd)=fail) then
Add(b2,bd);
Add(rep,i);
fi;
od;
return rep;
end);
#############################################################################
##
## IntermediateGroup(<G>,<U>) . . . . . . . . . subgroup of G containing U
##
## This routine tries to find a subgroup E of G, such that G>E>U. If U is
## maximal, it returns fail. This is done by using the maximal subgroups machinery or
## finding minimal blocks for
## the operation of G on the Right Cosets of U.
##
InstallGlobalFunction( IntermediateGroup, function(G,U)
local o,b,img,G1,c,m,hardlimit,gens,t,k,intersize;
if U=G then
return fail;
fi;
intersize:=Size(G);
m:=ValueOption("intersize");
if IsInt(m) and m<=intersize then
return fail; # avoid infinite recursion
fi;
# use maximals, use `Try` as we call with limiting options
IsNaturalAlternatingGroup(G);
IsNaturalSymmetricGroup(G);
if ValueOption("usemaximals")<>false then
m:=TryMaximalSubgroupClassReps(G:cheap,intersize:=intersize,nolattice);
if m<>fail and Length(m)>0 then
m:=Filtered(m,x->Size(x) mod Size(U)=0 and Size(x)>Size(U));
SortBy(m,x->Size(G)/Size(x));
gens:=SmallGeneratingSet(U);
for c in m do
if Index(G,c)<50000 then
t:=RightTransversal(G,c:noascendingchain); # conjugates
for k in t do
if ForAll(gens,x->k*x/k in c) then
Info(InfoCoset,2,"Found Size ",Size(c));
# U is contained in c^k
return c^k;
fi;
od;
else
t:=DoConjugateInto(G,c,U,true:intersize:=intersize,onlyone:=true);
if t<>fail and t<>[] then
Info(InfoCoset,2,"Found Size ",Size(c));
return c^(Inverse(t));
fi;
fi;
od;
Info(InfoCoset,2,"Found no intermediate subgroup ",Size(G)," ",Size(U));
return fail;
fi;
fi;
c:=ValueOption("refineChainActionLimit");
if IsInt(c) then
hardlimit:=c;
else
hardlimit:=1000000;
fi;
if Index(G,U)>hardlimit/10
and ValueOption("callinintermediategroup")<>true then
# try the `AscendingChain` mechanism
c:=AscendingChain(G,U:cheap,refineIndex:=QuoInt(IndexNC(G,U),2),
callinintermediategroup);
if Length(c)>2 then
return First(c,x->Size(x)>Size(U));
fi;
fi;
if ValueOption("cheap")=true then
return fail; # do not do hard work
fi;
if Index(G,U)>hardlimit then
Info(InfoWarning,1,
"will have to use permutation action of degree bigger than ", hardlimit);
fi;
# old code -- obsolete
if IsPermGroup(G) and Length(GeneratorsOfGroup(G))>3 then
G1:=Group(SmallGeneratingSet(G));
if HasSize(G) then
SetSize(G1,Size(G));
fi;
G:=G1;
fi;
o:=ActionHomomorphism(G,RightTransversal(G,U:noascendingchain),
OnRight,"surjective");
img:=Range(o);
b:=Blocks(img,MovedPoints(img));
if Length(b)=1 then
return fail;
else
b:=StabilizerOfBlockNC(img,First(b,i->1 in i));
b:=PreImage(o,b);
return b;
fi;
end );
#############################################################################
##
#F RefinedChain(<G>,<c>) . . . . . . . . . . . . . . . . refine chain links
##
InstallGlobalFunction(RefinedChain,function(G,cc)
local bound,a,b,c,cnt,r,i,j,bb,normalStep,gens,cheap,olda;
bound:=(10*LogInt(Size(G),10)+1)*Maximum(Factors(Size(G)));
bound:=Minimum(bound,20000);
cheap:=ValueOption("cheap")=true;
c:=ValueOption("refineIndex");
if IsInt(c) then
bound:=c;
fi;
c:=[];
for i in [2..Length(cc)] do
Add(c,cc[i-1]);
if Index(cc[i],cc[i-1]) > bound then
a:=AsSubgroup(Parent(cc[i]),cc[i-1]);
olda:=TrivialSubgroup(a);
while Index(cc[i],a)>bound and Size(a)>Size(olda) do
olda:=a;
# try extension via normalizer
b:=Normalizer(cc[i],a);
if Size(b)>Size(a) then
# extension by normalizer surely is a normal step
normalStep:=true;
bb:=b;
else
bb:=cc[i];
normalStep:=false;
b:=Centralizer(cc[i],Centre(a));
fi;
if Size(b)=Size(a) or Index(b,a)>bound then
cnt:=8+2^(LogInt(Index(bb,a),9));
#if cheap then cnt:=Minimum(cnt,50);fi;
cnt:=Minimum(cnt,40); # as we have better intermediate
repeat
if cnt<20 and not cheap then
# if random failed: do hard work
b:=IntermediateGroup(bb,a);
if b=fail then
b:=bb;
fi;
cnt:=0;
else
# larger indices may take more tests...
Info(InfoCoset,5,"Random");
repeat
r:=Random(bb);
until not(r in a);
if normalStep then
# NC is safe
b:=ClosureSubgroupNC(a,r);
else
# self normalizing subgroup: thus every element not in <a>
# will surely map one generator out
j:=0;
gens:=GeneratorsOfGroup(a);
repeat
j:=j+1;
until not(gens[j]^r in a);
r:=gens[j]^r;
# NC is safe
b:=ClosureSubgroupNC(a,r);
fi;
if Size(b)<Size(bb) then
Info(InfoCoset,1,"improvement found ",Size(bb)/Size(b));
bb:=b;
fi;
cnt:=cnt-1;
fi;
until Index(bb,a)<=bound or cnt<1;
fi;
if Index(b,a)>bound and Length(c)>1 then
bb:=IntermediateGroup(b,c[Length(c)-1]);
if bb<>fail and Size(bb)>Size(Last(c)) then
c:=Concatenation(c{[1..Length(c)-1]},[bb],Filtered(cc,x->Size(x)>=Size(b)));
return RefinedChain(G,c);
fi;
fi;
a:=b;
if a<>cc[i] then #not upper level
Add(c,a);
fi;
od;
fi;
od;
Add(c,Last(cc));
a:=Last(c);
for i in [Length(c)-1,Length(c)-2..1] do
#enforce parent relations
if not HasParent(c[i]) then
SetParent(c[i],a);
a:=c[i];
else
a:=AsSubgroup(a,c[i]);
c[i]:=a;
fi;
od;
return c;
end);
InstallMethod( AscendingChainOp, "generic", IsIdenticalObj, [IsGroup,IsGroup],0,
function(G,U)
return RefinedChain(G,[U,G]);
end);
InstallMethod(DoubleCoset,"generic",IsCollsElmsColls,
[IsGroup,IsObject,IsGroup],0,
function(U,g,V)
local d,fam;
fam:=FamilyObj(U);
if not IsBound(fam!.doubleCosetsDefaultType) then
fam!.doubleCosetsDefaultType:=NewType(fam,IsDoubleCosetDefaultRep
and HasLeftActingGroup and HasRightActingGroup
and HasRepresentative);
fi;
d:=rec();
ObjectifyWithAttributes(d,fam!.doubleCosetsDefaultType,
LeftActingGroup,U,RightActingGroup,V,Representative,g);
return d;
end);
InstallOtherMethod(DoubleCoset,"with size",true,
[IsGroup,IsObject,IsGroup,IsPosInt],0,
function(U,g,V,sz)
local d,fam,typ;
fam:=FamilyObj(U);
typ:=NewType(fam,IsDoubleCosetDefaultRep
and HasIsFinite and IsFinite
and HasLeftActingGroup and HasRightActingGroup
and HasRepresentative);
d:=rec();
ObjectifyWithAttributes(d,typ,
LeftActingGroup,U,RightActingGroup,V,Representative,g);
SetSize(d,sz); # Size has private setter which will cause problems with
# HasSize triggering an immediate method.
return d;
end);
InstallMethod(\=,"DoubleCosets",IsIdenticalObj,[IsDoubleCoset,IsDoubleCoset],0,
function(a,b)
if LeftActingGroup(a)<>LeftActingGroup(b) or
RightActingGroup(a)<>RightActingGroup(b) then
return false;
fi;
# avoid forcing RepresentativesContainedRightCosets on both if one has
if HasRepresentativesContainedRightCosets(b) then
if HasRepresentativesContainedRightCosets(a) then
return RepresentativesContainedRightCosets(a)
=RepresentativesContainedRightCosets(b);
else
return CanonicalRightCosetElement(LeftActingGroup(a),
Representative(a)) in
RepresentativesContainedRightCosets(b);
fi;
else
return CanonicalRightCosetElement(LeftActingGroup(b),
Representative(b)) in
RepresentativesContainedRightCosets(a);
fi;
end);
InstallMethod(ViewString,"DoubleCoset",true,[IsDoubleCoset],0,
function(d)
return(STRINGIFY("DoubleCoset(\<",
ViewString(LeftActingGroup(d)),",\>",
ViewString(Representative(d)),",\>",
ViewString(RightActingGroup(d)),"\<)"));
end);
InstallMethodWithRandomSource(Random,
"for a random source and a double coset",
[IsRandomSource, IsDoubleCoset],0,
function(rs, d)
return Random(rs,LeftActingGroup(d))*Representative(d)
*Random(rs,RightActingGroup(d));
end);
InstallMethod(PseudoRandom,"double coset",true,[IsDoubleCoset],0,
function(d)
return PseudoRandom(LeftActingGroup(d))*Representative(d)
*PseudoRandom(RightActingGroup(d));
end);
InstallMethod(RepresentativesContainedRightCosets,"generic",true,
[IsDoubleCoset],0,
function(c)
local u,v,o,i,j,img;
u:=LeftActingGroup(c);
v:=RightActingGroup(c);
o:=[CanonicalRightCosetElement(u,Representative(c))];
# orbit alg.
for i in o do
for j in GeneratorsOfGroup(v) do
img:=CanonicalRightCosetElement(u,i*j);
if not img in o then
Add(o,img);
fi;
od;
od;
return Set(o);
end);
InstallMethod(\in,"double coset",IsElmsColls,
[IsMultiplicativeElementWithInverse,IsDoubleCoset],0,
function(e,d)
return CanonicalRightCosetElement(LeftActingGroup(d),e)
in RepresentativesContainedRightCosets(d);
end);
InstallMethod(Size,"double coset",true,[IsDoubleCoset],0,
function(d)
return
Size(LeftActingGroup(d))*Length(RepresentativesContainedRightCosets(d));
end);
InstallMethod(AsList,"double coset",true,[IsDoubleCoset],0,
function(d)
local l;
l:=Union(List(RepresentativesContainedRightCosets(d),
i->RightCoset(LeftActingGroup(d),i)));
return l;
end);
#############################################################################
##
#M Enumerator
##
BindGlobal( "ElementNumber_DoubleCoset", function( enum, pos )
pos:= pos-1;
return enum!.leftgroupEnumerator[ ( pos mod enum!.leftsize )+1 ]
* enum!.rightCosetReps[ QuoInt( pos, enum!.leftsize )+1 ];
end );
BindGlobal( "NumberElement_DoubleCoset", function( enum, elm )
local p;
p:= First( [ 1 .. Length( enum!.rightCosetReps ) ],
i -> elm / enum!.rightCosetReps[i] in enum!.leftgroup );
p:= (p-1) * enum!.leftsize
+ Position( enum!.leftgroupEnumerator,
elm / enum!.rightCosetReps[p], 0 );
return p;
end );
InstallMethod( Enumerator,
"for a double coset",
[ IsDoubleCoset ],
d -> EnumeratorByFunctions( d, rec(
NumberElement := NumberElement_DoubleCoset,
ElementNumber := ElementNumber_DoubleCoset,
leftgroupEnumerator := Enumerator( LeftActingGroup( d ) ),
leftgroup := LeftActingGroup( d ),
leftsize := Size( LeftActingGroup( d ) ),
rightCosetReps := RepresentativesContainedRightCosets( d ) ) ) );
RightCosetCanonicalRepresentativeDeterminator :=
function(U,a)
return [CanonicalRightCosetElement(U,a)];
end;
InstallMethod(RightCoset,"generic",IsCollsElms,
[IsGroup,IsObject],0,
function(U,g)
local d,fam;
# noch tests...
fam:=FamilyObj(U);
if not IsBound(fam!.rightCosetsDefaultType) then
fam!.rightCosetsDefaultType:=NewType(fam,IsRightCosetDefaultRep and
HasActingDomain and HasFunctionAction and HasRepresentative and
HasCanonicalRepresentativeDeterminatorOfExternalSet);
fi;
d:=rec();
ObjectifyWithAttributes(d,fam!.rightCosetsDefaultType,
ActingDomain,U,FunctionAction,OnLeftInverse,Representative,g,
CanonicalRepresentativeDeterminatorOfExternalSet,
RightCosetCanonicalRepresentativeDeterminator);
if HasSize(U) then
# We cannot set the size in the previous ObjectifyWithAttributes as there
# is a custom setter method. In such a case ObjectifyWithAttributes just
# does `Objectify` and calls all setters separately which is what we want
# to avoid here.
SetSize(d,Size(U));
fi;
return d;
end);
InstallOtherMethod(\*,"group times element",IsCollsElms,
[IsGroup,IsMultiplicativeElementWithInverse],0,
function(s,a)
return RightCoset(s,a);
end);
InstallMethod(ViewString,"RightCoset",true,[IsRightCoset],0,
function(d)
return STRINGIFY("RightCoset(\<",
ViewString(ActingDomain(d)),",\>",
ViewString(Representative(d)),")");
end);
InstallMethod(PrintString,"RightCoset",true,[IsRightCoset],0,
function(d)
return STRINGIFY("RightCoset(\<",
PrintString(ActingDomain(d)),",\>",
PrintString(Representative(d)),")");
end);
InstallMethod(PrintObj,"RightCoset",true,[IsRightCoset],0,
function(d)
Print(PrintString(d));
end);
InstallMethod(ViewObj,"RightCoset",true,[IsRightCoset],0,
function(d)
Print(ViewString(d));
end);
InstallMethod(IsBiCoset,"test property",true,[IsRightCoset],0,
function(c)
local s,r;
s:=ActingDomain(c);
r:=Representative(c);
return ForAll(GeneratorsOfGroup(s),x->x^r in s);
end);
InstallMethodWithRandomSource(Random,
"for a random source and a RightCoset",
[IsRandomSource, IsRightCoset],0,
function(rs, d)
return Random(rs, ActingDomain(d))*Representative(d);
end);
InstallMethod(PseudoRandom,"RightCoset",true,[IsRightCoset],0,
function(d)
return PseudoRandom(ActingDomain(d))*Representative(d);
end);
InstallMethod(\=,"RightCosets",IsIdenticalObj,[IsRightCoset,IsRightCoset],0,
function(a,b)
return ActingDomain(a)=ActingDomain(b) and
Representative(a)/Representative(b) in ActingDomain(a);
end);
InstallOtherMethod(\*,"RightCoset with element",IsCollsElms,
[IsRightCoset,IsMultiplicativeElementWithInverse],0,
function(a,g)
return RightCoset( ActingDomain( a ), Representative( a ) * g );
end);
InstallOtherMethod(\*,"RightCosets",IsIdenticalObj,
[IsRightCoset,IsRightCoset],0,
function(a,b)
local c;
if ActingDomain(a)<>ActingDomain(b) then TryNextMethod();fi;
if not IsBiCoset(a) then # product does not require b to be bicoset
ErrorNoReturn("right cosets can only be multiplied if the left operand is a bicoset");
fi;
c:=RightCoset(ActingDomain(a), Representative(a) * Representative(b) );
if HasIsBiCoset(b) then
SetIsBiCoset(c,IsBiCoset(b));
fi;
return c;
end);
InstallOtherMethod(InverseOp,"Right cosets",true,
[IsRightCoset],0,
function(a)
local s,r;
s:=ActingDomain(a);
r:=Representative(a);
if not IsBiCoset(a) then
ErrorNoReturn("only right cosets which are bicosets can be inverted");
fi;
r:=RightCoset(s,Inverse(r));
SetIsBiCoset(r,true);
return r;
end);
InstallOtherMethod(OneOp,"Right cosets",true,
[IsRightCoset],0,
function(a)
return RightCoset(ActingDomain(a),One(Representative(a)));
end);
InstallMethod(IsGeneratorsOfMagmaWithInverses,"cosets",true,
[IsMultiplicativeElementWithInverseCollColl],0,
function(l)
local a,r;
if Length(l)>0 and ForAll(l,IsRightCoset) then
a:=ActingDomain(l[1]);
r:=List(l,Representative);
if ForAll(l,x->ActingDomain(x)=a) and
ForAll(r,x->ForAll(GeneratorsOfGroup(a),y->y^x in a)) then
return true;
fi;
fi;
TryNextMethod();
end);
InstallMethod(Intersection2, "general cosets", IsIdenticalObj,
[IsRightCoset,IsRightCoset],
function(cos1,cos2)
local swap, H1, H2, x1, x2, sigma, U, rho;
if Size(cos1) < 10 then
TryNextMethod();
elif Size(cos2) < 10 then
return Intersection2(cos2, cos1);
fi;
if Size(cos1) > Size(cos2) then
swap := cos1;
cos1 := cos2;
cos2 := swap;
fi;
H1:=ActingDomain(cos1);
H2:=ActingDomain(cos2);
x1:=Representative(cos1);
x2:=Representative(cos2);
sigma := x1 / x2;
if Size(H1) = Size(H2) and H1 = H2 then
if sigma in H1 then
return cos1;
else
return [];
fi;
fi;
# We want to compute the intersection of cos1 = H1*x1 with cos2 = H2*x2.
# This is equivalent to intersecting H1 with H2*x2/x1, which is either empty
# or equal to a coset U*rho, where U is the intersection of H1 and H2.
# In the non-empty case, the overall result then is U*rho*x1.
#
# To find U*rho, we iterate over all cosets of U in H1 and for each test
# if it is contained in H2*x2/x1, which is the case if and only if rho is
# in H2*x2/x1, if and only if rho/(x2/x1) = rho*x1/x2 is in H2
U:=Intersection(H1, H2);
for rho in RightTransversal(H1, U) do
if rho * sigma in H2 then
return RightCoset(U, rho * x1);
fi;
od;
return [];
end);
# disabled because of comparison incompatibilities
#InstallMethod(\<,"RightCosets",IsIdenticalObj,[IsRightCoset,IsRightCoset],0,
#function(a,b)
# # this comparison is *NOT* necessarily equivalent to a comparison of the
# # element lists!
# if ActingDomain(a)<>ActingDomain(b) then
# return ActingDomain(a)<ActingDomain(b);
# fi;
# return CanonicalRepresentativeOfExternalSet(a)
# <CanonicalRepresentativeOfExternalSet(b);
#end);
InstallGlobalFunction( DoubleCosets, function(G,U,V)
if not (IsSubset(G,U) and IsSubset(G,V)) then
Error("not contained");
fi;
return DoubleCosetsNC(G,U,V);
end );
InstallGlobalFunction( RightCosets, function(G,U)
if not IsSubset(G,U) then
Error("not contained");
fi;
return RightCosetsNC(G,U);
end );
InstallMethod(CanonicalRightCosetElement,"generic",IsCollsElms,
[IsGroup,IsObject],0,
function(U,e)
local l;
l:=List(AsList(U),i->i*e);
return Minimum(l);
end);
# TODO: In the long run this should become a more general operation,
# but for the moment it is specialized for the application at hand
BindGlobal("DCFuseSubgroupOrbits",function(P,G,reps,act,lim,count)
local live,orbs,orbset,done,nr,p,o,os,orbitextender,bahn,i,j,enum,dict,map,pam;
# return positive fuse number or negative position how far it got
orbitextender:=function(o,os,start,limit,this)
local i,gen,img,e;
i:=start;
while i<=Length(o) and Length(o)<limit do
for gen in GeneratorsOfGroup(G) do
img:=act(o[i],gen);
e:=Position(enum,img);
if not e in os then # duplicate? Still use os as we need to grow o
#p:=PositionProperty(orbset,x->e in x);
p:=LookupDictionary(dict,e);
if p<>fail and
# we could have found an element that we know (because of
# fusion) to be already in this orbit (but must store)
pam[map[p]]<>this then
p:=map[p]; #retrieved image position might have been fused away
return p;
fi;
Add(o,img);
AddSet(os,e);
if p=fail then AddDictionary(dict,e,this);fi;
fi;
od;
i:=i+1;
od;
#if i>Length(o) and Length(os)>Length(o) then Error("ran out of orbit");fi;
return -(i-1);
end;
bahn:=[];
enum:=Enumerator(P);
live:=[];
orbs:=[];
orbset:=[];
done:=[];
dict:=NewDictionary(1,true,rec(hashfun:=x->x));
map:=[];
pam:=[]; # reverse of map
for nr in [1..Length(reps)] do
#p:=PositionProperty(orbset,x->Position(enum,reps[nr]) in x);
p:=LookupDictionary(dict,Position(enum,reps[nr]));
if p=fail then
# start orbit algorithm
o:=[reps[nr]];
os:=[Position(enum,reps[nr])];
AddDictionary(dict,os[1],nr);
p:=orbitextender(o,os,1,lim,nr);
if p<0 then
# new orbit
Info(InfoCoset,4,nr," lives");
Add(live,nr);
Add(orbs,o);
Add(orbset,os);
Add(done,-p);
Add(bahn,[nr]);
map[nr]:=Length(orbs);
pam[Length(orbs)]:=nr;
i:=1;
while Length(orbs)>count do
# one orbit too many
if ForAll(orbs,x->Length(x)>=lim) then
if lim<20000 then
lim:=lim*2;
else
lim:=(QuoInt(lim,8000)+1)*8000;
fi;
fi;
Info(InfoCoset,4,"Redo ",i," ",lim);
p:=orbitextender(orbs[i],orbset[i],done[i],lim,pam[i]);
if p>0 then
Info(InfoCoset,4,"Join ",i," to ",p);
if p=i then Error("selfjoin cannot happen");fi;
bahn[p]:=Union(bahn[p],bahn[i]);
#UniteSet(orbset[p],orbset[i]);
for j in [1..Length(map)] do
if IsBound(map[j]) and map[j]=i then map[j]:=p; fi;
od;
# delete entry i, move higher ones one up
for j in [1..Length(map)] do
if IsBound(map[j]) and map[j]>i then map[j]:=map[j]-1; fi;
od;
# Remove entry i
Remove(orbs,i);
Remove(orbset,i);
Remove(done,i);
Remove(bahn,i);
Remove(pam,i);
else
done[i]:=-p;
fi;
i:=i+1; if i>Length(orbs) then i:=1;fi;
od;
else
Info(InfoCoset,4,nr," fuses into ",p," @",Length(os));
map[nr]:=p; # and indeed nr itself maps to p
AddSet(bahn[p],nr);
#UniteSet(orbset[p],os);
# not needed
#for j in os do AddDictionary(dict,j,p); od;
fi;
else
p:=map[p]; #retrieved image position might have been fused away
Info(InfoCoset,4,nr," lies in ",p);
map[nr]:=p;
AddSet(bahn[p],nr);
#AddDictionary(dict,Position(enum,reps[nr]),p);
fi;
od;
return bahn;
end);
#############################################################################
##
#F CalcDoubleCosets( <G>, <A>, <B> ) . . . . . . . . . double cosets: A\G/B
##
## DoubleCosets routine using an
## ascending chain of subgroups from A to G, using the fact, that a
## double coset is an union of right cosets
##
InstallGlobalFunction(CalcDoubleCosets,function(G,a,b)
local c, flip, maxidx, cano, tryfct, p, r, t,
stabs, dcs, homs, tra, a1, a2, indx, normal, hom, omi, omiz,c1,
unten, compst, s, nr, nstab, lst, sifa, pinv, blist, bsz, cnt,
ps, e, mop, mo, lstgens, lstgensop, rep, st, o, oi, i, img, ep,
siz, rt, j, canrep,step,nu,doneidx,orbcnt,posi,
sizes,cluster,sel,lr,lstabs,ssizes,num,actfun,mayflip,rs,
actlimit, uplimit, badlimit,avoidlimit,start,includestab,quot;
actlimit:=300000; # maximal degree on which we try blocks
uplimit:=500000; # maximal index for up step
avoidlimit:=200000; # beyond this index we want to get smaller
badlimit:=5000000; # beyond this index things might break down
mayflip:=true; # are we allowed to flip for better chain as well?
# Do we *want* stabilizers
includestab:=ValueOption("includestab")=true;
# if a is small and b large, compute cosets b\G/a and take inverses of the
# representatives: Since we compute stabilizers in b and a chain down to
# a, this is notably faster
if ValueOption("noflip")<>true and 3*Size(a)<2*Size(b) then
c:=b;
b:=a;
a:=c;
flip:=true;
Info(InfoCoset,1,"DoubleCosetFlip");
else
flip:=false;
fi;
if Index(G,a)=1 then
return [[One(G),Size(G)]];
fi;
# maximal index of a series
maxidx:=function(ser)
return Maximum(List([1..Length(ser)-1],x->Size(ser[x+1])/Size(ser[x])));
end;
# compute ascending chain and refine if necessarily (we anyhow need action
# on cosets).
#c:=AscendingChain(G,a:refineChainActionLimit:=Index(G,a));
c:=AscendingChain(G,a:refineChainActionLimit:=actlimit,indoublecoset);
# do we first go into a factor group?
quot:=ValueOption("usequotient");
PushOptions(rec(usequotient:=fail));# not to be used within itself
if not IsBool(quot) then
if IsMapping(quot) then
a1:=KernelOfMultiplicativeGeneralMapping(quot);
else
a1:=quot;
quot:=NaturalHomomorphismByNormalSubgroupNC(G,quot);
fi;
r:=RestrictedMapping(quot,b);
a2:=ClosureGroup(a1,a);
Size(a2);
start:=PositionProperty(c,
x->Size(x)=Size(a2) and ForAll(GeneratorsOfGroup(x),y->y in a2));
if start=fail then Error("closure not in chain");fi;
p:=Image(quot,G);
c1:=Image(quot,a);
tra:=Image(quot,b);
dcs:=CalcDoubleCosets(p,c1,tra:includestab,usequotient:=fail);
for i in dcs do
# add missing stabilizers (caused by flip)
if not IsBound(i[3]) then
i[3]:=Intersection(c1^i[1],tra);
fi;
od;
mayflip:=false;
Info(InfoCoset,1,"Factor returns ",Length(dcs)," double cosets");
# try kernel
a2:=Filtered(GeneratorsOfGroup(b),x->IsOne(ImagesRepresentative(quot,x)));
a2:=SubgroupNC(Parent(b),a2);
Assert(2,Size(a2)*Size(tra)=Size(b));
SetKernelOfMultiplicativeGeneralMapping(r,a2);
dcs:=List(dcs,x->[PreImagesRepresentative(quot,x[1]),Size(a1)*x[2],
PreImage(r,x[3])]);
r:=List(dcs,x->x[1]);
stabs:=List(dcs,x->x[3]);
else
start:=1;
r:=[One(G)];
stabs:=[b];
quot:=fail;
fi;
# cano indicates whether there is a final up step (and thus we need to
# form canonical representatives). ```Canonical'' means that on each
# transversal level the orbit representative is chosen to be minimal (in
# the transversal position).
cano:=false;
doneidx:=[]; # indices done already -- avoid duplicate
if maxidx(c)>avoidlimit and mayflip then
# try to do better
# what about flipping (back)?
c1:=AscendingChain(G,b:refineChainActionLimit:=actlimit,indoublecoset);
if maxidx(c1)<=avoidlimit then
Info(InfoCoset,1,"flip to get better chain");
c:=b;
b:=a;
a:=c;
flip:=not flip;
c:=c1;
stabs:=[b]; # make sure stabs also flips over
elif IsPermGroup(G) then
actlimit:=Maximum(actlimit,NrMovedPoints(G));
avoidlimit:=Maximum(avoidlimit,NrMovedPoints(G));
tryfct:=function(obj,act)
local G1,a1,c1;
if IsList(act) and Length(act)=2 then
G1:=act[1];
a1:=act[2];
else
#Print(maxidx(c),obj,Length(Orbit(G,obj,act))," ",
# Length(Orbit(a,obj,act)),"\n");
G1:=Stabilizer(G,obj,act);
if Index(G,G1)<maxidx(c) then
a1:=Stabilizer(a,obj,act);
else
a1:=G;
fi;
fi;
Info(InfoCoset,4,"attempt up step ",obj," index:",Size(a)/Size(a1));
if Index(G,G1)<maxidx(c) and Index(a,a1)<=uplimit and (
maxidx(c)>avoidlimit or Size(a1)>Size(c[1])) then
c1:=AscendingChain(G1,a1:refineIndex:=avoidlimit,
refineChainActionLimit:=actlimit,
indoublecoset);
if maxidx(c1)<maxidx(c) then
c:=Concatenation(c1,[G]);
cano:=true;
Info(InfoCoset,1,"improved chain with up step ",obj,
" index:",Size(a)/Size(a1)," maxidx=",maxidx(c));
fi;
fi;
end;
rs:=Filtered(TryMaximalSubgroupClassReps(G:cheap),
x->Index(G,x)<=5*avoidlimit);
SortBy(rs,a->-Size(a));
for i in rs do
if Index(G,i)<maxidx(c) then
p:=Intersection(a,i);
AddSet(doneidx,Index(a,p));
if Index(a,p)<=uplimit then
Info(InfoCoset,3,"Try maximal of Indices ",Index(G,i),":",
Index(a,p));
tryfct("max",[i,p]);
fi;
fi;
od;
p:=LargestMovedPoint(a);
tryfct(p,OnPoints);
for i in Orbits(Stabilizer(a,p),Difference(MovedPoints(a),[p])) do
tryfct(Set([i[1],p]),OnSets);
od;
fi;
if maxidx(c)>badlimit then
rs:=ShallowCopy(TryMaximalSubgroupClassReps(a:cheap));
rs:=Filtered(rs,x->Index(a,x)<uplimit and not Index(a,x) in doneidx);
SortBy(rs,a->-Size(a));
for j in rs do
#Print("j=",Size(j),"\n");
t:=AscendingChain(G,j:refineIndex:=avoidlimit,
refineChainActionLimit:=actlimit,indoublecoset);
Info(InfoCoset,4,"maxidx ",Index(a,j)," yields ",maxidx(t),": ",
List(t,Size));
if maxidx(t)<maxidx(c) and (maxidx(c)>badlimit or
# only increase up-step if index gets better by extra index
(maxidx(c)>maxidx(t)*Size(c[1])/Size(t[1])) ) then
c:=t;
cano:=true;
Info(InfoCoset,1,"improved chain with up step index:",
Size(a)/Size(j));
fi;
od;
fi;
elif ValueOption("sisyphus")=true then
# purely to allow for tests of up-step mechanism in smaller examples.
# This is creating unnecessary extra work and thus should never be used
# in practice, but will force some code to be run through.
c:=Concatenation([TrivialSubgroup(G)],c);
cano:=true;
fi;
dcs:=[];
# Do we want to keep result for a smaller group (as cheaper fuse is possible
# outside function at a later stage)?
if ValueOption("noupfuse")=true then cano:=false;fi;
Info(InfoCoset,1,"Chosen series is ",List(c,Size));
#if ValueOption("indoublecoset")<>true then Error("GNASH");fi;
# calculate setup for once
homs:=[];
tra:=[];
for step in [start..Length(c)-1] do
a1:=c[Length(c)-step+1];
a2:=c[Length(c)-step];
indx:=Index(a1,a2);
normal:=IsNormal(a1,a2);
# don't try to refine again for transversal, we've done so already.
t:=RightTransversal(a1,a2:noascendingchain);
tra[step]:=t;
# is it worth using a permutation representation?
if (step>1 or cano) and Length(t)<badlimit and IsPermGroup(G) and
not normal then
# in this case, we can beneficially compute the action once and then use
# homomorphism methods to obtain the permutation image
Info(InfoCoset,2,"using perm action on step ",step,": ",Length(t));
hom:=Subgroup(G,SmallGeneratingSet(a1));
hom:=ActionHomomorphism(hom,t,OnRight,"surjective");
else
hom:=fail;
fi;
homs[step]:=hom;
od;
omi:=[];
omiz:=[];
for step in [start..Length(c)-1] do
a1:=c[Length(c)-step+1];
a2:=c[Length(c)-step];
normal:=IsNormal(a1,a2);
indx:=Index(a1,a2);
if normal then
Info(InfoCoset,1,"Normal Step :",indx,": ",Length(r)," double cosets");
else
Info(InfoCoset,1,"Step :",indx,": ",Length(r)," double cosets");
fi;
# is this the last step?
unten:=step=Length(c)-1 and cano=false;
# shall we compute stabilizers?
compst:=(not unten) or normal or includestab;
t:=tra[step];
hom:=homs[step];
s:=[];
nr:=[];
nstab:=[];
for nu in [1..Length(r)] do
lst:=stabs[nu];
Info(InfoCoset,4,"number ",nu,", |stab|=",Size(lst));
sifa:=Size(a2)*Size(b)/Size(lst);
p:=r[nu];
pinv:=p^-1;
blist:=BlistList([1..indx],[]);
bsz:=indx;
orbcnt:=0;
# if a2 is normal in a1, the stabilizer is the same for all Orbits of
# right cosets. Thus we need to compute only one, and will receive all
# others by simple calculations afterwards
if normal then
cnt:=1;
else
cnt:=indx;
fi;
if cano=false and indx>20 and IsSolvableGroup(lst) then
lstgens:=Pcgs(lst);
else
lstgens:=GeneratorsOfGroup(lst);
if Length(lstgens)>2 and Length(t)>100 then
lstgens:=SmallGeneratingSet(lst);
fi;
fi;
lstgensop:=List(lstgens,i->i^pinv); # conjugate generators: operation
# is on cosets a.p; we keep original cosets: Ua.p.g/p, this
# corresponds to conjugate operation
if hom<>fail then
lstgensop:=List(lstgensop,i->Image(hom,i));
actfun:=OnPoints;
else
actfun:=function(num,gen)
return PositionCanonical(t,t[num]*gen);
end;
fi;
posi:=0;
while bsz>0 and cnt>0 do
cnt:=cnt-1;
# compute orbit and stabilizers for the next step
# own Orbitalgorithm and stabilizer computation
#while blist[posi] do posi:=posi+1;od;
posi:=Position(blist,false,posi);
ps:=posi;
blist[ps]:=true;
bsz:=bsz-1;
e:=t[ps];
mop:=1;
mo:=ps;
rep := [ One(b) ];
o:=[ps];
if cano or compst then
oi:=[];
oi[ps]:=1; # reverse index
fi;
orbcnt:=orbcnt+1;
if cano=false and IsPcgs(lstgens) then
if compst then
o:=OrbitStabilizer(lst,o[1],lstgens,lstgensop,actfun);
st:=o.stabilizer;
o:=o.orbit;
else
o:=Orbit(lst,o[1],lstgens,lstgensop,actfun);
fi;
for i in o do
blist[i]:=true;
od;
bsz:=bsz-Length(o)+1;
else
if compst then
# stabilizing generators
st:=Filtered(GeneratorsOfGroup(lst),
x->PositionCanonical(r,t[ps]*x)=ps);
if Length(st)=Length(GeneratorsOfGroup(lst)) then
st:=lst; # immediate end -- orbit 1
else
st := SubgroupNC(lst,st);
fi;
else
st:=TrivialSubgroup(lst);
fi;
i:=1;
while i<=Length(o)
# will not grab if nonreg,. orbit and stabilizer not computed,
# but comparatively low cost and huge help if hom=fail
and Size(st)*Length(o)<Size(lst) do
for j in [1..Length(lstgens)] do
if hom=fail then
img:=t[o[i]]*lstgensop[j];
ps:=PositionCanonical(t,img);
else
ps:=o[i]^lstgensop[j];
fi;
if blist[ps] then
if compst then
# known image
#NC is safe (initializing as TrivialSubgroup(G)
st := ClosureSubgroupNC(st,rep[i]*lstgens[j]/rep[oi[ps]]);
fi;
else
# new image
blist[ps]:=true;
bsz:=bsz-1;
Add(o,ps);
if cano or compst then
Add(rep,rep[i]*lstgens[j]);
if cano and ps<mo then
mo:=ps;
mop:=Length(rep);
fi;
oi[ps]:=Length(o);
fi;
fi;
od;
i:=i+1;
od;
fi;
Info(InfoCoset,5,"|o|=",Length(o));
ep:=e*rep[mop]*p;
Add(nr,ep);
if compst then
st:=st^rep[mop];
Add(nstab,st);
fi;
if cano and step=1 and not normal then
Add(omi,mo);
Add(omiz,Length(o));
fi;
siz:=sifa*Length(o); #order
if unten then
if includestab then
if flip then
Add(dcs,[ep^(-1),siz]);
else
Add(dcs,[ep,siz,st]);
fi;
else
if flip then
Add(dcs,[ep^(-1),siz]);
else
Add(dcs,[ep,siz]);
fi;
fi;
fi;
od;
Info(InfoCoset,4,"Get ",orbcnt," orbits");
if normal then
# in the normal case, we can obtain the other orbits easily via
# the orbit theorem (same stabilizer)
if Size(lst)/Size(st)<10 then
# if the group `st` is handled by a nice monomorphism, the
# identity might not be the canonical element for the subgroup.
rt:=Orbit(lst,CanonicalRightCosetElement(st,One(st)),
function(rep,g) return CanonicalRightCosetElement(st,rep*g);end);
else
rt:=RightTransversal(lst,st:noascendingchain);
fi;
Assert(1,Length(rt)=Length(o));
while bsz>0 do
ps:=Position(blist,false);
e:=t[ps];
blist[ps]:=true;
ep:=e*p;
mo:=ep;
mop:=ps;
# tick off the orbit
for i in rt do
#ps:=PositionCanonical(t,e*p*i/p);
j:=ep*i/p;
ps:=PositionCanonical(t,ep*i/p);
if cano then
if ps<mop then
mop:=ps;
mo:=j;
fi;
fi;
blist[ps]:=true;
od;
bsz:=bsz-Length(rt);
Add(nr,mo);
Add(nstab,st);
if unten then
if includestab then
if flip then
Add(dcs,[ep^(-1),siz]);
else
Add(dcs,[ep,siz,st]);
fi;
else
if flip then
Add(dcs,[ep^(-1),siz]);
else
Add(dcs,[ep,siz]);
fi;
fi;
fi;
od;
fi;
od;
stabs:=nstab;
r:=nr;
Info(InfoCoset,3,Length(r)," double cosets so far.");
od;
if cano then
# do the final up step
IsSSortedList(omi);
# canonization fct
canrep:=function(x)
local stb, p, pinv, t, hom,ps, mop, mo, o, oi, rep, st, lstgens, lstgensop,
i, img, step, j,calcs;
stb:=b;
p:=One(G);
for step in [1..Length(c)-1] do
calcs:=step<Length(c)-1;
pinv:=p^-1;
t:=tra[step];
hom:=homs[step];
# orbit-stabilizer algorithm
ps:=PositionCanonical(t,x);
mop:=1;
mo:=ps;
o:=[ps];
oi:=[];
oi[ps]:=1;
rep:=[One(stb)];
st:=TrivialSubgroup(b);
lstgens:=GeneratorsOfGroup(stb);
if Length(lstgens)>4 and
Length(lstgens)/(AbelianRank(stb)+1)*2>5 then
lstgens:=SmallGeneratingSet(stb);
fi;
lstgensop:=List(lstgens,i->i^pinv); # conjugate generators: operation
if hom<>fail then
lstgensop:=List(lstgensop,i->Image(hom,i));
fi;
i:=1;
while i<=Length(o) do
for j in [1..Length(lstgensop)] do
if hom=fail then
img:=t[o[i]]*lstgensop[j];
ps:=PositionCanonical(t,img);
else
ps:=o[i]^lstgensop[j];
fi;
if IsBound(oi[ps]) then
# known image
# if there is only one orbit on the top step, we know the
# stabilizer!
if calcs then
#NC is safe (initializing as TrivialSubgroup(G)
st := ClosureSubgroupNC(st,rep[i]*lstgens[j]/rep[oi[ps]]);
if Size(st)*Length(o)=Size(b) then i:=Length(o);fi;
fi;
#fi;
else
Add(o,ps);
Add(rep,rep[i]*lstgens[j]);
if ps<mo then
mo:=ps;
mop:=Length(rep);
if step=1 and mo in omi then
#Print("found\n");
if Size(st)*omiz[Position(omi,mo)]=Size(stb) then
# we have the minimum and the right stabilizer: break
#Print("|Orbit|=",Length(o),
#" of ",omiz[Position(omi,mo)]," min=",mo,"\n");
i:=Length(o);
fi;
fi;
fi;
oi[ps]:=Length(o);
if Size(st)*Length(o)=Size(b) then i:=Length(o);fi;
fi;
od;
i:=i+1;
od;
if calcs then
stb:=st^(rep[mop]);
fi;
#if HasSmallGeneratingSet(st) then
# SetSmallGeneratingSet(stb,List(SmallGeneratingSet(st),x->x^rep[mop]));
#fi;
#else
# stb:=omis;
#fi;
x:=x*(rep[mop]^pinv)/t[mo];
p:=t[mo]*p;
#Print("step ",step," |Orbit|=",Length(o),"nmin=",mo,"\n");
#if ForAny(GeneratorsOfGroup(stb),
# i->not x*p*i/p in t!.subgroup) then
# Error("RRR");
#fi;
od;
return p;
end;
# now fuse orbits under the left action of a
indx:=Index(a,a2);
Info(InfoCoset,2,"fusion index ",indx);
#t:=Filtered(RightTransversal(a,a2),x->not x in a2);
t:=RightTransversal(a,a2);
sifa:=Size(a2)*Size(b);
# cluster according to A-double coset sizes and C lengths
#sizes:=List(r,x->Size(a)*Size(b)/Size(Intersection(b,a^x)));
hom:=ActionHomomorphism(a,t,OnRight,"surjective");
sizes:=[];
for i in [1..Length(r)] do
lr:=Intersection(a,b^(r[i]^-1));
# size of double coset and
Add(sizes,[Size(a)*Size(b)/Size(lr),
Length(OrbitsDomain(Image(hom,lr),[1..Length(t)],OnPoints))]);
od;
ps:=ShallowCopy(sizes);
sizes:=Set(sizes); # sizes corresponding to clusters
cluster:=List(sizes,s->Filtered([1..Length(r)],x->ps[x]=s));
# now process per cluster
for i in [1..Length(sizes)] do
sel:=cluster[i];
lr:=r{sel};
lstabs:=stabs{sel};
SortParallel(lr,lstabs); # quick find
IsSSortedList(lr);
ssizes:=List(lstabs,x->sifa/Size(x));
num:=Sum(ssizes)/sizes[i][1]; # number of double cosets to be created
if num>1 and sizes[i][1]/Size(a)<=10*Index(a,a2)^2 then
# fuse orbits together
lr:=List(lr,x->CanonicalRightCosetElement(a,x));
o:=DCFuseSubgroupOrbits(G,b,lr,function(r,g)
return CanonicalRightCosetElement(a,r*g);
end,1000,num);
for j in o do
# record double coset
if flip then
Add(dcs,[lr[j[1]]^(-1),sizes[i][1]]);
else
Add(dcs,[lr[j[1]],sizes[i][1]]);
fi;
Info(InfoCoset,2,"orbit fusion ",Length(dcs)," orblen=",Length(j));
od;
lr:=[];lstabs:=[];
else
while num>1 do
# take first representative as rep for double coset
#stab:=Intersection(b,a^lr[1]);
# check how does its double coset a*lr[1]*b split up into a2-DC's
o:=OrbitsDomain(Image(hom,Intersection(a,b^(lr[1]^-1))),
[1..Length(t)],OnPoints);
# identify which of the a2-cosets they are they are (so we can
# remove them)
o:=List(o,x->Position(lr,canrep(t[x[1]]*lr[1])));
# record double coset
if flip then
Add(dcs,[lr[1]^(-1),sizes[i][1]]);
else
Add(dcs,[lr[1],sizes[i][1]]);
fi;
sel:=Difference([1..Length(lr)],o);
lr:=lr{sel};lstabs:=lstabs{sel};
Info(InfoCoset,2,"new fusion ",Length(dcs)," orblen=",Length(o),
" remainder ",Length(lr));
num:=num-1;
od;
# remainder must be a single double coset
if flip then
Add(dcs,[lr[1]^(-1),sizes[i][1]]);
else
Add(dcs,[lr[1],sizes[i][1]]);
fi;
Info(InfoCoset,2,"final fusion ",Length(dcs)," orblen=",Length(lr),
" remainder ",0);
fi;
od;
fi;
if AssertionLevel()>2 then
# test
bsz:=Size(G);
t:=[];
if flip then
# flip back
c:=a;
a:=b;
b:=c;
fi;
for i in dcs do
bsz:=bsz-i[2];
if AssertionLevel()>0 then
r:=CanonicalRightCosetElement(a,i[1]);
if ForAny(t,j->r in RepresentativesContainedRightCosets(j)) then
Error("duplicate!");
fi;
fi;
r:=DoubleCoset(a,i[1],b);
if AssertionLevel()>0 and Size(r)<>i[2] then
Error("size error!");
fi;
Add(t,r);
od;
if bsz<>0 then
Error("number");
fi;
fi;
PopOptions(); # the usequotient option
return dcs;
end);
InstallMethod(DoubleCosetsNC,"generic",true,
[IsGroup,IsGroup,IsGroup],0,
function(G,U,V)
return List(DoubleCosetRepsAndSizes(G,U,V),i->DoubleCoset(U,i[1],V,i[2]));
end);
InstallMethod(DoubleCosetRepsAndSizes,"generic",true,
[IsGroup,IsGroup,IsGroup],0,
CalcDoubleCosets);
#############################################################################
##
#M RightTransversal generic
##
DeclareRepresentation( "IsRightTransversalViaCosetsRep",
IsRightTransversalRep,
[ "group", "subgroup", "cosets" ] );
InstallMethod(RightTransversalOp, "generic, use RightCosets",
IsIdenticalObj,[IsGroup,IsGroup],0,
function(G,U)
return Objectify( NewType( FamilyObj( G ),
IsRightTransversalViaCosetsRep and IsList and
IsDuplicateFreeList and IsAttributeStoringRep ),
rec( group := G,
subgroup := U,
cosets:=RightCosets(G,U)));
end);
InstallMethod(Length, "for a right transversal in cosets representation",
[IsList and IsRightTransversalViaCosetsRep],
t->Length(t!.cosets));
InstallMethod( \[\], "rt via coset", true,
[ IsList and IsRightTransversalViaCosetsRep, IsPosInt ], 0,
function( cs, num )
return Representative(cs!.cosets[num]);
end );
InstallMethod( PositionCanonical,"rt via coset", IsCollsElms,
[ IsList and IsRightTransversalViaCosetsRep,
IsMultiplicativeElementWithInverse ], 0,
function( cs, elm )
return First([1..Index(cs!.group,cs!.subgroup)],i->elm in cs!.cosets[i]);
end );
InstallMethod(RightCosetsNC,"generic: orbit",IsIdenticalObj,
[IsGroup,IsGroup],0,
function(G,U)
return Orbit(G,RightCoset(U,One(U)),OnRight);
end);
# methods for groups which have a better 'RightTransversal' function
InstallMethod(RightCosetsNC,"perm groups, use RightTransversal",IsIdenticalObj,
[IsPermGroup,IsPermGroup],0,
function(G,U)
return List(RightTransversal(G,U),i->RightCoset(U,i));
end);
InstallMethod(RightCosetsNC,"pc groups, use RightTransversal",IsIdenticalObj,
[IsPcGroup,IsPcGroup],0,
function(G,U)
return List(RightTransversal(G,U),i->RightCoset(U,i));
end);
#############################################################################
##
#M RightTransversalOp( <G>, <U> ) . . . . . . . . . . . . . for trivial <U>
##
InstallMethod( RightTransversalOp,
"for trivial subgroup, call `EnumeratorSorted' for the big group",
IsIdenticalObj,
[ IsGroup, IsGroup and IsTrivial ],
100, # the method for pc groups has this offset but shall be avoided
# because the element enumerator is faster.
function( G, U )
if IsSubgroupFpGroup(G) then
TryNextMethod(); # this method is bad for the fp groups.
fi;
return Enumerator( G );
end );
#############################################################################
##
#R Length, \in functions for transversals via cosets rep
##
InstallMethod(Length, "for a right transversal in cosets representation",
[IsList and IsRightTransversalViaCosetsRep],
t->Length(t!.cosets));
InstallMethod(\in, "for a right coset with representative",
IsElmsColls, [IsObject,IsRightCosetDefaultRep and
HasActingDomain and HasFunctionAction and HasRepresentative],
function(x,C)
return x/Representative(C) in ActingDomain(C);
end);
#############################################################################
##
#R IsFactoredTransversalRep
##
## A transversal stored as product of several shorter transversals
DeclareRepresentation( "IsFactoredTransversalRep",
IsRightTransversalRep,
[ "transversals", "moduli" ] );
# group, subgroup, list of transversals (descending)
BindGlobal("FactoredTransversal",function(G,S,t)
local trans,m,i;
Assert(1,ForAll([1..Length(t)-1],i->t[i]!.subgroup=t[i+1]!.group));
m:=[1];
for i in [Length(t),Length(t)-1..2] do
Add(m,Last(m)*Length(t[i]));
od;
m:=Reversed(m);
trans:=Objectify(NewType(FamilyObj(G),
IsFactoredTransversalRep and IsList
and IsDuplicateFreeList and IsAttributeStoringRep),
rec(group:=G,
subgroup:=S,
transversals:=t,
moduli:=m) );
return trans;
end);
InstallMethod( \[\],"factored transversal",true,
[ IsList and IsFactoredTransversalRep, IsPosInt ], 0,
function( t, num )
local e, m, q, i;
num:=num-1; # indexing with 0 start
e:=One(t!.group);
m:=t!.moduli;
for i in [1..Length(m)] do
q:=QuoInt(num,m[i]);
e:=t!.transversals[i][q+1]*e;
num:=num mod m[i];
od;
return e;
end );
InstallMethod( PositionCanonical, "factored transversal", IsCollsElms,
[ IsList and IsFactoredTransversalRep,
IsMultiplicativeElementWithInverse ], 0,
function( t, elm )
local num, m, p, i;
num:=0;
m:=t!.moduli;
for i in [1..Length(m)] do
p:=PositionCanonical(t!.transversals[i],elm);
elm:=elm/t!.transversals[i][p];
num:=num+(p-1)*m[i];
od;
return num+1;
end );
[ Dauer der Verarbeitung: 0.64 Sekunden
(vorverarbeitet)
]
|
