(* ::Package:: *)

(************************************************************************)
(* This file was generated automatically by the Mathematica front end.  *)
(* It contains Initialization cells from a Notebook file, which         *)
(* typically will have the same name as this file except ending in      *)
(* ".nb" instead of ".m".                                               *)
(*                                                                      *)
(* This file is intended to be loaded into the Mathematica kernel using *)
(* the package loading commands Get or Needs.  Doing so is equivalent   *)
(* to using the Evaluate Initialization Cells menu command in the front *)
(* end.                                                                 *)
(*                                                                      *)
(* DO NOT EDIT THIS FILE.  This entire file is regenerated              *)
(* automatically each time the parent Notebook file is saved in the     *)
(* Mathematica front end.  Any changes you make to this file will be    *)
(* overwritten.                                                         *)
(************************************************************************)



SetAttributes[Diag, Orderless];
Draw[expr_] := expr /. diag_Diag :> Draw[diag];
Draw[diag_Diag] := Module[
{n=Length[diag]},
Graphics[{
Line[{{0,0}, {2n+1,0}}],
(List@@diag) /.  {
Ar[i_, j_] :> {
BezierCurve[{
{i,0}, {(i+j)/2, Abs[i-j]}, {j,0}
}],
Line[{{j-0.2+0.2Sign[i-j], 0.4}, {j,0}, {j+0.2+0.2Sign[i-j], 0.4}}]
}
}
}
]
];


Place[{Ar},{i_, j_}] := {Diag[Ar[i,j]], Diag[Ar[j,i]]};
Place[{Ar, objs__}, {i_, rest__}] := Flatten[Table[
Outer[Join,
Place[{Ar}, {i, {rest}[[k]]}],
Place[{objs}, Delete[{rest}, k]]
],
{k, Length[{rest}]}
]];
Diagrams[k_.*Ar]  := Place[Table[Ar, {k}], Range[2k]];


DiagramRotateLeft[d_Diag]:=Module[
{labels=Union@@(List@@Apply[List,d,{1}])},
d/.Thread[Rule[labels,RotateLeft[labels]]]
];
DiagramRotateToMinimal[expr_]:=expr/.d_Diag:>Module[
{bd=d,rd=DiagramRotateLeft[d]},
While[rd=!=d,
bd=First[Sort[{bd,rd}]];
rd=DiagramRotateLeft[rd]];
bd
];


