(* ::Package:: *)

(* This program computes the set of all homomorphisms between the commutative monoids X and Y, where X has to be of order 4.
	Note that the program does not check whether X and Y are indeed commutative monoids! *)
MonoidHom4[X_,Y_]:= Module[{Hom, x,y,z,Cand,good,v,w}, Hom = List[]; (* Module defines local variables; Hom = List of all Homomorphisms: will be the output *)
If[Dimensions[X]!={4,4},Print["First Module has to be of order 4!"];Abort[]]; (* Check if first module is of right order *)
If[Or[Length[Dimensions[Y]]!=2, Dimensions[Y][[1]]!= Dimensions[Y][[2]]],Print["Second entry no module!"];Abort[]]; (* Check if second entry is a square matrix *)
For[x=1, x<= Dimensions[Y][[2]], x++1,
For[y=1, y<= Dimensions[Y][[2]], y++1,
For[z=1, z<= Dimensions[Y][[2]], z++1, Cand = {1,x,y,z}; good = True; (* Cand = Candidate for a Hom; saved as {f(1),f(2),f(3),f(4)} *)
For[v=2,v<= 4,v++1, If[good == False,Break[]]; (* always check if inner loop already produced a False \[Rule] then we can disregard Cand already *)
For[w=v,w<=4,w++1, (* Due to symmetry only half of of the combinations have to be checked, hence v starts at w *)
If[Cand[[X[[v,w]]]] != Y[[Cand[[v]],Cand[[w]]]], good = False; Break[]]; (* Check for Hom, has to hold for all comb. of v and w *)
]];
If[good == True, AppendTo[Hom,Cand]] (* add Cand to Hom if good still True after v and w loops *)
]]];
Return[Hom]]



(* This program computes the set of all homomorphisms between the commutative monoids X and Y, where X has to be of order 3.
	Note that the program does not check whether X and Y are indeed commutative monoids! *)
MonoidHom3[X_,Y_]:= Module[{Hom, x,y,Cand,good,v,w}, Hom = List[]; (* for comments see MonoidHom4 *)
If[Dimensions[X]!={3,3},Print["First Module has to be of order 3!"];Abort[]];
If[Or[Length[Dimensions[Y]]!=2, Dimensions[Y][[1]]!= Dimensions[Y][[2]]],Print["Second entry no module!"];Abort[]];
For[x=1, x<= Dimensions[Y][[2]], x++1,
For[y=1, y<= Dimensions[Y][[2]], y++1, Cand = {1,x,y}; good = True;
For[v=2,v<= 3,v++1, If[good==False,Break[]];
For[w=v,w<=3,w++1,
If[Cand[[X[[v,w]]]] != Y[[Cand[[v]],Cand[[w]]]], good = False; Break[]];
]];
If[good == True, AppendTo[Hom,Cand]]
]];
Return[Hom]]




(* This program computes the set of all homomorphisms between the commutative monoids X and Y, where X has to be of order 2.
	Note that the program does not check whether X and Y are indeed commutative monoids! *)
MonoidHom2[X_,Y_]:= Module[{Hom, x,Cand,v,w}, Hom = List[]; (* for comments see MonoidHom4 *)
					If[Dimensions[X]!={2,2},Print["First Module has to be of order 2!"];Abort[]];
					If[Or[Length[Dimensions[Y]]!=2, Dimensions[Y][[1]]!= Dimensions[Y][[2]]],Print["Second entry no module!"];Abort[]];
					For[x=1, x<= Dimensions[Y][[2]], x++1, Cand = {1,x}; 
					If[Cand[[X[[2,2]]]] == Y[[Cand[[2]],Cand[[2]]]], AppendTo[Hom,Cand]]];
					Return[Hom]]


(* This program computes a table of the cardinality of Hom(Y,Z), where Y and Z have to be commutative monoids of order between 2 and 4.
	It requires a list of such commutative monoids as its input.
	Note that the program does not check whether the elements of its input are indeed commutative monoids! *)
