Изменения

==Что такое convex Convex hull trick==Convex hull {{-- выпуклая оболочка; Convex hull trick - }} один из методов оптимизации [[Динамическое_программирование | динамического программирования]], использующий идею [[Статические_выпуклые_оболочки:_Джарвис,_Грэхем,_Эндрю,_Чен,_QuickHull|выпуклой оболочки]]. Позволяет улучшить ассимптотику асимптотику решения некоторых задач, решемых методом динамического программирования , с <math>O(n^2)</math> до <mathtex>O(n*logn\cdot\log(n))</mathtex>. Техника впервые появилась в 1995 году (задачу на нее предложили в USACO {{- --}} национальной олимпиаде США по программированию). Массовую известность получила после IOI (международной олимпиады по программированию для школьников) 2002.

{{Задача|definition = Есть <math>n </math> деревьев с высотами <mathtex>a1a_1, a2a_2, \dots an, a_n</mathtex> (в метрах). Требуется спилить их все, потратив минимальное количество монет на заправку бензопилы. Но пила устроена так, что она может спиливать только по <math>1 </math> метру от дерева, к которому ее применили. Также после срубленного метра (любого дерева) пилу нужно заправлять, платя за бензин определенной определенное кол-во монет. Причем стоимость бензина зависит от срубленных (полностью) деревьев. Если сейчас максимальный индекс срубленного дерева равен <tex>i</tex>, то цена заправки равна ci<tex>c_i</tex>. Изначально пила заправлена. Также известны следующие ограничения : <mathtex>c[n] c_n = 0, a[1] a_1 = 1, a[i]a_i</mathtex> возрастают, <mathtex>c[i]c_i</mathtex> убывают. Изначально пила заправлена. (убывание и возрастание нестрогие)}}(Задача H отсюда : с Санкт-Петербургских сборов к РОИ 2016 <ref>[http://neerc.ifmo.ru/school/camp-2016/problems/20160318a.pdf Сайт с задачами Санкт-Петербургских сборов к РОИ 2016]</ref>)</noinclude><includeonly>{{#if: {{{neat|}}}|<div style="background-color: #fcfcfc; float:left;"><div style="background-color: #ddd;">'''Задача:'''</div><div style="border:1px dashed #2f6fab; padding: 8px; font-style: italic;">{{{definition}}}</div></div>|<table border="0" width="100%"><tr><td style="background-color: #ddd">'''Задача:'''</td></tr><tr><td style="border:1px dashed #2f6fab; padding: 8px; background-color: #fcfcfc; font-style: italic;">{{{definition}}}</td></tr></table>}}</includeonly>

Сначала заметим важный факт : т.к. <tex>c[i] </tex> убывают(нестрого) и <tex>c[n] = 0</tex>, то все <tex>c[i] </tex> неотрицательны.Понятно, что нужно затратив минимальную стоимость срубить последнее (<mathtex>n</mathtex>-е) дерево, т.к. после него все деревья можно будет рубить бесплатно (т.к. <mathtex>c[n] = 0</mathtex>). Посчитаем следующую динамику : <mathtex>dp[i]</mathtex> {{--- }} минимальная стоимость, заплатив которую можно добиться того, что дерево номер <mathtex>i.</mathtex> будет срублено. База динамики : <mathtex>dp[1] = 0</mathtex>, т.к. изначально пила заправлена и высота первого дерева равна <math>1</math>, по условию задачи.Переход динамики : понятно, что выгодно рубить сначала более дорогие и низкие деревья, а потом более высокие и дешевые (док-во этого факта оставляется читателям как несложное упражнение, т.к. эта идея относится скорее к теме [[Теорема_Радо-Эдмондса_(жадный_алгоритм)|жадных алгоритмновалгоритмов]], чем к теме данной статьи). Поэтому перед <mathtex>i</mathtex>-м деревом мы обязательно срубили какое-то <mathtex>j</mathtex>-е, причем <mathtex>j < \leqslant i- 1</mathtex>. Поэтому чтобы найти <mathtex>dp[i]</mathtex> нужно перебрать все <mathtex>1 <= \leqslant j < \leqslant i- 1</mathtex> и попытаться использовать ответ для дерева намер номер <mathtex>j</mathtex>. Итак, пусть перед <mathtex>i</mathtex>-м деревом мы полностью срубили <mathtex>j</mathtex>-е, причем высота <mathtex>i</mathtex>-го дерева составляет <mathtex>a[i]</mathtex>, а т.к. последнее дерево, которое мы срубили , имеет индекс <mathtex>j</mathtex>, то стоимость каждого метра <mathtex>i</mathtex>-го дерева составит <mathtex>c[j]</mathtex>. Поэтому на сруб <mathtex>i</mathtex>-го дерева мы потратим <mathtex>a[i] * \cdot c[j]</mathtex> монет. Также не стоит забывать, что ситуацию, когда <mathtex>j</mathtex>-е дерево полностью срублено, мы получили не бесплатно, а за <mathtex>dp[j]</mathtex> монет.Итогвая Итоговая формула пересчета : <mathtex>dp[i] = min_\min\limits_{j=01...i-1}(dp[j] + a[i] * \cdot c[j])</mathtex>. 
 Посмотрим на код выше описанного вышеописанного решения: '''int''' <tex>\mathtt{simpleDP}</tex>('''int''' a[n], '''int''' c[n]) dp[1] = 0 dp[2] = dp[3] = ... = dp[n] = <tex>\infty</tex> '''for'' ' i = 2 = 1..n-1 { dp[i] = +<tex>+\infty</tex> '''for''' j = 01..i-1 { '''if''' (dp[j] + a[i] * <tex>\cdot</tex> c[j] < dp[i]) dp[i] = dp[j] + a[i] * <tex>\cdot</tex> c[j] } } '''return''' dp[n]Нетрудно видеть, что такая динамика работает за <mathtex>O(n^2)</mathtex>.

1) Сделаем Для начала сделаем замену обозначений. Давайте обозначим <mathtex>dp[j]</mathtex> за <mathtex>b[j]</mathtex>, <mathtex>a[i]</mathtex> за <mathtex>x[i]</mathtex>, а <mathtex>c[j]</mathtex> за <mathtex>k[j]</mathtex>.

2) Теперь формула приняла вид <mathtex>dp[i] = min_\min\limits_{j=0...i-1}(k[j]*\cdot x[i] + b[j])</mathtex>. Выражение <mathtex>k[j]*\cdot x + b[j]</mathtex> {{- --}} это в точности уравнение прямой вида <mathtex>y = kx + b</mathtex>.

3) Сопоставим каждому <mathtex>j</mathtex>, обработанному ранее, прямую <mathtex>y[j](x) = k[j]*\cdot x + b[j]</mathtex>. Из условия «<mathtex>c[i]</mathtex> убывают <math><=tex> \Leftrightarrow k[j]</mathtex> уменьшаются с номером <mathtex>j</mathtex>» следует то, что прямые, полученные ранее отсортированы в порядке убывания углового коэффициент. Давайте нарисуем несколько таких прямых :

4) Выделим множество точек <mathtex>(x0x_0, y0y_0)</mathtex> , таких что все они принадлежат одной из прямых и при этом нету ни одной прямой <mathtex>y’(x)</mathtex>, такой что <mathtex>y’(x0x_0) < y0y_0</mathtex>. Иными словами возьмем «выпуклую (вверх) оболочку» нашего множества прямых (её еще называют нижней ошибающей огибающей множества прямых на плоскости). Назовем ее «<mathtex>y = convex(x)</mathtex>». Видно, что множество точек <math>(x, convex(x)) </math> представляет собой выпуклую вверх функцию.

==Для чего нужна нижняя ошибающая Цель нижней огибающей множества прямых== Пусть мы считаем динамику для <mathtex>i</mathtex>-го дерева. Его задает <mathtex>x[i]</mathtex>. Итак, нам нужно для данного <mathtex>x[i]</mathtex> найти <mathtex>min_\min\limits_{j=0..i-1}(k[j]*\cdot x[i] + b[ij]) = min_\min\limits_{j=0..i-1}(y[j](x[i]))</mathtex>. Это выражение есть <math>convex(x[i])</math>. Из монотонности угловых коэффицентов отрезков, задающих выпуклую оболочку, и их расположения по координаты координатам x следует то, что отрезок, который пересекает прямую <mathtex>x = x[i]</mathtex>, можно найти [[Целочисленный_двоичный_поиск|бинарным поиском]]. Это потребует <mathtex>O(logn\log(n))</mathtex> времени на поиск такого <tex>j</tex>, что <tex>dp[i] = k[j] * \cdot x[i] + b[j]</tex>. Теперь осталось научиться поддерживать множество прямых и быстро добавлять <mathtex>i</mathtex>-ю прямую после того, как мы посчитали <mathtex>b[i] = dp[i]</mathtex>.

Воспользуемся идеей алгоритма построения выпуклой оболочки множества точек. Заведем 2 стека <mathtex>k[]</mathtex> и <mathtex>b[]</mathtex>, которые задают прямые в отсортированном порядке их угловыми коэффицентами и свободными членами. Рассмотрим ситуацию, когда мы хотим добавить новую (<mathtex>i</mathtex>-тую) прямую в множество. Найдем точки Пусть сейчас в множестве лежит <tex>sz</tex> прямых (нумерация с 1). Пусть <tex>(x_L, y_L)</tex> {{---}} точка пересечения <tex>sz - 1</tex>-й прямой множества и <tex>sz</tex>-й, а <tex>(x_R, y_R)</tex> {{---}} точка пересечения с последними 2мя прямыми из стекановой прямой, которую мы хотим добавить в конец множества и <tex>sz</tex>-й. Нас будут интересовать только их <mathtex>x</mathtex>-овые координаты. Назовем их <mathtex>xLx_L</mathtex> и <mathtex>xRx_R</mathtex>, соответственно. Если оказалось, что новая прямая пересекает предпоследнюю <tex>sz</tex>-ю прямую выпуклой оболочки позже, чем последнюю <mathtex>sz</tex>-я <tex>sz - 1</tex>-ю, т.е. <tex>(xL >= xRx_L \geqslant x_R)</mathtex>, то последнюю <tex>sz</tex>-ю удалим из нашего множества, иначе - остановимся. Так будем делать, пока либо кол-во число прямых в стеке не станет равным 2, либо <mathtex>xLx_L</mathtex> не станет меньше <mathtex>xRx_R.</mathtex>

Асимптотика : аналогично обычному алгоритму построения выпуклой оболочки, каждая прямая ровно <math>1 </math> раз добавится в стек и максимум <math>1 </math> раз удалится. Значит время работы перестройки выпуклой оболочки займет <mathtex>O(n)</mathtex> суммарно.

Корректность : достаточно показать, что последнюю прямую нужно удалить из множества т.и т.т., когда она наша новая прямая пересекает ее в точке с координатой по оси X, меньшей, чем предпоследнюю. [[Файл:picture2convexhull.png]]

План доказательства : пусть новая прямая пересекает последнюю прямую {{Теорема|id=th1239. |statement=Алгоритм построения нижней огибающей множества позже, чем предпоследнюю (риспрямых корректен.2 - красная прямая новая, фиолетовая - предпоследняя, желтая - последняя), тогда найдется такой отрезок <math>[x1;x2]</math>|proof=Достаточно показать, что последняя (желтая) прямая при последнюю прямую нужно удалить из множества <tex>x \in [x1;x2]Leftrightarrow</tex> лежит ниже всех остальных и ее следует оставить в множестве.А если , когда наша новая прямая пересекает последнюю прямую множества раньшеее в точке с координатой по оси X, меньшей, чем предпоследнюю (рис.1), последняя прямая при любых x лежит выше какой-то другой прямой множества и значит ее нужно удалитьпредпоследнюю.

Более формальноПусть <tex>Y(x) = Kx + B</tex> {{---}} уравнение новой прямой, <tex>y[i](x) = K[i]x + B[i]</tex> {{---}} уравнения прямых множества. Тогда т.к. <tex>K < K[sz]</tex>, то при <tex>x \in [- \infty; x_R] : y[sz](x) <= Y(x)</tex>, а т.к. <tex> K[sz] < K[sz - 1]</tex>, то при <tex>x \in [x_L; + \infty] : y[sz - 1](x) \geqslant y[sz](x)</tex>. Пусть Если <tex>x_L < x_R</tex>, то при <tex>x \in [x_L; x_R] : y[sz - 1] \geqslant y[sz](x) и Y(x) \geqslant y[sz](x)</tex>, т.е. на отрезке <tex>[x_L; x_R]</tex> прямая номер sz лежит ниже остальных и её нужно оставить в множестве. Если же <tex>x_L > x_R<math/tex>xL , то она ниже всех на отрезке <tex>[x_L; x_R] = \varnothing </mathtex>, т.е. её можно удалить из множества.}}

Будем хранить 2 массива : <tex>front[]</tex> {{---}} <tex>x</tex>-координаты, начиная с которых прямые совпадают с выпуклой оболочкой (имитирующих стеки) : т.е. i-я прямая совпадает с выпуклой оболочкой текущего множества прямых при <tex>x</tex> <tex>\in</tex> <mathtex>[front[i]; front[i + 1])</mathtex> ) и <mathtex>st[]</mathtex> {{--- начало }} номера деревьев, соответствующих прямым (по x) соответствующей прямой выпуклой оболочки и т.е. <tex>i</tex>-я прямая множества, где <tex>i</tex> <tex>\in</tex> <tex>[1; sz]</tex> соответствует дереву номер этой прямой (в глобальной нумерации<tex>sz[i]</tex>). Также воспользуемся тем, что <mathtex>x[i] = a[i]</mathtex> возрастают (по условию задачи), а значит мы можем искать первое такое <mathtex>j</mathtex>, что <mathtex>x[i] >= \geqslant front[j]</mathtex> не бинарным поиском, а методом 2х двух указателей за <mathtex>O(n)</mathtex> операций суммарно. Также массив <math>front[] </math> можно хранить в целых числах, округляя х-координаты в сторону лежащих правее по оси x до ближайшего целого (*), т. Почему так? На к. на самом деле мы, считая динамику, подставляем в уравнения прямых только целые <mathtex>x[i]</mathtex>, а значит если <mathtex>k</mathtex>-я прямая пересекается с <mathtex>k+1</mathtex>-й в точке <mathtex>z +</mathtex> <tex>\alpha</tex> (<math>z</math>-целое, <tex>\alpha</tex> <tex>\in</tex> <mathtex>[0;1)</mathtex>), то мы будем подставлять в их уравнения <mathtex>z</mathtex> или <mathtex>z + 1</mathtex>. Поэтому можно считать, что новая прямая имеет "область действия"начинает совпадать с выпуклой оболочкой, начиная с <mathtex>x = z+1</mathtex>

'''voidint''' Convex-hull-trick<tex>\mathtt{ConvexHullTrick}</tex>('''int''' a[n], '''int''' c[n]) st[01] = 01 from front[01] = -<tex>\infty</tex><font color=green>// первая прямая покрывает все x-ы, начиная с -∞ </font> sz = 1 <font color=green>// текущий размер выпуклой оболочки </font> pos = 0 1 <font color=green>// текущая позиция первого такого j, что x[i] >= \geqslant front[st[j]] </font > '''for''' i = 12..n-1 { '''while''' (front[pos] < x[i]) <font color=green>// метод 1 указателя (ищем первое pos, такое что x[i] покрывается "областью действия" st[pos]-той прямой </font > pos = pos ++pos1 j = st[pos] dp[i] = K[j] * <math>\cdot</math> a[i] + B[j] '''if''' (i < n - 1) { <font color=green>// если у нас добавляется НЕ последняя прямая, то придется пересчитать выпуклую оболочку </font > K[i] = c[i] <font color=green>// наши переобозначения переменных </font > B[i] = dp[i] <font color=green>// наши переобозначения переменных </font > x = -<tex>\infty</tex> '''while''' (''true) {'' j = st[sz - 1] x = divide(B[j] - B[i], K[i] - K[j]) <font color=green>// x-координата пересечения с последней прямой оболочки, округленное в нужную сторону (*) </font > '''if ''' (x > from[sz - 1]) '''break''' <font color=green>// перестаем удалять последнюю прямую из множества, если новая прямая пересекает ее позже, чем начинается ее "область действия" </font > sz = sz --sz 1<font color=green>// удаляем последнюю прямую, если она лишняя </font > } st[sz+ 1] = i from front[sz++1] = x <font color=green>// добавили новую прямую </font > } sz = sz + 1 } '''return''' dp[n] (Здесь функция <tex>\mathtt{divide(a, b) }</tex> возвращает нужное(*) округление <tex>\frac{a }{b}</tex>. Приведем её код : '''int''' <tex>\mathtt{divide}</ tex>('''int''' a, '''int''' b) delta = 0 '''if''' (a '''mod''' b ≠ 0) delta = 1 '''if''' ((a > 0 '''and''' b > 0) '''or''' (a < 0 '''and''' b< 0)) '''return''' [a / b] + delta '''return''' -[|a| / |b|]  Такая реализация будет работать за <math>O(n)</math>.

Заметим, что условия на прямые, что возрастание/убывание <mathtex>k[i]</mathtex> возрастаетна убывание/убывает возрастание и <mathtex>x[i]</mathtex> убывает/возрастает выглядят не достаточно общимиредкими для большинства задач. Как же быть, если Пусть в задаче таких ограничений нет. Иногда можно Первый способ борьбы с этой проблемой - отсортировать прямые входные данные нужным образом, не испортив свойств задачи (пример : задача G отсюда c Санкт-Петербургских сборов к РОИ 2016 <ref>[http://neerc.ifmo.ru/school/camp-2016/problems/20160318a.pdfСайт с задачами Санкт-Петербургских сборов к РОИ 2016]</ref>).  Но рассмотрим общий случай. Наша задача поменялась следующим образом : поПо-прежнему у нас есть выпуклая оболочкапрямых, имея которую с помощью которой мы за <mathtex>O(logn\log(n))</mathtex> или быстрее можем найти <mathtex>dp[i]</mathtex>, но теперь вставку <tex>i</tex>-й прямой в оболочку уже нельзя выполнить старым описанным ранее способом за <mathtex>O(1)</mathtex> (в среднем). У нас есть выпуклая оболочка, наша прямая пересекает ее, возможно, «отрезая» «отсекая» несколько отрезков выпуклой оболочки в середине (рис. 4: красная прямая - та, которую мы хотим вставить в наше множество). Более формально : теперь наша новая прямая будет ниже остальных при <tex>x \in [x_1; x_2]</tex>, где <tex>x_1, x_2 \in R</tex> - точки пересечения с некоторыми прямыми, причем <tex>x_2</tex> не обязательно равно <tex>+ \infty</tex>

Т.е. нужно Чтобы уметь быстро (за <math>O(logn)</math>?) назодить, после какой прямой стоит пытаться вставить текущую (красную рис.4) вставлять прямую и удалять лишние справав множество будем хранить [[Красно-черное_дерево|двоичное дерево поиска]], начиная с нее, потом проводить аналогичные операции слева. Итак, давайте хранить <math>std::setв вершинах которого будут пары </mathtex> (или любой аналог в других языках) пар <math><k, st>)</mathtex> = <(коэффицент прямой, ее номер в глобальной нумерации>). Когда приходит новая прямая, делаем lower_bound - 1 в сете, т.е. ищем ближайшую последнюю прямую с меньшим углом наклонаугловым коэффицентом, чем у той прямой, и начиная которую мы хотим добавить в множество. Поиск такой прямой занимает <tex>O(\log(n))</tex>. Начиная с нее повторяем найденной прямой выполняем "старый " алгоритм (удаляем, пока текущая прямая бесполезнаямножества бесполезна). И симметричный алгоритм применяем ко всем прямым справа от нашей(удаляем правого соседа нашей прямой, пока она пересекает нас позже, чем своего правого соседа). Асимптотика решения составит <mathtex>O(logn\log(n))</mathtex> на каждый из <tex>n </tex> запросов «добавить прямую» + <mathtex>O(n\cdot\log(n))</mathtex> суммарно на удаление прямых, т.к. по-прежнему каждая прямая не более одного раза удалится из множества, а каждое удаление из std::set занимает <tex>O(\log(n))</tex> времени. Итого <math>O(nlognn\cdot\log(n))</math>.

Другой способ интерпретировать выражение <mathtex>dp[i] = min_\min\limits_{j=0...i-1}(c[j]*\cdot a[i] + dp[j])</mathtex> заключается в следующемтом, что мы будем пытаться свести задачу к стандартной выпуклой оболочке множества точек. Перепишем выражение средующим образом : давайте перепишем выражение <mathtex>dp[j] + a[i] * \cdot c[j] = (dp[j], c[j]) * \cdot (1, a[i])</mathtex>, т.е. запишем как скалярное произведение векторов <mathtex>v[j] = (dp[j], c[j])</mathtex> и <mathtex >u[i] = (1, a[i])</mathtex >. Вектора <mathtex >v[j] = (dp[j], c[j])</mathtex> хотелось бы организовать так, чтобы за <mathtex >O(logn\log(n))</mathtex> находить вектор, максимизирующий выражение <mathtex>v[j] * \cdot u[i]</mathtex>. Посмотрим на рис. 5. Заметим довольно интуитивно очевидный факт : красная точка(вектор) <mathtex>j</mathtex> не может давать более оптимальное значение <mathtex>v[j] * \cdot u[i]</mathtex> одновременно чем обе синие точки, т.к. <math>v[j] * u[i]</math> - это на самом деле проекция вектора <math>v[j]</math> на <math>u[i]</math>. По этой причине нам достаточно оставить выпуклую оболочку векторов <mathtex>v[j]</mathtex>, а ответ на запрос {{- --}} это поиск <mathtex>v[j]</mathtex>, максимизирующего проекцию на <mathtex>u[i]</mathtex>. Это задача поиска ближайшей точки выпуклого многоугольника (составленного из точек выпуклой оболочки) к заданной прямой (из <mathtex>(0, 0)</mathtex> в <mathtex>(1, a[i])</mathtex>). Ее можно решить за <tex>O(logn\log(n)) </tex> двумя бинарными или одним тернарным поискомАсимптотика алгоритма по-прежнему составит <mathtex>O(n*\cdot \log(n))</mathtex>