From 941ccd051fa0bdd11d76f5fce82c95a0750da7cd Mon Sep 17 00:00:00 2001 From: Daniel Rademacher Date: Thu, 13 Jun 2024 14:59:58 +0200 Subject: [PATCH 01/13] Include SL Code. --- .../constructive_recognition/SL/BaseCase.gi | 477 ++++ .../constructive_recognition/SL/GoingDown.gi | 697 ++++++ .../constructive_recognition/SL/GoingUp.gi | 547 +++++ .../constructive_recognition/SL/main.gi | 245 ++ .../SL/sl2_BlackBox.gi | 637 ++++++ .../constructive_recognition/main.gi | 31 + .../utils/achieve.gi} | 996 +-------- .../constructive_recognition/utils/utils.gi | 1986 +++++++++++++++++ read.g | 20 +- 9 files changed, 4679 insertions(+), 957 deletions(-) create mode 100644 gap/projective/constructive_recognition/SL/BaseCase.gi create mode 100644 gap/projective/constructive_recognition/SL/GoingDown.gi create mode 100644 gap/projective/constructive_recognition/SL/GoingUp.gi create mode 100644 gap/projective/constructive_recognition/SL/main.gi create mode 100644 gap/projective/constructive_recognition/SL/sl2_BlackBox.gi create mode 100644 gap/projective/constructive_recognition/main.gi rename gap/projective/{classicalnatural.gi => constructive_recognition/utils/achieve.gi} (71%) create mode 100644 gap/projective/constructive_recognition/utils/utils.gi diff --git a/gap/projective/constructive_recognition/SL/BaseCase.gi b/gap/projective/constructive_recognition/SL/BaseCase.gi new file mode 100644 index 00000000..48097c37 --- /dev/null +++ b/gap/projective/constructive_recognition/SL/BaseCase.gi @@ -0,0 +1,477 @@ +############################################################################# +## +## This file is part of recog, a package for the GAP computer algebra system +## which provides a collection of methods for the constructive recognition +## of groups. +## +## This files's authors include Daniel Rademacher. +## +## Copyright of recog belongs to its developers whose names are too numerous +## to list here. Please refer to the COPYRIGHT file for details. +## +## SPDX-License-Identifier: GPL-3.0-or-later +## +############################################################################# + + + +################################################################################################### +################################################################################################### +######## Constructive recognition SL(2,q) in natural representation ############################### +################################################################################################### +################################################################################################### + + + +# TODO: Add MSLPs to these algorithms +# This algorithm is an implementation of the paper +# "Fast Recognition of Classical Groups over Large Fields" +# by Marston Conder and Charles Leedham-Green +RECOG.ConstructiveRecognitionSL2NaturalRepresentation := function(G, q, epsilon) +local F, factors, counter, one, factor, foundEle, passed, max, eigenvalues, eigenvalues2, eigenvectors, eigenvectors2, A, B, rand, foundConjugate, test, BDiag, C, basechange, v, d, CC, check, i, j, D, DD, T, S, o, t1, t2, zero, slp; + + factors := PrimeDivisors(q-1); + counter := 1; + F := GF(q); + + if not(IsObjWithMemory(GeneratorsOfGroup(G)[1])) then + G := GroupWithMemory(G); + fi; + + one := One(G); + zero := Zero(F); + max := Int(1/epsilon); + while counter < max do + + # Step 1: Find random element of order q-1 + # foundEle := false; + # while not(foundEle) do + # A := PseudoRandom(G); + # if A^(q-1) = one then + # passed := true; + # for factor in factors do + # if A^((q-1)/factor) = one then + # passed := false; + # break; + # fi; + # od; + # if passed then + # foundEle := true; + # else + # counter := counter + 1; + # fi; + # else + # counter := counter +1; + # fi; + # if counter >= max then + # return fail; + # fi; + #od; + + # Step 1: Find random element of order q-1 (first improvement based on magma code) + foundEle := false; + while not(foundEle) do + A := PseudoRandom(G); + o := Order (A); + if (o mod (q - 1) = 0) then + A := A^(o/(q - 1)); + foundEle := true; + else + counter := counter +1; + fi; + if counter >= max then + return fail; + fi; + od; + + # Step 2: Eigenvectors and eigenvalues + eigenvalues := RootsOfPolynomial(CharacteristicPolynomial(A)); + # Note eigenvalues[1] = eigenvalues[2]^(-1) + eigenvectors := [1,2]; + eigenvectors[1] := RECOG.EigenspaceMat(StripMemory(A), eigenvalues[1]); + eigenvectors[2] := RECOG.EigenspaceMat(StripMemory(A), eigenvalues[2]); + + # Step 3: Conjugate of A that does not intersect with A's eigenspaces + foundConjugate := false; + eigenvectors2 := [1,2]; + while not(foundConjugate) do + rand := PseudoRandom(G); + B := A^rand; + eigenvectors2[1] := RECOG.EigenspaceMat(StripMemory(B), eigenvalues[1]); + eigenvectors2[2] := RECOG.EigenspaceMat(StripMemory(B), eigenvalues[2]); + test := []; + Append(test,eigenvectors[1]); + Append(test,eigenvectors2[1]); + if NullspaceMat(test) = [] then + test := []; + Append(test,eigenvectors[2]); + Append(test,eigenvectors2[2]); + if NullspaceMat(test) = [] then + basechange := [eigenvectors2[1,1],eigenvectors2[2,1]]; + t1 := eigenvectors[1,1] * basechange^(-1); + t2 := eigenvectors[2,1] * basechange^(-1); + if (t1[1] <> zero) and (t1[2] <> zero) and (t2[1] <> zero) and (t2[2] <> zero) then + foundConjugate := true; + else + counter := counter + 1; + fi; + else + counter := counter + 1; + fi; + else + counter := counter + 1; + fi; + + if counter >= max then + return fail; + fi; + od; + + # Step 4: Find suitable C + BDiag := IdentityMat(2,GF(q)); + BDiag[1,1] := eigenvalues[1]; + BDiag[2,2] := eigenvalues[2]; + v := eigenvectors[1,1] * basechange^(-1); + d := v[2] * v[1]^-1; + + S := IdentityMat(2,F); + while S = one do + foundEle := false; + while not(foundEle) do + C := PseudoRandom(G); + CC := basechange * C * basechange^(-1); + if ((CC[1,2]-d*CC[1,1])) <> zero and ((d^2*CC[2,1] - d * CC[2,2]) <> zero) then + check := (d^2*CC[2,1] - d * CC[2,2])/(CC[1,2]-d*CC[1,1]); + if not(Zero(F) = check) and IsSquareFFE(F,check) then + i := LogFFE(check,eigenvalues[1])/2; + foundEle := true; + else + counter := counter + 1; + fi; + else + counter := counter + 1; + fi; + + if counter >= max then + return fail; + fi; + od; + S := A^(-1) * (B^i*C)^(-1) * A * (B^i*C); + od; + + # Step 5: Find suitable D + v := eigenvectors[2,1] * basechange^(-1); + d := v[2] * v[1]^-1; + + T := IdentityMat(2,F); + while T = one do + foundEle := false; + while not(foundEle) do + D := PseudoRandom(G); + DD := basechange * D * basechange^(-1); + if ((DD[1,2]-d*DD[1,1]) <> zero) and ((d^2*DD[2,1] - d * DD[2,2]) <> zero) then + check := (d^2*DD[2,1] - d * DD[2,2])/(DD[1,2]-d*DD[1,1]); + if not(Zero(F) = check) and IsSquareFFE(F,check) then + j := LogFFE(check,eigenvalues[1])/2; + foundEle := true; + else + counter := counter + 1; + fi; + else + counter := counter + 1; + fi; + + if counter >= max then + return fail; + fi; + od; + T := A^(-1) * (B^j*D)^(-1) * A * (B^j*D); + od; + + basechange := [eigenvectors[1,1],eigenvectors[2,1]]; + basechange[2] := basechange[2] * Determinant(basechange)^(-1); + + slp := SLPOfElms([A,S,T]); + + A := basechange * A * basechange^(-1); + S := basechange * S * basechange^(-1); + T := basechange * T * basechange^(-1); + + return [[A,S,T],basechange,slp]; + od; + + return fail; + +end; + + + +# Note that we use the discrete logarithm to normalise the primitive element at position [1,1]. But this is not necessarly as the entry at position [1,1] is primitive. +# Hence, this function can be adapted to larger fields by avoiding the normalisation step +RECOG.ConstructiveRecognitionSL2NaturalRepresentationCompleteBasis := function(list,F,q,p,f) +local w, k, Diag, coeffs, coeff, cong, t, s, SC, i, res, res2, A, S, T, upper, lower, basis, slp; + + list := GeneratorsWithMemory(list); + + A := list[1]; + S := list[2]; + T := list[3]; + + # Normalisation step + w := PrimitiveElement(F); + k := LogFFE(w,A[1,1]); + Diag := A^k; + + t := T[1,2]; + basis := [1..f]; + for i in [0..f-1] do + basis[i+1] := w^(2*i); + od; + basis := Basis(F,basis); + coeffs := Coefficients(basis,t^(-1)); + + # standard basis element [1, 1, 0, 1] + res := A^0; + for i in [0..f-1] do + coeff := Int(coeffs[i+1]); + cong := Diag^(-i); + res := res * ((T^(cong))^coeff); + od; + upper := [1..f]; + upper[1] := res; + + # set up complete basis for upper triangular matrices + #UB := [GL(2, E) ! [1, x^(2 * i), 0, 1]: i in [0..e - 1]]; + #wUB := [wUS^(wD^-i): i in [0..e - 1]]; + for i in [1..f-1] do + upper[i+1] := res^(Diag^(-i)); + od; + + s := S[2,1]; + coeffs := Coefficients(basis,s^(-1)); + + # standard basis element [1, 0, 1, 1] + res2 := A^0; + for i in [0..f-1] do + coeff := Int(coeffs[i+1]); + cong := Diag^(i); + res2 := res2 * ((S^(cong))^coeff); + od; + lower := [1..f]; + lower[1] := res2; + + # set up complete basis for lower triangular matrices + for i in [1..f-1] do + lower[i+1] := res2^(Diag^(i)); + od; + + slp := SLPOfElms(Concatenation([Diag],upper,lower)); + + return [[Diag,upper,lower],slp]; + +end; + + +################################################################################################### +################################################################################################### +######## Rewriting SL(2,q) ######################################################################## +################################################################################################### +################################################################################################### + + +# Note that we use the discrete logarithm to normalise the primitive element at position [1,1]. But this is not necessarly as the entry at position [1,1] is primitive. +# Hence, this function can be adapted to larger fields by avoiding the normalisation step +RECOG.RewritingSL2 := function(gens,F,q,p,f,ele) +local list, re, ell, mat, base, i, coeff, wMat, l1, l2, stdgens; + + stdgens := StripMemory(gens); + stdgens := Concatenation([stdgens[1]],stdgens[2],stdgens[3]); + stdgens := GeneratorsWithMemory(stdgens); + l1 := stdgens{[2..f+1]}; + l2 := stdgens{[f+2..Size(stdgens)]}; + stdgens := [stdgens[1],l1,l2]; + + base := [1..f]; + for i in [1..f] do + base[i] := (stdgens[2,i])[1,2]; + od; + base := Basis(F,base); + re := stdgens[1]^0; + + if ele[1,2] = Zero(F) then + re := re*stdgens[2,1]; + mat := IdentityMat(2,F); + mat[1,2] := One(F); + ele := ele*mat; + fi; + + if not(ele[1,1] = One(F)) then + ell := (1-ele[1,1])*ele[1,2]^(-1); + mat := IdentityMat(2,F); + mat[2,1] := ell; + ele := ele*mat; + coeff := Coefficients(base,ell); + wMat := stdgens[1]^0; + for i in [1..f] do + wMat := wMat*(stdgens[3,i])^Int(coeff[i]); + od; + re := re*wMat; + fi; + + ell := -1*ele[1,2]; + mat := IdentityMat(2,F); + mat[1,2] := ell; + ele := ele*mat; + coeff := Coefficients(base,ell); + wMat := stdgens[1]^0; + for i in [1..f] do + wMat := wMat*(stdgens[2,i])^Int(coeff[i]); + od; + re := re*wMat; + + if not(ele[2,1] = Zero(F)) then + ell := -1 * ele[2,1]; + mat := IdentityMat(2,F); + mat[2,1] := ell; + ele := mat*ele; + coeff := Coefficients(base,ell); + wMat := stdgens[1]^0; + for i in [1..f] do + wMat := wMat*(stdgens[3,i])^Int(coeff[i]); + od; + re := re*wMat; + fi; + + return SLPOfElm(re^(-1)); +end; + + + +################################################################################################### +################################################################################################### +######## Constructive Recognition of SL(4,q) (Leedham-Green and O'Brien algorithm) ################ +################################################################################################### +################################################################################################### + + + +# Input: X where is isomorphic to SL(4,q), F where X are dxd matrices over F_q = F (q = p^f prime power) +RECOG.OneEvenSL4 := function(X,F) + local d, G, g, h, foundStrongInvoluation, order, gensCentraliser, EPlus; + + d := 4; + G := GroupByGenerators(X); + foundStrongInvoluation := false; + while not(foundStrongInvoluation) do + g := PseudoRandom(G); + order := Order(g); + if (order mod 2 = 0) then + h := g^(order/2); + EPlus := RECOG.FixspaceMat(h); + if Size(EPlus) = 2 then + foundStrongInvoluation := true; + fi; + fi; + od; + + gensCentraliser := RECOG.InvolutionCentraliser(G,h,100); + # TODO: Compute generating set for OmegaX(E) (see paper from LGO) + + #TODO: CONTINUE HERE + +end;; + + + +# Input: h an involuation in a BB group G, natural number N > 0 +RECOG.InvolutionCentraliser := function(G, h, N) + local C, i, g; + + C := [1..N]; + for i in [1..N] do + g := PseudoRandom(G); + C[i] := RECOG.ChFromg(h,g); + od; + + return DuplicateFreeList(C); +end;; + + + +# Input: h and g group elements of the same group. Returns an element as in Bray's Lemma +RECOG.ChFromg := function(h,g) + local order, com; + + com := h^(-1)*g^(-1)*h*g; + order := Order(com); + + if (order mod 2 = 0) then + return com^(order/2); + else + return com^((order+1)/2)*g^(-1); + fi; +end;; + + + +################################################################################################### +################################################################################################### +######## Older algorithms ######################################################################### +################################################################################################### +################################################################################################### + + + +RECOG.RecogniseSL2NaturalOddCharUsingBSGS := function(g,f) + local ext,p,q,res,slp,std; + p := Characteristic(f); + ext := DegreeOverPrimeField(f); + q := Size(f); + std := RECOG.MakeSL_StdGens(p,ext,2,2); + slp := RECOG.FindStdGensUsingBSGS(g,std.all,false,true); + if slp = fail then return fail; fi; + res := rec( g := g, one := One(f), One := One(g), f := f, q := q, + p := p, ext := ext, d := 2, bas := IdentityMat(2,f), + basi := IdentityMat(2,f) ); + res.all := ResultOfStraightLineProgram(slp,GeneratorsOfGroup(g)); + res.s := res.all{[1..ext]}; + res.t := res.all{[ext+1..2*ext]}; + res.a := res.all[2*ext+1]; + res.b := res.all[2*ext+2]; + return res; +end; + + + +RECOG.FindStdGensUsingBSGS := function(g,stdgens,projective,large) + # stdgens generators for the matrix group g + # returns an SLP expressing stdgens in the generators of g + # set projective to true for projective mode + # set large to true if we should not bother finding nice base points! + local S,dim,gens,gm,i,l,strong; + dim := DimensionOfMatrixGroup(g); + if IsObjWithMemory(GeneratorsOfGroup(g)[1]) then + gm := GroupWithMemory(StripMemory(GeneratorsOfGroup(g))); + else + gm := GroupWithMemory(g); + fi; + if HasSize(g) then SetSize(gm,Size(g)); fi; + if large then + S := StabilizerChain(gm,rec( Projective := projective, + Cand := rec( points := One(g), + ops := ListWithIdenticalEntries(dim, OnLines) ) ) ); + else + S := StabilizerChain(gm,rec( Projective := projective ) ); + fi; + strong := ShallowCopy(StrongGenerators(S)); + ForgetMemory(S); + l := List(stdgens,x->SiftGroupElementSLP(S,x)); + gens := EmptyPlist(Length(stdgens)); + for i in [1..Length(stdgens)] do + if not l[i].isone then + return fail; + fi; + Add(gens,ResultOfStraightLineProgram(l[i].slp,strong)); + od; + return SLPOfElms(gens); +end; diff --git a/gap/projective/constructive_recognition/SL/GoingDown.gi b/gap/projective/constructive_recognition/SL/GoingDown.gi new file mode 100644 index 00000000..b425fddc --- /dev/null +++ b/gap/projective/constructive_recognition/SL/GoingDown.gi @@ -0,0 +1,697 @@ +############################################################################# +## +## This file is part of recog, a package for the GAP computer algebra system +## which provides a collection of methods for the constructive recognition +## of groups. +## +## This files's authors include Daniel Rademacher. +## +## Copyright of recog belongs to its developers whose names are too numerous +## to list here. Please refer to the COPYRIGHT file for details. +## +## SPDX-License-Identifier: GPL-3.0-or-later +## +############################################################################# + + + +############################################################################# +############################################################################# +######## GoingDown method for special linear groups ######################### +############################################################################# +############################################################################# + + + +# TODO: Work on comments and documentation + + + +# Find random element s = r^PseudRandom(g) such that is isomorphic to SL(4,q) +# and check whether they are isomorphic +RECOG.SLn_constructsl4:=function(g,dim,q,r) + local s,h,count,readydim4,readydim3,ready,u,orderu, + nullr,nulls,nullspacer,nullspaces,int,intbasis,nullintbasis, + newu,newbasis,newbasisinv,newr,news,outputu,mat,i,shorts,shortr; + nullr:=RECOG.FixspaceMat(r); + nullspacer:=VectorSpace(GF(q),nullr); + mat:=One(r); + ready:=false; + repeat + s:=r^PseudoRandom(g); + nulls:=RECOG.FixspaceMat(s); + nullspaces:=VectorSpace(GF(q),nulls); + int:=Intersection(nullspacer,nullspaces); + intbasis:=Basis(int); + newbasis:=[]; + for i in [1..Length(intbasis)] do + Add(newbasis,intbasis[i]); + od; + i:=0; + repeat + i:=i+1; + if not mat[i] in int then + Add(newbasis,mat[i]); + int:=VectorSpace(GF(q),newbasis); + fi; + until Dimension(int)=dim; + ConvertToMatrixRep(newbasis); + newbasisinv:=newbasis^(-1); + newr:=newbasis*r*newbasisinv; + news:=newbasis*s*newbasisinv; + + #shortr, shorts do not need memory + #we shall throw away the computations in h + #check that we have SL(4,q), by non-constructive recognition + + # Remark D.R.: Tries to reduce matrix multiplications + # by working with 4 dimensional matrices + shortr:=newr{[dim-3..dim]}{[dim-3..dim]}; + shorts:=news{[dim-3..dim]}{[dim-3..dim]}; + h:=Group(shortr,shorts); + count:=0; + readydim4:=false; + readydim3:=false; + repeat + u:=PseudoRandom(h); + orderu:=Order(u); + if orderu mod ((q^4-1)/(q-1)) = 0 then + readydim4:=true; + elif Gcd(orderu,(q^2+q+1)/Gcd(3,q-1))>1 then + readydim3:=true; + fi; + if readydim4 = true and readydim3 = true then + ready:=true; + break; + fi; + count:=count+1; + until count=30; + until ready=true; + + return Group(r,s); +end; + + + +#g=SL(d,q), given as a subgroup of SL(dim,q) +#output: [SL(2,q), and a basis for the 2-dimensional subspace where it acts +RECOG.SLn_godownfromd:=function(g,q,d,dim) + local y,yy,ready,order,es,dims,subsp,z,x,a,b,c,h,vec,vec2, + pol,factors,degrees,comm1,comm2,comm3,image,basis,action,vs,readyqpl1, + readyqm1,count,u,orderu; + + repeat + ready:=false; + y:=PseudoRandom(g); + pol:=CharacteristicPolynomial(y); + factors:=Factors(pol); + degrees:=List(factors,Degree); + if d-1 in degrees then + order:=Order(y); + if order mod (q-1)=0 then + yy:=y^(order/(q-1)); + else + yy:=One(y); + fi; + if not IsOne(yy) then + es:= Eigenspaces(GF(q),yy); + dims:=List(es,Dimension); + if IsSubset(Set([1,d-1,dim-d]),Set(dims)) and + (1 in Set(dims)) then + es:=Filtered(es,x->Dimension(x)=1); + vec:=Basis(es[1])[1]; + if vec*yy=vec then + vec:=Basis(es[2])[1]; + fi; + repeat + z:=PseudoRandom(g); + x:=yy^z; + a:=Comm(x,yy); + b:=a^yy; + c:=a^x; + comm1:= Comm(a,c); + comm2:=Comm(a,b); + comm3:=Comm(b,c); + if comm1<>One(a) and comm2<>One(a) and + comm3<>One(a) and Comm(comm1,comm2)<>One(a) then + vec2:=vec*z; + vs:=VectorSpace(GF(q),[vec,vec2]); + basis:=Basis(vs); + #check that the action in 2 dimensions is SL(2,q) + #by non-constructive recognition, finding elements of + #order (q-1) and (q+1) + #we do not need memory in the group image + action:=List([a,b,c],x->RECOG.LinearAction(basis,q,x)); + image:=Group(action); + count:=0; + readyqpl1:=false; + readyqm1:=false; + repeat + u:=PseudoRandom(image); + orderu:=Order(u); + if orderu = q-1 then + readyqm1:=true; + elif orderu = q+1 then + readyqpl1:=true; + fi; + if readyqm1 = true and readyqpl1 = true then + ready:=true; + break; + fi; + count:=count+1; + until count=20; + fi; + until ready=true; + fi; + fi; + fi; + until ready; + + h:=Group(a,b,c); + subsp:=VectorSpace(GF(q),[vec,vec2]); + Info(InfoRecog,2,"New Dimension: 2"); + return [h,subsp]; + +end; + + + +#going down from 4 to 2 dimensions, when q=2,3,4,5,9 +#just construct the 4-dimensional invariant space and generators +#for the group acting on it +RECOG.SLn_exceptionalgodown:=function(h,q,dim) + local basis, v, vs, i, gen; + + vs:=VectorSpace(GF(q),One(h)); + basis:=[]; + repeat + if InfoLevel(InfoRecog) >= 3 then Print("C"); fi; + for i in [1..4] do + v:=PseudoRandom(vs); + for gen in GeneratorsOfGroup(h) do + Add(basis,v*gen-v); + od; + od; + basis:=ShallowCopy(SemiEchelonMat(basis).vectors); + until Length(basis)=4; + Info(InfoRecog,2,"New Dimension: 2"); + return [h,VectorSpace(GF(q),basis)]; +end; + + + +RECOG.SLn_constructsl2:=function(g,d,q) + local r,h; + + r := RECOG.constructppdTwoStingray(g,d,q,"SL",fail); + Info(InfoRecog,2,"Finished main GoingDown, i.e. we found a stringray element which operates irreducible on a 2 dimensional subspace. \n"); + Info(InfoRecog,2,"Next goal: Find an random element s.t. the two elements generate SL(4,q). \n"); + h := RECOG.SLn_constructsl4(g,d,q,r); + # Remark D.R.: at this point we know that h is isomorphic to SL(4,q) + Info(InfoRecog,2,"Succesful. "); + Info(InfoRecog,2,"Current Dimension: 4\n"); + Info(InfoRecog,2,"Next goal: Generate SL(2,q). \n"); + if not (q in [2,3,4,5,9]) then + return RECOG.SLn_godownfromd(h,q,4,d); + else + return RECOG.SLn_exceptionalgodown(h,q,d); + fi; +end; + + + +############################################################################## +# The going down method while constructing smaller matrices: +############################################################################## + + + +RECOG.CheckNewStingrayGroupSmallerMatrices := function(g1,dim1,base1,eigenspace1,g2,dim2,base2,eigenspace2,q) + local baseSum, b, action1, action2, fld, module, eigenspaceintersection; + + baseSum := ShallowCopy(base1); + Append(baseSum,base2); + + if NullspaceMat(baseSum) <> [] then + return [false,[]]; + fi; + + fld := GF(q); + b := Basis(VectorSpace(fld,baseSum),baseSum); + action1 := List(baseSum,v->Coefficients(b,v*g1)); + action2 := List(baseSum,v->Coefficients(b,v*g2)); + module := GModuleByMats( [action1,action2], fld ); + if MTX.IsIrreducible(module) then + eigenspaceintersection := SumIntersectionMat(eigenspace1,eigenspace2)[2]; + return [true,[action1,action2],BasisVectors(b),eigenspaceintersection]; + else + return [false,[]]; + fi; +end; + + + +RECOG.ConstructSL4SmallerMatrices := function(g1,g2,q) + local baseSum, b, action1, action2, fld, module, base, EleBase, fixbase, ele, eigenspace1, eigenspace2, eigenspaceintersection; + + base := []; + fld := GF(q); + for ele in [g1,g2] do + fixbase := RECOG.FixspaceMat(TransposedMat(ele)); + EleBase := NullspaceMat(TransposedMat(fixbase)); + Append(base,EleBase); + od; + + eigenspace1 := RECOG.FixspaceMat(StripMemory(g1)); + eigenspace2 := RECOG.FixspaceMat(StripMemory(g2)); + eigenspaceintersection := SumIntersectionMat(eigenspace1,eigenspace2)[2]; + + b := Basis(VectorSpace(fld,base),base); + action1 := List(base,v->Coefficients(b,v*g1)); + action2 := List(base,v->Coefficients(b,v*g2)); + return [GroupByGenerators([action1,action2]),BasisVectors(b),eigenspaceintersection]; +end; + + + +RECOG.SLn_godownStingrayWithSmallerMatrices:=function(list) +local d, first, q, p, g, i, r, pol, factors, degrees, newdim, power, rr, ss, + newgroup, colldegrees, exp, count, check, ocount, beta, NiList, Maxi, qFactors, irrfact, invbase, oneEigenspace, maxdim; + + first := function(list) + local i, j, goodElement; + for i in [1..Length(list)] do + if list[i]>1 then + if Gcd(list[i],Product(list)/list[i]) < list[i] then + return i; + else + goodElement := true; + for j in [1..Length(list)] do + if not(j = i) and Gcd(list[i],list[j]) = list[i] then + goodElement := false; + break; + fi; + od; + if goodElement then + return i; + fi; + fi; + fi; + od; + return fail; + end; + + g:=list[1]; + d:=list[2]; + q:=list[3]; + qFactors:=Factors(q); + p := qFactors[1]; + if d <= 700 then + maxdim := Maximum([Log2Int(d),3]); + else + # Test a heuristic + maxdim := Int(d/20); + fi; + + # Overall count. Replace by formula and unequality + ocount := 0; + while ocount < 100 do + + Info(InfoRecog,2,"Dimension: ",d); + #find an element with irreducible action of relative prime dimension to + #all other invariant subspaces + #count is just safety, if things go very bad + count:=0; + + repeat + count:=count+1; + if InfoLevel(InfoRecog) >= 3 then Print(".\c"); fi; + r:=PseudoRandom(g); + pol:=CharacteristicPolynomial(r); + factors:=Factors(pol); + degrees:=List(factors,Degree); + newdim:=first(degrees); + until (count>100) or (newdim <> fail and (degrees[newdim] < maxdim)); + # Be careful if Log2Int(d) = 2! In this case we search for stingray elements with k < 2. Hence use newdim < Maximum([Log2Int(d),3]) + + if count>100 then + return fail; + fi; + + # Split result from first: + irrfact := factors[newdim]; + newdim := degrees[newdim]; + + if newdim > 2 then + # Check whether the stingray candidate is a ppd-stingray element + check := RECOG.IsPpdStingrayElement(p,Length(qFactors),newdim,irrfact); + if check then + + # raise r to a power so that acting trivially outside one invariant irreducible subspace + NiList := Collected(degrees); + NiList := Filtered(NiList,x->not(x[1] = newdim)); + colldegrees := List(NiList,x->x[1]); + NiList := List(NiList,x->x[2]); + Maxi := Maximum(NiList); + beta := LogInt(Maxi,p); + if not(p^beta = Maxi) then + beta := beta + 1; + fi; + + # power further to cancel q-part of element order + power := Lcm(List(colldegrees, x->q^x-1))*p^beta; + rr:=r^power; + + invbase := NullspaceMat(TransposedMat(RECOG.FixspaceMat(TransposedMat(StripMemory(rr))))); + oneEigenspace := RECOG.FixspaceMat(StripMemory(rr)); + return [rr,newdim,invbase,oneEigenspace]; + + fi; + else + NiList := Collected(degrees); + NiList := Filtered(NiList,x->not(x[1] = newdim)); + colldegrees := List(NiList,x->x[1]); + NiList := List(NiList,x->x[2]); + Maxi := Maximum(NiList); + beta := LogInt(Maxi,p); + if not(p^beta = Maxi) then + beta := beta + 1; + fi; + + # power further to cancel q-part of element order + power := Lcm(List(colldegrees, x->q^x-1))*p^beta; + rr:=r^power; + + invbase := NullspaceMat(TransposedMat(RECOG.FixspaceMat(TransposedMat(StripMemory(rr))))); + if Size(invbase) = 2 then + oneEigenspace := RECOG.FixspaceMat(StripMemory(rr)); + return [rr,newdim,invbase,oneEigenspace]; + fi; + fi; + + ocount := ocount + 1; + od; + + return fail; + +end; + + + +RECOG.SLn_constructppdTwoStingraySmallerMatrices:=function(g,dim,q) +local out, list, out2, currentdim, check, slplist, slpToSmallerGroup, baselist, eigenspacelist; + + Info(InfoRecog,2,"Current Dimension: "); + Info(InfoRecog,2,dim); + Info(InfoRecog,2,"\n"); + + list:=[g,dim,q,g]; + slplist:=[]; + currentdim := dim; + baselist:=[]; + eigenspacelist := []; + repeat + out:=RECOG.SLn_godownStingrayWithSmallerMatrices(list); + if out=fail or out[1]*out[1]=One(out[1]) then + if InfoLevel(InfoRecog) >= 3 then Print("B\c"); fi; + Info(InfoRecog,2,"Restart. \n"); + Info(InfoRecog,2,"Current Dimension: "); + Info(InfoRecog,2,dim); + Info(InfoRecog,2,"\n"); + list:=[g,dim,q,g]; + slplist:=[]; + baselist:=[]; + eigenspacelist := []; + out:=fail; + else + if out[2]>2 then + repeat + out2:=RECOG.SLn_godownStingrayWithSmallerMatrices(list); + if out2=fail or out2[1]*out2[1]=One(out2[1]) then + if InfoLevel(InfoRecog) >= 3 then Print("B\c"); fi; + list:=[g,dim,q,g]; + slplist:=[]; + baselist:=[]; + eigenspacelist := []; + out2:=fail; + fi; + until out2<>fail and out2[2] > 2; + check := RECOG.CheckNewStingrayGroupSmallerMatrices(out[1],out[2],out[3],out[4],out2[1],out2[2],out2[3],out2[4],q); + if check[1] then + # At this point we can use the smaller matrices and use them during the next loop body execution + slpToSmallerGroup := SLPOfElms([out[1],out2[1]]); + Add(slplist,slpToSmallerGroup); + Add(baselist,check[3]); + Add(eigenspacelist,check[4]); + list:=[GroupWithMemory(check[2]),out[2]+out2[2],q,g]; + currentdim := list[2]; + + # We still have to compute the vector space on which the matrices act in the input group + + Info(InfoRecog,2,"Debugg Info:\n"); + Info(InfoRecog,2,"Dimension FirstElement: "); + Info(InfoRecog,2,out[2]); + Info(InfoRecog,2,"\n"); + Info(InfoRecog,2,"Dimension SecondElement: "); + Info(InfoRecog,2,out2[2]); + Info(InfoRecog,2,"\n"); + Info(InfoRecog,2,"End Debugg Info. \n"); + + Info(InfoRecog,2,"New Dimension: "); + Info(InfoRecog,2,out[2]+out2[2]); + Info(InfoRecog,2,"\n"); + else + if InfoLevel(InfoRecog) >= 3 then Print("B\c"); fi; + Info(InfoRecog,2,"Restart. \n"); + Info(InfoRecog,2,"Current Dimension: "); + Info(InfoRecog,2,dim); + Info(InfoRecog,2,"\n"); + list:=[g,dim,q,g]; + slplist:=[]; + baselist:=[]; + eigenspacelist := []; + out:=fail; + fi; + fi; + fi; + until out<>fail and out[2]=2; + + return [out[1],list[1],list[2],slplist,baselist,eigenspacelist]; + +end; + + + +RECOG.SLn_constructsl2WithSmallerMatrices:=function(g,d,q) +local r,h,slp,sl2,baselist,gens,b,sl2gens,sl2genss,f,eigenspacelist,subspaces,eigenspace1,eigenspace2,eigenspace3,eigenspaceintersection; + + r := RECOG.SLn_constructppdTwoStingraySmallerMatrices(g,d,q); + slp:=r[4]; + baselist:=r[5]; + eigenspacelist:=r[6]; + Info(InfoRecog,2,"Finished main GoingDown, i.e. we found a stringray element which operates irreducible on a 2 dimensional subspace. \n"); + Info(InfoRecog,2,"Next goal: Find an random element s.t. the two elements generate SL(4,q). \n"); + h := RECOG.SLn_constructsl4(r[2],r[3],q,r[1]); + Add(slp,SLPOfElms(GeneratorsOfGroup(h))); + #h := RECOG.LinearActionRepresentationSmallerMatrices(h); + #Add(baselist,h[3]); + h := GeneratorsOfGroup(h); h := RECOG.ConstructSL4SmallerMatrices(h[1],h[2],q); + Add(baselist,h[2]); + Add(eigenspacelist,h[3]); + h[1] := GroupWithMemory(h[1]); + #Error("here"); + # Remark D.R.: at this point we know that h is isomorphic to SL(4,q) + Info(InfoRecog,2,"Succesful. "); + Info(InfoRecog,2,"Current Dimension: 4\n"); + Info(InfoRecog,2,"Next goal: Generate SL(2,q). \n"); + if not (q in [2,3,4,5,9]) then + sl2 := RECOG.SLn_godownfromd(h[1],q,4,h[2]); + b := Basis(sl2[2]); + f := GF(q); + sl2gens := StripMemory(GeneratorsOfGroup(sl2[1])); + Add(eigenspacelist,RECOG.FixspaceMat(sl2gens[1])); + #eigenspace1 := RECOG.FixspaceMat(sl2gens[1]); + #eigenspace2 := RECOG.FixspaceMat(sl2gens[2]); + #eigenspace3 := RECOG.FixspaceMat(sl2gens[3]); + #eigenspaceintersection := SumIntersectionMat(eigenspace1,eigenspace2)[2]; + #eigenspaceintersection := SumIntersectionMat(eigenspaceintersection,eigenspace3)[2]; + #Add(eigenspacelist,eigenspaceintersection); + # Test by DR: + #sl2genss := List(sl2gens,x-> List(b,v->Coefficients(b,v*x))); + sl2genss := List(sl2gens,x->RECOG.LinearAction(b,f,x)); + Add(slp,SLPOfElms(GeneratorsOfGroup(sl2[1]))); + Add(baselist,BasisVectors(Basis(sl2[2]))); + # Add(eigenspacelist,RECOG.FixspaceMat(sl2gens[1])); + # Just for tests: Add(eigenspacelist,RECOG.FixspaceMat(TransposedMat(sl2gens[1]))); + Add(sl2,RECOG.ConnectSLPs(slp)); + Add(sl2,sl2genss); + subspaces := RECOG.Computesl2Subspace(baselist,eigenspacelist); + sl2[2] := subspaces[1]; + Add(sl2,subspaces[2]); + # Error("why"); + return sl2; + else + sl2 := RECOG.SLn_exceptionalgodown(h[1],q,h[2]); + b := Basis(sl2[2]); + f := GF(q); + sl2gens := StripMemory(GeneratorsOfGroup(sl2[1])); + sl2genss := List(sl2gens,x->RECOG.LinearAction(b,f,x)); + Add(slp,SLPOfElms(GeneratorsOfGroup(sl2[1]))); + Add(baselist,BasisVectors(Basis(sl2[2]))); + #Add(eigenspacelist,RECOG.FixspaceMat(sl2gens[1])); + Add(sl2,RECOG.ConnectSLPs(slp)); + Add(sl2,sl2genss); + subspaces := RECOG.Computesl2Subspace(baselist,eigenspacelist); + sl2[2] := subspaces[1]; + Add(sl2,subspaces[2]); + return sl2; + # return ["sorry only SL(4,q)",h]; + fi; +end; + + + +RECOG.ConnectSLPs:=function(slps) +local slp, currentslp, i; + + if Size(slps) = 0 then + Error("This should not happen."); + elif Size(slps) = 1 then + return slps[1]; + else + slp := slps[1]; + for i in [2..Size(slps)] do + slp := CompositionOfStraightLinePrograms(slps[i],slp); + od; + return slp; + fi; + +end; + + + +RECOG.Computesl2Subspace:=function(generators,eigenspaceGenerators) +local result, i, gens, j, combination, vec, comb, zerovec, sl2eigenspacebase, newsl2eigenspacevectors, ele; + + if Size(generators) = 1 then + # We startet with a SL(4,q) + # Just return the 2-dimensional subspace + + # TODO return eigenspacebase!!! See else case + + sl2eigenspacebase := eigenspaceGenerators[1]; + zerovec := Zero(result[1,1]) * result[1]; + for ele in eigenspaceGenerators[2] do + vec := zerovec; + for j in [1..Size(ele)] do + vec := vec + ele[j] * result[j]; + od; + Add(sl2eigenspacebase,vec); + od; + + return [generators[1],sl2eigenspacebase]; + else + # We start with the 1xd vectors + result := generators[1]; + sl2eigenspacebase := eigenspaceGenerators[1]; + zerovec := Zero(result[1,1]) * result[1]; + for ele in eigenspaceGenerators[2] do + vec := zerovec; + for j in [1..Size(ele)] do + vec := vec + ele[j] * result[j]; + od; + Add(sl2eigenspacebase,vec); + od; + + for i in [2..Size(generators)] do + combination := generators[i]; + gens := []; + for comb in combination do + vec := zerovec; + for j in [1..Size(comb)] do + vec := vec + comb[j] * result[j]; + od; + Add(gens,vec); + od; + if i+1 <= Size(eigenspaceGenerators) then + for ele in eigenspaceGenerators[i+1] do + vec := zerovec; + for j in [1..Size(ele)] do + vec := vec + ele[j] * result[j]; + od; + Add(sl2eigenspacebase,vec); + od; + fi; + result := gens; + od; + + return [result,sl2eigenspacebase]; + fi; + +end; + + +############################################################################## +# LGO approach for GoingDown to Dimension 2 +############################################################################## + + + +RECOG.SL_godownToDimension2WithInvolutions := function(h,q) +local counter, ele, ele2, x, x2, ord, invo, found, cent, product, productEle, fact1, fact2, eigenspace, result, +Minuseigenspace, newbasis, dimeigen, dimMinuseigen, gens, SL2, reco, SL2sub, pseudoorderlist, cord1, cord2, r1, r2; + + # First we construct an involution i in h + + found := false; + for counter in [1..100] do + ele := PseudoRandom(h); + x := RECOG.EstimateOrder(ele); + x2 := x[2]; + ord := x[3]; + if x2 <> One(h) then + invo := x2^(ord/2); + else + invo := One(h); + fi; + + if invo <> One(h) and invo^2 = One(h) then + eigenspace := Eigenspaces(GF(q),invo); + if Size(eigenspace) <> 1 then + Minuseigenspace := eigenspace[2]; + eigenspace := eigenspace[1]; + dimeigen := Dimension(eigenspace); + dimMinuseigen := Dimension(Minuseigenspace); + if dimeigen = 2 then + found := true; + break; + fi; + fi; + fi; + od; + + if not(found) then + Error("could not find an involution"); + fi; + + newbasis := MutableCopyMat(BasisVectors(Basis(eigenspace))); + Append(newbasis,BasisVectors(Basis(Minuseigenspace))); + + # Second we compute the two factors by computing the centralizer of the involution i + + cent := RECOG.CentraliserOfInvolution(h,invo,[],false,100,RECOG.CompletionCheck,PseudoRandom); + product := GroupByGenerators(cent[1]); + + # Third we continue as in "Constructive recognition of classical groups in odd characteristic" part 11 to find generator + + r1 := [1,2]; + r2 := [3,4]; + for counter in [1..100] do + result := RECOG.ConstructSmallSub(r1, r2, product, newbasis, g -> RECOG.IsThisSL2Natural(GeneratorsOfGroup(g),GF(q))); + if result <> fail then + break; + fi; + od; + + return fail; +end; diff --git a/gap/projective/constructive_recognition/SL/GoingUp.gi b/gap/projective/constructive_recognition/SL/GoingUp.gi new file mode 100644 index 00000000..ac89bf4e --- /dev/null +++ b/gap/projective/constructive_recognition/SL/GoingUp.gi @@ -0,0 +1,547 @@ +############################################################################# +## +## This file is part of recog, a package for the GAP computer algebra system +## which provides a collection of methods for the constructive recognition +## of groups. +## +## This files's authors include Daniel Rademacher. +## +## Copyright of recog belongs to its developers whose names are too numerous +## to list here. Please refer to the COPYRIGHT file for details. +## +## SPDX-License-Identifier: GPL-3.0-or-later +## +############################################################################# + + + +############################################################################# +############################################################################# +######## GoingUp method for special linear groups ########################### +############################################################################# +############################################################################# + + + +RECOG.SLn_UpStep := function(w) +# w has components: +# d : size of big SL +# n : size of small SL +# slnstdf : fakegens for SL_n standard generators +# bas : current base change, first n vectors are where SL_n acts +# rest of vecs are invariant under SL_n +# basi : current inverse of bas +# sld : original group with memory generators, PseudoRandom +# delivers random elements +# sldf : fake generators to keep track of what we are doing +# f : field +# The following are filled in automatically if not already there: +# p : characteristic +# ext : q=p^ext +# One : One(slnstdf[1]) +# can : CanonicalBasis(f) +# canb : BasisVectors(can) +# transh : fakegens for the "horizontal" transvections n,i for 1<=i<=n-1 +# entries can be unbound in which case they are made from slnstdf +# transv : fakegens for the "vertical" transvections i,n for 1<=i<=n-1 +# entries can be unbound in which case they are made from slnstdf +# +# We keep the following invariants (going from n -> n':=2n-1) +# bas, basi is a base change to the target base +# slnstdf are SLPs to reach standard generators of SL_n from the +# generators of sld +local DoColOp_n,DoRowOp_n,FixSLn,Fixc,MB,Vn,Vnc,aimdim,c,c1,c1f,cf,cfi, + ci,cii,coeffs,flag,i,id,int1,int3,j,k,lambda,list,mat,newbas,newbasf, + newbasfi,newbasi,newdim,newpart,perm,pivots,pivots2,pos,pow,s,sf, + slp,std,sum1,tf,trans,transd,transr,v,vals,zerovec,counter; + + Info(InfoRecog,3,"Going up: ",w.n," (",w.d,")..."); + + # Before we begin, we upgrade the data structure with a few internal + # things: + + if not(IsBound(w.can)) then w.can := CanonicalBasis(w.f); fi; + if not(IsBound(w.canb)) then w.canb := BasisVectors(w.can); fi; + if not(IsBound(w.One)) then w.One := One(w.slnstdf[1]); fi; + if not(IsBound(w.transh)) then w.transh := []; fi; + if not(IsBound(w.transv)) then w.transv := []; fi; + # Update our cache of *,n and n,* transvections because we need them + # all over the place: + std := RECOG.InitSLstd(w.f,w.n, + w.slnstdf{[1..w.ext]}, + w.slnstdf{[w.ext+1..2*w.ext]}, + w.slnstdf[2*w.ext+1], + w.slnstdf[2*w.ext+2]); + for i in [1..w.n-1] do + for k in [1..w.ext] do + pos := (i-1)*w.ext + k; + if not(IsBound(w.transh[pos])) then + RECOG.ResetSLstd(std); + RECOG.DoColOp_SL(false,w.n,i,w.canb[k],std); + w.transh[pos] := std.right; + fi; + if not(IsBound(w.transv[pos])) then + RECOG.ResetSLstd(std); + RECOG.DoRowOp_SL(false,i,w.n,w.canb[k],std); + w.transv[pos] := std.left; + fi; + od; + od; + + Unbind(std); + + # Now we can define two helper functions: + DoColOp_n := function(el,i,j,lambda,w) + # This adds lambda times the i-th column to the j-th column. + # Note that either i or j must be equal to n! + local coeffs,k; + coeffs := IntVecFFE(Coefficients(w.can,lambda)); + if i = w.n then + for k in [1..w.ext] do + if not(IsZero(coeffs[k])) then + if IsOne(coeffs[k]) then + el := el * w.transh[(j-1)*w.ext+k]; + elif not(IsZero(coeffs[k])) then + el := el * w.transh[(j-1)*w.ext+k]^coeffs[k]; + fi; + fi; + od; + elif j = w.n then + for k in [1..w.ext] do + if not(IsZero(coeffs[k])) then + if IsOne(coeffs[k]) then + el := el * w.transv[(i-1)*w.ext+k]; + else + el := el * w.transv[(i-1)*w.ext+k]^coeffs[k]; + fi; + fi; + od; + else + Error("either i or j must be equal to n"); + fi; + return el; + end; + DoRowOp_n := function(el,i,j,lambda,w) + # This adds lambda times the j-th row to the i-th row. + # Note that either i or j must be equal to n! + local coeffs,k; + coeffs := IntVecFFE(Coefficients(w.can,lambda)); + if j = w.n then + for k in [1..w.ext] do + if not(IsZero(coeffs[k])) then + if IsOne(coeffs[k]) then + el := w.transv[(i-1)*w.ext+k] * el; + else + el := w.transv[(i-1)*w.ext+k]^coeffs[k] * el; + fi; + fi; + od; + elif i = w.n then + for k in [1..w.ext] do + if not(IsZero(coeffs[k])) then + if IsOne(coeffs[k]) then + el := w.transh[(j-1)*w.ext+k] * el; + else + el := w.transh[(j-1)*w.ext+k]^coeffs[k] * el; + fi; + fi; + od; + else + Error("either i or j must be equal to n"); + fi; + return el; + end; + + # Here everything starts, some more preparations: + + # We compute exclusively in our basis, so we occasionally need an + # identity matrix: + id := IdentityMat(w.d,w.f); + FixSLn := VectorSpace(w.f,id{[w.n+1..w.d]}); + Vn := VectorSpace(w.f,id{[1..w.n]}); + + Info(InfoRecog,2,"Current dimension: " ); + Info(InfoRecog,2,w.n); + Info(InfoRecog,2,"\n"); + Info(InfoRecog,2,"New dimension: "); + Info(InfoRecog,2,Minimum(2*w.n-1,w.GoalDim)); + Info(InfoRecog,2,"\n"); + + Info(InfoRecog,2,"Preparation done."); + + ## + ## Step 1 + ## + + # First pick an element in SL_n with fixed space of dimension d-n+1: + # We already have an SLP for an n-1-cycle: it is one of the std gens. + # For n=2 we use a transvection for this purpose. + if w.n > 2 then + if IsOddInt(w.n) then + if w.p > 2 then + s := id{Concatenation([1,w.n],[2..w.n-1],[w.n+1..w.d])}; + ConvertToMatrixRepNC(s,w.f); + if IsOddInt(w.n) then s[2] := -s[2]; fi; + sf := w.slnstdf[2*w.ext+2]; + else # in even characteristic we take the n-cycle: + s := id{Concatenation([w.n],[1..w.n-1],[w.n+1..w.d])}; + ConvertToMatrixRepNC(s,w.f); + sf := w.slnstdf[2*w.ext+1]; + fi; + else + Error("this program only works for odd n or n=2"); + fi; + else + # In this case the n-1-cycle is the identity, so we take a transvection: + s := MutableCopyMat(id); + s[1,2] := One(w.f); + sf := w.slnstdf[1]; + fi; + + Info(InfoRecog,2,"Step 1 done."); + + # Find a good random element: + w.count := 0; + aimdim := Minimum(2*w.n-1,w.GoalDim); + newdim := aimdim - w.n; + counter := 0; + while true do # will be left by break + + ## + ## Step 2 + ## + v := fail; + repeat + counter := counter + 1; + if InfoLevel(InfoRecog) >= 3 then Print(".\c"); fi; + w.count := w.count + 1; + c1 := PseudoRandom(w.sld); + + # Do the base change into our basis: + #c1 := w.bas * c1 * w.basi; + c := s^(w.bas * c1 * w.basi); + + # Check how these elements look like. Where is the SLP and what elements do we really use + + # Now check that Vn + Vn*s^c1 has dimension 2n-1: + sum1 := SumIntersectionMat(id{[1..w.n]}, c{[1..w.n]}); + if Size(sum1[1]) = aimdim then + # intersect Fix(c) = Nullspace(c-id) with V_n in order to + # find a suitable vector which we can later to our basis + int1 := SumIntersectionMat(RECOG.FixspaceMat(c),id{[1..w.n]})[2]; + v := First(int1, v -> not IsZero(v[w.n])); + if v = fail then + Info(InfoRecog,2,"Ooops: Component n was zero!"); + fi; + fi; + until v <> fail; + + v := v / v[w.n]; # normalize to 1 in position n + Assert(1,v*c=v); + + # now that we have our c and c1, compute some associated + # values for later use + ci := c^-1; + slp := SLPOfElm(c1); + c1f := ResultOfStraightLineProgram(slp,w.sldf); + cf := sf^c1f; + cfi := cf^-1; + + Info(InfoRecog,2,"Step 2 done."); + + ## + ## Steps 3 and 4 + ## + + # Now we found our aimdim-dimensional space W. Since SL_n + # has a d-n-dimensional fixed space W_{d-n} and W contains a complement + # of that fixed space, the intersection of W and W_{d-n} has dimension + # newdim. + + # Change basis: + newpart := ExtractSubMatrix(c,[1..(w.n-1)],[1..(w.d)]); + # Clean out the first n entries to go to the fixed space of SL_n: + zerovec := Zero(newpart[1]); + for i in [1..(w.n-1)] do + CopySubVector(zerovec,newpart[i],[1..w.n],[1..w.n]); + od; + MB := MutableBasis(w.f,[],zerovec); + i := 1; + pivots := EmptyPlist(newdim); + while i <= Length(newpart) and NrBasisVectors(MB) < newdim do + if not(IsContainedInSpan(MB,newpart[i])) then + Add(pivots,i); + CloseMutableBasis(MB,newpart[i]); + fi; + i := i + 1; + od; + newpart := newpart{pivots}; + newbas := Concatenation(id{[1..w.n-1]},[v],newpart); + if 2*w.n-1 < w.d then + + # intersect Fix(c) with F_{d-n} + int3 := SumIntersectionMat(RECOG.FixspaceMat(c),id{[w.n+1..w.d]})[2]; + if Size(int3) <> w.d - aimdim then + Info(InfoRecog,2,"Ooops, FixSLn \cap Fixc wrong dimension"); + continue; + fi; + Append(newbas,int3); + fi; + ConvertToMatrixRep(newbas,Size(w.f)); + newbasi := newbas^-1; + if newbasi = fail then + Info(InfoRecog,2,"Ooops, Fixc intersected too much, we try again"); + continue; + fi; + + ci := newbas * ci * newbasi; + cii := ExtractSubMatrix(ci,[w.n+1..aimdim],[1..w.n-1]); + ConvertToMatrixRep(cii,Size(w.f)); + cii := TransposedMat(cii); + # The rows of cii are now what used to be the columns, + # their length is newdim, we need to span the full newdim-dimensional + # row space and need to remember how: + zerovec := Zero(cii[1]); + MB := MutableBasis(w.f,[],zerovec); + i := 1; + pivots2 := EmptyPlist(newdim); + while i <= Length(cii) and NrBasisVectors(MB) < newdim do + if not(IsContainedInSpan(MB,cii[i])) then + Add(pivots2,i); + CloseMutableBasis(MB,cii[i]); + fi; + i := i + 1; + od; + if Length(pivots2) = newdim then + break; + fi; + Info(InfoRecog,2,"Ooops, no nice bottom..."); + # Otherwise simply try again + od; + + cii := cii{pivots2}^-1; + ConvertToMatrixRep(cii,w.f); + c := newbas * c * newbasi; + w.bas := newbas * w.bas; + w.basi := w.basi * newbasi; + + + Info(InfoRecog,2," found c1 and c."); + # Now SL_n has to be repaired according to the base change newbas: + + # Now write this matrix newbas as an SLP in the standard generators + # of our SL_n. Then we know which generators to take for our new + # standard generators, namely newbas^-1 * std * newbas. + + newbasf := w.One; + for i in [1..w.n-1] do + if not(IsZero(v[i])) then + newbasf := DoColOp_n(newbasf,w.n,i,v[i],w); + fi; + od; + newbasfi := newbasf^-1; + w.slnstdf := List(w.slnstdf,x->newbasfi * x * newbasf); + # Now update caches: + w.transh := List(w.transh,x->newbasfi * x * newbasf); + w.transv := List(w.transv,x->newbasfi * x * newbasf); + + Info(InfoRecog,2,"Step 3 and 4 done"); + + ## + ## Step 5 + ## + + # Now consider the transvections t_i: + # t_i : w.bas[j] -> w.bas[j] for j <> i and + # t_i : w.bas[i] -> w.bas[i] + ww + # We want to modify (t_i)^c such that it fixes w.bas{[1..w.n]}: + trans := []; + for i in pivots2 do + # This does t_i + for lambda in w.canb do + # This does t_i : v_j -> v_j + lambda * v_n + tf := w.One; + tf := DoRowOp_n(tf,i,w.n,lambda,w); + # Now conjugate with c: + tf := cfi*tf*cf; + # Now cleanup in column n above row n, the entries there + # are lambda times the stuff in column i of ci: + for j in [1..w.n-1] do + tf := DoRowOp_n(tf,j,w.n,-ci[j,i]*lambda,w); + od; + Add(trans,tf); + od; + od; + + # Now put together the clean ones by our knowledge of c^-1: + transd := []; + for i in [1..Length(pivots2)] do + for lambda in w.canb do + tf := w.One; + vals := BlownUpVector(w.can,cii[i]*lambda); + for j in [1..w.ext * newdim] do + pow := IntFFE(vals[j]); + if not(IsZero(pow)) then + if IsOne(pow) then + tf := tf * trans[j]; + else + tf := tf * trans[j]^pow; + fi; + fi; + od; + Add(transd,tf); + od; + od; + Unbind(trans); + + Info(InfoRecog,2,"Step 5 done"); + + ## + ## Step 6 + ## + + # Now to the "horizontal" transvections, first create them as SLPs: + transr := []; + for i in pivots do + # This does u_i : v_i -> v_i + v_n + tf := w.One; + tf := DoColOp_n(tf,w.n,i,One(w.f),w); + # Now conjugate with c: + tf := cfi*tf*cf; + # Now cleanup in rows above row n: + for j in [1..w.n-1] do + tf := DoRowOp_n(tf,j,w.n,-ci[j,w.n],w); + od; + # Now cleanup in rows below row n: + for j in [1..newdim] do + coeffs := IntVecFFE(Coefficients(w.can,-ci[w.n+j,w.n])); + for k in [1..w.ext] do + if not(IsZero(coeffs[k])) then + if IsOne(coeffs[k]) then + tf := transd[(j-1)*w.ext + k] * tf; + else + tf := transd[(j-1)*w.ext + k]^coeffs[k] * tf; + fi; + fi; + od; + od; + + # Now cleanup column n above row n: + for j in [1..w.n-1] do + tf := DoColOp_n(tf,j,w.n,ci[j,w.n],w); + od; + + # Now cleanup row n left of column n: + for j in [1..w.n-1] do + tf := DoRowOp_n(tf,w.n,j,-c[i,j],w); + od; + + # Now cleanup column n below row n: + for j in [1..newdim] do + coeffs := IntVecFFE(Coefficients(w.can,ci[w.n+j,w.n])); + for k in [1..w.ext] do + if not(IsZero(coeffs[k])) then + if IsOne(coeffs[k]) then + tf := tf * transd[(j-1)*w.ext + k]; + else + tf := tf * transd[(j-1)*w.ext + k]^coeffs[k]; + fi; + fi; + od; + od; + Add(transr,tf); + od; + + Info(InfoRecog,2,"Step 6 done"); + + ## + ## Step 7 + ## + + # From here on we distinguish three cases: + # * w.n = 2 + # * we finish off the constructive recognition + # * we have to do another step as the next thing + if w.n = 2 then + w.slnstdf[2*w.ext+2] := transd[1]*transr[1]^-1*transd[1]; + w.slnstdf[2*w.ext+1] := w.transh[1]*w.transv[1]^-1*w.transh[1] + *w.slnstdf[2*w.ext+2]; + Unbind(w.transh); + Unbind(w.transv); + w.n := 3; + Info(InfoRecog,2,"Step 7 done"); + return w; + fi; + # We can finish off: + if aimdim = w.GoalDim then + # In this case we just finish off and do not bother with + # the transvections, we will only need the standard gens: + # Now put together the (newdim+1)-cycle: + # n+newdim -> n+newdim-1 -> ... -> n+1 -> n -> n+newdim + flag := false; + s := w.One; + for i in [1..newdim] do + if flag then + # Make [[0,-1],[1,0]] in coordinates w.n and w.n+i: + tf:=transd[(i-1)*w.ext+1]*transr[i]^-1*transd[(i-1)*w.ext+1]; + else + # Make [[0,1],[-1,0]] in coordinates w.n and w.n+i: + tf:=transd[(i-1)*w.ext+1]^-1*transr[i]*transd[(i-1)*w.ext+1]^-1; + fi; + s := s * tf; + flag := not(flag); + od; + + # Finally put together the new 2n-1-cycle and 2n-2-cycle: + s := s^-1; + w.slnstdf[2*w.ext+1] := w.slnstdf[2*w.ext+1] * s; + w.slnstdf[2*w.ext+2] := w.slnstdf[2*w.ext+2] * s; + Unbind(w.transv); + Unbind(w.transh); + w.n := aimdim; + Info(InfoRecog,2,"Step 7 done"); + return w; + fi; + + # Otherwise we do want to go on as the next thing, so we want to + # keep our transvections. This is easily done if we change the + # basis one more time. Note that we know that n is odd here! + + # Put together the n-cycle: + # 2n-1 -> 2n-2 -> ... -> n+1 -> n -> 2n-1 + flag := false; + s := w.One; + for i in [w.n-1,w.n-2..1] do + if flag then + # Make [[0,-1],[1,0]] in coordinates w.n and w.n+i: + tf := transd[(i-1)*w.ext+1]*transr[i]^-1*transd[(i-1)*w.ext+1]; + else + # Make [[0,1],[-1,0]] in coordinates w.n and w.n+i: + tf := transd[(i-1)*w.ext+1]^-1*transr[i]*transd[(i-1)*w.ext+1]^-1; + fi; + s := s * tf; + flag := not(flag); + od; + + # Finally put together the new 2n-1-cycle and 2n-2-cycle: + w.slnstdf[2*w.ext+1] := s * w.slnstdf[2*w.ext+1]; + w.slnstdf[2*w.ext+2] := s * w.slnstdf[2*w.ext+2]; + + list := Concatenation([1..w.n-1],[w.n+1..2*w.n-1],[w.n],[2*w.n..w.d]); + perm := PermList(list); + mat := PermutationMat(perm^-1,w.d,w.f); + ConvertToMatrixRep(mat,w.f); + w.bas := w.bas{list}; + ConvertToMatrixRep(w.bas,w.f); + w.basi := w.basi*mat; + + # Now add the new transvections: + for i in [1..w.n-1] do + w.transh[w.ext*(w.n-1)+w.ext*(i-1)+1] := transr[i]; + od; + Append(w.transv,transd); + w.n := 2*w.n-1; + + Info(InfoRecog,2,"Step 7 done"); + return w; +end; diff --git a/gap/projective/constructive_recognition/SL/main.gi b/gap/projective/constructive_recognition/SL/main.gi new file mode 100644 index 00000000..d5796ebf --- /dev/null +++ b/gap/projective/constructive_recognition/SL/main.gi @@ -0,0 +1,245 @@ +############################################################################# +## +## This file is part of recog, a package for the GAP computer algebra system +## which provides a collection of methods for the constructive recognition +## of groups. +## +## This files's authors include Daniel Rademacher. +## +## Copyright of recog belongs to its developers whose names are too numerous +## to list here. Please refer to the COPYRIGHT file for details. +## +## SPDX-License-Identifier: GPL-3.0-or-later +## +############################################################################# + + + +############################################################################# +############################################################################# +######## Main function for special linear groups ############################ +############################################################################# +############################################################################# + + + +RECOG.FindStdGens_SL := function(sld) + + return RECOG.FindStdGens_SL2(sld,DimensionOfMatrixGroup(sld)); + +end; + + + +RECOG.FindStdGens_SL2 := function(sld,IsoDim) + + # Group generated by input must be isomorphic SL(IsoDim,q) + + # gens of sld must be gens for SL(d,q) in its natural rep with memory + # This function calls RECOG.SLn_constructsl2 and then extends + # the basis to a basis of the full row space and calls + # RECOG.SLn_UpStep often enough. Finally it returns an slp such + # that the SL(d,q) standard generators with respect to this basis are + # expressed by the slp in terms of the original generators of sld. + local V,b,bas,basi,basit,d,data,ext,fakegens,id,nu,nu2,p,q,resl2,sl2,sl2gens, + sl2gensf,sl2genss,sl2stdf,slp,slpsl2std,slptosl2,st,std,stdgens,i,ex,f; + + # Some setup: + f := FieldOfMatrixGroup(sld); + p := Characteristic(f); + q := Size(f); + ext := DegreeOverPrimeField(f); + d := DimensionOfMatrixGroup(sld); + if not(IsObjWithMemory(GeneratorsOfGroup(sld)[1])) then + sld := GroupWithMemory(sld); + fi; + + # First find an SL2 with the space it acts on; + Info(InfoRecog,2,"Finding an SL2..."); + Info(InfoRecog,2,"-----"); + Info(InfoRecog,2,"Start of the GoingDown Algorithm."); + data := RECOG.SLn_constructsl2(sld,d,q); + Info(InfoRecog,2,"The GoingDown Algorithm was successful."); + Info(InfoRecog,2,"-----"); + + bas := ShallowCopy(BasisVectors(Basis(data[2]))); + sl2 := data[1]; + slptosl2 := SLPOfElms(GeneratorsOfGroup(sl2)); + sl2gens := StripMemory(GeneratorsOfGroup(sl2)); + V := data[2]; + b := Basis(V,bas); + sl2genss := List(sl2gens,x->RECOG.LinearAction(b,f,x)); + + Info(InfoRecog,2,"-----"); + Info(InfoRecog,2,"Solving the base case"); + if q in [2,3,4,5,9] then + Info(InfoRecog,2,"In fact found an SL4..."); + stdgens := RECOG.MakeSL_StdGens(p,ext,4,4).all; + slpsl2std := RECOG.FindStdGensUsingBSGS(Group(sl2genss),stdgens, + false,false); + nu := List(sl2gens,x->RECOG.FixspaceMat(x)); + ex := SumIntersectionMat(nu[1],nu[2])[2]; + for i in [3..Length(nu)] do + ex := SumIntersectionMat(nu[3],ex); + od; + Append(bas,ex); + ConvertToMatrixRep(bas,q); + basi := bas^-1; + else + # Now compute the natural SL2 action and run constructive recognition: + Info(InfoRecog,2, + "Recognising this SL2 constructively in 2 dimensions..."); + sl2genss := GeneratorsWithMemory(sl2genss); + if IsEvenInt(q) then + resl2 := RECOG.RecogniseSL2NaturalEvenChar(Group(sl2genss),f,false); + else + resl2 := RECOG.RecogniseSL2NaturalOddCharUsingBSGS(Group(sl2genss),f); + fi; + slpsl2std := SLPOfElms(resl2.all); + bas := resl2.bas * bas; + # We need the actual transvections: + slp := SLPOfElms([resl2.s[1],resl2.t[1]]); + st := ResultOfStraightLineProgram(slp, + StripMemory(GeneratorsOfGroup(sl2))); + + # Extend basis by something invariant under SL2: + id := IdentityMat(d,f); + nu := NullspaceMat(StripMemory(st[1]-id)); + nu2 := NullspaceMat(StripMemory(st[2]-id)); + Append(bas,SumIntersectionMat(nu,nu2)[2]); + ConvertToMatrixRep(bas,q); + basi := bas^-1; + fi; + Info(InfoRecog,2,"Finished the base case."); + Info(InfoRecog,2,"-----"); + + # Now set up fake generators for keeping track what we do: + fakegens := ListWithIdenticalEntries(Length(GeneratorsOfGroup(sld)),()); + fakegens := GeneratorsWithMemory(fakegens); + sl2gensf := ResultOfStraightLineProgram(slptosl2,fakegens); + sl2stdf := ResultOfStraightLineProgram(slpsl2std,sl2gensf); + std := rec( f := f, d := d, GoalDim := IsoDim, n := 2, bas := bas, basi := basi, + sld := sld, sldf := fakegens, slnstdf := sl2stdf, + p := p, ext := ext ); + Info(InfoRecog,2,"Going up to SL_d again..."); + Info(InfoRecog,2,"-----"); + Info(InfoRecog,2,"Start of the GoingUp Algorithm"); + while std.n < std.GoalDim do + RECOG.SLn_UpStep(std); + od; + Info(InfoRecog,2,"The GoingUp Algorithm was successful."); + Info(InfoRecog,2,"-----"); + return rec( slpstd := SLPOfElms(std.slnstdf), + bas := std.bas, basi := std.basi ); +end; + + + +RECOG.FindStdGensSmallerMatrices_SL := function(sld) + + return RECOG.FindStdGensSmallerMatrices_SL2(sld,DimensionOfMatrixGroup(sld)); + +end; + + + + +RECOG.FindStdGensSmallerMatrices_SL2 := function(sld,IsoDim) + +# Group generated by input must be isomorphic SL(IsoDim,q) + +# gens of sld must be gens for SL(d,q) in its natural rep with memory +# This function calls RECOG.SLn_constructsl2 and then extends +# the basis to a basis of the full row space and calls +# RECOG.SLn_UpStep often enough. Finally it returns an slp such +# that the SL(d,q) standard generators with respect to this basis are +# expressed by the slp in terms of the original generators of sld. +local V,b,bas,basi,basit,d,data,ext,fakegens,id,nu,nu2,p,q,resl2,sl2,sl2gens, + sl2gensf,sl2genss,sl2stdf,slp,slpsl2std,slptosl2,st,std,stdgens,i,ex,f; + + # Some setup: + f := FieldOfMatrixGroup(sld); + p := Characteristic(f); + q := Size(f); + ext := DegreeOverPrimeField(f); + d := DimensionOfMatrixGroup(sld); + if not(IsObjWithMemory(GeneratorsOfGroup(sld)[1])) then + sld := GroupWithMemory(sld); + #### Added by DR! Optimize this line! + ### second argument for length of list, third argument for number of shuffles + Group_InitPseudoRandom(sld,Size(GeneratorsOfGroup(sld))+2,5); + # Noch weniger initialisierungen testen + # TODO: Wenn initialisiert, nicht nochmal, if hinzufügen + fi; + + # First find an SL2 with the space it acts on; + Info(InfoRecog,2,"Finding an SL2..."); + Info(InfoRecog,2,"-----"); + Info(InfoRecog,2,"Start of the GoingDown Algorithm."); + data := RECOG.SLn_constructsl2WithSmallerMatrices(sld,d,q); + Info(InfoRecog,2,"The GoingDown Algorithm was successful."); + Info(InfoRecog,2,"-----"); + + bas := ShallowCopy(data[2]); + sl2 := data[1]; + slptosl2 := data[3]; + sl2gens := StripMemory(GeneratorsOfGroup(sl2)); + b := Basis(VectorSpace(GF(q),bas),bas); + sl2genss := data[4]; + + Info(InfoRecog,2,"-----"); + Info(InfoRecog,2,"Solving the base case"); + if q in [2,3,4,5,9] then + Info(InfoRecog,2,"In fact found an SL4..."); + stdgens := RECOG.MakeSL_StdGens(p,ext,4,4).all; + slpsl2std := RECOG.FindStdGensUsingBSGS(Group(sl2genss),stdgens, + false,false); + Append(bas,data[5]); + ConvertToMatrixRep(bas,q); + basi := bas^-1; + else + # Now compute the natural SL2 action and run constructive recognition: + Info(InfoRecog,2, + "Recognising this SL2 constructively in 2 dimensions..."); + sl2genss := GeneratorsWithMemory(sl2genss); + if IsEvenInt(q) then + resl2 := RECOG.RecogniseSL2NaturalEvenChar(Group(sl2genss),f,false); + else + resl2 := RECOG.RecogniseSL2NaturalOddCharUsingBSGS(Group(sl2genss),f); + fi; + slpsl2std := SLPOfElms(resl2.all); + if resl2.bas <> [[1,0],[0,1]]*One(f) then + Error("So i have to deal with this case..."); + # RECOG.FindStdGensSmallerMatrices_SL(SL(200,2^6),200); gives an example for this case. So we have + # to add at least one more base change + else + #bas := resl2.bas * bas; + # We need the actual transvections: + slp := SLPOfElms([resl2.s[1],resl2.t[1]]); + Append(bas,data[5]); + ConvertToMatrixRep(bas,q); + basi := bas^-1; + fi; + fi; + Info(InfoRecog,2,"Finished the base case."); + Info(InfoRecog,2,"-----"); + + # Now set up fake generators for keeping track what we do: + fakegens := ListWithIdenticalEntries(Length(GeneratorsOfGroup(sld)),()); + fakegens := GeneratorsWithMemory(fakegens); + sl2gensf := ResultOfStraightLineProgram(slptosl2,fakegens); + sl2stdf := ResultOfStraightLineProgram(slpsl2std,sl2gensf); + std := rec( f := f, d := d, GoalDim := IsoDim, n := 2, bas := bas, basi := basi, + sld := sld, sldf := fakegens, slnstdf := sl2stdf, + p := p, ext := ext ); + Info(InfoRecog,2,"Going up to SL_d again..."); + Info(InfoRecog,2,"-----"); + Info(InfoRecog,2,"Start of the GoingUp Algorithm"); + while std.n < std.GoalDim do + RECOG.SLn_UpStep(std); + od; + Info(InfoRecog,2,"The GoingUp Algorithm was successful."); + Info(InfoRecog,2,"-----"); + return rec( slpstd := SLPOfElms(std.slnstdf), + bas := std.bas, basi := std.basi ); +end; diff --git a/gap/projective/constructive_recognition/SL/sl2_BlackBox.gi b/gap/projective/constructive_recognition/SL/sl2_BlackBox.gi new file mode 100644 index 00000000..db11537e --- /dev/null +++ b/gap/projective/constructive_recognition/SL/sl2_BlackBox.gi @@ -0,0 +1,637 @@ +############################################################################# +## +## This file is part of recog, a package for the GAP computer algebra system +## which provides a collection of methods for the constructive recognition +## of groups. +## +## This files's authors include Daniel Rademacher. +## +## Copyright of recog belongs to its developers whose names are too numerous +## to list here. Please refer to the COPYRIGHT file for details. +## +## SPDX-License-Identifier: GPL-3.0-or-later +## +############################################################################# + + + +################################################################################################### +################################################################################################### +######## Solve Symmetric Powers (Ducs Code) ####################################################### +################################################################################################### +################################################################################################### + + + +# TODO: Use this for the constructive recognition of SL(2,q) + +# Code has been written by Duc Khuat during his Bachelors thesis +# This partly implements an algorithm based on the paper ”Constructive Recognition of SL(2, q)” by Leedham Green and E. A. O’Brien. +# For q being a p-power, the algorithm can only be applied to representations where the degree is smaller than p. + +## computes for an element of SL(2,q) its representation in the n-th symmetric power +## F is the field +## n-th symmetric power +## A element of SL(2,q) +## return : Matrix of dimension n+1 corresponding to the action of SL(2,q) on T_n represented in the basis ( x^n, ...., y^n) +RECOG.SymPowRepSL2 :=function(F,n, A) + local res,i,t,k,sum; + res := IdentityMat(n+1,F); + for i in [0..n] do + for t in [0..n] do + sum :=0; + for k in [0..i] do + if n-i-(t-k) > -1 and t-k > -1 then + sum:= sum + Binomial(n-i,t-k)*Binomial(i,k)* (A[1][1])^(t-k)*(A[2][1])^k*(A[1][2])^(n -i -(t-k))*(A[2][2])^(i-k); + fi; + od; + res[i+1][n-t+1]:= sum; + od; + od; + return res; +end;; + + + +## randomly looks for an element of order q-1. +##input: +## G the group where we look randomly for an element of order n. +## n the order of g. +## return g an element with order n, and the number of tries. +RECOG.RandomElementOfOrder:= function(G, n) + local nrTries ,g; nrTries := 0; + while nrTries < 1000 do + g := PseudoRandom(G); + if (Order(g) = n) then + return [g, nrTries]; + fi; + nrTries := nrTries +1; + od; + ####### Added by Daniel Rademacher ####### + return ["fail", "fail"]; + ########################################## + ErrorNoReturn ( " No element of order ", n, " has been found.\n"); +end ;; + + + +## z^r is the eigenvalue of g on natural module. +## find the unique (up to multiplication of -1) element in Z_q-1, to obtain the expected set of exponents, namely {-n, -n+2, ..., n-2, n }. +##input: +## g has order q-1 and q-1 eigenvectors. +## E is the set of exponent of eigenvalues in respect to a fixed primitve element of the field ## F the underlying field F of order q. +## r in the unit group of Z_{q-1} such that E*(r^{-1}) = { -s,-s+2, ..., s-2,s} +## return M matrix with eigenvectors as rows such that the i-th row is the eigenvector to -s + 2(i-1) for i =1 , ..., s+1. +RECOG.OrderEigenvectors := function(g , E ,F, r) + local + i, # every row of M. + s, # Eigenvector of g + M, # Output Matrix + EVs, # Eigenvectors of g. + z; + z:= PrimitiveElement(F); + M := []; + EVs := Eigenvectors(F,g); + for i in [ 1 .. DimensionsMat(g)[1]] do + for s in EVs do + if s*g = (z^(E[i]*r))*s then + M[i] := s; + fi; + od; + od; + return M; +end;; + + + +## find k such that List([0..n], x -> (n-2*x) mod (q-1)) * k = E. +## this k is unique, up to adding (q-1) multiples. +## input: +## E is the set of exponent of eigenvalues in respect to a fixed primitve element of the field +## n the n-th symmetric power. +## q the order of the underlying field. +## return k such that E = k * {-s,-s+2,..., s-2, s}. +## Note: k is unique up to sign for the considered cases. +RECOG.EVNatRep := function(E,n,q) + local k,l ,res, Exp; + res := 0; + Exp:=List([0..n], x -> (n-2*x) mod (q-1)); + for k in PrimeResidues(q-1) do + l := OrderMod(k,q-1)-1; # inverse of k in Z_q-1. + if Set(E*(k^l) mod (q-1)) = Set(Exp) then + res := k; + break; + fi; + od; + return res; +end;; + + + +## Symmetric Power Basis ### +## representation of SL(2,q) in GL(T_n) for n < p. +## find an element g having N distinct eigenvalues; +## if G is symmetric power of SL(2, q) n = G. And h conjugated to g. + con , # random conjugation element + i, # iterating through columns + mu_i , # coefficients of basis element + abZero , # if ab = 0 we take the last row, if ab not zero, we takte the first row of h. + ab, # coefficents to + mu_bet # coefficents to + ; + F := FieldOfMatrixGroup(G); + z:= PrimitiveElement(F); + n := DimensionsMat(PseudoRandom(G))[1] -1; + q := Size(F); + p := Characteristic(F); + k :=1; + if q < 6 or ((p mod 2 = 1) and q = p and p > 6 and n > (p-7)/2 and (not (p = 13 and n =4) )) then + ErrorNoReturn (" Exceptional Case , use another method "); + fi; + if (Size(PseudoRandom(G)) mod 2 = 0) then + g := RECOG.RandomElementOfOrder(G,q-1)[1]; + else + g:= RECOG.RandomElementOfOrder(G,(q-1)/2)[1]; + fi; + ####### Added by me ####### + if g = "fail" then + return ["fail"]; + fi; + ########################### + Ek:= List(Eigenvalues(F,g), x -> LogFFE(x,z)); + k := RECOG.EVNatRep(Ek,n,q); + M:= RECOG.OrderEigenvectors(g,List([0..n], x -> (n-2*x) mod (q-1) ), F, k); + ####### Added by me ####### + if M = [] then + return ["fail"]; + fi; + ########################### + #correct coefficients of ( x^(n-2)y^2,..,y^n) + con := IdentityMat(n+1,F); + while g*g^con = g^con*g do + con:= PseudoRandom(G); + od; + h := M*g^con*M^(-1); # = G. + if not (h[1] in Subspace(VectorSpace(F,g), [IdentityMat(n+1,F)[1]]) or h[1] in Subspace(VectorSpace(F,g), [IdentityMat(n+1,F)[n+1] ])) then + abZero :=1; + else + abZero := n+1; + fi; + ab := h[abZero][1] / (n^(-1)* h[abZero][2]); mu_bet := z^0; + for i in [2..n] do + mu_i := mu_bet*Binomial(n,i-1)^(-1)*ab^(-1)* h[abZero][i] /(Binomial(n,i)^(-1)*h[abZero][i+1]); + mu_bet := mu_i; + M[i+1]:= mu_i^(-1)*M[i+1]; + od; + return [M,g]; +end ;; + + + +## For a symmetric power G and elm of G construct image in PSL(2,q). +## input: G symmetric power of SL(2,q) of degree n < p. +## elm arbitrary matrix in G. +## Trafo the basis of the form (x^n, ..., y^n) for some element of order q-1 and eigenvectors x and y on the natrual module of SL(2,q). +## return: image of elm in PSL(2,q) for one possible homomorphism of +RECOG.HomToPSL := function (G, elm, Trafo, nOdd) + local + F, # field of matrix group + n, # degree of the symmetric power + z, # primitives element of the field + M, # the basis of the form (x^n,..., y^n) + h, # elm represented in M + k, # exponend of a^2 or d^2 + ba,ca,da,a2,a,cd,bd,d2,d,bc,c2,c, #quotients ba = b/a. + V; + F:= FieldOfMatrixGroup(G); z:= PrimitiveElement(F); + M :=Trafo; + n:= Size(M)-1; + h := M * elm * M^(-1); + V:= VectorSpace(F,M); # equals F^(n+1) + if not h[1] in Subspace(V,[IdentityMat(n+1, F)[n+1]]) then #check if a=0 + ba := (h[1][2]*n^(-1))/ h[1][1]; + ca := h[2][1]/ h[1][1]; + da := h[2][2]/h[1][1] - (n-1)*ba*ca; a2 :=1/ (da - ba*ca); k:=LogFFE(a2,z); + if nOdd then + a := h[1][1]/(a2)^(QuoInt(n,2)); + elif k mod 2 = 0 then + a := z^(k/2); + else + ErrorNoReturn("Something is wrong."); + fi; + return [[a,ba*a],[ca*a,da*a]]; + elif not h[n+1] in Subspace(V,[IdentityMat(n+1,F)[1] ]) then #check if d =0 + cd := (h[n+1][n]*n^(-1))/h[n+1][n+1]; bd := h[n][n+1] / h[n+1][n+1]; + d2 :=1 / ( - bd * cd); + k := LogFFE(d2,z); + if nOdd then + d := h[n+1][n+1] / (d2)^(QuoInt(n,2)); + elif k mod 2 = 0 then d := z^(k/2); + else + ErrorNoReturn("Something is wrong."); + fi; + return [[0,bd*d],[cd*d,d]]; + else + #if a=d=0 + bc := h[1][n+1]/ h[2][n]; c2 := -(bc)^(-1); + k := LogFFE(c2,z); + if nOdd then + c := h[n+1][1] / (c2)^(QuoInt(n,2)); + elif k mod 2 = 0 then + c := z^(k/2); + else + ErrorNoReturn("Something is wrong."); + fi; + return [[0,bc*c],[c,0]]; + fi; +end;; + + + +## MAIN FUNCTION +## G n-th symmetric power of SL(2,q) in GL(T_n) for n < p. +## return [homomorphism from G to PSL(2,q), homomorphism from SL(2,q) to G] +RECOG.RecTest := function(G) + local Trafo ,d, F; + Trafo:= RECOG.SymmetricPowerBasis(G)[1]; + if Trafo = "fail" then + return fail; + fi; + F:= FieldOfMatrixGroup(G); # underlying field + d:= Size(PseudoRandom(G))-1; # d-th symmetric power + return [x-> RECOG.HomToPSL(G,x,Trafo, d mod 2 = 1), x->Trafo^(-1)*RECOG.SymPowRepSL2(F,d,x)*Trafo]; +end;; + + + +################################################################################################### +################################################################################################### +######## Decomposition of Tensor Products ######################################################### +################################################################################################### +################################################################################################### + + + +# given sequence E of elements of finite field, construct certain +# arithmetic progressions; if MaxAP is true, give up as soon as +# we find one of length #E */ + +# If unset: MaxAP := false +RECOG.FindAPs := function (E, deg, p, MaxAP) +local AP, first, x, y, index, d, list, i, term; + + if Size(E) <= 1 then + return []; + fi; + + AP := []; + for first in [1..Size(E)] do + x := E[first]; + for index in [1..Size(E)] do + if index <> first then + y := E[index]; + d := y - x; + list := [x, y]; + if ((deg mod 2) = 0) then + Append (AP, list); + fi; + for i in [3..Size(E)] do + if MaxAP = false and Size(list) >= p then + break; + fi; + term := list[i - 1] + d; + if (term in E) and not(term in list) then + list[i] := term; + if ((deg mod i) = 0) then + Append (AP, list); + fi; + else + break; + fi; + od; + if MaxAP and Size(list) = Size(E) then + return [list]; + fi; + fi; + od; + od; + + return Set(AP); + #return SetToSequence (Set (AP)); + +end; + + + +RECOG.RandomElementOfOrderModuleCentre := function(G, required, MaxTries, Central) +local nrTries ,g, ord, centre; + + nrTries := 0; + centre := Centre(G.group); + while nrTries < MaxTries do + g := PseudoRandom(G.group); + ord := Order(g); + if (ord = 2*required) then + if (g^(required) in centre) then + return [true, g, nrTries]; + fi; + Display("nop"); + #return [true, g, nrTries]; + fi; + nrTries := nrTries +1; + od; + ####### Added by Daniel Rademacher ####### + return ["fail", "fail", nrTries]; + +end; + + + +# construct space determined by g and arithmetic progression a +# of its eigenvalues, and send space to FindPoint + +RECOG.ProcessSequence := function(G, g, a) +local F, w, ev, Space, CB, Status, gens, vec; + + F := G.fld; + w := PrimitiveElement(F); + ev := List(a, x-> w^(Int(x))); + # Maybe this line means to make one list with all the spaces? If yes, modify other lines like this too.... Space := &+[Eigenspace (g, e): e in ev]; + Space := List(ev, e -> RECOG.EigenspaceMat(g,e)); + gens := []; + for vec in Space do + Append(gens,vec); + od; + Space := VectorSpace(F,gens); + CB := "unknown"; + Status := false; + Space := RECOG.GeneralFindPoint(G, Space, Dimension(Space), Status, CB, true); + Status := Space[2]; + CB := Space[3]; + Space := Space[1]; + return [Status, CB]; +end; + + + +RECOG.StoreDetails := function(G, Result) +local U,W, CB; + CB := Result[1]; + RECOG.SetTensorProductFlag(G, true); + RECOG.SetTensorBasis(G, CB); + U := RECOG.ConstructTensorFactors(G, Result); + W := U[2]; + U := U[1]; + RECOG.SetTensorFactors (G, [U, W]); +end; + + + +# G is a tensor product of symmetric powers of SL(2, q) +# twisted under action of Frobenius maps; +# decompose one symmetric power + +RECOG.DecomposeTensor := function (G, F) +local d, p, f, q, i, factors, g, list, u, w, NmrTries, required, lambda, original, eigenvectors, multiplicity, nmr, E, least, flag, Result, index, first, step, a, x, y, term, lp, Zq, primitiveEle; + + # TODO: Add check whether G is already prepared + G := RECOG.PrepareTensor(G); + + d := G.d; + if d = 2 then + return [false,false]; + fi; + + q := Size(F); + f := Factors(q); + p := f[1]; + f := Size(f); + + if ( ( p mod 2 = 1) and ((f = 2 and (d in [(p - 1)^2, p * (p - 1)])) or (d = p^f))) or (p = 2 and f = 2 and d = 4) then + Print("sl2q: Need special version of DecomposeTensor for these cases \n"); + return fail; + + # TODO: NON-GENERIC CASES. DEAL WITH THEM LATER + + # if d mod p = 0 then + # factors := [[p, Int(d/p)]]; + # elif d mod (p - 1) = 0 then + # factors := [[p - 1, Int(d/(p - 1))]]; + # fi; + # # TODO: Need is Tensor + # flag := RECOG.IsTensor(G, factors); + # if flag then + # list := RECOG.TensorDimensions(G); + # u := list[1]; + # w := list[2]; + # #TODO: Write this in GAP + # #Result := ; + # return [true, Result]; + # else + # return [false,false]; + # fi; + fi; + + NmrTries := 100; + if p = 2 then + required := (q-1); + else + required := (q-1)/2; + fi; + flag := RECOG.RandomElementOfOrderModuleCentre(G, required, NmrTries, true); + g := flag[2]; + flag := flag[1]; + if not(flag) then + Print("sl2q: Failed to find element of order ", required); + return [false, false]; + fi; + + lambda := Eigenvalues(F,g); + #eigenvectors := Eigenvectors(F,g); + original := Size(lambda); + Print("sl2q: Number of eigenvalues is ", original, "\n"); + + lambda := List([1..original], i -> [lambda[i],Size(RECOG.EigenspaceMat(g,lambda[i]))]); + lambda := Filtered(lambda, x -> x[2] = 1); + + Print("sl2q: Number of eigenvalues of multiplicity 1 = ", Size(lambda), "\n"); + if Size(lambda) <= 1 then + Print("sl2q: Too few eigenvalues left \n"); + return [false, false]; + fi; + + # are there multiplicities in eigenvalues? + multiplicity := Size(lambda) < original; + + primitiveEle := PrimitiveElement(GF(q)); + E := List(lambda, x -> LogFFE(primitiveEle,x[1])); + Zq := ZmodnZ(q-1); + E := List(E, x-> ZmodnZObj(ElementsFamily(FamilyObj(Zq)),x)); + + nmr := 0; + + # largest prime dividing d + lp := Maximum(Factors(d)); + + # minimum length of AP if multiplicity among EVs + least := (p + 1)/2; + + # construct arithmetic progressions of length ell, + # where ell is at most p and ell is a multiple of lp; + # if multiplicity-free then ell is proper divisor of d; + # if not, then ell >= least + + for first in [1..Size(E)] do + x := E[first]; + for index in [1..Size(E)] do + if index = first then + continue; + fi; + y := E[index]; + step := y - x; + a := [x, y]; + if d mod Size(a) = 0 then + if (multiplicity and Size(a) >= least) or (not(multiplicity) and Size(a) mod lp = 0) then + flag := RECOG.ProcessSequence(G, g, a); + Result := flag[2]; + flag := flag[1]; + nmr := nmr + 1; + if flag then + RECOG.StoreDetails (G, Result); + Print("sl2q: Arithmetic progression is ", a); + Print("sl2q: Number of calls to FindPoint is", nmr); + return [true, Result]; + fi; + fi; + fi; + + # construct APs of length at most p; if the length + # of the AP properly divides the degree of G, + # then send space to FindPoint + for i in [3..p] do + term := a[i - 1] + step; + if not(term in E) then + break; + fi; + if (term in a) then + break; + fi; + a[i] := term; + if Size(a) > (d/2) then + break; + fi; + if (d mod Size(a)) = 0 then + if (multiplicity and Size(a) >= least) or (not(multiplicity) and (Size(a) mod lp) = 0) then + flag := RECOG.ProcessSequence(G, g, a); + Display("process 1"); + Result := flag[2]; + flag := flag[1]; + nmr := nmr + 1; + if flag then + Print("sl2q: Arithmetic progression is ", a); + Print("sl2q: Number of calls to FindPoint is", nmr); + RECOG.StoreDetails(G, Result); + return [true, Result]; + fi; + fi; + fi; + od; + od; + od; + + return [false,false]; + +end; + + + +################################################################################################### +################################################################################################### +######## Symmetric Powers and Twisted Tensor Products ############################################# +################################################################################################### +################################################################################################### + + + +RECOG.SymmetricPowerSL2 := function (q, n) +local G, fld, gens, res, g; + + gens := GeneratorsOfGroup(SL(2,q)); + res := []; + fld := GF(q); + for g in gens do + Add(res,RECOG.SymPowRepSL2(q,n,g)); + od; + return GroupByGenerators(res); + #return ActionGroup (A); +end; + +# given matrix x, return result of applying Frobenius automorphism +# x[i][j] --> x[i][j]^n to it +RECOG.GammaLMatrix:=function(x,d,n) +local G,i,j,y; + + if not(IsMutable(x)) then + x := MutableCopyMat(x); + fi; + + for i in [1..d] do + for j in [1..d] do + x[i,j] := x[i,j]^n; + od; + od; + return x; +end; + + + +RECOG.TwistedSymmetricPower:=function(q,s,f) +local F,Gens,S,Twisted,p,i,d; + F:=GF(q); + p:=Characteristic(F); + S:=RECOG.SymmetricPowerSL2(q,s); + Gens:=GeneratorsOfGroup(S); + d := NumberRows(Gens[1]); + Twisted:=List([1..Size(Gens)],i->RECOG.GammaLMatrix(Gens[i],d,p^f)); + return GroupByGenerators(Twisted); +end; + + + +# tensor product of twisted symmetric powers defined +# over GF (p^e); s lists the symmetric powers, +# f is the twisting via powers of the Frobenius +# automorphism to be applied to each symmetric power +RECOG.SymmetricPowerExample:=function(p,e,s,f) +local F,q,G,L,T2,i,gens,j,gens2; + F:=GF(p^e); + q:= p^e; + L:=SL(2,F); + G:=RECOG.TwistedSymmetricPower(q,s[1],f[1]); + if Size(s) = 1 then + return G; + fi; + for i in [2..Size(s)] do + T2:=RECOG.TwistedSymmetricPower(q,s[i],f[i]); + gens := GeneratorsOfGroup(G); + gens2 := [1..Size(gens)]; + for j in [1..Size(gens)] do + gens2[j] := KroneckerProduct(gens[j],GeneratorsOfGroup(T2)[j]); + od; + G:=GroupByGenerators(gens2); + od; + return G; +end; diff --git a/gap/projective/constructive_recognition/main.gi b/gap/projective/constructive_recognition/main.gi new file mode 100644 index 00000000..d605ebf7 --- /dev/null +++ b/gap/projective/constructive_recognition/main.gi @@ -0,0 +1,31 @@ +############################################################################# +## +## This file is part of recog, a package for the GAP computer algebra system +## which provides a collection of methods for the constructive recognition +## of groups. +## +## This files's authors include Daniel Rademacher. +## +## Copyright of recog belongs to its developers whose names are too numerous +## to list here. Please refer to the COPYRIGHT file for details. +## +## SPDX-License-Identifier: GPL-3.0-or-later +## +############################################################################# + + + +############################################################################# +############################################################################# +######## Main function for all classical groups in natural reps ############# +############################################################################# +############################################################################# + + + +# TODO: main function +RECOG.ConstructiveRecognitionClassicalGroupsNaturalRepresentation := function(G) + + # TODO: Check which classical group is G and apply the correct subfuntion + +end; diff --git a/gap/projective/classicalnatural.gi b/gap/projective/constructive_recognition/utils/achieve.gi similarity index 71% rename from gap/projective/classicalnatural.gi rename to gap/projective/constructive_recognition/utils/achieve.gi index c53dd2fe..8cabd4f9 100644 --- a/gap/projective/classicalnatural.gi +++ b/gap/projective/constructive_recognition/utils/achieve.gi @@ -4,7 +4,7 @@ ## which provides a collection of methods for the constructive recognition ## of groups. ## -## This files's authors include Max Neunhöffer, Ákos Seress. +## This files's authors include Daniel Rademacher, Max Neunhöffer, Ákos Seress. ## ## Copyright of recog belongs to its developers whose names are too numerous ## to list here. Please refer to the COPYRIGHT file for details. @@ -17,24 +17,6 @@ ## ############################################################################# -InstallMethod( CharacteristicPolynomial, "for a memory element matrix", - [ IsMatrix and IsObjWithMemory ], - function(m) - return CharacteristicPolynomial(m!.el); - end ); - -InstallOtherMethod( \-, "for two memory elements", - [ IsMatrix and IsObjWithMemory, IsMatrix and IsObjWithMemory ], - function(m,n) - return m!.el - n!.el; - end ); - -InstallMethod( Eigenspaces, "for a field and a memory element matrix", - [ IsField, IsMatrix and IsObjWithMemory ], - function( f, m ) - return Eigenspaces(f,m!.el); - end ); - # Obsolete stuff? # RECOG.RelativePrimeToqm1Part := function(q,n) @@ -435,39 +417,6 @@ InstallMethod( Eigenspaces, "for a field and a memory element matrix", # return subgplist; # end; -RECOG.FindStdGensUsingBSGS := function(g,stdgens,projective,large) - # stdgens generators for the matrix group g - # returns an SLP expressing stdgens in the generators of g - # set projective to true for projective mode - # set large to true if we should not bother finding nice base points! - local S,dim,gens,gm,i,l,strong; - dim := DimensionOfMatrixGroup(g); - if IsObjWithMemory(GeneratorsOfGroup(g)[1]) then - gm := GroupWithMemory(StripMemory(GeneratorsOfGroup(g))); - else - gm := GroupWithMemory(g); - fi; - if HasSize(g) then SetSize(gm,Size(g)); fi; - if large then - S := StabilizerChain(gm,rec( Projective := projective, - Cand := rec( points := One(g), - ops := ListWithIdenticalEntries(dim, OnLines) ) ) ); - else - S := StabilizerChain(gm,rec( Projective := projective ) ); - fi; - strong := ShallowCopy(StrongGenerators(S)); - ForgetMemory(S); - l := List(stdgens,x->SiftGroupElementSLP(S,x)); - gens := EmptyPlist(Length(stdgens)); - for i in [1..Length(stdgens)] do - if not l[i].isone then - return fail; - fi; - Add(gens,ResultOfStraightLineProgram(l[i].slp,strong)); - od; - return SLPOfElms(gens); -end; - RECOG.ResetSLstd := function(r) r.left := One(r.a); r.right := One(r.a); @@ -1100,113 +1049,7 @@ end; # return rec( slpstd := SLPOfElms(std.all), bas := bas, basi := basi ); # end; -# TODO: which algorithm is this? reference? -RECOG.FindStdGens_SL := function(sld,f) - # gens of sld must be gens for SL(d,q) in its natural rep with memory - # This function calls RECOG.SLn_constructsl2 and then extends - # the basis to a basis of the full row space and calls - # RECOG.SLn_UpStep often enough. Finally it returns an slp such - # that the SL(d,q) standard generators with respect to this basis are - # expressed by the slp in terms of the original generators of sld. - local V,b,bas,basi,basit,d,data,ext,fakegens,id,nu,nu2,p,q,resl2,sl2,sl2gens, - sl2gensf,sl2genss,sl2stdf,slp,slpsl2std,slptosl2,st,std,stdgens,i,ex; - - # Some setup: - p := Characteristic(f); - q := Size(f); - ext := DegreeOverPrimeField(f); - d := DimensionOfMatrixGroup(sld); - if not IsObjWithMemory(GeneratorsOfGroup(sld)[1]) then - sld := GroupWithMemory(sld); - fi; - - # First find an SL2 with the space it acts on: - Info(InfoRecog,2,"Finding an SL2..."); - data := RECOG.SLn_constructsl2(sld,d,q); - - bas := ShallowCopy(BasisVectors(Basis(data[2]))); - sl2 := data[1]; - slptosl2 := SLPOfElms(GeneratorsOfGroup(sl2)); - sl2gens := StripMemory(GeneratorsOfGroup(sl2)); - V := data[2]; - b := Basis(V,bas); - sl2genss := List(sl2gens,x->RECOG.LinearAction(b,f,x)); - - if q in [2,3,4,5,9] then - Info(InfoRecog,2,"In fact found an SL4..."); - stdgens := RECOG.MakeSL_StdGens(p,ext,4,4).all; - slpsl2std := RECOG.FindStdGensUsingBSGS(Group(sl2genss),stdgens, - false,false); - nu := List(sl2gens,x->NullspaceMat(x-One(x))); - ex := SumIntersectionMat(nu[1],nu[2])[2]; - for i in [3..Length(nu)] do - ex := SumIntersectionMat(nu[3],ex); - od; - Append(bas,ex); - ConvertToMatrixRep(bas,q); - basi := bas^-1; - else - # Now compute the natural SL2 action and run constructive recognition: - Info(InfoRecog,2, - "Recognising this SL2 constructively in 2 dimensions..."); - sl2genss := GeneratorsWithMemory(sl2genss); - if IsEvenInt(q) then - resl2 := RECOG.RecogniseSL2NaturalEvenChar(Group(sl2genss),f,false); - else - resl2 := RECOG.RecogniseSL2NaturalOddCharUsingBSGS(Group(sl2genss),f); - fi; - slpsl2std := SLPOfElms(resl2.all); - bas := resl2.bas * bas; - # We need the actual transvections: - slp := SLPOfElms([resl2.s[1],resl2.t[1]]); - st := ResultOfStraightLineProgram(slp, - StripMemory(GeneratorsOfGroup(sl2))); - - # Extend basis by something invariant under SL2: - id := IdentityMat(d,f); - nu := NullspaceMat(StripMemory(st[1]-id)); - nu2 := NullspaceMat(StripMemory(st[2]-id)); - Append(bas,SumIntersectionMat(nu,nu2)[2]); - ConvertToMatrixRep(bas,q); - basi := bas^-1; - fi; - # Now set up fake generators for keeping track what we do: - fakegens := ListWithIdenticalEntries(Length(GeneratorsOfGroup(sld)),1); - fakegens := GeneratorsWithMemory(fakegens); - sl2gensf := ResultOfStraightLineProgram(slptosl2,fakegens); - sl2stdf := ResultOfStraightLineProgram(slpsl2std,sl2gensf); - std := rec( f := f, d := d, n := 2, bas := bas, basi := basi, - sld := sld, sldf := fakegens, slnstdf := sl2stdf, - p := p, ext := ext ); - Info(InfoRecog,2,"Going up to SL_d again..."); - while std.n < std.d do - RECOG.SLn_UpStep(std); - od; - return rec( slpstd := SLPOfElms(std.slnstdf), - bas := std.bas, basi := std.basi ); -end; - -RECOG.RecogniseSL2NaturalOddCharUsingBSGS := function(g,f) - local ext,p,q,res,slp,std; - p := Characteristic(f); - ext := DegreeOverPrimeField(f); - q := Size(f); - std := RECOG.MakeSL_StdGens(p,ext,2,2); - slp := RECOG.FindStdGensUsingBSGS(g,std.all,false,true); - if slp = fail then - return fail; - fi; - res := rec( g := g, one := One(f), One := One(g), f := f, q := q, - p := p, ext := ext, d := 2, bas := IdentityMat(2,f), - basi := IdentityMat(2,f) ); - res.all := ResultOfStraightLineProgram(slp,GeneratorsOfGroup(g)); - res.s := res.all{[1..ext]}; - res.t := res.all{[ext+1..2*ext]}; - res.a := res.all[2*ext+1]; - res.b := res.all[2*ext+2]; - return res; -end; RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) # f a finite field, g equal to SL(2,Size(f)), t either an involution @@ -1229,7 +1072,8 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) fi; if torig = false then # if no involution t has been given, compute one, using Proposition 4 from - # [KK15]. + # "Black box groups isomorphic to PGL(2,2^e)" by Kantor & Kassabov, + # Journal of Algebra, 421 (2015) 16–26. repeat am:=PseudoRandom(g); until not IsOneProjective(am); @@ -1257,7 +1101,7 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) ch := Factors(CharacteristicPolynomial(f,f,t,1)); if Length(ch) <> 2 or ch[1] <> ch[2] then - ErrorNoReturn("matrix is not triagonalizable - this should never happen!"); + Error("matrix is not triagonalizable - this should never happen!"); fi; one := OneMutable(t); @@ -1294,12 +1138,12 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) xm := o[j]; j := j + 1; c := Comm(tm,xm); - until not IsOne(c^2); + until not(IsOne(c^2)); xm := xm * c^(((q-1)*(q+1)-1)/2); x := StripMemory(xm); xb := bas*x*bas^-1; co := Coefficients(can,xb[2,1]); - until not IsContainedInSpan(mb,co); + until not(IsContainedInSpan(mb,co)); CloseMutableBasis(mb,co); Add(tt,x); Add(ttm,xm); @@ -1333,20 +1177,20 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) xm := o[j]; j := j + 1; x := MutableCopyMat(bas*StripMemory(xm)*bas^-1); - until not IsZero(x[1,2]); + until not(IsZero(x[1,2])); - if not IsOne(x[2,2]) then + if not(IsOne(x[2,2])) then el := (One(f)-x[2,2])/x[1,2]; co := Coefficients(can,el) * mati; for i in [1..Length(co)] do - if not IsZero(co[i]) then + if not(IsZero(co[i])) then xm := ttm[i] * xm; fi; od; x[2] := x[2] + x[1] * el; if x <> bas*StripMemory(xm)*bas^-1 then # FIXME: sometimes triggered by RecognizeGroup(GL(2,16)); - ErrorNoReturn("!!!"); + Error("!!!"); fi; fi; # now x[2,2] is equal to One(f) @@ -1358,7 +1202,7 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) el := x[2,1]; co2 := Coefficients(can,el) * mati; for i in [1..Length(co2)] do - if not IsZero(co2[i]) then + if not(IsZero(co2[i])) then xm := xm * ttm[i]; fi; od; @@ -1424,7 +1268,7 @@ end; # j := p[2]-1; # while j >= 0 do # z := y^(s/p[1]^(p[2]-j)); -# if not IsOne(z) then break; fi; +# if not(IsOne(z)) then break; fi; # j := j - 1; # od; # o := o * p[1]^(j+1); @@ -1441,9 +1285,7 @@ RECOG.GuessProjSL2ElmOrder := function(x,f) fi; if p > 2 then y := x^p; - if IsOneProjective(y) then - return p; - fi; + if IsOneProjective(y) then return p; fi; fi; if IsOneProjective(x^(q-1)) then facts := Collected(FactInt(q-1:cheap)[1]); @@ -1461,7 +1303,7 @@ RECOG.GuessProjSL2ElmOrder := function(x,f) j := p[2]-1; while j >= 0 do z := y^(s/p[1]^(p[2]-j)); - if not IsOneProjective(z) then break; fi; + if not(IsOneProjective(z)) then break; fi; j := j - 1; od; o := o * p[1]^(j+1); @@ -1480,34 +1322,24 @@ RECOG.IsThisSL2Natural := function(gens,f) CheckElm := function(x) local o; o := RECOG.GuessProjSL2ElmOrder(x,f); - if o in [1,2] then - return false; - fi; + if o in [1,2] then return false; fi; if o > 5 then if notA5 = false then Info(InfoRecog,4,"SL2: Group is not A5"); fi; notA5 := true; - if seenqp1 and seenqm1 then - return true; - fi; - fi; - if o = p or o <= 5 then - return false; + if seenqp1 and seenqm1 then return true; fi; fi; + if o = p or o <= 5 then return false; fi; if (q+1) mod o = 0 then - if not seenqp1 then + if not(seenqp1) then Info(InfoRecog,4,"SL2: Found element of order dividing q+1."); seenqp1 := true; - if seenqm1 and notA5 then - return true; - fi; + if seenqm1 and notA5 then return true; fi; fi; else - if not seenqm1 then + if not(seenqm1) then Info(InfoRecog,4,"SL2: Found element of order dividing q-1."); seenqm1 := true; - if seenqp1 and notA5 then - return true; - fi; + if seenqp1 and notA5 then return true; fi; fi; fi; return false; @@ -1520,7 +1352,7 @@ RECOG.IsThisSL2Natural := function(gens,f) q := Size(f); p := Characteristic(f); - # For small q, comput the order of the group via a stabilizer chain. + # For small q, compute the order of the group via a stabilizer chain. # Note that at this point we are usually working projective, and thus # scalars are factored out "implicitly". Thus the generators we are # looking at may generate a group which only contains SL2 as a subgroup. @@ -1528,7 +1360,8 @@ RECOG.IsThisSL2Natural := function(gens,f) Info(InfoRecog,4,"SL2: Computing stabiliser chain."); S := StabilizerChain(Group(gens)); Info(InfoRecog,4,"SL2: size is ",Size(S)); - return Size(S) mod (q*(q-1)*(q+1)) = 0; + # return Size(S) mod (q*(q-1)*(q+1)) = 0; + return Size(S) = (q*(q-1)*(q+1)); fi; seenqp1 := false; @@ -1536,9 +1369,7 @@ RECOG.IsThisSL2Natural := function(gens,f) notA5 := false; for i in [1..Length(gens)] do - if CheckElm(gens[i]) then - return true; - fi; + if CheckElm(gens[i]) then return true; fi; od; CheckElm(gens[1]*gens[2]); if Length(gens) >= 3 then @@ -1554,9 +1385,7 @@ RECOG.IsThisSL2Natural := function(gens,f) for i in [1..l-1] do for j in [i+1..l] do x := Comm(gens[i],gens[j]); - if CheckElm(x) then - return true; - fi; + if CheckElm(x) then return true; fi; Add(coms,x); od; od; @@ -1566,24 +1395,21 @@ RECOG.IsThisSL2Natural := function(gens,f) a := RECOG.RandomSubproduct(gens,rec()); b := RECOG.RandomSubproduct(gens,rec()); x := Comm(a,b); - if CheckElm(x) then - return true; - fi; + if CheckElm(x) then return true; fi; Add(coms,x); od; fi; - if ForAll(coms, IsDiagonalMat) then + if ForAll(coms,c->RECOG.IsScalarMat(c) <> false) then Info(InfoRecog,4,"SL2: Group is soluble, commutators are central"); return false; fi; Info(InfoRecog,4,"SL2: Computing normal closure..."); clos := FastNormalClosure(gens,coms,5); for i in [Length(coms)+1..Length(clos)] do - if CheckElm(clos[i]) then - return true; - fi; + if CheckElm(clos[i]) then return true; fi; od; - if ForAll(clos{[Length(coms)+1..Length(clos)]}, IsDiagonalMat) then + if ForAll(clos{[Length(coms)+1..Length(clos)]}, + c->RECOG.IsScalarMat(c) <> false) then Info(InfoRecog,4,"SL2: Group is soluble, derived subgroup central"); return false; fi; @@ -1605,750 +1431,10 @@ RECOG.IsThisSL2Natural := function(gens,f) return false; end; -# The going down method: -#Version 1.2 - -# finds first element of a list that is relative prime to all others -# input: list=[SL(d,q), d, q, SL(n,q)] acting as a subgroup of some big SL(n,q) -# output: list=[rr, dd] for a ppd(2*dd;q)-element rr -RECOG.SLn_godown:=function(list) - local d, first, q, g, gg, i, r, pol, factors, degrees, newdim, power, rr, ss, - newgroup, colldegrees, exp, count; - - first:=function(list) - local i; - - for i in [1..Length(list)] do - if list[i]>1 and Gcd(list[i],Product(list)/list[i])=1 then - return list[i]; - fi; - od; - - return fail; - end; - - g:=list[1]; - d:=list[2]; - q:=list[3]; - gg:=list[4]; - - Info(InfoRecog,2,"Dimension: ",d); - #find an element with irreducible action of relative prime dimension to - #all other invariant subspaces - #count is just safety, if things go very bad - count:=0; - - repeat - count:=count+1; - if InfoLevel(InfoRecog) >= 3 then Print(".\c"); fi; - r:=PseudoRandom(g); - pol:=CharacteristicPolynomial(r); - factors:=Factors(pol); - degrees:=AsSortedList(List(factors,Degree)); - newdim:=first(degrees); - until (count>10) or (newdim <> fail and newdim<=Maximum(2,d/4)); - - if count>10 then - return fail; - fi; - - # raise r to a power so that acting trivially outside one invariant subspace - degrees:=Filtered(degrees, x->x<>newdim); - colldegrees:=Collected(degrees); - power:=Lcm(List(degrees, x->q^x-1))*q; - # power further to cancel q-part of element order - if degrees[1]=1 then - exp:=colldegrees[1][2]-(DimensionOfMatrixGroup(gg)-d); - if exp>0 then - power:=power*q^exp; - fi; - fi; - rr:=r^power; - - #conjugate rr to hopefully get a smaller dimensional SL - #ss:=rr^PseudoRandom(gg); - #newgroup:=Group(rr,ss); - - return [rr,newdim]; -end; - -# input is (group,dimension,q) -# output is a group element acting irreducibly in two dimensions, and fixing -# a (dimension-2)-dimensional subspace -RECOG.SLn_constructppd2:=function(g,dim,q) - local out, list ; - - list:=[g,dim,q,g]; - repeat - out:=RECOG.SLn_godown(list); - if out=fail or out[1]*out[1]=One(out[1]) then - if InfoLevel(InfoRecog) >= 3 then Print("B\c"); fi; - list:=[g,dim,q,g]; - out:=fail; - else - if out[2]>2 then - list:=[Group(out[1],out[1]^PseudoRandom(g)),2*out[2],q,g]; - fi; - fi; - until out<>fail and out[2]=2; - - return out[1]; - -end; - -RECOG.SLn_constructsl4:=function(g,dim,q,r) - local s,h,count,readydim4,readydim3,ready,u,orderu, - nullr,nulls,nullspacer,nullspaces,int,intbasis,nullintbasis, - newu,newbasis,newbasisinv,newr,news,outputu,mat,i,shorts,shortr; - nullr:=NullspaceMat(r-One(r)); - nullspacer:=VectorSpace(GF(q),nullr); - mat:=One(r); - ready:=false; - repeat - s:=r^PseudoRandom(g); - nulls:=NullspaceMat(s-One(s)); - nullspaces:=VectorSpace(GF(q),nulls); - int:=Intersection(nullspacer,nullspaces); - intbasis:=Basis(int); - newbasis:=[]; - for i in [1..Length(intbasis)] do - Add(newbasis,intbasis[i]); - od; - i:=0; - repeat - i:=i+1; - if not mat[i] in int then - Add(newbasis,mat[i]); - int:=VectorSpace(GF(q),newbasis); - fi; - until Dimension(int)=dim; - ConvertToMatrixRep(newbasis); - newbasisinv:=newbasis^(-1); - newr:=newbasis*r*newbasisinv; - news:=newbasis*s*newbasisinv; - - #shortr, shorts do not need memory - #we shall throw away the computations in h - #check that we have SL(4,q), by non-constructive recognition - shortr:=newr{[dim-3..dim]}{[dim-3..dim]}; - shorts:=news{[dim-3..dim]}{[dim-3..dim]}; - h:=Group(shortr,shorts); - count:=0; - readydim4:=false; - readydim3:=false; - repeat - u:=PseudoRandom(h); - orderu:=Order(u); - if orderu mod ((q^4-1)/(q-1)) = 0 then - readydim4:=true; - elif Gcd(orderu,(q^2+q+1)/Gcd(3,q-1))>1 then - readydim3:=true; - fi; - if readydim4 = true and readydim3 = true then - ready:=true; - break; - fi; - count:=count+1; - until count=30; - until ready=true; - - return Group(r,s); -end; - - -#g=SL(d,q), given as a subgroup of SL(dim,q) -#output: [SL(2,q), and a basis for the 2-dimensional subspace where it acts -RECOG.SLn_godownfromd:=function(g,q,d,dim) - local y,yy,ready,order,es,dims,subsp,z,x,a,b,c,h,vec,vec2, - pol,factors,degrees,comm1,comm2,comm3,image,basis,action,vs,readyqpl1, - readyqm1,count,u,orderu; - - repeat - ready:=false; - y:=PseudoRandom(g); - pol:=CharacteristicPolynomial(y); - factors:=Factors(pol); - degrees:=List(factors,Degree); - if d-1 in degrees then - order:=Order(y); - if order mod (q-1)=0 then - yy:=y^(order/(q-1)); - else - yy:=One(y); - fi; - if not IsOne(yy) then - es:= Eigenspaces(GF(q),yy); - dims:=List(es,Dimension); - if IsSubset(Set([1,d-1,dim-d]),Set(dims)) and - (1 in Set(dims)) then - es:=Filtered(es,x->Dimension(x)=1); - vec:=Basis(es[1])[1]; - if vec*yy=vec then - vec:=Basis(es[2])[1]; - fi; - repeat - z:=PseudoRandom(g); - x:=yy^z; - a:=Comm(x,yy); - b:=a^yy; - c:=a^x; - comm1:= Comm(a,c); - comm2:=Comm(a,b); - comm3:=Comm(b,c); - if comm1<>One(a) and comm2<>One(a) and - comm3<>One(a) and Comm(comm1,comm2)<>One(a) then - vec2:=vec*z; - vs:=VectorSpace(GF(q),[vec,vec2]); - basis:=Basis(vs); - #check that the action in 2 dimensions is SL(2,q) - #by non-constructive recognition, finding elements of - #order (q-1) and (q+1) - #we do not need memory in the group image - action:=List([a,b,c],x->RECOG.LinearAction(basis,q,x)); - image:=Group(action); - count:=0; - readyqpl1:=false; - readyqm1:=false; - repeat - u:=PseudoRandom(image); - orderu:=Order(u); - if orderu = q-1 then - readyqm1:=true; - elif orderu = q+1 then - readyqpl1:=true; - fi; - if readyqm1 = true and readyqpl1 = true then - ready:=true; - break; - fi; - count:=count+1; - until count=20; - fi; - until ready=true; - fi; - fi; - fi; - until ready; - - h:=Group(a,b,c); - subsp:=VectorSpace(GF(q),[vec,vec2]); - return [h,subsp]; - -end; - -#going down from 4 to 2 dimensions, when q=2,3,4,5,9 -#just construct the 4-dimensional invariant space and generators -#for the group acting on it -RECOG.SLn_exceptionalgodown:=function(h,q,dim) - local basis, v, vs, i, gen; - - vs:=VectorSpace(GF(q),One(h)); - basis:=[]; - repeat - if InfoLevel(InfoRecog) >= 3 then Print("C"); fi; - for i in [1..4] do - v:=PseudoRandom(vs); - for gen in GeneratorsOfGroup(h) do - Add(basis,v*gen-v); - od; - od; - basis:=ShallowCopy(SemiEchelonMat(basis).vectors); - until Length(basis)=4; - return [h,VectorSpace(GF(q),basis)]; -end; - - -RECOG.SLn_constructsl2:=function(g,d,q) - local r,h; - - r:=RECOG.SLn_constructppd2(g,d,q); - h:=RECOG.SLn_constructsl4(g,d,q,r); - if not (q in [2,3,4,5,9]) then - return RECOG.SLn_godownfromd(h,q,4,d); - else - return RECOG.SLn_exceptionalgodown(h,q,d); - # return ["sorry only SL(4,q)",h]; - fi; -end; - -# Now the going up code: - -RECOG.LinearAction := function(bas,field,el) - local mat,vecs; - if IsGroup(el) then - return Group(List(GeneratorsOfGroup(el), - x->RECOG.LinearAction(bas,field,x))); - fi; - if IsBasis(bas) then - vecs := BasisVectors(bas); - else - vecs := bas; - bas := Basis(VectorSpace(field,bas),bas); - fi; - mat := List(vecs,v->Coefficients(bas,v*el)); - ConvertToMatrixRep(mat,field); - return mat; -end; - -RECOG.SLn_UpStep := function(w) - # w has components: - # d : size of big SL - # n : size of small SL - # slnstdf : fakegens for SL_n standard generators - # bas : current base change, first n vectors are where SL_n acts - # rest of vecs are invariant under SL_n - # basi : current inverse of bas - # sld : original group with memory generators, PseudoRandom - # delivers random elements - # sldf : fake generators to keep track of what we are doing - # f : field - # The following are filled in automatically if not already there: - # p : characteristic - # ext : q=p^ext - # One : One(slnstdf[1]) - # can : CanonicalBasis(f) - # canb : BasisVectors(can) - # transh : fakegens for the "horizontal" transvections n,i for 1<=i<=n-1 - # entries can be unbound in which case they are made from slnstdf - # transv : fakegens for the "vertical" transvections i,n for 1<=i<=n-1 - # entries can be unbound in which case they are made from slnstdf - # - # We keep the following invariants (going from n -> n':=2n-1) - # bas, basi is a base change to the target base - # slnstdf are SLPs to reach standard generators of SL_n from the - # generators of sld - local DoColOp_n,DoRowOp_n,FixSLn,Fixc,MB,Vn,Vnc,aimdim,c,c1,c1f,cf,cfi, - ci,cii,coeffs,flag,i,id,int1,int3,j,k,lambda,list,mat,newbas,newbasf, - newbasfi,newbasi,newdim,newpart,perm,pivots,pivots2,pos,pow,s,sf, - slp,std,sum1,tf,trans,transd,transr,v,vals,zerovec; - - Info(InfoRecog,3,"Going up: ",w.n," (",w.d,")..."); - - # Before we begin, we upgrade the data structure with a few internal - # things: - - if not IsBound(w.can) then w.can := CanonicalBasis(w.f); fi; - if not IsBound(w.canb) then w.canb := BasisVectors(w.can); fi; - if not IsBound(w.One) then w.One := One(w.slnstdf[1]); fi; - if not IsBound(w.transh) then w.transh := []; fi; - if not IsBound(w.transv) then w.transv := []; fi; - # Update our cache of *,n and n,* transvections because we need them - # all over the place: - std := RECOG.InitSLstd(w.f,w.n, - w.slnstdf{[1..w.ext]}, - w.slnstdf{[w.ext+1..2*w.ext]}, - w.slnstdf[2*w.ext+1], - w.slnstdf[2*w.ext+2]); - for i in [1..w.n-1] do - for k in [1..w.ext] do - pos := (i-1)*w.ext + k; - if not IsBound(w.transh[pos]) then - RECOG.ResetSLstd(std); - RECOG.DoColOp_SL(false,w.n,i,w.canb[k],std); - w.transh[pos] := std.right; - fi; - if not IsBound(w.transv[pos]) then - RECOG.ResetSLstd(std); - RECOG.DoRowOp_SL(false,i,w.n,w.canb[k],std); - w.transv[pos] := std.left; - fi; - od; - od; - Unbind(std); - - # Now we can define two helper functions: - DoColOp_n := function(el,i,j,lambda,w) - # This adds lambda times the i-th column to the j-th column. - # Note that either i or j must be equal to n! - local coeffs,k; - coeffs := IntVecFFE(Coefficients(w.can,lambda)); - if i = w.n then - for k in [1..w.ext] do - if not IsZero(coeffs[k]) then - if IsOne(coeffs[k]) then - el := el * w.transh[(j-1)*w.ext+k]; - else - el := el * w.transh[(j-1)*w.ext+k]^coeffs[k]; - fi; - fi; - od; - elif j = w.n then - for k in [1..w.ext] do - if not IsZero(coeffs[k]) then - if IsOne(coeffs[k]) then - el := el * w.transv[(i-1)*w.ext+k]; - else - el := el * w.transv[(i-1)*w.ext+k]^coeffs[k]; - fi; - fi; - od; - else - ErrorNoReturn("either i or j must be equal to n"); - fi; - return el; - end; - DoRowOp_n := function(el,i,j,lambda,w) - # This adds lambda times the j-th row to the i-th row. - # Note that either i or j must be equal to n! - local coeffs,k; - coeffs := IntVecFFE(Coefficients(w.can,lambda)); - if j = w.n then - for k in [1..w.ext] do - if not IsZero(coeffs[k]) then - if IsOne(coeffs[k]) then - el := w.transv[(i-1)*w.ext+k] * el; - else - el := w.transv[(i-1)*w.ext+k]^coeffs[k] * el; - fi; - fi; - od; - elif i = w.n then - for k in [1..w.ext] do - if not IsZero(coeffs[k]) then - if IsOne(coeffs[k]) then - el := w.transh[(j-1)*w.ext+k] * el; - else - el := w.transh[(j-1)*w.ext+k]^coeffs[k] * el; - fi; - fi; - od; - else - ErrorNoReturn("either i or j must be equal to n"); - fi; - return el; - end; - - # Here everything starts, some more preparations: - - # We compute exclusively in our basis, so we occasionally need an - # identity matrix: - id := IdentityMat(w.d,w.f); - FixSLn := VectorSpace(w.f,id{[w.n+1..w.d]}); - Vn := VectorSpace(w.f,id{[1..w.n]}); - - # First pick an element in SL_n with fixed space of dimension d-n+1: - # We already have an SLP for an n-1-cycle: it is one of the std gens. - # For n=2 we use a transvection for this purpose. - if w.n > 2 then - if IsOddInt(w.n) then - if w.p > 2 then - s := id{Concatenation([1,w.n],[2..w.n-1],[w.n+1..w.d])}; - ConvertToMatrixRepNC(s,w.f); - if IsOddInt(w.n) then s[2] := -s[2]; fi; - sf := w.slnstdf[2*w.ext+2]; - else # in even characteristic we take the n-cycle: - s := id{Concatenation([w.n],[1..w.n-1],[w.n+1..w.d])}; - ConvertToMatrixRepNC(s,w.f); - sf := w.slnstdf[2*w.ext+1]; - fi; - else - ErrorNoReturn("this program only works for odd n or n=2"); - fi; - else - # In this case the n-1-cycle is the identity, so we take a transvection: - s := MutableCopyMat(id); - s[1,2] := One(w.f); - sf := w.slnstdf[1]; - fi; - - # Find a good random element: - w.count := 0; - aimdim := Minimum(2*w.n-1,w.d); - newdim := aimdim - w.n; - while true do # will be left by break - while true do # will be left by break - if InfoLevel(InfoRecog) >= 3 then Print(".\c"); fi; - w.count := w.count + 1; - c1 := PseudoRandom(w.sld); - slp := SLPOfElm(c1); - c1f := ResultOfStraightLineProgram(slp,w.sldf); - # Do the base change into our basis: - c1 := w.bas * c1 * w.basi; - c := s^c1; - cf := sf^c1f; - cfi := cf^-1; - # Now check that Vn + Vn*s^c1 has dimension 2n-1: - Vnc := VectorSpace(w.f,c{[1..w.n]}); - sum1 := ClosureLeftModule(Vn,Vnc); - if Dimension(sum1) = aimdim then - Fixc := VectorSpace(w.f,NullspaceMat(c-One(c))); - int1 := Intersection(Fixc,Vn); - for i in [1..Dimension(int1)] do - v := Basis(int1)[i]; - if not IsZero(v[w.n]) then break; fi; - od; - if IsZero(v[w.n]) then - Info(InfoRecog,2,"Ooops: Component n was zero!"); - continue; - fi; - v := v / v[w.n]; # normalize to 1 in position n - Assert(1,v*c=v); - ci := c^-1; - break; - fi; - od; - - # Now we found our aimdim-dimensional space W. Since SL_n - # has a d-n-dimensional fixed space W_{d-n} and W contains a complement - # of that fixed space, the intersection of W and W_{d-n} has dimension - # newdim. - - # Change basis: - newpart := ExtractSubMatrix(c,[1..w.n-1],[1..w.d]); - # Clean out the first n entries to go to the fixed space of SL_n: - zerovec := Zero(newpart[1]); - for i in [1..w.n-1] do - CopySubVector(zerovec,newpart[i],[1..w.n],[1..w.n]); - od; - MB := MutableBasis(w.f,[],zerovec); - i := 1; - pivots := EmptyPlist(newdim); - while i <= Length(newpart) and NrBasisVectors(MB) < newdim do - if not IsContainedInSpan(MB,newpart[i]) then - Add(pivots,i); - CloseMutableBasis(MB,newpart[i]); - fi; - i := i + 1; - od; - newpart := newpart{pivots}; - newbas := Concatenation(id{[1..w.n-1]},[v],newpart); - if 2*w.n-1 < w.d then - int3 := Intersection(FixSLn,Fixc); - if Dimension(int3) <> w.d-2*w.n+1 then - Info(InfoRecog,2,"Ooops, FixSLn \cap Fixc wrong dimension"); - continue; - fi; - Append(newbas,BasisVectors(Basis(int3))); - fi; - ConvertToMatrixRep(newbas,Size(w.f)); - newbasi := newbas^-1; - if newbasi = fail then - Info(InfoRecog,2,"Ooops, Fixc intersected too much, we try again"); - continue; - fi; - ci := newbas * ci * newbasi; - cii := ExtractSubMatrix(ci,[w.n+1..aimdim],[1..w.n-1]); - ConvertToMatrixRep(cii,Size(w.f)); - cii := TransposedMat(cii); - # The rows of cii are now what used to be the columns, - # their length is newdim, we need to span the full newdim-dimensional - # row space and need to remember how: - zerovec := Zero(cii[1]); - MB := MutableBasis(w.f,[],zerovec); - i := 1; - pivots2 := EmptyPlist(newdim); - while i <= Length(cii) and NrBasisVectors(MB) < newdim do - if not IsContainedInSpan(MB,cii[i]) then - Add(pivots2,i); - CloseMutableBasis(MB,cii[i]); - fi; - i := i + 1; - od; - if Length(pivots2) = newdim then - cii := cii{pivots2}^-1; - ConvertToMatrixRep(cii,w.f); - c := newbas * c * newbasi; - w.bas := newbas * w.bas; - w.basi := w.basi * newbasi; - break; - fi; - Info(InfoRecog,2,"Ooops, no nice bottom..."); - # Otherwise simply try again - od; - Info(InfoRecog,2," found c1 and c."); - # Now SL_n has to be repaired according to the base change newbas: - - # Now write this matrix newbas as an SLP in the standard generators - # of our SL_n. Then we know which generators to take for our new - # standard generators, namely newbas^-1 * std * newbas. - newbasf := w.One; - for i in [1..w.n-1] do - if not IsZero(v[i]) then - newbasf := DoColOp_n(newbasf,w.n,i,v[i],w); - fi; - od; - newbasfi := newbasf^-1; - w.slnstdf := List(w.slnstdf,x->newbasfi * x * newbasf); - # Now update caches: - w.transh := List(w.transh,x->newbasfi * x * newbasf); - w.transv := List(w.transv,x->newbasfi * x * newbasf); - - # Now consider the transvections t_i: - # t_i : w.bas[j] -> w.bas[j] for j <> i and - # t_i : w.bas[i] -> w.bas[i] + ww - # We want to modify (t_i)^c such that it fixes w.bas{[1..w.n]}: - trans := []; - for i in pivots2 do - # This does t_i - for lambda in w.canb do - # This does t_i : v_j -> v_j + lambda * v_n - tf := w.One; - tf := DoRowOp_n(tf,i,w.n,lambda,w); - # Now conjugate with c: - tf := cfi*tf*cf; - # Now cleanup in column n above row n, the entries there - # are lambda times the stuff in column i of ci: - for j in [1..w.n-1] do - tf := DoRowOp_n(tf,j,w.n,-ci[j,i]*lambda,w); - od; - Add(trans,tf); - od; - od; - - # Now put together the clean ones by our knowledge of c^-1: - transd := []; - for i in [1..Length(pivots2)] do - for lambda in w.canb do - tf := w.One; - vals := BlownUpVector(w.can,cii[i]*lambda); - for j in [1..w.ext * newdim] do - pow := IntFFE(vals[j]); - if not IsZero(pow) then - if IsOne(pow) then - tf := tf * trans[j]; - else - tf := tf * trans[j]^pow; - fi; - fi; - od; - Add(transd,tf); - od; - od; - Unbind(trans); - - # Now to the "horizontal" transvections, first create them as SLPs: - transr := []; - for i in pivots do - # This does u_i : v_i -> v_i + v_n - tf := w.One; - tf := DoColOp_n(tf,w.n,i,One(w.f),w); - # Now conjugate with c: - tf := cfi*tf*cf; - # Now cleanup in rows above row n: - for j in [1..w.n-1] do - tf := DoRowOp_n(tf,j,w.n,-ci[j,w.n],w); - od; - # Now cleanup in rows below row n: - for j in [1..newdim] do - coeffs := IntVecFFE(Coefficients(w.can,-ci[w.n+j,w.n])); - for k in [1..w.ext] do - if not IsZero(coeffs[k]) then - if IsOne(coeffs[k]) then - tf := transd[(j-1)*w.ext + k] * tf; - else - tf := transd[(j-1)*w.ext + k]^coeffs[k] * tf; - fi; - fi; - od; - od; - # Now cleanup column n above row n: - for j in [1..w.n-1] do - tf := DoColOp_n(tf,j,w.n,ci[j,w.n],w); - od; - # Now cleanup row n left of column n: - for j in [1..w.n-1] do - tf := DoRowOp_n(tf,w.n,j,-c[i,j],w); - od; - # Now cleanup column n below row n: - for j in [1..newdim] do - coeffs := IntVecFFE(Coefficients(w.can,ci[w.n+j,w.n])); - for k in [1..w.ext] do - if not IsZero(coeffs[k]) then - if IsOne(coeffs[k]) then - tf := tf * transd[(j-1)*w.ext + k]; - else - tf := tf * transd[(j-1)*w.ext + k]^coeffs[k]; - fi; - fi; - od; - od; - Add(transr,tf); - od; - - # From here on we distinguish three cases: - # * w.n = 2 - # * we finish off the constructive recognition - # * we have to do another step as the next thing - if w.n = 2 then - w.slnstdf[2*w.ext+2] := transd[1]*transr[1]^-1*transd[1]; - w.slnstdf[2*w.ext+1] := w.transh[1]*w.transv[1]^-1*w.transh[1] - *w.slnstdf[2*w.ext+2]; - Unbind(w.transh); - Unbind(w.transv); - w.n := 3; - return w; - fi; - # We can finish off: - if aimdim = w.d then - # In this case we just finish off and do not bother with - # the transvections, we will only need the standard gens: - # Now put together the (newdim+1)-cycle: - # n+newdim -> n+newdim-1 -> ... -> n+1 -> n -> n+newdim - flag := false; - s := w.One; - for i in [1..newdim] do - if flag then - # Make [[0,-1],[1,0]] in coordinates w.n and w.n+i: - tf:=transd[(i-1)*w.ext+1]*transr[i]^-1*transd[(i-1)*w.ext+1]; - else - # Make [[0,1],[-1,0]] in coordinates w.n and w.n+i: - tf:=transd[(i-1)*w.ext+1]^-1*transr[i]*transd[(i-1)*w.ext+1]^-1; - fi; - s := s * tf; - flag := not flag; - od; - - # Finally put together the new 2n-1-cycle and 2n-2-cycle: - s := s^-1; - w.slnstdf[2*w.ext+1] := w.slnstdf[2*w.ext+1] * s; - w.slnstdf[2*w.ext+2] := w.slnstdf[2*w.ext+2] * s; - Unbind(w.transv); - Unbind(w.transh); - w.n := aimdim; - return w; - fi; - - # Otherwise we do want to go on as the next thing, so we want to - # keep our transvections. This is easily done if we change the - # basis one more time. Note that we know that n is odd here! - - # Put together the n-cycle: - # 2n-1 -> 2n-2 -> ... -> n+1 -> n -> 2n-1 - flag := false; - s := w.One; - for i in [w.n-1,w.n-2..1] do - if flag then - # Make [[0,-1],[1,0]] in coordinates w.n and w.n+i: - tf := transd[(i-1)*w.ext+1]*transr[i]^-1*transd[(i-1)*w.ext+1]; - else - # Make [[0,1],[-1,0]] in coordinates w.n and w.n+i: - tf := transd[(i-1)*w.ext+1]^-1*transr[i]*transd[(i-1)*w.ext+1]^-1; - fi; - s := s * tf; - flag := not flag; - od; - - # Finally put together the new 2n-1-cycle and 2n-2-cycle: - w.slnstdf[2*w.ext+1] := s * w.slnstdf[2*w.ext+1]; - w.slnstdf[2*w.ext+2] := s * w.slnstdf[2*w.ext+2]; - - list := Concatenation([1..w.n-1],[w.n+1..2*w.n-1],[w.n],[2*w.n..w.d]); - perm := PermList(list); - mat := PermutationMat(perm^-1,w.d,w.f); - ConvertToMatrixRep(mat,w.f); - w.bas := w.bas{list}; - ConvertToMatrixRep(w.bas,w.f); - w.basi := w.basi*mat; - - # Now add the new transvections: - for i in [1..w.n-1] do - w.transh[w.ext*(w.n-1)+w.ext*(i-1)+1] := transr[i]; - od; - Append(w.transv,transd); - w.n := 2*w.n-1; - return w; -end; +################################################################################# +################################################################################# +################################################################################# # RECOG.MakeSLSituation := function(p,e,n,d) # local a,q,r; @@ -3184,7 +2270,7 @@ end; BindRecogMethod(FindHomMethodsProjective, "ClassicalNatural", "check whether it is a classical group in its natural representation", function(ri, g) - local changed,classical,d,det,ext,f,gcd,gens,gg,gm,i,p,pr,q,root,std,stdg,z; + local changed,classical,d,det,ext,f,gcd,gens,gm,i,p,pr,q,root,std,stdg,z; d := ri!.dimension; f := ri!.field; q := Size(f); @@ -3222,12 +2308,10 @@ function(ri, g) fi; od; if changed then - gg := GroupWithGenerators(gens); gm := GroupWithMemory(gens); pr := ProductReplacer(GeneratorsOfGroup(gm),rec(maxdepth := 500)); gm!.pseudorandomfunc := [rec( func := Next, args := [pr] )]; else - gg := g; gm := Group(ri!.gensHmem); gm!.pseudorandomfunc := [rec(func := function(ri,name,bool) return RandomElm(ri,name,bool).el; @@ -3246,10 +2330,10 @@ function(ri, g) Info(InfoRecog,2,"ClassicalNatural: this is PSL_2!"); if IsEvenInt(q) then std := RECOG.RecogniseSL2NaturalEvenChar(gm,f,false); - ri!.comment := "PSL2Even"; + ri!.comment := "_PSL2Even"; else std := RECOG.RecogniseSL2NaturalOddCharUsingBSGS(gm,f); - ri!.comment := "PSL2Odd"; + ri!.comment := "_PSL2Odd"; fi; Setslptonice(ri,SLPOfElms(std.all)); ri!.nicebas := std.bas; @@ -3270,12 +2354,12 @@ function(ri, g) # FIXME: Note d=3 currently has a problem in the SL2-finder. Info(InfoRecog,2,"Classical natural: SL(",d,",",q,"): small ", "case, handing over to Schreier-Sims."); - ri!.comment := Concatenation("SL(",String(d),",",String(q),")", + ri!.comment := Concatenation("_SL(",String(d),",",String(q),")", "_StabilizerChain"); return FindHomMethodsProjective.StabilizerChainProj(ri,g); fi; Info(InfoRecog,2,"ClassicalNatural: this is PSL_n!"); - std := RECOG.FindStdGens_SL(gm,f); + std := RECOG.FindStdGens_SL(gm); Setslptonice(ri,std.slpstd); ri!.nicebas := std.bas; ri!.nicebasi := std.basi; @@ -3285,7 +2369,7 @@ function(ri, g) x->std.basi*x*std.bas)); ri!.fakegens := RECOG.InitSLfake(f,d); ri!.fakegens.count := 0; - ri!.comment := "PSLd"; + ri!.comment := "_PSLd"; ri!.gcd := gcd; SetFilterObj(ri,IsLeaf); SetSize(ri,Product([0..d-1],i->(q^d-q^i))/((q-1)*gcd.gcd)); diff --git a/gap/projective/constructive_recognition/utils/utils.gi b/gap/projective/constructive_recognition/utils/utils.gi new file mode 100644 index 00000000..016f402a --- /dev/null +++ b/gap/projective/constructive_recognition/utils/utils.gi @@ -0,0 +1,1986 @@ +############################################################################# +## +## This file is part of recog, a package for the GAP computer algebra system +## which provides a collection of methods for the constructive recognition +## of groups. +## +## This files's authors include Daniel Rademacher, Max Neunhöffer, Ákos Seress. +## +## Copyright of recog belongs to its developers whose names are too numerous +## to list here. Please refer to the COPYRIGHT file for details. +## +## SPDX-License-Identifier: GPL-3.0-or-later +## +############################################################################# + + + +############################################################################# +############################################################################# +######## General utils ###################################################### +############################################################################# +############################################################################# + + + +InstallMethod( CharacteristicPolynomial, "for a memory element matrix", + [ IsMatrix and IsObjWithMemory ], + function(m) + return CharacteristicPolynomial(m!.el); + end ); + + + +InstallOtherMethod( \-, "for two memory elements", + [ IsMatrix and IsObjWithMemory, IsMatrix and IsObjWithMemory ], + function(m,n) + return m!.el - n!.el; + end ); + + + +InstallMethod( Eigenspaces, "for a field and a memory element matrix", + [ IsField, IsMatrix and IsObjWithMemory ], + function( f, m ) + return Eigenspaces(f,m!.el); + end ); + + +# compute the eigenspace of `mat` for the given eigenvalue lambda` +RECOG.EigenspaceMat := function(mat, lambda) + local i; + mat := MutableCopyMat( mat ); + # since mat is a copy, we can efficiently "subtract an identity matrix" + for i in [1..NrRows(mat)] do + mat[i,i] := mat[i,i] - lambda; + od; + # since mat is a copy we can use NullspaceMatDestructive instead of NullspaceMat + return NullspaceMatDestructive(mat); +end; + +# compute fixed space of mat, i.e. eigenspace for eigenvalue 1 +RECOG.FixspaceMat := function(mat) + return RECOG.EigenspaceMat(mat, 1); +end; + + + +############################################################################# +############################################################################# +######## CheckStingrayGroup ################################################# +############################################################################# +############################################################################# + + +RECOG.CheckNewStingrayGroup := function(g1,base1,g2,base2,q) +local baseSum, module; + + baseSum := Concatenation(base1,base2); + baseSum := TriangulizedMat(baseSum); + if IsZero(Last(baseSum)) then + return false; + fi; + + g1 := TransposedMat(StripMemory(g1)); + g2 := TransposedMat(StripMemory(g2)); + module := GModuleByMats( [g1,g2], GF(q) ); + module := MTX.InducedActionSubmoduleNB( module, baseSum ); + return MTX.IsIrreducible(module); +end; + + + +############################################################################# +############################################################################# +######## ConstructSmallSub ################################################## +############################################################################# +############################################################################# + + + +RECOG.ConstructSmallSub := function(r1, r2, product, newbasis, detectFun) + local gens, pseudoorderlist, Hsub, productEle, ele, ele2, H, cord1, cord2; + + gens := []; + pseudoorderlist := []; + Hsub := []; + repeat + productEle := PseudoRandom(product); + Add(Hsub, productEle); + ele := (productEle)^(newbasis^(-1)); + ele2 := ele{r2}{r2}; + ele := ele{r1}{r1}; + Add(pseudoorderlist, RECOG.EstimateOrder(ele2)[1]); + Add(gens,ele); + until Size(gens) = 2; + H := GroupByGenerators(gens); + if detectFun(H) = true then + cord1 := Order(gens[1]); + cord2 := Order(gens[2]); + if (Gcd(cord1,pseudoorderlist[1]) <> pseudoorderlist[1]) and (Gcd(cord2,pseudoorderlist[2]) <> pseudoorderlist[2]) then + gens[1] := gens[1]^pseudoorderlist[1]; + gens[2] := gens[2]^pseudoorderlist[2]; + H := GroupByGenerators(gens); + if detectFun(H) = true then + Hsub[1] := Hsub[1]^pseudoorderlist[1]; + Hsub[2] := Hsub[2]^pseudoorderlist[2]; + return [Hsub,H,newbasis]; + fi; + fi; + fi; + return fail; +end; + + + +############################################################################# +############################################################################# +######## constructppdTwoStingray ############################################ +############################################################################# +############################################################################# + + + +RECOG.constructppdTwoStingray := function(g,dim,q,type,form) + local out, list, out2, currentdim, aimdim, godown; + + if type = "SL" then + aimdim:=-1; + elif type = "O" then + aimdim:=8; + elif type = "Sp" then + aimdim:=8; + elif type = "SU" then + if IsEvenInt(q) then + aimdim := 10; + else + aimdim := 6; + fi; + else + Error("unsupported type ", type); + fi; + + Info(InfoRecog,2,"Current Dimension: ", dim, " for type ", type); + Info(InfoRecog,2,"\n"); + + list:=[g,dim,q,fail,form]; + currentdim := dim; + repeat + out:=RECOG.godownStingray(list,type); + if out=fail or IsOne(out[1]^2) then + Info(InfoRecog,2,"Restart. \n"); + Info(InfoRecog,2,"Current Dimension: "); + Info(InfoRecog,2,dim); + Info(InfoRecog,2,"\n"); + list:=[g,dim,q,fail,form]; + out:=fail; + else + if type = "SL" and out[2] = 2 then + return out[1]; + fi; + Assert(0, out[1] >= 2); + repeat + out2:=RECOG.godownStingray(list,type); + if out2=fail or out2[1]*out2[1]=One(out2[1]) then + if InfoLevel(InfoRecog) >= 3 then Print("B\c"); fi; + list:=[g,dim,q,fail,form]; + out2:=fail; + fi; + until out2<>fail and out2[2] >= 2; + if type = "SL" and out2[2] = 2 then + return out2[1]; + fi; + if RECOG.CheckNewStingrayGroup(out[1],out[3],out2[1],out2[3],q) then + list:=[Group(out[1],out2[1]),out[2]+out2[2],q,fail,form]; + currentdim := list[2]; + + Info(InfoRecog,2,"Debugg Info:\n"); + Info(InfoRecog,2,"Dimension FirstElement: "); + Info(InfoRecog,2,out[2]); + Info(InfoRecog,2,"\n"); + Info(InfoRecog,2,"Dimension SecondElement: "); + Info(InfoRecog,2,out2[2]); + Info(InfoRecog,2,"\n"); + Info(InfoRecog,2,"End Debugg Info. \n"); + + Info(InfoRecog,2,"New Dimension: "); + Info(InfoRecog,2,out[2]+out2[2]); + Info(InfoRecog,2,"\n"); + else + if InfoLevel(InfoRecog) >= 3 then Print("B\c"); fi; + Info(InfoRecog,2,"Restart. \n"); + Info(InfoRecog,2,"Current Dimension: "); + Info(InfoRecog,2,dim); + Info(InfoRecog,2,"\n"); + list:=[g,dim,q,fail,form]; + out:=fail; + fi; + fi; + until currentdim=aimdim; + + return list[1]; + +end; + + + +############################################################################# +############################################################################# +######## godownStingray ##################################################### +############################################################################# +############################################################################# + + + +# finds first element of a list that is relative prime to all others +# input: list=[Sp(d,q), d, q, Sp(n,q)] acting as a subgroup of some big Sp(n,q) +# output: list=[rr, dd] for a ppd(2*dd;q)-element rr +RECOG.godownStingray := function(list,type) +local d, firstSL, firstSU, q, p, g, i, r, pol, factors, degrees, newdim, power, rr, ss, max, +newgroup, colldegrees, exp, count, check, ocount, beta, NiList, Maxi, qFactors, +irrfact, invbase, form, CheckOtherFactors, CheckFactors, fld, restricted, b, j; + + CheckOtherFactors := function(i, deg, fact) + local j; + for j in [1..Length(deg)] do + if not(j = i) then + if RECOG.CheckPolynomialForSelfConjugate(fact[j]) then + if (deg[j] mod deg[i] = 0) then + return false; + fi; + else + if (deg[j] mod Int(deg[i]/2) = 0) then + return false; + fi; + fi; + fi; + od; + return true; + end; + + CheckFactors := function(deg, fact) + local i; + for i in [1..Length(deg)] do + if ((deg[i] mod 2) = 0) and RECOG.CheckPolynomialForSelfConjugate(fact[i]) and CheckOtherFactors(i,deg,fact) then + return i; + fi; + od; + return fail; + end; + + firstSU := function(list) + local i, j, goodElement; + for i in [1..Length(list)] do + if list[i]>1 and (list[i] mod 2 = 1) then + if Gcd(list[i],Product(list)/list[i]) < list[i] then + return i; + else + goodElement := true; + for j in [1..Length(list)] do + if not(j = i) and Gcd(list[i],list[j]) = list[i] then + goodElement := false; + break; + fi; + od; + if goodElement then + return i; + fi; + fi; + fi; + od; + return fail; + end; + + firstSL := function(list) + local i, j, goodElement; + for i in [1..Length(list)] do + if list[i]>1 then + if Gcd(list[i],Product(list)/list[i]) < list[i] then + return i; + else + goodElement := true; + for j in [1..Length(list)] do + if not(j = i) and Gcd(list[i],list[j]) = list[i] then + goodElement := false; + break; + fi; + od; + if goodElement then + return i; + fi; + fi; + fi; + od; + return fail; + end; + + g:=list[1]; + d:=list[2]; + q:=list[3]; + qFactors:=Factors(q); + p := qFactors[1]; + form := list[5]; + fld := GF(q); + + if type = "SL" then + max := Maximum([Log2Int(d),3]); + elif type = "Sp" then + max := Maximum([2*Log2Int(d),3]); + elif type = "SU" then + max := Maximum([2*Log2Int(d),3]); + elif type = "O" then + max := Maximum([2*Log2Int(d),3]); + else + Error("type not supported"); + fi; + + # Overall count. Replace by formula and unequality + ocount := 0; + while ocount < 100 do + + Info(InfoRecog,2,"Dimension: ",d); + #find an element with irreducible action of relative prime dimension to + #all other invariant subspaces + #count is just safety, if things go very bad + count:=0; + + repeat + count:=count+1; + r:=PseudoRandom(g); + pol:=CharacteristicPolynomial(r); + factors:=Factors(pol); + degrees:=List(factors,Degree); + if type = "SL" then + newdim:= firstSL(degrees); + elif type = "SU" then + newdim:= firstSU(degrees); + elif type = "O" or type = "Sp" then + newdim := CheckFactors(degrees, factors); + else + Error("type not supported"); + fi; + until (count>100) or (newdim <> fail and (degrees[newdim] < max)); + # Be careful if Log2Int(d) = 2! In this case we search for stingray elements with k < 2. Hence use newdim < Maximum([Log2Int(d),3]) + + if count>100 then + return fail; + fi; + + # Split result from first: + irrfact := factors[newdim]; + newdim := degrees[newdim]; + + if newdim = 2 and type = "SL" then + check := true; + else + # Check whether the stingray candidate is a ppd-stingray element + check := RECOG.IsPpdStingrayElement(p,Length(qFactors),newdim,irrfact); + fi; + + if check then + + # raise r to a power so that acting trivially outside one invariant irreducible subspace + NiList := Collected(degrees); + NiList := Filtered(NiList,x->not(x[1] = newdim)); + colldegrees := List(NiList,x->x[1]); + NiList := List(NiList,x->x[2]); + Maxi := Maximum(NiList); + beta := LogInt(Maxi,p); + if not(p^beta = Maxi) then + beta := beta + 1; + fi; + + # power further to cancel q-part of element order + power := Lcm(List(colldegrees, x->q^x-1))*p^beta; + rr:=r^power; + + invbase := NullspaceMat(TransposedMat(RECOG.FixspaceMat(StripMemory(rr)))); + + if newdim = 2 and type = "SL" then + if Size(invbase) = 2 then + return [rr,newdim,invbase]; + fi; + else + + #if (type = "SL") or (IsEvenInt(q) and type = "SU") then + # return [rr,newdim,invbase]; + #fi; + + #b := Basis(VectorSpace(fld,invbase),invbase); + #restricted := IdentityMat(newdim,fld); + #for i in [1..newdim] do + # for j in [1..newdim] do + # restricted[i,j] := b[i]*form*b[j]; + # od; + #od; + + #if IsEmpty(NullspaceMat(restricted)) then + return [rr,newdim,invbase]; + #else + # Error("here"); + #fi; + fi; + fi; + + ocount := ocount + 1; + od; + + return fail; + +end; + + + +############################################################################# +############################################################################# +######## PPD Check ########################################################## +############################################################################# +############################################################################# + + + +RECOG.CheckPPDdqe := function(m,d,q,e) +local factors, i, ele, good, ord; + + factors := Factors(q^e-1); + ord := Order(m); + for ele in factors do + good := true; + for i in [1..(e-1)] do + if ((q^i-1) mod ele) = 0 then + good := false; + break; + fi; + od; + if good then + if (ord mod ele) = 0 then + return true; + fi; + fi; + od; + + return false; +end; + + + +############################################################################# +############################################################################# +######## Coefficients Primitive Element ##################################### +############################################################################# +############################################################################# + + + +# The following function has been written by Thomas Breuer. +# It expresses an element alpha in a field fld as +# a linear combination of a Primitive Element. + +# Input: fld: A field, +# alpha : An element of fld + +# Output: Coefficients: A vector c sth for omega primitive element +# alpha = sum c[i] omega^(i-1) + +RECOG.CoefficientsPrimitiveElement := function ( fld, alpha ) + + if Size( fld ) <= MAXSIZE_GF_INTERNAL then + + return Coefficients( CanonicalBasis( fld ), alpha ); + else + + alpha := FFECONWAY.WriteOverLargerField( alpha, DegreeOverPrimeField( fld ) ); + + if IsCoeffsModConwayPolRep( alpha ) then + return alpha![1]; + elif IsModulusRep(alpha) then + return [alpha]; + else + Error( "this should not happen" ); + fi; + fi; + +end; + + + +############################################################################# +############################################################################# +######## Check PPD-Property and tests ####################################### +############################################################################# +############################################################################# + + + +## This function takes as input: +## +## field +## a characteristic polynomial +## degree of +##

