donman
donman

Reputation: 67

Constraints for Hitori puzzle

The Hitori puzzle consists of an nxn grid filled with numbers 1..n, some of them repeated in rows and/or columns. The goal is to cover some of the repeated numbers so that

i) there are no repeated numbers in any row or column,

ii) no covered square is vertically or horizontally adjacent to another covered square,

iii) and the uncovered squares are vertically or horizontally connected (in a single graph), ie there are no isolated islands of uncovered cells.

Below is a MiniZinc model and a small dataset.

When this model is run using the given dataset, it produces the following output.

Running hitori.mzn, hitori.dzn
79msec

hitorimask           hitorisoln 
  0  1  0  1  1     .  5  .  4  3
  1  1  1  1  0     4  2  5  1  .
  1  0  1  0  1     5  .  1  .  4
  0  1  0  1  1     .  4  .  5  2
  1  1  1  1  0     2  1  4  3  .
% time elapsed: 86msec
----------

hitorimask           hitorisoln 
  0  1  0  1  0     .  5  .  4  .
  1  1  1  1  1     4  2  5  1  3
  1  0  1  0  1     5  .  1  .  4
  0  1  0  1  1     .  4  .  5  2
  1  1  1  1  0     2  1  4  3  .
% time elapsed: 86msec
----------
==========
Finished in 79msec.

The first “solution” shows two regions of connected uncovered cells. The second solution is the correct one, where all the open cells are connected.

My last constraint in the MiniZinc model is a first stab at the connectedness condition. It simply avoids single cell islands. If one comments it out, 33 solutions are generated, many with single isolated cells.

Question: How can I model the connected condition in MiniZinc?

I have thought of modelling the connections of a candidate solution as an n^2 by n^2 Markov transition matrix and then riasing it to the n^2 power. I have not actually tried it, but it seems inefficient and awefully inelegant.

Suggestions will be appreciated.

Thank you.

% hitori.mzn  

% Given square a grid of sidelength n
% The grid is filled with numbers 1..n, some of them repeated in rows and columns
% Goal: Cover some of the repeated numbers so that
% i)    there are no duplicates in any row of column
% ii)   no covered square is vertically of horizontally adjacent to another coverde square
% iii)  the uncovered blocks are connected 

% data instance defined in *.dzn 
int: nsize ;
array[SDX,SDX] of SDX: hitoritabl ;

set of int: SDX = 1..nsize ;
set of int: SDX0 = 0..nsize ;  
set of int: SDX2 = 2..nsize ;

% primary unknowns, 0 or 1 whether whether a square is covered or open
array[SDX,SDX] of var 0..1 : hitorimask ;

% modelling variable to display the solution, 0 is used to indicate a covered cell, ie a blacked cell
array[SDX,SDX] of var SDX0: hitorisoln ;
constraint forall( ir in SDX, jc in SDX ) ( hitorisoln[ir,jc] = hitorimask[ir,jc] * hitoritabl[ir,jc] ) ;

% model constraints

% use primary variables to construct quantity on which nonrepeat constraint can be imposed
include "alldifferent_except_0.mzn" ; 
constraint forall ( ir in SDX) ( alldifferent_except_0 ( [ hitorisoln[ir,jc]| jc in SDX ] ) ) ;  % no repeats in a row
constraint forall ( jc in SDX) ( alldifferent_except_0 ( [ hitorisoln[ir,jc]| ir in SDX ] ) ) ;  % no repeats in a column

% impose condition on non-adjacent masked cells
% no black square directly below anothe black square
constraint forall (ir in SDX2, jc in SDX) ( 0 < hitorimask[ir-1,jc] + hitorimask[ir,jc] ) ;
% no black square directly to the right of anothe black square
constraint forall (ir in SDX, jc in SDX2) ( 0 < hitorimask[ir,jc-1] + hitorimask[ir,jc] ) ;

% this constraint simply excludes islands of singleton open cells
constraint forall (ir in SDX, jc in SDX) 
           ( if hitorimask[ir,jc]=1 then 
                ( if ir > 1     then hitorimask[ir-1,jc  ] else 0 endif 
                + if jc < nsize then hitorimask[ir  ,jc+1] else 0 endif
                + if ir < nsize then hitorimask[ir+1,jc  ] else 0 endif
                + if jc > 1     then hitorimask[ir  ,jc-1] else 0 endif
                 > 0 
                )
             endif ) ;

%
% still need a condition of connectedness of open squares
%

solve satisfy ;

output ["\nhitorimask           hitorisoln \n" ] ;
output [ 
         if jc <= nsize then "" ++ show_int(3, hitorimask[ir,jc] ) else "" endif
      ++ if jc = nsize then "   " else "" endif 
      ++ if nsize < jc /\ jc <= 2*nsize 
         then if fix(hitorisoln[ir,jc-nsize])>0 
              then show_int(3, hitorisoln[ir,jc-nsize] )  
              else "  ." 
              endif  
         else "" 
         endif 
      ++ if jc = 2*nsize then "\n" else "" endif
      | ir in SDX, jc in 1..2*nsize ] ;
 

data file

% hitori.dzn
nsize = 5 ;   % size of square grid

hitoritabl = [|
4, 5, 4, 4, 3 |
4, 2, 5, 1, 3 |
5, 1, 1, 3, 4 |
5, 4, 2, 5, 2 |
2, 1, 4, 3, 1 |] ;


Upvotes: 1

Views: 259

Answers (4)

Anonymous
Anonymous

Reputation: 111

This second approach to ensure connectivity among all unshaded cells is based on flow conservation:

2. flow conservation method

  • based on 2 directed edges (incoming and outgoing) between unshaded cells as nodes
  • "Starting a flow from a root cell, unshaded cells/nodes keep 1 unit while passing on the flow to their connected, adjacent neighbors.
  • i.e. we can avoid recursive propagation by formulating the sum of incoming flows and the sum of outgoing flows for each unshaded cell
  • the main constraints are:
    • sum of incoming flows == sum of outgoing flows - 1
    • outgoing flow of root == numbers of unshaded cells - 1

Edit:
adding PoC (Python CP-SAT):
performance is slower than method (3)


### FLOW CONSERVATION

from ortools.sat.python import cp_model
import numpy as np


