(*This Mathematica code computes the Lagrangian Floer potential of the orbifold sphere that has three cone points of order a,b,c. This code works under the assumption that a,b,c>=3 and that the computation is up to the J-th generation.*)

(*Modify the following parameter if needed.*)
a=3;b=4;c=5;J=4;

(*Global Variables*)
ord[0]=a;ord[1]=b;ord[2]=c;

(* The boundary word
x1 x2^{-e1} x3 x4^{-e2}... x_{6k-1}x_ {6k}^{-e_{3k}}
will be represented by {k, {e1,e2,...,e_{3k}}};
Here, k>=1, x1=alpha=let[0], x2=beta=let[1],x3=gamma=let[2],x4=alpha=let[0] etc.;*)

(*We will set;
p[0]=p[1]=p[2]=k : number of \alpha's, \beta's,\gamma's;
P[0]=e2+e5+e8+... : number of x-corners;
P[1]=e1+e4+e7+...: number of y-corners;
P[2]=e3+e6+e9+...: number of z-corners; *)

p[0]=0;p[1]=0;p[2]=0;P[0]=0;P[1]=0;P[2]=0;

(* Write down the exponents Mod a,b or c *)
StandardExp[input_List]:=
Module[{k,el,i,eln},
If[k<=0,Return[Null]];
k=input[[1]];
el=input[[2]];
eln={};
(* Output is {k, {e1 mod ord[1], e2 mod ord[0], e3 mod ord[2],...}} *)
For[i=1,i<=3k,i++,
eln=Append[eln,Mod[el[[i]],ord[Mod[2i-1,3]]]];
];
Return[{k,eln}];
];

(* Compute Area, p,q,r,P,Q,R values *)
(* Output form is {Area, {P, Q, R}} 
where Area is the ratio Area (polygon) / Area(minimal xyz triangle) *)

AreaOfWord[input_List]:=
Module[{k,els,i,P,Q,R},
k=input[[1]];
If[k<=0,Return[Null]];
els=StandardExp[input][[2]];
P=0;
Q=0;
R=0;
For[i=1,i<=k,i++,
P+=els[[3i-1]];
Q+=els[[3i-2]];
R+=els[[3i]];
];
Return[{3(P+Q+R)+8(P/ord[0]+Q/ord[1]+R/ord[2]-1)/(1-1/ord[0]-1/ord[1]-1/ord[2]),{P,Q,R}}];
];

(* Area of the word whose parameters are given by (i,P,Q,R), where i is the generation *)
Area[pl_List]:=Module[{P,Q,R},
P=pl[[2]];Q=pl[[3]];R=pl[[4]];
Return[{3(P+Q+R)+8(P/ord[0]+Q/ord[1]+R/ord[2]-1)/(1-1/ord[0]-1/ord[1]-1/ord[2]),{P,Q,R}}];
];

(* Results of all possible cut-glue operations on a given input *)
CutGlue[input_List]:=Module[{k,el,i,els,eln,eln0,j,output},
k=input[[1]];
If[k<=0,Return[Null]];
els=StandardExp[input][[2]];
(* To consider cyclic conjugations, we double the exponents *)
els=Join[els,els];
i=1;eln={};output={};
While[i<=3k,
If[els[[i]]>= 2,
For[j=0,j<=els[[i]]-2,j++,
eln=Join[Take[els,i-1],{j,-2,-2,els[[i]]-j-2},Take[els,i-3k]];
output=Append[output,StandardExp[{k+1,eln}]];
]];
If[els[[i]]>=1,
j=i+1;
While[els[[j]]==0,j++];
If[j<=3k,
eln=Join[Take[els,i-1],{els[[i]]-1,-2},Table[-3,{j-i}],{-2,els[[j]]-1},Take[els,j-3k]];
output=Append[output,StandardExp[{k+1,eln}]];
,
(*If i<= 3k < j, then we consider cyclic conjugation so that the first letter is still let[0]. *)
eln0=Join[{els[[i]]-1,-2},Table[-3,{j-i}],{-2,els[[j]]-1},Take[els,{j+1,i+3k-1}]];
eln=Join[Take[eln0,-(i+2)],Take[eln0,3k-i+1]];
output=Append[output,StandardExp[{k+1,eln}]];];
];
i++;
];
Return[output];
];
(*Compare two words in1 and in2*)
Compare[in1_List,in2_List]:=Module[{k,el1,el2,i},
k=in1[[1]];
If[k!=in2[[1]],Return[False]];
el1=in1[[2]];
el2=in2[[2]];
For[i=0,i<=k-1,i++,
(* Take cyclic conjugation, with 6 letters as one unit *)
el2=Join[Take[el2,-3k+3i],Take[el2,3i]];
If[el1==el2,Return[True]];
];
Return[False];
];
(*  Write w = u^n such that u is power-free, and output n. Note that we make sure each copy of u starts with the same letter, let[0]. Input form is {k, {e1,e2,...,e_{3k}}}*)
FindPower[input_List]:=Module[{m,l,dl,dll,i,p,j,Flag},
m=input[[1]];l=input[[2]];
dl=Divisors[m];dll=Length[dl];
For[i=1,i<=dll-1,i++,
p=3dl[[i]];
j=1;
Flag=True;
While[j<=3m-p&&l[[j]]==l[[j+p]],j++];
If[j<3m-p+1,Flag=False];
If[Flag,Break[]];
];
Return[If[Flag,3m/p,1]];
];
(* Print the Potential term given by {coefficient,{P,Q,R}} *)
PotTerm[pl_List]:=Module[{ar,P,Q,R,t},
t=pl[[1]];
P=pl[[2,1]];
Q=pl[[2,2]];
R=pl[[2,3]];
ar = 3(P+Q+R)+8(P/a+Q/b+R/c-1)/(1-1/a-1/b-1/c);
Return[StringJoin[
If[t==1,"",If[t==-1,"-",ToString[TeXForm[t]]]],
If[P==0,"",If[P==1,"x",StringJoin["x^{",ToString[P],"}"]]],
If[Q==0,"",If[Q==1,"y",StringJoin["y^{",ToString[Q],"}"]]],
If[R==0,"",If[R==1,"z",StringJoin["z^{",ToString[R],"}"]]],
"q^{",ToString[TeXForm[ar]],"}"]];
];

