Orthogonal Mate of a Latin Square

This code requires the Exact Cover package.

```% A Latin square of order N is an NxN array filled with N different symbols in
% such a way that each symbol occurs exactly once in each row and exactly once
% in each column.

% Two Latin squares of the same size are said to be orthogonal if no pair of
% corresponding elements occurs more than once.

% Given two orthogonal Latin squares, each is said to be an orthogonal mate of
% the other.

% A transversal of a Latin square of order N is a set of N cells, one from each
% row and one from each column, and each cell containing a different symbol.

% A Latin square of order N has an orthogonal mate if and only if it possesses
% a set of N transversals which partition the N^2 cells of the square.

/*
* Finds an orthogonal mate of a given Latin square using Exact Cover twice,
*  first of all to find all transversals of the Latin square, then to use
*  these transversals to construct the orthogonal mate.
*/

% Constructs an orthogonal mate from a given Latin square.  It is required, but
%   not tested, that the symbols used in a Latin square of order N be the
%   integers 1 to N.  On backtracking, all solutions will be found.
% e.g. findall(M, mate([1,2,3,4,4,3,2,1,2,1,4,3,3,4,1,2], M), Mates) gives
%   Mates=[[1,2,3,4,2,1,4,3,3,4,1,2,4,3,2,1],[1,2,3,4,3,4,1,2,4,3,2,1,2,1,4,3]]
mate(LatinSquare, Mate):-
findall(T, transversal(LatinSquare, T), Transversals),
mate_1(Transversals, LatinSquare, Mate).

% Constructs an orthogonal mate from a given Latin square and its transversals.
%   On backtracking, all solutions will be found.
% e.g. L=[1,2,3,3,1,2,2,3,1], findall(T,transversal(L,T),Ts), mate_1(Ts,L,Mate)
%   gives Mate=[1,2,3,2,3,1,3,1,2]
mate_1(Transversals, LatinSquare, Mate):-
length(LatinSquare, N2),
N is sqrt(N2),
N =:= ip(N),
TwoN2 is 2 * N2,
findall(X, for(1, TwoN2, X), Set),
subsets(Transversals, LatinSquare, N, -1, Subsets),
excov(Set, Subsets, Cover),
transversals_in_cover(Cover, Transversals, CoverTransversals),
cells_and_symbols(CoverTransversals, N, [], CellsAndSymbols),
sort(CellsAndSymbols, SortedCellsAndSymbols),
symbols(SortedCellsAndSymbols, Mate).

% Returns the sorted subsets for all of the transversals. A negative index into
%   the transversals is put as a prefix for each subset so as to retrieve the
%   transversals used in the solution.
% e.g. subsets([[1,3,2],[3,2,1],[2,1,3]], [1,2,3,3,1,2,2,3,1], 3, -1, Yss)
%   gives Yss=[[-1,1,5,9,10,15,17], [-2,2,6,7,11,13,18], [-3,3,4,8,12,14,16]]
subsets([], _, _, _, []).
subsets([Transversal|Transversals], LatinSquare, N, Index, [[Index|Ys]|Yss]):-
subsets_1(Transversal, LatinSquare, 1, N, Ys0),
sort(Ys0, Ys),
Index1 is Index - 1,
subsets(Transversals, LatinSquare, N, Index1, Yss).

% Returns the unsorted subset for one of the transversals.
% e.g. subsets_1([3,2,1], [1,2,3,3,1,2,2,3,1], 1, 3, Ys) gives
%   Ys=[7,18,2,11,6,13]
subsets_1([], _, _, _, []).
subsets_1([Column|Transversal], LatinSquare, Row, N, [Y1,Y2|Ys]):-
Cell is (Row - 1) * N + Column,
nth_member(LatinSquare, Cell, Symbol), !,
Y1 is (Symbol - 1) * N + Row,
Y2 is N * N + (Symbol - 1) * N + Column,
Row1 is Row + 1,
subsets_1(Transversal, LatinSquare, Row1, N, Ys).

% Retreives the transversals corresponding to the exact cover solution.
transversals_in_cover([], _, []).
transversals_in_cover([[Index|_]|Cover], Transversals0, [T|Transversals]):-
Index1 is -Index,
nth_member(Transversals0, Index1, T), !,
transversals_in_cover(Cover, Transversals0, Transversals).

% Retrieves the (unordered) cells and symbols, for all transversals, which will
%   be in the mate.
% e.g. cells_and_symbols([[1,3,2],[2,1,3],[3,2,1]], 3, [], Xss)
%   gives Xss=[[3,3],[5,3],[7,3],[9,2],[2,2],[4,2],[6,1],[8,1],[1,1]]
cells_and_symbols([], _, Xss, Xss).
cells_and_symbols([Transversal|Transversals], N, Xss0, Xss):-
% The symbol put in this transversal is the index of 1 in the transversal,
%   that is, the column in which this transversal occurs in the top row.
nth_member(Transversal, Symbol, 1), !,
cells_and_symbols_1(Transversal, 1, N, Symbol, Xss0, Xss1),
cells_and_symbols(Transversals, N, Xss1, Xss).

% Retrieves the (unordered) cells and symbols, for one transversal, which will
%   be in the mate.
% e.g. cells_and_symbols_1([2,1,3], 1, 3, 1, [], Xs)
%   gives Xs=[[9,1],[2,1],[4,1]]
cells_and_symbols_1([], _, _, _, Xss, Xss).
cells_and_symbols_1([Row|Transversal], Column, N, Symbol, Xss0, Xss):-
Cell is (Row - 1) * N + Column,
Column1 is Column + 1,
Xss1=[[Cell,Symbol]|Xss0],
cells_and_symbols_1(Transversal, Column1, N, Symbol, Xss1, Xss).

% Retrieves the symbols which will be in the cells of the mate.
% e.g. symbols([[1,1],[2,2],[3,3],[4,2],[5,3],[6,1],[7,3],[8,1],[9,2]], Mate)
%   gives Mate=[1,2,3,2,3,1,3,1,2]
symbols([], []).
symbols([[_,Symbol]|CellsAndSymbols], [Symbol|Mate]):-
symbols(CellsAndSymbols, Mate).

/*
* Finds a transversal of a given Latin square using Exact Cover.
*/

% Finds a transversal of a given Latin square.  On backtracking, all solutions
%   will be found.  It is required, but not tested, that the symbols used in a
%   Latin square of order N be the integers 1 to N.
% e.g. transversal([1,2,3,4,5,2,1,4,5,3,3,5,1,2,4,4,3,5,1,2,5,4,2,3,1], T)
%   gives T=[1,3,2,5,4], T=[1,4,5,2,3] and T=s=[1,5,4,3,2]
%   where T=[1,3,2,5,4] indicates that the cells forming this transversal are
%   (row,column) = (1,1), (3,2), (2,3), (5,4) and (4,5).
transversal(LatinSquare, Transversal):-
length(LatinSquare, N2),
N is sqrt(N2),
N =:= ip(N),
findall(X, transversal_1(N, X), Set),
transversal_2(LatinSquare, 1, 1, N, Subsets),
excov(Set, Subsets, Cover),
sort(Cover, Cover1),
transversal_3(Cover1, N, Transversal).

% Generates an element of the set.  On backtracking, generates all elements.
transversal_1(N, H):-for(1, N, C), H is C.     % Column
transversal_1(N, H):-for(1, N, R), H is N+R.   % Row
transversal_1(N, H):-for(1, N, S), H is N+N+S. % Symbol

% Generates all subsets.
transversal_2([], _, _, _, []):-!.
transversal_2(LatinSquare, I, J, N, Subsets):-
J > N, !,         % End of row
I1 is I + 1,      % Next row
transversal_2(LatinSquare, I1, 1, N, Subsets).
transversal_2([X|LatinSquare], I, J, N, [[J,NI,NNX]|Subsets]):-
NI is N + I,      % Row
NNX is N + N + X, % Symbol
J1 is J + 1,      % Next column
transversal_2(LatinSquare, I, J1, N, Subsets).

% Interprets the output from excov/3.
transversal_3([], _, []).
transversal_3([[_,Z,_]|Cover], N, [T|Transversal]):-
T is Z - N,
transversal_3(Cover, N, Transversal).

/*
* Finds an orthogonal mate of a given Latin square directly using Exact Cover.
*/

% naive(LatinSquare, Mate) is true if Mate is an orthogonal mate of the
%   LatinSquare.  It is required, but not tested, that the symbols used in a
%   Latin square of order N be the integers 1 to N.  On backtracking, all
%   solutions are found.
%   This method is more straightforward than mate/2 as it doesn't use
%   transversals, but is much slower.
% e.g. findall(M, naive([1,2,3,4,4,3,2,1,2,1,4,3,3,4,1,2], M), Mates) gives
%   Mates=[[1,2,3,4,2,1,4,3,3,4,1,2,4,3,2,1],[1,2,3,4,3,4,1,2,4,3,2,1,2,1,4,3]]
naive(LatinSquare, Mate):-
length(LatinSquare, N2),
N is sqrt(N2),
N =:= ip(N),
N_1 is N - 1,
findall(X, naive_1(N, N_1, N2, X), Set),
findall(Ys, naive_2(LatinSquare, N, N_1, N2, Ys), Subsets),
excov(Set, Subsets, Cover),
sort(Cover, Cover1),
naive_4(Cover1, N, Mate).

% Generates an element of the set.  On backtracking, generates all elements.
naive_1(N, N_1, _,  H):-for(0, N_1, R), for(0, N_1, C), H is        N*R + C.
naive_1(N, N_1, N2, H):-for(0, N_1, V), for(0, N_1, R), H is   N2 + N*V + R.
naive_1(N, N_1, N2, H):-for(0, N_1, V), for(0, N_1, C), H is 2*N2 + N*V + C.
naive_1(N, N_1, N2, H):-for(0, N_1, V), for(0, N_1, W), H is 3*N2 + N*V + W.

% Generates a subset for one value in one cell.  On backtracking generates the
%   subsets for all values in all cells.
naive_2(LatinSquare, N, N_1, N2, [RC,VR,VC,VW]):- % Subsets for the top row
% For the top row we want R=0 and V=C
for(0, N_1, C),
C1 is C + 1,
nth_member(LatinSquare, C1, W),
RC is C,
VR is N2 + N*C,
VC is 2*N2 + N*C + C,
VW is 3*N2 + N*C + W-1.
naive_2(LatinSquare, N, N_1, N2, [RC,VR,VC,VW]):- % Subsets for the other rows
for(1, N_1, R),
for(0, N_1, C),
I is R*N + C + 1,
nth_member(LatinSquare, I, W),
for(0, N_1, V),
RC is N*R + C,
VR is N2+ N*V + R,
VC is 2*N2 + N*V + C,
VW is 3*N2 + N*V + W-1.

% Interprets the output from excov/3.
naive_4([], _, []).
naive_4([[_,VR|_]|Cover], N, [X|Mate]):-
X is (VR - N*N) // N + 1,
naive_4(Cover, N, Mate).

% for(I, J, K) is true if K is an integer between I and J inclusive.
for(I, J, I):-I =< J.
for(I, J, K):-I < J, I1 is I + 1, for(I1, J, K).

% length(Xs, L) is true if L is the number of elements in the list Xs.
%length(Xs, L):-length_1(Xs, 0, L).

%length_1([], L, L).
%length_1([_|Xs], L0, L):-L1 is L0 + 1, length_1(Xs, L1, L).

% member(X, Xs) is true if the element X is contained in the list Xs.
%member(X, [X|_]).
%member(X, [_|Xs]):-member(X, Xs).

% nth_member(+Xs, ?N, ?X) is true if X is the N-th (base 1) element of the
%   list Xs.
nth_member(Xs, N, X):-nth_member_1(Xs, X, 1, N).

nth_member_1([X|_], X, I, I).
nth_member_1([_|Xs], X, I0, I):-
I1 is I0 + 1,
nth_member_1(Xs, X, I1, I).
```