def check_connectivity_flow(matrix):
   model = cp_model.CpModel()
   rows, cols = matrix.shape

   # Identify all True cells and count them.
   true_cells = [(r, c) for r in range(rows) for c in range(cols) if matrix[r, c]]
   num_true = len(true_cells)

   # Pick first True cell as  a root
   root = true_cells[0]

   # Allowed edges: adjacent in same row (col diff 1) or same column (row diff 1)
   allowed_edges = []
   for (r, c) in true_cells:
       for (dr, dc) in [(0, 1), (1, 0), (0, -1), (-1, 0)]:
           nr, nc = r + dr, c + dc
           if (nr, nc) in true_cells:
               allowed_edges.append(((r, c), (nr, nc)))

   # Flow variables for each allowed directed edge
   flow_vars = {}
   for (u, v) in allowed_edges:
       flow_vars[(u, v)] = model.new_int_var(0, num_true - 1, f"flow_{u}_{v}")
       flow_vars[(v, u)] = model.new_int_var(0, num_true - 1, f"flow_{v}_{u}")
       flow_direct = model. new_bool_var(f"flow_direct_{u}_{v}")

       # Prevent bidirectional flow
       model.add(flow_vars[(u, v)] > 0).only_enforce_if(flow_direct)
       model.add(flow_vars[(v, u)] == 0).only_enforce_if(flow_direct)
       model.add(flow_vars[(v, u)] > 0).only_enforce_if(flow_direct.Not())
       model.add(flow_vars[(u, v)] == 0).only_enforce_if(flow_direct.Not())

   # Flow conservation 
   for cell in true_cells:
       out_flow = sum(flow_vars[(cell, v)] for (u, v) in flow_vars if u == cell)
       in_flow = sum(flow_vars[(u, cell)] for (u, v) in flow_vars if v == cell)

       if cell == root:
           # Root must send out exactly num_true - 1 units of flow
           model.add(out_flow == num_true - 1)
           # Optionally, force in_flow == 0
           model.add(in_flow == 0)
       else:
           # Every other cell must have net inflow 1
           model.add(in_flow == out_flow + 1)

   # Solve the model
   solver = cp_model.CpSolver()
   status = solver.Solve(model)

   print(solver.ResponseStats())
   print(f"Input mask:\n{matrix}")

   # If the model is infeasible, connectivity cannot be established.
   if status != cp_model.OPTIMAL and status != cp_model.FEASIBLE:
       return False

   # Flow difference matrix for each cell
   flow_diff_matrix = np.zeros((rows, cols), dtype=int)
   in_flow_matrix = np.zeros((rows, cols), dtype=int)
   out_flow_matrix = np.zeros((rows, cols), dtype=int)

   for cell in true_cells:
       in_flow = sum(solver.value(flow_vars[(u, cell)]) for (u, v) in flow_vars if v == cell)
       out_flow = sum(solver.value(flow_vars[(cell, v)]) for (u, v) in flow_vars if u == cell)
       flow_diff_matrix[cell[0], cell[1]] = in_flow - out_flow

       in_flow_matrix[cell[0], cell[1]] = in_flow
       out_flow_matrix[cell[0], cell[1]] = out_flow

   print(f"\nflow conservation:\n{flow_diff_matrix}")
   print(f"\nin_flow:\n{in_flow_matrix}")
   print(f"\nout_flow:\n{out_flow_matrix}")

   return True


# Examples
small_mask = np.array([
   [True,  False, False],
   [True,  False, True ],
   [True,  False, True ]
])

true_mask = np.array([
   [0, 1, 0, 1, 0],
   [1, 1, 1, 1, 1],
   [1, 0, 1, 0, 1],
   [0, 1, 0, 1, 1],
   [1, 1, 1, 1, 0]
])

false_mask = np.array([
   [0, 1, 0, 1, 1],
   [1, 1, 1, 1, 0],
   [1, 0, 1, 0, 1],
   [0, 1, 0, 1, 1],
   [1, 1, 1, 1, 0]
])

bigger_mask = np.array([
   [1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 1],
   [0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0],
   [1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1],
   [1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0],
   [0, 1, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1],
   [1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0],
   [0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0],
   [1, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0],
   [1, 0, 1, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 1],
   [1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0],
   [1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1],
   [0, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1],
   [1, 1, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1],
   [1, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1],
   [0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1],
   [1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1],
   [1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1]
])


print(f"\n> Is Connected: {check_connectivity_flow(false_mask)}\n***************\n")

print(f"\n> Is Connected: {check_connectivity_flow(bigger_mask)}\n***************\n")

Output:

CpSolverResponse summary:
status: INFEASIBLE
objective: NA
best_bound: NA
integers: 0
booleans: 0
conflicts: 0
branches: 0
propagations: 0
integer_propagations: 0
restarts: 0
lp_iterations: 0
walltime: 0.000847509
usertime: 0.000847564
deterministic_time: 0
gap_integral: 0

Input mask:
[[0 1 0 1 1]
 [1 1 1 1 0]
 [1 0 1 0 1]
 [0 1 0 1 1]
 [1 1 1 1 0]]

> Is Connected: False
***************

CpSolverResponse summary:
status: OPTIMAL
objective: 0
best_bound: 0
integers: 263
booleans: 131
conflicts: 0
branches: 262
propagations: 512
integer_propagations: 2161
restarts: 262
lp_iterations: 126
walltime: 0.0671798
usertime: 0.0671798
deterministic_time: 0.00821366
gap_integral: 0
solution_fingerprint: 0x611a630ef34b9365

Input mask:
[[1 1 1 1 0 1 1 0 1 1 1 0 1 1 0 1 1]
 [0 1 1 0 1 0 1 1 0 1 1 1 1 0 1 1 0]
 [1 0 1 1 1 1 1 0 1 1 0 1 0 1 0 1 1]
 [1 1 1 0 1 0 1 1 0 1 1 1 1 1 1 1 0]
 [0 1 0 1 0 1 0 1 1 1 0 1 1 1 0 1 1]
 [1 1 1 1 1 1 1 1 0 1 0 1 1 0 1 1 0]
 [0 1 1 0 1 1 0 1 1 0 1 0 1 1 0 1 0]
 [1 1 0 1 0 1 1 1 0 1 1 1 1 0 1 1 0]
 [1 0 1 1 1 0 1 0 1 0 1 0 1 1 0 1 1]
 [1 1 0 1 0 1 1 1 1 1 0 1 0 1 1 1 0]
 [1 0 1 1 1 1 0 1 0 1 1 1 1 0 1 0 1]
 [0 1 0 1 0 1 1 1 1 0 1 0 1 1 0 1 1]
 [1 1 1 0 1 0 1 0 1 1 0 1 1 1 1 0 1]
 [1 1 0 1 1 1 0 1 0 1 1 1 0 1 0 1 1]
 [0 1 1 1 0 1 1 1 1 1 0 1 1 1 1 0 1]
 [1 1 0 1 1 0 1 1 1 0 1 1 0 1 1 1 1]
 [1 0 1 1 0 1 1 0 1 1 0 1 1 0 1 0 1]]

