1#############################################################################
2##
3##  This file is part of GAP, a system for computational discrete algebra.
4##  This file's authors include Alexander Hulpke.
5##
6##  Copyright of GAP belongs to its developers, whose names are too numerous
7##  to list here. Please refer to the COPYRIGHT file for details.
8##
9##  SPDX-License-Identifier: GPL-2.0-or-later
10##
11##  This  file contains declarations for the subgroup lattice functions for
12##  pc groups.
13##
14
15#############################################################################
16##
17#F  InvariantElementaryAbelianSeries( <G>, <morph>, [ <N> ] )
18##           find <morph> invariant EAS of G (through N)
19##
20InstallGlobalFunction(InvariantElementaryAbelianSeries,function(arg)
21local G,morph,N,s,p,e,i,j,k,ise,fine,cor;
22  G:=arg[1];
23  morph:=arg[2];
24  fine:=false;
25  if Length(arg)>2 then
26    N:=arg[3];
27    e:=[G,N];
28    if Length(arg)>3 then
29      fine:=arg[4];
30    fi;
31    if fine then
32      e:=ElementaryAbelianSeries(e);
33    else
34      e:=ElementaryAbelianSeriesLargeSteps(e);
35    fi;
36  else
37    N:=TrivialSubgroup(G);
38    e:=DerivedSeriesOfGroup(G);
39    e:=ElementaryAbelianSeriesLargeSteps(e);
40  fi;
41  s:=[G];
42  i:=2;
43  while i<=Length(e) do
44    # intersect all images of normal subgroup to obtain invariant one
45    # as G is invariant, we dont have to deal with special cases
46    ise:=[e[i]];
47    cor:=e[i];
48    for j in ise do
49      for k in morph do
50	p:=Image(k,j);
51	if not IsSubset(p,cor) then
52	  Add(ise,p);
53	  cor:=Intersection(cor,p);
54        fi;
55      od;
56    od;
57    Assert(1,HasElementaryAbelianFactorGroup(s[Length(s)],cor));
58    ise:=cor;
59    Add(s,ise);
60    p:=Position(e,ise);
61    if p<>fail then
62      i:=p+1;
63    elif fine then
64      e:=ElementaryAbelianSeries([G,ise,TrivialSubgroup(G)]);
65      i:=Position(e,ise)+1;
66    else
67      e:=ElementaryAbelianSeriesLargeSteps([G,ise,TrivialSubgroup(G)]);
68      i:=Position(e,ise)+1;
69    fi;
70    Assert(1,ise in e);
71  od;
72  return s;
73end);
74
75#############################################################################
76##
77#F  InducedAutomorphism(<epi>,<aut>)
78##
79InstallGlobalFunction(InducedAutomorphism,function(epi,aut)
80local f;
81  f:=Range(epi);
82  if HasIsConjugatorAutomorphism( aut ) and IsConjugatorAutomorphism( aut )
83     and ConjugatorOfConjugatorIsomorphism( aut ) in Source( epi ) then
84    aut:= ConjugatorAutomorphismNC( f,
85              Image( epi, ConjugatorOfConjugatorIsomorphism( aut ) ) );
86  else
87    aut:= GroupHomomorphismByImagesNC(f,f,GeneratorsOfGroup(f),
88				   List(GeneratorsOfGroup(f),
89	       i->Image(epi,Image(aut,PreImagesRepresentative(epi,i)))));
90    SetIsInjective(aut,true);
91    SetIsSurjective(aut,true);
92  fi;
93  return aut;
94end);
95
96#############################################################################
97##
98#F  InvariantSubgroupsElementaryAbelianGroup(<G>,<homs>[,<dims])  submodules
99#F    find all subgroups of el. ab. <G>, which are invariant under all <homs>
100#F    which have dimension in dims
101##
102InstallGlobalFunction(InvariantSubgroupsElementaryAbelianGroup,function(arg)
103local g,op,a,pcgs,ma,mat,d,f,i,j,new,newmat,id,p,dodim,compldim,compl,dims,nm;
104  g:=arg[1];
105  op:=arg[2];
106  if not IsElementaryAbelian(g) then
107    Error("<g> must be a vector space");
108  fi;
109  if IsTrivial(g) then
110    return [g];
111  fi;
112  pcgs:=Pcgs(g);
113  d:=Length(pcgs);
114  p:=RelativeOrderOfPcElement(pcgs,pcgs[1]);
115  f:=GF(p);
116  if Length(arg)=2 then
117    dims:=[0..d];
118  else
119    dims:=arg[3];
120  fi;
121
122  if Length(dims)=0 then
123    return [];
124  fi;
125
126  if Length(op)=0 then
127
128    # trivial operation: enumerate subspaces
129    # check which dimensions we'll need
130    ma:=QuoInt(d,2);
131    dodim:=[];
132    compldim:=[];
133    for i in dims do
134      if i<=ma then
135        AddSet(dodim,i);
136      else
137        AddSet(dodim,d-i);
138	AddSet(compldim,d-i);
139      fi;
140    od;
141    if d<3 then compldim:=[]; fi;
142    dodim:=Maximum(dodim);
143
144    # enumerate spaces
145    id:= Immutable( IdentityMat(d, One(f)) );
146    ma:=[[],[ShallowCopy(id[1])]];
147    ConvertToMatrixRep(ma[2],f);
148    # the complements to ma
149    if d>1 then
150      compl:=[ShallowCopy(id)];
151    else
152      compl:=[];
153    fi;
154    if d>2 then
155      nm:=TriangulizedNullspaceMat(TransposedMat(id{[1]}));
156      ConvertToMatrixRep(nm,f);
157      Add(compl,nm);
158    fi;
159    for i in [2..d] do
160      new:=[];
161      for mat in ma do
162	# subspaces of equal dimension
163	for j in [0..p^Length(mat)-1] do
164	  if j=0 then
165	    # special case for subspace of higher dimension
166	    if Length(mat)<dodim then
167	      newmat:=Concatenation(mat,[id[i]]);
168	      ConvertToMatrixRep(newmat,f);
169	    else
170	      newmat:=false;
171	    fi;
172	  else
173	    # possible extension number d
174	    a:=CoefficientsQadic(j,p)*One(f);
175	    newmat:=List(mat,ShallowCopy);
176	    for j in [1..Length(a)] do
177		newmat[j][i]:=a[j];
178	    od;
179	    ConvertToMatrixRep(newmat,f);
180	  fi;
181	  if newmat<>false then
182	    # we will need the space for the next level
183	    Add(new,newmat);
184
185	    # note complements if necc.
186	    if Length(newmat) in compldim then
187	      nm:=NullspaceMat(TransposedMat(newmat));
188	      ConvertToMatrixRep(nm,f);
189	      Add(compl,nm);
190	      #Add(compl,List(NullspaceMat(TransposedMat(newmat*One(f))),
191	      #               i->List(i,IntFFE)));
192	    fi;
193	  fi;
194        od;
195      od;
196      ma:=Concatenation(ma,new);
197    od;
198
199    ma:=Concatenation(ma,compl);
200
201    # take only those of right dim
202    ma:=Filtered(ma,i->Length(i) in dims);
203
204    # convert to grps (noting also the triv. one)
205    new:=ma;
206    for i in [1..Length(new)] do
207      #a:=SubgroupNC(Parent(g),List(i,j->Product([1..d],k->pcgs[k]^j[k])));
208      ma:=new[i];
209      a:=SubgroupNC(Parent(g),List(ma,
210	                  j->PcElementByExponentsNC(pcgs,List(j,IntFFE))));
211#      a:=MySubgroupNC(Parent(g),List(i,j->PcElementByExponentsNC(pcgs,j)),
212#                      IsFinite and IsSubsetLocallyFiniteGroup and
213#		      IsSupersolvableGroup and IsNilpotentGroup and
214#		      IsCommutative and IsElementaryAbelian);
215
216      SetSize(a,p^Length(ma));
217      new[i]:=a;
218    od;
219    ma:=new;
220
221  else
222
223    # compute representation
224    ma:=[];
225    for i in op do
226      mat:=[];
227      for j in pcgs do
228	Add(mat,ExponentsOfPcElement(pcgs,Image(i,j))*One(f));
229      od;
230      mat:=ImmutableMatrix(f,mat);
231      Add(ma,mat);
232    od;
233
234    ma:=GModuleByMats(ma,f);
235    mat:=MTX.BasesSubmodules(ma);
236
237    ma:=[];
238    for i in mat do
239      Add(ma,SubgroupNC(Parent(g),
240		      List(i,j->PcElementByExponentsNC(pcgs,j))));
241		      #List(i,j->Product([1..d],k->pcgs[k]^IntFFE(j[k])))));
242    od;
243  fi;
244  return ma;
245end);
246
247#############################################################################
248##
249#F  ActionSubspacesElementaryAbelianGroup(<P>,<G>[,<dims>])
250##
251##  compute the permutation action of <P> on the subspaces of the
252##  elementary abelian subgroup <G> of <P>. Returns
253##  a list [<subspaces>,<action>], where <subspaces> is a list of all the
254##  subspaces (as groups) and <action> a homomorphism from <P> in a
255##  permutation group, which is equal to the action homomrophism for the
256##  action of <P> on <subspaces>. If <dims> is given, only subspaces of
257##  dimension <dims> are considered.  Instead of <G> also a (modulo) pcgs
258##  may be given, in this case <subspaces> are pre-images of the subspaces.
259##
260InstallGlobalFunction(ActionSubspacesElementaryAbelianGroup,function(arg)
261local P,g,op,act,a,pcgs,ma,mat,d,f,i,j,new,newmat,id,p,dodim,compldim,compl,
262      dims,Pgens,Pcgens,Pu,Pc,perms,one,par,ker,kersz,pcelm,pccache,asz;
263
264  P:=arg[1];
265  if IsModuloPcgs(arg[2]) then
266    pcgs:=arg[2];
267    g:=Group(NumeratorOfModuloPcgs(pcgs));
268    if not IsSubset(Parent(P),g) then # for matrix groups we need a parent here.
269      par:=ClosureGroup(Parent(P),g);
270    else
271      par:=P;
272    fi;
273    Pu:=AsSubgroup(par,g);
274    ker:=SubgroupNC(par,DenominatorOfModuloPcgs(pcgs));
275    kersz:=Size(ker);
276  else
277    kersz:=1;
278    g:=arg[2];
279    par:=Parent(g);
280    Pu:=Centralizer(P,g);
281    if not IsElementaryAbelian(g) then
282      Error("<g> must be a vector space");
283    fi;
284    if IsTrivial(g) then
285      return [g];
286    fi;
287
288    pcgs:=Pcgs(g);
289  fi;
290
291  d:=Length(pcgs);
292  p:=RelativeOrderOfPcElement(pcgs,pcgs[1]);
293  f:=GF(p);
294  one:=One(f);
295  if Length(arg)=2 then
296    dims:=[0..d];
297  else
298    dims:=arg[3];
299  fi;
300
301  if Length(dims)=0 then
302    return [];
303  fi;
304
305  # find representatives generating the acting factor
306  Pgens:=[];
307  Pc:=Pu;
308  Pcgens:=GeneratorsOfGroup(Pu);
309  while Size(Pu)<Size(P) do
310    repeat
311      i:=Random(P);
312    until not i in Pu;
313    Add(Pgens,i);
314    Pu:=ClosureGroup(Pu,i);
315  od;
316  if Length(Pgens)>2 and Length(Pgens)>Length(SmallGeneratingSet(P)) then
317    Pgens:=SmallGeneratingSet(P);
318  fi;
319
320  # compute representation
321  op:=[];
322  for i in Pgens do
323    mat:=[];
324    for j in pcgs do
325      Add(mat,ExponentsConjugateLayer(pcgs,j,i)*One(f));
326    od;
327    mat:=ImmutableMatrix(f,mat);
328    Add(op,mat);
329  od;
330
331  # and action on canonical bases
332  #act:=function(bas,m)
333  #  bas:=bas*m;
334  #  bas:=List(bas,ShallowCopy);
335  #  TriangulizeMat(bas);
336  #  bas:=List(bas,IntVecFFE);
337  #  return bas;
338  #end;
339  if p=2 then
340    act:=OnSubspacesByCanonicalBasisGF2;
341  else
342    act:=OnSubspacesByCanonicalBasis;
343  fi;
344
345  # enumerate subspaces
346  # check which dimensions we'll need
347  ma:=QuoInt(d,2);
348  dodim:=[];
349  compldim:=[];
350  for i in dims do
351    if i<=ma then
352      AddSet(dodim,i);
353    else
354      AddSet(dodim,d-i);
355      AddSet(compldim,d-i);
356    fi;
357  od;
358  if d<3 then compldim:=[]; fi;
359  dodim:=Maximum(dodim);
360
361  # enumerate spaces
362  id:= Immutable( IdentityMat(d, one) );
363  ma:=[[],[id[1]]];
364  # the complements to ma
365  if d>1 then
366    compl:=[ShallowCopy(id)];
367  else
368    compl:=[];
369  fi;
370  if d>2 then
371    Add(compl,List(TriangulizedNullspaceMat(TransposedMat(id{[1]})),
372                   ShallowCopy));
373  fi;
374  for i in [2..d] do
375    new:=[];
376    for mat in ma do
377      # subspaces of equal dimension
378      for j in [0..p^Length(mat)-1] do
379	if j=0 then
380	  # special case for subspace of higher dimension
381	  if Length(mat)<dodim then
382	    newmat:=Concatenation(mat,[id[i]]);
383	  else
384	    newmat:=false;
385	  fi;
386	else
387	  # possible extension number d
388	  a:=CoefficientsQadic(j,p)*one;
389	  newmat:=List(mat,ShallowCopy);
390	  for j in [1..Length(a)] do
391	      newmat[j][i]:=a[j];
392	  od;
393	fi;
394	if newmat<>false then
395	  # we will need the space for the next level
396	  Add(new,Immutable(newmat));
397
398	  # note complements if necc.
399	  if Length(newmat) in compldim then
400	    a:=List(TriangulizedNullspaceMat(MutableTransposedMat(newmat)),
401	            ShallowCopy);
402	    Add(compl,Immutable(a));
403	  fi;
404	fi;
405      od;
406    od;
407
408    ma:=Concatenation(ma,new);
409  od;
410
411  ma:=Concatenation(ma,compl);
412
413  # take only those of right dim
414  ma:=Filtered(ma,i->Length(i) in dims);
415
416  perms:=List(Pgens,i->());
417  new:=[];
418  for i in dims do
419    mat:=Immutable(Set(Filtered(ma,j->Length(j)=i)));
420    # compute action on mat
421    if i>0 and i<d then
422      for j in [1..Length(Pgens)] do
423	#a:=Permutation(op[j],mat,act);
424	a:=List([1..Length(mat)],k->PositionSorted(mat,act(mat[k],op[j])));
425	a:=PermList(a);
426	perms[j]:=perms[j]*a^MappingPermListList([1..Length(mat)],
427				[Length(new)+1..Length(new)+Length(mat)]);
428      od;
429    fi;
430    Append(new,mat);
431  od;
432  ma:=new;
433
434  # convert to grps
435  pccache:=[]; # avoid recerating different copies of same element
436  pcelm:=function(vec)
437  local e,p;
438    e:=Immutable([vec]);
439    p:=PositionSorted(pccache,e);
440    if IsBound(pccache[p]) and pccache[p][1]=vec then
441      return pccache[p][2];
442    else
443      e:=Immutable([vec,PcElementByExponentsNC(pcgs,vec)]);
444      AddSet(pccache,e);
445      return e[2];
446    fi;
447  end;
448
449  new:=[];
450  for i in ma do
451    #a:=SubgroupNC(Parent(g),List(i,j->Product([1..d],k->pcgs[k]^j[k])));
452    asz:=kersz*p^Length(i);
453    if kersz=1 then
454      a:=SubgroupNC(par,List(i,j->pcelm(j)));
455    else
456      #a:=ClosureSubgroup(ker,List(i,j->pcelm(j)):knownClosureSize:=asz);
457      a:=SubgroupNC(par,Concatenation(GeneratorsOfGroup(ker),
458        List(i,j->pcelm(j))));
459    fi;
460    SetSize(a,asz);
461    Add(new,a);
462  od;
463
464  ma:= GroupByGenerators( perms, () );
465  #Assert(1,Group(perms)=Action(P,new));
466
467  op:=GroupHomomorphismByImagesNC(P,ma,Concatenation(Pcgens,Pgens),
468    Concatenation(List(Pcgens,i->()),perms));
469#  Assert(1,Size(P)=Size(KernelOfMultiplicativeGeneralMapping(op))
470#                   *Size(Image(op)));
471  return [new,op];
472
473end);
474
475# test whether the c-conjugate of g is h-invariant, internal
476HasInvariantConjugateSubgroup:=function(g,c,h)
477  # This should be done better!
478  g:=ConjugateSubgroup(g,c);
479  return ForAll(h,i->Image(i,g)=g);
480end;
481
482#############################################################################
483##
484#F  SubgroupsSolvableGroup(<G>[,<opt>]) . classreps of subgrps of <G>,
485##   				             <homs>-inv. with options.
486##    Options are:
487##                  actions:  list of automorphisms: search for invariants
488##		    funcnorm: N_G(actions) (speeds up calculation)
489##                  normal:   just search for normal subgroups
490##                  consider: function(A,N,B,M) indicator function, whether
491##			      complements of this type would be needed
492##                  retnorm:  return normalizers
493##
494InstallGlobalFunction(SubgroupsSolvableGroup,function(arg)
495local g,	# group
496      isom,	# isomorphism onto AgSeries group
497      func,	# automorphisms to be invariant under
498      funcs,    # <func>
499      funcnorm, # N_G(funcs)
500      efunc,	# induced automs on factor
501      efnorm,	# funcnorm^epi
502      e,	# EAS
503      len,	# Length(e)
504      start,	# last index with EA factor
505      i,j,k,l,
506      m,kp,	# loop
507      kgens,	# generators of k
508      kconh,	# complemnt conjugacy storage
509      opt,	# options record
510      normal,	# flag for 'normal' option
511      consider,	# optional 'consider' function
512      retnorm,	# option: return all normalizers
513      f,	# g/e[i]
514      home,	# HomePcgs(f)
515      epi,	# g -> f
516      lastepi,  # epi of last step
517      n,	# e[i-1]^epi
518      fa,	# f/n = g/e[i-1]
519      hom,	# f -> fa
520      B,	# subgroups of n
521      ophom,	# perm action of f on B (or false if not computed)
522      a,	# preimg. of group over n
523      no,	# N_f(a)
524#      aop,	# a^ophom
525#      nohom,	# ophom\rest no
526      oppcgs,	# acting pcgs
527      oppcgsimg,# images under ophom
528      ex,	# external set/orbits
529      bs,	# b\in B normal under a, reps
530      bsp,	# bs index
531      bsnorms,	# respective normalizers
532      b,	# in bs
533      bpos,	# position in bs
534      hom2,	# N_f(b) -> N_f(b)/b
535      nag,	# AgGroup(n^hom2)
536      fghom,	# assoc. epi
537      t,s,	# dnk-transversals
538      z,	# Cocycles
539      coboundbas,# Basis(OneCobounds)
540      field,	# GF(Exponent(n))
541      com,	# complements
542      comnorms,	# normalizers supergroups
543      isTrueComnorm, # is comnorms the true normalizer or a supergroup
544      comproj,	# projection onto complement
545      kgn,
546      kgim,	# stored decompositions, translated to matrix language
547      kgnr,	# assoc index
548      ncom,	# dito, tested
549      idmat,	# 1-matrix
550      mat,	# matrix action
551      mats,	# list of mats
552      conj,	# matrix action
553      chom,	# homom onto <conj>
554      shom,	# by s induced autom
555      shoms,	# list of these
556      smats,	# dito, matrices
557      conjnr,	# assoc. index
558      glsyl,
559      glsyr,	# left and right side of eqn system
560      found,	# indicator for success
561      grps,	# list of subgroups
562      ngrps,	# dito, new level
563      gj,	# grps[j]
564      grpsnorms,# normalizers of grps
565      ngrpsnorms,# dito, new level
566      bgids,    # generators of b many 1's (used for copro)
567      opr,	# operation on complements
568      xo;	# xternal orbits
569
570  g:=arg[1];
571  if Size(g)=1 then
572    return [g];
573  fi;
574  if Length(arg)>1 and IsRecord(arg[Length(arg)]) then
575    opt:=arg[Length(arg)];
576  else
577    opt:=rec();
578  fi;
579
580  # parse options
581  normal:=IsBound(opt.normal) and opt.normal=true;
582  if IsBound(opt.consider) then
583    consider:=opt.consider;
584  else
585    consider:=false;
586  fi;
587
588  retnorm:=IsBound(opt.retnorm) and opt.retnorm;
589
590  isom:=fail;
591
592  # get automorphisms and compute their normalizer, if applicable
593  if IsBound(opt.actions) then
594    func:=opt.actions;
595    hom2:= Filtered( func,     HasIsConjugatorAutomorphism
596			   and IsInnerAutomorphism );
597    hom2:= List( hom2, ConjugatorOfConjugatorIsomorphism );
598
599    if IsBound(opt.funcnorm) then
600      # get the func. normalizer
601      funcnorm:=opt.funcnorm;
602      b:=g;
603    else
604      funcs:= GroupByGenerators( Filtered( func,
605                  i -> not ( HasIsConjugatorAutomorphism( i ) and
606                             IsInnerAutomorphism( i ) ) ),
607		   IdentityMapping(g));
608      IsGroupOfAutomorphismsFiniteGroup(funcs); # set filter
609      if IsTrivial( funcs ) then
610	b:=ClosureGroup(Parent(g),List(func,x->ConjugatorOfConjugatorIsomorphism(x)));
611	func:=hom2;
612      else
613        if IsSolvableGroup(funcs) then
614	  a:=IsomorphismPcGroup(funcs);
615	else
616	  a:=IsomorphismPermGroup(funcs);
617	fi;
618	hom:=InverseGeneralMapping(a);
619	IsTotal(hom); IsSingleValued(hom); # to be sure (should be set anyway)
620	b:=SemidirectProduct(Image(a),hom,g);
621	hom:=Embedding(b,1);
622	funcs:=List(GeneratorsOfGroup(funcs),i->Image(hom,Image(a,i)));
623	isom:=Embedding(b,2);
624	hom2:=List(hom2,i->Image(isom,i));
625	func:=Concatenation(funcs,hom2);
626	g:=Image(isom,g);
627      fi;
628
629      # get the normalizer of <func>
630      funcnorm:=Normalizer(g,SubgroupNC(b,func));
631      func:=List(func,i->ConjugatorAutomorphism(b,i));
632    fi;
633
634    Assert(1,IsSubgroup(g,funcnorm));
635
636    # compute <func> characteristic series
637    e:=InvariantElementaryAbelianSeries(g,func);
638  else
639    func:=[];
640    funcnorm:=g;
641    e:=ElementaryAbelianSeriesLargeSteps(g);
642  fi;
643
644  if IsBound(opt.series) then
645    e:=opt.series;
646  else
647    f:=DerivedSeriesOfGroup(g);
648    if Length(e)>Length(f) and
649      ForAll([1..Length(f)-1],i->IsElementaryAbelian(f[i]/f[i+1])) then
650      Info(InfoPcSubgroup,1,"  Preferring Derived Series");
651      e:=f;
652    fi;
653  fi;
654
655#  # check, if the series is compatible with the AgSeries and if g is a
656#  # parent group. If not, enforce this
657#  if not(IsParent(g) and ForAll(e,i->IsElementAgSeries(i))) then
658#    Info(InfoPcSubgroup,1,"  computing better series");
659#    isom:=IsomorphismAgGroup(e);
660#    g:=Image(isom,g);
661#    e:=List(e,i->Image(isom,i));
662#    funcnorm:=Image(isom,funcnorm);
663#
664#    #func:=List(func,i->isom^-1*i*isom);
665#    hom:=[];
666#    for i in func do
667#      hom2:=GroupHomomorphismByImagesNC(g,g,g.generators,List(g.generators,
668#                 j->Image(isom,Image(i,PreImagesRepresentative(isom,j)))));
669#      hom2.isMapping:=true;
670#      Add(hom,hom2);
671#    od;
672#    func:=hom;
673#  else
674#    isom:=false;
675#  fi;
676
677  len:=Length(e);
678
679  if IsBound(opt.groups) then
680    start:=0;
681    while start+1<=Length(e) and ForAll(opt.groups,i->IsSubgroup(e[start+1],i)) do
682      start:=start+1;
683    od;
684    Info(InfoPcSubgroup,1,"starting index ",start);
685    epi:=NaturalHomomorphismByNormalSubgroup(g,e[start]);
686    lastepi:=epi;
687    f:=Image(epi,g);
688    grps:=List(opt.groups,i->Image(epi,i));
689    if not IsBound(opt.grpsnorms) then
690      opt:=ShallowCopy(opt);
691      opt.grpsnorms:=List(opt.groups,i->Normalizer(e[1],i));
692    fi;
693    grpsnorms:=List(opt.grpsnorms,i->Image(epi,i));
694  else
695    # search the largest elementary abelian quotient
696    start:=2;
697    while start<len and IsElementaryAbelian(g/e[start+1]) do
698      start:=start+1;
699    od;
700
701    # compute all subgroups there
702    if start<len then
703      # form only factor groups if necessary
704      epi:=NaturalHomomorphismByNormalSubgroup(g,e[start]);
705      LockNaturalHomomorphismsPool(g,e[start]);
706      f:=Image(epi,g);
707    else
708      f:=g;
709      epi:=IdentityMapping(f);
710    fi;
711    lastepi:=epi;
712    efunc:=List(func,i->InducedAutomorphism(epi,i));
713    grps:=InvariantSubgroupsElementaryAbelianGroup(f,efunc);
714    Assert(1,ForAll(grps,i->ForAll(efunc,j->Image(j,i)=i)));
715    grpsnorms:=List(grps,i->f);
716    Info(InfoPcSubgroup,5,List(grps,Size),List(grpsnorms,Size));
717
718  fi;
719
720  for i in [start+1..len] do
721    Info(InfoPcSubgroup,1," step ",i,": ",Index(e[i-1],e[i]),", ",
722                    Length(grps)," groups");
723    # compute modulo e[i]
724    if i<len then
725      # form only factor groups if necessary
726      epi:=NaturalHomomorphismByNormalSubgroup(g,e[i]);
727      f:=Image(epi,g);
728    else
729      f:=g;
730      epi:=IdentityMapping(g);
731    fi;
732    home:=HomePcgs(f); # we want to compute wrt. this pcgs
733    n:=Image(epi,e[i-1]);
734
735    # the induced factor automs
736    efunc:=List(func,i->InducedAutomorphism(epi,i));
737    # filter the non-trivial ones
738    efunc:=Filtered(efunc,i->ForAny(GeneratorsOfGroup(f),j->Image(i,j)<>j));
739
740    if Length(efunc)>0 then
741      efnorm:=Image(epi,funcnorm);
742    fi;
743
744    if Length(efunc)=0 then
745      ophom:=ActionSubspacesElementaryAbelianGroup(f,n);
746      B:=ophom[1];
747      Info(InfoPcSubgroup,2,"  ",Length(B)," normal subgroups");
748      ophom:=ophom[2];
749
750      ngrps:=[];
751      ngrpsnorms:=[];
752      oppcgs:=Pcgs(Source(ophom));
753      oppcgsimg:=List(oppcgs,i->Image(ophom,i));
754      ex:=[1..Length(B)];
755      IsSSortedList(ex);
756      ex:=ExternalSet(Source(ophom),ex,oppcgs,oppcgsimg,OnPoints);
757      ex:=ExternalOrbitsStabilizers(ex);
758
759      for j in ex do
760        Add(ngrps,B[Representative(j)]);
761	Add(ngrpsnorms,StabilizerOfExternalSet(j));
762#	Assert(1,Normalizer(f,B[j[1]])=ngrpsnorms[Length(ngrps)]);
763      od;
764
765    else
766      B:=InvariantSubgroupsElementaryAbelianGroup(n,efunc);
767      ophom:=false;
768      Info(InfoPcSubgroup,2,"  ",Length(B)," normal subgroups");
769
770      # note the groups in B
771      ngrps:=SubgroupsOrbitsAndNormalizers(f,B,false);
772      ngrpsnorms:=List(ngrps,i->i.normalizer);
773      ngrps:=List(ngrps,i->i.representative);
774    fi;
775
776    # Get epi to the old factor group
777    # as hom:=NaturalHomomorphism(f,fa); does not work, we have to play tricks
778    hom:=lastepi;
779    lastepi:=epi;
780    fa:=Image(hom,g);
781
782    hom:= GroupHomomorphismByImagesNC(f,fa,GeneratorsOfGroup(f),
783           List(GeneratorsOfGroup(f),i->
784	     Image(hom,PreImagesRepresentative(epi,i))));
785    Assert(2,KernelOfMultiplicativeGeneralMapping(hom)=n);
786
787    # lift the known groups
788    for j in [1..Length(grps)] do
789
790      gj:=grps[j];
791      if Size(gj)>1 then
792	a:=PreImage(hom,gj);
793	Assert(1,Size(a)=Size(gj)*Size(n));
794	Add(ngrps,a);
795	no:=PreImage(hom,grpsnorms[j]);
796
797	Add(ngrpsnorms,no);
798
799	if Length(efunc)>0 then
800	  # get the double cosets
801	  t:=List(DoubleCosets(f,no,efnorm),Representative);
802	  Info(InfoPcSubgroup,2,"  |t|=",Length(t));
803	  t:=Filtered(t,i->HasInvariantConjugateSubgroup(a,i,efunc));
804	  Info(InfoPcSubgroup,2,"invar:",Length(t));
805        fi;
806
807	# we have to extend with those b in B, that are normal in a
808	if ophom<>false then
809	  #aop:=Image(ophom,a);
810	  #SetIsSolvableGroup(aop,true);
811
812	  if Length(GeneratorsOfGroup(a))>2 then
813	    bs:=SmallGeneratingSet(a);
814	  else
815	    bs:=GeneratorsOfGroup(a);
816	  fi;
817	  bs:=List(bs,i->Image(ophom,i));
818
819	  bsp:=Filtered([1..Length(B)],i->ForAll(bs,j->i^j=i)
820	                                 and Size(B[i])<Size(n));
821	  bs:=B{bsp};
822	else
823	  bsp:=false;
824	  bs:=Filtered(B,i->IsNormal(a,i) and Size(i)<Size(n));
825	fi;
826
827        if Length(efunc)>0 and Length(t)>1 then
828	  # compute also the invariant ones under the conjugates:
829	  # equivalently: Take all equivalent ones and take those, whose
830	  # conjugates lie in a and are normal under a
831	  for k in Filtered(t,i->not i in no) do
832	    bs:=Union(bs,Filtered(List(B,i->ConjugateSubgroup(i,k^(-1))),
833		  i->IsSubset(a,i) and IsNormal(a,i) and Size(i)<Size(n) ));
834	  od;
835	fi;
836
837	# take only those bs which are valid
838	if consider<>false then
839	  Info(InfoPcSubgroup,2,"  ",Length(bs)," subgroups lead to ");
840	  if bsp<>false then
841	    bsp:=Filtered(bsp,j->consider(no,a,n,B[j],e[i])<>false);
842	    IsSSortedList(bsp);
843	    bs:=bsp; # to get the 'Info' right
844	  else
845	    bs:=Filtered(bs,j->consider(no,a,n,j,e[i])<>false);
846	  fi;
847	  Info(InfoPcSubgroup,2,Length(bs)," valid ones");
848	fi;
849
850	if ophom<>false then
851	  #nohom:=List(GeneratorsOfGroup(no),i->Image(ophom,i));
852	  #aop:=SubgroupNC(Image(ophom),nohom);
853	  #nohom:=GroupHomomorphismByImagesNC(no,aop,
854	  #                                   GeneratorsOfGroup(no),nohom);
855
856	  if Length(bsp)>0 then
857	    oppcgs:=Pcgs(no);
858	    oppcgsimg:=List(oppcgs,i->Image(ophom,i));
859	    ex:=ExternalSet(no,bsp,oppcgs,oppcgsimg,OnPoints);
860	    ex:=ExternalOrbitsStabilizers(ex);
861
862	    bs:=[];
863	    bsnorms:=[];
864	    for bpos in ex do
865	      Add(bs,B[Representative(bpos)]);
866	      Add(bsnorms,StabilizerOfExternalSet(bpos));
867#	    Assert(1,Normalizer(no,B[bpos[1]])=bsnorms[Length(bsnorms)]);
868	    od;
869          fi;
870
871	else
872	  # fuse under the action of no and compute the local normalizers
873	  bs:=SubgroupsOrbitsAndNormalizers(no,bs,true);
874	  bsnorms:=List(bs,i->i.normalizer);
875	  bs:=List(bs,i->i.representative);
876        fi;
877
878Assert(1,ForAll(bs,i->ForAll(efunc,j->Image(j,i)=i)));
879
880	# now run through the b in bs
881	for bpos in [1..Length(bs)] do
882	  b:=bs[bpos];
883	  Assert(2,IsNormal(a,b));
884	  # test, whether we'll have to consider this case
885
886# this test has basically be done before the orbit calculation already
887#	  if consider<>false and consider(a,n,b,e[i])=false then
888#	    Info(InfoPcSubgroup,2,"  Ignoring case");
889#	    s:=[];
890
891	  # test, whether b is invariant
892	  if Length(efunc)>0 then
893	    # extend to dcs of bnormalizer
894	    s:=RightTransversal(no,bsnorms[bpos]);
895	    nag:=Length(s);
896	    s:=Concatenation(List(s,i->List(t,j->i*j)));
897	    z:=Length(s);
898	    #NOCH: Fusion
899	    # test, which ones are usable at all
900	    s:=Filtered(s,i->HasInvariantConjugateSubgroup(b,i,efunc));
901	    Info(InfoPcSubgroup,2,"  |s|=",nag,"-(m)>",z,"-(i)>",Length(s));
902	  else
903	    s:=[()];
904	  fi;
905
906          if Length(s)>0 then
907	    nag:=InducedPcgs(home,n);
908	    nag:=nag mod InducedPcgs(nag,b);
909#	    if Index(Parent(a),a.normalizer)>1 then
910#	      Info(InfoPcSubgroup,2,"  normalizer index ",
911#	                      Index(Parent(a),a.normalizer));
912#	    fi;
913
914	    z:=rec(group:=a,
915	        generators:=InducedPcgs(home,a) mod NumeratorOfModuloPcgs(nag),
916	        modulePcgs:=nag);
917	    OCOneCocycles(z,true);
918	    if IsBound(z.complement) and
919	      # normal complements exist, iff the coboundaries are trivial
920	      (normal=false or Dimension(z.oneCoboundaries)=0)
921	      then
922	      # now fetch the complements
923
924	      z.factorGens:=z.generators;
925	      coboundbas:=Basis(z.oneCoboundaries);
926	      com:=BaseSteinitzVectors(BasisVectors(Basis(z.oneCocycles)),
927	                               BasisVectors(coboundbas));
928	      field:=LeftActingDomain(z.oneCocycles);
929	      if Size(field)^Length(com.factorspace)>100000 then
930		Info(InfoWarning,1, "Many (",
931		  Size(field)^Length(com.factorspace),") complements!");
932	      fi;
933	      com:=Enumerator(VectorSpace(field,com.factorspace,
934	                                       Zero(z.oneCocycles)));
935	      Info(InfoPcSubgroup,3,"  ",Length(com),
936	           " local complement classes");
937
938	      # compute fusion
939	      kconh:=List([1..Length(com)],i->[i]);
940	      if i<len or retnorm then
941		# we need to compute normalizers
942		comnorms:=[];
943	      else
944		comnorms:=fail;
945	      fi;
946
947	      if Length(com)>1 and Size(a)<Size(bsnorms[bpos]) then
948
949	        opr:=function(cyc,elm)
950		      local l,i;
951			l:=z.cocycleToList(cyc);
952			for i in [1..Length(l)] do
953			  l[i]:=(z.complementGens[i]*l[i])^elm;
954			od;
955			l:=CorrespondingGeneratorsByModuloPcgs(z.origgens,l);
956			for i in [1..Length(l)] do
957			  l[i]:=LeftQuotient(z.complementGens[i],l[i]);
958			od;
959			l:=z.listToCocycle(l);
960			return SiftedVector(coboundbas,l);
961		      end;
962
963		xo:=ExternalOrbitsStabilizers(
964		     ExternalSet(bsnorms[bpos],com,opr));
965
966                for k in xo do
967		  l:=List(k,i->Position(com,i));
968		  if comnorms<>fail then
969		    comnorms[l[1]]:=StabilizerOfExternalSet(k);
970		    isTrueComnorm:=false;
971		  fi;
972		  l:=Set(l);
973		  for kp in l do
974		    kconh[kp]:=l;
975		  od;
976		od;
977
978	      elif comnorms<>fail then
979		if Size(a)=Size(bsnorms[bpos]) then
980		  comnorms:=List(com,i->z.cocycleToComplement(i));
981		  isTrueComnorm:=true;
982		  comnorms:=List(comnorms,
983			      i->ClosureSubgroup(CentralizerModulo(n,b,i),i));
984	        else
985		  isTrueComnorm:=false;
986		  comnorms:=List(com,i->bsnorms[bpos]);
987		fi;
988	      fi;
989
990
991              if Length(efunc)>0 then
992		ncom:=[];
993
994	        #search for invariant ones
995
996		# force exponents corresponding to vector space
997
998                # get matrices for the inner automorphisms
999#		conj:=[];
1000#		for k in GeneratorsOfGroup(a) do
1001#		  mat:=[];
1002#		  for l in nag do
1003#		    Add(mat,One(field)*ExponentsOfPcElement(nag,l^k));
1004#		  od;
1005#		  Add(conj,mat);
1006#		od;
1007                conj:=LinearOperationLayer(a,GeneratorsOfGroup(a),nag);
1008
1009                idmat:=conj[1]^0;
1010		mat:= GroupByGenerators( conj, idmat );
1011		chom:= GroupHomomorphismByImagesNC(a,mat,
1012		        GeneratorsOfGroup(a),conj);
1013
1014		smats:=[];
1015		shoms:=[];
1016
1017                fghom:=Concatenation(z.factorGens,GeneratorsOfGroup(n));
1018		bgids:=List(GeneratorsOfGroup(n),i->One(b));
1019
1020		# now run through the complements
1021		for kp in [1..Length(com)] do
1022
1023		  if kconh[kp]=fail then
1024		    Info(InfoPcSubgroup,3,"already conjugate");
1025		  else
1026
1027		    l:=z.cocycleToComplement(com[kp]);
1028		    # the projection on the complement
1029		    k:=ClosureSubgroup(b,l);
1030		    if Length(s)=1 and IsOne(s[1]) then
1031		      # special case -- no conjugates
1032		      if ForAll(efunc,x->ForAll(GeneratorsOfGroup(l),
1033			   y->ImagesRepresentative(x,y) in k)) then
1034			l:=rec(representative:=k);
1035			if comnorms<>fail then
1036			  if IsBound(comnorms[kp]) then
1037			    l.normalizer:=comnorms[kp];
1038			  else
1039			    l.normalizer:=Normalizer(bsnorms[bpos],
1040				    ClosureSubgroup(b,k));
1041			  fi;
1042			fi;
1043			Add(ncom,l);
1044
1045			# tag all conjugates
1046			for l in kconh[kp] do
1047			  kconh[l]:=fail;
1048			od;
1049		      fi;
1050
1051		    else
1052		      # generic case
1053
1054		      comproj:= GroupHomomorphismByImagesNC(a,a,fghom,
1055				Concatenation(GeneratorsOfGroup(l),bgids));
1056
1057		      # now run through the conjugating elements
1058		      conjnr:=1;
1059		      found:=false;
1060		      while conjnr<=Length(s) and found=false do
1061			if not IsBound(smats[conjnr]) then
1062			  # compute the matrix action for the induced, jugated
1063			  # morphisms
1064			  m:=s[conjnr];
1065			  smats[conjnr]:=[];
1066			  shoms[conjnr]:=[];
1067			  for l in efunc do
1068			    # the induced, jugated morphism
1069			    shom:= GroupHomomorphismByImagesNC(a,a,
1070				    GeneratorsOfGroup(a),
1071				    List(GeneratorsOfGroup(a),
1072				    i->Image(l,i^m)^Inverse(m)));
1073
1074			    mat:=List(nag,
1075				  i->One(field)*ExponentsOfPcElement(nag,
1076				  Image(shom,i)));
1077			    Add(smats[conjnr],mat);
1078			    Add(shoms[conjnr],shom);
1079			  od;
1080			fi;
1081
1082			mats:=smats[conjnr];
1083			# now test whether the complement k can be conjugated to
1084			# be invariant under the morphisms to mats
1085			glsyl:=List(nag,i->[]);
1086			glsyr:=[];
1087			for l in [1..Length(efunc)] do
1088			  kgens:=GeneratorsOfGroup(k);
1089			  for kgnr in [1..Length(kgens)] do
1090
1091			    kgn:=Image(shoms[conjnr][l],kgens[kgnr]);
1092			    kgim:=Image(comproj,kgn);
1093			    Assert(2,kgim^-1*kgn in n);
1094			    # nt part
1095			    kgn:=kgim^-1*kgn;
1096
1097			    # translate into matrix terms
1098			    kgim:=Image(chom,kgim);
1099			    kgn:=One(field)*ExponentsOfPcElement(nag,kgn);
1100
1101			    # the matrix action
1102			    mat:=idmat+(mats[l]-idmat)*kgim-mats[l];
1103
1104			    # store action and vector
1105			    for m in [1..Length(glsyl)] do
1106			      glsyl[m]:=Concatenation(glsyl[m],mat[m]);
1107			    od;
1108			    glsyr:=Concatenation(glsyr,kgn);
1109
1110			  od;
1111			od;
1112
1113			# a possible conjugating element is a solution of the
1114			# large LGS
1115			l:= SolutionMat(glsyl,glsyr);
1116			if l <> fail then
1117			  m:=Product([1..Length(l)],
1118				    i->nag[i]^IntFFE(l[i]));
1119			  # note that we found one!
1120			  found:=[s[conjnr],m];
1121			fi;
1122
1123			conjnr:=conjnr+1;
1124		      od;
1125
1126		      # there is an invariant complement?
1127		      if found<>false then
1128			found:=found[2]*found[1];
1129			l:=ConjugateSubgroup(ClosureSubgroup(b,k),found);
1130			Assert(1,ForAll(efunc,i->Image(i,l)=l));
1131			l:=rec(representative:=l);
1132			if comnorms<>fail then
1133			  if IsBound(comnorms[kp]) then
1134			    l.normalizer:=ConjugateSubgroup(comnorms[kp],found);
1135			  else
1136			    l.normalizer:=ConjugateSubgroup(
1137					    Normalizer(bsnorms[bpos],
1138				    ClosureSubgroup(b,k)), found);
1139			  fi;
1140			fi;
1141			Add(ncom,l);
1142
1143			# tag all conjugates
1144			for l in kconh[kp] do
1145			  kconh[l]:=fail;
1146			od;
1147
1148		      fi;
1149
1150		    fi;
1151
1152                  fi; # if not already a conjugate
1153
1154		od;
1155
1156		# if invariance test needed
1157	      else
1158		# get representatives of the fused complement classes
1159		l:=Filtered([1..Length(com)],i->kconh[i][1]=i);
1160
1161		ncom:=[];
1162		for kp in l do
1163		  m:=rec(representative:=
1164			  ClosureSubgroup(b,z.cocycleToComplement(com[kp])));
1165		  if comnorms<>fail then
1166		    m.normalizer:=comnorms[kp];
1167		  fi;
1168		  Add(ncom,m);
1169		od;
1170	      fi;
1171	      com:=ncom;
1172
1173	      # take the preimages
1174	      for k in com do
1175
1176		Assert(1,ForAll(efunc,i->Image(i,k.representative)
1177		                         =k.representative));
1178		Add(ngrps,k.representative);
1179		if IsBound(k.normalizer) then
1180		  if isTrueComnorm then
1181		    Add(ngrpsnorms,k.normalizer);
1182		  else
1183		    Add(ngrpsnorms,Normalizer(k.normalizer,k.representative));
1184		  fi;
1185		fi;
1186	      od;
1187	    fi;
1188	  fi;
1189	od;
1190
1191      fi;
1192    od;
1193
1194    grps:=ngrps;
1195    grpsnorms:=ngrpsnorms;
1196    Info(InfoPcSubgroup,5,List(grps,Size),List(grpsnorms,Size));
1197  od;
1198
1199  if isom<>fail then
1200    grps:=List(grps,j->PreImage(isom,j));
1201    if retnorm then
1202      grpsnorms:=List(grpsnorms,j->PreImage(isom,j));
1203    fi;
1204  fi;
1205
1206  if retnorm then
1207    return [grps,grpsnorms];
1208  else
1209    return grps;
1210  fi;
1211end);
1212
1213
1214#############################################################################
1215##
1216#M  CharacteristicSubgroups( <G> )
1217##
1218InstallMethod(CharacteristicSubgroups,"solvable, automorphisms",true,
1219  [IsGroup and IsSolvableGroup],0,
1220function(G)
1221local A,s;
1222  if AbelianRank(G)<5 then
1223    TryNextMethod();
1224  fi;
1225  A:=AutomorphismGroup(G);
1226  s:=SubgroupsSolvableGroup(G,rec(normal:=true,actions:=GeneratorsOfGroup(A)));
1227  return Filtered(s,x->IsCharacteristicSubgroup(G,x));
1228end);
1229
1230#############################################################################
1231##
1232#M  LatticeSubgroups(<G>)  . . . . . . . . . .  lattice of subgroups
1233##
1234InstallMethod(LatticeSubgroups,"elementary abelian extension",true,
1235  [IsGroup and IsFinite and CanComputeFittingFree],
1236  # want to be better than cyclic extension.
1237  1,
1238function(G)
1239local s,i,c,classes, lattice,map,GI;
1240
1241  if Size(G)=1 or not IsSolvableGroup(G) then #or not CanEasilyComputePcgs(G) then
1242    TryNextMethod();
1243  fi;
1244  if not IsPcGroup(G) or IsPermGroup(G) then
1245    map:=IsomorphismPcGroup(G);
1246    GI:=Image(map,G);
1247  else
1248    map:=fail;
1249    GI:=G;
1250  fi;
1251  s:=SubgroupsSolvableGroup(GI,rec(retnorm:=true));
1252  classes:=[];
1253  for i in [1..Length(s[1])] do
1254    if map=fail then
1255      c:=ConjugacyClassSubgroups(G,s[1][i]);
1256      SetStabilizerOfExternalSet(c,s[2][i]);
1257    else
1258      c:=ConjugacyClassSubgroups(G,PreImage(map,s[1][i]));
1259      SetStabilizerOfExternalSet(c,PreImage(map,s[2][i]));
1260    fi;
1261    Add(classes,c);
1262  od;
1263  Sort(classes,function(a,b)
1264                 return Size(Representative(a))<Size(Representative(b));
1265	       end);
1266
1267  # create the lattice
1268  lattice:=Objectify(NewType(FamilyObj(classes),IsLatticeSubgroupsRep),
1269		     rec());
1270  lattice!.conjugacyClassesSubgroups:=classes;
1271  lattice!.group     :=G;
1272
1273  # return the lattice
1274  return lattice;
1275
1276end);
1277
1278# #############################################################################
1279# ##
1280# #M  NormalSubgroups(<G>)  . . . . . . . . . .  list of normal subgroups
1281# ##
1282# InstallMethod(NormalSubgroups,"elementary abelian extension",true,
1283#   [CanEasilyComputePcgs],0,
1284# function(G)
1285# local n;
1286#   n:=SubgroupsSolvableGroup(G,rec(
1287#        actions:=List(GeneratorsOfGroup(G),i->InnerAutomorphism(G,i)),
1288#        normal:=true));
1289#
1290#   # sort the normal subgroups according to their size
1291#   Sort(n,function(a,b) return Size(a) < Size(b); end);
1292#
1293#   return n;
1294# end);
1295
1296#############################################################################
1297##
1298#F  SizeConsiderFunction(<size>)  returns auxiliary function for
1299##  'SubgroupsSolvableGroup' that allows one to discard all subgroups whose
1300##  size is not divisible by <size>
1301##
1302InstallGlobalFunction(SizeConsiderFunction,function(size)
1303  return function(c,a,n,b,m)
1304	   return IsInt(Size(a)/Size(n)*Size(b)*Size(m)/size);
1305         end;
1306end);
1307
1308#############################################################################
1309##
1310#F  ExactSizeConsiderFunction(<size>)  returns auxiliary function for
1311##  'SubgroupsSolvableGroup' that allows one to discard all subgroups whose
1312##  size is not <size>
1313##
1314InstallGlobalFunction(ExactSizeConsiderFunction,function(size)
1315  return function(c,a,n,b,m)
1316           local result;
1317	   result:=IsInt(Size(a)/Size(n)*Size(b)*Size(m)/size)
1318	      and not (Size(a)/Size(n)*Size(b))>size;
1319	   return result;
1320         end;
1321end);
1322
1323BindGlobal("ElementaryAbelianConsider",
1324function(c,a,n,b,m)
1325  return HasElementaryAbelianFactorGroup(a,n) and
1326    ForAll(GeneratorsOfGroup(a),x->ForAll(GeneratorsOfGroup(b),y->Comm(x,y)
1327    in m));
1328end);
1329