В прошлых статьях я разбирал линейное и целочисленное программирование: PuLP, OR-Tools, pyomo, задачи о назначениях, коммивояжёра, раскрой, генерацию столбцов. Здесь рассмотрим немного иной подход: выполнимость булевой формулы, она же SAT. С одной стороны, это фундамент теории NP-полноты. С другой, вполне рабочий инструмент, который на чисто комбинаторных задачах нередко обходит классические MIP-солверы.

Дальше по плану: что такое SAT, как устроены алгоритмы, позволяющие солверам тянуть сотни тысяч переменных, как с этим работать из PySAT и какие солверы он даёт. В финале решим прикладную оптимизационную задачу через MaxSAT.

Кому пригодится: тем, кто занимается математической оптимизацией или ML-инженерам, как ещё один инструмент в своем наборе. Если вы отвечаете за процессы, то в примерах, возможно, узнаете знакомые сюжеты: ручные расписания, подбор совместимых комплектаций, дележ дефицитного ресурса. Всё это неплохо ложится на SAT.

1. Задача SAT: о чём вообще речь

Формулировка задачи о выполнимости (Boolean satisfiability problem, SAT) короткая:

Дана формула логики высказываний. Существует ли набор значений переменных (True/False), при котором она истинна?

Если такой набор есть, формула выполнима (SAT), и солвер обязан вернуть конкретные значения. Если набора нет ни одного, формула невыполнима (UNSAT).

С SAT связан исторический сюжет: в теореме Кука и Левина (1971) именно она стала первой задачей с доказанной NP-полнотой. Любую задачу из NP можно за полиномиальное время свести к SAT. Отсюда и универсальность, и неприятная оборотная сторона в виде экспоненты в худшем случае.

Конъюнктивная нормальная форма (КНФ)

Подавляющее большинство солверов принимают формулу в конъюнктивной нормальной форме. КНФ это конъюнкция клауз (логическое «И»), где каждая клауза есть дизъюнкция литералов (логическое «ИЛИ»), а литерал это переменная либо её отрицание.

В программистской записи переменные нумеруются целыми, а отрицание задаётся минусом. Та же формула:

[[1, -2], [-1, 3], [-2, -3]]

Клауза истинна, если истинен хотя бы один её литерал. Формула истинна, если истинны все клаузы. Любую булеву формулу можно привести к КНФ (например, преобразованием Цейтина), так что требование КНФ ничего не отнимает по части выразимости.

Для нас, «оптимизаторов», ценность вот в чём. Дискретные ограничения естественно записываются клаузами, а булевы переменные это те же 0/1, что и в целочисленном программировании. По сути перед нами язык решений «да/нет»: ставить сотрудника на смену или нет, совместима ли опция в заказе, разрешён ли слот. Каждое такое правило задаёт клаузу, а вопрос «складывается ли план разом при всех правилах» и есть SAT.

2. Основные алгоритмы

2.1. DPLL

Алгоритм Дэвиса, Патнема, Логемана и Лавленда (DPLL, 1962) это перебор с возвратом (backtracking), усиленный двумя приёмами, которые сильно обрезают дерево поиска:

• Распространение единичных клауз (unit propagation). Если в клаузе остался один неопределённый литерал, а все прочие уже ложны, этот литерал вынужденно истинен. Присваиваем его и каскадом смотрим, не стали ли единичными другие клаузы.

• Удаление чистых литералов (pure literal elimination). Если переменная встречается во всей формуле в одной полярности (везде без отрицания либо везде с ним), её сразу фиксируют так, чтобы удовлетворить все содержащие её клаузы.

Когда оба приёма исчерпаны, а ответа нет, алгоритм ветвится: берёт переменную, пробует значение, спускается рекурсивно, на конфликте откатывается и пробует второе.

Псевдокод:

function DPLL(formula):
    formula = unit_propagate(formula)
    formula = pure_literal_assign(formula)
    if formula пуста:            
      return SAT
    if есть пустая клауза:       
      return UNSAT   # конфликт
    x = выбрать переменную
    return DPLL(formula ∧ x) or DPLL(formula ∧ ¬x)

2.2. CDCL

В промышленных солверах работает не голый DPLL, а CDCL (Conflict-Driven Clause Learning). Отличается он отношением к конфликтам: конфликт здесь повод не просто откатиться, а извлечь из него новую клаузу и больше в эту яму не падать.

• Анализ конфликта и обучение клауз (clause learning). По графу импликаций солвер находит причину конфликта и дописывает в формулу «выученную» клаузу, которая запрещает повтор той же комбинации присваиваний.

• Нехронологический откат (backjumping). Откатываемся не на шаг назад, а сразу на тот уровень решений, который и спровоцировал конфликт.

• Two-watched literals. В каждой клаузе отслеживаются лишь два литерала, и за счёт этого unit propagation остаётся дешёвым даже на миллионах клауз.

• VSIDS (Variable State Independent Decaying Sum). Эвристика выбора переменной: приоритет у тех, что часто всплывают в недавних конфликтах.

• Рестарты. Периодически сбрасываем дерево поиска (выученные клаузы сохраняем), чтобы не залипать в неудачной области пространства.

Связка «обучение клауз + watched literals + VSIDS + рестарты» и сделала из SAT практический инструмент: задачи с сотнями тысяч переменных решаются за секунды.

2.3. Локальный поиск

Есть и отдельная ветвь, стохастический локальный поиск (GSAT, WalkSAT): стартуем со случайного назначения и переворачиваем значения переменных, уменьшая число невыполненных клауз. Методы неполные, UNSAT они не доказывают, зато выполнимую формулу иногда находят очень быстро. Для бизнес-задач я по умолчанию беру полный CDCL: предсказуемость тут важнее редких удачных забегов.

3. От разрешимости к оптимизации: MaxSAT

Обычный SAT отвечает только «да/нет». Нам же чаще нужно лучшее решение, а не любое допустимое. Для этого есть MaxSAT.

В базовом виде MaxSAT это «выполнить как можно больше клауз». На практике полезнее взвешенный частичный вариант (Weighted Partial MaxSAT), где клаузы двух сортов:

• жёсткие (hard), их нужно выполнять всегда, это ограничения модели;

• мягкие (soft), у каждой есть вес, и солвер старается набрать максимально «дорогой» выполнимый набор.

Эквивалентно солвер минимизирует суммарный вес нарушенных мягких клауз; эту величину называют cost. Дальше аналогия с MIP читается напрямую: жёсткие клаузы играют роль ограничений, веса мягких задают коэффициенты целевой функции. Хранится такая модель в формате WCNF (взвешенная КНФ).

Cardinality и псевдобулевы ограничения

Неудобный момент при кодировании это ограничения вида «не более k из набора» или «ровно один». В лоб «не более одного» через попарные запреты даёт клауз, и на больших наборах это бьёт по производительности. Поэтому используют экономные кодировки:

• sequential counter и totalizer для cardinality () имеют линейный или размер;

• псевдобулевы ограничения (pseudo-Boolean, ) пригодны для взвешенных сумм.

В PySAT всё это уже есть (pysat.card, pysat.pb), вручную кодировать не придётся.

4. PySAT и солверы

PySAT это Python-обёртка над набором конкурсных SAT-солверов плюс утилиты для кодирования. Внутри доступны, среди прочих, Glucose, MiniSat, CaDiCaL, Lingeling, MapleSAT, MergeSat, MiniCard. Сменить движок можно одним аргументом.

Модули, которые понадобятся дальше:

«Hello, SAT»

Решим формулу из первого раздела:

# Установка
!pip install python-sat[pblib,aiger]

from pysat.solvers import Solver