flow conservation:
[[-201   1   1   1   0   1   1   0   1   1   1   0   1   1   0   1   1]
 [   0   1   1   0   1   0   1   1   0   1   1   1   1   0   1   1   0]
 [   1   0   1   1   1   1   1   0   1   1   0   1   0   1   0   1   1]
 [   1   1   1   0   1   0   1   1   0   1   1   1   1   1   1   1   0]
 [   0   1   0   1   0   1   0   1   1   1   0   1   1   1   0   1   1]
 [   1   1   1   1   1   1   1   1   0   1   0   1   1   0   1   1   0]
 [   0   1   1   0   1   1   0   1   1   0   1   0   1   1   0   1   0]
 [   1   1   0   1   0   1   1   1   0   1   1   1   1   0   1   1   0]
 [   1   0   1   1   1   0   1   0   1   0   1   0   1   1   0   1   1]
 [   1   1   0   1   0   1   1   1   1   1   0   1   0   1   1   1   0]
 [   1   0   1   1   1   1   0   1   0   1   1   1   1   0   1   0   1]
 [   0   1   0   1   0   1   1   1   1   0   1   0   1   1   0   1   1]
 [   1   1   1   0   1   0   1   0   1   1   0   1   1   1   1   0   1]
 [   1   1   0   1   1   1   0   1   0   1   1   1   0   1   0   1   1]
 [   0   1   1   1   0   1   1   1   1   1   0   1   1   1   1   0   1]
 [   1   1   0   1   1   0   1   1   1   0   1   1   0   1   1   1   1]
 [   1   0   1   1   0   1   1   0   1   1   0   1   1   0   1   0   1]]

in_flow:
[[  0 201 131   1   0   1   2   0   1   6   4   0   2   1   0   2   1]
 [  0  69 197   0   1   0   4   1   0   7   3   4   3   0   1   4   0]
 [  1   0 196  91  90  87  86   0   1   8   0   4   0   1   0   6   1]
 [  2 103 104   0   1   0  81  80   0  54  45  44  16  18  16  15   0]
 [  0 100   0   1   0   1   0  79  57  56   0  23  17  10   0   8   1]
 [  1  99  89  88  86  24   9  21   0   1   0  11  16   0   1   6   0]
 [  0  16   9   0  69  90   0  11   1   0   1   0  15   1   0   4   0]
 [  5   6   0   1   0  89  96   9   0   1   4   5  13   0   1   3   0]
 [  4   0   1   4   1   0  95   0   1   0   1   0   7   6   0   2   1]
 [  3   1   0   5   0  45  94  48  18  16   0   1   0   5   4   2   0]
 [  1   0   1   8   9  44   0  29   0  15  14  12  10   0   1   0   1]
 [  0   1   0   1   0  34  33  59  58   0   1   0  10   9   0   1   3]
 [  3   6   1   0   1   0   1   0  57  56   0   2   2   8   1   0   4]
 [  2   7   0  17  19  20   0   1   0  55  21  20   0   5   0   1   6]
 [  0  10  11  16   0  21  26  28  32  33   0  17  12  22   8   0   7]
 [  2   3   0   4   1   0   4   2   4   0   1   4   0  21  20  10   9]
 [  1   0   1   2   0   1   2   0   2   1   0   2   1   0   1   0   1]]

out_flow:
[[201 200 130   0   0   0   1   0   0   5   3   0   1   0   0   1   0]
 [  0  68 196   0   0   0   3   0   0   6   2   3   2   0   0   3   0]
 [  0   0 195  90  89  86  85   0   0   7   0   3   0   0   0   5   0]
 [  1 102 103   0   0   0  80  79   0  53  44  43  15  17  15  14   0]
 [  0  99   0   0   0   0   0  78  56  55   0  22  16   9   0   7   0]
 [  0  98  88  87  85  23   8  20   0   0   0  10  15   0   0   5   0]
 [  0  15   8   0  68  89   0  10   0   0   0   0  14   0   0   3   0]
 [  4   5   0   0   0  88  95   8   0   0   3   4  12   0   0   2   0]
 [  3   0   0   3   0   0  94   0   0   0   0   0   6   5   0   1   0]
 [  2   0   0   4   0  44  93  47  17  15   0   0   0   4   3   1   0]
 [  0   0   0   7   8  43   0  28   0  14  13  11   9   0   0   0   0]
 [  0   0   0   0   0  33  32  58  57   0   0   0   9   8   0   0   2]
 [  2   5   0   0   0   0   0   0  56  55   0   1   1   7   0   0   3]
 [  1   6   0  16  18  19   0   0   0  54  20  19   0   4   0   0   5]
 [  0   9  10  15   0  20  25  27  31  32   0  16  11  21   7   0   6]
 [  1   2   0   3   0   0   3   1   3   0   0   3   0  20  19   9   8]
 [  0   0   0   1   0   0   1   0   1   0   0   1   0   0   0   0   0]]

> Is Connected: True
***************

Addendum
This is the result for a 15x15 Hitori on Colab:

CpSolverResponse summary:
status: OPTIMAL
objective: 0
best_bound: 0
integers: 356
booleans: 46
conflicts: 0
branches: 92
propagations: 1058
integer_propagations: 2398
restarts: 92
lp_iterations: 350
walltime: 0.10211
usertime: 0.102111
deterministic_time: 0.00997731
gap_integral: 0
solution_fingerprint: 0xa5547dd0b45d7440


Input:
[ 15   3   1   4  15  13   1   6  15   5   7   2   1   7  11 ]
[  8  13   6  13   4  14  13   7  11  12  10   1  15   5   1 ]
[  6   1  15   8   6  15   4   6  13  14  15   3   5   6   9 ]
[ 11  12   3  12  12   1  12  15  14  12   5  12   8   2   4 ]
[  6   8   3   1  13   6   5   9   3  15   3   4  10  11   9 ]
[ 10   6  11  10   7   9   5   5   3   6   2  15  10   4   8 ]
[ 12  11   7   2  12   6  15   4   7   9  13   7   3  12  10 ]
[  3  11   4  11   5  11   9  11   6   2  11  14   2  10  11 ]
[  6   5   6   9   6  12   6  14   6  10  15   6  11   6   3 ]
[  4  14   2  13  15   5  11   3   7   8   1  10   5   9   3 ]
[ 11   2  14   3   8  10  11   1  11   7   3  11  13   3   6 ]
[  2   4  12  15   1  12   6  12   8  12  11   5  12  14  12 ]
[  5  13  12  14  13   8   7  11  13   4   9  13   2  13  15 ]
[ 14  15   5  12   9   4  12  10   2  11   1  13   1  12   7 ]
[ 13   7  13   6   4   3   2   4  10  13  14   3   4  15   4 ]

Mask (unshaded == 1):
[  1   1   1   1   0   1   0   1   0   1   0   1   0   1   1 ]
[  1   0   1   0   1   1   1   1   1   1   1   1   1   1   0 ]
[  0   1   1   1   1   0   1   0   1   1   0   1   1   0   1 ]
[  1   0   1   0   1   1   0   1   1   0   1   0   1   1   1 ]
[  1   1   0   1   1   0   1   1   0   1   1   1   1   1   0 ]
[  0   1   1   1   1   1   0   1   1   0   1   1   0   1   1 ]
[  1   1   0   1   0   1   1   1   0   1   1   1   1   0   1 ]
[  1   0   1   1   1   0   1   0   1   1   0   1   0   1   0 ]
[  0   1   0   1   0   1   0   1   0   1   1   0   1   1   1 ]
[  1   1   1   1   1   1   1   1   1   1   1   1   0   1   0 ]
[  0   1   1   0   1   1   0   1   0   1   0   1   1   1   1 ]
[  1   1   0   1   1   0   1   1   1   0   1   1   0   1   0 ]
[  1   0   1   1   0   1   1   1   0   1   1   0   1   1   1 ]
[  1   1   1   1   1   1   0   1   1   1   0   1   1   0   1 ]
[  0   1   0   1   0   1   1   0   1   1   1   0   1   1   0 ]

