glzmodmz.gi

#############################################################################
##
#W  glzmodmz.gi                    GAP library                    Stefan Kohl
#W                                                           Alexander Hulpke
##
##
#Y  Copyright (C) 2011 The GAP Group
##
##  This file contains the functionality for constructing clasical groups over
##  residue class rings.
##

#############################################################################
##
#F  SizeOfGLdZmodmZ( <d>, <m> ) . . . . . .  Size of the group GL(<d>,Z/<m>Z)
##
##  Computes the order of the group `GL( <d>, Integers mod <m> )' for
##  positive integers <d> and <m> > 1.
##
InstallGlobalFunction( SizeOfGLdZmodmZ,

  function ( d, m )

    local  size, pow, p, q, k, i;

    if   not (IsPosInt(d) and IsInt(m) and m > 1)
    then Error("GL(",d,",Integers mod ",m,") is not a well-defined group, ",
               "resp. not supported.\n");
    fi;
    size := 1;
    for pow in Collected(Factors(m)) do
      p := pow[1]; k := pow[2]; q := p^k;
      size := size * Product([d*k - d .. d*k - 1], i -> q^d - p^i);
    od;
    return size;
  end );

#############################################################################
##
#M  SpecialLinearGroupCons( IsNaturalSL, <d>, Integers mod <m> )
##
InstallMethod( SpecialLinearGroupCons,
               "natural SL for dimension and residue class ring",
               [ IsMatrixGroup and IsFinite, IsPosInt,
                 IsRing and IsFinite and IsZmodnZObjNonprimeCollection ],

  function ( filter, d, R )

    local  G, gens, g, m, T;

    m := Size(R);
    if R <> Integers mod m or m = 1 then TryNextMethod(); fi;
    if IsPrime(m) then return SpecialLinearGroupCons(IsMatrixGroup,d,m); fi;
    if   d = 1
    then gens := [IdentityMat(d,R)];
    else gens := List(GeneratorsOfGroup(SymmetricGroup(d)),
                      g -> PermutationMat(g,d) * One(R));
         for g in gens do
           if DeterminantMat(g) <> One(R) then g[1] := -g[1]; fi;
         od;
         T := IdentityMat(d,R); T[1][2] := One(R); Add(gens,T);
    fi;
    G := GroupByGenerators(gens);
    SetName(G,Concatenation("SL(",String(d),",Z/",String(m),"Z)"));
    SetIsNaturalSL(G,true);
    SetDimensionOfMatrixGroup(G,d);
    SetIsFinite(G,true);
    SetSize(G,SizeOfGLdZmodmZ(d,m)/Phi(m));
    return G;
  end );

#############################################################################
##
#M  GeneralLinearGroupCons( IsNaturalGL, <d>, Integers mod <m> )
##
InstallMethod( GeneralLinearGroupCons,
               "natural GL for dimension and residue class ring",
               [ IsMatrixGroup and IsFinite, IsPosInt,
                 IsRing and IsFinite and IsZmodnZObjNonprimeCollection ],

  function ( filter, d, R )

    local  G, gens, g, m, T, D;

    m := Size(R);
    if R <> Integers mod m or m = 1 then TryNextMethod(); fi;
    if IsPrime(m) then return GeneralLinearGroupCons(IsMatrixGroup,d,m); fi;
    if   d = 1
    then gens := List(GeneratorsOfGroup(Units(R)), g -> [[g]]);
    else gens := List(GeneratorsOfGroup(SymmetricGroup(d)),
                      g -> PermutationMat(g,d) * One(R));
         T := IdentityMat(d,R); T[1][2] := One(R); Add(gens,T);
         for g in GeneratorsOfGroup(Units(R)) do
           D := IdentityMat(d,R); D[1][1] := g; Add(gens,D);
         od;
    fi;
    G := GroupByGenerators(gens);
    SetName(G,Concatenation("GL(",String(d),",Z/",String(m),"Z)"));
    SetIsNaturalGL(G,true);
    SetDimensionOfMatrixGroup(G,d);
    SetIsFinite(G,true);
    SetSize(G,SizeOfGLdZmodmZ(d,m));
    return G;
  end );

