Well, I haven't yet thought of a P-class solution, but it did occur to me that this problem may be a good candidate for stochastic solutions.
Notably, there's an easily-defined feasible starting point: just set all of the cover rectangles to the extents of the bounding box of the target squares.
From this initial state, new valid states can be generated by reducing one of the bounds of the cover rectangles and checking to see that all of the target squares are still covered.
Further, the path between any two states is likely to be short (each rectangle can be reduced to its appropriate dimension in O(√n) time, where n is the number of squares in the bounding box), meaning that it's easy to move around the search space. Though this comes with the caveat that some possible solutions are separated by a narrow path back through the initial state, which means that rerunning the algorithm we're about to develop a few times is probably good.
Given the foregoing, simulated annealing is a possible means to address the problem. The following Python script implements it:
#!/usr/bin/env python3
import random
import numpy as np
import copy
import math
import scipy
import scipy.optimize
#Generate a grid
class Grid:
def __init__(self,grid_array):
self.grid = np.array(grid_array)
self.width = len(self.grid[0]) #Use inclusive coordinates
self.height = len(self.grid) #Use inclusive coordinates
#Convert into a list of cells
self.cells = {}
for y in range(len(self.grid)):
for x in range(len(self.grid[y])):
self.cells[(x,y)] = self.grid[y][x]
#Find all cells which are border cells (the ones we need covered)
self.borders = []
for c in self.cells:
for dx in [-1,0,1]: #Loop through neighbors
for dy in [-1,0,1]:
n = (c[0]+dx,c[1]+dy) #This is the neighbor
if self.cells[c]==1 and self.cells.get(n, 1)==0: #See if this cell has a neighbor with value 0. Use default return to simplify code
self.borders.append(c)
#Ensure grid contains only valid target cells
self.grid = np.zeros((self.height,self.width))
for b in self.borders:
self.grid[b[1],b[0]] = 1
self.ntarget = np.sum(self.grid)
def copy(self):
return self.grid.copy()
#A state is valid if the bounds of each rectangle are inside the bounding box of
#the target squares and all the target squares are covered.
def ValidState(rects):
#Check bounds
if not (np.all(0<=rects[0::4]) and np.all(rects[0::4]<g.width)): #x
return False
if not (np.all(0<=rects[1::4]) and np.all(rects[1::4]<g.height)): #y
return False
if not (np.all(0<=rects[2::4]) and np.all(rects[2::4]<=g.width)): #w
return False
if not (np.all(0<=rects[3::4]) and np.all(rects[3::4]<=g.height)): #h
return False
fullmask = np.zeros((g.height,g.width))
for r in range(0,len(rects),4):
fullmask[rects[r+1]:rects[r+3],rects[r+0]:rects[r+2]] = 1
return np.sum(fullmask * g.grid)==g.ntarget
#Mutate a randomly chosen bound of a rectangle. Keep trying this until we find a
#mutation that leads to a valid state.
def MutateRects(rects):
current_state = rects.copy()
while True:
rects = current_state.copy()
c = random.randint(0,len(rects)-1)
rects[c] += random.randint(-1,1)
if ValidState(rects):
return rects
#Determine the score of a state. The score is the sum of the number of times
#each empty space is covered by a rectangle. The best solutions will minimize
#this count.
def EvaluateState(rects):
score = 0
invgrid = -(g.grid-1) #Turn zeros into ones, and ones into zeros
for r in range(0,len(rects),4):
mask = np.zeros((g.height,g.width))
mask[rects[r+1]:rects[r+3],rects[r+0]:rects[r+2]] = 1
score += np.sum(mask * invgrid)
return score
#Print the list of rectangles (useful for showing output)
def PrintRects(rects):
for r in range(0,len(rects),4):
mask = np.zeros((g.height,g.width))
mask[rects[r+1]:rects[r+3],rects[r+0]:rects[r+2]] = 1
print(mask)
#Input grid is here
gridi = [[0,0,1,0,0],
[0,1,1,1,0],
[1,1,0,1,1],
[0,1,1,1,0],
[0,1,0,1,0]]
g = Grid(gridi)
#Number of rectangles we wish to solve with
rect_count = 2
#A rectangle is defined as going from (x,y)-(w,h) where (w,h) is an upper bound
#on the array coordinates. This allows efficient manipulation of rectangles as
#numpy arrays
rects = []
for r in range(rect_count):
rects += [0,0,g.width,g.height]
rects = np.array(rects)
#Might want to run a few times since the initial state is something of a
#bottleneck on moving around the search space
sols = []
for i in range(10):
#Use simulated annealing to solve the problem
sols.append(scipy.optimize.basinhopping(
func = EvaluateState,
take_step = MutateRects,
x0 = rects,
disp = True,
niter = 3000
))
#Get a minimum solution and display it
PrintRects(min(sols, key=lambda x: x['lowest_optimization_result']['fun'])['x'])
Here's a display of the algorithm's progress for the the ten runs I specify in my example code above as a function of the number of iterations (I've added some jitter so you can see all of the lines):
You'll note that most (8/10) of the runs find the minima at 8 early on. Likewise, of the 6/10 runs that find the minima at 5, most of them do so early on. This suggests that it may be better to run many shorter searches rather than a few long searches. Choosing appropriate lengths and numbers of runs will be a matter of experimentation.
Note that EvaluateState adds points for each time an empty square is covered by a rectangle. This disincentivizes redundant coverage which may be necessary to find a solution or may result in getting to a solution faster. It's pretty common for cost functions to include this sort of thing. Experimenting with a cost function that directly asks for what you want is easy - just replace EvaluateState as follows:
#Determine the score of a state. The score is the sum of the number of times
#each empty space is covered by a rectangle. The best solutions will minimize
#this count.
def EvaluateState(rects):
score = 0
invgrid = -(g.grid-1) #Turn zeros into ones, and ones into zeros
mask = np.zeros((g.height,g.width))
for r in range(0,len(rects),4):
mask[rects[r+1]:rects[r+3],rects[r+0]:rects[r+2]] = 1
score += np.sum(mask * invgrid)
return score
Using this cost function does seem to produce better results in this case:
That may be because it provides more transition paths for rectangles between feasible states. But I'd keep the other function in mind if you hit difficulties.