Solution:
[ 15   3   1   4   _  13   _   6   _   5   _   2   _   7  11 ]
[  8   _   6   _   4  14  13   7  11  12  10   1  15   5   _ ]
[  _   1  15   8   6   _   4   _  13  14   _   3   5   _   9 ]
[ 11   _   3   _  12   1   _  15  14   _   5   _   8   2   4 ]
[  6   8   _   1  13   _   5   9   _  15   3   4  10  11   _ ]
[  _   6  11  10   7   9   _   5   3   _   2  15   _   4   8 ]
[ 12  11   _   2   _   6  15   4   _   9  13   7   3   _  10 ]
[  3   _   4  11   5   _   9   _   6   2   _  14   _  10   _ ]
[  _   5   _   9   _  12   _  14   _  10  15   _  11   6   3 ]
[  4  14   2  13  15   5  11   3   7   8   1  10   _   9   _ ]
[  _   2  14   _   8  10   _   1   _   7   _  11  13   3   6 ]
[  2   4   _  15   1   _   6  12   8   _  11   5   _  14   _ ]
[  5   _  12  14   _   8   7  11   _   4   9   _   2  13  15 ]
[ 14  15   5  12   9   4   _  10   2  11   _  13   1   _   7 ]
[  _   7   _   6   _   3   2   _  10  13  14   _   4  15   _ ]

If we don’t need to see the inflows and outflows of the unshaded cells separately, we can simplify the model.
By directly aggregating all flows of a cell at once, we can reduce the numbers of decision variables and constraints.

This simplified model matches more with method (3), both in terms of complexity and performance.

Stats for the same 15x15 Hitori (simplified model):

CpSolverResponse summary:
status: OPTIMAL
objective: 0
best_bound: 0
integers: 36
booleans: 0
conflicts: 0
branches: 0
propagations: 0
integer_propagations: 0
restarts: 0
lp_iterations: 39
walltime: 0.0686117
usertime: 0.0686119
deterministic_time: 0.001577
gap_integral: 0
solution_fingerprint: 0x3c10622297f28acb

Stats for the 20x20 Hitori (simplified model):

CpSolverResponse summary:
status: OPTIMAL
objective: 0
best_bound: 0
integers: 33
booleans: 0
conflicts: 0
branches: 0
propagations: 0
integer_propagations: 0
restarts: 0
lp_iterations: 28
walltime: 0.11172
usertime: 0.111721
deterministic_time: 0.00226062
gap_integral: 0
solution_fingerprint: 0x2274178dfe825f1e

Upvotes: 1

Anonymous
Anonymous

Reputation: 111

I figured out 2 approaches to formulate the connectedness constraint for all unshaded (not covered) cells.

The first approach is based on transitive relation among the unshaded cells:

1. transitivity method

  • based on undirected edges between unshaded cells as nodes
  • "If A is related to B, and B is related to C, then A is also related to C"
  • basically, this approach matches in complexity with OP’s description of a “...candidate solution as an n^2 by n^2 Markov transition matrix and then raising it to the n^2 power”

Edit:
The (1) transitivity method struggles to scale effectively … therefore, replacing it with
the (3) min. distance method, which outperforms the (2) flow conservation method as well:

3. minimum distance method

  • based on 1 directed edge between cells as nodes
  • "Distance of an unshaded cell from a root cell is defined as the min. distance of its neighbors + 1."
  • the main constraints are:
    • distance of shaded cells == max_distance == n * n
    • distance of unshaded cells to root == 1 + min distance of its adjacent neighbors < max_distance
    • root distance == 0

PoC (Python CP-SAT) to check a given matrix whether the True/1-cells form a connected area:

note: it was easier for me to setup colab with ortools than with minizinc and additional solvers, and I am also more familiar with CP-SAT

### MIN DISTANCE

from ortools.sat.python import cp_model
import numpy as np
import math


def check_connectivity_min_dist(matrix):
   model = cp_model.CpModel()
   rows, cols = matrix.shape

   # Max distance is n * n
   max_distance = rows * cols

   # Identify all True cells and count them.
   true_cells = [(r, c) for r in range(rows) for c in range(cols) if matrix[r, c]]

   # Create a distance variable for each cell (including non-colored cells)
   distance_vars = {}
   for (r, c) in true_cells:
       distance_vars[(r, c)] = model.new_int_var(0, max_distance, f"distance_{r}_{c}")

   # Set distance of root cell (first true cell) to 0
   root = true_cells[0]

   # Function to get neighbors for any cell
   def get_neighbors(r, c):
       neighbors = []
       for dr, dc in [(0, 1), (1, 0), (0, -1), (-1, 0)]:
           nr, nc = r + dr, c + dc
           if 0 <= nr < rows and 0 <= nc < cols:
               neighbors.append((nr, nc))
       return neighbors

   # Distance Constraints:
   for (r, c) in true_cells:
       if (r, c) == root:
           model.Add(distance_vars[(r, c)] == 0)
       else:
           min_neighbor_distance = model.new_int_var(0, max_distance, f"min_neighbor_distance_{r}_{c}")
           neighbors = [distance_vars[nr, nc] for nr, nc in get_neighbors(r, c) if matrix[nr, nc]]
           model.add_min_equality(target=min_neighbor_distance, exprs=neighbors)
           model.add(distance_vars[(r, c)] ==  1 + min_neighbor_distance)
           model.add(distance_vars[(r, c)] < max_distance)

   # Solve the model
   solver = cp_model.CpSolver()
   status = solver.Solve(model)

   print(solver.ResponseStats())
   print(f"Input mask:\n{matrix}")

   # If the model is infeasible, connectivity cannot be established.
   if status != cp_model.OPTIMAL and status != cp_model.FEASIBLE:
       return False

   # Distance matrix
   distance_matrix = np.full((rows, cols), -1, dtype=int)
   for (r, c) in true_cells:
           distance_matrix[r, c] = solver.Value(distance_vars[(r, c)])
   print(f"\nDistances:\n{distance_matrix}")

   return True


# Examples
small_mask = np.array([
   [True,  False, False],
   [True,  False, True ],
   [True,  False, True ]
])

true_mask = np.array([
   [0, 1, 0, 1, 0],
   [1, 1, 1, 1, 1],
   [1, 0, 1, 0, 1],
   [0, 1, 0, 1, 1],
   [1, 1, 1, 1, 0]
])

false_mask = np.array([
   [0, 1, 0, 1, 1],
   [1, 1, 1, 1, 0],
   [1, 0, 1, 0, 1],
   [0, 1, 0, 1, 1],
   [1, 1, 1, 1, 0]
])

bigger_mask = np.array([
   [1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 1],
   [0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0],
   [1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1],
   [1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0],
   [0, 1, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1],
   [1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0],
   [0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0],
   [1, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0],
   [1, 0, 1, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 1],
   [1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0],
   [1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1],
   [0, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1],
   [1, 1, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1],
   [1, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1],
   [0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1],
   [1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1],
   [1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1]
])

