Оптимизации динамики: когда лучший выбор сдвигается монотонно

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

Как думаешь, где в динамике лишняя работа

Задача: разрезать тысячу чисел на восемь кусков подряд так, чтобы сумма КВАДРАТОВ сумм кусков была наименьшей (так распределяют нагрузку, когда перегрузка обходится квадратично дороже). Динамика по слоям: dp[k][j] = min по i (dp[k−1][i] + (S_j − S_i)²). Восемь слоёв по n²/2 - четыре миллиона действий.

Прежде чем читать дальше, ответь себе: если для j = 500 лучшим началом последнего куска оказалось i = 430, может ли для j = 501 лучшим быть i = 100? Интуиция говорит «нет»: последний кусок сдвинулся вправо, и его начало вряд ли прыгнет далеко влево.

На что это похоже

Поиск подходящих брюк по длине ног в ряду людей, построенных по росту. Если для человека ростом 175 подошёл размер 32, то для 180 смотреть размеры меньше 32 бессмысленно. Разделяй и властвуй - это подобрать брюки сначала человеку из середины ряда, а потом каждой половине ряда искать только по свою сторону найденного размера.

Выпуклая оболочка прямых - это табло «какой тариф выгоднее при таком-то объёме». Тарифы - прямые: абонентская плата плюс цена за единицу. При любом объёме выгоднее всех один тариф, и при росте объёма выгодный тариф меняется в одну сторону - тарифы, которые не выгоднее никогда, можно выбросить сразу.

Механизм: три сокращения

РАЗДЕЛЯЙ И ВЛАСТВУЙ          условие: opt(j) ≤ opt(j + 1)
решить(jl..jr, il..ir):
    j = середина; перебрать i в [il, min(j−1, ir)] → opt(j)
    решить(jl..j−1, il..opt(j));  решить(j+1..jr, opt(j)..ir)
на каждом уровне рекурсии перебор суммарно O(n), уровней log n

КНУТ (динамика по отрезкам)
    условие: root[i][j−1] ≤ root[i][j] ≤ root[i+1][j]
перебирать корень только в этих границах; по диагонали длины d
отрезки «съедают» границы друг друга, и сумма O(n) на диагональ - O(n²)

ВЫПУКЛАЯ ОБОЛОЧКА ПРЯМЫХ
    условие: cost = a_i · x_j + b_i + (часть только от j)
каждое i - прямая y = a_i · x + b_i;  dp[j] = min прямых в x_j
прямые добавляются с убывающими наклонами, x_j растут:
    держим стек прямых нижней огибающей, лишнюю (никогда не
    минимальную) снимаем при добавлении; указатель на лучшую
    прямую только сдвигается вперёд → O(n) на слой

Откуда условия. Монотонность opt для разрезания с выпуклой ценой следует из неравенства четырёхугольника: cost(a, c) + cost(b, d) ≤ cost(a, d) + cost(b, c) при a ≤ b ≤ c ≤ d - «два перекрывающихся куска не дороже вложенных». Для квадрата суммы оно выполняется. Для прямых в нашей задаче: (S_j − S_i)² = S_i² − 2·S_i·S_j + S_j², и при фиксированном i это прямая по x = S_j с наклоном −2·S_i (наклоны убывают, потому что S_i растёт) и сдвигом prev[i] + S_i².

Правило

Динамику вида dp[j] = min по i (prev[i] + cost(i, j)) ускоряют, если у задачи есть структура. Лучшее i монотонно по j (неравенство четырёхугольника) - разделяй и властвуй, O(n log n) на слой. Лучший корень отрезка между корнями соседних отрезков - Кнут, O(n²) вместо O(n³). Цена раскладывается в a_i · x_j + b_i - минимум прямых, выпуклая оболочка, O(n) на слой при монотонных наклонах и точках. Условие проверяют до применения: без него ускорение даёт неверный ответ.

Разбор: огибающая из трёх прямых

прямые: y1 = −4x + 20,  y2 = −2x + 13,  y3 = −1x + 9
пересечения: y1 = y2 при x = 3.5;  y2 = y3 при x = 4
             y1 = y3 при x = 11/3

добавляем y3 к стеку [y1, y2]:
    y2 нужна, если она ниже обеих соседок хоть где-то: ниже y3
    она при x > 4, ниже y1 - при x < 3.5. Таких x нет →
    y2 никогда не минимальна, снимаем
стек [y1, y3]: при x < 11/3 минимальна y3, при x > 11/3 - y1

запросы x = 1, 3, 5:   y3(1) = 8,  y3(3) = 6,  y1(5) = 0
указатель: начал с y3, на x = 5 сдвинулся к y1 и назад не вернётся

Проверка «лишней прямой» в программе - то же сравнение точек пересечения, переписанное без деления: (b − b1)·(m1 − m2) ≤ (b2 − b1)·(m1 − m).

Разбор: три способа на одной задаче и Кнут на дереве поиска

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

import random

rng = random.Random(28)
n, K = 1000, 8
a = [rng.randint(1, 100) for _ in range(n)]
S = [0]
for x in a:
    S.append(S[-1] + x)
INF = float("inf")
schet = {}
# задача: разрезать a на K кусков подряд, цена куска - квадрат его суммы