BindGlobal("OrderMatrixIntegerResidue",function(p,a,M)
local f,M2,o,e,MM,i;
  MM:=M;
  f:=GF(p);
  M2:=ImmutableMatrix(f,List(M,x->List(x,y->Int(y)*One(f))));
  o:=Order(M2);
  M:=M^o;
  e:=p;
  i:=1;
  while i<a do
    i:=i+1;
    e:=e*p;
    M2:=M-M^0;
    if ForAny(M2,x->ForAny(x,x->Int(x) mod e<>0)) then
      o:=o*p;
      M:=M^p;
    fi;
  od;
  
  Assert(1,IsOne(M));
  return o;
end);

BindGlobal("SPRingGeneric",function(n,ring)
local t,geni,m,slmats,gens,f,rels,i,j,k,l,mat,mat1,mats,id,nh;
  nh:=n;
  n:=2*n;
  t:=function(i,j)
    return slmats[geni[i][j]];
  end;
  geni:=List([1..n],x->[]);
  mats:=[];
  slmats:=[];
  id:=IdentityMat(n,One(ring));
  m:=0;
  for i in [1..nh] do
    #t_{i,n+i}
    mat:=List(id,ShallowCopy);
    mat[i][nh+i]:=One(ring);
    Add(slmats,mat);
    #t_{n+i,i}
    mat:=List(id,ShallowCopy);
    mat[nh+i][i]:=One(ring);
    Add(slmats,mat);
  od;

  for i in [1..nh] do
    for j in [i+1..nh] do
      # t_{i,n+j}
      mat:=List(id,ShallowCopy);
      mat[i][nh+j]:=One(ring);
      mat1:=mat;
      # t_{j,n+i}
      mat:=List(id,ShallowCopy);
      mat[j][nh+i]:=One(ring);
      Add(slmats,mat1*mat);
      # t_{n+i,j}
      mat:=List(id,ShallowCopy);
      mat[nh+i][j]:=One(ring);
      mat1:=mat;
      # t_{n+j,i}
      mat:=List(id,ShallowCopy);
      mat[nh+j][i]:=One(ring);
      Add(slmats,mat1*mat);
    od;
  od;

  return Group(slmats);
end);

