Megaron-Cube — это головоломка, в которой шесть разных блоков должны быть объединены в куб 3x3x3. Я пытался решить пару раз, но так и не понял. Я знаю, что использовать программу — обман…

… ну, это можно рассматривать как еще один вызов, и, возможно, я лучше справлюсь с этим.

Что такое Мегарон-Куб?

Megaron-Cube состоит из следующих шести блоков, которые необходимо объединить, чтобы построить куб 3x3x3.

Как показано, все шесть блоков разные. Для нашей визуализации мы будем использовать эти разные цвета. Обычно блоки изготавливаются из дерева.

Как нарисовать эти блоки с помощью Python?

Я не хочу решать эту головоломку, не визуализируя решение. Проще проверить правильность шагов решения, если их можно визуализировать.

В первую очередь нужны эти пакеты:

import numpy as np 
import matplotlib.pyplot as plt 
from mpl_toolkits.mplot3d import axes3d 
from itertools import product, combinations 
import math3d as m3d 
import math

Вероятно, вам нужно установить math3d с помощью pip install math3d. Он используется для вращения наших блоков.

Мы начнем с рисования одного куба в определенной позиции в сетке 3x3x3. Поэтому определяется трехмерная фигура matplotlib.

fig = plt.figure() 
fig.suptitle('Megaron-Cube', fontsize=14, fontweight='bold')
ax = fig.gca(projection='3d')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')

Чтобы нарисовать один куб, мы будем использовать

drawCube(ax,0,0,1,"blue")
ax.set_xlim3d(0, 3)
ax.set_ylim3d(0, 3)
ax.set_zlim3d(0, 3)
plt.show()

Здесь мы вызываем функцию drawCube. Синий куб должен быть нарисован в позиции [0,0,1]. Ряды set_ определяют, что мы видим всю сетку 3x3x3. Теперь необходимо определить функцию drawCube.

def drawCube(ax,x,y,z,color):
    r = [0,1]
    X, Y = np.meshgrid(r, r)
    ax.plot_surface(x+X,y+Y,z+1, alpha=1, color=color) # top
    ax.plot_surface(x+X,y+Y,z+0, alpha=1, color=color) # bottom
    ax.plot_surface(x+X,y+0,z+Y, alpha=1, color=color) # front
    ax.plot_surface(x+X,y+1,z+Y, alpha=1, color=color) # back
    ax.plot_surface(x+0,y+X,z+Y, alpha=1, color=color) # left
    ax.plot_surface(x+1,y+X,z+Y, alpha=1, color=color) # right

Функция meshgrid возвращает матрицы координат X и Y для вектора координат [0,1]. Вывод:

X: 
[[0 1]
 [0 1]]
Y: 
[[0 0]
 [1 1]]

Единую поверхность для нашего куба можно нарисовать с помощью функции plot_surface. Этой функции нужны эти матрицы координат или одно значение (здесь z+1).

ax.plot_surface(x+X,y+Y,z+1, alpha=1, color=color) # top

Весь куб для нашего вызова функции drawCube(ax,0,0,1,"blue") будет выглядеть так:

Давайте определим наши шесть блоков, задав положение этих кубов. На этом этапе мы также определяем цвета блоков.

color = {}
color["a"] = "red"
color["b"] = "green"
color["c"] = "blue"
color["d"] = "orange"
color["e"] = "yellow"
color["f"] = "purple"
"""
    define all blocks
"""
blocks = {}
blocks["a"] = np.array([[0,0,0],[1,0,0],[2,0,0],[1,0,1]])
blocks["b"] = np.array([[0,0,0],[1,0,0],[2,0,0],[1,0,1],[0,1,0]])
blocks["c"] = np.array([[0,0,0],[1,0,0],[1,0,1],[1,1,0],[2,1,0]])
blocks["d"] = np.array([[0,0,0],[0,0,1],[0,1,0],[1,1,0],[1,2,0]])
blocks["e"] = np.array([[0,0,0],[1,0,0],[0,1,0],[0,0,1]])
blocks["f"] = np.array([[0,0,0],[1,0,0],[0,1,0],[1,0,1]])

В первом блоке есть кубики в четырех позициях [[0,0,0],[1,0,0],[2,0,0],[1,0,1]].

