Алгоритм Форчуна: заметающая прямая и береговая линия из парабол

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

Как думаешь, почему это трудно

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

Ответь себе до чтения: какая часть плоскости над прямой точно не изменится, что бы ни лежало ниже? Ответ урока: точки, до которых от какого-нибудь пройденного центра не дальше, чем до самой прямой. Любой будущий центр лежит ниже прямой, то есть заведомо дальше. Граница этой области - береговая линия.

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

Отлив. Вода уходит вниз, и над ней обнажается песок - он уже точно сухой. Граница воды не прямая, а волнистая: около каждого камня песок высыхает раньше. Камни - пройденные центры, линия воды - береговая линия, сухой песок - готовая часть диаграммы. Когда вода отступает ниже нового камня, около него появляется новая выпуклость берега.

Механизм: параболы, изломы и два вида событий

Парабола. Прямая стоит на высоте l, центр p выше неё. Точки, равноудалённые от p и от прямой, - парабола: y = ((x − px)² + py² − l²) / (2 · (py − l)). Она получается из равенства квадратов расстояний (x − px)² + (y − py)² = (y − l)²: квадраты y² сокращаются, и остаётся линейное уравнение относительно y. Точка над параболой ближе к p, чем к прямой, - а значит, ближе к p, чем к любому будущему центру. Береговая линия - нижняя огибающая всех таких парабол: над ней каждая точка уже знает своего ближайшего центра.

Изломы. Две соседние дуги береговой линии встречаются в точке, равноудалённой от двух центров и от прямой. Когда прямая опускается, этот излом движется - и движется он по серединному перпендикуляру двух центров, то есть по ребру диаграммы. Береговая линия, опускаясь, прочерчивает рёбра Вороного.

Событие «центр». Прямая доходит до нового центра. Его парабола в этот момент - вертикальный луч вверх, и он разрезает дугу, над которой стоит, на две: последовательность … X … становится … X Y X …. Отсюда граница: каждое такое событие добавляет не больше двух дуг, всего дуг не больше 2n − 1.

Событие «окружность». Три соседние дуги X Y Z, у которых изломы по краям средней дуги сходятся (так бывает не у каждой тройки - иногда изломы расходятся, и средняя дуга только растёт; проверка - поворот центров X, Y, Z слева направо: изломы сходятся, когда он направо, по часовой стрелке). Средняя сжимается и исчезает в точке, которая равноудалена от X, Y, Z и от прямой, - это центр окружности через три центра, а прямая в этот момент касается окружности снизу. Значит, событие наступает, когда прямая доходит до нижней точки окружности. Центр окружности - вершина Вороного. Если раньше появится новый центр внутри этой окружности - событие отменяется: окружность перестала быть пустой.

Структуры и цена. Статус - дуги береговой линии в сбалансированном дереве. Очередь событий: все центры по y сразу и события-окружности, которые добавляются и отменяются по ходу - при каждом изменении тройки соседних дуг. Событий O(n), каждое - логарифм, итого O(n · log n). Быстрее нельзя (в той же модели алгебраических решений): по диаграмме выпуклая оболочка читается за линию, а она не дешевле сортировки.

Правило

Алгоритм Форчуна ведёт прямую сверху вниз и держит береговую линию - нижнюю огибающую парабол пройденных центров; над ней диаграмма окончательна, изломы чертят рёбра. Событие «центр» разрезает дугу, событие «окружность» - в нижней точке пустой окружности через три соседние дуги - убирает дугу и рождает вершину. Дуг ≤ 2n − 1, время O(n · log n).

Разбор: береговая линия на семи уровнях

Пять центров из урока о диаграмме Вороного. Для каждого уровня прямой находим, какая парабола ниже всех над каждой точкой x от −40 до 80 с шагом 1/4, и записываем дуги по порядку. Высоты парабол - дроби, сравниваются в целых.

centry = [[6, 3], [30, 4], [18, 9], [8, 13], [33, 14]]
bukvy = "ABCDE"

def bereg(l):
    # береговая линия при заметающей прямой y = l, идущей сверху вниз:
    # над каждой точкой x - ближайшая из парабол пройденных центров.
    # x от -40 до 80 с шагом 1/4: x = k / 4. Высота параболы центра p -
    # дробь ((x - px)² + py² - l²) / (2 (py - l)): крест-накрест
    duga = []
    for k in range(-160, 321):
        luchshij = -1
        for i in range(len(centry)):
            p = centry[i]
            if p[1] <= l:
                continue
            ch = (k - 4 * p[0]) ** 2 + 16 * (p[1] * p[1] - l * l)
            zn = 32 * (p[1] - l)
            if luchshij == -1 or ch * zl < cl * zn:
                luchshij = i
                cl = ch
                zl = zn
        if luchshij >= 0 and (len(duga) == 0 or duga[-1] != luchshij):
            duga.append(luchshij)
    stroka = ""
    for i in duga:
        stroka = stroka + bukvy[i]
    return stroka