InstallGlobalFunction("ConstructFormPreservingGroup",function(arg)
local oper,n,R,o,nrit,
  q,p,field,zero,one,oner,a,f,pp,b,d,fb,btf,eq,r,i,j,e,k,ogens,gens,gensi,
  bp,sol,
  g,prev,proper,fp,ho,evrels,hom,bas,basm,em,ngens,addmat,sub,transpose;

  oper:=arg[1];
  R:=arg[Length(arg)];
  n:=arg[Length(arg)-1];
  q:=Size(R);
  if not IsPrimePowerInt(q) then
    TryNextMethod();
  fi;
  p:=Factors(q)[1];
  if p=2 then 
    if oper=SP then
      return SPRingGeneric(n/2,R);
    else
      return fail;
    fi;
  fi;
  field:=GF(p);
  zero:=Zero(field);
  one:=One(field);
  if Length(arg)=3 then
    g:=oper(n,p);
  else
    g:=oper(arg[2],n,p);
  fi;
    
  # get the form and get the correct -1's
  f:=InvariantBilinearForm(g).matrix;

  transpose:=not ForAll(GeneratorsOfGroup(g),
    x->TransposedMat(x)*f*x=f);
  if transpose then
    Info(InfoGroup,1,"transpose!");
    if HasSize(g) then
      e:=Size(g);
    else
      e:=fail;
    fi;
    g:=Group(List(GeneratorsOfGroup(g),TransposedMat));
    if e<>fail then
      SetSize(g,e);
    fi;
  fi;
  #IsomorphismFpGroup(g); # force hom for next steps
  f:=List(f,r->List(r,Int));
  for i in [1..n] do

    for j in [1..n] do
      if f[i][j]=p-1 then
	f[i][j]:=-1;
      fi;
    od;
  od;

  nrit:=0;
  pp:=p; # previous p
  while pp<q do

    nrit:=nrit+1;
    prev:=g;

    if HasIsomorphismFpGroup(prev) then
      hom:=IsomorphismFpGroup(prev);
      fp:=Range(hom);
      ogens:=List(GeneratorsOfGroup(fp),
	      x->List(PreImagesRepresentative(hom,x)));
    else
      fp:=fail;
      ogens:=GeneratorsOfGroup(prev);
    fi;
    ogens:=List(ogens,x->List(x,r->List(r,Int)));
    gens:=[];

    for bp in [1..Length(ogens)+1] do
      if bp<=Length(ogens) then
	b:=ogens[bp];
      else
	b:=One(ogens[1]);
      fi;
      d:=(TransposedMat(b)*f*b-f)*1/pp;
      # solve  D+E^T*F*B+B^T*F*E=0
      fb:=f*b;
      btf:=TransposedMat(b)*f;
      eq:=[];
      r:=[];
      for i in [1..n] do
	for j in [1..n] do
	  # eq for entry i,j
	  e:=ListWithIdenticalEntries(n^2,zero);
	  for k in [1..n] do
	    e[(k-1)*n+i]:=e[(k-1)*n+i]+fb[k][j];
	    e[(k-1)*n+j]:=e[(k-1)*n+j]+btf[i][k];
	  od;
	  Add(eq,e);

	  #RHS is -d entry
	  Add(r,-d[i][j]*one);
	od;
      od;
      eq:=TransposedMat(eq); # columns were corresponding to variables

      if bp<=Length(ogens) then
	# lift generator
	sol:=SolutionMat(eq,r);

	# matrix from it
	sol:=List([1..n],x->sol{[(x-1)*n+1..x*n]});
	sol:=List(sol,x->List(x,Int));
	Add(gens,b+pp*sol);
      else
	# we know all gens

	oner:=One(Integers mod (pp*p));
	gens:=List(gens,x->x*oner);

	g:=Group(gens);

	# d will be zero, so homogeneous

	sol:=NullspaceMat(eq);
	#Info(InfoGroup,1,"extend by dim",Length(sol));

	proper:=p^Length(sol)*Size(prev); # proper order of group

	if ValueOption("avoidkerneltest")<>true then
	  # vector space in kernel that is generated
	  bas:=[];
	  basm:=[];
	  sub:=VectorSpace(field,bas,Zero(e));

	  addmat:=function(em)
	  local c;
	    e:=List(em,r->List(r,Int))-b;
	    e:=1/pp*e;
	    e:=Concatenation(e)*one;
	    if p<257 then
	      ConvertToVectorRep(e,p);
	    fi;
	    if not e in sub then
	      Add(bas,e);
	      Add(basm,em);
	      sub:=VectorSpace(field,bas);
	    fi;
	  end;

	  if fp<>fail then
	    # evaluate relators
	    evrels:=RelatorsOfFpGroup(fp);

	    i:=1;
	    while i<=Length(evrels) and Length(bas)<Length(sol) do
	      em:=MappedWord(evrels[i],FreeGeneratorsOfFpGroup(fp),gens);
	      addmat(em);
	      i:=i+1;
	    od;
	  else
	    evrels:=Source(EpimorphismFromFreeGroup(prev));
	    repeat
	      j:=PseudoRandom(evrels:radius:=10);
	      k:=MappedWord(j,GeneratorsOfGroup(evrels),GeneratorsOfGroup(prev));
	      o:=OrderMatrixIntegerResidue(p,nrit,k);
	      k:=MappedWord(j,GeneratorsOfGroup(evrels),gens)^o;
	    until not IsOne(k);
	    addmat(k);
	      
	  fi;

	  # close under action
	  gensi:=List(gens,Inverse);
	  i:=1;
	  while i<=Length(basm) and Length(bas)<Length(sol) do
	    for j in [1..Length(gens)] do
	      #em:=basm[i]^j;
	      em:=gensi[j]*basm[i]*gens[j];
	      addmat(em);
	    od;
	    i:=i+1;
	  od;

	  if Length(bas)=Length(sol) then
	    Info(InfoGroup,1,"kernel generated ",Length(bas));
	  else
	    Info(InfoGroup,1,"kernel partially generated ",Length(bas));
	    ngens:=ShallowCopy(gens);
	    i:=Iterator(sol); # just run through basis as linear
	    while Length(bas)<Length(sol) do
	      e:=NextIterator(i);
	      e:=List(e,Int);
	      e:=b+pp*List([1..n],x->e{[(x-1)*n+1..x*n]});
	      addmat(e);
	      if e=basm[Length(basm)] then
		# was added
		Add(ngens,e);
		g:=Group(ngens);
		Info(InfoGroup,1,"added generator");
	      fi;
	    od;
	  fi;

	  if fp <>fail then
	    # extend presentation
	    bas:=Basis(sub,bas);
	    RUN_IN_GGMBI:=true;
	    hom:=GroupGeneralMappingByImagesNC(g,fp,gens,GeneratorsOfGroup(fp));
	    hom:=LiftFactorFpHom(hom,g,"M",SubgroupNC(g,basm),rec(
		  pcgs:=basm,
		  prime:=p,
		  decomp:=function(em)
		  local e;
		    e:=List(em,r->List(r,Int))-b;
		    e:=1/pp*e;
		    e:=Concatenation(e)*one;
		    return List(Coefficients(bas,e),Int);
		  end
		  ));
	    RUN_IN_GGMBI:=false;
	    #simplify Image to avoid explosion of generator number
	    fp:=Range(hom);
	    if true then
	      # remove redundant generators
	      e:=PresentationFpGroup(fp);
	      TzOptions(e).printLevel:=0;
	      j:=Filtered(Reversed([1..Length(e!.generators)]),
		x->not MappingGeneratorsImages(hom)[1][x] in ngens);
	      j:=e!.generators{j};

	      TzInitGeneratorImages(e);
	      for i in j do
		TzEliminate(e,i);
	      od;
	      fp:=FpGroupPresentation(e);
	      j:=MappingGeneratorsImages(hom);
	      k:=TzPreImagesNewGens(e);
	      k:=List(k,x->j[1][Position(OldGeneratorsOfPresentation(e),x)]);

	      RUN_IN_GGMBI:=true;
	      hom:=GroupHomomorphismByImagesNC(g,fp,
		    k,
		    GeneratorsOfGroup(fp));
	      RUN_IN_GGMBI:=false;
	    fi;

	    SetIsomorphismFpGroup(g,hom);
	  fi;
	fi;

	SetSize(g,Size(prev)*Size(field)^Length(sol));
      fi;

    od;

    pp:=pp*p;
  od;

  if transpose then
    e:=Size(g);
    g:=Group(List(GeneratorsOfGroup(g),TransposedMat));
    SetSize(g,e);
  fi;
  SetInvariantBilinearForm(g,rec(matrix:=f*oner));
  
  return g;
end);

