🧩 Constraint Solving POTD:Latin Squares & Orthogonality — Classic CSP Meets Combinatorial Design

[github-actions[bot]]

Today's problem sits at the intersection of classical combinatorics and modern constraint solving: Latin Squares. Elegant in structure, surprisingly deep in difficulty, and ubiquitous in applications from experimental design to cryptography.

---

Problem Statement

A Latin square of order n is an n × n grid filled with n distinct symbols such that each symbol appears exactly once in every row and exactly once in every column.

Small Concrete Instance (n = 4)

1  2  3  4
2  1  4  3
3  4  1  2
4  3  2  1
```

Every row and every column contains each of `{1, 2, 3, 4}` exactly once. ✓

**Input:** Integer `n` (order of the square), optionally a set of pre-assigned cells (a *partial* Latin square).

**Output:** A completed `n × n` grid satisfying the Latin square property, or `UNSAT` if no completion exists.

#### Extension: Mutually Orthogonal Latin Squares (MOLS)

Two Latin squares `L1` and `L2` of order `n` are **orthogonal** if, when superimposed, every ordered pair `(L1[i,j], L2[i,j])` appears exactly once across all `n²` cells.

```
L1:          L2:          Superimposed:
1  2  3  4   1  2  3  4   (1,1) (2,2) (3,3) (4,4)
2  1  4  3   3  4  1  2   (2,3) (1,4) (4,1) (3,2)
3  4  1  2   4  3  2  1   (3,4) (4,3) (1,2) (2,1)
4  3  2  1   2  1  4  3   (4,2) (3,1) (2,4) (1,3)
```

All 16 pairs `(a,b)` with `a,b ∈ {1..4}` appear — orthogonality holds. ✓

Finding a **maximum set of MOLS** of order `n` is an open problem for several values of `n`.

---

### Why It Matters

**Experimental design.** A Latin square of order `n` encodes a balanced incomplete block design: `n` treatments, `n` blocks, `n` time periods, with every treatment appearing once per block and once per period — eliminating two sources of bias simultaneously. Orthogonal Latin squares (called *Graeco-Latin squares*) eliminate a third factor.

**Error-correcting codes.** Sets of MOLS construct optimal *orthogonal arrays*, which are the combinatorial backbone of Reed-Solomon codes and other algebraic codes used in QR codes, CDs, and satellite communications.

**Cryptography.** Latin squares are the mathematical abstraction behind the *S-boxes* of many symmetric ciphers. Highly non-linear Latin squares resist differential and linear cryptanalysis.

**Sudoku puzzles.** A 9×9 Sudoku is a Latin square with the additional constraint that each of the nine 3×3 boxes also contains each digit exactly once.

---

### Modeling Approaches

#### Approach 1: Constraint Programming (CP)

**Decision variables:** `x[i,j] ∈ {1..n}` for each cell `(i,j)`.

**Constraints:**
```
for each row i:    AllDifferent(x[i,1], x[i,2], ..., x[i,n])
for each column j: AllDifferent(x[1,j], x[2,j], ..., x[n,j])
```

**Symmetry breaking** (optional, to prune the search space):
```
x[1,j] = j  for j = 1..n   (fix first row to 1,2,...,n)
x[i,1] = i  for i = 1..n   (fix first column to 1,2,...,n → reduced form)
```

This yields a **normalized / reduced Latin square**. Every equivalence class has exactly one representative — reducing the search space by a factor of `n! × (n−1)!`.

**Trade-offs:** `AllDifferent` propagates strongly via Hall's theorem (arc consistency in O(n^1.5)). CP handles partial Latin squares (pre-assigned cells) natively and scales well to `n ≈ 30` for completion problems.

---

#### Approach 2: Mixed Integer Programming (MIP)

Introduce binary variables `y[i,j,k] = 1` iff cell `(i,j)` contains symbol `k`.

**Constraints:**
```
∑_k y[i,j,k] = 1       ∀ i,j    (each cell has exactly one symbol)
∑_j y[i,j,k] = 1       ∀ i,k    (each symbol once per row)
∑_i y[i,j,k] = 1       ∀ j,k    (each symbol once per column)
y[i,j,k] ∈ {0,1}

Variables: n³ binary variables. Constraints: n² + 2n² equalities.

Trade-offs: The LP relaxation of this formulation is integral — the constraint matrix is totally unimodular — so MIP solvers often solve it without branching. However, n³ variables become expensive for n > 50. MIP shines when combined with side constraints (e.g., minimize some objective over a set of Latin squares).

---

Approach 3: SAT Encoding

