## Tuesday, December 22, 2020

### Numbrix puzzle via MIP

The puzzle of numbrix [1,2] is simple to state.

• We have a 9x9 board with 81 cells.
• We need to fill each cell with a unique integer value between 1 and 81.
• There are some prefilled cells.
• We should be able to start with the cell containing value 1 and then following a path (only going one up, down, left or right) to the next values and thus eventually arriving at the cell with value 81.

#### Example

Below is an example of an input (the "given" cells) and a solution. The colors may help to see that we have a path from the cell with value 1 to the cell with value 81.

Of course, we want to see if we can create a MIP model for this.

Like in Sudoku models we use a binary variable $$\color{darkred}x_{i,j,k} \in \{0,1\}$$ defined by $\color{darkred}x_{i,j,k} = \begin{cases} 1 & \text{if cell (i,j) has value k} \\ 0 & \text{otherwise}\end{cases}$ This introduces a large number of binary variables ($$9 \times 9 \times 81=6,561$$), but it makes the problem much simpler when formulating as a MIP model. This is really the big decision: once you have the main variables defined, things follow from there.

We need to formulate just three constraints:

1. Each cell can contain just one value. $\sum_k \color{darkred}x_{i,j,k} = 1 \>\> \forall i,j$
2. Each value should occur exactly once. $\sum_{i,j}\color{darkred}x_{i,j,k} = 1 \>\> \forall k$
3. We need: if $$\color{darkred}x_{i,j,k} = 1$$ then one of its neighbors should contain value $$k+1$$. A neighbor is defined as going one up, down, left or right. This can be formulated as an implication $\color{darkred}x_{i,j,k} = 1 \Rightarrow \color{darkred}x_{i-1,j,k+1}+\color{darkred}x_{i+1,j,k+1}+\color{darkred}x_{i,j-1,k+1}+\color{darkred}x_{i,j+1,k+1} = 1$ Or, as inequality: $\color{darkred}x_{i-1,j,k+1}+\color{darkred}x_{i+1,j,k+1}+\color{darkred}x_{i,j-1,k+1}+\color{darkred}x_{i,j+1,k+1} \ge \color{darkred}x_{i,j,k}$

This is all!

Note that I don't explicitly form a path in the model. Just the constraint about "next k should be a neighbor" is enough. So, our complete model looks like:

MIP Model
\begin{align} & \sum_k \color{darkred}x_{i,j,k} = 1 && \forall i,j \\ & \sum_{i,j}\color{darkred}x_{i,j,k} = 1 && \forall k \\ & \color{darkred}x_{i-1,j,k+1}+\color{darkred}x_{i+1,j,k+1}+\color{darkred}x_{i,j-1,k+1}+\color{darkred}x_{i,j+1,k+1} \ge \color{darkred}x_{i,j,k} && \forall i,j, \forall k\lt 81 \\ & \color{darkred}x_{i,j,k}=1 && \text{for given cells}\\ & \color{darkred}x_{i,j,k} \in \{0,1\} \end{align}

