Vtyjkyuk

I am very interested in using mutually orthogonal latin squares (MOLS) to reduce the number of test cases but I struggle to find a way how to implement the algorithm. In an article published in a foreign language, I did find a formula including the source code.

This is said to be a formula to construct a series of squares (if $n$ is a prime number):

For $k, i, j in 1,2, cdots, n$
$$A_k(i, j) = [k cdot (i-1) + (j-1)] mod n$$

Is this correct? I tried to find other sources but failed.

If this is correct, what is meant by $k$, $i$ and $j$? In the code, they pass "size" of the latin square to all of these variables.

The info given in the OP is almost correct. As noted in Jyrki Lahtonen's answer $k$ cannot be congruent to $n$.

In $A_k(i,j)$, $k$ is the index of the square: the algorithm generates $n-1$ MOLS, one for each $k$ in $1,cdots,n-1$. Using normal matrix conventions, $i$ is the row index and $j$ is the column index, but the formula still works if those are swapped.

Here's some rudimentary Python code that uses that formula to generate a set of mutually orthogonal latin squares and then tests each combination for orthogonality. This code was developed using Python 2.6 but it should run correctly on Python 3.

In this code $i$ and $j$ run from 0 to $n-1$ and $k$ runs from 1 to $n-1$; note that $n$ must be prime.

from __future__ import print_function

def show_grid(g):
 for row in g:
 print(row)
 print()

def test_mols(g):
 ''' check that all entries in g are unique '''
 size = len(g) ** 2
 a = set()
 for row in g:
 a.update(row)
 return len(a) == size

def mols(n):
 ''' Generate a set of mutually orthogonal latin squares 
 n must be prime
 ''' 
 r = range(n)

 #Generate each Latin square
 allgrids = 
 for k in range(1, n):
 grid = 
 for i in r:
 row = 
 for j in r:
 a = (k*i + j) % n
 row.append(a)
 grid.append(row)
 allgrids.append(grid)

 for g in allgrids:
 show_grid(g)

 print('- ' * 20 + 'n')

 #Combine the squares to show their orthoganility
 m = len(allgrids)
 for i in range(m):
 g0 = allgrids[i]
 for j in range(i+1, m):
 g1 = allgrids[j]
 newgrid = 
 for r0, r1 in zip(g0, g1):
 newgrid.append(zip(r0, r1))
 print(test_mols(newgrid))
 show_grid(newgrid)

mols(3)

output for n=3

[0, 1, 2]
[1, 2, 0]
[2, 0, 1]

[0, 1, 2]
[2, 0, 1]
[1, 2, 0]

- - - - - - - - - - - - - - - - - - - - 

True
[(0, 0), (1, 1), (2, 2)]
[(1, 2), (2, 0), (0, 1)]
[(2, 1), (0, 2), (1, 0)]

output for n=5

[0, 1, 2, 3, 4]
[1, 2, 3, 4, 0]
[2, 3, 4, 0, 1]
[3, 4, 0, 1, 2]
[4, 0, 1, 2, 3]

[0, 1, 2, 3, 4]
[2, 3, 4, 0, 1]
[4, 0, 1, 2, 3]
[1, 2, 3, 4, 0]
[3, 4, 0, 1, 2]

[0, 1, 2, 3, 4]
[3, 4, 0, 1, 2]
[1, 2, 3, 4, 0]
[4, 0, 1, 2, 3]
[2, 3, 4, 0, 1]

[0, 1, 2, 3, 4]
[4, 0, 1, 2, 3]
[3, 4, 0, 1, 2]
[2, 3, 4, 0, 1]
[1, 2, 3, 4, 0]

- - - - - - - - - - - - - - - - - - - - 

True
[(0, 0), (1, 1), (2, 2), (3, 3), (4, 4)]
[(1, 2), (2, 3), (3, 4), (4, 0), (0, 1)]
[(2, 4), (3, 0), (4, 1), (0, 2), (1, 3)]
[(3, 1), (4, 2), (0, 3), (1, 4), (2, 0)]
[(4, 3), (0, 4), (1, 0), (2, 1), (3, 2)]

True
[(0, 0), (1, 1), (2, 2), (3, 3), (4, 4)]
[(1, 3), (2, 4), (3, 0), (4, 1), (0, 2)]
[(2, 1), (3, 2), (4, 3), (0, 4), (1, 0)]
[(3, 4), (4, 0), (0, 1), (1, 2), (2, 3)]
[(4, 2), (0, 3), (1, 4), (2, 0), (3, 1)]

True
[(0, 0), (1, 1), (2, 2), (3, 3), (4, 4)]
[(1, 4), (2, 0), (3, 1), (4, 2), (0, 3)]
[(2, 3), (3, 4), (4, 0), (0, 1), (1, 2)]
[(3, 2), (4, 3), (0, 4), (1, 0), (2, 1)]
[(4, 1), (0, 2), (1, 3), (2, 4), (3, 0)]

True
[(0, 0), (1, 1), (2, 2), (3, 3), (4, 4)]
[(2, 3), (3, 4), (4, 0), (0, 1), (1, 2)]
[(4, 1), (0, 2), (1, 3), (2, 4), (3, 0)]
[(1, 4), (2, 0), (3, 1), (4, 2), (0, 3)]
[(3, 2), (4, 3), (0, 4), (1, 0), (2, 1)]

True
[(0, 0), (1, 1), (2, 2), (3, 3), (4, 4)]
[(2, 4), (3, 0), (4, 1), (0, 2), (1, 3)]
[(4, 3), (0, 4), (1, 0), (2, 1), (3, 2)]
[(1, 2), (2, 3), (3, 4), (4, 0), (0, 1)]
[(3, 1), (4, 2), (0, 3), (1, 4), (2, 0)]

True
[(0, 0), (1, 1), (2, 2), (3, 3), (4, 4)]
[(3, 4), (4, 0), (0, 1), (1, 2), (2, 3)]
[(1, 3), (2, 4), (3, 0), (4, 1), (0, 2)]
[(4, 2), (0, 3), (1, 4), (2, 0), (3, 1)]
[(2, 1), (3, 2), (4, 3), (0, 4), (1, 0)]

Yes, that formula works. If $n=p$ is a prime it gives you a set of $p-1$ MOLS - one for each $k=1,2,ldots,p-1$. The variables $i$ and $j$ stand for the row and column index respectively.

If $n=p^m$ is a power of a prime $p$, then there exists a similar formula for $n-1$ MOLS, but that requires you to do arithmetic in an extension field of integers modulo $p$.