As far as I know, this is still an open problem. It was posed by Keedwell in "Some problems concerning complete Latin squares" (1974).
Observation: If L and M are orthogonal complete Latin squares, then we can permute the symbols in L (or M) without affecting the orthogonality property, nor the completeness property. Hence, we can assume that the first row of L and M are in order (i.e. "normalised"). I wrote some code to check up to $8 \times 8$ squares, and found:
\[
\begin{array}{|r|rrrrrrrr|}
\hline
n & 2 & 3 & 4 & 5 & 6 & 7 & 8 & 9 \\\\
\text{nr. normalised complete Latin squares} & 1 & 0 & 2 & 0 & 8 & 0 & 192 & 0 \\\\
\hline
\end{array}
\]
None of the Latin squares encountered were orthogonal.
My guess is they exist, as there doesn't seem to be any reason for them not to exist. I suspect they don't exist for small orders simply because there's too many conditions placed on them.
Here's my C code below:
#include
/* change this to the order of your latin square */
#define ORDER_LATIN_SQUARE 9
#define MAX_NR_LATIN_SQUARES 100000
int n;
int latin[ORDER_LATIN_SQUARE][ORDER_LATIN_SQUARE];
int latin2[ORDER_LATIN_SQUARE][ORDER_LATIN_SQUARE];
int complete_latin_squares[MAX_NR_LATIN_SQUARES][ORDER_LATIN_SQUARE][ORDER_LATIN_SQUARE];
int nr_complete_latin_squares=0;
int nr_orthogonal_pairs_of_complete_latin_squares=0;
void print_latin_square() {
for(int i=0;i0 && i==1 && j==0) { printf("%f percent complete.\n",(float) (s-1)/(n-2)*100); }
for(int c=0;c0) {
int t=latin[i][j-1];
for(int r=0;r0) {
int t=latin[i-1][j];
for(int r=0;r
Here's a GAP version of the same program (it's slower, but it's nice to check the results are correct):
SetOfCompleteLatinSquares:=[];
GenerateCompleteLatinSquaresBacktrackingAlgorithm:=function(L,i,j,n)
local s,k,S1,S2,try_symbol;
# S1 is the set of extant pairs (l[r][c],l[r][c+1])
S1:=List([1..i-1],r->List([1..n-1],c->[L[r][c],L[r][c+1]]));
S1:=Concatenation(S1);
S1:=Concatenation([S1,List([1..j-2],c->[L[i][c],L[i][c+1]])]);
# S2 is the set of extant pairs (l[r+1][c],l[r][c])
S2:=List([2..i-1],r->List([1..n],c->[L[r-1][c],L[r][c]]));
S2:=Concatenation(S2);
if(i>1) then S2:=Concatenation([S2,List([1..j-1],c->[L[i-1][c],L[i][c]])]); fi;
for s in [1..n] do
try_symbol:=true;
for k in [1..j-1] do if s=L[i][k] then try_symbol:=false; fi; od; # row clash
for k in [1..i-1] do if s=L[k][j] then try_symbol:=false; fi; od; # column clash
if(j>1) then if([L[i][j-1],s] in S1) then try_symbol:=false; fi; fi; # row pair clash
if(i>1) then if([L[i-1][j],s] in S2) then try_symbol:=false; fi; fi; # column pair clash
if(try_symbol) then
L[i][j]:=s;
if(i=n and j=n) then
# Print(L,"\n");
Append(SetOfCompleteLatinSquares,[Immutable(L)]);
elif(j=n) then
GenerateCompleteLatinSquaresBacktrackingAlgorithm(L,i+1,1,n);
else
GenerateCompleteLatinSquaresBacktrackingAlgorithm(L,i,j+1,n);
fi;
L[i][j]:=0;
fi;
od;
end;;
GenerateCompleteLatinSquares:=function(n)
local L,j;
L:=List([1..n],i->List([1..n],j->0));
for j in [1..n] do L[1][j]:=j; od; # assumes the first row is 1,2,..,n
SetOfCompleteLatinSquares:=[];
GenerateCompleteLatinSquaresBacktrackingAlgorithm(L,2,1,n);
return SetOfCompleteLatinSquares;
end;;