Решаем судоку с помощью Алгоритма X

В этой статье рассмотрим «Алгоритм X» Кнута и его применение для решения судоку. Прелесть алгоритма в том, что судоку при этом решается быстро без программирования каких-то продвинутых техник решения.

Началось всё, собственно, с задачки из Project Euler, где, чтобы получить ответ, нужно решить 50 судоку. И вроде ответ на неё получил, написав программку для решения довольно тупым перебором, но как-то осталась неудовлетворённость скоростью решения. Посмотрев, как решают судоку нормальные люди, я обнаружил, что сейчас для этого используется Алгоритм X, придуманный тем самым Дональдом Кнутом.

Алгоритм X

Этот алгоритм решает задачу точного покрытия множества. Или, если хотите, собирает «дискретный паззл», имея в наличии информацию о форме доступных кусочков. Более формально:

• Есть множество S
• Есть набор подмножеств Y этого множества
• Задача состоит в том, чтобы выбрать из Y такие элементы Y*, что каждый элемент из S содержится только в одном из множеств, входящих в Y*
• То есть объединение всех множеств в Y* и составляет (покрывает) множество S, и при этом в Y* нет попарно пересекающихся множеств

Например, рассмотрим множества
S = {1, 2, 3, 4, 5} и
Y = { A = {1, 2},
         B = {2, 3},
         C = {1, 5},
         D = {1, 4},
         E = {5} }
Множества B, D и E формируют точное покрытие множества S.

Для алгоритма X Кнута множество Y представляется в виде двоичной матрицы, где строки соответствуют элементам Y, и A_i,j = 1, если S_j находится в Y_i. Т.е. для примера выше:

Алгоритм поиска точного покрытия следующий:

• Входные данные: множества S и Y; стэк Stack множеств, потенциально входящих в покрытие (может изначально быть пустой или уже иметь какие-то элементы)
  
  1. Если множество S пустое — на стэке лежит искомое покрытие.
  2. Если множество Y пустое — покрытие не найдено.
  3. Ищем в множестве S элемент s, входящий в минимальное число множеств из Y.
  4. Выбираем из Y все строчки X, содержащие s.
  5. Для каждого множества X повторяем 6-9.
  6. Добавляем X в стэк Stack как потенциальный элемент покрытия.
  7. Формируем множества S’ и Y’: S’ — это S, из которого удалены все элементы, содержащиеся в X, Y’ — множества из Y, не пересекающиеся с X.
  8. Вызываем алгоритм X для S’, Y’ и Stack.
  9. Если на шаге 7 получено, что покрытие невозможно — снимаем с вершины Stackа элемент и переходим к следующему X. Если решение найдено — возвращаем его.
  10. Если ни для какого X решения нет — для этих входных данных задача не решается.

В общем, ничего особо сложного. По существу — обычный поиск в глубину. Заметим, кстати, что если изначально задать стэк непустым, то задачу можно сформулировать как «найти точное покрытие, в которое входят элементы, уже лежащие на стэке».

Тонкость в том, что на практике этот алгоритм применяется для задач, где множества в Y — «маленькие», т.е. матрица весьма разреженная, из-за чего, например, поиск пересечений между столбцами при стандартном хранении в виде матрицы занимает непозволительно много времени.
Поэтому Кнут дополняет этот алгоритм механизмом «пляшущих ссылок». Матрица представляется в виде двумерного двусвязного списка: для каждой строки в списке хранятся только номера столбцов, где в этой строке содержатся единицы. Также в списке хранятся ссылки на следующий и предыдущий элемент в строке и столбце. Такая организация позволяет удалять из разреженной матрицы столбцы и строки за время O(1) по сравнению с O(m * n) при хранении в двумерном массиве.

Интересно, что Ali Assaf предлагает альтернативу механизму пляшущих ссылок с использованием ассоциативных списков, что позволяет на высокоуровневых языках реализовывать алгоритм X буквально в несколько десятков строчек.
Идея в том, чтобы хранить как столбцы, так и строки матрицы в ассоциативных списках. В столбцах храним индексы строк, на пересечении с которыми находятся ненулевые элементы, в строках — соответственно, индексы столбцов. Причём в строках будем индексы хранить упорядоченно, в массиве — заметим, что в алгоритме Кнута модифицировать строки, по существу, не требуется, поэтому оптимизация под быстрое удаление элемента из строки не нужна. А вот столбцы будут задаваться в виде множеств, т.к. при удалении строки из матрицы нужно удалить её идентификатор из всех столбцов (и при удалении его из всех столбцов — строка исчезает «сама собой»).

Рассмотрим реализацию алгоритма на Julia.
Матрица из примера будет выглядеть теперь так:

Y = Dict(     'A' => [1, 2],     'B' => [2, 3],     'C' => [1, 5],     'D' => [1, 4],     'E' => [5] )  S = Dict(     1 => Set(['A', 'C', 'D']),     2 => Set(['A', 'B']),     3 => Set(['B']),     4 => Set(['D']),     5 => Set(['C', 'E']) )

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