Place[{r:(TC|R3b), objs__}, {i_, rest__}] := Flatten[Table[
Outer[Join,
Place[{r}, {i, {rest}[[j]], {rest}[[k]]}],
Place[{objs}, Delete[{rest}, {{j}, {k}}]]
],
{k, 2, Length[{rest}]}, {j,1, k-1}
]];
Place[{R3b},{i_, j_, k_}] := Permutations[{i,j,k}] /. {i1_,j1_,k1_} :> Diag[R3b[i1,j1,k1]];
Diagrams[R3b] := Place[{R3b}, {1,2,3}];
Diagrams[R3b+k_.*Ar] /; k>0 := Flatten[
Place[#, Range[2k+3]]& /@ Permutations[Table[Ar, {k}]~Append~R3b]
];
Place[{TC},{i_, j_, k_}] := Diag /@ {TC[i,j,k], TC[j,k,i], TC[k,i,j]};
Diagrams[TC] := Place[{TC}, {1,2,3}];
Diagrams[TC+k_.*Ar] /; k>0 := Flatten[
Place[#, Range[2k+3]]& /@ Permutations[Table[Ar, {k}]~Append~TC]
];
Diagrams[TC+k_.*Ar] /; k<0 := {};
Place[{R2c},{i_, j_}] := {Diag[R2c[i,j]], Diag[R2c[j,i]]};
Place[{r:(R2c), objs__}, {i_, rest__}] := Flatten[Table[
Outer[Join,
Place[{r}, {i, {rest}[[j]]}],
Place[{objs}, Delete[{rest}, j]]
],
{j, 1, Length[{rest}]}
]];
Diagrams[R2c] := Place[{R2c}, {1,2}];
Diagrams[R2c+k_.*Ar] /; k>0 := Flatten[
Place[#, Range[2k+2]]& /@ Permutations[Table[Ar, {k}]~Append~R2c]
];
Diagrams[R2c+k_.*Ar] /; k<0 := {};
Place[{R1d},{i_}] := {Diag[R1d[i]]};
Place[{R1a},{i_}] := {Diag[R1a[i]]};
Place[{r:(R1d|R1a), objs__}, {i_, rest__}] := Flatten[
Outer[Join,
Place[{r}, {i}],
Place[{objs},{rest}]
]];
Diagrams[R1d] := Place[{R1d}, {1}];
Diagrams[R1a] := Place[{R1a}, {1}];
Diagrams[R1d+k_.*Ar] /; k>0 := Flatten[
Place[#, Range[2k+1]]& /@ Permutations[Table[Ar, {k}]~Append~R1d]
];
Diagrams[R1a+k_.*Ar] /; k>0 := Flatten[
Place[#, Range[2k+1]]& /@ Permutations[Table[Ar, {k}]~Append~R1a]
];
Diagrams[R1d+k_.*Ar] /; k<0 := {};
Diagrams[R1a+k_.*Ar] /; k<0 := {};


NormalizeDiag[diag_Diag] := Module[
{indices = Union@@ (List @@ diag /. Ar -> List)},
diag /. Thread[indices -> Range[Length[indices]]]
];
R[Diag[lft___, R3b[i_,j_,k_], rgt___]] := (
(* RHS of (9) *)
NormalizeDiag[Diag[lft, Ar[i,j], Ar[i+0.5, k], Ar[j+0.5, k+0.5], rgt]]
+NormalizeDiag[Diag[lft, Ar[i,j], Ar[i+0.5, k], rgt]]
+NormalizeDiag[Diag[lft, Ar[i,j], Ar[j+0.5, k], rgt]]
+NormalizeDiag[Diag[lft, Ar[i,k], Ar[j, k+0.5], rgt]]
(* LHS of (9) *)
-NormalizeDiag[Diag[lft, Ar[j, k], Ar[i,k+0.5], Ar[i+0.5, j+0.5], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,k], Ar[i+0.5, j], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,j+0.5], Ar[j, k], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,k+0.5], Ar[j, k], rgt]]
);
RWrong1[Diag[lft___, R3b[i_,j_,k_], rgt___]] := (
(* RHS of (9) *)
NormalizeDiag[Diag[lft, Ar[i,j], Ar[i+0.5, k], Ar[j+0.5, k+0.5], rgt]]
+NormalizeDiag[Diag[lft, Ar[i,j], Ar[i+0.5, k], rgt]]
+NormalizeDiag[Diag[lft, Ar[i,j], Ar[j+0.5, k], rgt]]
+NormalizeDiag[Diag[lft, Ar[i,k], Ar[j, k+0.5], rgt]]
(* LHS of (9) *)
+NormalizeDiag[Diag[lft, Ar[j, k], Ar[i,k+0.5], Ar[i+0.5, j+0.5], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,k], Ar[i+0.5, j], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,j+0.5], Ar[j, k], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,k+0.5], Ar[j, k], rgt]]
);
RWrong2[Diag[lft___, R3b[i_,j_,k_], rgt___]] := (
(* RHS of (9) *)
NormalizeDiag[Diag[lft, Ar[i,j], Ar[i+0.5, k], Ar[j+0.5, k+0.5], rgt]]
+NormalizeDiag[Diag[lft, Ar[i,j], Ar[i+0.5, k], rgt]]
+NormalizeDiag[Diag[lft, Ar[i,j], Ar[j+0.5, k], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,k], Ar[j, k+0.5], rgt]]
(* LHS of (9) *)
-NormalizeDiag[Diag[lft, Ar[j, k], Ar[i,k+0.5], Ar[i+0.5, j+0.5], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,k], Ar[i+0.5, j], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,j+0.5], Ar[j, k], rgt]]
-NormalizeDiag[Diag[lft, Ar[i,k+0.5], Ar[j, k], rgt]]
);
RMinus[n_, Diag[lft___, R3b[i_,j_,k_], rgt___]] := Module[
{Expand2, Expand3, n1},
n1=n-Length[{lft, rgt}];
Expand2[Diag[Ar[i1_,j1_], Ar[i2_,j2_]]] := Sum[
(-1)^(k1+k2) * NormalizeDiag[Flatten[Diag[
lft,
Diag@@Table[Ar[i1+0.01q, j1+0.01q],{q, k1}],
Diag@@Table[Ar[i2+0.01q, j2+0.01q],{q, k2}],
rgt
]]],
{k1, n1-1}, {k2, n1-k1}
];
Expand3[Diag[Ar[i1_,j1_], Ar[i2_,j2_], Ar[i3_,j3_]]] := Sum[
(-1)^(k1+k2+k3) * NormalizeDiag[Flatten[Diag[
lft,
Diag@@Table[Ar[i1+0.01q, j1+0.01q],{q, k1}],
Diag@@Table[Ar[i2+0.01q, j2+0.01q],{q, k2}],
Diag@@Table[Ar[i3+0.01q, j3+0.01q],{q, k3}],
rgt
]]],
{k1, n1-2}, {k2, n1-k1-1}, {k3,n1-k1-k2}
];
Expand[Plus[
(
Diag[Ar[i,j], Ar[i+0.5, k]]+Diag[Ar[i,j], Ar[j+0.5, k]]+Diag[Ar[i,k], Ar[j, k+0.5]]
-Diag[Ar[i,k], Ar[i+0.5, j]]-Diag[Ar[i,j+0.5], Ar[j, k]]-Diag[Ar[i,k+0.5], Ar[j, k]]
) /. d_Diag :> Expand2[d],
(Diag[Ar[i,j], Ar[i+0.5, k], Ar[j+0.5, k+0.5]] - Diag[Ar[j, k], Ar[i,k+0.5], Ar[i+0.5, j+0.5]]) /. d_Diag:> Expand3[d]
]]
];
R[n_, Diag[lft___, R2c[i_,j_], rgt___]] := Module[
{n1},
n1=n-Length[{lft, rgt}];
Sum[
Expand[(-1)^k * s* NormalizeDiag[Flatten[Diag[
lft,
Diag@@Table[Ar[i+0.01q, j+s*0.01q],{q, k}],
rgt
]]]],
{k,2,n1},{s,-1,1,2}
]
];
R[n_, Diag[lft___, R1d[i_], rgt___]] := Module[
{n1},
n1=n-Length[{lft, rgt}];
Sum[
Expand[(-1)^k * NormalizeDiag[Flatten[Diag[
lft,
Diag@@Table[Ar[i+0.01q, i+0.5+0.01q],{q, k}],
rgt
]]]],
{k,2,n1}
]
];
R[n_, Diag[lft___, R1a[i_], rgt___]] := Module[
{n1},
n1=n-Length[{lft, rgt}];
Sum[
Expand[(-1)^k * NormalizeDiag[Flatten[Diag[
lft,
Diag@@Table[Ar[i+0.5+0.01q, i+0.01q],{q, k}],
rgt
]]]],
{k,2,n1}
]
];
R[Diag[lft___, TC[i_,j_,k_], rgt___]] := (
+NormalizeDiag[Diag[lft, Ar[k,i], Ar[k+0.5, j], rgt]]
-NormalizeDiag[Diag[lft, Ar[k+0.5, i], Ar[k, j], rgt]]
);
ContainsShortArrow[diag_Diag] := (1==Min[List@@diag /. Ar[i_,j_] :> Abs[i-j]]);
DescendingQ[diag_Diag] := And @@ (OrderedQ /@ diag);