for l in [13, 11, 8, 6, 2, -4, -20]:
    print(l, bereg(l))

13 E
11 DE
8 DCE
6 DCE
2 DADCBE
-4 DACBE
-20 DACBE

Прямая на высоте 13 прошла только центр E (y = 14) - одна дуга. Ниже 13 появилась D, ниже 9 - C. Ниже 4 и 3 пришли B и A, и на уровне 2 видно событие «центр» в действии: A разрезала дугу D на две, отсюда D A D. Между уровнями 2 и −4 средняя D исчезла - это событие «окружность». К уровню −20 больше ничего не изменилось: новых центров нет, а следующие события-окружности - если они есть - наступают ниже.

Строка DE на уровне 11 показывает предел окна: когда D разрезала дугу E, кусок E слева от D ушёл левее x = −40, и в окне его нет. Береговая линия бесконечна в обе стороны.

Разбор: когда исчезает средняя D

Средняя дуга D зажата между A и C. Она исчезнет, когда прямая коснётся снизу окружности через A, D и C. Посчитаем эту окружность и посмотрим на береговую линию на 0,3 выше и ниже события. Уровень прямой здесь дробный, поэтому высоты парабол сравниваются во float, - для того чтобы посмотреть на дуги, этого хватает.

centry = [[6, 3], [30, 4], [18, 9], [8, 13], [33, 14]]
bukvy = "ABCDE"

def bereg(l):
    # береговая линия при заметающей прямой y = l, идущей сверху вниз:
    # над каждой точкой x - ближайшая из парабол пройденных центров.
    # x от -40 до 80 с шагом 1/4: x = k / 4. Высота параболы центра p -
    # дробь ((x - px)² + py² - l²) / (2 (py - l)): крест-накрест
    duga = []
    for k in range(-160, 321):
        luchshij = -1
        for i in range(len(centry)):
            p = centry[i]
            if p[1] <= l:
                continue
            ch = (k - 4 * p[0]) ** 2 + 16 * (p[1] * p[1] - l * l)
            zn = 32 * (p[1] - l)
            if luchshij == -1 or ch * zl < cl * zn:
                luchshij = i
                cl = ch
                zl = zn
        if luchshij >= 0 and (len(duga) == 0 or duga[-1] != luchshij):
            duga.append(luchshij)
    stroka = ""
    for i in duga:
        stroka = stroka + bukvy[i]
    return stroka

import math

def okruzhnost(a, b, c):
    # центр и радиус окружности через три точки
    sa = a[0] ** 2 + a[1] ** 2
    sb = b[0] ** 2 + b[1] ** 2
    sc = c[0] ** 2 + c[1] ** 2
    d = a[0] * (b[1] - c[1]) + b[0] * (c[1] - a[1])
    d = d + c[0] * (a[1] - b[1])
    ux = sa * (b[1] - c[1]) + sb * (c[1] - a[1]) + sc * (a[1] - b[1])
    uy = sa * (c[0] - b[0]) + sb * (a[0] - c[0]) + sc * (b[0] - a[0])
    ux = ux / (2 * d)
    uy = uy / (2 * d)
    return [ux, uy, math.sqrt((a[0] - ux) ** 2 + (a[1] - uy) ** 2)]

ux, uy, r = okruzhnost(centry[0], centry[3], centry[2])
sobytie = uy - r
print(round(ux, 3), round(uy, 3), round(r, 3), round(sobytie, 3))
print(bereg(sobytie + 0.3))
print(bereg(sobytie - 0.3))

11.444 7.111 6.822 0.289
DADCBE
DACBE

Центр окружности (11,444; 7,111), радиус 6,822, нижняя точка - на высоте 0,289. На 0,3 выше дуга ещё есть, на 0,3 ниже - уже нет. Точка (11,444; 7,111) - вершина диаграммы, где сходятся ячейки A, C и D; та же вершина была на рисунке ячеек в прошлом уроке. Совсем рядом с событием дуга уже уже шага 1/4, и выборка её не видит, - поэтому здесь отступ 0,3, а не 0,01.

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

События-окружности для несоседних дуг. Тройка центров даёт событие, только если их дуги соседние на береговой линии в этот момент. Перебор всех троек - это снова n³ и ложные вершины.

Неотменённые ложные тревоги. Если между тройкой появился новый центр или дуги перестали быть соседями, событие надо удалить из очереди. Иначе алгоритм «стянет» дугу, которой уже нет, и рождает вершину внутри чужой ячейки.

Изломы во float. Положение излома - корень квадратного уравнения. Если по нему решать, какую дугу разрезает новый центр, ошибка округления отправит центр не в ту дугу. Решение принимают по точным предикатам - тем же определителям, что в прошлом уроке.

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

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

Теперь сам

1. Центры (0, 10) и (6, 4). Прямая на высоте 2. Сколько дуг у береговой линии и в каком порядке?
2. Центры (0, 0), (4, 0), (2, 2). На какой высоте прямой наступит событие-окружность и где родится вершина?
3. Почему дуг не больше 2n − 1?