Convex hull trick
Что такое convex hull trick
Convex hull - выпуклая оболочка; Convex hull trick - один из методов оптимизации динамического программирования, использующий идею выпуклой оболочки. Позволяет улучшить ассимптотику решения некоторых задач, решемых методом динамического программирования с
до . Техника впервые появилась в 1995 году (задачу на нее предложили в USACO - национальной олимпиаде США по программированию). Массовую известность получила после IOI (международной олимпиады по программированию для школьников) 2002Постановка примера задачи
Рассмотрим задачу на ДП:
Есть n деревьев с высотами(в метрах). Требуется спилить их все, потратив минимальное количество монет на заправку бензопилы. Но пила устроена так, что она может спиливать только по 1 метру от дерева, к которому ее применили. Также после срубленного метра (любого дерева) пилу нужно заправлять, платя за бензин определенной кол-во монет. Причем стоимость бензина зависит от срубленных (полностью) деревьев. Если сейчас максимальный индекс срубленного дерева равен i, то цена заправки равна ci. Изначально пила заправлена. Также известны следующие ограничения : возрастают, убывают. Изначально пила заправлена. (убывание и возрастание нестрогие)
(Задача H отсюда : http://neerc.ifmo.ru/school/camp-2016/problems/20160318a.pdf)
Наивное решение
Перед началом решения заметим важный факт, напрямую следующий из условия : т.к. c[i] убывают(нестрого) и c[n] = 0, то все c[i] не отрицательны.
Понятно, что нужно затратив минимальную стоимость срубить последнее (
-е) дерево, т.к. после него все деревья можно будет пилить бесплатно (т.к. ). Посчитаем следующую динамику : - минимальная стоимость, заплатив которую будет срублено дерево номер Тогда . То есть понятно, что выгодно рубить сначала более дорогие и низкие деревья, а потом более высокие и дешевые (док-во этого факта оставляется читателям как несложное упражнение). Тогда переберем - индекс предыдущего срубленного дерева. Пусть мы его срубили оптимальным (в смысле денег) способом. Тогда просто раз уменьшим высоту дерева i на 1. Каждый такой раз будем платить за последующую заправку пилы. Итак, на сруб -го дерева мы заплатили . Нетрудно видеть, что такая динамика работает за . Посмотрим на код наивного решения:for i = 1..n-1 {
dp[i] = +
for j = 0..i-1 {
if (dp[j] + a[i] * c[j] < dp[i])
dp[i] = dp[j] + a[i] * c[j]
}
}
О-Оптимизация
Давайте обозначим
за , за , а за . Теперь формула приняла вид . Выражение напоминает уравнение прямой . Сопоставим каждому , обработанному ранее прямую . Из условия « убывают уменьшаются с номером » следует то, что прямые, полученные ранее отсортированы в порядке убывания углового коэффициент. Давайте нарисуем несколько таких прямых :Итак, давайте выделим множество точек
, таких что все они принадлежат одной из прямых и при этом нету ни одной прямой , такой что . Иными словами возьмем «выпуклую (вверх) оболочку» нашего множества прямых. На картинке множество точек выделено жирным оранжевым цветом и представляет собой выпуклую вверх функцию. Назовем ее « »Для чего нам нужна эта выпуклая оболочка прямых?
Действительно, пусть новая прямая пересекает последнюю прямую множества позже, чем предпоследнюю (рис.2 - красная прямая новая, фиолетовая - предпоследняя, желтая - последняя), то найдется такой отрезок
, что последняя(желтая) прямая при этих x-ах лежит ниже всех остальных и ее следует оставить в множестве. Теперь пусть новая прямая пересекает последнюю прямую множества раньше, чем предпоследнюю (рис.1), последняя прямая при любых x лежит выше какой-то другой прямой множества и значит ее можно удалить, чтд.Детали реализации:
Будем хранить 2 массива (имитирующих стеки) :
и - начало (по x) соответствующей прямой выпуклой оболочки и номер этой прямой (в глобальной нумерации). Также воспользуемся тем, что возрастают (по условию задачи), а значит мы можем искать первое такое , что не бинарным поиском, а методом 2х указателей за суммарно. Также массив front[] можно хранить в целых числах, округляя х-координаты в сторону лежащих правее по оси x до ближайшего целого (*). Почему так? На самом деле мы, считая динамику, подставляем в уравнения прямых только целые , а значит если -я прямая пересекается с -й в точке (z-целое, ), то мы будем подставлять или . Поэтому можно считать, что новая прямая имеет "область действия", начиная сР.Реализация
void Convex-hull-trick st[0] = 0 from[0] = -// первая прямая покрывает все x-ы, начиная с -∞ sz = 1 // текущий размер выпуклой оболочки pos = 0 // текущая позиция первого такого j, что x[i] >= front[st[j]] for i = 1..n-1 { while (front[pos] < x[i]) // метод 1 указателя (ищем первое pos, такое что x[i] покрывается "областью действия" st[pos]-той прямой ++pos j = st[pos] dp[i] = K[j] * a[i] + B[j] if (i < n - 1) { // если у нас добавляется НЕ последняя прямая, то придется пересчитать выпуклую оболочку K[i] = c[i] // наши переобозначения переменных B[i] = dp[i] // наши переобозначения переменных x = - while (true) { j = st[sz - 1] x = divide(B[j] - B[i], K[i] - K[j]) // x-координата пересечения с последней прямой оболочки, округленное в нужную сторону (*) if (x > from[sz - 1]) break // перестаем удалять последнюю прямую из множества, если новая прямая пересекает ее позже, чем начинается ее "область действия" --sz // удаляем последнюю прямую, если она лишняя } st[sz] = i from[sz++] = x // добавили новую прямую } }
(Здесь функция divide(a, b) возвращает нужное(*) округление a / b) Такая реализация будет работать за O(n).
Динамический convex hull trick
Заметим, что условия на прямые, что http://neerc.ifmo.ru/school/camp-2016/problems/20160318a.pdf).
Но рассмотрим общий случай. Наша задача поменялась следующим образом : по-прежнему у нас есть выпуклая оболочка, имея которую мы за или быстрее можем найти , но теперь вставку i-й прямой в оболочку уже нельзя выполнить старым способом за (в среднем). У нас есть выпуклая оболочка, наша прямая пересекает ее, возможно, «отрезая» несколько отрезков выпуклой оболочки в середине (рис. 4).
Т.е. нужно уметь быстро (за
?) назодить, после какой прямой стоит пытаться вставить текущую (красную рис.4) прямую и удалять лишние справа, начиная с нее, потом проводить аналогичные операции слева. Итак, давайте хранить (или любой аналог в других языках) пар = <коэффицент прямой, ее номер в глобальной нумерации>. Когда приходит новая прямая, делаем lower_bound - 1 в сете, т.е. ищем ближайшую прямую с меньшим углом наклона, и начиная с нее повторяем старый алгоритм (удаляем, пока прямая бесполезная). И симметричный алгоритм применяем ко всем прямым справа от нашей. Асимптотика решения составит на каждый из n запросов «добавить прямую» + суммарно на удаление прямых. Итого .Альтернативный подход
Другой способ интерпретировать выражение
заключается в следующем: давайте перепишем выражение , т.е. запишем как скалярное произведение векторов и . Вектора хотелось бы организовать так, чтобы за находить максимизирующий выражение . Посмотрим на рис. 5. Заметим довольно очевидный факт : красная точка(вектор) не может давать более оптимальное значение одновременно чем обе синие точки, т.к. - это на самом деле проекция вектора на . По этой причине нам достаточно оставить выпуклую оболочку векторов , а ответ на запрос - это поиск , максимизирующего проекцию на . Это задача поиска ближайшей точки выпуклого многоугольника (составленного из точек выпуклой оболочки) к заданной прямой (из в ). Ее можно решить за O(logn) двумя бинарными или одним тернарным поиском Асимптотика алгоритма по-прежнему составит