In the neighbor constraint, I am careful about $$k$$ being possibly outside the board, but less so about $$i$$ and $$j$$. The reason is, that with GAMS, when addressing using leads/lags outside of the domain, a value of zero is returned. That is exactly what we need here, for $$i$$ and $$j$$. The $$k$$ index needs protection. The GAMS model can look like:

 sets   i 'rows'    /i1*i9/   j 'columns' /j1*j9/   k 'values'  /1*81/ ; table board(i,j)         j1 j2 j3 j4 j5 j6 j7 j8 j9   i1   i2       25 20 13 14 15  8 41   i3       26                40   i4        1                39   i5       30                38   i6       67                37   i7       70                50   i8       73 72 79 60 61 52 53   i9 ; binary variable x(i,j,k); variable dummy; equations    onevalue(i,j) 'each cell has one value (between 1 and 81)'    unique(k)     'each value must be in exactly one cell'    next(i,j,k)   'x(i,j,k)=1 => a neighbor is k+1'    obj           'dummy' ; obj.. dummy =e= 0; onevalue(i,j).. sum(k, x(i,j,k)) =e= 1; unique(k)..     sum((i,j),x(i,j,k)) =e= 1; next(i,j,k+1).. x(i-1,j,k+1)+x(i+1,j,k+1)+x(i,j-1,k+1)+x(i,j+1,k+1) =g= x(i,j,k); x.fx(i,j,k)$(board(i,j)=k.val) = 1; model m /all/; option mip=cplex; solve m minimizing dummy using mip; display x.l; parameter result(i,j); result(i,j) = sum(k, k.val*x.l(i,j,k)); option result:0; display result; This model is solved completely in preprocessing. Cplex needs 0 iterations and 0 nodes. The output looks like: ---- 46 PARAMETER result j1 j2 j3 j4 j5 j6 j7 j8 j9 i1 23 22 21 12 11 10 9 42 43 i2 24 25 20 13 14 15 8 41 44 i3 27 26 19 18 17 16 7 40 45 i4 28 1 2 3 4 5 6 39 46 i5 29 30 31 32 33 34 35 38 47 i6 68 67 66 65 64 63 36 37 48 i7 69 70 71 80 81 62 51 50 49 i8 74 73 72 79 60 61 52 53 54 i9 75 76 77 78 59 58 57 56 55  Notes: • We don't need an optcr (gap) setting: we are looking for a feasible solution. • GAMS has a rule: "return zero when a lead/lag is outside the domain". We need to know this when looking at equation next. • The k+1 in next(i,j,k+1) means that we skip the last k when generating this equation. This protects us against generating $$0 \ge \color{darkred}x_{i,j,81}$$. That would make the model infeasible. • We were lucky with addressing $$i\pm1, j\pm 1$$. GAMS returning zero there when outside the domain is exactly what we want. • We could have defined equation next in terms of $$k-1$$ and $$k$$. Then you would not need the k+1 in next(i,j,k+1). I.e., next(i,j,k).. x(i-1,j,k)+x(i+1,j,k)+x(i,j-1,k)+x(i,j+1,k) =g= x(i,j,k-1); This formulation does not need any protection, as the GAMS rule: "return zero when a lead/lag is outside the domain" is good in all cases. I still like the original formulation a little bit better as it is closer to what we stated in English. • Interestingly when solving as an RMIP (Relaxed MIP), Cplex will need to do a few iterations. It will still find the integer solution. I suspect this is just by accident. • A well-formed Numbrix puzzle has one unique solution. To verify this, we can add a cut to the model that forbids the current solution $$x^*$$: $\sum_{i,j,k} \color{darkblue}x^*_{i,j,k}\cdot \color{darkred} x_{i,j,k}\le 80$ and resolve. This must be infeasible. This is indeed the case for this input. In GAMS this test can be stated as:  equation cut 'cut off current solution x.l'; cut.. sum((i,j,k)$(x.l(i,j,k)>0.5), x(i,j,k)) =l= card(k)-1; model m2 /all/; solve m2 minimizing dummy using mip; abort\$(m2.modelstat <> %modelStat.infeasible% and        m2.modelstat <> %modelStat.integerInfeasible%) "Should be infeasible";
• Q: can we drop some given cells while making sure there is still one unique solution? A: yes
• Q: what is the minimum of given cells that we need to keep so the solution is unique? A: don't know. Looks like a difficult question to answer. Maybe write a model that tries to find two different solutions with given cells represented by binary variables.
• I think this model looks quite nice!

#### References

1. Numbrix puzzles by Marilyn vos Savant, https://parade.com/numbrix/
2. Sara Jensen, What is a Numbrix Puzzle? https://www.youtube.com/watch?v=FypFESO6rVo
3. Numbrix puzzle via TSP, http://yetanothermathprogrammingconsultant.blogspot.com/2020/12/numbrix-puzzle-via-tsp.html. Follow up post. Here we model the problem as a TSP problem. The MTZ formulation will give us the path as a byproduct. Although much smaller in size, the TSP problem is not as cheap to solve.