print(f"\n> Is Connected: {check_connectivity_min_dist(false_mask)}\n***************\n")

print(f"\n> Is Connected: {check_connectivity_min_dist(bigger_mask)}\n***************\n")

Output:

CpSolverResponse summary:
status: INFEASIBLE
objective: NA
best_bound: NA
integers: 0
booleans: 0
conflicts: 0
branches: 0
propagations: 0
integer_propagations: 0
restarts: 0
lp_iterations: 0
walltime: 0.000979077
usertime: 0.000979146
deterministic_time: 0
gap_integral: 0

Input mask:
[[0 1 0 1 1]
 [1 1 1 1 0]
 [1 0 1 0 1]
 [0 1 0 1 1]
 [1 1 1 1 0]]

> Is Connected: False
***************

CpSolverResponse summary:
status: OPTIMAL
objective: 0
best_bound: 0
integers: 0
booleans: 0
conflicts: 0
branches: 0
propagations: 0
integer_propagations: 0
restarts: 0
lp_iterations: 0
walltime: 0.0105172
usertime: 0.0105172
deterministic_time: 6.5e-07
gap_integral: 0
solution_fingerprint: 0xd22cc995cb38921f

Input mask:
[[1 1 1 1 0 1 1 0 1 1 1 0 1 1 0 1 1]
 [0 1 1 0 1 0 1 1 0 1 1 1 1 0 1 1 0]
 [1 0 1 1 1 1 1 0 1 1 0 1 0 1 0 1 1]
 [1 1 1 0 1 0 1 1 0 1 1 1 1 1 1 1 0]
 [0 1 0 1 0 1 0 1 1 1 0 1 1 1 0 1 1]
 [1 1 1 1 1 1 1 1 0 1 0 1 1 0 1 1 0]
 [0 1 1 0 1 1 0 1 1 0 1 0 1 1 0 1 0]
 [1 1 0 1 0 1 1 1 0 1 1 1 1 0 1 1 0]
 [1 0 1 1 1 0 1 0 1 0 1 0 1 1 0 1 1]
 [1 1 0 1 0 1 1 1 1 1 0 1 0 1 1 1 0]
 [1 0 1 1 1 1 0 1 0 1 1 1 1 0 1 0 1]
 [0 1 0 1 0 1 1 1 1 0 1 0 1 1 0 1 1]
 [1 1 1 0 1 0 1 0 1 1 0 1 1 1 1 0 1]
 [1 1 0 1 1 1 0 1 0 1 1 1 0 1 0 1 1]
 [0 1 1 1 0 1 1 1 1 1 0 1 1 1 1 0 1]
 [1 1 0 1 1 0 1 1 1 0 1 1 0 1 1 1 1]
 [1 0 1 1 0 1 1 0 1 1 0 1 1 0 1 0 1]]

Distances:
[[ 0  1  2  3 -1 11 10 -1 18 17 18 -1 20 21 -1 23 24]
 [-1  2  3 -1  7 -1  9 10 -1 16 17 18 19 -1 23 22 -1]
 [ 8 -1  4  5  6  7  8 -1 16 15 -1 17 -1 19 -1 21 22]
 [ 7  6  5 -1  7 -1  9 10 -1 14 15 16 17 18 19 20 -1]
 [-1  7 -1 11 -1 13 -1 11 12 13 -1 17 18 19 -1 21 22]
 [ 9  8  9 10 11 12 13 12 -1 14 -1 18 19 -1 23 22 -1]
 [-1  9 10 -1 12 13 -1 13 14 -1 24 -1 20 21 -1 23 -1]
 [11 10 -1 24 -1 14 15 14 -1 24 23 22 21 -1 25 24 -1]
 [12 -1 24 23 24 -1 16 -1 20 -1 24 -1 22 23 -1 25 26]
 [13 14 -1 22 -1 18 17 18 19 20 -1 24 -1 24 25 26 -1]
 [14 -1 22 21 20 19 -1 19 -1 21 22 23 24 -1 26 -1 38]
 [-1 38 -1 22 -1 20 21 20 21 -1 23 -1 25 26 -1 38 37]
 [38 37 38 -1 32 -1 22 -1 22 23 -1 27 26 27 28 -1 36]
 [37 36 -1 32 31 30 -1 28 -1 24 25 26 -1 28 -1 36 35]
 [-1 35 34 33 -1 29 28 27 26 25 -1 27 28 29 30 -1 34]
 [37 36 -1 34 35 -1 29 28 27 -1 29 28 -1 30 31 32 33]
 [38 -1 36 35 -1 31 30 -1 28 29 -1 29 30 -1 32 -1 34]]

> Is Connected: True
***************

Addendum
This is the result for a 20x20 Hitori on Colab:

CpSolverResponse summary:
status: OPTIMAL
objective: 0
best_bound: 0
integers: 207
booleans: 1
conflicts: 0
branches: 1
propagations: 0
integer_propagations: 420
restarts: 0
lp_iterations: 0
walltime: 0.0888325
usertime: 0.0888326
deterministic_time: 0.00802969
gap_integral: 0
solution_fingerprint: 0x88d601be3265190d


Input:
[  7  10  16  11   2   9   2   3  14  14  17  19   4   4  13  17  13  15   5  11 ]
[  1  18  13  19   9   8   2   1   7   5  15   7   3   4   9   3  11  11  16  16 ]
[  9  19   6  11   7  16  14   2   3  20  11   5  12  11  15   4  18   3  14  10 ]
[  5   6   9  11  12   8  17  17   4   8  18   8  15  14  12  19  13  18   2   1 ]
[ 20   4  11   9   7  17  13   5   1  15  19  15   7  10  12  11   2   8  18  11 ]
[ 18  12   4  17  15  20  11   7  14   1  11   9  18  13   2   6   3  18  19  17 ]
[ 11  13   3  16  18   2  20  15   7   3   2  13   8  19   2   1  17  11   2   9 ]
[  4   8  10   3  15  14   1   4   9  13  12   2  13  18  19   2   6  18  20   6 ]
[ 10  16   5  20   1  15   9   9  18  19   6  10   6  11   7  13  12   3  17   4 ]
[  5  14   8   1  13   4  12   7  11  16   9  11  17  19   3   5  16  19   1  12 ]
[  4   7   7  15  11   6  19  12  11   2  19   4  16  10  18  20   9  11   8   5 ]
[  3  17   5   6   5  11  15  19   8  18   1  12   3  15  10  10   3  13   9  13 ]
[ 10   5   9  18  20   4  10  18  16  10   7  12   6  12   1   2   8   7   3  19 ]
[ 14   3  18   2  16   5  12  13   1   9  16  11  20  11  10   8  17   7   4  17 ]
[  1  20  19   4  17   8  16  19  19   6  17   8  16   7   4  11   1  12  13  18 ]
[ 11  12  19   2   8  18   9  14  19  17  11  19   1   5  20   5  10   8   6  14 ]
[ 19  19  17  13   4  12   3  12  15  10  14   1  18  15   8   2   4   5  16   7 ]
[ 16  11   7  18   4   3  14   7   7  15   1   3  10   9   6  17  19   2  12   2 ]
[  7  15  12  14   3  19   4  10  17  17  20  16   1  18  14   9   1   2   7   8 ]
[ 13   3  16   8  10   2  10   4  10  16   3  18  19   8  17  12  12   6   2  20 ]