Мы хотим рисовать блоки в том же порядке, поэтому мы определяем список:

list_blocks = ["a","b","c","d","e","f"]

Цель состоит в том, чтобы нарисовать каждый блок во всех возможных позициях. Эти две линии дают нам все позиции внутри куба 3x3x3.

r3 = range(3) poss = list(product(r3,r3,r3))

Чтобы получить все вращения, используются следующие векторы вращения:

deg = 90
rad = deg*math.pi/180
rz = m3d.Orientation.new_axis_angle([0,0,1], rad)
ry = m3d.Orientation.new_axis_angle([0,1,0], rad)
rx = m3d.Orientation.new_axis_angle([1,0,0], rad)

Как видно из названия, rz определяет вращение вокруг оси Z.

Чтобы получить все вращения, мы снова используем декартово произведение.

protations = [False,rz,ry,rx]
protations = list(product(protations,protations,protations,protations))

False означает, что мы не вращаемся. Список некоторых protations (возможных поворотов), где некоторые из них могут привести к одному и тому же объекту.

[(False, False, False, False),
 (False, rx, ry, False),
 (ry, ry, ry, False), ... ]

Последняя строка описывает комбинацию трех поворотов вокруг оси Y. Каждое значение False будет игнорироваться, поэтому (False, rx, ry, False) совпадает с (rx, False, ry, False).

Этот массив протаций должен быть преобразован в объект вращения.

rotations = [False]
for r in protations:
    cr = m3d.Orientation.new_axis_angle([1,0,0], 0)
    for t in r:
        if t is not False:
            cr *= t
    if cr not in rotations:
        rotations.append(cr)

Здесь инициализируется массив rotations, в котором хранится первое вращение (без вращения). Для всех кортежей в protations отдельные повороты объединяются с использованием cr *= t, если t не равно False. Текущее вращение cr сначала определяется как вращение вокруг оси x с углом 0 рад (без вращения). В конце этого блока мы проверяем, является ли вращение уже частью rotations, чтобы удалить дубликаты, подобные упомянутому выше.

Цель состоит в том, чтобы показать все возможные позиции блоков для каждого блока.

for cblock in list_blocks:
    for rotation in rotations:
        for pos in poss:
            boolCheck, cubes = checkCubes(cblock,blocks[cblock],pos,rotation)
            if boolCheck:    
                ax.clear()
                drawCubes(ax,cubes,color[cblock])
ax.set_xlim3d(0, 3)
                ax.set_ylim3d(0, 3)
                ax.set_zlim3d(0, 3)
                plt.pause(0.5)
while True:
    plt.pause(50)

Для каждого определенного блока в list_blocks мы перебираем все повороты и позиции. Затем мы проверяем, являются ли кубики разумными. Здесь разумный означает, что блок находится внутри сетки 3x3x3 и ранее не рисовался. В противном случае некоторые блоки будут перерисованы несколько раз, потому что они симметричны. Если куб является разумным, мы нарисуем его и сделаем паузу на 0,5 секунды, прежде чем будет нарисован следующий разумный блок. В конце мы делаем паузу, чтобы показать последний блок, пока...

Давайте сначала посмотрим на простую функцию drawCubes.

def drawCubes(ax,cubes,color):
    for cube in cubes: 
        drawCube(ax,cube[0],cube[1],cube[2],color)

Каждый блок определяется как список кубов, которые можно нарисовать как отдельные кубы.

Теперь мы проверим, является ли блок разумным. Используются параметры name, block, pos и rot. name определяет имя для текущего блока, который является одним из a-f. block — это список кубов для блока. pos определяет позицию первого куба в блоке, а rot — параметр поворота.

def checkCubes(name,block,pos,rot):    
    global obj_pos
    global obj_fields
    
    cubes = []
    for cube in block: 
        cube = getCube(cube,pos,rot)
        cubes.append(cube)
        for val in cube:
            if val < 0 or val > 2:
                return False, cubes
    cubes.sort()
    if cubes in obj_pos[name]:
        return False, cubes
    obj_pos[name].append(cubes)  
    field = np.zeros((3,3,3))
    for cube in cubes:
        field[cube] += 1
    obj_fields[name].append(field)    
    return True, cubes

Эта функция использует две глобальные переменные, которые определены в верхней части программы:

