ILP Part 65 – Fancy sudoku – Random IT Utensils

This is the sixtieth fifth part of the ILP series. For your convenience you can find other parts in the table of contents in Part 1 – Boolean algebra

Today we’re going to solve a modification of classical Sudoku riddle shown in Precz z cyframi.

We need to fill the board with digits as in Sudoku and meet some additional constraints:

Field with X indicates that number in that field is equal to the number of different values in corner neighbours of the field (so top-left field, top-right, bottom-left and bottom-right)
Field with + indicates that among four straight neighbours (field above, below, to the left, to the right) there are two which added values are equal to the value in the field
Field with circle indicates that field is equal to the number of different values in all neighbours

Additionally, the board is filled with all crosses, pluses and circles completely. This means that if given field is empty then no cross/field/circle can be placed in it.

Solution:

var board = new int[][]{
	new int []{0,0,1,2,3,0,2,0,0},
	new int []{0,2,0,0,4,3,0,2,2},
	new int []{0,2,0,1,2,0,0,0,0},
	new int []{1,2,0,0,2,1,0,0,0},
	new int []{2,0,0,2,2,0,2,1,3},
	new int []{0,1,1,2,0,2,0,2,0},
	new int []{4,0,2,1,0,0,2,4,0},
	new int []{0,2,0,0,1,2,1,2,1},
	new int []{4,1,0,0,2,0,3,0,0}
};

var width = board[0].Length;
var height = board.Length;

var solver = new CplexMilpSolver(new CplexMilpSolverSettings
{
	IntegerWidth = 15,
	Epsilon = 0.1,
	StoreDebugExpressions = false,
	CacheConstants = false,
	StoreDebugConstraints = false
});

Console.WriteLine("Create variables");
var fields = Enumerable.Range(0, height).Select(r => Enumerable.Range(0, width).Select(c => solver.CreateAnonymous(Domain.PositiveOrZeroInteger)).ToArray()).ToArray();

Console.WriteLine("Set ranges");
for(int row = 0; row < height;++row){
	for(int column = 0; column < width; ++ column){
		solver.Set<LessOrEqual>(fields[row][column], solver.FromConstant(9));
		solver.Set<GreaterOrEqual>(fields[row][column], solver.FromConstant(1));
	}
}

Console.WriteLine("Rows");
for(int row = 0; row < height;++row){
	solver.Set<AllDifferent>(fields[row][0], fields[row].Skip(1).ToArray());
}

Console.WriteLine("Columns");
for(int column = 0; column < width; ++ column){
	solver.Set<AllDifferent>(fields[0][column], Enumerable.Range(1, height-1).Select(row => fields[row][column]).ToArray());
}

Console.WriteLine("Quadrants");
for(int row = 0; row < height / 3;++row){
	for(int column = 0; column < width / 3; ++ column){
		var quandrant = Enumerable.Range(0, 3).SelectMany(i => Enumerable.Range(0, 3).Select(j => fields[row * 3 + i][column * 3 + j])).ToArray();
		solver.Set<AllDifferent>(quandrant[0], quandrant.Skip(1).ToArray());
	}
}

Console.WriteLine("Shapes");
for(int row = 0; row < height;++row){
	for(int column = 0; column < width; ++ column){
		var corners = Enumerable.Range(-1, 3).SelectMany(i => Enumerable.Range(-1, 3)
		.Where(j => i % 2 != 0 && j % 2 != 0)
		.Where(j => row + i >= 0 && row + i < height)
		.Where(j => column + j >= 0 && column + j < width)
		.Select(j => 
			fields[row+i][column+j]
		)).ToArray();
		
		var straights = Enumerable.Range(-1, 3).SelectMany(i => Enumerable.Range(-1, 3)
		.Where(j => (i % 2 == 0 && j % 2 != 0) || (i % 2 != 0 && j % 2 == 0))
		.Where(j => row + i >= 0 && row + i < height)
		.Where(j => column + j >= 0 && column + j < width)
		.Select(j => 
			fields[row+i][column+j]
		)).ToArray();
		
		var neighbours = corners.Concat(straights).ToArray();
		
		var allSums = Enumerable.Range(0, straights.Length).SelectMany(i => Enumerable.Range(0, straights.Length).Where(j => j > i).Select(j => solver.Operation<Addition>(straights[i], straights[j]))).ToArray();
		var isEqual = solver.Operation<Addition>(allSums.Select(s => solver.Operation<IsEqual>(fields[row][column], s)).ToArray());
		
		switch(board[row][column]){
			case 0:
				solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);
				solver.Set<Equal>(isEqual, solver.FromConstant(0));
				solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);
			break;
			case 1:
				solver.Set<Equal>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);
				solver.Set<Equal>(isEqual, solver.FromConstant(0));
				solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);
			break;
			case 2:
				solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);
				solver.Set<GreaterOrEqual>(isEqual, solver.FromConstant(1));
				solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);
			break;
			case 3:
				solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);
				solver.Set<Equal>(isEqual, solver.FromConstant(0));
				solver.Set<Equal>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);
			break;
			case 4:
				solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);
				solver.Set<GreaterOrEqual>(isEqual, solver.FromConstant(1));
				solver.Set<Equal>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);
			break;
		}
	}
}