#############################################################################
##
#M  SymplecticGroupCons( <IsMatrixGroup>, <d>, Integers mod <q> )
##
InstallOtherMethod( SymplecticGroupCons,
  "symplectic group for dimension and residue class ring for prime powers",
  [ IsMatrixGroup and IsFinite, IsPosInt,
    IsRing and IsFinite and IsZmodnZObjNonprimeCollection ],
function ( filter, n, R )
local g;
  g:=ConstructFormPreservingGroup(SP,n,R);
  SetName(g,Concatenation("Sp(",String(n),",Z/",String(Size(R)),"Z)"));
  return g;
end);

#############################################################################
##
#M  GeneralOrthogonalGroupCons ( <IsMatrixGroup>, <d>, Integers mod <q> )
##
InstallOtherMethod( GeneralOrthogonalGroupCons,
  "GO for dimension and residue class ring for prime powers",
  [ IsMatrixGroup and IsFinite, IsInt,IsPosInt,
    IsRing and IsFinite and IsZmodnZObjNonprimeCollection ],
function ( filter, sign,n, R )
local g;
  if sign=0 then
    g:=ConstructFormPreservingGroup(GO,n,R);
    SetName(g,Concatenation("GO(",String(n),",Z/",String(Size(R)),"Z)"));
  else
    g:=ConstructFormPreservingGroup(GO,sign,n,R);
    SetName(g,Concatenation("GO(",String(sign),",",String(n),
      ",Z/",String(Size(R)),"Z)"));
  fi;
  return g;
end);

#############################################################################
##
#M  SpecialOrthogonalGroupCons( <IsMatrixGroup>, <d>, Integers mod <q> )
##
InstallOtherMethod( SpecialOrthogonalGroupCons,
  "GO for dimension and residue class ring for prime powers",
  [ IsMatrixGroup and IsFinite, IsInt,IsPosInt,
    IsRing and IsFinite and IsZmodnZObjNonprimeCollection ],
function ( filter, sign,n, R )
local g;
  if sign=0 then
    g:=ConstructFormPreservingGroup(SO,n,R);
    if g=fail then TryNextMethod();fi;
    SetName(g,Concatenation("SO(",String(n),",Z/",String(Size(R)),"Z)"));
  else
    g:=ConstructFormPreservingGroup(SO,sign,n,R);
    if g=fail then TryNextMethod();fi;
    SetName(g,Concatenation("SO(",String(sign),",",String(n),
      ",Z/",String(Size(R)),"Z)"));
  fi;
  return g;
end);

#############################################################################
##
#E  glzmodmz.gi . . . . . . . . . . . . . . . . . . . . . . . . . . ends here