# (x1 ∨ x2) ∧ (¬x1 ∨ x3) ∧ (¬x2 ∨ ¬x3)
clauses = [[1, 2], [-1, 3], [-2, -3]]

with Solver(name='glucose3', bootstrap_with=clauses) as s:
    if s.solve():
        print('SAT:', s.get_model())   # напр. [1, -2, 3]  ->  x1=True, x2=False, x3=True
    else:
        print('UNSAT')

Менеджер переменных IDPool

В реальной модели адресовать переменные голыми числами неудобно. IDPool сопоставляет произвольному Python-объекту уникальный целочисленный id и умеет ходить обратно:

from pysat.formula import IDPool

vpool = IDPool()
v = vpool.id(('x', 2, 3))   # вернёт целое; повторный вызов вернёт то же число
print(vpool.obj(v))         # -> ('x', 2, 3)

Cardinality-ограничения «из коробки»

«Ровно один из списка» можно закодировать масштабируемо. Для крупных наборов вместо EncType.pairwise (те самые клауз) берут EncType.seqcounter или EncType.totalizer:

from pysat.card import CardEnc, EncType
from pysat.formula import CNF, IDPool

vpool = IDPool()
lits = [vpool.id(('x', c)) for c in range(5)]

cnf = CNF()
cnf.append(lits)                                          # хотя бы один
am1 = CardEnc.atmost(lits=lits, bound=1,
                     vpool=vpool, encoding=EncType.pairwise)
cnf.extend(am1.clauses)                                   # не более одного

5. Практика I: раскраска графа как SAT (разрешимость)

Начнём с decision-задачи, так проще почувствовать, как вообще выглядит кодирование. Возьмём раскраску графа (её я уже упоминал в Типовых моделях). За абстракцией прячутся обычные прикладные истории: распределение радиочастот, составление расписаний и экзаменов (вершины это предметы или смены, ребро это запрет «нельзя в один слот»), рассадка несовместимых процессов по серверам, распределение регистров в компиляторах. Везде один вопрос: влезут ли сущности в ограниченное число слотов без конфликтов.

Можно ли покрасить вершины графа в k цветов так, чтобы концы любого ребра были разного цвета?

Переменные кодируем по схеме one-hot: булева означает «вершина покрашена в цвет ».

Ограничения:

1. У каждой вершины ровно один цвет. Раскладываем на два: «хотя бы один» (клауза-дизъюнкция) и «не более одного» (попарные запреты).

2. Для каждого ребра и цвета концы не делят цвет: .

from pysat.solvers import Glucose3
from pysat.formula import IDPool

def color_graph(edges, n_vertices, k_colors):
    vpool = IDPool()
    x = lambda v, c: vpool.id(('x', v, c))   # вершине v назначен цвет c

    solver = Glucose3()
    for v in range(n_vertices):
        solver.add_clause([x(v, c) for c in range(k_colors)])        # >= 1 цвет
        for c1 in range(k_colors):                                   # <= 1 цвет
            for c2 in range(c1 + 1, k_colors):
                solver.add_clause([-x(v, c1), -x(v, c2)])
    for (u, v) in edges:
        for c in range(k_colors):
            solver.add_clause([-x(u, c), -x(v, c)])                  # разные цвета

    if not solver.solve():
        solver.delete()
        return None

    model = set(solver.get_model())
    coloring = {v: c for v in range(n_vertices)
                     for c in range(k_colors) if x(v, c) in model}
    solver.delete()
    return coloring

# Цикл из 5 вершин C5 (нечётный цикл)
edges = [(0, 1), (1, 2), (2, 3), (3, 4), (4, 0)]

print(color_graph(edges, 5, 2))   # -> None  (UNSAT)
print(color_graph(edges, 5, 3))   # -> {0: 0, 1: 1, 2: 0, 3: 1, 4: 2}  (SAT)