def naivno():
    prev = [0] + [INF] * n
    for k in range(K):
        cur = [INF] * (n + 1)
        for j in range(1, n + 1):
            for i in range(j):        # начало последнего куска
                schet["наивно"] = schet.get("наивно", 0) + 1
                c = prev[i] + (S[j] - S[i]) ** 2
                if c < cur[j]:
                    cur[j] = c
        prev = cur
    return prev[n]


def razdelyaj():                      # лучшее i не убывает с ростом j
    prev = [0] + [INF] * n
    for k in range(K):
        cur = [INF] * (n + 1)

        def reshit(jl, jr, il, ir):   # j в [jl, jr], i в [il, ir]
            if jl > jr:
                return
            j = (jl + jr) // 2
            luchshee = il
            for i in range(il, min(j, ir + 1)):
                schet["разделяй"] = schet.get("разделяй", 0) + 1
                c = prev[i] + (S[j] - S[i]) ** 2
                if c < cur[j]:
                    cur[j], luchshee = c, i
            reshit(jl, j - 1, il, luchshee)
            reshit(j + 1, jr, luchshee, ir)

        reshit(1, n, 0, n - 1)
        prev = cur
    return prev[n]


def pryamye():                        # prev[i] + (S_j − S_i)² =
    prev = [0] + [INF] * n            # (−2·S_i)·S_j + (prev[i] + S_i²)
    for k in range(K):                # + S_j²: прямая с наклоном −2·S_i
        cur = [INF] * (n + 1)
        obolochka, p = [], 0          # наклоны убывают, S_j растёт
        for j in range(1, n + 1):
            i = j - 1
            if prev[i] < INF:
                m, b = -2 * S[i], prev[i] + S[i] ** 2
                while len(obolochka) >= 2:   # лишняя прямая - убрать
                    (m1, b1), (m2, b2) = obolochka[-2], obolochka[-1]
                    if (b - b1) * (m1 - m2) <= (b2 - b1) * (m1 - m):
                        obolochka.pop()
                    else:
                        break
                obolochka.append((m, b))
                p = min(p, len(obolochka) - 1)
            if not obolochka:         # ещё нельзя набрать k + 1 кусок
                continue
            x = S[j]
            while p + 1 < len(obolochka) and \
                    obolochka[p + 1][0] * x + obolochka[p + 1][1] <= \
                    obolochka[p][0] * x + obolochka[p][1]:
                p += 1
                schet["прямые"] = schet.get("прямые", 0) + 1
            schet["прямые"] = schet.get("прямые", 0) + 1
            m, b = obolochka[p]
            cur[j] = m * x + b + x * x
        prev = cur
    return prev[n]


def derevo_poiska(f, knut):           # оптимальное дерево поиска
    m = len(f)
    pref = [0]
    for x in f:
        pref.append(pref[-1] + x)
    e = [[0] * (m + 1) for _ in range(m + 2)]
    kor = [[0] * (m + 1) for _ in range(m + 2)]
    for i in range(1, m + 1):
        e[i][i], kor[i][i] = f[i - 1], i
    for d in range(2, m + 1):
        for i in range(1, m - d + 2):
            j = i + d - 1
            ot, do = (kor[i][j - 1], kor[i + 1][j]) if knut else (i, j)
            e[i][j] = INF
            for r in range(ot, do + 1):   # Кнут: корень между соседними
                schet[knut] = schet.get(knut, 0) + 1
                c = e[i][r - 1] + e[r + 1][j]
                if c < e[i][j]:
                    e[i][j], kor[i][j] = c, r
            e[i][j] += pref[j] - pref[i - 1]
    return e[1][m]


otvety = [naivno(), razdelyaj(), pryamye()]
assert len(set(otvety)) == 1
print("наименьшая цена:", otvety[0])
for imya in ("наивно", "разделяй", "прямые"):
    print(f"{imya:9} действий {schet[imya]}")
f = [rng.randint(1, 100) for _ in range(300)]
assert derevo_poiska(f, False) == derevo_poiska(f, True)
print("дерево поиска на 300 ключах: перебор корней", schet[False],
      " с Кнутом", schet[True])

наименьшая цена: 306954246
наивно    действий 4004000
разделяй  действий 77519
прямые    действий 13103
дерево поиска на 300 ключах: перебор корней 4544800  с Кнутом 88642

Все три способа дали одну цену. Разделяй и властвуй сократил работу в пятьдесят раз, выпуклая оболочка - в триста. Кнут на дереве поиска перебрал корни в пятьдесят раз реже и пришёл к той же цене.

Что ломается без этого

Главная ошибка - применить ускорение без условия: разделяй и властвуй на задаче, где лучшее i не монотонно, тихо пропустит оптимум - код работает, ответ хуже. Проверка простая: прогнать наивную и ускоренную версии на тысяче маленьких случайных входов. Вторая - выпуклая оболочка с немонотонными наклонами или точками: стек и указатель перестают работать, нужна оболочка с двоичным поиском или дерево Ли Чао. Третья - деление при сравнении пересечений: дробные точки и округление выкинут нужную прямую; сравнение произведений целых точно. Четвёртая - в Кнуте забыть, что границы берутся из УЖЕ посчитанных более коротких отрезков: порядок по длине обязателен.

Зачем это в работе

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

Теперь сам

dp[j] = min по i < j (prev[i] + (S_j − S_i)²). Запиши эту цену как прямую для фиксированного i: что будет наклоном, что - сдвигом и что останется общим для всех i?