solver.AddGoal("goal", solver.FromConstant(0));
solver.Solve();

for(int row = 0; row < height;++row){
	if(row % 3 == 0)Console.WriteLine("-------------");
	for(int column = 0; column < width; ++ column){
		if(column % 3 == 0)Console.Write("|");
		Console.Write(Math.Round(solver.GetValue(fields[row][column])));
	}
	Console.Write("|");
	Console.WriteLine();
}
Console.WriteLine("-------------");

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

var board = new int[][]{

new int []{0,0,1,2,3,0,2,0,0},

new int []{0,2,0,0,4,3,0,2,2},

new int []{0,2,0,1,2,0,0,0,0},

new int []{1,2,0,0,2,1,0,0,0},

new int []{2,0,0,2,2,0,2,1,3},

new int []{0,1,1,2,0,2,0,2,0},

new int []{4,0,2,1,0,0,2,4,0},

new int []{0,2,0,0,1,2,1,2,1},

new int []{4,1,0,0,2,0,3,0,0}

};

var width = board[0].Length;

var height = board.Length;

var solver = new CplexMilpSolver(new CplexMilpSolverSettings

{

IntegerWidth = 15,

Epsilon = 0.1,

StoreDebugExpressions = false,

CacheConstants = false,

StoreDebugConstraints = false

});

Console.WriteLine("Create variables");

var fields = Enumerable.Range(0, height).Select(r => Enumerable.Range(0, width).Select(c => solver.CreateAnonymous(Domain.PositiveOrZeroInteger)).ToArray()).ToArray();

Console.WriteLine("Set ranges");

for(int row = 0; row < height;++row){

for(int column = 0; column < width; ++ column){

solver.Set<LessOrEqual>(fields[row][column], solver.FromConstant(9));

solver.Set<GreaterOrEqual>(fields[row][column], solver.FromConstant(1));

}

Console.WriteLine("Rows");

for(int row = 0; row < height;++row){

solver.Set<AllDifferent>(fields[row][0], fields[row].Skip(1).ToArray());

}

Console.WriteLine("Columns");

for(int column = 0; column < width; ++ column){

solver.Set<AllDifferent>(fields[0][column], Enumerable.Range(1, height-1).Select(row => fields[row][column]).ToArray());

}

Console.WriteLine("Quadrants");

for(int row = 0; row < height / 3;++row){

for(int column = 0; column < width / 3; ++ column){

var quandrant = Enumerable.Range(0, 3).SelectMany(i => Enumerable.Range(0, 3).Select(j => fields[row * 3 + i][column * 3 + j])).ToArray();

solver.Set<AllDifferent>(quandrant[0], quandrant.Skip(1).ToArray());

}

Console.WriteLine("Shapes");

for(int row = 0; row < height;++row){

for(int column = 0; column < width; ++ column){

var corners = Enumerable.Range(-1, 3).SelectMany(i => Enumerable.Range(-1, 3)

.Where(j => i % 2 != 0 && j % 2 != 0)

.Where(j => row + i >= 0 && row + i < height)

.Where(j => column + j >= 0 && column + j < width)

.Select(j =>

fields[row+i][column+j]

)).ToArray();

var straights = Enumerable.Range(-1, 3).SelectMany(i => Enumerable.Range(-1, 3)

.Where(j => (i % 2 == 0 && j % 2 != 0) || (i % 2 != 0 && j % 2 == 0))

.Where(j => row + i >= 0 && row + i < height)

.Where(j => column + j >= 0 && column + j < width)

.Select(j =>

fields[row+i][column+j]

)).ToArray();

var neighbours = corners.Concat(straights).ToArray();

var allSums = Enumerable.Range(0, straights.Length).SelectMany(i => Enumerable.Range(0, straights.Length).Where(j => j > i).Select(j => solver.Operation<Addition>(straights[i], straights[j]))).ToArray();

var isEqual = solver.Operation<Addition>(allSums.Select(s => solver.Operation<IsEqual>(fields[row][column], s)).ToArray());

switch(board[row][column]){

case 0:

solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);

solver.Set<Equal>(isEqual, solver.FromConstant(0));

solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);

break;

case 1:

solver.Set<Equal>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);

solver.Set<Equal>(isEqual, solver.FromConstant(0));

solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);

break;

case 2:

solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);