MonoidHomCardTable[X_]:=Module[{n, myTable, FailCount, k, l}, (* Module defines local variables *)
						If[Or[Length[Dimensions[X]]==1, Dimensions[X][[2]]==Dimensions[X][[3]]], n = Dimensions[X][[1]], Print["Please enter a list of monoids!"];Abort[]]; 
						(* check whether all entries of X are square matrices (only works if all monoids in X are of the same order) *)
						myTable = ConstantArray[0, {n, n}];
						For[k=1, k<=n, k++1,
						For[l=1, l<=n, l++1,
						FailCount = 0;
						If[Dimensions[X[[k]]]=={2,2},myTable[[k,l]] = Length[MonoidHom2[X[[k]],X[[l]]]],FailCount = FailCount +1]; (* fill table using above functions *)
						If[Dimensions[X[[k]]]=={3,3},myTable[[k,l]] = Length[MonoidHom3[X[[k]],X[[l]]]],FailCount = FailCount +1];
						If[Dimensions[X[[k]]]=={4,4},myTable[[k,l]] = Length[MonoidHom4[X[[k]],X[[l]]]],FailCount = FailCount +1];
						If[FailCount == 3, Print[StringForm["Element number `` is not of the right order!",k]];Abort[]] (* Second check for wrong inputs *)
						]];
						Return[myTable//MatrixForm]]


(* This program represents Hom(X,Y) as a monoid via its Cayley table, where X has to be a commutative monoid of order between 2 and 4 and Y can be an arbitrary commutative monoid.
	Note that the program does not check whether X and Y are indeed commutative monoids! *)
HomASTable[X_,Y_]:= Module[{FailCount, Hom, M, el, n, f, g, k}, FailCount = 0;
					If[Dimensions[X]=={2,2}, Hom = MonoidHom2[X,Y],FailCount = FailCount +1]; (* compute Hom in dependence of order of X *)
					If[Dimensions[X]=={3,3}, Hom = MonoidHom3[X,Y],FailCount = FailCount +1];
					If[Dimensions[X]=={4,4}, Hom = MonoidHom4[X,Y],FailCount = FailCount +1];
					If[FailCount == 3, Print["First monoid not of order between 2 and 4!"];Abort[]]; (* check if X is of order between 2 and 4 *)
					n = Length[Hom]; M = ConstantArray[0, {n, n}]; el = ConstantArray[0, {Length[Hom[[1]]]}]; (* M is the table that will be the output;
					el will be the "multiplication" of Hom[[f]] and Hom[[g]] with the operation of Y  *)
					For[f=1, f<=n, f++1, (* Hom[[f]] and Hom[[g]] are elements of Hom, thus homomorphisms between X and Y *)
					For[g=1, g<=n, g++1,
					For[k=1, k<=Length[Hom[[1]]], k++1, el[[k]] = Y[[Hom[[f]][[k]],Hom[[g]][[k]]]]]; (* computation of el= Hom[[f]] Hom[[g]] with the operation of Y *)
					M[[f,g]] = Position[Hom,el][[1]][[1]] ]]; (* el is matched to one of the elements of Hom (there has to exist exactly one) *)
					Return[M]]


(* This program computes a table stating to which of its input elements Hom(X,Y) is isomorphic to, where X and Y are elements of its input list.
	It requires a list of commutative monoids as its input. Moreover, the neutral element has to be the first element!
	Note that the program does not check whether the elements of its input are indeed commutative monoids or whether the neutral element is the first one! *)
MonoidHomTable[X_]:= Module[{n, myTable, Hom, k, l, f, Hom2, Hom3, i, j, F, G, m}, (* Module defines local variables *)
					 If[Or[Length[Dimensions[X]]==1, Dimensions[X][[2]]==Dimensions[X][[3]]], n = Dimensions[X][[1]], Print["Please enter a list of monoids!"];Abort[]]; 
					 (* check whether all entries of X are square matrices (only works if all monoids in X are of the same order) *)
					 myTable = ConstantArray[0, {n, n}]; (* final table that will be the output *)
					 f={1,3,2}; (* the only permutation of {1,2,3} leaving 1 first not equal to the identity *)
					 F = {{1,2,4,3},{1,3,2,4},{1,3,4,2},{1,4,2,3},{1,4,3,2}}; (* list of all permutations of {1,2,3,4} leaving 1 first not equal to the identity *)
					 G = {{1,2,4,3},{1,3,2,4},{1,4,2,3},{1,3,4,2},{1,4,3,2}}; (* inverse Elements of elements in F *)
					 For[k=1, k<=n, k++1,
					 For[l=1, l<=n, l++1, Hom = HomASTable[X[[k]],X[[l]]]; (* computes the Cayley table of Hom(X,Y) for all combinations of the input elements *)
					 If[Hom == {{1}}, myTable[[k,l]] = "M0", (* check whether Hom contains only the trivial element *)
					 (* Note that all checks are nested in each other and are only performed if the previous ones returned "False" *)
					 If[Position[X,Hom]!={}, myTable[[k,l]] = Position[X,Hom][[1]][[1]], (* checks whether Hom is already in input list and returns the position if that is the case *)
					 If[Dimensions[Hom][[1]] > 4, myTable[[k,l]] = StringForm["`` El.",Dimensions[Hom][[1]]], (* If Hom has >4 elements, the cardinality is returned as a string *)
					 If[Dimensions[Hom][[1]] == 3, Hom2 = ConstantArray[0, {3, 3}]; Hom3 = ConstantArray[0, {3, 3}]; (* isomorphism check for |Hom| = 3 *)
					 For[i=1, i<=3, i++1,
					 For[j=1, j<=3, j++1, Hom2[[i,j]]=Hom[[ f[[i]],f[[j]] ]]; Hom3[[i,j]] = f[[Hom2[[i,j]]]];  (* exchanges the roles of "1" and "2" *)
					 ]];
					 If[Position[X,Hom3]!={}, myTable[[k,l]] = Position[X,Hom3][[1]][[1]], myTable[[k,l]] = "3 El, \[NotElement] input"], (* check if Hom3 is in input; if not some monoid of order 3 has to be missing in input! *)
					 If[Dimensions[Hom][[1]] == 4, Hom2 = ConstantArray[0, {4, 4}]; Hom3 = ConstantArray[0, {4, 4}]; (* isomorphism check for |Hom|= 4 *)
					 For[m=1, m<=5, m++1, (* same code as for |Hom| = 3 but with additional for loop running over all permutations *)
					 For[i=1, i<=4, i++1,
					 For[j=1, j<=4, j++1, Hom2[[i,j]]=Hom[[ F[[m]][[i]],F[[m]][[j]] ]]; Hom3[[i,j]] = G[[m]][[Hom2[[i,j]]]];
					 ]];
					 If[Position[X,Hom3]!={}, myTable[[k,l]] = Position[X,Hom3][[1]][[1]]]
					 ], myTable[[k,l]] = StringForm["`` El, \[NotElement] input!",Dimensions[Hom][[1]]] (* else condition to check whether some monoid of order 2 or 4 is missing in input *)
					 ]
					 ]]]]]];
					 Return[myTable//MatrixForm]]


(* This program checks whether the natural homomorphism x \mapsto L_x, where L_x(f):=f(x) for f \in Hom[X,Y] is injective.
	X and Y have to be commutative monoids and X has to be of order between 2 and 4.
	Note that this is not checked to full extend by the program!
	The output is either True or False. *)
NatHomInjCheck[X_,Y_]:= Module[{Hom, FailCount, L, k, m}, FailCount = 0;
						    If[Dimensions[X]=={2,2}, Hom = MonoidHom2[X,Y],FailCount = FailCount +1]; (* compute Hom in dependence of order of X *)
							If[Dimensions[X]=={3,3}, Hom = MonoidHom3[X,Y],FailCount = FailCount +1];
							If[Dimensions[X]=={4,4}, Hom = MonoidHom4[X,Y],FailCount = FailCount +1];
							If[FailCount == 3, Print["First monoid not of order between 2 and 4!"];Abort[]]; (* check if X is of order between 2 and 4 *)
							L = ConstantArray[0, {Length[X], Length[Hom]}]; 
							For[k=1, k<=Length[X], k++1,
							For[m=1, m<=Length[Hom], m++1, L[[k,m]] = Hom[[m]][[k]]; (* L lists L_x in its rows for x \in X *)
							]];
							Return[L == DeleteDuplicates[L]]]						


(* This program checks for a list of commutative monoids X which monoids on X are reflexive with respect to one monoid on X.
	X has to consist of commutative monoids of order between 2 and 4. Note that this is not checked to full extend by the program!
	The program prints a string for each reflexive combination its finds and gives out a count of all combinations it has found in the end.*)
ReflexivityCheck[X_]:= Module[{myTable, Dims, Count, k, m, d},
						myTable = MonoidHomTable[X][[1]]; Dims = Dimensions[myTable]; Count = 0;
						For[k=1, k<=Dims[[1]], k++1,
						For[m=1, m<=Dims[[2]], m++1,
						If[Position[Range[Dims[[1]] ], myTable[[k,m]] ] == {},Continue[], (* if an entry is not a natural number, it is skipped *)
						d = myTable[[myTable[[k,m]], m]]; (* d is the dual of the dual of X w.r.t Y. *)
						If[Position[Range[Dims[[1]] ], d] == {},Continue[], (* if d is not a natural number, it is skipped *)
						If[ d == k  &&  NatHomInjCheck[X[[k]],X[[m]] ], Print[StringForm["El. `` is reflexive w.r.t. El ``.",k,m]]; Count = Count + 1, 
						If[d == k, Print[StringForm["El. `` is NOT reflexive w.r.t. El ``.",k,m]]]; (* checks if a monoid M is isomorphic to M'' without the x \mapsto L_x being an iso. *)
						];
						]] ]]; Print[StringForm["In total, `` reflexive combinations were found.",Count]];
						]


(* This program checks for a list of commutative monoids of order 3, denoted by L, and a list of monoids with absorbing el. of order 3, denoted by M, which combinations yield a semiring, where the el. from M
	denotes the addition and the el. from M denotes the multiplication.
	All monoids have to have their neutral el. denoted as 1 and their absorbing el. as 3.
	The output is a table marking the combinations that yield a semiring. *)
SemiringCheckOrder3[L_,M_]:= Module[{n,m,myTable,Count,k,l,f,g,H,PLUS,MAL,r,i,j,HELP1,HELP2,good,s,t},
					n = Dimensions[L][[1]]; m = Dimensions[M][[1]]; myTable = ConstantArray[0, {n, m}];
					For[k=1, k<=n, k++1, (* checks for the first list *)
					If[Dimensions[L[[k]]]!={3,3}, Print[StringForm["El. `` of 1. list is not a 3x3 matrix.",k]];Abort[]];
					For[l=1, l<=3,l++1,
					If[L[[k,1,l]]!=l, Print[StringForm["El. `` of 1. list has 1 not as its neutral el.",k]];Abort[]];
					If[L[[k,l,1]]!=l, Print[StringForm["El. `` of 1. list has 1 not as its neutral el.",k]];Abort[]];
					]];
					For[k=1, k<=m, k++1, (* checks for the second list *)
					If[Dimensions[L[[k]]]!={3,3}, Print[StringForm["El. `` of 2. list is not a 3x3 matrix.",k]];Abort[]];
					For[l=1, l<=3,l++1,
					If[M[[k,1,l]]!=l, Print[StringForm["El. `` of 2. list has 1 not as its neutral el.",k]];Abort[]];
					If[M[[k,l,1]]!=l, Print[StringForm["El. `` of 2. list has 1 not as its neutral el.",k]];Abort[]];
					If[M[[k,3,l]]!=3, Print[StringForm["El. `` of 2. list has 3 not as its abs. el.",k]];Abort[]];
					If[M[[k,l,3]]!=3, Print[StringForm["El. `` of 2. list has 3 not as its abs. el.",k]];Abort[]];
					]];
					f = {3,1,2}; g = {2,3,1}; (* f is the permutation to write the el. of M in "product form" with the 0 first; g is its inverse map *)
					H = {{1,2,3},{1,3,2}}; (* list of permutations that can be applied to el. of L leaving the neutral el. in the first position *)
					For[k=1, k<=n, k++1,
					For[l=1, l<=m, l++1,
					PLUS = ConstantArray[0, {3, 3}]; MAL = ConstantArray[0, {3, 3}]; HELP1 = ConstantArray[0, {3, 3}]; HELP2 = ConstantArray[0, {3, 3}];
					For[r=1, r<=2, r++1, (* r runs through all possible "versions" of the add. monoid, i.e. through the el. of H *)
					For[i=1, i<=3, i++1,
					For[j=1, j<=3, j++1, HELP2[[i,j]]=M[[l,f[[i]],f[[j]] ]]; MAL[[i,j]] = g[[HELP2[[i,j]] ]];  (* exchanges the roles of "1" and "3" *)
					HELP1[[i,j]]=L[[k,H[[r,i]],H[[r,j]] ]]; PLUS[[i,j]] = H[[r,HELP1[[i,j]] ]]; (* applies el. of H to el. of L *)
					]];
					good = True; (* check parameter for distributivity *)
					For[s=2, s<=3, s++1,
					For[t=2, t<=3, t++1,
					If[Or[MAL[[3,PLUS[[s,t]] ]]!= PLUS[[MAL[[3,s]],MAL[[3,t]] ]], MAL[[PLUS[[s,t]],3 ]]!= PLUS[[MAL[[s,3]],MAL[[t,3]] ]] ], good = False; Break[]]; (*  distributivity check *)
					];
					If[good == False, Break[] ]]; (* if check fails break out of s-loop, i.e. end checks *)
					If[good == True, myTable[[k,l]] = "SR";Print[StringForm["+: `` with x: `` semiring.",PLUS,MAL]](*; Break[]*)]; (* if all checks are true, break out of r-loop *)
					];
					]]; 
					Return[myTable//MatrixForm]]


(* This program checks for a list of commutative monoids of order 4, denoted by L, and a list of monoids with absorbing el. of order 4, denoted by M, which combinations yield a semiring, where the el. from M
	denotes the addition and the el. from M denotes the multiplication.
	All monoids have to have their neutral el. denoted as 1 and their absorbing el. as 4.
	The output is a table marking the combinations that yield a semiring. *)
SemiringCheckOrder4[L_,M_]:= Module[{n,m,myTable,Count,k,l,F,G,H,J,PLUS,MAL,r,i,j,HELP1,HELP2,good,s,t,u,v},
					n = Dimensions[L][[1]]; m = Dimensions[M][[1]]; myTable = ConstantArray[0, {n, m}];
					For[k=1, k<=n, k++1, (* checks for the first list *)
					If[Dimensions[L[[k]]]!={4,4}, Print[StringForm["El. `` of 1. list is not a 4x4 matrix.",k]];Abort[]];
					For[l=1, l<=4,l++1,
					If[L[[k,1,l]]!=l, Print[StringForm["El. `` of 1. list has 1 not as its neutral el.",k]];Abort[]];
					If[L[[k,l,1]]!=l, Print[StringForm["El. `` of 1. list has 1 not as its neutral el.",k]];Abort[]];
					]];
					For[k=1, k<=m, k++1, (* checks for the second list *)
					If[Dimensions[L[[k]]]!={4,4}, Print[StringForm["El. `` of 2. list  is not a 3x3 matrix.",k]];Abort[]];
					For[l=1, l<=4,l++1,
					If[M[[k,1,l]]!=l, Print[StringForm["El. `` of 2. list has 1 not as its neutral el.",k]];Abort[]];
					If[M[[k,l,1]]!=l, Print[StringForm["El. `` of 2. list has 1 not as its neutral el.",k]];Abort[]];
					If[M[[k,4,l]]!=4, Print[StringForm["El. `` of 2. list has 4 not as its abs. el.",k]];Abort[]];
					If[M[[k,l,4]]!=4, Print[StringForm["El. `` of 2. list has 4 not as its abs. el.",k]];Abort[]];
					]];
					F = {{4,1,2,3},{4,1,3,2}}; (* list of permutation to write the el. of M in "product form" with the 0 first *)
					G = {{2,3,4,1},{2,4,3,1}}; (* inverse maps of F *)
					H = {{1,2,3,4},{1,2,4,3},{1,3,2,4},{1,3,4,2},{1,4,2,3},{1,4,3,2}};  (* list of permutations that can be applied to el. of L leaving the neutral el. in the first position *)
					J = {{1,2,3,4},{1,2,4,3},{1,3,2,4},{1,4,2,3},{1,3,4,2},{1,4,3,2}};  (* inverse maps of H *)
					For[k=1, k<=n, k++1,
					For[l=1, l<=m, l++1,
					PLUS = ConstantArray[0, {4, 4}]; MAL = ConstantArray[0, {4, 4}]; HELP1 = ConstantArray[0, {4, 4}]; HELP2 = ConstantArray[0, {4, 4}];
					For[u=1, u<=2, u++1,
					For[r=1, r<=5, r++1, (* r runs through all possible "versions" of the add. monoid, i.e. through the el. of H *)
					For[i=1, i<=4, i++1,
					For[j=1, j<=4, j++1, HELP2[[i,j]]=M[[l,F[[u,i]],F[[u,j]] ]]; MAL[[i,j]] = G[[u,HELP2[[i,j]] ]];  (* exchanges the roles of "1" and "3" *)
					HELP1[[i,j]]=L[[k,H[[r,i]],H[[r,j]] ]]; PLUS[[i,j]] = J[[r,HELP1[[i,j]] ]]; (* applies el. of H to el. of L *)
					]];
					good = True; (* check parameter for distributivity *)
					For[v=3, v<=4, v++1,
					For[s=2, s<=4, s++1,
					For[t=2, t<=4, t++1,
					If[Or[MAL[[v,PLUS[[s,t]] ]]!= PLUS[[MAL[[v,s]],MAL[[v,t]] ]], MAL[[PLUS[[s,t]],v ]]!= PLUS[[MAL[[s,v]],MAL[[t,v]] ]] ], good = False; Break[]]; (*  distributivity check *)
					];
					If[good == False, Break[] ] ]; (* if check fails break out of s-loop *)
					If[good == False, Break[] ] ]; (* if check fails break also out of v-loop, i.e. end checks *)
					If[good == True, myTable[[k,l]] = "SR";Print[StringForm["+: `` with x: `` semiring.",PLUS,MAL]]; (*Break[]*)]; (* if all checks are true, break out of r-loop *)
					]; (* end of r-loop *)
					(*If[myTable[[k,l]] \[Equal] "SR", Break[]] *)]; (* end of u-loop *)
					]]; (* end of k- and l-loop *)
					Return[myTable//MatrixForm]]


DistributivityCheck[A_,M_]:= Module[{good,n,v,s,t},
							n = Dimensions[A][[1]]; good = True;
							For[v=1, v<=n, v++1,
							For[s=1, s<=n, s++1,
							For[t=1, t<=n, t++1,
							If[Or[M[[v,A[[s,t]] ]]!= A[[M[[v,s]],M[[v,t]] ]], M[[A[[s,t]],v ]]!= A[[M[[s,v]],M[[t,v]] ]] ], good = False; Break[]]; (*  distributivity check *)
							]]];
							Return[good] ]