Mask (unshaded == 1):
[  1   1   1   0   1   1   0   1   0   1   0   1   1   0   1   1   0   1   1   1 ]
[  0   1   1   1   0   1   1   1   1   1   1   0   1   1   1   0   1   0   1   0 ]
[  1   1   1   0   1   1   0   1   0   1   0   1   1   0   1   1   1   1   1   1 ]
[  0   1   0   1   1   0   1   0   1   1   1   0   1   1   0   1   1   0   1   1 ]
[  1   1   1   1   0   1   1   1   1   0   1   1   1   0   1   0   1   1   1   0 ]
[  1   0   1   0   1   1   1   0   1   1   0   1   0   1   1   1   1   0   1   1 ]
[  0   1   0   1   1   0   1   1   0   1   1   0   1   1   0   1   0   1   0   1 ]
[  1   1   1   1   0   1   1   0   1   0   1   1   1   0   1   0   1   1   1   0 ]
[  0   1   0   1   1   1   0   1   1   1   1   1   0   1   1   1   1   0   1   1 ]
[  1   1   1   0   1   0   1   0   1   0   1   0   1   0   1   0   1   1   1   0 ]
[  0   1   0   1   1   1   1   1   0   1   0   1   1   1   1   1   1   0   1   1 ]
[  1   1   1   1   0   1   0   1   1   1   1   1   0   1   0   1   0   1   1   0 ]
[  0   1   1   0   1   1   1   1   1   0   1   0   1   1   1   0   1   0   1   1 ]
[  1   0   1   1   0   1   0   1   0   1   1   1   1   0   1   1   1   1   1   0 ]
[  1   1   0   1   1   0   1   0   1   1   0   1   0   1   0   1   0   1   1   1 ]
[  0   1   1   0   1   1   1   1   0   1   1   0   1   1   1   0   1   0   1   0 ]
[  1   0   1   1   0   1   1   0   1   1   1   1   1   0   1   1   1   1   0   1 ]
[  1   1   0   1   1   0   1   1   0   1   0   1   1   1   1   0   1   0   1   1 ]
[  0   1   1   1   1   1   1   1   1   0   1   1   0   1   0   1   1   1   1   1 ]
[  1   1   0   1   0   1   0   1   1   1   0   1   1   0   1   1   0   1   0   1 ]

Solution:
[  7  10  16   _   2   9   _   3   _  14   _  19   4   _  13  17   _  15   5  11 ]
[  _  18  13  19   _   8   2   1   7   5  15   _   3   4   9   _  11   _  16   _ ]
[  9  19   6   _   7  16   _   2   _  20   _   5  12   _  15   4  18   3  14  10 ]
[  _   6   _  11  12   _  17   _   4   8  18   _  15  14   _  19  13   _   2   1 ]
[ 20   4  11   9   _  17  13   5   1   _  19  15   7   _  12   _   2   8  18   _ ]
[ 18   _   4   _  15  20  11   _  14   1   _   9   _  13   2   6   3   _  19  17 ]
[  _  13   _  16  18   _  20  15   _   3   2   _   8  19   _   1   _  11   _   9 ]
[  4   8  10   3   _  14   1   _   9   _  12   2  13   _  19   _   6  18  20   _ ]
[  _  16   _  20   1  15   _   9  18  19   6  10   _  11   7  13  12   _  17   4 ]
[  5  14   8   _  13   _  12   _  11   _   9   _  17   _   3   _  16  19   1   _ ]
[  _   7   _  15  11   6  19  12   _   2   _   4  16  10  18  20   9   _   8   5 ]
[  3  17   5   6   _  11   _  19   8  18   1  12   _  15   _  10   _  13   9   _ ]
[  _   5   9   _  20   4  10  18  16   _   7   _   6  12   1   _   8   _   3  19 ]
[ 14   _  18   2   _   5   _  13   _   9  16  11  20   _  10   8  17   7   4   _ ]
[  1  20   _   4  17   _  16   _  19   6   _   8   _   7   _  11   _  12  13  18 ]
[  _  12  19   _   8  18   9  14   _  17  11   _   1   5  20   _  10   _   6   _ ]
[ 19   _  17  13   _  12   3   _  15  10  14   1  18   _   8   2   4   5   _   7 ]
[ 16  11   _  18   4   _  14   7   _  15   _   3  10   9   6   _  19   _  12   2 ]
[  _  15  12  14   3  19   4  10  17   _  20  16   _  18   _   9   1   2   7   8 ]
[ 13   3   _   8   _   2   _   4  10  16   _  18  19   _  17  12   _   6   _  20 ]

Distances:
[  0   1   2   _  14  13   _  15   _  17   _  27  26   _  28  29   _  35  34  35 ]
[  _   2   3   4   _  12  13  14  15  16  17   _  25  26  27   _  31   _  33   _ ]
[  4   3   4   _  10  11   _  15   _  17   _  25  24   _  28  29  30  31  32  33 ]
[  _   4   _   8   9   _  23   _  19  18  19   _  23  24   _  30  31   _  33  34 ]
[  6   5   6   7   _  23  22  21  20   _  20  21  22   _  32   _  32  33  34   _ ]
[  7   _   7   _  25  24  23   _  21  22   _  22   _  30  31  32  33   _  35  36 ]
[  _  31   _  27  26   _  24  25   _  23  24   _  28  29   _  33   _  49   _  37 ]
[ 31  30  29  28   _  26  25   _  29   _  25  26  27   _  45   _  47  48  49   _ ]
[  _  31   _  29  28  27   _  29  28  27  26  27   _  45  44  45  46   _  48  49 ]
[ 33  32  33   _  29   _  33   _  29   _  27   _  41   _  43   _  45  46  47   _ ]
[  _  33   _  31  30  31  32  33   _  37   _  39  40  41  42  43  44   _  48  49 ]
[ 35  34  33  32   _  32   _  34  35  36  37  38   _  42   _  44   _  50  49   _ ]
[  _  35  34   _  34  33  34  35  36   _  38   _  42  43  44   _  48   _  50  51 ]
[ 55   _  35  36   _  34   _  36   _  40  39  40  41   _  45  46  47  48  49   _ ]
[ 54  53   _  37  38   _  42   _  42  41   _  41   _  49   _  47   _  49  50  51 ]
[  _  52  51   _  39  40  41  42   _  42  43   _  47  48  49   _  53   _  51   _ ]
[ 52   _  50  49   _  41  42   _  44  43  44  45  46   _  50  51  52  53   _  59 ]
[ 51  50   _  48  47   _  43  44   _  44   _  46  47  48  49   _  53   _  57  58 ]
[  _  49  48  47  46  45  44  45  46   _  48  47   _  49   _  55  54  55  56  57 ]
[ 51  50   _  48   _  46   _  46  47  48   _  48  49   _  57  56   _  56   _  58 ]