Options[PSpace] = {Debug -> False, Gr -> False, ReduceEquivalences -> True};
PSpace[m_] /; m<2 := QuotientSpace[Diagrams[m Ar], {}];
PSpace[m_, reltypes_List, opts___] /; m>=2 := Module[
{
debug = Debug /. {opts} /. Options[PSpace],
gr = Gr /. {opts}  /. Options[PSpace],
reduceequivalences = ReduceEquivalences /. {opts}  /. Options[PSpace],
diags, rels, i, filter, equalities, equivclasses, pos, redrule
},
If[debug, Print[Date[], ": Starting work..."]];
diags=Join@@Table[Diagrams[k Ar], {k,If[gr, m, 1], m}];
filter=Identity;
reltypes /. {
"R1" :> (
diags = Select[diags, !ContainsShortArrow[#]&]
),
"Descending" :> (
filter = Select[#, DescendingQ]&;
diags = Select[diags, DescendingQ]
)
};
If[debug, Print[Date[], ": Computed ",Length[diags]," diagrams..."]];
rels = Join@@(reltypes /. {
"R1" -> {},
"Descending" -> {},
"Round" :> (diags /.
d_Diag :> (NormalizeDiag[d /. {1 -> 2Length[d], i_Integer :> i-1}] -d)
),
"R3b" :> (
(R /@ Join@@Table[Diagrams[R3b+k Ar] // filter, {k, If[gr, m-2, 0], m-2}]) /.diag_Diag /; Length[diag]>m :> 0
),
"RWrong1" :> (
(RWrong1 /@ Join@@Table[Diagrams[R3b+k Ar], {k, If[gr, m-2, 0], m-2}]) /.diag_Diag /; Length[diag]>m :> 0
),
"RWrong2" :> (
(RWrong2 /@ Join@@Table[Diagrams[R3b+k Ar], {k, If[gr, m-2, 0], m-2}]) /.diag_Diag /; Length[diag]>m :> 0
),
"R3bMinus" :> (
RMinus[m, #]& /@ Join@@Table[Diagrams[R3b+k Ar] // filter, {k, If[gr, m-2, 0], m-2}]
),
"R2c" :> (
R[m, #]& /@ Join@@Table[Diagrams[R2c+k Ar] // filter, {k, If[gr, m-2, 0], m-2}]
),
"R1d" :> (
R[m, #]& /@ Join@@Table[Diagrams[R1d+k Ar] // filter, {k, If[gr, m-2, 0], m-2}]
),
"R1a" :> (
R[m, #]& /@ Join@@Table[Diagrams[R1a+k Ar] // filter, {k, If[gr, m-2, 0], m-2}]
),
"TC" :> (
R /@ Join@@Table[Diagrams[TC+k Ar], {k, If[gr, m-2, 0], m-2}]
)
});
reltypes /. {
"R1" :> (
rels = DeleteCases[rels /. {d_Diag /; ContainsShortArrow[d] :> 0}, 0]
)
};
If[debug, Print[Date[], ": Computed ",Length[rels]," relations..."]];
If[reduceequivalences,
equalities = Cases[rels, d1_Diag-d2_Diag :> {d1,d2}];
equivclasses = Fold[
(
pos=First/@Position[#1,#2];
Append[Delete[#1,List/@pos],Union@@(#1[[pos]])]
)&,
equalities,
Union@@equalities
];
redrule = Dispatch[Flatten[Table[
class=equivclasses[[i]];
(# -> First[class])& /@ Rest[class],
{i,Length[equivclasses]}
]]];
diags = Union[diags /. redrule];
rels = Union[rels /. redrule];
];
QuotientSpace[diags, rels]
]


Options[Dim] = {Debug -> False, UseLinBox -> False};
Dim[QuotientSpace[diags_List, rels_List], opts___] := Module[
{
debug = Debug /. {opts} /. Options[Dim],
uselinbox = UseLinBox /. {opts} /. Options[Dim],
diagtoindex, mat, rel, matfile, rank
},
diagtoindex = Dispatch[Thread[Rule[diags, Range[Length[diags]]]]];
If[uselinbox,
matfile = OpenWrite["matfile.matrix"];
WriteString[matfile, Length[rels], " ", Length[diags], " S\n"];
Do[
rel=rels[[i]];
WriteString[matfile, StringJoin @@ Flatten[{
ToString[Count[{rel}, _Diag, Infinity]], " ",
{
ToString[(# /. diagtoindex)-1], " ",
ToString[Coefficient[rel, #]], " "
}& /@ Cases[{rel}, diag_Diag, Infinity],
"\n"
}]],
{i,Length[rels]}
];
matfile = Close[matfile];
rank = RunThrough["./ComputeRank "<>matfile, Random[]][[2]];
Length[diags]-rank,
(* else *) mat = SparseArray[
Join @@ Table[
rel=rels[[i]];
{i,# /. diagtoindex} -> Coefficient[rel, #]& /@ Cases[{rel}, diag_Diag, Infinity],
{i,Length[rels]}
],
{Length[rels], Length[diags]}
];
If[debug, Print[Date[], ": Computed mat..."]];
Length[diags]-MatrixRank[mat]
]
];
DimP[m_, reltypes_List, opts___] := Dim[PSpace[m,reltypes, opts], opts]


Analyze[m_Integer, reltypes_List, opts___Rule] := Module[
{qs},
Table[
{
qs=PSpace[n, reltypes, Gr -> True, opts];
(Length /@ qs) -> Dim[qs, opts],
qs=PSpace[n, reltypes, Gr -> False, opts];
(Length /@ qs) -> Dim[qs, opts]
},
{n,2,m}
]
];
Analyze[reltypes_List, opts___Rule] := Analyze[4, reltypes, opts];