from collections import defaultdict
obj_pos = defaultdict(list)
obj_fields = defaultdict(list)

Они используются для проверки того, нарисован ли уже блок.

Все кубы для блоков вычисляются, и проверяется, находится ли одна из этих точек куба за пределами куба 3x3x3.

cubes = []
for cube in block: 
    cube = getCube(cube,pos,rot)
    cubes.append(cube)
    for val in cube:
        if val < 0 or val > 2:
            return False, cubes

Он использует функцию getCube, которая вращает куб и переводит его.

def getCube(cube,pos,rot):
    if rot is not False:
        cube_vec = m3d.Vector(*cube)
        cube = list(rot*cube_vec)
        cube = [ int(round(x)) for x in cube ]
    return tuple([x + y for x, y in zip(cube, pos)])

После этого список кубов сортируется. Этот массив:

[(0, 0, 0), (1, 0, 0), (2, 0, 0), (1, 0, 1)]

сортируется по:

[(0, 0, 0), (1, 0, 0), (1, 0, 1), (2, 0, 0)]

Это определяет тот же блок. Повороты, которые снова получают тот же блок, как и поворот на 360 градусов, имеют такое же отсортированное представление. Используется в следующих рядах.

cubes.sort()
if cubes in obj_pos[name]:
    return False, cubes
obj_pos[name].append(cubes)

Позиции сохраняются, чтобы проверить, нарисован ли уже блок, и чтобы быстрее рисовать кубы для следующих шагов. Затем сохраняется представление куба 3x3x3 (field), где отображается только текущий блок.

field = np.zeros((3,3,3))
for cube in cubes:
    field[cube] = 1
obj_fields[name].append(field)    
return True, cubes

Представительство для

сопоставляется с:

[[[ 1.  0.  0.]
  [ 0.  0.  0.]
  [ 0.  0.  0.]]
[[ 1.  1.  0.]
  [ 0.  0.  0.]
  [ 0.  0.  0.]]
[[ 1.  0.  0.]
  [ 0.  0.  0.]
  [ 0.  0.  0.]]]

Используя эти функции, мы можем рисовать все возможные позиции каждого блока. Следующим шагом является объединение этих блоков. Важно ускорить сохранение всех возможных позиций для каждого блока, что делается в checkCubes. Поэтому мы вызываем эти функции для всех блоков, позиций и вращений.

for cblock in list_blocks:
    for rotation in rotations:
        for pos in poss:
            checkCubes(cblock,blocks[cblock],pos,rotation)

Теперь мы можем получить количество возможных позиций для каждого блока.

pos_lengths = [range(len(obj_pos[x])) for x in list_blocks] print(pos_lengths)

Сохраняется не сама длина, а range, что упрощает дальнейшие шаги.

Получим все комбинации для первых двух блоков.

possibleCombinations = list(product(pos_lengths[0],pos_lengths[1]))

Эти possibleCombinations могут быть невозможны, но я не знаю, как их назвать ;)

Идея следующей части нашей программы состоит в том, чтобы получить все реальные возможности двух блоков. Если комбинация возможна, мы будем добавлять к ней следующий блок, пока не будут возможны все шесть блоков. Тогда мы нашли решение.

Мы добавим возможности, поэтому используем цикл while, который останавливается, если не осталось возможных комбинаций, которые следует проверить.

p = 0 
while p < len(possibleCombinations):

Затем мы проверяем, возможно ли текущее положение блоков:

combination = possibleCombinations[p]
checkBool = checkListOfCubes(combination)

checkBool равно True, если текущая комбинация действительно возможна, и False в противном случае. Комбинированный кортеж выглядит примерно так: (1,4) или (1,4,8). Первое целое число определяет положение блока «a», второе — положение блока «b» и так далее.

Блоки не могут объединяться, если позиция используется более одного раза. Мы определяем пустой field, который представляет куб 3x3x3, и добавляем блоки, определенные в combination. Если позиция используется более одного раза, возвращается False, а в противном случае — True.

def checkListOfCubes(combination):    
    global list_blocks
    global obj_fields
    
    field = np.zeros((3,3,3))
    for p_idx in range(len(combination)):
        p = combination[p_idx]
        field += obj_fields[list_blocks[p_idx]][p]
if np.max(field) > 1:
        return False
    return True

