diff --git a/lib/grppcaut.gi b/lib/grppcaut.gi index 2840694725..8f776a142d 100644 --- a/lib/grppcaut.gi +++ b/lib/grppcaut.gi @@ -1306,7 +1306,9 @@ InstallGlobalFunction(AutomorphismGroupSolvableGroup,function( G ) B := NormalizingReducedGL( spec, s, n, M,B ); # A and B will not be used later, so it is no problem to # replace them by other groups with fewer generators - B:=SubgroupNC(B,SmallGeneratingSet(B)); + if Length(GeneratorsOfGroup(B))>4 then + B:=SubgroupNC(B,SmallGeneratingSet(B)); + fi; if weights[s][2] = 1 then #Info( InfoAutGrp, 2,"compute reduced gl "); @@ -1515,8 +1517,14 @@ InstallGlobalFunction(AutomorphismGroupFrattFreeGroup,function( G ) # go over to perm rep Info( InfoAutGrp, 2, "compute perm rep "); - iso := IsomorphismPermGroup( A ); + if Length(AbelianInvariants(F))>4 then + # likely the element action is best anyhow + iso := ActionHomomorphism(A,Elements(Group(gensF)),"surjective");; + else + iso := IsomorphismPermGroup( A ); + fi; P := Image( iso ); + if HasSize(A) then SetSize(P,Size(A));fi; # compute subgroup Info( InfoAutGrp, 2, "compute subgroup "); @@ -1525,7 +1533,7 @@ InstallGlobalFunction(AutomorphismGroupFrattFreeGroup,function( G ) imgs := List( gensF, y -> y ^ k ); aut := GroupHomomorphismByImagesNC( F, F, gensF, imgs ); # CheckAuto( aut ); - Add( gensU, Image( iso, aut ) ); + Add( gensU, ImagesRepresentative( iso, aut ) ); od; U := SubgroupNC( P, gensU ); hom := GroupHomomorphismByImagesNC( K, U, gensK, gensU ); diff --git a/lib/grppcext.gi b/lib/grppcext.gi index 524f2c5fd1..4f23bd2f1e 100644 --- a/lib/grppcext.gi +++ b/lib/grppcext.gi @@ -361,8 +361,12 @@ local ag, p1iso, agp, p2iso, DP, p1, p2, gens, genimgs, triso,s,i,u,opt, SetSize(ag,s); fi; IsGroupOfAutomorphismsFiniteGroup(ag); - p1iso:=IsomorphismPermGroup(ag); - agp:=Image(p1iso); + # go nice route to avoid an abelian method getting in the way + p1iso:=NiceMonomorphism(ag); + if not IsPermGroup(Image(p1iso)) then + p1iso:=IsomorphismPermGroup(ag); + fi; + agp:=Image(p1iso,ag); # are both groups solvable? p2iso:=IsomorphismPermGroup(DirectProductInfo(D).groups[2]); @@ -502,7 +506,8 @@ end); InstallGlobalFunction( CompatiblePairs, function( arg ) local G, M, Mgrp, oper, A, B, D, translate, gens, genimgs, triso, K, K1, K2, f, tmp, Ggens, pcgs, l, idx, u, tup,Dos,elmlist,preimlist,pows, - baspt,newimgs,i,j,basicact,neu,K1nontriv,epi,hf,pool,modulehom,test; + baspt,newimgs,i,j,basicact,neu,K1nontriv,epi,hf,pool,modulehom,test, + direct; # catch arguments if Length(arg)>2 and IsGroupOfAutomorphismsFiniteGroup(arg[1]) and @@ -517,9 +522,10 @@ local G, M, Mgrp, oper, A, B, D, translate, gens, genimgs, triso, K, K1, Mgrp := GroupByGenerators( M.generators ); Ggens:=Pcgs(G); oper:=fail; - if IsPcgs(Ggens) and Length(Ggens)=Length(M.generators) then - oper := GroupHomomorphismByImagesNC( G, Mgrp, Ggens, M.generators ); - elif Length(arg)=2 then + + direct:=false; + + if A=fail and Length(arg)=2 then # search through automorphism group for projection image and reps, # then add module automorphisms gens:=GeneratorsOfGroup(G); @@ -527,10 +533,37 @@ local G, M, Mgrp, oper, A, B, D, translate, gens, genimgs, triso, K, K1, Info( InfoCompPairs, 1, " CompP: compute aut group"); A:=AutomorphismGroup(G); fi; - triso:=IsomorphismPermGroup(A); + u:=fail; + if Size(A)>1 then + direct:=true; + fi; + + elif A=fail and Length(arg)=3 and HasDirectProductInfo(arg[3]) + and IsGroupOfAutomorphismsFiniteGroup(Image(Projection(arg[3],1))) then + + A:=Image(Projection(arg[3],1)); + B:=Image(Projection(arg[3],2)); + u:=B; + gens:=GeneratorsOfGroup(G); + if Size(A)>1 and Size(B)>1 and Size(A)*Size(B)>1000 then + direct:=true; + fi; + + elif IsPcgs(Ggens) and Length(Ggens)=Length(M.generators) then + oper := GroupHomomorphismByImagesNC( G, Mgrp, Ggens, M.generators ); + fi; + + if direct then + triso:=NiceMonomorphism(A:autactbase:=fail); + if not IsPermGroup(Image(triso)) then + triso:=IsomorphismPermGroup(A:autactbase:=fail); + fi; + pool:=[]; modulehom:=GroupHomomorphismByImages(G,Group(M.generators), gens,M.generators); + M.isMTXModule:=true; + M.IsOverFiniteField:=true; test:=function(perm) local aut,imgs,mat; aut:=PreImagesRepresentative(triso,perm); @@ -569,7 +602,20 @@ local G, M, Mgrp, oper, A, B, D, translate, gens, genimgs, triso, K, K1, A:=GroupWithGenerators(pool); SetSize(A,Size(K)*Size(B)); fi; - return A; + + # make sure the B-part is OK + if u<>fail then + B:=List(GeneratorsOfGroup(A),x->x[2]); + B:=Group(B); + if not IsSubset(u,B) then + Info(InfoCompPairs,1," projection wrong, go normal"); + else + return A; + fi; + else + return A; + fi; + fi; if oper=fail then @@ -577,6 +623,7 @@ local G, M, Mgrp, oper, A, B, D, translate, gens, genimgs, triso, K, K1, oper := GroupHomomorphismByImagesNC( G, Mgrp, Ggens, M.generators ); fi; + # automorphism groups of G and M if Length( arg ) = 2 then if A=fail then @@ -839,7 +886,8 @@ local G, M, Mgrp, oper, A, B, D, translate, gens, genimgs, triso, K, K1, tmp := StabilizerOp( D, rec(hashfun:= lst->lst*pows),tup, gens,newimgs, f ); else - tmp := Stabilizer( D, tup,gens,genimgs, f ); + #tmp := Stabilizer( D, tup,gens,genimgs, f ); + tmp := SubnormalOrbitExtension( D, false, tup,gens,genimgs, f ); fi; else tmp := Stabilizer( D, tup,gens,genimgs, f ); diff --git a/lib/morpheus.gi b/lib/morpheus.gi index 3569b1b529..1dec6ddc9f 100644 --- a/lib/morpheus.gi +++ b/lib/morpheus.gi @@ -475,6 +475,11 @@ local hom, gens, c, ran, r, cen, img, u, orbs, Error("illegal actbase given!"); fi; baddegree:=RootInt(Sum(actbase,Size)^2,3); + if IsPcGroup(g) then + # for pc groups action and storage is cheaper, and there are more + # classes + baddegree:=Maximum(List(actbase,Size)); + fi; else baddegree:=RootInt(Size(g)^3,4); fi; @@ -1926,6 +1931,12 @@ local combi,Gr,Gcl,Ggc,Hr,Hcl,bg,bpri,x,dat, fi; result.inner:=inns; else + dat:=ValueOption("costlimit"); + if IsInt(dat) and Product(List(combi,x->Sum(x,Size)))>dat then + Info(InfoMorph,2,"Morpheus seems to be to costly: ", + Product(List(combi,x->Sum(x,Size)))," vs ",dat); + return -1; # not fail, as this is valid + fi; result:=MorClassLoop(H,combi,result,7); fi; @@ -1933,6 +1944,55 @@ local combi,Gr,Gcl,Ggc,Hr,Hcl,bg,bpri,x,dat, end); + +############################################################################# +## +#F AutOrderAbelian( ) . . order of Aut(A), A the abelian group with +## cyclic factors given by the prime powers +## +## function written by claude, based on formula in HillarRhea07, +## DOI 10.1080/00029890.2007.11920485 +## +AutOrderAbelian := function( L ) +local total, primes, p, q, exps, e, f, power, num, den, c, t; + + for q in L do + if not IsInt( q ) or q <= 1 or not IsPrimePowerInt( q ) then + Error( " must be a list of prime powers" ); + fi; + od; + + total := 1; + primes := Set( List( L, q -> FactorsInt( q )[1] ) ); + + for p in primes do + + # exponents of the p-primary part + exps := List( Filtered( L, q -> q mod p = 0 ), q -> LogInt( q, p ) ); + + # |End(A_p)| = p^( sum_{i,j} min(e_i,e_j) ) + power := 0; + for e in exps do + for f in exps do + power := power + Minimum( e, f ); + od; + od; + + # times prod_k |GL_{m_k}(F_p)| / p^{m_k^2}, kept integral throughout + num := p^power; + den := 1; + for c in Collected( exps ) do # c = [ exponent, multiplicity ] + for t in [ 1 .. c[2] ] do + num := num * ( p^t - 1 ); + den := den * p^t; + od; + od; + total := total * ( num / den ); + od; + + return total; +end; + ############################################################################# ## #F AutomorphismGroupAbelianGroup() @@ -2050,6 +2110,7 @@ local i,j,k,l,m,o,nl,nj,max,r,e,au,p,gens,offs; SetIsGroupOfAutomorphismsFiniteGroup(au,true); fi; + SetSize(au,AutOrderAbelian(List(gens,Order))); return au; end); @@ -2932,6 +2993,19 @@ local m; or Size(SolvableRadical(G))^2>Size(G) or ValueOption("forcetest")=true) and ValueOption("forcetest")<>"old" then + + # catch 2-generator groups with few images + if Length(SmallGeneratingSet(G))=2 then + # 2^28 is an experimental limit from some examples + m:=Morphium(G,H,false:costlimit:=Minimum(Size(G)^2,2^28)); + if IsList(m) and Length(m)=0 then + return fail; + elif m<>-1 then + # otherwise it just failed + return m; + fi; + fi; + # In place until a proper implementation of Cannon/Holt isomorphism is # done return PatheticIsomorphism(G,H); diff --git a/lib/oprtperm.gi b/lib/oprtperm.gi index 0678e7eafb..6f6263a0e4 100644 --- a/lib/oprtperm.gi +++ b/lib/oprtperm.gi @@ -59,6 +59,8 @@ InstallOtherMethod( OrbitStabilizerOp, "permgroup", true, end ); + + ############################################################################# ## #M Orbits( , , , , ) . . . . . . . on integers @@ -177,7 +179,7 @@ InstallMethod( BlocksOp, "permgroup on integers", if Length( orbit ) <> Length( D ) then Error(" must operate transitively on "); fi; - Info( InfoAction, 1, "BlocksNoSeed transversal computed" ); + Info( InfoAction, 4, "BlocksNoSeed transversal computed" ); nrorbs := Length( orbit ); # since $i \in k^{G_1}$ implies $\beta(i)=\beta(k)$, we initialize @@ -251,7 +253,7 @@ InstallMethod( BlocksOp, "permgroup on integers", od; od; - Info( InfoAction, 1, "BlocksNoSeed ", + Info( InfoAction, 4, "BlocksNoSeed ", "number of orbits of < _1 is ",nrorbs ); # take arbitrary point , and an element taking 1 to @@ -328,18 +330,18 @@ InstallMethod( BlocksOp, "permgroup on integers", od; block := Set( block ); blocks := [ block ]; - Info( InfoAction, 1, "BlocksNoSeed ", + Info( InfoAction, 4, "BlocksNoSeed ", "length of alleged block is ",Length(block) ); # quick test to see if the group is primitive if Length( block ) = Length( orbit ) then - Info( InfoAction, 1, "BlocksNoSeed is primitive" ); + Info( InfoAction, 4, "BlocksNoSeed is primitive" ); return Immutable( [ D ] ); fi; # quick test to see if the orbit can be a block if Length( orbit ) mod Length( block ) <> 0 then - Info( InfoAction, 1, "BlocksNoSeed ", + Info( InfoAction, 4, "BlocksNoSeed ", "alleged block is clearly not a block" ); changed := -1000; fi; @@ -372,7 +374,7 @@ InstallMethod( BlocksOp, "permgroup on integers", # check that all points in the image are new for pnt in img do if rep[pnt] <> 0 then - Info( InfoAction, 1, "BlocksNoSeed ", + Info( InfoAction, 4, "BlocksNoSeed ", "alleged block is not a block" ); changed := -1000; fi; @@ -385,7 +387,7 @@ InstallMethod( BlocksOp, "permgroup on integers", # check that all points in the image lie in the block for pnt in img do if rep[pnt] <> rep[img[1]] then - Info( InfoAction, 1, "BlocksNoSeed ", + Info( InfoAction, 4, "BlocksNoSeed ", "alleged block is not a block" ); changed := -1000; fi; @@ -597,7 +599,7 @@ local blocks, # block system of , result if Length( orbit ) <> Length( D ) then Error( " must act transitively on " ); fi; - Info(InfoAction,1,"RepresentativesMinimalBlocks transversal computed"); + Info(InfoAction,4,"RepresentativesMinimalBlocks transversal computed"); nrorbs := Length( orbit ); # since $i \in k^{G_1}$ implies $\beta(i)=\beta(k)$, we initialize @@ -680,7 +682,7 @@ local blocks, # block system of , result od; od; - Info(InfoAction,1,"RepresentativesMinimalBlocks ", + Info(InfoAction,4,"RepresentativesMinimalBlocks ", "number of orbits of < _1 is ",nrorbs); # take arbitrary point , and an element taking 1 to @@ -786,18 +788,18 @@ local blocks, # block system of , result od; block := Set( block ); blocks := [ block ]; - Info(InfoAction,1,"RepresentativesMinimalBlocks ", + Info(InfoAction,4,"RepresentativesMinimalBlocks ", "length of alleged block is ",Length(block)); # quick test to see if the group is primitive if Length( block ) = Length( orbit ) then - Info(InfoAction,1,"RepresentativesMinimalBlocks is primitive"); + Info(InfoAction,4,"RepresentativesMinimalBlocks is primitive"); return Immutable([ D ]); fi; # quick test to see if the orbit can be a block if Length( orbit ) mod Length( block ) <> 0 then - Info(InfoAction,1,"RepresentativesMinimalBlocks ", + Info(InfoAction,4,"RepresentativesMinimalBlocks ", "alleged block is clearly not a block"); changed := -1000; fi; @@ -830,7 +832,7 @@ local blocks, # block system of , result # check that all points in the image are new for pnt in img do if rep[pnt] <> 0 then - Info(InfoAction,1, + Info(InfoAction,4, "RepresentativesMinimalBlocks, alleged block is not a block"); changed := -1000; fi; @@ -843,7 +845,7 @@ local blocks, # block system of , result # check that all points in the image lie in the block for pnt in img do if rep[pnt] <> rep[img[1]] then - Info(InfoAction,1, + Info(InfoAction,4, "RepresentativesMinimalBlocks , alleged block is not a block"); changed := -1000; fi; @@ -1614,6 +1616,304 @@ local S,j; end ); +############################################################################ +## +#F SubnormalOrbitExtension( , , , , , ) +## +## computes the stabilizer of the point +## under the action of the finite group , where the action is +## given in a nonstandard way: the group generators of act as +## the corresponding elements of the list (of the same length), i.e. +## the actual acting objects are the images under the homomorphism +## ->:. If = no such translation takes place. +## is a domain (or a sample list of points) that is used only to create a +## suitable dictionary for the points in the orbit; must be an +## element of the same kind. +## +## The orbit and the stabilizer are built up along a composition series of +## , from the bottom up: for each composition factor the orbit of the +## subgroup reached so far is extended by the new generator(s) and the +## Schreier generators arising from collisions are added to the stabilizer. +## Because every subgroup in the series is normal in the next one, the +## orbit of the smaller group decomposes the orbit of the larger one into +## blocks of equal length; this block structure is used to extend the orbit +## in blocks (avoiding one application of per point and generator) +## and to keep, for every level, only one word (a `representative') per +## block. These block representatives are stored as words in free groups, +## one per level, so that they can be evaluated both in (to build the +## stabilizer) and in the acting group (to trace points back), see the +## subfunction `invrep'. +## +## Only the part of the series below the subgroup generated by the radical +## and the perfect residuum is treated this way. If the orbit becomes large +## in comparison to the index of the current subgroup, the subfunction +## `tryup' attempts to determine the stabilizer instead by extending it +## upwards by normalizing elements of prime power order (`zuppos'), which +## can be substantially cheaper. +## +## The function returns the full stabilizer of in . +## +BindGlobal("SubnormalOrbitExtension",function(G,dom,pnt,gens,acts,act) +local invrep, # function: word mapping a point of level back to + cs, # composition series of G, in increasing order + d, # dictionary: point -> its position in `orb' + orb, # orbit of computed so far + lens, # lens[i]: length of the orbit after level i + rlens, # rlens[i]: orbit growth factor (number of blocks) at level i + levgens, # levgens[i]: new generators of cs[i] (modulo cs[i-1]) + levacts, # levacts[i]: the images of levgens[i] under `hom' (actors) + stb, # stabilizer of in the subgroup treated so far + ablev, # ablev[i]: is the i-th composition factor abelian? + frees, # frees[i]: free group on the generators levgens[i] + i, # loop variable: index of the current level in `cs' + j, # loop variable (generator index/orbit position) + k, # loop variable, also Schreier generator being constructed + a, # new generator(s) of the current level, also scratch/group + b, # image(s) of `a' under `hom', also scratch/group + relo, # relative order (prime index) of an abelian comp. factor + free, # free group for the current level + img, # image of an orbit point under an acting generator + pos, # position of `img' in `orb' (or fail if new) + reps, # reps[i]: words in frees[i] mapping block 1 to each block + lreps, # list of block representative words for the current level + hom, # homomorphism from G to the acting group + leno, # orbit length of the previous level, i.e. the block length + p, # position of the block currently processed in `orb' + fgens, # generators of `free' + fcache, # cache of free groups -- reuse + myfree, # free group fetcher + tryup; # function: try to find the stabilizer by stepping upwards + + # Given the position of a point in the orbit of level , return + # an element of G that maps this point back to , by peeling off one + # level at a time. + invrep:=function(pos,lev) + local elm, # the group element constructed so far + p, # number of the block of level containing the point + rep, # word mapping block 1 to block p (to be inverted) + a, # value of `rep^-1' in G + b; # value of `rep^-1' in the acting group + elm:=One(G); + # go down through levels + while lev>1 do + if rlens[lev]>1 then # otherwise there was no growth + p:=QuoInt(pos-1,lens[lev-1])+1; # number of block image + rep:=reps[lev][p]^-1; # word for mapper back + a:=MappedWord(rep,GeneratorsOfGroup(frees[lev]),levgens[lev]); + b:=MappedWord(rep,GeneratorsOfGroup(frees[lev]),levacts[lev]); + elm:=elm*a; + pos:=act(orb[pos],b); + pos:=LookupDictionary(d,pos); # mapped back + fi; + lev:=lev-1; + od; + Assert(1,pos=1); # the point must have been mapped back to + return elm; + end; + + # Try to enlarge the stabilizer within cs[i] by adjoining normalizing + # elements of prime power order. Returns the enlarged group, or `fail' if + # there are too many candidates for this to pay off. + tryup:=function(s,i) + local n, # normalizer of in cs[i] + map, # natural homomorphism n -> n/s + q, # the factor group n/s + prev, # image of the part of n below cs[i-1] (already treated) + c, # class representatives of q of prime power order + zup, # the relevant zuppos (prime power order elements) of q + o, # order of the class representative `a' + cn, # normalizer of in q (to form the class of `a') + t, # loop variable over a transversal of cn in q + u, # the stabilizer as extended so far + a; # loop variable over class representatives resp. zuppos + n:=Normalizer(cs[i],s); + map:=NaturalHomomorphismByNormalSubgroup(n,s); + q:=Image(map,n); + prev:=Image(map,Intersection(n,cs[i-1])); + # relevant Zuppos + c:=List(ConjugacyClasses(q),Representative); + c:=Filtered(c,x->IsPrimePowerInt(Order(x)) and not x in prev); + zup:=[]; + for a in c do + o:=Order(a); + if ForAll([2..o],x->Gcd(o,x)<>1 or not a^x in zup) then + cn:=Normalizer(q,Subgroup(q,[a])); + for t in RightTransversal(q,cn) do + Add(zup,a^t); + od; + fi; + if Length(zup)>Length(orb) then + return fail; # too many + fi; + od; + Info(InfoAction,1, + "try upstep ",i,": ",Size(q),"; ",Length(zup)," zuppos"); + + # now try expanding with all possible normalizing zuppos + u:=s; + for a in zup do + a:=PreImagesRepresentative(map,a); + if not a in u then + if act(pnt,ImagesRepresentative(hom,a))=pnt then + u:=ClosureSubgroup(u,a); + fi; + fi; + od; + Info(InfoAction,1,"Found extra ",Size(u)/Size(s)); + return u; + end; + + # reuse the same free groups again and again + fcache:=[]; + myfree:=function(n) + if not IsBound(fcache[n]) then + fcache[n]:=FreeGroup(n); + fi; + return fcache[n]; + end; + + b:=ClosureSubgroup(RadicalGroup(G),PerfectResiduum(G)); + cs:=Reversed(CompositionSeriesThrough(G,[RadicalGroup(G),b])); + + d:=NewDictionary(pnt,true,dom); + if gens=acts then + hom:=IdentityMapping(G); + else + hom:=GroupGeneralMappingByImagesNC(G,Group(acts),gens,acts); + fi; + orb:=[pnt]; + AddDictionary(d,pnt,1); + lens:=[1]; + rlens:=[1]; + levgens:=[[]]; + levacts:=[[]]; + stb:=TrivialSubgroup(G); + ablev:=[fail]; + frees:=[0]; + reps:=[0]; + #for i in [2..Length(cs)] do + for i in [2..Position(cs,b)] do + Info(InfoAction,1,"Step ",i,": ",Length(orb)," Stb: ",Size(stb), + " Group ",Index(cs[i],cs[i-1])); + leno:=lens[i-1]; + + # do we want to explore through normalizers up? + if Length(orb)>10^4 and Index(G,cs[i])fail do + a:=tryup(a,j); + j:=j+1; + od; + if a<>fail then + return a; + fi; + fi; + + if HasAbelianFactorGroup(cs[i],cs[i-1]) then + Add(ablev,true); + # standard pc stabilizer + # the new generator + a:=First(GeneratorsOfGroup(cs[i]),x->not x in cs[i-1]); + b:=ImagesRepresentative(hom,a); + relo:=Size(cs[i])/Size(cs[i-1]); + Add(levgens,[a]); + Add(levacts,[b]); + free:=myfree(1); # reuse free groups + Add(frees,free); + img:=act(pnt,b); + pos:=LookupDictionary(d,img); + if pos=fail then + Add(reps,List([0..relo-1],x->free.1^x)); + Add(rlens,relo); + # move orbit as block + Add(orb,img); + AddDictionary(d,img,Length(orb)); + for j in [2..leno] do + img:=act(orb[j],b); + Add(orb,img); + AddDictionary(d,img,Length(orb)); + od; + for k in [3..relo] do + for j in [(k-2)*leno+1..(k-1)*leno] do + img:=act(orb[j],b); + Add(orb,img); + AddDictionary(d,img,Length(orb)); + od; + od; + else + # generator fixes orbit, add stabilizing element + Add(rlens,1); + Add(reps,0); + b:=a*invrep(pos,i-1); + Assert(1,act(pnt,ImagesRepresentative(hom,b))=pnt); # b stabilizes + stb:=ClosureSubgroup(stb,b); + fi; + + else + # nonabelian part + Add(ablev,false); + a:=Filtered(GeneratorsOfGroup(cs[i]),x->not x in cs[i-1]); + if ClosureSubgroup(cs[i-1],a{[1,2]})<>cs[i] then + # pull back small gens + b:=NaturalHomomorphismByNormalSubgroup(cs[i],cs[i-1]); + a:=SmallGeneratingSet(Image(b,cs[i])); + a:=List(a,x->PreImagesRepresentative(b,x)); + fi; + b:=List(a,x->ImagesRepresentative(hom,x)); + Add(levgens,a); + Add(levacts,b); + free:=myfree(Length(a)); + fgens:=GeneratorsOfGroup(free); + Add(frees,free); + # now run an orbit algorithm on orbits + lreps:=[One(free)]; + Add(reps,lreps); + p:=1; + while p<=Length(orb) do + for j in [1..Length(a)] do + img:=act(orb[p],b[j]); + pos:=LookupDictionary(d,img); + if pos=fail then + # new orbit image + Add(lreps,lreps[QuoInt(p,leno)+1]*fgens[j]); + Add(orb,img); + AddDictionary(d,img,Length(orb)); + for k in [1..leno-1] do + img:=act(orb[p+k],b[j]); + Add(orb,img); + AddDictionary(d,img,Length(orb)); + od; + else + # stabilize + k:=lreps[QuoInt(p,leno)+1]*fgens[j]/lreps[QuoInt(pos,leno)+1]; + k:=MappedWord(k,fgens,a); # stabilizes orbit + pos:=LookupDictionary(d,act(pnt,ImagesRepresentative(hom,k))); + k:=k*invrep(pos,i-1); + stb:=ClosureSubgroup(stb,k); + fi; + od; + p:=p+leno; + if Length(orb)*Size(stb)=Size(cs[i]) then + Info(InfoAction,2,"early break ",QuoInt(p,leno)+1," of ",Length(lreps)); + p:=Length(orb)+1; + fi; + od; + Add(rlens,Length(lreps)); + + fi; + # orbit-stabilizer theorem must hold for the group treated so far + Assert(2,Length(orb)*Size(stb)=Size(cs[i])); + + Add(lens,Length(orb)); + + od; + + return stb; + +end); + + ############################################################################# ## #F OnSetsSets( , )