(* Importing GAP data *)
ClearAll["Global`*"];
SetDirectory[NotebookDirectory[]];
R=Import["dataRmath"]/.{a->2,b->3,c->5,d->7}//Simplify;
S=ToExpression[Import["dataSmath"]]/.{a->2,b->3,c->5,d->7}//Simplify;
basis=ToExpression[Import["basisMath"]];
CharIrr=Import["dataCharIrrMath"]//Simplify;
Hyp = ToExpression[Import["dataHypMath"]]; (* 1:s, 2:e_s, 3:classno, 4:N_s, 5:det_s*)
Params = ToExpression[Import["dataParamsMath"]]/.{a->2,b->3,c->5,d->7}//Simplify;
L=Length[R];

(* Part of the first function *)
dotR[word_,l_]:=dotR[word,l]=If[word=={},R[[l]][[1]],
If[Length[word]==1,MatrixPower[R[[l]][[Abs[word[[1]]]]],Sign[word[[1]]]],
dotR[word[[1;;Floor[Length[word]/2]]],l].dotR[word[[Floor[Length[word]/2]+1;;Length[word]]],l]
]]//Simplify;

(* Calculation on the basis elements *)
AbsoluteTiming[basisDotR=ParallelTable[dotR[basis[[i]],j],{i,1,Length[basis]},{j,1,L}]][[1]]
Clear[dotR]

(* Function for calculating tr(w1 w2) *)
FbasisProd[w1_,w2_]:=Simplify[
Which[w1=={}&&w2=={},Total[Table[Length[basisDotR[[1,l]]]/S[[l]],{l,1,L}]],
w1=={},Total[Table[Tr[basisDotR[[Flatten[Position[basis,w2]][[1]],l]]]/S[[l]],{l,1,L}]],
w2=={},Total[Table[Tr[basisDotR[[Flatten[Position[basis,w1]][[1]],l]]]/S[[l]],{l,1,L}]],
True,Total[Table[Tr[basisDotR[[Flatten[Position[basis,w1]][[1]],l]].basisDotR[[Flatten[Position[basis,w2]][[1]],l]]/S[[l]]],{l,1,L}]]]];

(* Symmetrize the matrix tr(w w') *)
ToSym[M_]:=Module[{i,j,temp},
temp=M;
For[i=1,i≤Length[M],i++,
For[j=i+1,j≤Length[M],j++,
temp[[j,i]]=temp[[i,j]]
]
];temp];

(* Functions for the calculation of the third function *)
m[i_,s_]:=Module[{o,classnos},
classnos=Prepend[Hyp[[s,3]],1];
o=Hyp[[s,2]];
FullSimplify[Table[Sum[CharIrr[[i,classnos[[k+1]]]]*Exp[-2*Pi*I*k*j/o]/o,{k,0,o-1}],{j,0,o-1}] ]
];

M=Table[m[i,s],{i,1,Length[CharIrr]},{s,1,Length[Hyp]}];

omegapiFun[i_]:=Module[{d,res,s,t,C,o,m},
d=CharIrr[[i,1]];
res=1;
For[j=1, j≤Length[Hyp],j++,
s=Hyp[[j,1]];
t=Params[[s]];
C=Hyp[[s]];
o=C[[2]];
m=M[[i,s]];
res*=Product[(Exp[2*Pi*I/o]^(-k+1)*t[[k]])^(m[[k]]*C[[4]]*C[[2]]*d^(-1)),{k,1,o}]
];
res
];

omegapi=Table[omegapiFun[i],{i,1,Length[CharIrr]}]//Simplify;

G[word_]:=Module[{ginv},
If[word=={},
Sum[omegapi[[l]]*Length[R[[l,1]]]/S[[l]],{l,1,L}],
ginv=(-1)*Reverse[word];
Sum[omegapi[[l]]*Tr[Inverse[basisDotR[[Flatten[Position[basis,word]][[1]],l]]]]/S[[l]],{l,1,L}]//Simplify
]];

(* First function: tr(w), w a basis element *)
trMat=ToSym[PadLeft[ParallelTable[FbasisProd[basis[[i]],basis[[j]]],{i,1,Length[basis]},{j,i,Length[basis]}]]];
Export["trMat.txt",trMat,"TSV"];

(* Third function: tr(w), w a basis element *)
trMatG=ParallelMap[G,basis]//Simplify
Export["trMatG.txt",trMatG,"TSV"];

Det[trMat]//Simplify//FactorInteger

G[{}]//Simplify//FactorInteger