Variables: p[i,j,k] = "cell (i,j) contains symbol k" (Boolean).

Constraints:

• At-least-one per cell: p[i,j,1] ∨ p[i,j,2] ∨ … ∨ p[i,j,n]
• At-most-one per cell: ¬p[i,j,k] ∨ ¬p[i,j,k'] for all k ≠ k'
• Each symbol once per row: pairwise ¬p[i,j,k] ∨ ¬p[i,j',k] for j ≠ j'
• Each symbol once per column: pairwise ¬p[i,j,k] ∨ ¬p[i',j,k] for i ≠ i'

At-most-one encoding grows as O(n²) clauses per cell; commander-variable or product encodings keep it O(n log n).

Trade-offs: SAT excels at deciding satisfiability (e.g., "does a completion exist?") and scales to moderate n with CDCL solvers. Unit propagation naturally performs row/column arc consistency. Less natural for optimization variants.

---

Example Model (MiniZinc CP)

include "alldifferent.mzn";

int: n = 5;
array[1..n, 1..n] of var 1..n: x;

% Latin square constraints
constraint forall(i in 1..n)(alldifferent(row(x, i)));
constraint forall(j in 1..n)(alldifferent(col(x, j)));

% Symmetry breaking: fix first row and first column (reduced form)
constraint forall(j in 1..n)(x[1,j] = j);
constraint forall(i in 1..n)(x[i,1] = i);

solve satisfy;

output [ show(x[i,j]) ++ if j = n then "\n" else "  " endif
       | i in 1..n, j in 1..n ];

---

Key Techniques

1. AllDifferent Propagation via Hall's Theorem

The AllDifferent global constraint is the workhorse here. Arc-consistency for AllDifferent is achieved by detecting Hall sets: subsets S of variables whose union of domains |∪D(x_i, i∈S)| = |S| — those values cannot appear elsewhere. Regin's algorithm achieves full arc-consistency in O(n^1.5) using bipartite matching.

Bounds consistency is cheaper but weaker. For Latin squares with small domains, full arc-consistency pays off.

2. Symmetry Breaking and Isotopy Classes

Latin squares have massive symmetry: permuting rows, columns, or symbols yields another valid square. The autotopy group of a Latin square of order n can be enormous. Effective symmetry breaking strategies:

• Static breaking: fix first row and/or column (reduced form), reducing search by n! × (n−1)!
• Dynamic breaking: lex_lesseq constraints between rows or columns
• NAUTY / symmetry detection: compute canonical forms to detect isotopy equivalence

Without symmetry breaking, solvers waste enormous effort on redundant solutions.

3. Conflict-Driven Backjumping for Partial Completions

Deciding whether a partial Latin square can be completed is NP-complete (proved by Colbourn, 1984). The hardest instances are partially filled squares near 50–60% density. Conflict-driven clause learning (CDCL) in SAT solvers, or Nogood recording in CP solvers, dramatically prunes backtracking by learning from failed assignments — turning exponential search into something tractable in practice.

---

Challenge Corner

🔍 Open question: The maximum number of MOLS of order n is denoted N(n). It is known that:

• N(n) ≤ n − 1 for all n (the projective plane bound)
• N(p^k) = p^k − 1 when n is a prime power
• N(6) = 1 (Euler's 36-officers problem, no pair of MOLS exists!)
• N(10) is still unknown — believed to be 2, but unproven

Can you write a constraint model to search for three MOLS of order 10 — or prove none exist? What symmetry-breaking constraints would make the search tractable? (Hint: the space has 10^100 assignments — you'll need very strong pruning!)

As a warm-up: verify that N(4) = 3 by finding three mutually orthogonal Latin squares of order 4.

---

References

1. J. Dénes & A. D. Keedwell — Latin Squares and Their Applications (1974). The definitive reference on Latin square theory.

2. C. J. Colbourn — "The complexity of completing partial Latin squares," Discrete Applied Mathematics 8(1):25–30, 1984. Proves NP-completeness of the completion problem.

3. I. P. Gent, C. Jefferson, I. Miguel — "Minion: A Fast Scalable Constraint Solver," ECAI 2006. Demonstrates efficient AllDifferent propagation in the context of Latin squares and Sudoku.

4. Handbook of Constraint Programming, Ch. 4 (Global Constraints) — Rossi, van Beek & Walsh (eds.), Elsevier 2006. Covers Hall's theorem and Regin's algorithm in depth.

Generated by Constraint Solving — Problem of the Day · ● 178.2K · ◷

• expires on Apr 15, 2026, 12:56 PM UTC