TopElement := function(perm,level)
  local blocksize,ims,nbar;
  nbar := List( [1..level], x -> 2 );
  blocksize := Product(nbar{[2..level]});
  ims := Concatenation(List( [1,2],
    i -> [(i^perm-1)*blocksize+1..i^perm*blocksize] ));
  return PermList(ims);
end;
SequenceElement := function(seq,level)
  local blocksize,i,ims,n,nbar,offset,perm,perms;
  offset := 0;
  ims := [];
  nbar := List( [1..level], i -> 2 );
  blocksize := Product(nbar{[2..level]});
  for i in [1..level-1] do
    n := nbar[i+1]; blocksize := blocksize/n;
    perm := seq[(i -1) mod Length(seq) +1];
    Append(ims,Concatenation(List( [1..n],
      i-> offset + [(i^perm-1)*blocksize+1..i^perm*blocksize] ) ));
    offset := offset + n*blocksize;
  od;
  Add(ims,offset+1);
  return PermList(ims);
end;
a := level -> TopElement((1,2),level);
b := level -> SequenceElement([(1,2),(1,2),()],level);
c := level -> SequenceElement([(1,2),(),(1,2)],level);
d := level -> SequenceElement([(),(1,2),(1,2)],level);
List( [1..10], i -> Size(Group(b(i),c(i))) );
List( [1..10], i -> Group(b(i),c(i)) = Group(b(i),c(i),d(i)) );
G := level -> Group( a(level), b(level), c(level));
List( [1..8], i -> Collected(Factors(Size(G(i)))) );
List( [1..8], i -> NilpotencyClassOfGroup(G(i)) );
List( [1..8], i -> DerivedLength(G(i)) );
List( [1..5], i -> Exponent(G(i)) );
List( [1..5], i -> Length(ConjugacyClasses(G(i))) );
classes :=  ConjugacyClasses(G(5));;   
Collected(List( classes, c -> [Order(Representative(c)),Size(c)] ));
lower := LowerCentralSeries(G(6));;
List( lower, h -> Size(h) );
List( [1..Length(lower)-1],
i -> AbelianInvariants(FactorGroup(lower[i],lower[i+1])) );
upper := UpperCentralSeries(G(6));;                    
List( upper, h -> Size(h) );                              
List( [1..Length(upper)-1],                               
i -> AbelianInvariants(FactorGroup(upper[i],upper[i+1])) );
List(List( [1..6], i ->    
LowerCentralSeries(G(i)) = UpperCentralSeries(G(i)) ));
grigs := List( [1..6], i -> G(i) );;
list := List( grigs, g -> NormalSubgroups(g) );;
List( list, n -> Length(n) );
for normals in list do
Print(Collected(List( normals, h -> Size(h) )),"\n");
od;
LoadPackage("autpgrp");
for g in grigs do
ConvertHybridAutGroup(AutomorphismGroupPGroup(g));
Print(Collected(List( NormalSubgroups(g),
h -> IsCharacteristicSubgroup(g,h)) ), "\n");
od;
for i in [2..5] do
Print( AbelianInvariantsMultiplier(grigs[i]), "\n");
od;
quit;