Вывод показателен: нечётный цикл в два цвета не раскрашивается (солвер возвращает UNSAT), а в три уже да. Если перебирать k от меньшего к большему и брать первый SAT, получим хроматическое число. По сути это уже оптимизация поверх серии вызовов SAT.

6. Практика II: распределение каналов как Weighted Partial MaxSAT (оптимизация)

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

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

Это взвешенная раскраска в условиях нехватки цветов. Обычной раскраской (decision) тут не обойтись: допустимого решения может вовсе не быть. Зато на Weighted Partial MaxSAT задача ложится ровно:

• жёсткие клаузы: у каждого передатчика ровно один канал (структура, нарушать нельзя);

• мягкие клаузы: для каждой мешающей пары и канала берём с весом этой пары. Если пара всё же села на общий канал, ровно одна такая клауза нарушится и добавит её вес в cost.

Пару слов про RC2 (Relaxable Cardinality Constraints), точный MaxSAT-солвер из PySAT и победитель основного трека MaxSAT Evaluation 2018. Перебором решений он не занимается. Работает снизу вверх: гоняет обычный SAT-солвер, на невыполнимом запросе достаёт конфликтующее ядро (unsat core) из мягких клауз, которые нельзя удовлетворить вместе, и «отпускает» из них ровно столько, сколько неизбежно, поднимая нижнюю оценку cost. Как только запрос стал выполнимым, накопленный cost уже доказанно минимален: первое допустимое решение и есть оптимум. То есть найденный минимум глобальный.

from pysat.formula import WCNF, IDPool
from pysat.examples.rc2 import RC2

def assign_channels(n_tx, k_ch, interference):
    vpool = IDPool()
    x = lambda t, c: vpool.id(('x', t, c))   # передатчику t назначен канал c

    wcnf = WCNF()
    # HARD: ровно один канал на передатчик
    for t in range(n_tx):
        wcnf.append([x(t, c) for c in range(k_ch)])              # >= 1
        for c1 in range(k_ch):
            for c2 in range(c1 + 1, k_ch):
                wcnf.append([-x(t, c1), -x(t, c2)])              # <= 1
    # SOFT: интерферирующая пара не должна делить канал; штраф = вес пары
    for (u, v, w) in interference:
        for c in range(k_ch):
            wcnf.append([-x(u, c), -x(v, c)], weight=w)

    with RC2(wcnf) as rc2:
        model = set(rc2.compute())
        cost = rc2.cost

    assignment = {t: c for t in range(n_tx)
                       for c in range(k_ch) if x(t, c) in model}
    return assignment, "cost: " + str(cost)


# 6 передатчиков; две 'тесные' группы 0-1-2 и 3-4-5 связаны слабым ребром (1,4)
interference = [
    (0, 1, 4), (1, 2, 3), (2, 0, 2),     # треугольник 0-1-2
    (2, 3, 5),                            # перемычка
    (3, 4, 4), (4, 5, 3), (5, 3, 2)     # треугольник 3-4-5
]

print(assign_channels(6, 2, interference))   # -> ({...}, 4)   при 2 каналах
print(assign_channels(6, 3, interference))   # -> ({...}, 0)   при 3 каналах

Что получилось. При k = 2 минимум суммарной интерференции равен 4: два нечётных треугольника двумя каналами полностью не разводятся, и оптимум жертвует двумя самыми дешёвыми конфликтами (2 + 2). При k = 3 каналов уже достаточно, cost = 0, помех нет.

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

7. Чем SAT и MaxSAT полезны на практике

За академической постановкой стоит довольно приземлённый класс задач. Чаще всего отдача видна на расписаниях и графиках смен: жёсткие правила вроде трудового кодекса, квалификаций и запрещённых сочетаний обязанностей кодируются как hard-клаузы. В свою очередь, пожелания сотрудников и приоритеты представляются как soft-клаузы с весами.