a prime power +## an integer +## an irreducible factor of and of degree a + +RECOG.IsPpdStingrayElement := function( p, f, k, irrfact ) + local c, e, R, pm, g, islarge, F; + + F := GF(p^f); + c := irrfact; + R := PolynomialRing(F); + + e := k; + ## find the noppds and ppds parts + pm := PrimitivePrimeDivisors( f*e, p ); + ## pm contains two fields, noppds and ppds. + ## ppds is the product of all ppds of p^(ae)-1 + ## and noppds is p^(ae)-1/ppds. + + ## get rid of the non-ppd part + ## g will be x^noppds in F[x]/ + g := PowerMod( Indeterminate(F), pm.noppds, c ); + + ## now we know that is a ppd-element + + ## if g is one there is no ppd involved + if IsOne(g) then + return false; + else + return true; + fi; + + #if IsOne(g) then + # ## (e+1) is the only ppd dividing || and only once + # islarge := false; + # return [ e, islarge ]; + #else + # ## Either another ppd also divides || and this one is large or + # ## (e+1)^2 divides || and hence still large + # islarge := true; + # return [ e, islarge ]; + #fi; + + +end; + + + +RECOG.ComparePPDFunction := function(m,d,q,e,irrfact) + local f,p,factors; + + factors := Factors(q); + p := factors[1]; + f := Length(factors); + + if not(RECOG.CheckPPDdqe(m,d,q,e) = RECOG.IsPpdStingrayElement(p,f,e,irrfact[1])) then + Error("PPD error"); + fi; + +end; + + + +############################################################################# +############################################################################# +######## Linear action representation ####################################### +############################################################################# +############################################################################# + + + +RECOG.LinearAction := function(bas,field,el) + local mat,vecs; + if IsGroup(el) then + return Group(List(GeneratorsOfGroup(el), + x->RECOG.LinearAction(bas,field,x))); + fi; + if IsBasis(bas) then + vecs := BasisVectors(bas); + else + vecs := bas; + bas := Basis(VectorSpace(field,bas),bas); + fi; + mat := List(vecs,v->Coefficients(bas,v*el)); + ConvertToMatrixRep(mat,field); + return mat; +end; + + + +RECOG.LinearActionRepresentation := function(G) +local OldGens, newGens, i, base, fld, d, EleBase, fixbase, B, action, ele, V; + + OldGens := ShallowCopy(GeneratorsOfGroup(G)); + for i in [1..Length(OldGens)] do + if IsObjWithMemory(OldGens[i]) then + OldGens[i] := StripMemory(OldGens[i]); + fi; + od; + + fld := FieldOfMatrixList(OldGens); + d := Size(OldGens[1]); + base := []; + for i in [1..Length(OldGens)] do + ele := OldGens[i]; + fixbase := RECOG.FixspaceMat(TransposedMat(ele)); + if fixbase = [] then + return fail; + fi; + EleBase := NullspaceMat(TransposedMat(fixbase)); + Append(base,EleBase); + od; + + V := VectorSpace(fld,base); + B := Basis(V); + base := BasisVectors(B); + newGens := []; + for i in [1..Length(OldGens)] do + ele := OldGens[i]; + action := List(base,v->Coefficients(B,v*ele)); + Add(newGens,action); + od; + + return GroupByGenerators(newGens); +end; + + + +RECOG.LinearActionRepresentationSmallerMatrices := function(G) +local OldGens, newGens, i, base, fld, d, EleBase, fixbase, B, action, ele, V, baseCopy; + + OldGens := ShallowCopy(GeneratorsOfGroup(G)); + for i in [1..Length(OldGens)] do + if IsObjWithMemory(OldGens[i]) then + OldGens[i] := StripMemory(OldGens[i]); + fi; + od; + # Einfach StripMemory OldGens := StripMemory(GeneratorsOfGroup(G)) + + fld := FieldOfMatrixList(OldGens); + d := NumberRows(OldGens[1]); + base := []; + for i in [1..Length(OldGens)] do + ele := OldGens[i]; + fixbase := RECOG.FixspaceMat(TransposedMat(ele)); + EleBase := NullspaceMat(TransposedMat(fixbase)); + Append(base,EleBase); + od; + baseCopy := base; + + V := VectorSpace(fld,base); + B := Basis(V); + base := BasisVectors(B); + newGens := []; + for i in [1..Length(OldGens)] do + ele := OldGens[i]; + action := List(base,v->Coefficients(B,v*ele)); + + # DR: Change here so that we still operate from the same side + Add(newGens,action); + od; + + return [GroupByGenerators(newGens),Size(B),baseCopy]; +end; + + + +RECOG.LinearActionRepresentationSmallerMatricesVersion2 := function(G) +local OldGens, newGens, i, base, fld, d, EleBase, fixbase, B, action, ele, V; + + OldGens := ShallowCopy(GeneratorsOfGroup(G)); + for i in [1..Length(OldGens)] do + if IsObjWithMemory(OldGens[i]) then + OldGens[i] := StripMemory(OldGens[i]); + fi; + od; + # Einfach StripMemory OldGens := StripMemory(GeneratorsOfGroup(G)) + + fld := FieldOfMatrixList(OldGens); + d := NumberRows(OldGens[1]); + base := []; + for i in [1..Length(OldGens)] do + ele := OldGens[i]; + fixbase := RECOG.FixspaceMat(TransposedMat(ele)); + EleBase := NullspaceMat(TransposedMat(fixbase)); + Append(base,EleBase); + od; + + V := VectorSpace(fld,base); + B := Basis(V); + base := BasisVectors(B); + newGens := []; + for i in [1..Length(OldGens)] do + ele := OldGens[i]; + action := List(base,v->Coefficients(B,v*ele)); + + # DR: Change here so that we still operate from the same side + Add(newGens,action); + od; + + return [GroupByGenerators(newGens),Size(B),BasisVectors(B)]; +end; + + + +############################################################################# +############################################################################# +######## Self-conjugate polynomial check #################################### +############################################################################# +############################################################################# + + + +RECOG.CheckPolynomialForSelfConjugate := function (f) +local ind, coeff, aZero, i, fld, deg, pol; + + fld := FieldOfPolynomial(f); + ind := IndeterminateOfLaurentPolynomial(f); + coeff := CoefficientsOfLaurentPolynomial(f)[1]; + deg := Length(coeff); + aZero := coeff[1]; + + pol := ind^0 * Zero(fld); + for i in [1..deg] do + pol := pol + ind^(deg-i)*coeff[i]; + od; + + pol := aZero * pol; + + if pol = f then + return true; + else + return false; + fi; +end; + + + +############################################################################# +############################################################################# +######## Extract and rescale block matrices ################################# +############################################################################# +############################################################################# + + + +RECOG.ComputeBlockBaseChangeMatrix := function(list,d,q) +local fixbase, elebase, basis, matrix, fix, moved, currentmove, currentfix, k, newbase, OldGens, i; + + OldGens := ShallowCopy(list); + for i in [1..Length(OldGens)] do + if IsObjWithMemory(OldGens[i]) then + OldGens[i] := StripMemory(OldGens[i]); + fi; + od; + list := OldGens; + + fix := []; + moved := []; + + for matrix in list do; + fixbase := RECOG.FixspaceMat(TransposedMat(matrix)); + elebase := NullspaceMat(TransposedMat(fixbase)); + Add(moved, elebase); + + fixbase := RECOG.FixspaceMat(matrix); + Add(fix,fixbase); + od; + + if Size(moved) = 1 then + newbase := MutableCopyMat(moved[1]); + Append(newbase,fix[1]); + return newbase; + else + currentmove := MutableCopyMat(moved[1]); + currentfix := MutableCopyMat(fix[1]); + k := 1; + while k < Size(moved) do + currentmove := SumIntersectionMat(currentmove,moved[k+1])[1]; + currentfix := SumIntersectionMat(currentfix,fix[k+1])[2]; + k := k + 1; + od; + Append(currentmove,currentfix); + return currentmove; + fi; + +end; + + + +RECOG.ExtractSmallerGroup := function(list,basechange,size) +local gens, ele, block, OldGens, i; + + OldGens := ShallowCopy(list); + for i in [1..Length(OldGens)] do + if IsObjWithMemory(OldGens[i]) then + OldGens[i] := StripMemory(OldGens[i]); + fi; + od; + list := OldGens; + + gens := []; + for ele in list do + block := (ele^(basechange^(-1))); + block := block{[1..size]}{[1..size]}; + Add(gens,block); + od; + + return [GroupByGenerators(gens),gens]; +end; + + + +RECOG.LiftGroup := function(list,size,q,d) +local gens, ele, block, OldGens, i; + + OldGens := ShallowCopy(list); + for i in [1..Length(OldGens)] do + if IsObjWithMemory(OldGens[i]) then + OldGens[i] := StripMemory(OldGens[i]); + fi; + od; + list := OldGens; + + gens := []; + for ele in list do + block := IdentityMat(d,GF(q)); + block{[1..size]}{[1..size]} := ele; + Add(gens,block); + od; + + return [GroupByGenerators(gens),gens]; +end; + + + +############################################################################# +############################################################################# +######## Membership test in groups preserving a form ######################## +############################################################################# +############################################################################# + + +# given a matrix `mat`, test if it is contained in G, which must be Omega(e,n,fld) +# +# TODO: add unit tests +# +# e:=+1; d:=8; q:=8; +# G:=Omega(e,d,q); +# H:=SO(e,d,q,InvariantQuadraticForm(G).matrix); +# ForAll(GeneratorsOfGroup(G), g -> g in H); +# ForAll(GeneratorsOfGroup(G), g -> IsInOmega(G, g)); +# ForAll(GeneratorsOfGroup(H), g -> IsInOmega(G, g)); +# +RECOG.IsInOmega:=function(G,mat) + local d, Q, form, fld; + d := DimensionOfMatrixGroup(G); + fld := FieldOfMatrixGroup(G); + Assert(0, NrRows(mat) = d); + + # first verify the quadratic form is preserved + Q := InvariantQuadraticForm(G).matrix; + if not RespectsQuadraticForm(Q, mat) then + return false; + fi; + + if Characteristic(fld) <> 2 then + form := InvariantBilinearForm(G).matrix; + return IsOne(SpinorNorm(form, fld, mat)); + elif IsOddInt(d) then + # Omega(0,2n+1,2^k) = SO(0,2n+1,2^k) = GO(0,2n+1,2^k) + return true; + else + # the following is based on Lemma 3.5(2) in Holt, Roney-Dougal: + # "Constructing maximal subgroups of orthogonal groups" + return IsEvenInt(RankMat(mat + One(G))); + fi; +end; + + + +############################################################################# +############################################################################# +# Code from ClassicalMaximal to check BilinearForm ########################## +# https://github.com/gap-packages/ClassicalMaximals/blob/main/gap/Forms.gi ## +############################################################################# +############################################################################# + + + +RECOG.MatrixByEntries := function(field, nrRows, nrCols, entries) + local i, m, o; + o := One(field); + m := NullMat(nrRows, nrCols, field); + for i in entries do + m[i[1],i[2]] := i[3]*o; + od; + return ImmutableMatrix(field, m); +end; + + + +RECOG.AntidiagonalMat := function(entries, field) + local d, m, i; + if IsInt(entries) then + d := entries; + entries := ListWithIdenticalEntries(d, One(field)); + else + d := Length(entries); + fi; + m := NullMat(d, d, field); + for i in [1..d] do + m[i, d - i + 1] := entries[i]; + od; + return ImmutableMatrix(field, m); +end; + + + +# Solving the congruence a ^ 2 + b ^ 2 = c in F_q by trial and error. +# +# A solution always exists by a simple counting argument using the pigeonhole +# principle and the fact that there are (q + 1) / 2 > q / 2 squares in F_q (for +# q odd; the case q even is trivial). The trial and error approach taken here +# should in principle almost always terminate quickly: Assuming that - 1 - a ^ 2 +# is evenly distributed in GF(q), the chance to hit a quadratic residue is about +# 1 / 2 in each trial. +RECOG.SolveQuadraticCongruence := function(c, q) + local F, a, b; + F := GF(q); + for a in F do + b := RootFFE(F, (c - a ^ 2) * Z(q) ^ 0, 2); + if not b = fail then + break; + fi; + od; + return rec(a := a, b := b); +end; + + + +RECOG.ApplyFunctionToEntries := function(M, func) + local numberRows, numberColumns, i, j, result; + if not IsMatrix(M) or Length(M) = 0 then + ErrorNoReturn(" must be a matrix but = ", M); + fi; + + numberRows := NrRows(M); + numberColumns := NrCols(M); + result := NullMat(numberRows, numberColumns, DefaultFieldOfMatrix(M)); + for i in [1..numberRows] do + for j in [1..numberColumns] do + result[i, j] := func(M[i, j]); + od; + od; + + return result; +end; + + + +RECOG.HermitianConjugate := function(M, q) + return TransposedMat(RECOG.ApplyFunctionToEntries(M, x -> x ^ q)); +end; + + + +# If type = "S" then find a beta in GF(q ^ 2) with beta + beta ^ q = alpha. +# If type = "P" then find a beta in GF(q ^ 2) with gamma * gamma ^ q = alpha. +# In both cases, alpha is an element of GF(q). +# Construction as in Lemma 2.2 of [HR05] +RECOG.SolveFrobeniusEquation := function(type, alpha, q) + local F, R, S, x, delta, polynomial, result; + + F := GF(q); + if not alpha in F then + ErrorNoReturn(" must be an element of GF() but = ", + alpha, " and = ", q); + fi; + if not type in ["S", "P"] then + ErrorNoReturn(" must be one of 'S' or 'P' but = ", type); + fi; + # We have to make an exception for this case since the construction below + # does not work here: x ^ 2 + delta is never irreducible over GF(q) since + # all elements of GF(q) are squares for q even. + if type = "S" and alpha = 0 and IsEvenInt(q) then + return Z(q) ^ 0; + fi; + + R := PolynomialRing(F, ["x"]); + S := PolynomialRing(GF(q ^ 2), ["x"]); + x := Indeterminate(F, "x"); + + # A quick argument using the quadratic formula for q odd or some + # algebraic manipulations and the non-surjectivity of the Artin-Schreier + # map x -> x ^ 2 + x for q odd and alpha <> 0 shows that the construction + # below always works. + if type = "S" then + for delta in F do + polynomial := x ^ 2 - alpha * x + delta; + if IsIrreducibleRingElement(R, polynomial) then + result := -CoefficientsOfUnivariatePolynomial(Factors(S, polynomial)[1])[1]; + return result; + fi; + od; + # A similar argument to the one used for type "S" works here as well. Note, + # however, that the argument for q odd via the quadratic formula now + # additionally needs that the squares in GF(q) do not form an arithmetic + # progression (which is "closed", i.e. not only a_i+1 = a_i + d, but also + # a_n + d = a_1), which can be proved in the following way: If they did, + # then they would be given by -kd, ..., -d, 0, d, 2d, ..., ((q - 1) / 2 - k) * d + # for some 0 <= k <= (q - 1) / 2; since they form a multiplicative + # subgroup, we can divide by -d or d, respectively, and see that + # -+k, ..., -+1, 0, +-1, +-2, ..., +-((q - 1) / 2 - k) are also all the + # squares in GF(q). Most notably they all are in GF(p) and thus there are + # at most p squares in GF(q), which is < (q + 1) / 2 if q >= p ^ 2 - a + # contradiction. Now we can restrict ourselves to q = p and reach a + # contradiction for p >= 7 (we leave out the details); this leaves p = 3 + # and p = 5, which can easily be checked manually. + elif type = "P" then + for delta in F do + polynomial := x ^ 2 + delta * x + alpha; + if IsIrreducibleRingElement(R, polynomial) then + result := -CoefficientsOfUnivariatePolynomial(Factors(S, polynomial)[1])[1]; + return result; + fi; + od; + fi; +end; + + + +# An n x n - matrix of zeroes with a 1 in position (row, column) +RECOG.SquareSingleEntryMatrix := function(field, n, row, column) + return RECOG.MatrixByEntries(field, n, n, [[row, column, 1]]); +end; + + + +# Compute Ceil(m / n) for two integers m, n +RECOG.QuoCeil := function(m, n) + if m mod n = 0 then + return QuoInt(m, n); + else + return QuoInt(m, n) + 1; + fi; +end; + + + +# Compute the size of Sp(n, q) according to Theorem 1.6.22 in [BHR13] +RECOG.SizeSp := function(n, q) + local m, result, powerOfq, i; + if IsOddInt(n) then + ErrorNoReturn("Dimension must be even but ", n, " was given."); + fi; + m := QuoInt(n, 2); + result := q ^ (m * m); + powerOfq := 1; + for i in [1..m] do + powerOfq := powerOfq * q * q; + result := result * (powerOfq - 1); + od; + return result; +end; + + + +# Compute the size of PSp(n, q) according to Table 1.3 in [BHR13], +RECOG.SizePSp := function(n, q) + return QuoInt(RECOG.SizeSp(n, q), Gcd(2, q - 1)); +end; + + + +# Compute the size of SU(n, q) according to Theorem 1.6.22 in [BHR13] +# using the formula for GU(n, q) but starting with i = 2 +# because Table 1.3 gives [GU(n, q):SU(n, q)] = q + 1. +RECOG.SizeSU := function(n, q) + local result, powerOfq, sign, i; + result := q ^ QuoInt(n * (n - 1), 2); + powerOfq := q; + sign := 1; + for i in [2..n] do + powerOfq := powerOfq * q; + sign := -sign; + result := result * (powerOfq + sign); + od; + return result; +end; + + + +# Compute the size of PSU(n, q) according to Table 1.3 in [BHR13] +# Namely, we have | G / Z(G) : S / Z(S) | = | G : S | * |Z(S)| / |Z(G)| so due +# to | G : S | = q + 1, |Z(G)| = q + 1 and | G / Z(G) : S / Z(S) | = (q + 1, n), +# which are given in said table, this gives |Z(S)| = (q + 1, n). +RECOG.SizePSU := function(n, q) + return RECOG.SizeSU(n, q) / Gcd(n, q + 1); +end; + + + +# Compute the size of GU(n, q) according to Table 1.3 in [BHR13] +RECOG.SizeGU := function(n, q) + return (q + 1) * RECOG.SizeSU(n, q); +end; + + + +# Compute the size of GO(epsilon, n, q) according to Theorem 1.6.22 in [BHR13] +RECOG.SizeGO := function(epsilon, n, q) + local m, result, powerOfq, i; + if epsilon = 0 then + + if IsEvenInt(n) then + ErrorNoReturn("for = ", epsilon, " the dimension must be odd but ", n, " was given."); + fi; + + if IsEvenInt(q) then + return RECOG.SizeSp(n - 1, q); + fi; + + m := QuoInt(n - 1, 2); + result := 2 * q ^ (m * m); + + elif epsilon in [-1, 1] then + + if IsOddInt(n) then + ErrorNoReturn("for = ", epsilon, " the dimension must be even but ", n, " was given."); + fi; + + m := QuoInt(n, 2); + result := 2 * q ^ (m * (m - 1)) * (q ^ m - epsilon); + m := m - 1; + else + ErrorNoReturn(" must be in [-1, 0, 1] but ", epsilon, " was given."); + fi; + + powerOfq := 1; + for i in [1..m] do + powerOfq := powerOfq * q * q; + result := result * (powerOfq - 1); + od; + + return result; +end; + + + +# Compute the size of SO(epsilon, n, q) according to Table 1.3 in [BHR13] +RECOG.SizeSO := function(epsilon, n, q) + return QuoInt(RECOG.SizeGO(epsilon, n, q), Gcd(2, q - 1)); +end; + + + +RECOG.ReflectionMatrix := function(n, q, gramMatrix, v) + local F, reflectionMatrix, i, basisVector, reflectBasisVector, beta; + F := GF(q); + reflectionMatrix := NullMat(n, n, F); + beta := BilinearFormByMatrix(gramMatrix); + if IsZero(EvaluateForm(beta, v, v)) then + ErrorNoReturn("The vector must have non-zero norm with respect to", + " the bilinear form given by "); + fi; + for i in [1..n] do + basisVector := List([1..n], j -> Zero(F)); + basisVector[i] := Z(q) ^ 0; + reflectBasisVector := basisVector + - 2 * EvaluateForm(beta, v, basisVector) + / EvaluateForm(beta, v, v) * v; + reflectionMatrix[i]{[1..n]} := reflectBasisVector; + od; + return reflectionMatrix; +end; + + + +# Construct generators for the orthogonal groups with the properties listed in +# Lemma 2.4 of [HR05]. +# Construction as in: C. M. Roney-Dougal. "Conjugacy of Subgroups of the +# General Linear Group." Experimental Mathematics, vol. 13 no. 2, 2004, pp. +# 151-163. Lemma 2.4. +# We take the notation from [HR05]. +RECOG.GeneratorsOfOrthogonalGroup := function(epsilon, n, q) + local F, gramMatrix, generatorsOfSO, vectorOfSquareNorm, D, E, zeta, a, b, + solutionOfQuadraticCongruence; + if IsEvenInt(q) then + ErrorNoReturn("This function was only designed for odd but = ", + n, "and = ", q); + fi; + + F := GF(q); + zeta := PrimitiveElement(F); + if IsOddInt(n) then + gramMatrix := IdentityMat(n, F); + generatorsOfSO := GeneratorsOfGroup(RECOG.ConjugateToSesquilinearForm(SO(epsilon, n, q), + "O", + gramMatrix)); + D := - IdentityMat(n, F); + E := zeta * IdentityMat(n, F); + else + if epsilon = 1 then + gramMatrix := RECOG.AntidiagonalMat(n, F); + generatorsOfSO := GeneratorsOfGroup(RECOG.ConjugateToSesquilinearForm(SO(epsilon, n, q), + "O", + gramMatrix)); + # Our standard bilinear form is given by the Gram matrix + # Antidiag(1, ..., 1). The norm of [1, 0, ..., 0, 2] under this + # bilinear form is 4, i.e. a square. (Recall q odd, so this is not zero!) + vectorOfSquareNorm := zeta ^ 0 * Concatenation([1], + List([1..n - 2], i -> 0), + [2]); + D := RECOG.ReflectionMatrix(n, q, gramMatrix, vectorOfSquareNorm); + E := DiagonalMat(Concatenation(List([1..n / 2], i -> zeta), + List([1..n / 2], i -> zeta ^ 0))); + elif epsilon = -1 then + + # Get a, b in GF(q) with a ^ 2 + b ^ 2 = zeta + solutionOfQuadraticCongruence := RECOG.SolveQuadraticCongruence(zeta, q); + a := solutionOfQuadraticCongruence.a; + b := solutionOfQuadraticCongruence.b; + + if IsOddInt(n * (q - 1) / 4) then + gramMatrix := IdentityMat(n, F); + generatorsOfSO := GeneratorsOfGroup(RECOG.ConjugateToSesquilinearForm(SO(epsilon, n, q), + "O", + gramMatrix)); + # Our standard bilinear form is given by the Gram matrix + # Diag(1, ..., 1). The norm of [1, 0, ..., 0] under this bilinear + # form is 1, i.e. a square. + vectorOfSquareNorm := zeta ^ 0 * Concatenation([1], + List([1..n - 1], i -> 0)); + D := RECOG.ReflectionMatrix(n, q, gramMatrix, vectorOfSquareNorm); + # Block diagonal matrix consisting of n / 2 blocks of the form + # [[a, b], [b, -a]]. + E := RECOG.MatrixByEntries(F, n, n, + Concatenation(List([1..n], i -> [i, i, (-1) ^ (i + 1) * a]), + List([1..n], i -> [i, i + (-1) ^ (i + 1), b]))); + else + gramMatrix := Z(q) ^ 0 * DiagonalMat(Concatenation([zeta], + List([1..n - 1], i -> 1))); + generatorsOfSO := GeneratorsOfGroup(RECOG.ConjugateToSesquilinearForm(SO(epsilon, n, q), + "O", + gramMatrix)); + # Our standard bilinear form is given by the Gram matrix + # Diag(zeta, 1, ..., 1). The norm of [0, ..., 0, 1] under this + # bilinear form is 1, i.e. a square. + vectorOfSquareNorm := zeta ^ 0 * Concatenation(List([1..n - 1], i -> 0), + [1]); + D := RECOG.ReflectionMatrix(n, q, gramMatrix, vectorOfSquareNorm); + # Block diagonal matrix consisting of one block [[0, zeta], [1, 0]] + # and n / 2 - 1 blocks of the form [[a, b], [b, -a]]. + E := RECOG.MatrixByEntries(F, n, n, + Concatenation([[1, 2, zeta], [2, 1, zeta ^ 0]], + List([3..n], i -> [i, i, (-1) ^ (i + 1) * a]), + List([3..n], i -> [i, i + (-1) ^ (i + 1), b]))); + fi; + fi; + fi; + + return rec(generatorsOfSO := generatorsOfSO, D := D, E := E); +end; + + + +RECOG.MatrixGroup := function(F, gens) + if IsEmpty(gens) then + ErrorNoReturn(" cannot be empty"); + elif not IsField(F) then + ErrorNoReturn(" must be a field"); + fi; + return Group(List(gens, g -> ImmutableMatrix(F, g))); +end; + + + +RECOG.MatrixGroupWithSize := function(F, gens, size) + local result; + result := RECOG.MatrixGroup(F, gens); + SetSize(result, size); + return result; +end; + + + +# Assuming that the group G acts absolutely irreducibly, try to find a +# * symplectic form (if = S) or a +# * symmetric bilinear form (if = O) +# which is G-invariant or prove that no such form exists. +# +# We use this function instead of PreservedBilinearForms form the Forms package +# since PreservedBilinearForms seems to be buggy and unreliable (see also +# comment above UnitaryForm). +# +# In general, this function should only be used if one can be sure that +# preserves a bilinear form (but one does not know which one). +RECOG.BilinearForm := function(G, type) + local F, M, inverseTransposeM, counter, formMatrix, condition; + + if not type in ["S", "O"] then + ErrorNoReturn(" must be one of 'S', 'O'"); + fi; + # Set the condition the Gram matrix needs to satisfy for each of the + # possible types. + if type = "S" then + condition := x -> (x = - TransposedMat(x)); + elif type = "O" then + condition := x -> (x = TransposedMat(x)); + fi; + + F := DefaultFieldOfMatrixGroup(G); + + # Return stored bilinear form if it exists and is symplectic / symmetric + if HasInvariantBilinearForm(G) then + formMatrix := InvariantBilinearForm(G).matrix; + if condition(formMatrix) then + return ImmutableMatrix(F, formMatrix); + fi; + fi; + + M := GModuleByMats(GeneratorsOfGroup(G), F); + + if not MTX.IsIrreducible(M) then + ErrorNoReturn("BilinearForm failed - group is not irreducible"); + fi; + + # An element A of G acts as A ^ (-T) in MTX.DualModule(M) + inverseTransposeM := MTX.DualModule(M); + + counter := 0; + # As the MeatAxe is randomised, we might have to make some more trials to + # find a preserved symplectic / symmetric bilinear form if there is one; + # breaking after 1000 trials is just a "safety net" in case a group + # that does not preserve a symplectic / symmetric bilinear form is input. + while counter < 1000 do + counter := counter + 1; + + # If f: M -> inverseTransposeM is an isomorphism, it must respect + # multiplication by group elements, i.e. for A in G + # f(x * A) = f(x) * A ^ (-T) + # Let f be given by the matrix F, i.e. f: x -> x * F. Then we have + # (x * A) * F = x * F * A ^ (-T) + # Putting these results together for all vectors x gives + # A * F = F * A ^ (-T) + # <==> A * F * A ^ T = F, + # which is what we need. + formMatrix := MTX.IsomorphismModules(M, inverseTransposeM); + + if formMatrix <> fail then + # check if formMatrix is antisymmetric + if condition(formMatrix) then + return ImmutableMatrix(F, formMatrix); + fi; + if not MTX.IsAbsolutelyIrreducible(M) then + ErrorNoReturn("BilinearForm failed - group is not absolutely irreducible"); + fi; + fi; + od; + + return fail; +end; + + + +RECOG.SymplecticForm := function(G) + return RECOG.BilinearForm(G, "S"); +end; + + + +RECOG.SymmetricBilinearForm := function(G) + return RECOG.BilinearForm(G, "O"); +end; + + + +RECOG.QuadraticForm := function(G) + local d, F, formMatrix, polarFormMatrix, i, g, RightSideForLinSys, + MatrixForLinSys, x; + + d := DimensionOfMatrixGroup(G); + F := DefaultFieldOfMatrixGroup(G); + if not IsFinite(F) then + ErrorNoReturn("The base field of must be finite"); + fi; + + if HasInvariantQuadraticForm(G) then + formMatrix := InvariantQuadraticForm(G).matrix; + return ImmutableMatrix(F, formMatrix); + fi; + + # We first look for an invariant symmetric bilinear form of G, which will + # be the polar form of the desired quadratic form + polarFormMatrix := RECOG.SymmetricBilinearForm(G); + # The Gram matrix formMatrix of the quadratic form is upper triangular and + # polarFormMatrix = formMatrix + formMatrix ^ T, so the entries above the + # main diagonal of polarFormMatrix and formMatrix must be the same + formMatrix := List([1..d], i -> Concatenation(ListWithIdenticalEntries(i, Zero(F)), + polarFormMatrix[i]{[i + 1..d]})); + if Characteristic(F) <> 2 then + # In this case, the polar form determines the quadratic form completely + formMatrix := formMatrix + 1 / 2 * DiagonalMat(DiagonalOfMatrix(polarFormMatrix)); + else + # We are left to determine the diagonal entries of formMatrix. Let them + # be x1, ..., xd and X = diag(x1, ..., xd); furthermore, let U be the + # part of polarFormMatrix above the main diagonal (i.e. the current + # value of formMatrix). Then for the quadratic form X + U to be + # preserved, we need g * (X + U) * g ^ T to have the same diagonal + # entries as X + U, i.e. as X, for each generator g of G. + # + # Hence, we need xi = (g * U * g ^ T)_ii + (x1 * g[i, 1] ^ 2 + ... + xd * g[i, d] ^ 2) + # This leads to a linear system for the xi, which we solve. + + RightSideForLinSys := Concatenation(List(GeneratorsOfGroup(G), + g -> DiagonalOfMatrix(g * formMatrix * TransposedMat(g)))); + MatrixForLinSys := Concatenation(List(GeneratorsOfGroup(G), + g -> List([1..d], + i -> DiagonalOfMatrix(TransposedMat([g[i]{[1..d]}]) * [g[i]{[1..d]}])) + + IdentityMat(d, F))); + x := SolutionMat(TransposedMat(MatrixForLinSys), RightSideForLinSys); + + if x = fail then + return fail; + fi; + + formMatrix := formMatrix + DiagonalMat(x); + fi; + + return formMatrix; +end; + + + +# Compute the Gram matrix of the quadratic form corresponding to the bilinear +# form described by in odd characteristic. +RECOG.BilinearToQuadraticForm := function(gramMatrix) + local F, n, i, Q; + + F := DefaultFieldOfMatrix(gramMatrix); + + if Characteristic(F) = 2 then + ErrorNoReturn("Characteristic must be odd"); + fi; + + n := NrRows(gramMatrix); + Q := List(gramMatrix, ShallowCopy); + for i in [1..n] do + Q{[i + 1..n]}{[i]} := NullMat(n - i, 1, F); + Q[i, i] := gramMatrix[i, i] / 2; + od; + + return Q; +end; + +# One needs to ensure that the attribute DefaultFieldOfMatrixGroup is set +# correctly for ; this can be done, for example, by making the +# generators used during construction of the group immutable matrices over the +# appropriate field. +RECOG.ConjugateToSesquilinearForm := function(group, type, gramMatrix) + local gapForm, newForm, gapToCanonical, canonicalToNew, field, formMatrix, + result, d, q, broadType; + if not type in ["S", "O-B", "O-Q", "U"] then + ErrorNoReturn(" must be one of 'S', 'U', 'O-B', 'O-Q'"); + fi; + d := DimensionOfMatrixGroup(group); + field := DefaultFieldOfMatrixGroup(group); + if type = "S" or type = "O-B" then + if type = "S" then + broadType := type; + else + broadType := "O"; + fi; + formMatrix := RECOG.BilinearForm(group, broadType); + if formMatrix = fail then + if type = "S" then + ErrorNoReturn("No preserved symplectic form found for "); + else + ErrorNoReturn("No preserved symmetric bilinear form found for", + " "); + fi; + fi; + gapForm := BilinearFormByMatrix(formMatrix, field); + newForm := BilinearFormByMatrix(gramMatrix, field); + elif type = "U" then + if IsOddInt(DegreeOverPrimeField(field)) then + q := Size(field); + field := GF(q ^ 2); + fi; + formMatrix := RECOG.UnitaryForm(group); + if formMatrix = fail then + ErrorNoReturn("No preserved unitary form found for "); + fi; + gapForm := HermitianFormByMatrix(formMatrix, field); + newForm := HermitianFormByMatrix(gramMatrix, field); + else + # This is the case type = "O-Q" + formMatrix := RECOG.QuadraticForm(group); + if formMatrix = fail then + ErrorNoReturn("No preserved quadratic form found for "); + fi; + gapForm := QuadraticFormByMatrix(formMatrix, field); + newForm := QuadraticFormByMatrix(gramMatrix, field); + fi; + if gapForm = newForm then + # nothing to be done + result := group; + # The Forms package has a bug for d = 1 so we need to make this exception + elif d <> 1 then + # the following if condition can only ever be fulfilled if is an + # orthogonal group; there the case of even dimension is problematic since, + # in that case, there are two similarity classes of bilinear forms + if not WittIndex(gapForm) = WittIndex(newForm) then + ErrorNoReturn("The form preserved by must be similar to the form ", + "described by the Gram matrix ."); + fi; + gapToCanonical := BaseChangeHomomorphism(BaseChangeToCanonical(gapForm), + field); + canonicalToNew := BaseChangeHomomorphism(BaseChangeToCanonical(newForm) ^ (-1), + field); + result := RECOG.MatrixGroup(field, canonicalToNew(gapToCanonical(GeneratorsOfGroup(group)))); + + # Set useful attributes + UseIsomorphismRelation(group, result); + else + # replaces the Witt index check above + if IsZero(gramMatrix) <> IsZero(formMatrix) then + ErrorNoReturn("The form preserved by must be similar to the", + " form described by the Gram matrix ."); + fi; + result := group; + fi; + + if type = "S" then + SetInvariantBilinearForm(result, rec(matrix := gramMatrix)); + elif type = "O-B" then + SetInvariantQuadraticFormFromMatrix(result, RECOG.BilinearToQuadraticForm(gramMatrix)); + elif type = "U" then + SetInvariantSesquilinearForm(result, rec(matrix := gramMatrix)); + else + SetInvariantQuadraticFormFromMatrix(result, gramMatrix); + fi; + + return result; +end; + +# If preserves a sesquilinear form of type (one of "S", "U", "O" +# (in odd dimension), "O+" or "O-" (both in even dimension), return a group +# conjugate to preserving the standard form of that type. +# +# Also, one need to ensure that the attribute DefaultFieldOfMatrixGroup is set +# correctly for ; this can be done, for example, by making the +# generators used during construction of the group immutable matrices over the +# appropriate field. +RECOG.ConjugateToStandardForm := function(group, type) + local d, F, q, gapForm, broadType; + + # determining d (dimension of matrix group), F (base field) and q (order of + # F) plus some sanity checks + if not type in ["S", "O+", "O-", "O", "U"] then + ErrorNoReturn(" must be one of 'S', 'U', 'O+', 'O-', 'O'"); + fi; + F := DefaultFieldOfMatrixGroup(group); + d := DimensionOfMatrixGroup(group); + if type = "O" and IsEvenInt(d) then + ErrorNoReturn(" cannot be 'O' if the dimension of is even"); + elif type in ["O+", "O-"] and IsOddInt(d) then + ErrorNoReturn(" cannot be 'O+' or 'O-' if the dimension of", + " is odd"); + elif IsEvenInt(Size(F)) and IsOddInt(d) and type in ["O+", "O-", "O"] then + ErrorNoReturn("If is 'O+', 'O-' or 'O' and the size of is", + " even, must be even"); + fi; + if type in ["S", "O", "O+", "O-"] then + q := Size(F); + else + if IsSquareInt(Size(F)) then + q := RootInt(Size(F)); + else + # It might be that G is to be understood as a matrix group over + # GF(q ^ 2), but the matrices can actually be represented over a + # smaller field, which causes DefaultFieldOfMatrixGroup to return GF(q) + # instead of GF(q ^ 2) - we have to remedy this somehow ... + q := Size(F); + fi; + fi; + + # get standard GAP form + if type = "S" then + gapForm := InvariantBilinearForm(Sp(d, q)).matrix; + elif type = "U" then + gapForm := InvariantSesquilinearForm(GU(d, q)).matrix; + elif type = "O" then + gapForm := InvariantBilinearForm(Omega(d, q)).matrix; + elif type = "O+" then + if Characteristic(F) = 2 then + gapForm := InvariantQuadraticForm(Omega(1, d, q)).matrix; + else + gapForm := InvariantBilinearForm(Omega(1, d, q)).matrix; + fi; + elif type = "O-" then + if Characteristic(F) = 2 then + gapForm := InvariantQuadraticForm(Omega(-1, d, q)).matrix; + else + gapForm := InvariantBilinearForm(Omega(-1, d, q)).matrix; + fi; + fi; + + if type in ["O", "O+", "O-"] then + if Characteristic(F) = 2 then + broadType := "O-Q"; + else + broadType := "O-B"; + fi; + else + broadType := type; + fi; + + return RECOG.ConjugateToSesquilinearForm(group, broadType, gapForm); +end; + +# Let = [f1, f2, ..., ft] be a list of sesquilinear forms on the vector +# spaces F ^ d1, F ^ d2, ..., F ^ dt. Then we can lift these to a sesquilinear +# form f on the tensor product F ^ d1 x F ^ d2 x ... x F ^ dt by defining +# f(v1 x v2 x ... x vt, w1 x w2 x ... x wt) = f1(v1, w1)f2(v2, w2)...ft(vt, wt). +# +# Return the Gram matrix of f; the forms f1, f2, ..., ft must be given as their +# respective Gram matrices. +RECOG.LiftFormsToTensorProduct := function(forms, F) + local dims, d, t, newForm, i, j, iteri, iterj, indicesi, indicesj; + + dims := List(forms, NrRows); + d := Product(dims); + t := Length(dims); + newForm := NullMat(d, d, F); + dims := List(dims,d->[1..d]); + + iteri := IteratorOfCartesianProduct(dims); + for i in [1..d] do + indicesi := NextIterator(iteri); + iterj := IteratorOfCartesianProduct(dims); + for j in [1..d] do + indicesj := NextIterator(iterj); + newForm[i, j] := Product([1..t], k -> (forms[k])[indicesi[k], indicesj[k]]); + od; + od; + + return newForm; +end; + +# Return the standard orthogonal and corresponding bilinear form as used for +# constructions in [HR10], cf. section 3.1 loc. cit. +RECOG.StandardOrthogonalForm := function(epsilon, d, q) + local field, m, F, Q, gamma, xi; + + if not epsilon in [-1, 0, 1] then + ErrorNoReturn(" must be one of -1, 0, 1"); + fi; + if epsilon = 0 and IsEvenInt(d) then + ErrorNoReturn(" must be one of -1 or 1 if is even"); + fi; + if epsilon <> 0 and IsOddInt(d) then + ErrorNoReturn(" must be 0 if is odd"); + fi; + if IsEvenInt(q) and IsOddInt(d) then + ErrorNoReturn(" must be even if is even"); + fi; + + field := GF(q); + m := QuoInt(d, 2); + F := RECOG.AntidiagonalMat(d, field); + + if IsOddInt(d * q) then + Q := RECOG.AntidiagonalMat(One(field) * Concatenation(ListWithIdenticalEntries(m, 1), + [1 / 2], + ListWithIdenticalEntries(m, 0)), + field); + else + Q := RECOG.AntidiagonalMat(One(field) * Concatenation(ListWithIdenticalEntries(m, 1), + ListWithIdenticalEntries(m, 0)), + field); + if epsilon = -1 then + if IsEvenInt(q) then + gamma := RECOG.FindGamma(q); + else + xi := PrimitiveElement(GF(q ^ 2)); + gamma := xi ^ (q + 1) * (xi + xi ^ q) ^ (-2); + fi; + + F[m, m] := 2 * gamma ^ 0; + F[m + 1, m + 1] := 2 * gamma; + Q[m, m] := gamma ^ 0; + Q[m + 1, m + 1] := gamma; + fi; + fi; + + return rec(Q := Q, F := F); +end; + +RECOG.ConjugateModule := function(M, q) + return GModuleByMats(List(MTX.Generators(M), A -> RECOG.ApplyFunctionToEntries(A, x -> x ^ q)), + MTX.Field(M)); +end; + +# Assuming that the group G acts absolutely irreducibly, try to find a unitary +# form which is G-invariant or prove that no such form exists. +# +# We use this function instead of PreservedSesquilinearForms from the Forms +# package since PreservedSesquilinearForms seems to be buggy and unreliable. +# As an example, take the group generated by +# [ 1 0 0 ] [ z^39 z^9 z^24 ] +# [ z^33 z^14 z^26 ] and [ z^25 z^16 z^6 ] +# [ z^19 z^31 z^5 ] [ z^7 z^32 z^28 ] +# where z = Z(49), which does preserve a unitary form, but this is not +# recognised by PreservedSesquilinearForms, even after some 1000 calls of the +# function. +# +# In general, this function should only be used if one can be sure that +# preserves a unitary form (but one does not know which one). +RECOG.UnitaryForm := function(G) + local d, F, q, M, inverseHermitianConjugateM, formMatrix, row, col, x, + scalar, counter; + + d := DimensionOfMatrixGroup(G); + F := DefaultFieldOfMatrixGroup(G); + if not IsFinite(F) then + ErrorNoReturn("The base field of must be finite"); + fi; + if not IsEvenInt(DegreeOverPrimeField(F)) then + # It might be that G is to be understood as a matrix group over + # GF(q ^ 2), but the matrices can actually be represented over a + # smaller field, which causes DefaultFieldOfMatrixGroup to return GF(q) + # instead of GF(q ^ 2) - we have to remedy this somehow ... + q := Size(F); + else + q := RootInt(Size(F)); + fi; + + # Return stored sesquilinear form if it exists and is hermitian + if HasInvariantSesquilinearForm(G) then + formMatrix := InvariantSesquilinearForm(G).matrix; + if formMatrix = RECOG.HermitianConjugate(formMatrix, q) then + return ImmutableMatrix(F, formMatrix); + fi; + fi; + + M := GModuleByMats(GeneratorsOfGroup(G), F); + # An element A of G acts as A ^ (-T) in MTX.DualModule(M) and hence as + # HermitianConjugate(A, q) ^ (-1) in inverseHermitianConjugateM + inverseHermitianConjugateM := RECOG.ConjugateModule(MTX.DualModule(M), q); + + counter := 0; + scalar := fail; + # As the MeatAxe is randomised, we might have to make some more trials to + # find a preserved unitary form if there is one; breaking after 1000 trials + # is just a "safety net" in case a group that does not preserve a + # unitary form is input. + while scalar = fail and counter < 1000 do + counter := counter + 1; + + # If f: M -> inverseHermitianConjugateM is an isomorphism, it must respect + # multiplication by group elements, i.e. for A in G + # f(x * A) = f(x) * HermitianConjugate(A, q) ^ (-1). + # Let f be given by the matrix F, i.e. f: x -> x * F. Then we have + # (x * A) * F = x * F * HermitianConjugate(A, q) ^ (-1). + # Putting these results together for all vectors x gives + # A * F = F * HermitianConjugate(A, q) ^ (-1) + # <==> A * F * HermitianConjugate(A, q) = F, + # which is what we need. + formMatrix := MTX.IsomorphismModules(M, inverseHermitianConjugateM); + + # We now need to ensure that formMatrix is actually the matrix of a + # unitary form, which can be achieved by multiplying it by a scalar + if formMatrix <> fail then + if formMatrix <> RECOG.HermitianConjugate(formMatrix, q) then + # find a non-zero entry of formMatrix + row := First([1..d], x -> not IsZero(formMatrix[x])); + col := First([1..d], x -> not IsZero(formMatrix[row][x])); + if not IsZero(formMatrix[col, row]) then + # this must be 1 for formMatrix to be hermitian + x := formMatrix[row, col] * formMatrix[col, row] ^ (-q); + # multiplying formMatrix by scalar will ensure that x = 1, i.e. that + # formMatrix is hermitian + scalar := RootFFE(x, q - 1); + fi; + + if IsZero(formMatrix[col, row]) or scalar = fail then + if not MTX.IsAbsolutelyIrreducible(M) then + ErrorNoReturn("UnitaryForm failed - group is not absolutely irreducible"); + fi; + continue; + fi; + + # make formMatrix hermitian + formMatrix := scalar * formMatrix; + fi; + + if formMatrix <> RECOG.HermitianConjugate(formMatrix, q) and not MTX.IsAbsolutelyIrreducible(M) then + ErrorNoReturn("UnitaryForm failed - group is not absolutely irreducible"); + fi; + + return ImmutableMatrix(F, formMatrix); + fi; + od; + + return fail; +end; + +# Assuming that the group G acts absolutely irreducibly, try to find a +# * symplectic form (if = S) or a +# * symmetric bilinear form (if = O) +# which is G-invariant or prove that no such form exists. +# +# We use this function instead of PreservedBilinearForms form the Forms package +# since PreservedBilinearForms seems to be buggy and unreliable (see also +# comment above UnitaryForm). +# +# In general, this function should only be used if one can be sure that +# preserves a bilinear form (but one does not know which one). +RECOG.BilinearForm := function(G, type) + local F, M, inverseTransposeM, counter, formMatrix, condition; + + if not type in ["S", "O"] then + ErrorNoReturn(" must be one of 'S', 'O'"); + fi; + # Set the condition the Gram matrix needs to satisfy for each of the + # possible types. + if type = "S" then + condition := x -> (x = - TransposedMat(x)); + elif type = "O" then + condition := x -> (x = TransposedMat(x)); + fi; + + F := DefaultFieldOfMatrixGroup(G); + + # Return stored bilinear form if it exists and is symplectic / symmetric + if HasInvariantBilinearForm(G) then + formMatrix := InvariantBilinearForm(G).matrix; + if condition(formMatrix) then + return ImmutableMatrix(F, formMatrix); + fi; + fi; + + M := GModuleByMats(GeneratorsOfGroup(G), F); + + if not MTX.IsIrreducible(M) then + ErrorNoReturn("BilinearForm failed - group is not irreducible"); + fi; + + # An element A of G acts as A ^ (-T) in MTX.DualModule(M) + inverseTransposeM := MTX.DualModule(M); + + counter := 0; + # As the MeatAxe is randomised, we might have to make some more trials to + # find a preserved symplectic / symmetric bilinear form if there is one; + # breaking after 1000 trials is just a "safety net" in case a group + # that does not preserve a symplectic / symmetric bilinear form is input. + while counter < 1000 do + counter := counter + 1; + + # If f: M -> inverseTransposeM is an isomorphism, it must respect + # multiplication by group elements, i.e. for A in G + # f(x * A) = f(x) * A ^ (-T) + # Let f be given by the matrix F, i.e. f: x -> x * F. Then we have + # (x * A) * F = x * F * A ^ (-T) + # Putting these results together for all vectors x gives + # A * F = F * A ^ (-T) + # <==> A * F * A ^ T = F, + # which is what we need. + formMatrix := MTX.IsomorphismModules(M, inverseTransposeM); + + if formMatrix <> fail then + # check if formMatrix is antisymmetric + if condition(formMatrix) then + return ImmutableMatrix(F, formMatrix); + fi; + if not MTX.IsAbsolutelyIrreducible(M) then + ErrorNoReturn("BilinearForm failed - group is not absolutely irreducible"); + fi; + fi; + od; + + return fail; +end; + +RECOG.SymplecticForm := function(G) + return RECOG.BilinearForm(G, "S"); +end; + +RECOG.SymmetricBilinearForm := function(G) + return RECOG.BilinearForm(G, "O"); +end; + +RECOG.QuadraticForm := function(G) + local d, F, formMatrix, polarFormMatrix, i, g, RightSideForLinSys, + MatrixForLinSys, x; + + d := DimensionOfMatrixGroup(G); + F := DefaultFieldOfMatrixGroup(G); + if not IsFinite(F) then + ErrorNoReturn("The base field of must be finite"); + fi; + + if HasInvariantQuadraticForm(G) then + formMatrix := InvariantQuadraticForm(G).matrix; + return ImmutableMatrix(F, formMatrix); + fi; + + # We first look for an invariant symmetric bilinear form of G, which will + # be the polar form of the desired quadratic form + polarFormMatrix := RECOG.SymmetricBilinearForm(G); + # The Gram matrix formMatrix of the quadratic form is upper triangular and + # polarFormMatrix = formMatrix + formMatrix ^ T, so the entries above the + # main diagonal of polarFormMatrix and formMatrix must be the same + formMatrix := List([1..d], i -> Concatenation(ListWithIdenticalEntries(i, Zero(F)), + polarFormMatrix[i]{[i + 1..d]})); + if Characteristic(F) <> 2 then + # In this case, the polar form determines the quadratic form completely + formMatrix := formMatrix + 1 / 2 * DiagonalMat(DiagonalOfMatrix(polarFormMatrix)); + else + # We are left to determine the diagonal entries of formMatrix. Let them + # be x1, ..., xd and X = diag(x1, ..., xd); furthermore, let U be the + # part of polarFormMatrix above the main diagonal (i.e. the current + # value of formMatrix). Then for the quadratic form X + U to be + # preserved, we need g * (X + U) * g ^ T to have the same diagonal + # entries as X + U, i.e. as X, for each generator g of G. + # + # Hence, we need xi = (g * U * g ^ T)_ii + (x1 * g[i, 1] ^ 2 + ... + xd * g[i, d] ^ 2) + # This leads to a linear system for the xi, which we solve. + + RightSideForLinSys := Concatenation(List(GeneratorsOfGroup(G), + g -> DiagonalOfMatrix(g * formMatrix * TransposedMat(g)))); + MatrixForLinSys := Concatenation(List(GeneratorsOfGroup(G), + g -> List([1..d], + i -> DiagonalOfMatrix(TransposedMat([g[i]{[1..d]}]) * [g[i]{[1..d]}])) + + IdentityMat(d, F))); + x := SolutionMat(TransposedMat(MatrixForLinSys), RightSideForLinSys); + + if x = fail then + return fail; + fi; + + formMatrix := formMatrix + DiagonalMat(x); + fi; + + return formMatrix; +end; + + +############################################################################# +############################################################################# +################## Old function from RECOG package ########################## +############################################################################# +############################################################################# + + +RECOG.DerivedSubgroupMonteCarlo := function(g, NumberGenerators) + local gens,gens2,i,x,y; + gens := []; + for i in [1..Maximum([NumberGenerators, Size(GeneratorsOfGroup(g)) * 2 + 10])] do + x := PseudoRandom(g); + y := PseudoRandom(g); + Add(gens,Comm(x,y)); + od; + gens2 := FastNormalClosure(GeneratorsOfGroup(g),gens,10); + return GroupWithGenerators(gens2); +end; diff --git a/read.g b/read.g index 76e33bcc..8bcf59b4 100644 --- a/read.g +++ b/read.g @@ -55,7 +55,25 @@ ReadPackage("recog","gap/projective/almostsimple/threeelorders.gi"); ReadPackage("recog","gap/projective/almostsimple.gi"); ReadPackage("recog","gap/projective/almostsimple/lietype.gi"); ReadPackage("recog","gap/projective/almostsimple/hints.gi"); -ReadPackage("recog","gap/projective/classicalnatural.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/utils/achieve.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/utils/utils.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/SL/BaseCase.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/SL/GoingDown.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/SL/GoingUp.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/SL/main.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/SL/sl2_BlackBox.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/Sp/BaseCase.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/Sp/GoingDown.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/Sp/GoingUp.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/Sp/main.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/SU/BaseCase.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/SU/GoingDown.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/SU/GoingUp.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/SU/main.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/O/BaseCase.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/O/GoingDown.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/O/GoingUp.gi"); +ReadPackage("recog","gap/projective/constructive_recognition/O/main.gi"); ReadPackage("recog","gap/projective/AnSnOnFDPM.gi"); # All the method installations are now here: From cf89f1492403396bae245a2bb76f97026ecfd4e8 Mon Sep 17 00:00:00 2001 From: Daniel Rademacher Date: Fri, 14 Jun 2024 11:08:53 +0200 Subject: [PATCH 02/13] Remove some deleted files from read.g --- read.g | 12 ------------ 1 file changed, 12 deletions(-) diff --git a/read.g b/read.g index 8bcf59b4..690d2066 100644 --- a/read.g +++ b/read.g @@ -62,18 +62,6 @@ ReadPackage("recog","gap/projective/constructive_recognition/SL/GoingDown.gi"); ReadPackage("recog","gap/projective/constructive_recognition/SL/GoingUp.gi"); ReadPackage("recog","gap/projective/constructive_recognition/SL/main.gi"); ReadPackage("recog","gap/projective/constructive_recognition/SL/sl2_BlackBox.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/Sp/BaseCase.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/Sp/GoingDown.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/Sp/GoingUp.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/Sp/main.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/SU/BaseCase.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/SU/GoingDown.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/SU/GoingUp.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/SU/main.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/O/BaseCase.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/O/GoingDown.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/O/GoingUp.gi"); -ReadPackage("recog","gap/projective/constructive_recognition/O/main.gi"); ReadPackage("recog","gap/projective/AnSnOnFDPM.gi"); # All the method installations are now here: From df70f5e9aeccb8b998962e5077aff4cd604f0556 Mon Sep 17 00:00:00 2001 From: Daniel Rademacher Date: Tue, 18 Jun 2024 15:21:50 +0200 Subject: [PATCH 03/13] Add fix for manual crash. --- gap/projective/classicalnatural.gi | 3 +++ 1 file changed, 3 insertions(+) create mode 100644 gap/projective/classicalnatural.gi diff --git a/gap/projective/classicalnatural.gi b/gap/projective/classicalnatural.gi new file mode 100644 index 00000000..b6e732b2 --- /dev/null +++ b/gap/projective/classicalnatural.gi @@ -0,0 +1,3 @@ +#! @BeginChunk ClassicalNatural +#! TODO +#! @EndChunk From 7fbadd0318da3b45f3d1c2fee46a824f0fe67210 Mon Sep 17 00:00:00 2001 From: Daniel Rademacher Date: Tue, 18 Jun 2024 15:41:56 +0200 Subject: [PATCH 04/13] Implement Max remark --- gap/projective/classicalnatural.gi | 3 --- makedoc.g | 3 +++ 2 files changed, 3 insertions(+), 3 deletions(-) delete mode 100644 gap/projective/classicalnatural.gi diff --git a/gap/projective/classicalnatural.gi b/gap/projective/classicalnatural.gi deleted file mode 100644 index b6e732b2..00000000 --- a/gap/projective/classicalnatural.gi +++ /dev/null @@ -1,3 +0,0 @@ -#! @BeginChunk ClassicalNatural -#! TODO -#! @EndChunk diff --git a/makedoc.g b/makedoc.g index 1647cd9a..73607924 100644 --- a/makedoc.g +++ b/makedoc.g @@ -32,6 +32,9 @@ scan_dirs := [ "gap/perm", "gap/projective", "gap/projective/almostsimple", + "gap/projective/constructive_recognition", + "gap/projective/constructive_recognition/SL", + "gap/projective/constructive_recognition/utils", ]; AutoDoc(rec( From e5e2bd7b13c3a960f3b9c59e9f4246487b097778 Mon Sep 17 00:00:00 2001 From: ThomasBreuer Date: Wed, 19 Jun 2024 14:02:29 +0200 Subject: [PATCH 05/13] preliminary fix in `makedoc.g` Make sure that the correct path of `regen_doc.g` is used when PackageManager builds the documentation. --- makedoc.g | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/makedoc.g b/makedoc.g index 73607924..66626f76 100644 --- a/makedoc.g +++ b/makedoc.g @@ -21,7 +21,7 @@ if fail = LoadPackage("AutoDoc", ">= 2019.07.03") then ErrorNoReturn("AutoDoc 2019.07.03 or newer is required"); fi; -Read("regen_doc.g"); +Read(Filename(DirectoryCurrent(), "regen_doc.g")); scan_dirs := [ "doc", From f2b49f55983b35d7c4da7ef39a1eb9ced0192ebe Mon Sep 17 00:00:00 2001 From: Daniel Rademacher Date: Wed, 19 Jun 2024 14:19:52 +0200 Subject: [PATCH 06/13] New version and add dependencies --- PackageInfo.g | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/PackageInfo.g b/PackageInfo.g index 335d9702..b8a05f04 100644 --- a/PackageInfo.g +++ b/PackageInfo.g @@ -18,7 +18,7 @@ SetPackageInfo( rec( PackageName := "recog", Subtitle := "A package for constructive recognition of permutation and matrix groups", -Version := "1.4.2", +Version := "1.4.3", Date := "27/09/2022", # dd/mm/yyyy format License := "GPL-3.0-or-later", @@ -268,6 +268,7 @@ Dependencies := rec( GAP := ">=4.11", NeededOtherPackages := [ ["AtlasRep", ">= 1.4.0"], + ["Alnuth", ">= 3.2.1"], ["FactInt", ">= 1.5.2"], ["Forms", ">= 1.2"], ["genss", ">= 1.3"], From 6502e3b2ce1370bf737116d82c5693dd372ea90c Mon Sep 17 00:00:00 2001 From: Max Horn Date: Mon, 17 Jun 2024 16:05:45 +0200 Subject: [PATCH 07/13] typo --- gap/projective/constructive_recognition/SL/GoingDown.gi | 4 ++-- gap/projective/constructive_recognition/utils/utils.gi | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/gap/projective/constructive_recognition/SL/GoingDown.gi b/gap/projective/constructive_recognition/SL/GoingDown.gi index b425fddc..c81ee180 100644 --- a/gap/projective/constructive_recognition/SL/GoingDown.gi +++ b/gap/projective/constructive_recognition/SL/GoingDown.gi @@ -448,14 +448,14 @@ local out, list, out2, currentdim, check, slplist, slpToSmallerGroup, baselist, # We still have to compute the vector space on which the matrices act in the input group - Info(InfoRecog,2,"Debugg Info:\n"); + Info(InfoRecog,2,"Debug Info:\n"); Info(InfoRecog,2,"Dimension FirstElement: "); Info(InfoRecog,2,out[2]); Info(InfoRecog,2,"\n"); Info(InfoRecog,2,"Dimension SecondElement: "); Info(InfoRecog,2,out2[2]); Info(InfoRecog,2,"\n"); - Info(InfoRecog,2,"End Debugg Info. \n"); + Info(InfoRecog,2,"End Debug Info. \n"); Info(InfoRecog,2,"New Dimension: "); Info(InfoRecog,2,out[2]+out2[2]); diff --git a/gap/projective/constructive_recognition/utils/utils.gi b/gap/projective/constructive_recognition/utils/utils.gi index 016f402a..08dd2fd4 100644 --- a/gap/projective/constructive_recognition/utils/utils.gi +++ b/gap/projective/constructive_recognition/utils/utils.gi @@ -194,14 +194,14 @@ RECOG.constructppdTwoStingray := function(g,dim,q,type,form) list:=[Group(out[1],out2[1]),out[2]+out2[2],q,fail,form]; currentdim := list[2]; - Info(InfoRecog,2,"Debugg Info:\n"); + Info(InfoRecog,2,"Debug Info:\n"); Info(InfoRecog,2,"Dimension FirstElement: "); Info(InfoRecog,2,out[2]); Info(InfoRecog,2,"\n"); Info(InfoRecog,2,"Dimension SecondElement: "); Info(InfoRecog,2,out2[2]); Info(InfoRecog,2,"\n"); - Info(InfoRecog,2,"End Debugg Info. \n"); + Info(InfoRecog,2,"End Debug Info. \n"); Info(InfoRecog,2,"New Dimension: "); Info(InfoRecog,2,out[2]+out2[2]); From 935344024d8676641e13fcb6bafae622abc806bf Mon Sep 17 00:00:00 2001 From: Max Horn Date: Thu, 20 Jun 2024 12:00:51 +0200 Subject: [PATCH 08/13] simplify --- .../constructive_recognition/SL/main.gi | 43 ++++++------------- 1 file changed, 13 insertions(+), 30 deletions(-) diff --git a/gap/projective/constructive_recognition/SL/main.gi b/gap/projective/constructive_recognition/SL/main.gi index d5796ebf..96cb3968 100644 --- a/gap/projective/constructive_recognition/SL/main.gi +++ b/gap/projective/constructive_recognition/SL/main.gi @@ -25,15 +25,7 @@ RECOG.FindStdGens_SL := function(sld) - return RECOG.FindStdGens_SL2(sld,DimensionOfMatrixGroup(sld)); - -end; - - - -RECOG.FindStdGens_SL2 := function(sld,IsoDim) - - # Group generated by input must be isomorphic SL(IsoDim,q) + # Group generated by input must be isomorphic SL(d,q) # gens of sld must be gens for SL(d,q) in its natural rep with memory # This function calls RECOG.SLn_constructsl2 and then extends @@ -118,7 +110,7 @@ RECOG.FindStdGens_SL2 := function(sld,IsoDim) fakegens := GeneratorsWithMemory(fakegens); sl2gensf := ResultOfStraightLineProgram(slptosl2,fakegens); sl2stdf := ResultOfStraightLineProgram(slpsl2std,sl2gensf); - std := rec( f := f, d := d, GoalDim := IsoDim, n := 2, bas := bas, basi := basi, + std := rec( f := f, d := d, GoalDim := d, n := 2, bas := bas, basi := basi, sld := sld, sldf := fakegens, slnstdf := sl2stdf, p := p, ext := ext ); Info(InfoRecog,2,"Going up to SL_d again..."); @@ -137,25 +129,16 @@ end; RECOG.FindStdGensSmallerMatrices_SL := function(sld) - return RECOG.FindStdGensSmallerMatrices_SL2(sld,DimensionOfMatrixGroup(sld)); - -end; - - - - -RECOG.FindStdGensSmallerMatrices_SL2 := function(sld,IsoDim) - -# Group generated by input must be isomorphic SL(IsoDim,q) - -# gens of sld must be gens for SL(d,q) in its natural rep with memory -# This function calls RECOG.SLn_constructsl2 and then extends -# the basis to a basis of the full row space and calls -# RECOG.SLn_UpStep often enough. Finally it returns an slp such -# that the SL(d,q) standard generators with respect to this basis are -# expressed by the slp in terms of the original generators of sld. -local V,b,bas,basi,basit,d,data,ext,fakegens,id,nu,nu2,p,q,resl2,sl2,sl2gens, - sl2gensf,sl2genss,sl2stdf,slp,slpsl2std,slptosl2,st,std,stdgens,i,ex,f; + # Group generated by input must be isomorphic SL(d,q) + + # gens of sld must be gens for SL(d,q) in its natural rep with memory + # This function calls RECOG.SLn_constructsl2 and then extends + # the basis to a basis of the full row space and calls + # RECOG.SLn_UpStep often enough. Finally it returns an slp such + # that the SL(d,q) standard generators with respect to this basis are + # expressed by the slp in terms of the original generators of sld. + local V,b,bas,basi,basit,d,data,ext,fakegens,id,nu,nu2,p,q,resl2,sl2,sl2gens, + sl2gensf,sl2genss,sl2stdf,slp,slpsl2std,slptosl2,st,std,stdgens,i,ex,f; # Some setup: f := FieldOfMatrixGroup(sld); @@ -229,7 +212,7 @@ local V,b,bas,basi,basit,d,data,ext,fakegens,id,nu,nu2,p,q,resl2,sl2,sl2gens, fakegens := GeneratorsWithMemory(fakegens); sl2gensf := ResultOfStraightLineProgram(slptosl2,fakegens); sl2stdf := ResultOfStraightLineProgram(slpsl2std,sl2gensf); - std := rec( f := f, d := d, GoalDim := IsoDim, n := 2, bas := bas, basi := basi, + std := rec( f := f, d := d, GoalDim := d, n := 2, bas := bas, basi := basi, sld := sld, sldf := fakegens, slnstdf := sl2stdf, p := p, ext := ext ); Info(InfoRecog,2,"Going up to SL_d again..."); From bbb4f706e791b36eb98b81bf1707d3f5958ca01b Mon Sep 17 00:00:00 2001 From: Max Horn Date: Thu, 20 Jun 2024 12:01:32 +0200 Subject: [PATCH 09/13] cleanup --- gap/projective/constructive_recognition/SL/GoingDown.gi | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/gap/projective/constructive_recognition/SL/GoingDown.gi b/gap/projective/constructive_recognition/SL/GoingDown.gi index c81ee180..9507c126 100644 --- a/gap/projective/constructive_recognition/SL/GoingDown.gi +++ b/gap/projective/constructive_recognition/SL/GoingDown.gi @@ -495,7 +495,8 @@ local r,h,slp,sl2,baselist,gens,b,sl2gens,sl2genss,f,eigenspacelist,subspaces,ei Add(slp,SLPOfElms(GeneratorsOfGroup(h))); #h := RECOG.LinearActionRepresentationSmallerMatrices(h); #Add(baselist,h[3]); - h := GeneratorsOfGroup(h); h := RECOG.ConstructSL4SmallerMatrices(h[1],h[2],q); + h := GeneratorsOfGroup(h); + h := RECOG.ConstructSL4SmallerMatrices(h[1],h[2],q); Add(baselist,h[2]); Add(eigenspacelist,h[3]); h[1] := GroupWithMemory(h[1]); @@ -574,13 +575,13 @@ RECOG.Computesl2Subspace:=function(generators,eigenspaceGenerators) local result, i, gens, j, combination, vec, comb, zerovec, sl2eigenspacebase, newsl2eigenspacevectors, ele; if Size(generators) = 1 then - # We startet with a SL(4,q) + # We started with a SL(4,q) # Just return the 2-dimensional subspace # TODO return eigenspacebase!!! See else case sl2eigenspacebase := eigenspaceGenerators[1]; - zerovec := Zero(result[1,1]) * result[1]; + zerovec := ZeroOfBaseDomain(result) * result[1]; for ele in eigenspaceGenerators[2] do vec := zerovec; for j in [1..Size(ele)] do @@ -594,7 +595,7 @@ local result, i, gens, j, combination, vec, comb, zerovec, sl2eigenspacebase, ne # We start with the 1xd vectors result := generators[1]; sl2eigenspacebase := eigenspaceGenerators[1]; - zerovec := Zero(result[1,1]) * result[1]; + zerovec := ZeroOfBaseDomain(result) * result[1]; for ele in eigenspaceGenerators[2] do vec := zerovec; for j in [1..Size(ele)] do From 261634462123c380e52f16b25dc7b09878ab0062 Mon Sep 17 00:00:00 2001 From: Max Horn Date: Thu, 20 Jun 2024 12:01:48 +0200 Subject: [PATCH 10/13] Debug GenerateRandomKernelElementsAndOptionallyVerifyThem --- gap/base/kernel.gi | 3 +++ 1 file changed, 3 insertions(+) diff --git a/gap/base/kernel.gi b/gap/base/kernel.gi index b5ad9258..bb389d21 100644 --- a/gap/base/kernel.gi +++ b/gap/base/kernel.gi @@ -36,7 +36,10 @@ BindGlobal( "GenerateRandomKernelElementsAndOptionallyVerifyThem", verificationSuccess := true; # We generate a random element of the kernel as the quotient of a random # element and the preimage of its image under the homomorphism. + Info(InfoRecog,3, "GenerateRandomKernelElementsAndOptionallyVerifyThem, n = ", n); + for i in [1 .. n] do + Info(InfoRecog,3, "GenerateRandomKernelElementsAndOptionallyVerifyThem: ", i, " / ", n); # Finding kernel generators and immediate verification must use # different random elements! This is ensured by using the same stamp # in both situations. From 01c27032bf23422a7432e4e086df7b4bdfcd1855 Mon Sep 17 00:00:00 2001 From: Max Horn Date: Thu, 20 Jun 2024 12:37:08 +0200 Subject: [PATCH 11/13] Remove some dead (commented out) code --- .../constructive_recognition/utils/achieve.gi | 1479 ----------------- 1 file changed, 1479 deletions(-) diff --git a/gap/projective/constructive_recognition/utils/achieve.gi b/gap/projective/constructive_recognition/utils/achieve.gi index 8cabd4f9..6a48d7c3 100644 --- a/gap/projective/constructive_recognition/utils/achieve.gi +++ b/gap/projective/constructive_recognition/utils/achieve.gi @@ -17,405 +17,6 @@ ## ############################################################################# -# Obsolete stuff? - -# RECOG.RelativePrimeToqm1Part := function(q,n) -# local x,y; -# x := (q^n-1)/(q-1); -# repeat -# y := x/(q-1); -# x := NumeratorRat(y); -# until DenominatorRat(y) = q-1; -# return x; -# end; -# -# RECOG.SearchForElByCharPolFacts := function(g,f,degs,limit) -# local count,degrees,factors,pol,y; -# count := 0; -# while true do # will be left by return -# if InfoLevel(InfoRecog) >= 3 then Print(".\c"); fi; -# y:=PseudoRandom(g); -# pol:=CharacteristicPolynomial(f,f,StripMemory(y),1); -# factors:=Factors(PolynomialRing(f),pol); -# degrees:=List(factors,Degree); -# SortParallel(degrees,factors); -# if degrees = degs then -# if InfoLevel(InfoRecog) >= 3 then Print("\n"); fi; -# return rec( el := y, factors := factors, degrees := degrees ); -# fi; -# count := count + 1; -# if count >= limit then return fail; fi; -# od; -# end; -# -# RECOG.SL_Even_godownone:=function(g,subspg,q,d) -# local n,y,yy,yyy,ready,order,es,null,subsph,z,x,a,b,c,h,r,divisors,cent,i, -# pol,factors,degrees; -# -# n:=DimensionOfMatrixGroup(g); -# #d:=Dimension(subspg); -# repeat -# ready:=false; -# y:=PseudoRandom(g); -# pol:=CharacteristicPolynomial(GF(q),GF(q),StripMemory(y),1); -# factors:=Factors(pol); -# degrees:=List(factors,Degree); -# if d-1 in degrees then -# order:=Order(y); -# yy:=y^(order/Gcd(order,q-1)); -# if not IsOne(yy) then -# es:= Eigenspaces(GF(q),StripMemory(yy)); -# es:=Filtered(es,x->Dimension(x)=d-1 and IsSubspace(subspg,x)); -# if Length(es)>0 then -# subsph:=es[1]; -# ready:=true; -# fi; -# yyy:=y^(Gcd(order,q-1)); -# fi; -# fi; -# until ready; -# -# cent:=[yyy]; -# for i in [1..4] do -# z:=PseudoRandom(g); -# x:=yy^z; -# a := x; -# b := x^yy; -# c := x^(yy^2); -# h := Group(a,b,c); -# ready:=false; -# repeat -# r:=PseudoRandom(h); -# r:=r^(q*(q+1)); -# if not IsOne(r) and r*yy=yy*r then -# Add(cent,yyy^r); -# ready:=true; -# fi; -# until ready=true; -# od; -# return [Group(cent), subsph]; -# end; -# -# RECOG.SL_FindSL2 := function(g,f) -# local V,a,bas,c,count,ev,gens,genss,genssm,gl4,h,i,j,n,ns,o,pos,pow,pr,q,r, -# res,sl2gens,sl3,slp,std,v,w,y,z,zz; -# q := Size(f); -# n := DimensionOfMatrixGroup(g); -# if q = 2 then -# # We look for a transvection: -# while true do # will be left by break -# r := RECOG.SearchForElByCharPolFacts(g,f,[1,1,n-2],3*n+20); -# if r = fail then return fail; fi; -# y := r.el^(q^(n-2)-1); -# if not IsOne(y) and IsOne(y^2) then break; fi; -# od; -# # Find a good random conjugate: -# repeat -# z := y^PseudoRandom(g); -# until Order(z*y) = 3; -# gens := [y,z]; -# o := IdentityMat(n,f); -# w := []; -# for i in [1..2] do -# for j in [1..n] do -# w[i] := o[j]*gens[i]-o[j]; -# if not IsZero(w[i]) then break; fi; -# od; -# od; -# return [Group(gens),VectorSpace(GF(q),w)]; -# fi; -# if q = 3 and n = 3 then -# std := RECOG.MakeSL_StdGens(3,1,2,3); -# slp := RECOG.FindStdGensUsingBSGS(g,Concatenation(std.s,std.t), -# false,true); -# if slp = fail then return fail; fi; -# h := Group(ResultOfStraightLineProgram(slp,GeneratorsOfGroup(g))); -# o := IdentityMat(3,f); -# return [h,VectorSpace(f,o{[1..2]})]; -# fi; -# if q = 3 and n = 4 then -# std := RECOG.MakeSL_StdGens(3,1,2,4); -# slp := RECOG.FindStdGensUsingBSGS(g,Concatenation(std.s,std.t), -# false,true); -# if slp = fail then return fail; fi; -# h := Group(ResultOfStraightLineProgram(slp,GeneratorsOfGroup(g))); -# o := IdentityMat(4,f); -# return [h,VectorSpace(f,o{[1..2]})]; -# fi; -# if q = 3 then -# # We look for a transvection: -# while true do # will be left by break -# r := RECOG.SearchForElByCharPolFacts(g,f,[1,1,n-2],3*n+20); -# if r = fail then return fail; fi; -# y := r.el^(q^(n-2)-1); -# if not IsOne(y) and IsOne(y^3) then break; fi; -# od; -# # Find a two good random conjugates: -# while true do # will be left by return -# z := y^PseudoRandom(g); -# zz := y^PseudoRandom(g); -# gens := [y,z,zz]; -# o := IdentityMat(n,f); -# ns := []; -# for j in [1..3] do -# for i in [1..n] do -# w := o[i]*gens[j]-o[i]; -# if not IsZero(w) then break; fi; -# od; -# # Since y has order y at least one basis vector is moved. -# ns[j] := w; -# od; -# V := VectorSpace(f,ns); -# bas := Basis(V,ns); -# genss := List(StripMemory(gens), -# x->List(ns,i->Coefficients(bas,i*x))); -# genssm := GeneratorsWithMemory(genss); -# sl3 := Group(genssm); -# pr := ProductReplacer(genssm,rec( maxdepth := 400, scramble := 0, -# scramblefactor := 0 ) ); -# sl3!.pseudorandomfunc := [rec(func := Next,args := [pr])]; -# res := RECOG.SL_FindSL2(sl3,f); -# if res = fail then -# if InfoLevel(InfoRecog) >= 3 then Print("#\c"); fi; -# continue; -# fi; -# slp := SLPOfElms(GeneratorsOfGroup(res[1])); -# sl2gens := ResultOfStraightLineProgram(slp,gens); -# ns := BasisVectors(Basis(res[2])) * ns; -# ConvertToMatrixRep(ns,q); -# return [Group(sl2gens),VectorSpace(f,ns)]; -# od; -# fi; -# if q = 4 and n = 3 then -# std := RECOG.MakeSL_StdGens(2,2,2,3); -# slp := RECOG.FindStdGensUsingBSGS(g,Concatenation(std.s,std.t), -# false,true); -# if slp = fail then return fail; fi; -# h := Group(ResultOfStraightLineProgram(slp,GeneratorsOfGroup(g))); -# o := IdentityMat(3,f); -# return [h,VectorSpace(f,o{[1..2]})]; -# fi; -# if q = 5 and n = 4 then -# std := RECOG.MakeSL_StdGens(5,1,2,4); -# slp := RECOG.FindStdGensUsingBSGS(g,Concatenation(std.s,std.t), -# false,true); -# if slp = fail then return fail; fi; -# h := Group(ResultOfStraightLineProgram(slp,GeneratorsOfGroup(g))); -# o := IdentityMat(4,f); -# return [h,VectorSpace(f,o{[1..2]})]; -# fi; -# if n mod (q-1) <> 0 and q <> 3 then # The generic case: -# # We look for an element with n-1 dimensional eigenspace: -# count := 0; -# while true do # will be left by break -# count := count + 1; -# if count > 20 then return fail; fi; -# r := RECOG.SearchForElByCharPolFacts(g,f,[1,n-1],3*n+20); -# if r = fail then return fail; fi; -# pow := RECOG.RelativePrimeToqm1Part(q,n-1); -# y := r.el^pow; -# o := Order(y); -# if o mod (q-1) = 0 then -# y := y^(o/(q-1)); -# break; -# fi; -# od; -# # Now y has order q-1 and and n-1 dimensional eigenspace -# ev := -Value(r.factors[1],0*Z(q)); -# ns := NullspaceMat(StripMemory(r.el)-ev*StripMemory(One(y))); -# # this is a 1xn matrix now -# ns := ns[1]; -# pos := PositionNonZero(ns); -# ns := (ns[pos]^-1) * ns; -# count := 0; -# while true do # will be left by break -# count := count + 1; -# if count > 20 then return fail; fi; -# a := PseudoRandom(g); -# v := OnLines(ns,a); -# z := y^a; -# if OnLines(v,y) <> v and OnLines(ns,z) <> ns then -# # Now y and z most probably generate a GL(2,q), we need -# # the derived subgroup and then check: -# c := Comm(y,z); -# sl2gens := FastNormalClosure([y,z],[c],1); -# V := VectorSpace(f,[ns,v]); -# bas := Basis(V,[ns,v]); -# genss := List(sl2gens,x->List([ns,v],i->Coefficients(bas,i*x))); -# if RECOG.IsThisSL2Natural(genss,f) then break; fi; -# if InfoLevel(InfoRecog) >= 3 then Print("$\c"); fi; -# else -# if InfoLevel(InfoRecog) >= 3 then Print("-\c"); fi; -# fi; -# od; -# if InfoLevel(InfoRecog) >= 3 then Print("\n"); fi; -# return [Group(sl2gens),VectorSpace(f,[ns,v])]; -# fi; -# # Now q-1 does divide n, we have to do something else: -# # We look for an element with n-2 dimensional eigenspace: -# while true do # will be left by break -# r := RECOG.SearchForElByCharPolFacts(g,f,[1,1,n-2],5*n+20); -# if r = fail then return fail; fi; -# pow := RECOG.RelativePrimeToqm1Part(q,n-2); -# y := r.el^pow; -# o := Order(y); -# if o mod (q-1) = 0 then -# y := y^(o/(q-1)); -# if RECOG.IsScalarMat(y) = false then break; fi; -# fi; -# od; -# # Now y has order q-1 and n-2 dimensional eigenspace -# if r.factors[1] <> r.factors[2] then -# ev := -Value(r.factors[1],0*Z(q)); -# ns := NullspaceMat(StripMemory(r.el)-ev*StripMemory(One(y))); -# if not IsMutable(ns) then ns := MutableCopyMat(ns); fi; -# # this is a 1xn matrix now -# ev := -Value(r.factors[2],0*Z(q)); -# Append(ns,NullspaceMat(StripMemory(r.el)-ev*StripMemory(One(y)))); -# # ns now is a 2xn matrix -# else -# ev := -Value(r.factors[1],0*Z(q)); -# ns := NullspaceMat((StripMemory(r.el) -# -ev*StripMemory(One(y)))^2); -# if not IsMutable(ns) then ns := MutableCopyMat(ns); fi; -# fi; -# -# count := 0; -# while true do # will be left by break -# count := count + 1; -# if count > 20 then return fail; fi; -# if Length(ns) > 2 then ns := ns{[1..2]}; fi; -# a := PseudoRandom(g); -# Append(ns,ns * a); -# if RankMat(ns) < 4 then -# if InfoLevel(InfoRecog) >= 3 then Print("+\c"); fi; -# continue; -# fi; -# z := y^a; -# # Now y and z most probably generate a GL(4,q), we need -# # the derived subgroup and then check: -# V := VectorSpace(f,ns); -# bas := Basis(V,ns); -# genss := List([y,z],x->List(ns,i->Coefficients(bas,i*x))); -# genssm := GeneratorsWithMemory(genss); -# gl4 := Group(genssm); -# pr := ProductReplacer(genssm,rec( maxdepth := 400, scramble := 0, -# scramblefactor := 0 ) ); -# gl4!.pseudorandomfunc := [rec(func := Next,args := [pr])]; -# res := RECOG.SL_FindSL2(gl4,f); -# if res = fail then -# if InfoLevel(InfoRecog) >= 3 then Print("#\c"); fi; -# continue; -# fi; -# slp := SLPOfElms(GeneratorsOfGroup(res[1])); -# sl2gens := ResultOfStraightLineProgram(slp,[y,z]); -# ns := BasisVectors(Basis(res[2])) * ns; -# return [Group(sl2gens),VectorSpace(f,ns)]; -# od; -# return fail; -# end; -# -# -# RECOG.SL_Even_constructdata:=function(g,q) -# local S,a,b,degrees,eva,factors,gens,h,i,n,ns,o,pol,pos,ready,ready2, -# ready3,subgplist,w,ww,y,yy,z; -# -# n:=DimensionOfMatrixGroup(g); -# -# if q=2 then -# repeat -# ready:=false; -# y:=PseudoRandom(g); -# pol:=CharacteristicPolynomial(GF(q),GF(q),StripMemory(y),1); -# factors:=Factors(pol); -# degrees:=List(factors,Degree); -# if SortedList(degrees)=[1,1,n-2] then -# yy := y^(q^(n-2)-1); -# if not IsOne(yy) and IsOne(yy^2) then ready := true; fi; -# fi; -# until ready; -# repeat -# z := yy^PseudoRandom(g); -# until Order(z*yy) = 3; -# o := OneMutable(z); -# i := 1; -# while i <= n do -# w := o[i]*yy-o[i]; -# if not IsZero(w) then break; fi; -# i := i + 1; -# od; -# i := 1; -# while i <= n do -# ww := o[i]*z-o[i]; -# if not IsZero(ww) then break; fi; -# i := i + 1; -# od; -# return [Group(z,yy),VectorSpace(GF(2),[w,ww])]; -# else -# #case of q <> 2 -# repeat -# ready:=false; -# y:=PseudoRandom(g); -# if InfoLevel(InfoRecog) >= 3 then Print(".\c"); fi; -# pol:=CharacteristicPolynomial(GF(q),GF(q),StripMemory(y),1); -# factors:=Factors(pol); -# degrees:=List(factors,Degree); -# if n-1 in degrees then -# yy := y^(RECOG.RelativePrimeToqm1Part(q,n-1)); -# o := Order(yy); -# if o mod (q-1) = 0 then -# yy := yy^(o/(q-1)); -# ready := true; -# fi; -# fi; -# until ready; -# if InfoLevel(InfoRecog) >= 3 then Print("\n"); fi; -# -# ready2:=false; -# ready3:=false; -# repeat -# gens:=[yy]; -# a := PseudoRandom(g); -# b := PseudoRandom(g); -# Add(gens,yy^a); -# Add(gens,yy^b); -# h:=Group(gens); -# if q = 4 then -# S := StabilizerChain(h); -# if not Size(S) in [60480,3*60480] then continue; fi; -# pos := Position(degrees,1); -# eva := -Value(factors[pos],0*Z(q)); -# ns := NullspaceMat(StripMemory(y)-eva*One(StripMemory(y))); -# return [h, -# VectorSpace(GF(q),[ns[1],ns[1]*StripMemory(a),ns[1]*StripMemory(b)])]; -# fi; -# -# # Now check using ppd-elements: -# for i in [1..10] do -# z:=PseudoRandom(h); -# pol:=CharacteristicPolynomial(GF(q),GF(q),StripMemory(z),1); -# factors:=Factors(pol); -# degrees:=List(factors,Degree); -# if Maximum(degrees)=2 then -# ready2:=true; -# elif Maximum(degrees)=3 then -# ready3:=true; -# fi; -# if ready2 and ready3 then -# break; -# fi; -# od; -# if not (ready2 and ready3) then -# ready2:=false; -# ready3:=false; -# fi; -# until ready2 and ready3; -# -# subgplist:=RECOG.SL_Even_godownone(h,VectorSpace(GF(q),One(g)),q,3); -# fi; -# -# return subgplist; -# end; RECOG.ResetSLstd := function(r) r.left := One(r.a); @@ -713,344 +314,6 @@ end; -# BindGlobal("FunnyProductObjsFamily",NewFamily("FunnyProductObjsFamily")); -# DeclareCategory("IsFunnyProductObject", -# IsPositionalObjectRep and IsMultiplicativeElement and -# IsMultiplicativeElementWithInverse ); -# BindGlobal("FunnyProductObjsType", -# NewType(FunnyProductObjsFamily,IsFunnyProductObject)); -# DeclareOperation("FunnyProductObj",[IsObject,IsObject]); -# -# -# InstallOtherMethod( \*, "for two funny product objects", -# [ IsFunnyProductObject, IsFunnyProductObject ], -# function(a,b) -# return Objectify(FunnyProductObjsType,[a![1]+a![2]*b![1],a![2]*b![2]]); -# end ); -# -# InstallOtherMethod( InverseSameMutability, "for a funny product object", -# [ IsFunnyProductObject ], -# function(a) -# local i; -# i := a![2]^-1; -# return Objectify(FunnyProductObjsType,[-i*a![1],i]); -# end ); -# -# InstallOtherMethod( OneMutable, "for a funny product object", -# [ IsFunnyProductObject ], -# function(a) -# return Objectify(FunnyProductObjsType,[Zero(a![1]),OneMutable(a![2])]); -# end ); -# -# InstallMethod( FunnyProductObj, "for two arbitrary objects", -# [ IsObject, IsObject ], -# function(a,b) -# return Objectify(FunnyProductObjsType,[a,b]); -# end ); -# -# FIXME: unused? but see misc/work/DOWORK. -# Perhaps this was / is meant as a replacement for RECOG.FindStdGens_SL -# in even characteristic. -# But in a quick test based on misc/work/DOWORK, the code there -# seems to be faster. -# RECOG.FindStdGens_SL_EvenChar := function(sld,f) -# # gens of sld must be gens for SL(d,q) in its natural rep with memory -# # This function calls RECOG.SL_Even_constructdata and then extends -# # the basis to a basis of the full row space and returns an slp such -# # that the SL(d,q) standard generators with respect to this basis are -# # expressed by the slp in terms of the original generators of sld. -# local V,b,bas,basi,basit,d,data,diffv,diffw,el,ext,fakegens,gens,i,id, -# lambda,mu,n,notinv,nu,nu2,oldyf,oldyy,p,pos,q,resl2,sl2,sl2gens, -# sl2gensf,sl2genss,sl2stdf,slp,slpsl2std,slptosl2,st,std,stdsl2, -# w,x,xf,y,y2f,y3f,yf,yy,yy2,yy3,yyy,yyy2,yyy3,z,zf,zzz,goodguy; -# -# # Some setup: -# p := Characteristic(f); -# q := Size(f); -# ext := DegreeOverPrimeField(f); -# d := DimensionOfMatrixGroup(sld); -# if not IsObjWithMemory(GeneratorsOfGroup(sld)[1]) then -# sld := GroupWithMemory(sld); -# fi; -# -# # First find an SL2 with the space it acts on: -# Info(InfoRecog,2,"Finding an SL2..."); -# #data := RECOG.SL_Even_constructdata(sld,q); -# repeat -# data := RECOG.SL_FindSL2(sld,f); -# until data <> fail; -# bas := ShallowCopy(BasisVectors(Basis(data[2]))); -# sl2 := data[1]; -# slptosl2 := SLPOfElms(GeneratorsOfGroup(sl2)); -# -# # Now compute the natural SL2 action and run constructive recognition: -# sl2gens := StripMemory(GeneratorsOfGroup(sl2)); -# V := VectorSpace(f,bas); -# b := Basis(V,bas); -# sl2genss := List(sl2gens,x->List(BasisVectors(b),v->Coefficients(b,v*x))); -# for i in sl2genss do -# ConvertToMatrixRep(i,q); -# od; -# Info(InfoRecog,2, -# "Recognising this SL2 constructively in 2 dimensions..."); -# sl2genss := GeneratorsWithMemory(sl2genss); -# if IsEvenInt(q) then -# resl2 := RECOG.RecogniseSL2NaturalEvenChar(Group(sl2genss),f,false); -# else -# resl2 := RECOG.RecogniseSL2NaturalOddCharUsingBSGS(Group(sl2genss),f); -# fi; -# slpsl2std := SLPOfElms(resl2.all); -# bas := resl2.bas * bas; -# # We need the actual transvections: -# slp := SLPOfElms([resl2.s[1],resl2.t[1]]); -# st := ResultOfStraightLineProgram(slp,StripMemory(GeneratorsOfGroup(sl2))); -# -# # Extend basis by something invariant under SL2: -# id := IdentityMat(d,f); -# nu := NullspaceMat(StripMemory(st[1]-id)); -# nu2 := NullspaceMat(StripMemory(st[2]-id)); -# Append(bas,SumIntersectionMat(nu,nu2)[2]); -# ConvertToMatrixRep(bas,q); -# basi := bas^-1; -# basit := TransposedMatMutable(basi); -# -# # Now set up fake generators for keeping track what we do: -# fakegens := ListWithIdenticalEntries(Length(GeneratorsOfGroup(sld)),()); -# fakegens := GeneratorsWithMemory(fakegens); -# sl2gensf := ResultOfStraightLineProgram(slptosl2,fakegens); -# sl2stdf := ResultOfStraightLineProgram(slpsl2std,sl2gensf); -# std := RECOG.InitSLstd(f,d,sl2stdf{[1..ext]},sl2stdf{[ext+1..2*ext]}, -# sl2stdf[2*ext+1],sl2stdf[2*ext+2]); -# -# # workrec := rec( n := 2, slnstdf := sl2stdf, bas := bas, basi := basi, -# # std := std, sld := sld, sldf := fakegens, f := f ); -# # -# #Error("... now go on with alternative going up..."); -# -# Info(InfoRecog,2,"Going up to SL_d again..."); -# for n in [Dimension(data[2])..d-1] do -# if InfoLevel(InfoRecog) >= 3 then Print(n," \c"); fi; -# while true do # will be left by break at the end -# x := PseudoRandom(sld); -# slp := SLPOfElm(x); -# xf := ResultOfStraightLineProgram(slp,fakegens); -# # From now on plain matrices, we have to keep track with the -# # fake ones! -# x := StripMemory(x); -# -# # Find a new basis vector: -# y := st[1]^x; -# notinv := First([1..n],i->bas[i]*y<>bas[i]); -# if notinv = fail then continue; fi; # try next x -# w := bas[notinv]*y-bas[notinv]; -# if ForAll(basit{[n+1..d]},v->IsZero(ScalarProduct(v,w))) then -# continue; # try next x -# fi; -# # Now make it so that w is invariant under SL_n by modifying -# # it by something in the span of bas{[1..n]}: -# for i in [1..n] do -# w := w - bas[i] * ScalarProduct(w,basit[i]); -# od; -# if w*y=w then -# if InfoLevel(InfoRecog) >= 3 then Print("!\c"); fi; -# continue; -# fi; -# -# # w is supposed to become the next basis vector number n+1. -# # So we need to throw away one of bas{[n+1..d]}: -# i := First([n+1..d],i->not IsZero(ScalarProduct(w,basit[i]))); -# Remove(bas,i); -# Add(bas,w,n+1); -# # However, we want that the rest of them bas{[n+2..d]} is invariant -# # under y which we can achieve by adding a multiple of w: -# diffw := w*y-w; -# pos := PositionNonZero(diffw); -# for i in [n+2..d] do -# diffv := bas[i]*y-bas[i]; -# if not IsZero(diffv) then -# bas[i] := bas[i] - (diffv[pos]/diffw[pos]) * w; -# fi; -# od; -# basi := bas^-1; -# basit := TransposedMat(basi); -# -# # Compute the action of y-One(y) on Span(bas{[1..n+1]}) -# yy := EmptyPlist(n+1); -# for i in [1..n+1] do -# Add(yy,(bas[i]*y-bas[i])*basi); -# yy[i] := yy[i]{[1..n+1]}; -# od; -# if q > 2 and IsOne(yy[n+1,n+1]) then -# if InfoLevel(InfoRecog) >= 3 then Print("#\c"); fi; -# continue; -# fi; -# ConvertToMatrixRep(yy,q); -# break; -# od; -# yf := xf^-1*std.s[1]*xf; -# -# # make sure that rows n-1 and n are non-zero: -# std.left := std.One; -# std.right := std.One; -# if IsZero(yy[n-1]) then -# RECOG.DoRowOp_SL(yy,n-1,notinv,std.one,std); -# RECOG.DoColOp_SL(yy,n-1,notinv,-std.one,std); -# fi; -# if IsZero(yy[n]) then -# RECOG.DoRowOp_SL(yy,n,notinv,std.one,std); -# RECOG.DoColOp_SL(yy,n,notinv,-std.one,std); -# fi; -# yf := std.left * yf * std.right; -# -# oldyy := MutableCopyMat(yy); -# oldyf := yf; -# -# if q = 2 then -# # In this case y is already good after cleaning out! -# # (remember that y+One(y) has rank 1 and does not fix bas[notinv]) -# std.left := std.One; -# std.right := std.One; -# for i in [1..n-1] do -# lambda := -yy[i,n+1]/yy[n,n+1]; -# RECOG.DoRowOp_SL(yy,i,n,lambda,std); -# RECOG.DoColOp_SL(yy,i,n,-lambda,std); -# od; -# yf := std.left * yf * std.right; -# z := yy+One(yy); -# zf := yf; -# if not IsZero(z[n,n]) or not IsOne(z[n,n+1]) or -# not IsZero(z[n+1,n+1]) or not IsOne(z[n+1,n]) then -# ErrorNoReturn("How on earth could this happen???"); -# fi; -# else # q > 2 -# while true do # will be left by break when we had success! -# # Note that by construction yy[n,n+1] is not zero! -# yy2 := MutableCopyMat(yy); -# std.left := std.One; -# std.right := std.One; -# # We want to be careful not to kill row n: -# repeat -# lambda := PrimitiveRoot(f)^Random(0,q-1); -# until lambda <> -yy2[n,n+1]/yy2[n-1,n+1]; -# RECOG.DoRowOp_SL(yy2,n,n-1,lambda,std); -# RECOG.DoColOp_SL(yy2,n,n-1,-lambda,std); -# mu := lambda; -# y2f := std.left * yf * std.right; -# -# yy3 := MutableCopyMat(yy); -# std.left := std.One; -# std.right := std.One; -# # We want to be careful not to kill row n: -# repeat -# lambda := PrimitiveRoot(f)^Random(0,q-1); -# until (lambda <> -yy3[n,n+1]/yy3[n-1,n+1]) and -# (lambda <> mu or q = 3); -# # in GF(3) there are not enough values! -# RECOG.DoRowOp_SL(yy3,n,n-1,lambda,std); -# RECOG.DoColOp_SL(yy3,n,n-1,-lambda,std); -# y3f := std.left * yf * std.right; -# -# # We now perform conjugations such that the ys leave -# # bas{[1..n-1]} fixed: -# -# # (remember that y-One(y) has rank 1 and does not fix bas[notinv]) -# std.left := std.One; -# std.right := std.One; -# for i in [1..n-1] do -# lambda := -yy[i,n+1]/yy[n,n+1]; -# RECOG.DoRowOp_SL(yy,i,n,lambda,std); -# RECOG.DoColOp_SL(yy,i,n,-lambda,std); -# od; -# yf := std.left * yf * std.right; -# -# std.left := std.One; -# std.right := std.One; -# for i in [1..n-1] do -# lambda := -yy2[i,n+1]/yy2[n,n+1]; -# RECOG.DoRowOp_SL(yy2,i,n,lambda,std); -# RECOG.DoColOp_SL(yy2,i,n,-lambda,std); -# od; -# y2f := std.left * y2f * std.right; -# -# std.left := std.One; -# std.right := std.One; -# for i in [1..n-1] do -# lambda := -yy3[i,n+1]/yy3[n,n+1]; -# RECOG.DoRowOp_SL(yy3,i,n,lambda,std); -# RECOG.DoColOp_SL(yy3,i,n,-lambda,std); -# od; -# y3f := std.left * y3f * std.right; -# -# gens :=[ExtractSubMatrix(yy,[n,n+1],[n,n+1])+IdentityMat(2,f), -# ExtractSubMatrix(yy2,[n,n+1],[n,n+1])+IdentityMat(2,f), -# ExtractSubMatrix(yy3,[n,n+1],[n,n+1])+IdentityMat(2,f)]; -# if RECOG.IsThisSL2Natural(gens,f) = true then break; fi; -# if InfoLevel(InfoRecog) >= 3 then Print("$\c"); fi; -# yy := MutableCopyMat(oldyy); -# yf := oldyf; -# od; -# -# # Now perform a constructive recognition in the SL2 in the lower -# # right corner: -# gens := GeneratorsWithMemory(gens); -# if IsEvenInt(q) then -# resl2 := RECOG.RecogniseSL2NaturalEvenChar(Group(gens),f,gens[1]); -# else -# resl2 := RECOG.RecogniseSL2NaturalOddCharUsingBSGS(Group(gens),f); -# fi; -# stdsl2 := RECOG.InitSLfake(f,2); -# goodguy := Reversed(IdentityMat(2,f)); -# goodguy[1,2] := - goodguy[1,2]; -# slp := RECOG.ExpressInStd_SL2(resl2.bas*goodguy*resl2.basi,stdsl2); -# el := ResultOfStraightLineProgram(slp,resl2.all); -# slp := SLPOfElm(el); -# -# yy := yy+One(yy); -# yy2 := yy2+One(yy2); -# yy3 := yy3+One(yy3); -# yyy := FunnyProductObj(ExtractSubMatrix(yy,[n,n+1],[1..n-1]), -# ExtractSubMatrix(yy,[n,n+1],[n,n+1])); -# yyy2 := FunnyProductObj(ExtractSubMatrix(yy2,[n,n+1],[1..n-1]), -# ExtractSubMatrix(yy2,[n,n+1],[n,n+1])); -# yyy3 := FunnyProductObj(ExtractSubMatrix(yy3,[n,n+1],[1..n-1]), -# ExtractSubMatrix(yy3,[n,n+1],[n,n+1])); -# zzz := ResultOfStraightLineProgram(slp,[yyy,yyy2,yyy3]); -# z := OneMutable(yy); -# CopySubMatrix(zzz![1],z,[1..2],[n,n+1],[1..n-1],[1..n-1]); -# CopySubMatrix(zzz![2],z,[1..2],[n,n+1],[1..2],[n,n+1]); -# zf := ResultOfStraightLineProgram(slp,[yf,y2f,y3f]); -# fi; -# -# std.left := std.One; -# std.right := std.One; -# # Now we clean out the last row of z: -# for i in [1..n-1] do -# if not IsZero(z[n+1,i]) then -# RECOG.DoColOp_SL(z,n,i,-z[n+1,i],std); -# fi; -# od; -# # Now we clean out the second last row of z: -# for i in [1..n-1] do -# if not IsZero(z[n,i]) then -# RECOG.DoRowOp_SL(z,n,i,-z[n,i],std); -# fi; -# od; -# zf := std.left * zf * std.right; -# -# # Now change the standard generators in the fakes: -# std.a := std.a * zf; -# std.b := std.b * zf; -# std.all[std.ext*2+1] := std.a; -# std.all[std.ext*2+2] := std.b; -# RECOG.ResetSLstd(std); -# -# od; -# if InfoLevel(InfoRecog) >= 3 then Print(".\n"); fi; -# return rec( slpstd := SLPOfElms(std.all), bas := bas, basi := basi ); -# end; - - - RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) # f a finite field, g equal to SL(2,Size(f)), t either an involution # or false. @@ -1238,44 +501,6 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) return res; end; -# RECOG.GuessSL2ElmOrder := function(x,f) -# local facts,i,j,o,p,q,r,s,y,z; -# p := Characteristic(f); -# q := Size(f); -# if IsOne(x) then return 1; -# elif IsOne(x^2) then return 2; -# fi; -# if p > 2 then -# y := x^p; -# if IsOne(y) then return p; -# elif IsOddInt(p) and IsOne(y^2) then -# return 2*p; -# fi; -# fi; -# if IsOne(x^(q-1)) then -# facts := Collected(FactInt(q-1:cheap)[1]); -# s := Product(facts,x->x[1]^x[2]); -# r := (q-1)/s; -# else -# facts := Collected(FactInt(q+1:cheap)[1]); -# s := Product(facts,x->x[1]^x[2]); -# r := (q+1)/s; -# fi; -# y := x^r; -# o := r; -# for i in [1..Length(facts)] do -# p := facts[i]; -# j := p[2]-1; -# while j >= 0 do -# z := y^(s/p[1]^(p[2]-j)); -# if not(IsOne(z)) then break; fi; -# j := j - 1; -# od; -# o := o * p[1]^(j+1); -# od; -# return o; -# end; - RECOG.GuessProjSL2ElmOrder := function(x,f) local facts,i,j,o,p,q,r,s,y,z; p := Characteristic(f); @@ -1432,710 +657,6 @@ RECOG.IsThisSL2Natural := function(gens,f) end; -################################################################################# -################################################################################# -################################################################################# - -# RECOG.MakeSLSituation := function(p,e,n,d) -# local a,q,r; -# q := p^e; -# a := RECOG.MakeSL_StdGens(p,e,n,d).all; -# Append(a,GeneratorsOfGroup(SL(d,q))); -# a := GeneratorsWithMemory(a); -# r := rec( f := GF(q), d := d, n := n, bas := IdentityMat(d,GF(q)), -# basi := IdentityMat(d,GF(q)), sld := Group(a), -# sldf := a, slnstdf := a{[1..2*e+2]}, p := p, ext := e ); -# return r; -# end; -# -# RECOG.MakeSLTest := function(p,e,n,d) -# local a,fake,q,r; -# q := p^e; -# a := RECOG.MakeSL_StdGens(p,e,n,d).all; -# Append(a,GeneratorsOfGroup(SL(d,q))); -# a := GeneratorsWithMemory(a); -# fake := GeneratorsWithMemory(List([1..Length(a)],i->())); -# r := rec( f := GF(q), d := d, n := n, bas := IdentityMat(d,GF(q)), -# basi := IdentityMat(d,GF(q)), sld := Group(a), -# sldf := fake, slnstdf := fake{[1..2*e+2]}, p := p, ext := e ); -# return r; -# end; -# -# RECOG.MakeSp2n := function(n,p,e) -# # n must be even -# local bas,basch,basi,form,g,gens,gg,i; -# g := Sp(2*n,p^e); -# form := InvariantBilinearForm(g).matrix; -# basch := EmptyPlist(2*n); -# for i in [1..n] do -# basch[i] := 2*i-1; -# basch[2*n+1-i] := 2*i; -# od; -# basi := PermutationMat(PermList(basch),2*n,GF(p,e)); -# bas := basi^-1; -# gens := List(GeneratorsOfGroup(g),x->bas*x*basi); -# form := bas * form * basi; -# gg := Group(gens); -# SetSize(gg,Size(g)); -# SetInvariantBilinearForm(gg,rec(matrix := form)); -# return [gg,form]; -# end; -# -# RECOG.MakeSpnTransvection := function(g,type,i,lambda) -# # g must be Sp(2n,q) as made by RECOG.MakeSpn, this defines n -# # type is either "e" or "f", i is in [0..2n-2] -# # Our basis is (b_1, ..., b_{2n}) = (e_1,f_1,...,e_n,f_n) -# # For type="e", this makes the following transvection: -# # x -> x + lambda * (x,e_n + b) * (e_n + b) -# # where b = b_i for i <> 0 and b = f_n for i = 0 -# # For type="f", this makes the following transvection: -# # x -> x + lambda * (x,f_n + b) * (f_n + b) -# # where b = b_i for i <> 0 and b = 0 for i = 0 -# local f,form,id,j,n,o,v; -# n := DimensionOfMatrixGroup(g)/2; -# f := FieldOfMatrixGroup(g); -# o := One(f); -# id := OneMutable(One(g)); -# v := ZeroMutable(id[1]); -# if type = "e" then -# v[2*n-1] := -o; -# else -# v[2*n] := o; -# fi; -# if i <> 0 then -# v[i] := o; -# fi; -# form := InvariantBilinearForm(g).matrix; -# for j in [1..2*n] do -# id[j] := id[j] + (lambda * (id[j]*form)*v) * v; -# od; -# return id; -# end; -# -# RECOG.ComputeGramSymplecticStandardForm := function(vecs) -# # vecs a matrix of vectors of length 2*n interpreted as written in -# # the standard symplectic form below. -# # This computes the Gram matrix of the vectors vecs using the -# # standard symplectic form, which is defined for the standard -# # basis e_1, f_1, ... e_n, f_n to be -# # (e_i|e_j) = 0, (f_i, f_j) = 0, (e_i,f_j) = \delta_{i,j} -# local f,gram,i,j,k,l,n,one,v,zero; -# l := Length(vecs); -# f := BaseDomain(vecs); -# zero := Zero(f); -# one := One(f); -# gram := ZeroMatrix(l,l,vecs); -# n := RowLength(vecs)/2; -# Assert(1,IsInt(n),ErrorNoReturn("RowLength must be even")); -# for i in [1..l] do -# for j in [i+1..l] do -# v := zero; -# for k in [1,3..2*n-1] do -# v := v + vecs[i,k]*vecs[j,k+1] - vecs[i,k+1]*vecs[j,k]; -# od; -# gram[i,j] := v; -# gram[j,i] := -v; -# od; -# od; -# return gram; -# end; -# -# RECOG.FindSymplecticPairBasis := function(vecs) -# local bas,d,dummy,gram,i,j,k,s; -# d := Length(vecs); -# if IsOddInt(d) then -# return [fail,"odd dimension"]; -# fi; -# gram := RECOG.ComputeGramSymplecticStandardForm(vecs); -# bas := IdentityMatrix(d,vecs); -# for i in [1,3..d-1] do -# j := i+1; -# while j <= d do -# if not IsZero(gram[i,j]) then -# s := gram[i,j]^-1; -# MultRowVector(bas[j],s); -# MultRowVector(gram[j],s); -# for k in [1..d] do -# gram[k,j] := gram[k,j]*s; -# od; -# Assert(1,gram = RECOG.ComputeGramSymplecticStandardForm(bas*vecs)); -# # Now exchange vectors i+1 and j: -# if i+1 <> j then -# bas{[i+1,j]} := bas{[j,i+1]}; -# gram{[i+1,j]} := gram{[j,i+1]}; -# for k in [1..d] do -# dummy := gram[k,i+1]; -# gram[k,i+1] := gram[k,j]; -# gram[k,j] := dummy; -# od; -# Assert(1,gram = RECOG.ComputeGramSymplecticStandardForm(bas*vecs)); -# fi; -# break; -# fi; -# j := j + 1; -# od; -# if j > d then return [fail,"degenerate"]; fi; -# # Now i,i+1 is a symplectic pair, clean out the rest: -# for j in [i+2..d] do -# if not IsZero(gram[i,j]) then -# s := gram[i,j]; -# AddRowVector(bas[j],bas[i+1],-s); -# AddRowVector(gram[j],gram[i+1],-s); -# for k in [1..d] do -# gram[k,j] := gram[k,j] - s*gram[k,i+1]; -# od; -# Assert(1,gram = RECOG.ComputeGramSymplecticStandardForm(bas*vecs)); -# fi; -# if not IsZero(gram[i+1,j]) then -# s := gram[i+1,j]; -# AddRowVector(bas[j],bas[i],s); -# AddRowVector(gram[j],gram[i],s); -# for k in [1..d] do -# gram[k,j] := gram[k,j] + s*gram[k,i]; -# od; -# Assert(1,gram = RECOG.ComputeGramSymplecticStandardForm(bas*vecs)); -# fi; -# od; -# # Now all further vectors are perpendicular to vecs i and i+1 -# od; -# return bas; -# end; -# -# RECOG.SetupSpExperiment := function(n,d,f) -# local em,formg,formh,g,h,ncycle; -# Assert(1,n < d); -# g := RECOG.MakeSp2n(d,Characteristic(f),DegreeOverPrimeField(f)); -# formg := g[2]; -# g := g[1]; -# h := RECOG.MakeSp2n(n,Characteristic(f),DegreeOverPrimeField(f)); -# formh := h[2]; -# h := h[1]; -# em := GroupHomomorphismByFunction(g,h, -# function(x) -# local i; -# i := IdentityMatrix(2*d,formg); -# CopySubMatrix(x,i,[1..2*n],[1..2*n],[1..2*n],[1..2*n]); -# return i; -# end); -# ncycle := PermutationMat(PermList(Concatenation([3,4..2*n],[1,2])),2*d,f); -# return rec(g := g,formg := formg,h := h,formh := formh,em := em, -# ncycle := ncycle, p := Characteristic(f), -# ext := DegreeOverPrimeField(f), d := d, n := n, f := f, -# q := Size(f), id := IdentityMat(2*d,f)); -# end; -# -# # Standard generators of Sp(2n,q) are given by a record with: -# # n n -# # q q=p^e -# # p p -# # ext e -# # f GF(q) -# # can CanonicalBasis(GF(q)) -# # s the element [[0,1],[-1,0]] on -# # delta the element [[zeta,0],[0,zeta^-1]] on -# # v the double-n-cycle (e_1,e_2,...,e_n)(f_1,f_2,...,f_n) -# # ten transvections t_{e_n} -# # a list of ext elements -# # tfn transvections t_{f_n} -# # a list of ext elements -# # tfnei transvections t_{f_n+e_i} (1 <= i <= n-1) -# # each entry is a list of ext elements -# # tfnfi transvections t_{f_n+e_i} (1 <= i <= n-1) -# # each entry is a list of ext elements -# -# RECOG.MakeSpnTfn := function(n,d,f,lambda) -# local t; -# t := IdentityMat(2*d,f); -# t[2*n-1,2*n] := lambda; -# return t; -# end; -# -# RECOG.MakeSpnTfnei := function(n,d,f,i,lambda) -# local t; -# t := IdentityMat(2*d,f); -# t[2*i,2*i-1] := -lambda; -# t[2*i,2*n] := -lambda; -# t[2*n-1,2*i-1] := lambda; -# t[2*n-1,2*n] := lambda; -# return t; -# end; -# -# RECOG.MakeSpnTfnfi := function(n,d,f,i,lambda) -# local t; -# t := IdentityMat(2*d,f); -# t[2*i-1,2*i] := lambda; -# t[2*i-1,2*n] := lambda; -# t[2*n-1,2*i] := lambda; -# t[2*n-1,2*n] := lambda; -# return t; -# end; -# -# RECOG.MakeSpnTen := function(n,d,f,lambda) -# local t; -# t := IdentityMat(2*d,f); -# t[2*n,2*n-1] := -lambda; -# return t; -# end; -# -# RECOG.MakeSpnTenei := function(n,d,f,i,lambda) -# local t; -# t := IdentityMat(2*d,f); -# t[2*i,2*i-1] := -lambda; -# t[2*i,2*n-1] := -lambda; -# t[2*n,2*i-1] := -lambda; -# t[2*n,2*n-1] := -lambda; -# return t; -# end; -# -# RECOG.MakeSpnTenfi := function(n,d,f,i,lambda) -# local t; -# t := IdentityMat(2*d,f); -# t[2*i-1,2*i] := lambda; -# t[2*i-1,2*n-1] := lambda; -# t[2*n,2*i] := -lambda; -# t[2*n,2*n-1] := -lambda; -# return t; -# end; -# -# RECOG.MakeSp_StdGens := function(p,ext,n,d) -# local f,g,id,l,one,q,res,zero,zeta; -# q := p^ext; -# f := GF(p,ext); -# res := rec( q := q, p := p, ext := ext, f := f, n := n, -# can := CanonicalBasis(f) ); -# zero := Zero(f); -# one := One(f); -# zeta := PrimitiveRoot(f); -# id := IdentityMat(2*d,f); -# res.s := MutableCopyMat(id); -# res.s[2*n-1,2*n-1] := zero; -# res.s[2*n-1,2*n] := one; -# res.s[2*n,2*n-1] := -one; -# res.s[2*n,2*n] := zero; -# res.delta := MutableCopyMat(id); -# res.delta[2*n-1,2*n-1] := zeta; -# res.delta[2*n,2*n] := zeta^-1; -# l := Concatenation([3..2*n],[1,2]); -# res.v := PermutationMat(PermList(l),2*d,f); -# res.ten := List([0..ext-1], -# k->RECOG.MakeSpnTen(n,d,f,zeta^k)); -# res.tfn := List([0..ext-1], -# k->RECOG.MakeSpnTfn(n,d,f,zeta^k)); -# res.tfnei := List([1..n-1],i-> -# List([0..ext-1], -# k->RECOG.MakeSpnTfnei(n,d,f,i,zeta^k))); -# res.tfnfi := List([1..n-1],i-> -# List([0..ext-1], -# k->RECOG.MakeSpnTfnfi(n,d,f,i,zeta^k))); -# res.all := Concatenation([res.s,res.delta,res.v], -# res.ten,res.tfn, -# Concatenation(res.tfnei), -# Concatenation(res.tfnfi)); -# return res; -# end; -# -# RECOG.MakeSp_FakeGens := function(p,ext,n) -# local count,f,fake,i,q,res; -# q := p^ext; -# f := GF(p,ext); -# res := rec( q := q, p := p, ext := ext, f := f, n := n, -# can := CanonicalBasis(f) ); -# fake := GeneratorsWithMemory( -# ListWithIdenticalEntries(3+(2*n+2)*ext,1)); -# res.s := fake[1]; -# res.delta := fake[2]; -# res.v := fake[3]; -# count := 3; -# res.tfn := fake{[count+1..count+ext]}; -# count := count + ext; -# res.ten := fake{[count+1..count+ext]}; -# count := count + ext; -# res.tfnei := EmptyPlist(n-1); -# for i in [1..n-1] do -# Add(res.tfnei,fake{[count+1..count+ext]}); -# count := count + ext; -# od; -# res.tfnfi := EmptyPlist(n-1); -# for i in [1..n-1] do -# Add(res.tfnfi,fake{[count+1..count+ext]}); -# count := count + ext; -# od; -# res.all := fake; -# return res; -# end; -# -# RECOG.SpMakeImage_en := -# function(v,s,M,usencycle) -# # v is a vector over F_q of length at least 2n and v[2n-1]=1. -# # s is a set of standard generators of Sp(2n,q) (see above). -# # This func. makes an element t of Sp(2n,q) that maps v to e_n and fixes -# # f_n. The result t is expressed as a product in the standard generators -# # of Sp(2n,q) in s (see above). If M is not equal to fail then it must -# # be a matrix of mutable vectors over F_q of at least length 2n and it is -# # modified as if it were multiplied by t. This means that if M is -# # a mutable identity matrix of size at least 2n x 2n, then it will -# # contain the matrix of t after the operation in its upper left corner. -# # usencycle must be either true or false. If it is set to true, -# # the n-cycle amongst the standard generators is used resulting -# # in shorter products. If usencycle is false, then the n-cycle is -# # not used, note that this does not work for q=2. -# # The function returns t and changes M if not equal to fail. -# local Morig,coeff,ei,ext,fI,i,k,l,n,one,sc,sc2,si,t,vorig,zero,zeta; -# -# # We want to put together an element that maps v to e_n and fixes f_n: -# # At the same time we map M under the result whilst building it up. -# # We start with (v,M) and apply transvections... -# t := s.tfn[1]^0; # start here -# n := s.n; -# zero := Zero(s.f); -# one := One(s.f); -# zeta := PrimitiveRoot(s.f); -# ext := s.ext; -# Assert(1,IsOne(v[2*n-1])); -# vorig := ShallowCopy(v); -# if M <> fail then Morig := MutableCopyMat(M); fi; -# for i in [1..s.n-1] do -# ei := 2*i-1; # these are the coordinates to modify -# fI := 2*i; -# coeff := one; -# if IsZero(one+v[ei]) and IsZero(one-v[fI]) then -# if usencycle then -# t := t * s.tfn[1]^(s.v^i); -# v[fI] := v[fI] + v[ei]; -# if M <> fail then -# for l in [1..Length(M)] do -# M[l,fI] := M[l,fI] + M[l,ei]; -# od; -# fi; -# else -# if Size(s.f) = 2 then -# ErrorNoReturn("This does not work for GF(2)."); -# fi; -# t := t * s.delta; -# v[2*n-1] := v[2*n-1] * zeta; -# v[2*n] := v[2*n] * zeta^-1; -# if M <> fail then -# for l in [1..Length(M)] do -# M[l,2*n-1] := M[l,2*n-1] * zeta; -# M[l,2*n] := M[l,2*n] * zeta^-1; -# od; -# fi; -# coeff := zeta; -# fi; -# Assert(1,v=vorig*t and (M = fail or Morig*t=M),ErrorNoReturn("Hallo 0")); -# fi; -# if IsZero(v[ei]) or not IsZero(coeff-v[fI]) then -# # The first easy case: -# # First kill v[ei] if need be: -# if not IsZero(v[ei]) then -# sc := -v[ei]/(coeff-v[fI]); -# si := IntVecFFE(Coefficients(s.can,sc)); -# for k in [1..ext] do -# t := t * s.tfnei[i,k]^si[k]; -# od; -# v[2*n] := v[2*n] - v[ei]; -# v[ei] := zero; -# if M <> fail then -# for l in [1..Length(M)] do -# sc2 := sc * (M[l,2*n-1]-M[l,fI]); -# M[l,ei] := M[l,ei] + sc2; -# M[l,2*n] := M[l,2*n] + sc2; -# od; -# fi; -# Assert(1,v=vorig*t and (M = fail or Morig*t=M),ErrorNoReturn("Hallo 1")); -# fi; -# # Now kill v[fI] if need be: -# if not IsZero(v[fI]) then -# sc := -v[fI]/coeff; -# si := IntVecFFE(Coefficients(s.can,sc)); -# for k in [1..ext] do -# t := t * s.tfnfi[i,k]^si[k]; -# od; -# v[2*n] := v[2*n] - v[fI]; -# v[fI] := zero; -# if M <> fail then -# for l in [1..Length(M)] do -# sc2 := sc * (M[l,2*n-1]+M[l,ei]); -# M[l,fI] := M[l,fI] + sc2; -# M[l,2*n] := M[l,2*n] + sc2; -# od; -# fi; -# Assert(1,v=vorig*t and (M = fail or Morig*t=M),ErrorNoReturn("Hallo 2")); -# fi; -# elif not IsZero(one+v[ei]) then -# # The second easy case: -# # Here v[fI] = 1 and v[ei] <> 0: -# # First kill v[fI]: -# sc := -v[fI]/(coeff+v[ei]); -# si := IntVecFFE(Coefficients(s.can,sc)); -# for k in [1..ext] do -# t := t * s.tfnfi[i,k]^si[k]; -# od; -# v[2*n] := v[2*n] - v[fI]; -# v[fI] := zero; -# if M <> fail then -# for l in [1..Length(M)] do -# sc2 := sc * (M[l,2*n-1]+M[l,ei]); -# M[l,fI] := M[l,fI] + sc2; -# M[l,2*n] := M[l,2*n] + sc2; -# od; -# fi; -# Assert(1,v=vorig*t and (M = fail or Morig*t=M),ErrorNoReturn("Hallo 3")); -# # Now kill v[ei] if need be: -# sc := -v[ei]/coeff; -# si := IntVecFFE(Coefficients(s.can,sc)); -# for k in [1..ext] do -# t := t * s.tfnei[i,k]^si[k]; -# od; -# v[2*n] := v[2*n] - v[ei]; -# v[ei] := zero; -# if M <> fail then -# for l in [1..Length(M)] do -# sc2 := sc * (M[l,2*n-1]-M[l,fI]); -# M[l,ei] := M[l,ei] + sc2; -# M[l,2*n] := M[l,2*n] + sc2; -# od; -# fi; -# Assert(1,v=vorig*t and (M = fail or Morig*t=M),ErrorNoReturn("Hallo 4")); -# fi; -# if coeff = zeta then -# # Fix the e_n coefficient again: -# t := t * s.delta^-1; -# v[2*n-1] := v[2*n-1] * zeta^-1; -# v[2*n] := v[2*n] * zeta; -# if M <> fail then -# for l in [1..Length(M)] do -# M[l,2*n-1] := M[l,2*n-1] * zeta^-1; -# M[l,2*n] := M[l,2*n] * zeta; -# od; -# fi; -# Assert(1,v=vorig*t and (M = fail or Morig*t=M),ErrorNoReturn("Hallo 5")); -# fi; -# od; -# # Finally arrange fn component to 0: -# if not IsZero(v[2*n]) then -# sc := -v[2*n]; -# si := IntVecFFE(Coefficients(s.can,sc)); -# for k in [1..ext] do -# t := t * s.tfn[k]^si[k]; -# od; -# v[2*n] := zero; -# if M <> fail then -# for l in [1..Length(M)] do -# M[l,2*n] := M[l,2*n] + M[l,2*n-1] * sc; -# od; -# fi; -# Assert(1,v=vorig*t and (M = fail or Morig*t=M),ErrorNoReturn("Hallo 6")); -# fi; -# return t; -# end; -# -# RECOG.SpMakeImage_enfn := function(v,w,s,usencycle) -# local t,ttt; -# # This produces an element of Sp(2n,q) mapping v to e_n and w to f_n -# # as a product of the standard generators. Obviously, the pair (v,w) -# # must be a symplectic pair, furthermore, the e_n-component of v -# # must be one. -# # This function destroys v and w, it uses the ncycle if and only if -# # usencycle is true. -# t := RECOG.SpMakeImage_en(v,s,[w],usencycle); -# # We have achieved that t maps v to e_n and w is changed according -# # to the action to t on it. -# # Now we want to find a tt that maps w to f_n and fixes e_n, since -# # we have s.s mapping e_n to f_n and f_n to -e_n, we can use a ttt -# # mapping w*s.s^-1 to e_n and fixing f_n, and set tt := s.s^-1 * ttt * s.s. -# # Recall that (e_n,w) is a symplectic pair since (int[1],int[2]) was. -# Assert(1,IsOne(w[2*s.n])); -# # Compute w*s.s^-1: -# w[2*s.n] := -w[2*s.n-1]; -# w[2*s.n-1] := One(s.f); -# ttt := RECOG.SpMakeImage_en(w,s,fail,usencycle); -# t := t * s.s^-1 * ttt * s.s; -# return t; -# end; -# -# RECOG.DoSpExperiment := function(r) -# local Vn,Vnc,bas,bigbas,bigbasi,c,c1,fixc,i,int,int2,int3,perp,s,sum,suminter,suminter2,suminter3,t,tt,ttt,u,v,vecs,w,zeta; -# c1 := PseudoRandom(r.g); -# c := r.ncycle^c1; -# Vn := ExtractSubMatrix(r.id,[1..2*r.n],[1..2*r.d]); -# Vnc := ExtractSubMatrix(c,[1..2*r.n],[1..2*r.d]); -# suminter := SumIntersectionMat(Vn,Vnc); -# sum := suminter[1]; -# ConvertToMatrixRep(sum,r.f); -# vecs := suminter[2]; -# ConvertToMatrixRep(vecs,r.f); -# if Length(vecs) <> 2 then -# return ["Vn \cap Vnc not 2-dim",c1]; -# fi; -# if RankMat(ExtractSubMatrix(vecs,[1,2],[2*r.n-1,2*r.n])) < 2 then -# return ["Vn \cap Vnc cannot replace ",c1]; -# fi; -# if IsZero(vecs[1,2*r.n-1]) then -# vecs[1] := vecs[1]+vecs[2]; -# fi; -# MultRowVector(vecs[1],vecs[1,2*r.n-1]^-1); -# bas := RECOG.FindSymplecticPairBasis(vecs); -# if bas[1] = fail then -# return ["Vn \cap Vnc degenerate",c1]; -# fi; -# int := bas*vecs; -# perp := ExtractSubMatrix(r.id,[2*r.n+1..2*r.d],[1..2*r.d]); -# suminter2 := SumIntersectionMat(sum,perp); -# vecs := suminter2[2]; -# ConvertToMatrixRep(vecs,r.f); -# if Length(vecs) <> 2*r.n-2 then -# return ["(Vn + Vnc) \cap Vnperp not 2*n-2-dim",c1]; -# fi; -# bas := RECOG.FindSymplecticPairBasis(vecs); -# if bas[1] = fail then -# return ["(Vn + Vnc) \cap Vnperp degenerate",c1]; -# fi; -# int2 := bas * vecs; -# if 2*r.n-1 < r.d then -# fixc := NullspaceMat(c-One(c)); -# suminter3 := SumIntersectionMat(fixc,perp); -# vecs := suminter3[2]; -# ConvertToMatrixRep(vecs); -# if Length(vecs) <> 2*r.d - 4*r.n + 2 then -# return ["Fixc \cap Vnperp not 2*d-4*n+2-dim",c1]; -# fi; -# bas := RECOG.FindSymplecticPairBasis(vecs); -# if bas[1] = fail then -# return ["Fixc \cap Vnperp degenerate",c1]; -# fi; -# int3 := bas*vecs; -# else -# int3 := []; -# fi; -# # Now we find a product of transvections mapping -# # int[1] to e_n and fixing f_n, we keep track where int[2] is going. -# s := RECOG.MakeSp_StdGens(Characteristic(r.f), -# DegreeOverPrimeField(r.f),r.n,r.d); -# zeta := PrimitiveRoot(r.f); -# v := ShallowCopy(int[1]); -# w := ShallowCopy(int[2]); -# t := RECOG.SpMakeImage_enfn(v,w,s,true); -# # We have achieved that t^-1*c*t fixes e_n and f_n. -# c := t^-1 * c * t; -# -# # Now we need to find the new nice basis vectors n_1, ..., n_2*n-2 -# # they ought to be a symplectic basis when mapped with c and then -# # truncated to coordinates 2*n+1..2*d -# vecs := ExtractSubMatrix(c,[1..2*r.n-2],[2*r.n+1..2*r.d]); -# bas := RECOG.FindSymplecticPairBasis(vecs); -# vecs := bas * ExtractSubMatrix(c,[1..2*r.n-2],[1..2*r.d]); -# # We shall clean out the first 2*n entries of these vectors later on, -# # however, for the time being we keep them for cleaning purposes: -# u := EmptyPlist(2*r.n-2); -# for i in [1..2*r.n-2] do -# v := ZeroMutable(v); -# CopySubVector(bas[i],v,[1..2*r.n-2],[1..2*r.n-2]); -# v[2*r.n-1] := One(r.f); -# ttt := RECOG.SpMakeImage_en(v,s,fail,true); -# # Now clean it in the upper right and lower left corner: -# w := ZeroMutable(w); -# CopySubVector(vecs[i],w,[1..2*r.n-2],[1..2*r.n-2]); -# w[2*r.n-1] := One(r.f); -# tt := RECOG.SpMakeImage_en(w,s,fail,true); -# u[i] := List(s.ten,t->t^(ttt^-1*c*tt)); -# od; -# CopySubMatrix(ZeroMatrix(2*r.n-2,2*r.n,vecs),vecs, -# [1..2*r.n-2],[1..2*r.n-2],[1..2*r.n],[1..2*r.n]); -# bigbas := Concatenation(ExtractSubMatrix(r.id,[1..2*r.n],[1..2*r.d]), -# vecs, -# int3); -# bigbasi := bigbas^-1; -# if bigbasi = fail then -# return ["bigbas is singular",c1]; -# fi; -# return rec( bigbas := bigbas, bigbasi := bigbasi, c := c, c1 := c1, -# t := t, std := s, u := u ); -# end; -# -# RECOG.FindOrder3Element := function(g) -# local a,f,fa,m,o,p,pp,ppp,q,x,y; -# f := FieldOfMatrixGroup(g); -# q := Size(f); -# p := Characteristic(f); -# while true do -# Print(":\c"); -# x := PseudoRandom(g); -# m := MinimalPolynomial(x); -# fa := Collected(Factors(PolynomialRing(f),m)); -# o := Lcm(List(fa,p->q^Degree(p[1])-1)); -# pp := Maximum(List(fa,x->x[2])); -# ppp := p; -# while ppp < pp do -# ppp := ppp * p; -# od; -# while true do -# Print("-\c"); -# a := QuotientRemainder(Integers,o,3); -# if a[2] <> 0 then break; fi; -# o := a[1]; -# od; -# x := x^(o*ppp); -# if IsOne(x) then continue; fi; -# while true do -# Print("+\c"); -# y := x^3; -# if IsOne(y) then break; fi; -# x := y; -# od; -# break; -# od; -# Print("!\n"); -# # Now x is an element of Order 3 -# return x; -# end; -# -# RECOG.MovedSpace := function(g) -# local gens,sp; -# gens := GeneratorsOfGroup(g); -# sp := SemiEchelonMat(Concatenation(List(gens,x->x-One(x)))).vectors; -# return sp; -# end; -# -# RECOG.FixedSpace := function(g) -# local gens,i,inter,sp; -# gens := GeneratorsOfGroup(g); -# sp := List(gens,x->NullspaceMat(x-One(x))); -# if Length(sp) = 1 then -# ConvertToMatrixRep(sp[1],FieldOfMatrixGroup(g)); -# return sp[1]; -# fi; -# inter := SumIntersectionMat(sp[1],sp[2])[2]; -# for i in [3..Length(sp)] do -# inter := SumIntersectionMat(inter,sp[i])[2]; -# od; -# ConvertToMatrixRep(inter,FieldOfMatrixGroup(g)); -# return inter; -# end; -# -# RECOG.guck := function ( w ) -# local i; -# for i in w.slnstdf do -# Display( w.bas * i * w.basi ); -# od; -# if IsBound( w.transh ) then -# for i in [ 1 .. Length( w.transh ) ] do -# Print( i, "\n" ); -# if IsBound(w.transh[i]) then -# Display( w.bas * w.transh[i] * w.basi ); -# fi; -# od; -# fi; -# if IsBound( w.transv ) then -# for i in [ 1 .. Length( w.transv ) ] do -# Print( i, "\n" ); -# if IsBound(w.transv[i]) then -# Display( w.bas * w.transv[i] * w.basi ); -# fi; -# od; -# fi; -# return; -# end; # Now the code for writing SLPs: From 3847eb73e92b8d175259b1b3b60289c473cc6c2d Mon Sep 17 00:00:00 2001 From: Max Horn Date: Thu, 20 Jun 2024 12:46:30 +0200 Subject: [PATCH 12/13] Restore some code --- .../constructive_recognition/utils/achieve.gi | 70 ++++++++++++------- 1 file changed, 45 insertions(+), 25 deletions(-) diff --git a/gap/projective/constructive_recognition/utils/achieve.gi b/gap/projective/constructive_recognition/utils/achieve.gi index 6a48d7c3..c0075724 100644 --- a/gap/projective/constructive_recognition/utils/achieve.gi +++ b/gap/projective/constructive_recognition/utils/achieve.gi @@ -364,7 +364,7 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) ch := Factors(CharacteristicPolynomial(f,f,t,1)); if Length(ch) <> 2 or ch[1] <> ch[2] then - Error("matrix is not triagonalizable - this should never happen!"); + ErrorNoReturn("matrix is not triagonalizable - this should never happen!"); fi; one := OneMutable(t); @@ -401,12 +401,12 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) xm := o[j]; j := j + 1; c := Comm(tm,xm); - until not(IsOne(c^2)); + until not IsOne(c^2); xm := xm * c^(((q-1)*(q+1)-1)/2); x := StripMemory(xm); xb := bas*x*bas^-1; co := Coefficients(can,xb[2,1]); - until not(IsContainedInSpan(mb,co)); + until not IsContainedInSpan(mb,co); CloseMutableBasis(mb,co); Add(tt,x); Add(ttm,xm); @@ -440,20 +440,20 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) xm := o[j]; j := j + 1; x := MutableCopyMat(bas*StripMemory(xm)*bas^-1); - until not(IsZero(x[1,2])); + until not IsZero(x[1,2]); - if not(IsOne(x[2,2])) then + if not IsOne(x[2,2]) then el := (One(f)-x[2,2])/x[1,2]; co := Coefficients(can,el) * mati; for i in [1..Length(co)] do - if not(IsZero(co[i])) then + if not IsZero(co[i]) then xm := ttm[i] * xm; fi; od; x[2] := x[2] + x[1] * el; if x <> bas*StripMemory(xm)*bas^-1 then # FIXME: sometimes triggered by RecognizeGroup(GL(2,16)); - Error("!!!"); + ErrorNoReturn("!!!"); fi; fi; # now x[2,2] is equal to One(f) @@ -465,7 +465,7 @@ RECOG.RecogniseSL2NaturalEvenChar := function(g,f,torig) el := x[2,1]; co2 := Coefficients(can,el) * mati; for i in [1..Length(co2)] do - if not(IsZero(co2[i])) then + if not IsZero(co2[i]) then xm := xm * ttm[i]; fi; od; @@ -510,7 +510,9 @@ RECOG.GuessProjSL2ElmOrder := function(x,f) fi; if p > 2 then y := x^p; - if IsOneProjective(y) then return p; fi; + if IsOneProjective(y) then + return p; + fi; fi; if IsOneProjective(x^(q-1)) then facts := Collected(FactInt(q-1:cheap)[1]); @@ -528,7 +530,7 @@ RECOG.GuessProjSL2ElmOrder := function(x,f) j := p[2]-1; while j >= 0 do z := y^(s/p[1]^(p[2]-j)); - if not(IsOneProjective(z)) then break; fi; + if not IsOneProjective(z) then break; fi; j := j - 1; od; o := o * p[1]^(j+1); @@ -547,24 +549,34 @@ RECOG.IsThisSL2Natural := function(gens,f) CheckElm := function(x) local o; o := RECOG.GuessProjSL2ElmOrder(x,f); - if o in [1,2] then return false; fi; + if o in [1,2] then + return false; + fi; if o > 5 then if notA5 = false then Info(InfoRecog,4,"SL2: Group is not A5"); fi; notA5 := true; - if seenqp1 and seenqm1 then return true; fi; + if seenqp1 and seenqm1 then + return true; + fi; + fi; + if o = p or o <= 5 then + return false; fi; - if o = p or o <= 5 then return false; fi; if (q+1) mod o = 0 then - if not(seenqp1) then + if not seenqp1 then Info(InfoRecog,4,"SL2: Found element of order dividing q+1."); seenqp1 := true; - if seenqm1 and notA5 then return true; fi; + if seenqm1 and notA5 then + return true; + fi; fi; else - if not(seenqm1) then + if not seenqm1 then Info(InfoRecog,4,"SL2: Found element of order dividing q-1."); seenqm1 := true; - if seenqp1 and notA5 then return true; fi; + if seenqp1 and notA5 then + return true; + fi; fi; fi; return false; @@ -594,7 +606,9 @@ RECOG.IsThisSL2Natural := function(gens,f) notA5 := false; for i in [1..Length(gens)] do - if CheckElm(gens[i]) then return true; fi; + if CheckElm(gens[i]) then + return true; + fi; od; CheckElm(gens[1]*gens[2]); if Length(gens) >= 3 then @@ -610,7 +624,9 @@ RECOG.IsThisSL2Natural := function(gens,f) for i in [1..l-1] do for j in [i+1..l] do x := Comm(gens[i],gens[j]); - if CheckElm(x) then return true; fi; + if CheckElm(x) then + return true; + fi; Add(coms,x); od; od; @@ -620,7 +636,9 @@ RECOG.IsThisSL2Natural := function(gens,f) a := RECOG.RandomSubproduct(gens,rec()); b := RECOG.RandomSubproduct(gens,rec()); x := Comm(a,b); - if CheckElm(x) then return true; fi; + if CheckElm(x) then + return true; + fi; Add(coms,x); od; fi; @@ -631,7 +649,9 @@ RECOG.IsThisSL2Natural := function(gens,f) Info(InfoRecog,4,"SL2: Computing normal closure..."); clos := FastNormalClosure(gens,coms,5); for i in [Length(coms)+1..Length(clos)] do - if CheckElm(clos[i]) then return true; fi; + if CheckElm(clos[i]) then + return true; + fi; od; if ForAll(clos{[Length(coms)+1..Length(clos)]}, c->RECOG.IsScalarMat(c) <> false) then @@ -851,10 +871,10 @@ function(ri, g) Info(InfoRecog,2,"ClassicalNatural: this is PSL_2!"); if IsEvenInt(q) then std := RECOG.RecogniseSL2NaturalEvenChar(gm,f,false); - ri!.comment := "_PSL2Even"; + ri!.comment := "PSL2Even"; else std := RECOG.RecogniseSL2NaturalOddCharUsingBSGS(gm,f); - ri!.comment := "_PSL2Odd"; + ri!.comment := "PSL2Odd"; fi; Setslptonice(ri,SLPOfElms(std.all)); ri!.nicebas := std.bas; @@ -875,7 +895,7 @@ function(ri, g) # FIXME: Note d=3 currently has a problem in the SL2-finder. Info(InfoRecog,2,"Classical natural: SL(",d,",",q,"): small ", "case, handing over to Schreier-Sims."); - ri!.comment := Concatenation("_SL(",String(d),",",String(q),")", + ri!.comment := Concatenation("SL(",String(d),",",String(q),")", "_StabilizerChain"); return FindHomMethodsProjective.StabilizerChainProj(ri,g); fi; @@ -890,7 +910,7 @@ function(ri, g) x->std.basi*x*std.bas)); ri!.fakegens := RECOG.InitSLfake(f,d); ri!.fakegens.count := 0; - ri!.comment := "_PSLd"; + ri!.comment := "PSLd"; ri!.gcd := gcd; SetFilterObj(ri,IsLeaf); SetSize(ri,Product([0..d-1],i->(q^d-q^i))/((q-1)*gcd.gcd)); From 36bbef13ee9313c6e23e9c02771752a6c4ee1a24 Mon Sep 17 00:00:00 2001 From: ThomasBreuer Date: Thu, 20 Jun 2024 12:46:10 +0200 Subject: [PATCH 13/13] fix more places where the current directory is needed If one wants to build the documentation with version 1.4.3 of PackageManager then also the calls of `OutputTextFile` require paths relative to `DirectoryCurrent()`. --- regen_doc.g | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/regen_doc.g b/regen_doc.g index a1a2630e..8ab3c633 100644 --- a/regen_doc.g +++ b/regen_doc.g @@ -69,9 +69,10 @@ ListOfUnusedMethods := function() end; GenerateMethodsTableXML := function(shortname, desc, db) -local xmlfile, meth; +local dir, xmlfile, meth; - xmlfile := Concatenation("doc/_methods_", shortname, "_table.xml"); + dir := Filename(DirectoryCurrent(), ""); + xmlfile := Concatenation(dir, "doc/_methods_", shortname, "_table.xml"); xmlfile := OutputTextFile(xmlfile, false); SetPrintFormattingStatus(xmlfile, false); @@ -99,7 +100,7 @@ end; GenerateUnusedMethodsTableXML := function() local xmlfile, meth; - xmlfile := "doc/_methods_unused_table.xml"; + xmlfile := Filename(DirectoryCurrent(), "doc/_methods_unused_table.xml"); xmlfile := OutputTextFile(xmlfile, false); SetPrintFormattingStatus(xmlfile, false); @@ -123,9 +124,10 @@ GenerateUnusedMethodsTableXML := function() end; GenerateMethodsListXML := function(shortname, db) - local xmlfile, dbsWhichUseMethod, nrDbsWhichUseMethod, s, meth; + local dir, xmlfile, dbsWhichUseMethod, nrDbsWhichUseMethod, s, meth; - xmlfile := Concatenation("doc/_methods_", shortname, "_list.xml"); + dir := Filename(DirectoryCurrent(), ""); + xmlfile := Concatenation(dir, "doc/_methods_", shortname, "_list.xml"); xmlfile := OutputTextFile(xmlfile, false); SetPrintFormattingStatus(xmlfile, false);