Upvotes: 1

donman
donman

Reputation: 67

Not often that somebody cares about my code. So here is the code of my MiniZinc model with two data files. I also include a snippet of the output for the smaller data file run with Gecode.

Solvers and Running Times: I typically start with Gecode, because it is the default. It also gives useful output statistics. However, when Gecode tests my patience, I try Chuffed, COIN, HiGHS and OR Tools. For the 7x7 instance these five solvers all did the job in a few 100 ms. HiGHS solved the 15x15 problem in less than half a second, COIN took 6 seconds and Chuffed needed 10 seconds. I stopped Gecode after 1 minute as it had not given any output yet. OR Tools produced some output after a few seconds, but it made very slow progress and I stopped it after 2 minutes.

----------

sum(flowremain): 34   obj: 47 

  hitorimask              hitorisoln              flowremain   
  0  1  1  1  0  1  0     .  1  5  7  .  3  .     0  2  1  1  0  1  0
  1  1  0  1  1  1  1     1  3  .  2  7  5  6     1  1  0  1  1  1  1
  0  1  1  0  1  1  0     .  2  3  .  6  7  .     0  1  1  0  1  1  0
  1  1  0  1  0  1  1     5  7  .  6  .  2  4     1  1  0  1  0  1  1
  1  0  1  1  1  1  0     2  .  7  4  1  6  .     1  0  1  1  1  1  0
  0  1  0  1  0  1  1     .  6  .  1  .  4  7     0  1  0  1  0  1  1
  1  1  1  1  1  0  1     7  4  6  3  2  .  5     1  1  1  1  1  0  0
% time elapsed: 132msec
----------

sum(flowremain): 34   obj: 49 

  hitorimask              hitorisoln              flowremain   
  0  1  1  1  0  1  0     .  1  5  7  .  3  .     0  1  1  1  0  1  0
  1  1  0  1  1  1  1     1  3  .  2  7  5  6     1  1  0  1  1  1  1
  0  1  1  0  1  1  0     .  2  3  .  6  7  .     0  1  1  0  1  1  0
  1  1  0  1  0  1  1     5  7  .  6  .  2  4     1  1  0  1  0  1  1
  1  0  1  1  1  1  0     2  .  7  4  1  6  .     1  0  1  1  1  1  0
  0  1  0  1  0  1  1     .  6  .  1  .  4  7     0  1  0  1  0  1  1
  1  1  1  1  1  0  1     7  4  6  3  2  .  5     1  1  1  1  1  0  1
% time elapsed: 133msec
----------
==========
%%%mzn-stat: failures=1716
%%%mzn-stat: initTime=0.00951658
%%%mzn-stat: nodes=3504
%%%mzn-stat: peakDepth=69
%%%mzn-stat: propagations=382613
%%%mzn-stat: propagators=917
%%%mzn-stat: restarts=0
%%%mzn-stat: solutions=34
%%%mzn-stat: solveTime=0.0152159
%%%mzn-stat: variables=731
%%%mzn-stat-end
%%%mzn-stat: nSolutions=34
%%%mzn-stat-end
Finished in 159msec.

Data files

% hitori_7x7.dzn

nsize = 7 ;  

hitoritabl = [|
5, 1, 5, 7, 7, 3, 3 |
1, 3, 3, 2, 7, 5, 6 |
2, 2, 3, 3, 6, 7, 6 |
5, 7, 7, 6, 6, 2, 4 |
2, 7, 7, 4, 1, 6, 4 |
1, 6, 6, 1, 4, 4, 7 |
7, 4, 6, 3, 2, 4, 5 |] ;


% hitori_15x15.dzn

nsize = 15 ;

hitoritabl = [|
15, 3, 1, 4, 15, 13, 1, 6, 15, 5, 7, 2, 1, 7, 11 |
8, 13, 6, 13, 4, 14, 13, 7, 11, 12, 10, 1, 15, 5, 1 |
6, 1, 15, 8, 6, 15, 4, 6, 13, 14, 15, 3, 5, 6, 9 |
11, 12, 3, 12, 12, 1, 12, 15, 14, 12, 5, 12, 8, 2, 4 |
6, 8, 3, 1, 13, 6, 5, 9, 3, 15, 3, 4, 10, 11, 9 |
10, 6, 11, 10, 7, 9, 5, 5, 3, 6, 2, 15, 10, 4, 8 |
12, 11, 7, 2, 12, 6, 15, 4, 7, 9, 13, 7, 3, 12, 10 |
3, 11, 4, 11, 5, 11, 9, 11, 6, 2, 11, 14, 2, 10, 11 |
6, 5, 6, 9, 6, 12, 6, 14, 6, 10, 15, 6, 11, 6, 3 |
4, 14, 2, 13, 15, 5, 11, 3, 7, 8, 1, 10, 5, 9, 3 |
11, 2, 14, 3, 8, 10, 11, 1, 11, 7, 3, 11, 13, 3, 6 |
2, 4, 12, 15, 1, 12, 6, 12, 8, 12, 11, 5, 12, 14, 12 |
5, 13, 12, 14, 13, 8, 7, 11, 13, 4, 9, 13, 2, 13, 15 |
14, 15, 5, 12, 9, 4, 12, 10, 2, 11, 1, 13, 1, 12, 7 |
13, 7, 13, 6, 4, 3, 2, 4, 10, 13, 14, 3, 4, 15, 4 |] ;

and finally my MiniZinc model


%% hitori.mzn


% Given square a grid of sidelength n
% The grid is filled with numbers 1..n, some of them repeated in rows and columns
% Goal: Cover some of the repeated numbers so that
% i)    there are no duplicates in any row of column
% ii)   no covered square is vertically of horizontally adjacent to another coverde square
% iii)  the uncovered blocks are connected 


% data instance defined in *.dzn 
int: nsize ;
array[SDX,SDX] of SDX: hitoritabl ;

int : maxinputflow = nsize*nsize ;
var 1..2 : toprowstartcol ;
var int : sumflowremain ;
var int: inputflow = sum(hitorimask) ;


set of int: SDX = 1..nsize ;
set of int: SDXmin1 = 1..nsize-1  ;              % index set over pairs of adjacent cells  
set of int: SDX0 = 0..nsize ;
set of int: SDX2 = 2..nsize ;


% primary unknowns, 0 or 1 whether a square is covered or open 
array[SDX,SDX] of var 0..1 : hitorimask ;

% modelling variable to display the solution, 0 is used for a covered/black cell
array[SDX,SDX] of var SDX0: hitorisoln ;
constraint forall( ir in SDX, jc in SDX) ( hitorisoln[ir,jc] = hitorimask[ir,jc] * hitoritabl[ir,jc] ) ;