Задачи конфигурирования продуктов и комплектаций: вопрос «можно ли вообще собрать такой заказ из совместимых опций» это чистый SAT, а «если нельзя, то какая ближайшая допустимая сборка» уже требует MaxSAT. Тем же способом решается распределение дефицитного ресурса без конфликтов (частоты, слоты, причалы, лицензии, переговорные, как в разделе 6) и назначение исполнителей на заявки и маршруты, когда между ними есть правила совместимости и приоритеты.

Что из этого видно не в логах солвера, а в бюджете и в нервах того, кто отвечает за план? Здесь наследуются все свойства MIP/CP решений: полное соблюдение жестких правил, доказанное оптимальное решение, возможность провдения what-if расчётов («а что, если добавить смену, станок или канал»).

Оговорюсь честно: SAT не панацея. Он выстреливает там, где решения дискретны, правил много и они переплетены, а ручками это уже не собрать. Если задача про непрерывные объёмы и потоки, это уже к линейному и целочисленному программированию (про него я пишу в цикле Make Optimization Simple).

8. Заметки на полях: SAT/MaxSAT или MIP?

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

• Размер кодировки решает многое. Попарное «ровно один» это клауз, и на сотнях значений модель раздувается. Для масштаба берите seqcounter/totalizer из pysat.card. Дисциплина та же, о которой я писал в статье про эффективные мат. модели: лишние переменные и клаузы напрямую съедают время поиска.

• SAT хорош на чисто комбинаторных ограничениях: выполнимость, конфигурирование, верификация, головоломки, расписания с жёсткой логикой. Тут CDCL с обучением клауз нередко обгоняет ветвление MIP.

• MIP хорош там, где много линейной арифметики над непрерывными величинами и плотная числовая цель. Дробные переменные чистый SAT не умеет: либо дискретизация, либо гибрид (SMT, ленивое добавление ограничений).

• MaxSAT это полноценная дискретная оптимизация. Если задача целиком 0/1, а «арифметика» сводится к взвешенным суммам, связка WCNF + RC2 становится серьёзным, а часто и более быстрым конкурентом целочисленному солверу.

• Сменить движок дешево. Одну и ту же КНФ можно прогнать через glucose3, cadical153, minisat22 и выбрать лучший под свой профиль задач; мини-бенчмарк пишется в десяток строк.

• Дебаггинг. При возникновении противоречий в ограничениях задачи, SAT-солвер вернет минимальное конфликтующее ядро (unsat core) и прямо показывает, какие правила несовместимы. Конечно, у MIP (Irreducible Infeasible Subset) и CP (Minimal Unsatisfiable Subset) солверов есть аналоги. В случае SAT задач обнаруженные конфликты, выглядят, более информативными в силу бинарности переменных.

Заключение

SAT удобно держать в наборе рядом с линейным и целочисленным программированием. Сама задача аскетична донельзя (булевы переменные и клаузы), но CDCL-солверы и надстройка MaxSAT превращают её в рабочий инструмент дискретной оптимизации. PySAT снимает большую часть порога входа: кодирование через IDPool и CardEnc, готовый RC2 для MaxSAT и десяток сменных движков под капотом.

Путь получился такой: определение и КНФ, затем DPLL и CDCL, затем рабочий код. Сначала проверили выполнимость (раскраска графа), потом решили настоящую оптимизацию (распределение каналов) с жёсткими и взвешенными мягкими ограничениями.

Если тема будет интересна, в следующих частях Make Optimization Simple могу разобрать кодировки cardinality «под капотом», сравнить SAT и CP-SAT из OR-Tools на одной задаче расписаний или показать гибрид через SMT для смешанных моделей.

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

Что почитать дальше

• Biere, Heule, van Maaren, Walsh, «Handbook of Satisfiability», фундаментальный справочник.

• Документация PySAT: pysathq.github.io.

• Marques-Silva, Lynce, Malik, обзорная глава про CDCL.

• Описание MaxSAT-солвера RC2 (MaxSAT Evaluation).