#Given x a type ("F","I", or "E") and y the number of roots,
#finds the associated coxeter group.  For type "I", creates the group as
#a subgroup of E_8
#and for types "F" and "E" uses the Lie Algebra package. W is the group.

weyl:=function(x,y)
local x1,x2,x3,x4,f,r,W,gen,gens,L,R,G, p, orb,pgens,b1,b2,b3,b4,w8;

if  not x  in ["E","F","G","I"] then 
Print("Error ",x," is not a type of a Weyl group\n");
return false;
fi;

if x="I" then
  if not y in [3,4] then
    Print("Error this is not a Weylgroup of type I\n");
   return false;
  fi;

# Create E_8. We will make I_3,I_4 as a subgroup of E_8.
  L:=SimpleLieAlgebra("E", 8, Rationals);
  R:=RootSystem(L);
  G:=WeylGroup(R);
#Redefine the group as a permutation group.
  gens:=GeneratorsOfGroup(G);
  p:=PositiveRoots(R);
  orb:=Orbit(G,p[1]);
  w8:=List(gens,y->PermList(List(orb,x->Position(orb,x^y))));
# now make the gens for I(4);
 b1:=w8[2]*w8[5];
 b2:=w8[4]*w8[6];
 b3:=w8[3]*w8[7];    
 b4:=w8[1]*w8[8];    
 if y = 3 then
  return Group(b1,b2,b3); #Gives I_3 as subgp of E_8
 else
  return Group(b1,b2,b3,b4); #Gives I_4 as subgp of E_8
 fi; 




##Creates our coxeter group given x a type "F" or "E" and y the number of roots
# e.g. x:="F"; y:=4;
else
  L:=SimpleLieAlgebra(x, y, Rationals);
  R:=RootSystem(L);
  G:=WeylGroup(R);
  gens:=GeneratorsOfGroup(G);

#Redefine the group as a permutation group.

  p:=PositiveRoots(R);
  orb:=Orbit(G,p[1]);
  pgens:=List(gens,y->PermList(List(orb,x->Position(orb,x^y))));
  W:=Group(pgens);

fi;

return W;
end;

checkstronglyrealofnormalizersofparabolics:=function(weyltype,weylnumber)
local i,l,P,m,j,k,lfilt,genforP,theanswer,N,ord2,conj,ct,
rep,existinv,sum,c,W,gen;
W:=weyl(weyltype,weylnumber);
if W = false then
Print("Sorry no Weyl group\n");
return fail;
fi;
gen:=GeneratorsOfGroup(W);
l:=List([1..Length(gen)], c->0);
P:=[1..2^Length(gen)-1]; 


#here we create a list of all possible parabolic 
#subgroups of W

for m in [1..2^Length(gen)-1]
do
    
	#add one in binary to a list of size 
	#Length(gen) = number of generating reflections

   for i in [1..Length(gen)]
   do
     if l[i] = 0 then
         l[i]:=1; break;
       else 
         l[i]:=0;
     fi;
   od;

   lfilt:=Filtered([1..Length(gen)], x->l[x]=1);    
   genforP:=List(lfilt, x->gen[x]);    
   P[m]:=Subgroup(W, genforP);  
od;


theanswer:=true;

#for each parabolic subgroup, we create the normalizer, create a list
#of conjugacy class representatives, and filter the list to find the
#conjugacy classes of involutions 

for m in [1..2^Length(gen)-1]
  do
    N:=Normalizer(W, P[m]);
    conj:=ConjugacyClasses(N);
    rep:=List(conj, x->Representative(x));
    ord2:=Filtered([1..Length(rep)], x->Order(rep[x])=2);
    ct:=CharacterTable(N);

#Create a list of the sums 
#of Class Multiplication Coefficients

    existinv:=true;
    sum:=List([1..Length(conj)], x->0);
    for k in [1..Length(rep)] 
    do
      for i in ord2 
      do  
        for j in ord2 
        do
          c:=ClassMultiplicationCoefficient(ct, i, j, k);
	  sum[k] := sum[k]+c;
        od;  
      od;    
    od;

#for each conjugacy class, we find whether
#there are two involutions multiplying into that class
#we define the boolean existinv to be "false" if a conjugacy class rep. of order >2
#does not have 2 such involutions.

    for k in [1..Length(rep)]
    do
      if sum[k] = 0 then
        if Order(rep[k]) > 2 then
	    existinv := false; break;
        fi;  
      fi;
    od;

    #If existinv is fals for any conjugacy class, 
    #define the boolean theanswer to be false and stop the loop 

    if existinv = false then
        theanswer:=false; break;
      else
        theanswer:=true;
    fi;
  od;
return theanswer;
end;