Вернемся к нашему циклу while:

p = 0
while p < len(possibleCombinations):
    combination = possibleCombinations[p]
    checkBool = checkListOfCubes(combination)

Теперь мы хотим добавить следующий блок, если текущая комбинация возможна.

if checkBool:
    block_no = len(combination)
    if block_no < len(list_blocks):
        for q in pos_lengths[block_no]:
            lcombination = list(combination)
            lcombination.append(q)
            possibleCombinations.append(tuple(lcombination))
    if len(combination) == len(list_blocks):
        solutions.append(combination)

block_no определяет номер следующего блока, который мы хотим добавить, если он существует. Затем мы добавляем все позиции, которые возможны для этого нового блока, в список possibleCombinations. Если мы использовали все шесть блоков, мы нашли решение и добавили его к решениям.

После цикла while мы сохраняем решения и позиции блоков с помощью pickle.

save_obj("solutions",solutions)
save_obj("obj_pos",obj_pos)

Функция определяется как:

def save_obj(filename, obj):
    with open(filename+'.pickle', 'wb') as handle:
        pickle.dump(obj, handle)

Визуализация решения

Давайте визуализируем решение в другой программе под названием drawSolution.py. Он использует некоторые переменные и функции, описанные ранее.

В первую очередь загружаем решения и позиции блоков.

def load_obj(filename):
    with open(filename+'.pickle', 'rb') as handle:
        return pickle.load(handle)
solutions = load_obj("solutions")
obj_pos = load_obj("obj_pos")

Теперь мы можем нарисовать решения:

for combination in solutions: 
    ax.clear()
    list_of_cubes = []
    list_of_colors = [] 
    for block_idx in range(len(combination)):
        block = list_blocks[block_idx]
        cubes = obj_pos[block][combination[block_idx]]
        list_of_cubes.append(cubes)
        list_of_colors.append(color[block])
for cube_idx in range(len(list_of_cubes)):
        cubes = list_of_cubes[cube_idx]
        c = list_of_colors[cube_idx]
        drawCubes(ax,cubes,c)
ax.set_xlim3d(0, 3)
    ax.set_ylim3d(0, 3)
    ax.set_zlim3d(0, 3)
    plt.pause(5)

Сначала все блоки добавляются к list_of_cubes, а цвета к list_of_colors. После этого они рисуются с помощью drawCubes.

Хорошо, но сколько существует решений и как быстро мы можем их решить?

Есть 24 решения, и их поиск занимает около минуты.

Если бы вы попытались решить эту головоломку, то могли бы подумать… Ух ты, 24 решения, и я не нашел ни одного из них? Успокойтесь :) Это не так уж и плохо. Существует только одно реальное решение, а остальные являются зеркальными.

Улучшения

Если мы сможем удалить эти зеркальные решения, это будет быстрее, потому что нужно будет проверять меньше комбинаций. Чтобы удалить их, можно использовать первый блок, чтобы уменьшить пространство поиска.

На данный момент есть 72 возможных позиции для первого блока.

Давайте посмотрим на некоторые из них.

Первые два изображения будут одинаковыми, если мы повернем куб 3x3x3. То же самое верно и для двух других визуализаций.

Всего реальных позиций для первого блока всего четыре.

Мы использовали следующее, чтобы получить все позиции для каждого блока.

for cblock in list_blocks:
    for rotation in rotations:
        for pos in poss:
            checkCubes(cblock,blocks[cblock],pos,rotation)

Теперь нам нужны только эти четыре позиции первого блока, поэтому вместо этого используется это:

cblock = list_blocks[0]
for pos in poss:
    if pos[1] != 2:
        checkCubes(cblock,blocks[cblock],pos,False)
        
for cblock in list_blocks[1:]:
    for rotation in rotations:
        for pos in poss:
            checkCubes(cblock,blocks[cblock],pos,rotation)

Теперь на моей машине требуется всего 2,5 секунды, и генерируется одно решение.

Возможно как минимум еще одно улучшение. Некоторые комбинации блоков возможны, но видно, что они не могут быть частью общего решения, подобного следующему:

Это неразумно, потому что есть единственный блок (внизу справа), который не может быть заполнен. В каждом блоке есть как минимум четыре кубика. Следовательно, мы можем удалить все комбинации, в которых наименьшее отверстие меньше четырех.

