用正方形优化填充网格图形 [英] Optimal filling of grid figure with squares

查看:106
本文介绍了用正方形优化填充网格图形的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

最近,我为孩子们设计了一个拼图.但是,我现在想找到最佳解决方案.

recently I have designed a puzzle for children to solve. However I would like to now the optimal solution.

问题如下:您已经用小方块组成了这个数字

The problem is as follows: You have this figure made up off small squares

您必须用较大的正方形填充它,并在下表中对其进行评分:

You have to fill it in with larger squares and it is scored with the following table:

| Square Size | 1x1 | 2x2 | 3x3 | 4x4 | 5x5 | 6x6 | 7x7 | 8x8 |
|-------------|-----|-----|-----|-----|-----|-----|-----|-----|
| Points      | 0   | 4   | 10  | 20  | 35  | 60  | 84  | 120 |

有很多可能的解决方案可以对它们全部进行检查.其他人建议使用动态编程,但是我不知道如何将图形划分为较小的图形,这些较小的图形具有相同的最佳解决方案.

There are simply to many possible solutions to check them all. Some other people suggested dynamic programming, but I don't know how to divide the figure in smaller ones which put together have the same optimal solution.

我想找到一种在合理的时间内(例如,在常规台式机上最多几天)为这些问题找到最佳解决方案的方法.到目前为止,通过猜测算法和一些手工工作得出的最高分是1112.

I would like to find a way to find the optimal solutions to these kinds of problems in reasonable time (like a couple of days max on a regular desktop). The highest score found so far with a guessing algorithm and some manual work is 1112.

通过组合子问题来解决类似问题的解决方案也受到赞赏. 我不需要写出所有代码.一个算法的提纲或想法就足够了.

Solutions to similar problems with combining sub-problems are also appreciated. I don't need all the code written out. An outline or idea for an algorithm would be enough.

注意:可以容纳的最大正方形是8x8,因此不包括较大正方形的分数.

Note: The biggest square that can fit is 8x8 so scores for bigger squares are not included.

[[1,1,0,0,0,1,0,0,0,0,0,0,1,1,1,1,1,1,0,0,1,1,1,1,1,0,0,1,1,1],
 [1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,1,1,0,0,0,1,1],
 [1,0,0,0,0,0,0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,1],
 [0,0,0,1,1,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0],
 [0,0,0,0,1,1,0,0,0,0,1,0,0,0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0],
 [0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,0,0,0,0,0,0,1,1,1],
 [0,0,0,0,0,0,0,0,1,1,0,0,0,0,1,1,1,1,1,1,1,1,0,0,0,0,1,1,1,1],
 [1,0,0,0,0,0,0,1,1,1,1,0,0,0,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,1],
 [1,1,0,0,0,0,0,1,1,1,1,0,0,0,1,1,1,1,1,1,1,1,0,0,1,0,0,0,0,1],
 [1,1,1,0,0,0,0,1,1,1,1,1,0,0,1,1,1,1,1,1,1,0,0,0,1,1,1,0,0,0],
 [0,1,1,1,0,0,0,1,1,1,1,1,0,0,0,0,1,1,1,0,0,0,0,1,1,1,1,0,0,0],
 [0,0,1,1,1,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0],
 [0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0],
 [0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1],
 [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1],
 [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,1],
 [0,0,0,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,1],
 [0,0,1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0],
 [0,1,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0],
 [0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0],
 [0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1],
 [0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1],
 [0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1],
 [0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1],
 [0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1],
 [0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1],
 [0,0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1],
 [0,0,0,0,0,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1],
 [0,0,0,0,0,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1],
 [1,1,1,0,0,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1],
 [1,1,1,0,0,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1],
 [1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1],
 [1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1],
 [1,1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1],
 [1,1,0,0,0,1,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,1,1,1,1],
 [1,0,0,0,0,1,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,1,1,1,1],
 [1,0,0,0,0,1,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,1,1,0,0,1,1,1,1,1],
 [1,0,0,0,0,1,0,0,0,1,1,1,0,0,0,0,0,0,0,0,1,1,1,0,0,1,1,1,1,1],
 [0,0,0,0,0,1,0,0,0,1,1,1,1,1,1,0,0,0,0,1,1,1,1,0,0,0,1,1,1,1],
 [0,0,0,0,0,1,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,1,1,1,1],
 [0,0,0,0,0,1,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,1,1,1,1]];

推荐答案

这是一个使用 Mixed-integer-programming 的通用原型,它可以优化解决您的实例(i就像推论自己一样获得了1112的值),也可能会解决其他问题.

Here is a quite general prototype using Mixed-integer-programming which solves your instance optimally (i obtained the value of 1112 like you deduced yourself) and might solve others too.

通常,您的问题是 np-complete ,这使您很难(在某些情况下会遇到麻烦).

