##
using Plots, LinearAlgebra, Random, StatsBase, SparseArrays
##
using Gurobi, JuMP
function sudoku_gurobi(X)
# SUDOKU_GUROBI Solve a sedoku puzzle using Gurobi for an Binary LP
#
# S=sudoku_gurobi(X)
# X has a set of initial conditions where anything ~= 1 ... 9 is an
# unknown
# S is the solved sudoku puzzle
@assert(all(size(X) == (9, 9)))
# create variables
nvars = 9*9*9
# index the variables for easier access
vars = zeros(Int64, 9, 9, 9)
ind = 1
for i = 1:9
for j = 1:9
for v = 1:9
vars[i, j, v] = ind
ind = ind + 1
end
end
end
#obj = ones(Int64, 1, nvars) # no objective, just constraints
# add known cells as constraints
A = spzeros(nvars, 0)
b = ones(1, 0)
for i = 1:9
for j = 1:9
if X[i, j] < 1 || X[i, j] > 9
# this is a free variable, ignore
else
# this is a constraint, add it
cons = zeros(nvars)
cons[vars[i, j, X[i, j]]] = 1
A = [A sparse(cons)] # append a row
b = [b 1]
end
end
end
#
# add sudoku constraints
#
# each cell only has one value
for i = 1:9
for j = 1:9
cons = zeros(Int64, nvars)
for v = 1:9
cons[vars[i, j, v]] = 1
end
A = [A sparse(cons)]
b = [b 1]
end
end
# each row has 1-9
for i = 1:9
for v = 1:9
cons = zeros(nvars)
# the value v only appears in each row once
for j = 1:9
cons[vars[i, j, v]] = 1
end
A = [A sparse(cons)]
b = [b 1]
end
end
# each column has 1-9
for j = 1:9
for v = 1:9
cons = zeros(nvars)
# the value v only appears in each row once
for i = 1:9
cons[vars[i, j, v]] = 1
end
A = [A sparse(cons)]
b = [b 1]
end
end
# subgrid has 1-9
for v = 1:9
for i0 = 0:2
for j0 = 0:2
cons = zeros(nvars)
for i = 1:3
for j = 1:3
cons[vars[i+3*i0, j+3*j0, v]] = 1
end
end
A = [A sparse(cons)]
b = [b 1]
end
end
end
# setup Gurobi call
env = Gurobi.Env()
m = Model(Gurobi.Optimizer)
# add binary variables
@variable(m, x[1:nvars], binary=true)
#add_bvars!(model, zeros(nvars))
# have the variables incorporated into the model
#update_model!(model)
#add_constrs_t!(model, A, '=', vec(b))
@constraint(m, A'*x .== vec(b))
# run optimization
optimize!(m)
# get results
x = JuMP.value.(x)
# assign solution
S = zeros(Int64, 9 ,9)
for i = 1:9
for j = 1:9
for v = 1:9
if x[vars[i, j, v]] > 0
S[i, j] = v
end
end
end
end
return S
end
X = [0 0 8 0 0 0 6 0 0;
0 2 0 0 0 0 0 1 0;
5 0 0 6 0 0 0 0 9;
0 0 7 0 4 0 0 0 3;
0 4 0 2 0 9 0 6 0;
6 0 0 0 1 0 8 0 0;
9 0 0 0 0 7 0 0 4;
0 5 0 0 0 0 0 7 0;
0 0 1 0 0 0 5 0 0];
S = sudoku_gurobi(X)
Set parameter Username Academic license - for non-commercial use only - expires 2024-01-15 Set parameter Username Academic license - for non-commercial use only - expires 2024-01-15 Gurobi Optimizer version 10.0.0 build v10.0.0rc2 (mac64[arm]) CPU model: Apple M1 Thread count: 8 physical cores, 8 logical processors, using up to 8 threads Optimize a model with 348 rows, 729 columns and 2940 nonzeros Model fingerprint: 0xe4cb3567 Variable types: 0 continuous, 729 integer (729 binary) Coefficient statistics: Matrix range [1e+00, 1e+00] Objective range [0e+00, 0e+00] Bounds range [0e+00, 0e+00] RHS range [1e+00, 1e+00] Presolve removed 348 rows and 729 columns Presolve time: 0.00s Presolve: All rows and columns removed Explored 0 nodes (0 simplex iterations) in 0.00 seconds (0.00 work units) Thread count was 1 (of 8 available processors) Solution count 1: 0 Optimal solution found (tolerance 1.00e-04) Best objective 0.000000000000e+00, best bound 0.000000000000e+00, gap 0.0000% User-callback calls 221, time in user-callback 0.00 sec
9×9 Matrix{Int64}: 1 9 8 4 7 2 6 3 5 4 2 6 5 9 3 7 1 8 5 7 3 6 8 1 4 2 9 2 1 7 8 4 6 9 5 3 8 4 5 2 3 9 1 6 7 6 3 9 7 1 5 8 4 2 9 6 2 1 5 7 3 8 4 3 5 4 9 6 8 2 7 1 7 8 1 3 2 4 5 9 6
## Get the NYTimes hard problem
# from 2023-01-20
X = [0 0 7 0 9 0 6 0 0
9 0 2 0 0 0 4 0 3
5 0 0 0 0 0 8 1 0
0 2 0 0 0 6 0 0 0
3 0 0 0 7 5 9 0 0
0 0 1 0 0 0 0 0 0
0 0 0 0 6 0 0 4 2
0 0 0 8 3 0 5 0 0
0 0 0 0 0 1 0 0 0]
S = sudoku_gurobi(X)
Set parameter Username Academic license - for non-commercial use only - expires 2024-01-15 Set parameter Username Academic license - for non-commercial use only - expires 2024-01-15 Gurobi Optimizer version 10.0.0 build v10.0.0rc2 (mac64[arm]) CPU model: Apple M1 Thread count: 8 physical cores, 8 logical processors, using up to 8 threads Optimize a model with 348 rows, 729 columns and 2940 nonzeros Model fingerprint: 0xe4655fd9 Variable types: 0 continuous, 729 integer (729 binary) Coefficient statistics: Matrix range [1e+00, 1e+00] Objective range [0e+00, 0e+00] Bounds range [0e+00, 0e+00] RHS range [1e+00, 1e+00] Presolve removed 348 rows and 729 columns Presolve time: 0.00s Presolve: All rows and columns removed Explored 0 nodes (0 simplex iterations) in 0.00 seconds (0.00 work units) Thread count was 1 (of 8 available processors) Solution count 1: 0 Optimal solution found (tolerance 1.00e-04) Best objective 0.000000000000e+00, best bound 0.000000000000e+00, gap 0.0000% User-callback calls 440, time in user-callback 0.00 sec
9×9 Matrix{Int64}: 4 8 7 1 9 3 6 2 5 9 1 2 6 5 8 4 7 3 5 3 6 4 2 7 8 1 9 7 2 5 9 1 6 3 8 4 3 4 8 2 7 5 9 6 1 6 9 1 3 8 4 2 5 7 8 5 3 7 6 9 1 4 2 1 7 4 8 3 2 5 9 6 2 6 9 5 4 1 7 3 8
## Get the NYTimes hard problem
# from 2023-01-20
X = [0 0 0 0 0 0 0 0 0
9 0 2 0 0 0 4 0 3
5 0 0 0 0 0 0 1 0
0 2 0 0 0 6 0 0 0
3 0 0 0 7 5 9 0 0
0 0 1 0 0 0 0 0 0
0 0 0 0 0 0 0 4 2
0 0 0 8 3 0 5 0 0
0 0 0 0 0 1 0 0 0]
S1 = sudoku_gurobi(X)
Set parameter Username Academic license - for non-commercial use only - expires 2024-01-15 Set parameter Username Academic license - for non-commercial use only - expires 2024-01-15 Gurobi Optimizer version 10.0.0 build v10.0.0rc2 (mac64[arm]) CPU model: Apple M1 Thread count: 8 physical cores, 8 logical processors, using up to 8 threads Optimize a model with 343 rows, 729 columns and 2935 nonzeros Model fingerprint: 0x2fbc93da Variable types: 0 continuous, 729 integer (729 binary) Coefficient statistics: Matrix range [1e+00, 1e+00] Objective range [0e+00, 0e+00] Bounds range [0e+00, 0e+00] RHS range [1e+00, 1e+00] Presolve removed 116 rows and 466 columns Presolve time: 0.01s Presolved: 227 rows, 263 columns, 1063 nonzeros Variable types: 0 continuous, 263 integer (263 binary) Found heuristic solution: objective 0.0000000 Explored 0 nodes (0 simplex iterations) in 0.01 seconds (0.01 work units) Thread count was 8 (of 8 available processors) Solution count 1: 0 Optimal solution found (tolerance 1.00e-04) Best objective 0.000000000000e+00, best bound 0.000000000000e+00, gap 0.0000% User-callback calls 397, time in user-callback 0.00 sec
9×9 Matrix{Int64}: 1 4 8 9 6 3 2 7 5 9 6 2 5 1 7 4 8 3 5 3 7 2 4 8 6 1 9 7 2 5 4 9 6 8 3 1 3 8 6 1 7 5 9 2 4 4 9 1 3 8 2 7 5 6 8 7 3 6 5 9 1 4 2 2 1 9 8 3 4 5 6 7 6 5 4 7 2 1 3 9 8
S1 == S
false
S1
9×9 Matrix{Int64}: 1 4 8 9 6 3 2 7 5 9 6 2 5 1 7 4 8 3 5 3 7 2 4 8 6 1 9 7 2 5 4 9 6 8 3 1 3 8 6 1 7 5 9 2 4 4 9 1 3 8 2 7 5 6 8 7 3 6 5 9 1 4 2 2 1 9 8 3 4 5 6 7 6 5 4 7 2 1 3 9 8
S
9×9 Matrix{Int64}: 4 8 7 1 9 3 6 2 5 9 1 2 6 5 8 4 7 3 5 3 6 4 2 7 8 1 9 7 2 5 9 1 6 3 8 4 3 4 8 2 7 5 9 6 1 6 9 1 3 8 4 2 5 7 8 5 3 7 6 9 1 4 2 1 7 4 8 3 2 5 9 6 2 6 9 5 4 1 7 3 8