load GVoronoi;

TestConjugacy:=function(F1,F2)
 //test conjugacy of automorphism groups of two minimal
 //classes (from the list REPRESENTATIVES)

 G1:=StabilizerOfMinimalClass(F1);
 G2:=StabilizerOfMinimalClass(F2);

 if not IdentifyGroup(G1) eq IdentifyGroup(G2) then
  return false;
 end if;

 L1:=[ReynoldsProjection(G1,x) : x in BasHermNorm];
 L1:=matbas3(L1);
 V:=sub<KMatrixSpace(Rationals(),2*n,2*n)|L1>;
 BG1:=Basis(V);
 BG1:=matbas(BG1);
 dim1:=Dimension(V);

 L2:=[ReynoldsProjection(G2,x) : x in BasHermNorm];
 L2:=matbas3(L2);
 V:=sub<KMatrixSpace(Rationals(),2*n,2*n)|L2>;
 BG2:=Basis(V);
 BG2:=matbas(BG2);
 dim2:=Dimension(V);

 FG1:=ReynoldsProjection(G1,F1);
 FG1:=(1/hermitianmin(FG1))*FG1;
 FG2:=ReynoldsProjection(G2,F2);
 FG2:=(1/hermitianmin(FG2))*FG2;

 if dim1 eq dim2 and dim1 eq 1 then
  if not isisom(FG1,FG2) then
   return false;
  end if;
 end if;
 
 if prank(FG1) eq n and prank(FG2) eq n then
  //In this case the Voronoi domain has only dead ends
  if Gprank(G1,FG1) eq dim1 and Gprank(G2,FG2) eq dim2 then
   if not isisom(FG1,FG2) then
    return false;
   end if;
  end if;
 end if;
 return "Unknown.";
end function;

IsMaximalFinite:=function(A)
 //tests maximal finiteness of Stabilizer for Representatives
 if not HasOneDimensionalIntersection(A) then
  return false;
 end if;

 G1:=StabilizerOfMinimalClass(A);
 FG:=ReynoldsProjection(G1,A);
 FG:=(1/hermitianmin(FG))*FG;
 G2:=StabilizerOfMinimalClass(FG);
 if not G1 subset G2 then
  return "Error in IsMaximalFinite: G1 is not contained in G2.";
 end if;
 if #G1 lt #G2 then
  return "G2 is a proper supergroup.";
 end if;
 if IdentifyGroup(G1) eq <2,1> then
  return false;
 end if;

 L1:=[ReynoldsProjection(G1,x) : x in BasHermNorm];
 L1:=matbas3(L1);
 V:=sub<KMatrixSpace(Rationals(),2*n,2*n)|L1>;
 BG1:=Basis(V);
 BG1:=matbas(BG1);
 dim1:=Dimension(V);

 if dim1 eq 1 then
  return true;
 end if;

 if Gprank(G1,FG) lt dim1 then
  //F not G-perfect
  return "Error in IsMaximalFinite: F is not G-perfect.";
 end if;

 if prank(FG) eq n then
  //In this case F is the only G-perfect form and
  //there are no other well-rounded G-minimal classes
  return true;
 end if;

 if Gprank(G1,FG) ge dim1 and prank(FG) ge n then
  //In this case there may be well rounded minimal classes
  FF:=GMinimalProjections(G1,FG);
  if #FF gt 2 then
   return "Unknown.";
  end if;
  //Determine G-facets and G-facet vectors
  F1:=[FF[1]]; F2:=[FF[2]];
  perp1:=findperpmatrix(F1,FF[2]); perp2:=findperpmatrix(F2,FF[1]);
  R1:=ReynoldsProjection(G1,perp1); R2:=ReynoldsProjection(G1,perp2);
  R1:=Sign(Rationals()!(Trace(R1*FF[2])))*R1;
  R2:=Sign(Rationals()!(Trace(R2*FF[1])))*R2;
  //Check for dead ends
  dead1:=false; dead2:=false;
  L1:=[x : x in minvecs(FG) | matbas3([GProjection(G1,x)])[1] in sub<KMatrixSpace(Rationals(),2*n,2*n)|matbas3(F1)>];
  L2:=[x : x in minvecs(FG) | matbas3([GProjection(G1,x)])[1] in sub<KMatrixSpace(Rationals(),2*n,2*n)|matbas3(F1)>];
  if Dimension(sub<KMatrixSpace(K,1,n)|L1>) lt n then
   //F1 is a dead end
   dead1:=true;
   print "dead1";
   FG1:=FG+R1;
   GG1:=StabilizerOfMinimalClass(FG1);
   if #GG1 gt #G1 and G1 subset GG1 then
    return false;
   end if;
  end if;
  if Dimension(sub<KMatrixSpace(K,1,n)|L2>) lt n then
   //F2 is a dead end
   dead2:=true;
   print "dead2";
   FG2:=FG+R2;
   GG2:=StabilizerOfMinimalClass(FG2);
   if #GG2 gt #G1 and G1 subset GG2 then
    return false;
   end if;
  end if;
  if dead1 and dead2 then
   //in this case FG is the only G-perfect form but
   //there were no supergroups found in the previous steps
   return true;
  end if;

  //Now: consider the case in which at least one of the facets is not a dead end.
  t:=1; counter:=0;
  while not (IsPositiveDefinite(spurform(FG+t*R1)) and hermitianmin(FG+t*R1) eq 1) do
   t:=t/2;
   if counter gt 100 then
    return "Error, counter.";
   end if;
   counter:=counter+1;
  end while;
  print "Found something.";
  GG1:=StabilizerOfMinimalClass(FG+t*R1);
  if #GG1 gt #G1 and G1 subset GG1 then
   return false;
  end if;
  print "#GG1=" cat IntegerToString(#GG1);

  t:=1; counter:=0;
  while not (IsPositiveDefinite(spurform(FG+t*R2)) and hermitianmin(FG+t*R2) eq 1) do
   t:=t/2;
   if counter gt 100 then
    return "Error, counter.";
   end if;
   counter:=counter+1;
  end while;
  print "Found something 2.";
  GG2:=StabilizerOfMinimalClass(FG+t*R2);
  if #GG2 gt #G1 and G1 subset GG2 then
   return false;
  end if;
  print "#GG2=" cat IntegerToString(#GG2);
 end if;

 return "Unknown.";
end function;