pythoncombinationsor-toolsconstraint-programmingcp-sat
Finding all the combinations of free polyominoes within a specific area with a SAT-solver (Python)
I am new to the world of SAT solvers and would need some guidance regarding the following problem.
Considering that:
❶ I have a selection of 14 adjacent cells in a 4*4 grid
❷ I have 5 polyominoes (A, B, C, D, E) of sizes 4, 2, 5, 2 and 1
❸ these polyominoes are free, i.e. their shape is not fixed and can form different patterns
How can I compute all the possible combinations of these 5 free polyominoes inside the selected area (cells in grey) with a SAT-solver ?
Borrowing both from @spinkus's insightful answer and the OR-tools documentation I could make the following example code (runs in a Jupyter Notebook):
from ortools.sat.python import cp_model
import numpy as np
import more_itertools as mit
import matplotlib.pyplot as plt
%matplotlib inline
W, H = 4, 4 #Dimensions of grid
sizes = (4, 2, 5, 2, 1) #Size of each polyomino
labels = np.arange(len(sizes)) #Label of each polyomino
colors = ('#FA5454', '#21D3B6', '#3384FA', '#FFD256', '#62ECFA')
cdict = dict(zip(labels, colors)) #Color dictionary for plotting
inactiveCells = (0, 1) #Indices of disabled cells (in 1D)
activeCells = set(np.arange(W*H)).difference(inactiveCells) #Cells where polyominoes can be fitted
ranges = [(next(g), list(g)[-1]) for g in mit.consecutive_groups(activeCells)] #All intervals in the stack of active cells
def main():
model = cp_model.CpModel()
#Create an Int var for each cell of each polyomino constrained to be within Width and Height of grid.
pminos = [[] for s in sizes]
for idx, s in enumerate(sizes):
for i in range(s):
pminos[idx].append([model.NewIntVar(0, W-1, 'p%i'%idx + 'c%i'%i + 'x'), model.NewIntVar(0, H-1, 'p%i'%idx + 'c%i'%i + 'y')])
#Define the shapes by constraining the cells relative to each other
## 1st polyomino -> tetromino ##
# #
# #
# # #
# ### #
# #
################################
p0 = pminos[0]
model.Add(p0[1][0] == p0[0][0] + 1) #'x' of 2nd cell == 'x' of 1st cell + 1
model.Add(p0[2][0] == p0[1][0] + 1) #'x' of 3rd cell == 'x' of 2nd cell + 1
model.Add(p0[3][0] == p0[0][0] + 1) #'x' of 4th cell == 'x' of 1st cell + 1
model.Add(p0[1][1] == p0[0][1]) #'y' of 2nd cell = 'y' of 1st cell
model.Add(p0[2][1] == p0[1][1]) #'y' of 3rd cell = 'y' of 2nd cell
model.Add(p0[3][1] == p0[1][1] - 1) #'y' of 3rd cell = 'y' of 2nd cell - 1
## 2nd polyomino -> domino ##
# #
# #
# # #
# # #
# #
#############################
p1 = pminos[1]
model.Add(p1[1][0] == p1[0][0])
model.Add(p1[1][1] == p1[0][1] + 1)
## 3rd polyomino -> pentomino ##
# #
# ## #
# ## #
# # #
# #
################################
p2 = pminos[2]
model.Add(p2[1][0] == p2[0][0] + 1)
model.Add(p2[2][0] == p2[0][0])
model.Add(p2[3][0] == p2[0][0] + 1)
model.Add(p2[4][0] == p2[0][0])
model.Add(p2[1][1] == p2[0][1])
model.Add(p2[2][1] == p2[0][1] + 1)
model.Add(p2[3][1] == p2[0][1] + 1)
model.Add(p2[4][1] == p2[0][1] + 2)
## 4th polyomino -> domino ##
# #
# #
# # #
# # #
# #
#############################
p3 = pminos[3]
model.Add(p3[1][0] == p3[0][0])
model.Add(p3[1][1] == p3[0][1] + 1)
## 5th polyomino -> monomino ##
# #
# #
# # #
# #
# #
###############################
#No constraints because 1 cell only
#No blocks can overlap:
block_addresses = []
n = 0
for p in pminos:
for c in p:
n += 1
block_address = model.NewIntVarFromDomain(cp_model.Domain.FromIntervals(ranges),'%i' % n)
model.Add(c[0] + c[1] * W == block_address)
block_addresses.append(block_address)
model.AddAllDifferent(block_addresses)
#Solve and print solutions as we find them
solver = cp_model.CpSolver()
solution_printer = SolutionPrinter(pminos)
status = solver.SearchForAllSolutions(model, solution_printer)
print('Status = %s' % solver.StatusName(status))
print('Number of solutions found: %i' % solution_printer.count)
class SolutionPrinter(cp_model.CpSolverSolutionCallback):
''' Print a solution. '''
def __init__(self, variables):
cp_model.CpSolverSolutionCallback.__init__(self)
self.variables = variables
self.count = 0
def on_solution_callback(self):
self.count += 1
plt.figure(figsize = (2, 2))
plt.grid(True)
plt.axis([0,W,H,0])
plt.yticks(np.arange(0, H, 1.0))
plt.xticks(np.arange(0, W, 1.0))
for i, p in enumerate(self.variables):
for c in p:
x = self.Value(c[0])
y = self.Value(c[1])
rect = plt.Rectangle((x, y), 1, 1, fc = cdict[i])
plt.gca().add_patch(rect)
for i in inactiveCells:
x = i%W
y = i//W
rect = plt.Rectangle((x, y), 1, 1, fc = 'None', hatch = '///')
plt.gca().add_patch(rect)
The problem is that I have hard-coded 5 unique/fixed polyominoes and I don't know to how define the constraints so as each possible pattern for each polyomino is taken into account (provided it is possible).
Solution
One relatively straightforward way to constrain a simply connected region in OR-Tools is to constrain its border to be a circuit. If all your polyominos are to have size less than 8, we don’t need to worry about non-simply connected ones.
This code finds all 3884 solutions:
from ortools.sat.python import cp_model
cells = {(x, y) for x in range(4) for y in range(4) if x > 1 or y > 0}
sizes = [4, 2, 5, 2, 1]
num_polyominos = len(sizes)
model = cp_model.CpModel()
# Each cell is a member of one polyomino
member = {
(cell, p): model.NewBoolVar(f"member{cell, p}")
for cell in cells
for p in range(num_polyominos)
}
for cell in cells:
model.Add(sum(member[cell, p] for p in range(num_polyominos)) == 1)
# Each polyomino contains the given number of cells
for p, size in enumerate(sizes):
model.Add(sum(member[cell, p] for cell in cells) == size)
# Find the border of each polyomino
vertices = {
v: i
for i, v in enumerate(
{(x + i, y + j) for x, y in cells for i in [0, 1] for j in [0, 1]}
)
}
edges = [
edge
for x, y in cells
for edge in [
((x, y), (x + 1, y)),
((x + 1, y), (x + 1, y + 1)),
((x + 1, y + 1), (x, y + 1)),
((x, y + 1), (x, y)),
]
]
border = {
(edge, p): model.NewBoolVar(f"border{edge, p}")
for edge in edges
for p in range(num_polyominos)
}
for (((x0, y0), (x1, y1)), p), border_var in border.items():
left_cell = ((x0 + x1 + y0 - y1) // 2, (y0 + y1 - x0 + x1) // 2)
right_cell = ((x0 + x1 - y0 + y1) // 2, (y0 + y1 + x0 - x1) // 2)
left_var = member[left_cell, p]
model.AddBoolOr([border_var.Not(), left_var])
if (right_cell, p) in member:
right_var = member[right_cell, p]
model.AddBoolOr([border_var.Not(), right_var.Not()])
model.AddBoolOr([border_var, left_var.Not(), right_var])
else:
model.AddBoolOr([border_var, left_var.Not()])
# Each border is a circuit
for p in range(num_polyominos):
model.AddCircuit(
[(vertices[v0], vertices[v1], border[(v0, v1), p]) for v0, v1 in edges]
+ [(i, i, model.NewBoolVar(f"vertex_loop{v, p}")) for v, i in vertices.items()]
)
# Print all solutions
x_range = range(min(x for x, y in cells), max(x for x, y in cells) + 1)
y_range = range(min(y for x, y in cells), max(y for x, y in cells) + 1)
solutions = 0
class SolutionPrinter(cp_model.CpSolverSolutionCallback):
def OnSolutionCallback(self):
global solutions
solutions += 1
for y in y_range:
print(
*(
next(
p
for p in range(num_polyominos)
if self.Value(member[(x, y), p])
)
if (x, y) in cells
else "-"
for x in x_range
)
)
print()
solver = cp_model.CpSolver()
solver.SearchForAllSolutions(model, SolutionPrinter())
print("Number of solutions found:", solutions)
- What does cls() function do inside a class method?
- Load Registered Component in Azure ML for Pipeline using Python sdk v2
- Can I add message to the tqdm progressbar?
- How can I preserve the previous value to find the row that is greater than it?
- What's a correct and good way to implement __hash__()?
- Python Harvesters Image Acquisition GigeCam
- Why do I get a recursion error when the depth of the expected recursion should be way less than 999?
- How to (intermittently) skip certain cells when running IPython notebook?
- How to give jupyter cell standard input in python?
- Connect EtherCAT Device with pysoem
- Regex to substitute the next two words after a matching point
- What is the time complexity of heapq.nlargest?
- How to check if a given number is a power of two?
- How to move Jupyter notebook cells up/down using keyboard shortcut?
- cannot override sys.excepthook
- Getting Python's unittest results in a tearDown() method
- Is redis in Python asynchronous?
- How to tell Click to always show option defaults
- flask sqlalchemy with circular imports with db models
- .join with a query
- Python's argparse to show program's version with prog and version string formatting
- classifiers in scikit-learn that handle nan/null
- How does improve detection of line in images by open cv?
- Performance Advantages to Iterators?
- Problem with selecting a specific web element with Playwright in Python
- Thread long running python task, doesn't work inside Docker
- How to add a column with JSON representation of rows in Polars DataFrame?
- Pandas insert empty row at 0th position
- What should HTTP 201 response body be when responding to a POST request with large data?
- Can Python's unittest test in parallel, like nose can?