In general, your problem is np-complete and this makes it hard (there are some instances which will be trouble).

虽然我怀疑基于 SAT-solver CP-solver 的方法可能更强大(由于具有组合性,我甚至对MIP在这里工作感到惊讶) ,MIP方法也有一些优点:

While i suspect that SAT-solver and CP-solver based approaches might be more powerful (because of the combinatoric nature; i even was surprised that MIP works here), the MIP-approach has also some advantages:

  • MIP求解器是完整(如SAT和CP;但是许多随机的-不是基于启发式的):
  • 如有需要,有许多商业级求解器可用
  • 公式相当简单(尤其是与SAT相比; SAT公式将需要n次公式中最多k个高级(用于计分公式),而这些公式正在逐步发展(二次方程式)方法成倍增长)!它们确实存在,但并非无关紧要)
  • 优化目标是自然的(SAT和CP需要迭代精炼=使用一些下界进行求解;增量约束并重新求解)
  • MIP求解器还可以非常强大地获得最佳解的近似值,并提供一些经过验证的界线(例如,最佳解小于x)
  • MIP-solvers are complete (as SAT and CP; but many random-based heuristics are not):
  • There are many commercial-grade solvers available if needed
  • The formulation is quite easy (especially compared to SAT; SAT-formulations will need advanced at most k out of n-formulations (for scoring-formulations) which are growing sub-quadratic (the naive approach grows exponentially)! They do exist, but are non-trivial)
  • The optimization-objective is just natural (SAT and CP would need iterative-refining = solve with some lower-bound; increment bound and re-solve)
  • MIP-solvers can also be quite powerful to obtain approximations of the optimal solution and also provide some proven bounds (e.g. optimum lower than x)

以下代码是使用可用的常见科学工具(全部为开源)在 python 中实现的.它允许设置图块范围(例如添加9x9图块)和不同的成本函数.评论应该足以理解这些想法.它将使用一些(可能是最好的)开源MIP求解器,但也可以使用商业性的(注释过的行显示用法).

The following code is implemented in python using common scientific tools available (all of these are open-source). It allows setting the tile-range (e.g. adding 9x9 tiles) and different cost-functions. The comments should be enough to understand the ideas. It will use some (probably the best) open-source MIP-solver, but can also use commercial ones (outcommented line shows usage).

import numpy as np
import itertools
from collections import defaultdict
import matplotlib.pyplot as plt     # visualization only
import seaborn as sns               # ""
from pulp import *                  # MIP-modelling & solver

""" INSTANCE """
instance = np.asarray([[1,1,0,0,0,1,0,0,0,0,0,0,1,1,1,1,1,1,0,0,1,1,1,1,1,0,0,1,1,1],
 [1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,1,1,0,0,0,1,1],
 [1,0,0,0,0,0,0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,1],
 [0,0,0,1,1,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0],
 [0,0,0,0,1,1,0,0,0,0,1,0,0,0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0],
 [0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,0,0,0,0,0,0,1,1,1],
 [0,0,0,0,0,0,0,0,1,1,0,0,0,0,1,1,1,1,1,1,1,1,0,0,0,0,1,1,1,1],
 [1,0,0,0,0,0,0,1,1,1,1,0,0,0,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,1],
 [1,1,0,0,0,0,0,1,1,1,1,0,0,0,1,1,1,1,1,1,1,1,0,0,1,0,0,0,0,1],
 [1,1,1,0,0,0,0,1,1,1,1,1,0,0,1,1,1,1,1,1,1,0,0,0,1,1,1,0,0,0],
 [0,1,1,1,0,0,0,1,1,1,1,1,0,0,0,0,1,1,1,0,0,0,0,1,1,1,1,0,0,0],
 [0,0,1,1,1,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0],
 [0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0],
 [0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1],
 [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1],
 [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,1],
 [0,0,0,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,1],
 [0,0,1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0],
 [0,1,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0],
 [1,1,1,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0],
 [0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0],
 [0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1],
 [0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1],
 [0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1],
 [0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1],
 [0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1],
 [0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1],
 [0,0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1],
 [0,0,0,0,0,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1],
 [0,0,0,0,0,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1],
 [1,1,1,0,0,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1],
 [1,1,1,0,0,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1],
 [1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1],
 [1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1],
 [1,1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1],
 [1,1,0,0,0,1,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,1,1,1,1],
 [1,0,0,0,0,1,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,1,1,1,1],
 [1,0,0,0,0,1,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,1,1,0,0,1,1,1,1,1],
 [1,0,0,0,0,1,0,0,0,1,1,1,0,0,0,0,0,0,0,0,1,1,1,0,0,1,1,1,1,1],
 [0,0,0,0,0,1,0,0,0,1,1,1,1,1,1,0,0,0,0,1,1,1,1,0,0,0,1,1,1,1],
 [0,0,0,0,0,1,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,1,1,1,1],
 [0,0,0,0,0,1,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,1,1,1,1]], dtype=bool)