solver.Set<GreaterOrEqual>(isEqual, solver.FromConstant(1));

solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);

break;

case 3:

solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);

solver.Set<Equal>(isEqual, solver.FromConstant(0));

solver.Set<Equal>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);

break;

case 4:

solver.Set<NotEqual>(solver.Operation<DifferentValuesCount>(corners), fields[row][column]);

solver.Set<GreaterOrEqual>(isEqual, solver.FromConstant(1));

solver.Set<Equal>(solver.Operation<DifferentValuesCount>(neighbours), fields[row][column]);

break;

}

solver.AddGoal("goal", solver.FromConstant(0));

solver.Solve();

for(int row = 0; row < height;++row){

if(row % 3 == 0)Console.WriteLine("-------------");

for(int column = 0; column < width; ++ column){

if(column % 3 == 0)Console.Write("|");

Console.Write(Math.Round(solver.GetValue(fields[row][column])));

}

Console.Write("|");

Console.WriteLine();

}

Console.WriteLine("-------------");

Lines 25-52 are regular Sudoku rules.

Next, in lines 54-106 we encode additional rules. First, we extract corner neighbours (lines 57-63), straight neighbours (65-71), and all neighbours together (line 73). Next, we prepare sums of pairs of straight neighbours (line 75) and then compare them with current field value.

Finally, in lines 78-103 we encode the rules. If field is of type 0 then it means field is empty so number of different values for corner neighbours must not be equal to the current value, number of straight neighbours summing to the current field value is zero (so no sum matches), and similarly for all neighbours – number of different values among them does not match current value.

Type 1 indicates X, type 2 is +, type 3 is circle, type 4 is + with circle. Notice that X with circle is not possible, nor is X with + with circle.

Output:

Tried aggregator 2 times.
MIP Presolve eliminated 5032 rows and 2922 columns.
MIP Presolve modified 32939 coefficients.
Aggregator did 915 substitutions.
Reduced MIP has 22057 rows, 10412 columns, and 68119 nonzeros.
Reduced MIP has 10331 binaries, 81 generals, 0 SOSs, and 0 indicators.
Presolve time = 0.06 sec. (63.34 ticks)
Probing fixed 2358 vars, tightened 75 bounds.
Probing changed sense of 858 constraints.
Probing time = 0.41 sec. (93.68 ticks)
Cover probing fixed 0 vars, tightened 803 bounds.
Tried aggregator 2 times.
MIP Presolve eliminated 6449 rows and 2938 columns.
MIP Presolve modified 8576 coefficients.
Aggregator did 854 substitutions.
Reduced MIP has 14754 rows, 6620 columns, and 46066 nonzeros.
Reduced MIP has 6553 binaries, 67 generals, 0 SOSs, and 0 indicators.
Presolve time = 0.13 sec. (100.36 ticks)
Probing fixed 148 vars, tightened 7 bounds.
Probing changed sense of 23 constraints.
Probing time = 0.16 sec. (42.71 ticks)
Cover probing fixed 0 vars, tightened 41 bounds.
Tried aggregator 2 times.
MIP Presolve eliminated 532 rows and 257 columns.
MIP Presolve modified 1534 coefficients.
Aggregator did 64 substitutions.
Reduced MIP has 14158 rows, 6299 columns, and 44200 nonzeros.
Reduced MIP has 6232 binaries, 67 generals, 0 SOSs, and 0 indicators.
Presolve time = 0.09 sec. (57.17 ticks)
Probing fixed 6232 vars, tightened 112 bounds.
Probing changed sense of 39 constraints.
Probing time = 0.75 sec. (187.24 ticks)
Clique table members: 18414.
MIP emphasis: balance optimality and feasibility.
MIP search method: dynamic search.
Parallel mode: deterministic, using up to 8 threads.
Root relaxation solution time = 0.02 sec. (5.76 ticks)

        Nodes                                         Cuts/
   Node  Left     Objective  IInf  Best Integer    Best Bound    ItCnt     Gap