Эта проверка будет выполнена в checkListOfCubes в конце вместо return True.

return checkZeroNeighborField(field)

Идея состоит в том, чтобы проверить, насколько большим может быть самый большой блок, который заполняет пустой куб.

Нам нужно несколько функций, чтобы получить это представление. Прежде всего, мы получаем всех соседей для каждой позиции и сохраняем их в глобальном объекте.

neighborsObj = {}
def getNeighborsPos(pos):
    global neighborsObj
    
    neighbors = []
    for t in [[0,0,1],[0,1,0],[1,0,0],[0,0,-1],[0,-1,0],[-1,0,0]]:
        test = [pos[0]+t[0],pos[1] + t[1],pos[2]+t[2]]
        if max(test) <= 2 and min(test) >= 0:
            neighbors.append(test)
    neighborsObj[pos] = neighbors
    return neighbors

Каждый куб имеет не более шести возможных соседей. По крайней мере один из них не находится внутри куба 3x3x3, если pos не является средней позицией (1,1,1).

Они сохраняются в neighborsObj где-то в начале нашей программы с помощью:

for pos in poss:
    getNeighborsPos(pos)

Нас интересуют все пустые позиции в нашем поле, поэтому релевантны только пустые соседи.

def getZeroNeighbors(field,pos):
    neighbors = neighborsObj[tuple(pos)]
    zero_neighbors = []
    for neighbor in neighbors:
        if field[neighbor[0],neighbor[1],neighbor[2]] == 0:
            zero_neighbors.append(neighbor)
    return zero_neighbors

Следующий шаг — получить самый длинный путь, используя только пустых соседей из определенной позиции. Этот путь вычисляется с помощью следующей функции:

def getLongestZeroPath(field, pos):
    path = []
    pathObj = {}
    neighbors2check = getZeroNeighbors(field,pos)
    for neighbor in neighbors2check:
        path.append(neighbor)
        pathObj[tuple(neighbor)] = 1
        
    c = 0
    while len(neighbors2check) > 0:
        c+=1
        neighbor = neighbors2check.pop(0)
        for n in getZeroNeighbors(field,neighbor):
            if tuple(n) not in pathObj:
                path.append(n)
                pathObj[tuple(n)] = 1
                neighbors2check.append(n)
    return path

Мы определяем путь, используя массив и более быстрый объект, чтобы проверить, находится ли позиция уже внутри пути. Сначала к пути добавляются все непосредственные соседи.

Затем для каждого из этих соседей добавляются пустые соседи, которые не являются частью пути, и так далее.

На последнем шаге, чтобы получить представление самых больших блоков для каждой позиции, используется следующая функция.

def checkZeroNeighborField(field):
    r3 = range(3)
    neighborsMap = np.zeros((3,3,3))
    for pos in product(r3,r3,r3):
        if field[pos[0],pos[1],pos[2]] == 0:
            if neighborsMap[pos[0],pos[1],pos[2]] == 0:
                zeroPath = getLongestZeroPath(field,pos)
                len_zeroPath = len(zeroPath)
                if len_zeroPath < 4:
                    return False
                neighborsMap[pos[0],pos[1],pos[2]] = len_zeroPath
                for point in zeroPath:
                    neighborsMap[point[0],point[1],point[2]] = len_zeroPath            
        else:
            neighborsMap[pos[0],pos[1],pos[2]] = -1
    return True

Для каждой позиции мы проверяем, является ли она пустой. Если не пусто, мы устанавливаем значение в нашем neighborsMap на -1. В противном случае мы вычисляем самый длинный нулевой путь и возвращаем False, если путь короче четырех, потому что наши наименьшие блоки состоят из четырех кубов. Нет необходимости вычислять getLongestZeroPath для каждой позиции. Функция разумна только в том случае, если текущая позиция не является частью пути ранее.

Это улучшение решает головоломку примерно за 0,6 секунды. Что составляет одну сотую от нашего первого времени выполнения.

Надеюсь, вам понравился этот необычный пост в блоге.

Если вы хотите решить головоломку самостоятельно, вам следует сначала купить этот амазонский мегарон-куб.

Вы можете скачать весь код в моем репозитории OpenSourcES GitHub.

Первоначально опубликовано на opensourc.es.