def plot_compare(instance, solution, subgrids):
    f, (ax1, ax2) = plt.subplots(2, sharex=True, sharey=True)
    sns.heatmap(instance, ax=ax1, cbar=False, annot=True)
    sns.heatmap(solution, ax=ax2, cbar=False, annot=True)
    plt.show()

""" PARAMETERS """
SUBGRIDS = 8  # 1x1 - 8x8
SUGBRID_SCORES = {1:0, 2:4, 3:10, 4:20, 5:35, 6:60, 7:84, 8:120}
N, M = instance.shape  # free / to-fill = zeros!

""" HELPER FUNCTIONS """
def get_square_covered_indices(instance, pos_x, pos_y, sg):
    """ Calculate all covered tiles when given a top-left position & size
            -> returns the base-index too! """
    N, M = instance.shape
    neighbor_indices = []
    valid = True
    for sX in range(sg):
        for sY in range(sg):
            if pos_x + sX < N:
                if pos_y + sY < M:
                    if instance[pos_x + sX, pos_y + sY] == 0:
                        neighbor_indices.append((pos_x + sX, pos_y + sY))
                    else:
                        valid = False
                        break
                else:
                    valid = False
                    break
            else:
                valid = False
                break
    return valid, neighbor_indices

def preprocessing(instance, SUBGRIDS):
    """ Calculate all valid placement / tile-selection combinations """
    placements = {}
    index2placement = {}
    placement2index = {}
    placement2type = {}
    type2placement = defaultdict(list)
    cover2index = defaultdict(list)  # cell covered by placement-index

    index_gen = itertools.count()
    for sg in range(1, SUBGRIDS+1):  # sg = subgrid size
        for pos_x in range(N):
            for pos_y in range(M):
                if instance[pos_x, pos_y] == 0:  # free
                    feasible, covering = get_square_covered_indices(instance, pos_x, pos_y, sg)
                    if feasible:
                        new_index = next(index_gen)
                        placements[(sg, pos_x, pos_y)] = covering
                        index2placement[new_index] = (sg, pos_x, pos_y)
                        placement2index[(sg, pos_x, pos_y)] = new_index
                        placement2type[new_index] = sg
                        type2placement[sg].append(new_index)
                        cover2index[(pos_x, pos_y)].append(new_index)

    return placements, index2placement, placement2index, placement2type, type2placement, cover2index

def calculate_collisions(placements, index2placement):
    """ Calculate collisions between tile-placements (position + tile-selection)
        -> only upper triangle is used: a < b! """
    n_p = len(placements)
    coll_mat = np.zeros((n_p, n_p), dtype=bool)  # only upper triangle is used

    for pA in range(n_p):
        for pB in range(n_p):
            if pA < pB:
                covered_A = placements[index2placement[pA]]
                covered_B = placements[index2placement[pB]]
                if len(set(covered_A).intersection(set(covered_B))) > 0:
                    coll_mat[pA, pB] = True

    return coll_mat

""" PREPROCESSING """
placements, index2placement, placement2index, placement2type, type2placement, cover2index = preprocessing(instance, SUBGRIDS)
N_P = len(placements)
coll_mat = calculate_collisions(placements, index2placement)

""" MIP-MODEL """
prob = LpProblem("GridFill", LpMaximize)

# Variables
X = np.empty(N_P, dtype=object)
for x in range(N_P):
    X[x] = LpVariable('x'+str(x), 0, 1, cat='Binary')

# Objective
placement_scores = [SUGBRID_SCORES[index2placement[p][0]] for p in range(N_P)]
prob += lpDot(placement_scores, X), "Score"

# Constraints
# C1: Forbid collisions of placements
for a in range(N_P):
    for b in range(N_P):
        if a < b:  # symmetry-reduction
            if coll_mat[a, b]:
                prob += X[a] + X[b] <= 1  # not both!

""" SOLVE """
print('solve')
#prob.solve(GUROBI())  # much faster commercial solver; if available
prob.solve(PULP_CBC_CMD(msg=1, presolve=True, cuts=True))
print("Status:", LpStatus[prob.status])

""" INTERPRET AND COMPLETE SOLUTION """
solution = np.zeros((N, M), dtype=int)
for x in range(N_P):
    if X[x].value() > 0.99:
        sg, pos_x, pos_y = index2placement[x]
        _, positions = get_square_covered_indices(instance, pos_x, pos_y, sg)
        for pos in positions:
            solution[pos[0], pos[1]] = sg

fill_with_ones = np.logical_and((solution == 0), (instance == 0))
solution[fill_with_ones] = 1

""" VISUALIZE """
plot_compare(instance, solution, SUBGRIDS)

假设/算法性质

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