*     0     0      integral     0        0.0000        0.0000        0    0.00%
Elapsed time = 1.88 sec. (690.31 ticks, tree = 0.00 MB, solutions = 1)

Root node processing (before b&c):
  Real time             =    1.88 sec. (690.79 ticks)
Parallel b&c, 8 threads:
  Real time             =    0.00 sec. (0.00 ticks)
  Sync time (average)   =    0.00 sec.
  Wait time (average)   =    0.00 sec.
                          ------------
Total (root+branch&cut) =    1.88 sec. (690.79 ticks)
-------------
|962|753|814|
|451|286|739|
|783|491|526|
-------------
|295|674|183|
|817|932|645|
|634|518|297|
-------------
|548|327|961|
|176|849|352|
|329|165|478|
-------------

Tried aggregator 2 times.

MIP Presolve eliminated 5032 rows and 2922 columns.

MIP Presolve modified 32939 coefficients.

Aggregator did 915 substitutions.

Reduced MIP has 22057 rows, 10412 columns, and 68119 nonzeros.

Reduced MIP has 10331 binaries, 81 generals, 0 SOSs, and 0 indicators.

Presolve time = 0.06 sec. (63.34 ticks)

Probing fixed 2358 vars, tightened 75 bounds.

Probing changed sense of 858 constraints.

Probing time = 0.41 sec. (93.68 ticks)

Cover probing fixed 0 vars, tightened 803 bounds.

Tried aggregator 2 times.

MIP Presolve eliminated 6449 rows and 2938 columns.

MIP Presolve modified 8576 coefficients.

Aggregator did 854 substitutions.

Reduced MIP has 14754 rows, 6620 columns, and 46066 nonzeros.

Reduced MIP has 6553 binaries, 67 generals, 0 SOSs, and 0 indicators.

Presolve time = 0.13 sec. (100.36 ticks)

Probing fixed 148 vars, tightened 7 bounds.

Probing changed sense of 23 constraints.

Probing time = 0.16 sec. (42.71 ticks)

Cover probing fixed 0 vars, tightened 41 bounds.

Tried aggregator 2 times.

MIP Presolve eliminated 532 rows and 257 columns.

MIP Presolve modified 1534 coefficients.

Aggregator did 64 substitutions.

Reduced MIP has 14158 rows, 6299 columns, and 44200 nonzeros.

Reduced MIP has 6232 binaries, 67 generals, 0 SOSs, and 0 indicators.

Presolve time = 0.09 sec. (57.17 ticks)

Probing fixed 6232 vars, tightened 112 bounds.

Probing changed sense of 39 constraints.

Probing time = 0.75 sec. (187.24 ticks)

Clique table members: 18414.

MIP emphasis: balance optimality and feasibility.

MIP search method: dynamic search.

Parallel mode: deterministic, using up to 8 threads.

Root relaxation solution time = 0.02 sec. (5.76 ticks)

Nodes Cuts/

Node Left Objective IInf Best Integer Best Bound ItCnt Gap

* 0 0 integral 0 0.0000 0.0000 0 0.00%

Elapsed time = 1.88 sec. (690.31 ticks, tree = 0.00 MB, solutions = 1)

Root node processing (before b&c):

Real time = 1.88 sec. (690.79 ticks)

Parallel b&c, 8 threads:

Real time = 0.00 sec. (0.00 ticks)

Sync time (average) = 0.00 sec.

Wait time (average) = 0.00 sec.

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

Total (root+branch&cut) = 1.88 sec. (690.79 ticks)

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

|962|753|814|

|451|286|739|

|783|491|526|

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

|295|674|183|

|817|932|645|

|634|518|297|

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

|548|327|961|

|176|849|352|

|329|165|478|

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