function extract_intersects!(rows, columns, base_row)     buf = []     for elt in rows[base_row]         # вынимаем текущий столбец из таблицы в буфер         push!(buf, pop!(columns, elt))         # удаляем все пересекающиеся строки из всех оставшихся столбцов         for intersecting_row in buf[end]             for other_elt in rows[intersecting_row]                 if other_elt != elt                     delete!(columns[other_elt], intersecting_row)                 end             end         end     end     return buf end  function restore_intersects!(rows, columns, base_row, buf)     # удаляли столбцы от первого пересечения с base_row к последнему, восстанавливать надо в обратном порядке     for elt in Iterators.reverse(rows[base_row])         columns[elt] = pop!(buf)         for added_row in columns[elt]             for col in rows[added_row]                 push!(columns[col], added_row)             end         end     end end

Чтобы эти две функции работали как надо, как раз и требовалось упорядоченное хранение элементов в строках матрицы. В функции extract_intersects!() на каждой итерации внешнего цикла из матрицы убираются те строки, которые пересекаются с base_row, но не содержат элементов, просмотренных на предыдущих итерациях. Это гарантирует, что, когда мы в restore_intersects!() вставляем столбцы в обратном порядке, в самом внутреннем цикле на момент вызова push!(columns[col], added_row) столбец columns[col] в матрицу уже будет возвращён, и все удалённые в extract_intersects!() элементы из столбцов будут возвращены на прежнее место.

Теперь сам алгоритм X:

function algorithm_x(rows, columns, cover = [])     if isempty(columns)         return cover     else         # ищем столбец с минимальным числом элементов         m, c = findmin(Dict(k => length(v) for (k, v) in columns))         for subset in collect(columns[c])             push!(cover, subset)             # удаляем пересекающиеся подмножества и             # содержащиеся в subset элементы             buf_cols = extract_intersects!(rows, columns, subset)             s = algorithm_x(rows, columns, cover)             # если нашлось непустое решение - готово, выходим             s == nothing || return s             restore_intersects!(rows, columns, subset, buf_cols)             pop!(cover)         end         # сюда дойдём либо если в columns[c] пусто,         # либо когда рекурсивный поиск не нашёл решения         return nothing     end end

Судоку

Алгоритм есть, дело за малым — представить судоку как задачу поиска точного покрытия.
Сформулируем требования, которым должно удовлетворять решённое судоку:

1. В каждой клетке стоит цифра от 1 до 9 (или до n², если решаются квадраты другого размера).
2. В каждой строке каждое число встречается по разу.
3. В каждом столбце каждое число встречается по разу.
4. В каждом квадранте каждое число встречается по разу.

Каждое из этих требований должно выполняться ровно по 1 разу, т.е. они и формируют множество, которое надо покрыть. В нём ровно 4n² элементов (столбцов в матрице).
Подмножества, которые рассматриваем, формируются подстановкой конкретного числа в конкретную клетку. Например, число 9 на пересечении 1 строки и 4 столбца «накрывает» подмножество «в клетке (1,4) есть число, в 1 строке есть число 9, в 4 столбце есть число 9, во 2 квадранте есть число 9» (подразумевая обычное судоку 9×9).
После этого алгоритм решения пишется тривиально.

# судоку задаётся матрицей 9×9, на месте неизвестных чисел нули # идентификаторы строк - кортежи вида (row, col, num) # идентификаторы столбцов: # (0, row, col) - на пересечении row и col стоит число # (1, row, num) - в строке row есть число num # (2, col, num) - в столбце col есть число num # (3, q, num) - в квадранте q есть число num function xsudoku(puzzle::AbstractMatrix{Int})     rows = Dict()     cols = Dict()     # заполняем строки     for row in 1:9, col in 1:9, num in 1:9         r = []         quad = ((row-1)÷3)*3 + (col-1)÷3 + 1         push!(r, (0, row, col), (1, row, num), (2, col, num), (3, quad, num))         rows[(row, col, num)] = r     end     # заполняем столбцы     for type in 0:3, n1 in 1:9, n2 in 1:9         cols[(type, n1, n2)] = Set()     end     for (rk, rv) in rows         for v in rv             push!(cols[v], rk)         end     end      # s - заготовка для ответа     # для начала, туда надо внести те цифры, которые уже заполнены     s = []     for i in 1:9, j in 1:9         if puzzle[i, j] > 0             elt = (i, j, puzzle[i,j])             push!(s, elt)             # добавив клетку в решение, удаляем из матрицы все несовместимые элементы             extract_intersects!(rows, cols, elt)         end     end      # всё, что осталось - найти покрытие     success = algorithm_x(rows, cols, s)     success != nothing || return nothing     # ответ выдадим в виде матрицы     ret = similar(puzzle)     for (i, j, num) in s         ret[i,j] = num     end     return ret end

Проверим на каком-нибудь примере:

 julia> @time xsudoku([0 0 0 0 0 0 4 0 0;        3 0 6 0 0 0 0 0 0;        0 0 0 1 9 6 0 3 0;        0 7 0 0 0 0 0 1 0;        8 0 0 2 5 0 0 9 0;        0 4 0 0 0 0 8 0 0;        0 6 0 4 0 9 0 0 8;        0 0 5 0 0 0 0 2 0;        0 0 0 5 0 0 0 0 7])   0.006326 seconds (54.91 k allocations: 3.321 MiB) 9×9 Array{Int64,2}:  1  5  7  8  3  2  4  6  9  3  9  6  7  4  5  2  8  1  2  8  4  1  9  6  7  3  5  6  7  2  9  8  4  5  1  3  8  3  1  2  5  7  6  9  4  5  4  9  6  1  3  8  7  2  7  6  3  4  2  9  1  5  8  4  1  5  3  7  8  9  2  6  9  2  8  5  6  1  3  4  7

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