% model constraint
% use primary variables to construct quantity on which nonrepeat constraint is imposed
include "alldifferent_except_0.mzn" ; 
constraint forall ( ir in SDX) ( alldifferent_except_0 ( [ hitorisoln[ir,jc]| jc in SDX ] ) ) ;
constraint forall ( jc in SDX) ( alldifferent_except_0 ( [ hitorisoln[ir,jc]| ir in SDX ] ) ) ;

% impose condition on non-adjacent masked cells
% no black squares directly below another black square
constraint forall (ir in SDX2, jc in SDX) ( 0 < hitorimask[ir-1,jc] + hitorimask[ir,jc] ) ;
% no black square directly to the right of another black square
constraint forall (ir in SDX, jc in SDX2) ( 0 < hitorimask[ir,jc-1] + hitorimask[ir,jc] ) ;


% this constraint weeds out numerous sibgle isolated open cells
constraint forall (ir in SDX, jc in SDX) 
           ( if hitorimask[ir,jc]=1 then 
                ( if ir > 1     then hitorimask[ir-1,jc  ] else 0 endif 
                + if jc < nsize then hitorimask[ir  ,jc+1] else 0 endif
                + if ir < nsize then hitorimask[ir+1,jc  ] else 0 endif
                + if jc > 1     then hitorimask[ir  ,jc-1] else 0 endif
                 > 0 
                )
             endif 
           ) ;


% the following implements the connectedness condition of the open cells
array[SDX,SDX] of var -maxinputflow..maxinputflow : xtop ;           % inflows from top     v
array[SDX,SDX] of var -maxinputflow..maxinputflow : xright ;         % inflows from right  <-
array[SDX,SDX] of var -maxinputflow..maxinputflow : xbottom ;        % inflows from bottom  ^
array[SDX,SDX] of var -maxinputflow..maxinputflow : xleft ;          % inflows from left   ->

array[SDX,SDX] of var 0..maxinputflow: flowremain ;

% initialize top, right, bottom and left rows and columns of xl, xr, xb and xl
constraint if hitorimask[1,1] = 0                                    % then necessarily hitorimask[1,2]=1
           then toprowstartcol = 2 /\ xtop[1,2]=inputflow /\ xtop[1,1] = 0 
           else toprowstartcol = 1 /\ xtop[1,1]=inputflow 
           endif ;           
constraint forall (jc in 3..nsize) (xtop[1,jc] = 0) ;   % make zero the inflow from the top for the rest of the top row
constraint forall (ir in SDX) (xright[ir,nsize] = 0) ; 
constraint forall (jc in SDX) (xbottom[nsize,jc] = 0) ; 
constraint forall (ir in SDX) (xleft[ir,1] = 0) ; 

% now make zero the inflow from all directions for all black cells
constraint forall (ir in SDX, jc in SDX) 
           ( if hitorimask[ir, jc] = 0
             then (xtop[ir,jc]=0) /\ (xright[ir,jc]=0) /\ (xbottom[ir,jc]=0) /\ (xleft[ir,jc]=0)
             else true
             endif
           ) ;

% flow balance across cell boundaries, first horizontally then vertically
constraint forall (ir in SDXmin1, jc in SDX) ( xbottom[ir,jc] = -xtop[ir+1,jc]  ) ;
constraint forall (ir in SDX, jc in SDXmin1) ( xright[ir,jc]  = -xleft[ir,jc+1] ) ;
                             
% sum of flows into a cell  
constraint forall (ir in SDX, jc in SDX) 
          (flowremain[ir,jc] = xtop[ir,jc] + xright[ir,jc] + xbottom[ir,jc] + xleft[ir,jc] ) ;
          
% objective function
var int: obj ;
constraint obj = sum (ir in SDX, jc in SDX) (hitorimask[ir,jc]=flowremain[ir,jc]) ;


solve maximize obj ;

output["\nsum(flowremain): \(sum(flowremain))   obj: \(obj) " ]  ++ ["\n"] ;


output ["\n  hitorimask              hitorisoln              flowremain   \n" ] ;
output [ 
         if jc <= nsize then "" ++ show_int(3, hitorimask[ir,jc] ) else "" endif
      ++ if jc = nsize then "   " else "" endif 
      ++ if nsize < jc /\ jc <= 2*nsize 
         then if fix(hitorisoln[ir,jc-nsize])>0 
              then show_int(3, hitorisoln[ir,jc-nsize] )  
              else "  ." 
              endif  
         else "" 
         endif 
      ++ if jc = 2*nsize then "   " else "" endif 
      ++ if 2*nsize < jc /\ jc <= 3*nsize 
         then show_int(3, flowremain[ir,jc-2*nsize] ) 
         else "" 
         endif 
      ++ if jc = 3*nsize then "\n" else "" endif
      | ir in SDX, jc in 1..3*nsize ] ;


% the following displays first the four matrices xt, xr, xb and xl and then their sum flowremain
% set of int: SD5n = 1..5*nsize ;                  % index set along columns of printout

%output ["row          xtop            x <-          xbottom           x ->        flowremain  "  ] ++ ["\n"] ;
%
%output [if jc=1 then show_int(3,ir) ++ "  " else "" endif
%        ++ if jc <= nsize then show_int(3,xtop[ir,jc]) else "" endif
%
%        ++ if jc = nsize then  "   " else "" endif 
%        ++ if nsize < jc /\ jc <= 2*nsize then show_int(3, xright[ir,jc-nsize]) else "" endif 
%        
%        ++ if jc = 2*nsize then  "   " else "" endif 
%        ++ if 2*nsize < jc /\ jc <= 3*nsize then show_int(3,xbottom[ir,jc-2*nsize]) else "" endif    
%        
%        ++ if jc = 3*nsize then  "   " else "" endif 
%        ++ if 3*nsize < jc /\ jc <= 4*nsize then show_int(3,xleft[ir,jc-3*nsize]) else "" endif    
%        
%        ++ if jc = 4*nsize then  "   " else "" endif 
%        ++ if 4*nsize < jc /\ jc <= 5*nsize then show_int(3,flowremain[ir,jc-4*nsize]) else "" endif    
%        
%        ++ if jc = 5*nsize then  "\n" else "" endif 
%           | ir in SDX, jc in SD5n ] ++ ["\n"] ;
          

Upvotes: 0

donman
donman

Reputation: 67

Your suggestion of flow into the network from the top-left corner works great!

I implemented it using four tables for the flows into a cell, from the top, the right, the bottom and from the cell to the left. These flows are constrained by zero flow through the outside boundaries, except for the single flow at the top left or second from left, by continuity across interior boundaries and zero flows through the boundaries of the black cells. I then “discovered” that the condition of a single unit of flow remaining in each uncovered cell and the zero flow condition for the covered cells correspond exactly to the mask I used to identify covered and uncovered cells. Then, by maximizing the number of cells where the remaining flow corresponds to whether the cells is covered or not, one finds the solution.

I tested my MiniZinc model on several instances, including a 15x15 hard one. Most of the solvers find the unique solution within half a second.

BTW, you must have quite some experience with network flow models to have realised that the connectedness condition could be modelled as network flow.

So, with great appreciation, thank you for your suggestions.

Upvotes: 0

Related Questions