(* Determine the sign of a string. Useful when printing terms separated by + or - *)
SignString[l_String]:=If[Characters[l][[1]]=="-","","+"];
ConCatStr[l_List]:=
Module[{res,i},
res=l[[1]];
For[i=2,i<=Length[l],i++,
res=StringJoin[res,SignString[l[[i]]],l[[i]]]];
Return[res]];

(* s[j] is the set of cyclic words in j-th generation.
ss[j] is s[j] after removing redundancies.
as[j] is the parameter set (i,P,Q,R,eta) of ss[j] *)
s[0]=Table[{0,i,{ord[i]}},{i,0,2}];
s[1]={StandardExp[{1,{-2,-2,-2}}]};

For[j=2,j<=J,j++,
s[j]={};
For[i=1,i<=Length[s[j-1]],i++,
s[j]=Join[s[j],CutGlue[s[j-1][[i]]]];
];
ss[j]={};
For[i=1,i<=Length[s[j]],i++,
flag=False;
For[l=1,l<=Length[ss[j]]&&!flag,l++,
If[Compare[ss[j][[l]],s[j][[i]]],flag=True]];
If[!flag,ss[j]=Append[ss[j],s[j][[i]]]];
];
as[j]={};
For[i=1,i<=Length[ss[j]],i++,
as[j]=Append[as[j],Join[{j},AreaOfWord[ss[j][[i]]][[2]],{FindPower[ss[j][[i]]]}]];
];
];

(* -1st generation word : \alpha\beta\gamma . p=-1, P=Q=R=1, eta = 1*)
as[-1]={{-1,1,1,1,1}};

(* 0th generation words : (\beta\gamma)^a, (\gamma\alpha)^b,(\alpha\beta)^c.
(p,P,Q,R,eta)=(0,a,0,0,a),(0,0,b,0,b),(0,0,0,c,c) *)
as[0]={{0,a,0,0,a},{0,0,b,0,b},{0,0,0,c,c}};
as[1]={{1,a-2,b-2,c-2,1}};

(* Compute potential coefficients.
Output is for the list of {coefficient, {P,Q,R}}.
Duplicate polynomial terms will be merged. *)
potential={{-1,{1,1,1}}};
For[j=0,j<=J,j++,
For[i=1,i<=Length[as[j]],i++,
duplicate=False;
p=as[j][[i,1]];
P=as[j][[i,2]];
Q=as[j][[i,3]];
R=as[j][[i,4]];
eta=as[j][[i,5]];
For[k=1,k<=Length[potential]&&!duplicate,k++,
If[{P,Q,R}== potential[[k,2]],
potential=ReplacePart[potential,{k,1}->potential[[k,1]]+(-1)^(p+Q+R) (2p+P+Q+R)/eta];
duplicate=True;];
];
If[!duplicate,potential=Append[potential,{(-1)^(p+Q+R) (2p+P+Q+R)/eta,{P,Q,R}}]
];
];
];
wl={};
For[i=1,i<=Length[potential],i++,
wl=Append[wl,PotTerm[potential[[i]]]];
];

(* The Lagrangian-Floer potential for (a,b,c) orbifold sphere, up to the J-th generation *)
ConCatStr[wl]
(* Print boundary words for each generation *)
(* The letter a represents \alpha, b for \beta and c for \gamma.*)
let[0]="a";let[1]="b";let[2]="c";
PrintWord[input_List]:=
Module[{k,i,str,els},
If[k<=0,Return[Null]];
k=input[[1]];
els=StandardExp[input][[2]];
str={};
For[i=1,i<=3k,i++,
(* x_{2i-1} is let[2i-2 mod 3]*)
str=StringJoin[str,let[Mod[2i-2,3]]];
If[els[[i]]!=0,str=StringJoin[str,let[Mod[2i-1,3]],"^{",ToString[-els[[i]]],"}"]];
];
Return[str];
];

For[j=2,j<=J,j++,
Print["Words in S^",j,":"];
For[i=1,i<=Length[ss[j]],i++,Print["(",PrintWord[ss[j][[i]]],")"];
]
]