Дорожная сеть: иерархии сокращений вместо поиска по всей карте

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

Как думаешь, сколько карты нужно, чтобы проложить маршрут

Из Москвы в Казань. Дейкстра расходится от старта кругами и, пока не дойдёт до Казани, просмотрит всё в радиусе 800 километров - каждый двор и каждый проезд. A* сжимает круг в вытянутый эллипс, но всё равно идёт по дворам.

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

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

Ветки и магистрали метро. Чтобы попасть с окраины на окраину, едешь сначала к центру, где линии крупнее, а потом расходишься. Маленькие станции между пересадками для маршрута не важны: проезжаешь насквозь. Иерархия сокращений строит такое метро для любой дорожной сети: «проезжать насквозь маленькие станции» - это и есть сокращения, а «сначала к центру, потом от центра» - запрос вверх по важности с двух сторон.

Механизм: сжатие со свидетелем и запрос вверх

ПРЕДОБРАБОТКА
важность(v) ≈ (сколько сокращений понадобится) − (степень v)
              + (сколько соседей уже сжато)
пока есть вершины: взять наименее важную v (пересчитав - важность
    меняется по ходу), и для каждой пары её соседей u, w:
        путь через v: c(u,v) + c(v,w)
        «свидетель» - Дейкстра от u до w БЕЗ v, не дальше этой длины
        свидетеля нет → сокращение u–w с весом c(u,v) + c(v,w)
    v получает очередной ранг и убирается из рабочего графа

ЗАПРОС s → t
    поиск от s только по рёбрам к вершинам БОЛЬШЕГО ранга
    поиск от t - так же
    ответ = min по общим вершинам x: dist_s(x) + dist_t(x)

Почему запрос верен. Возьмём кратчайший путь s → t и самую важную вершину на нём x. Любая менее важная вершина v пути была сжата раньше соседей по пути - и в этот момент либо добавилось сокращение в обход v, либо нашёлся свидетель не длиннее. Убирая так все «впадины», получаем путь той же длины, у которого ранги растут от s до x и падают от x до t. Подъём от s найдёт x, подъём от t - тоже, и сумма в x будет кратчайшей.

Правило

Иерархия сокращений сжимает вершины от неважных к важным; при сжатии v для каждой пары соседей ищут свидетеля - путь в обход v не длиннее пути через v - и, если его нет, добавляют сокращение. Важность - сколько сокращений придётся добавить минус степень плюс число сжатых соседей; её пересчитывают лениво. Запрос - две Дейкстры только к вершинам большего ранга и минимум суммы в общей вершине. Предобработка - один раз на карту, запрос - доли от обычного поиска.

Разбор: квадрат и сокращение

    a ---2--- b
    |         |
    1         3
    |         |
    d ---1--- c

сжимаем b (соседи a и c): через b путь a–c стоит 2 + 3 = 5
    свидетель: a → d → c = 1 + 1 = 2 ≤ 5   → сокращение НЕ нужно
сжимаем d (соседи a и c): через d путь a–c = 2
    свидетель без d: a → b → c = 5 > 2, но b уже сжата и убрана
    → сокращение a–c с весом 2

осталось a и c с ребром 2 - кратчайший путь a–c = 2 сохранён

Сжатие b прошло даром: дорога через неё длиннее объезда. Сжатие d потребовало сокращения: она лежала на кратчайшем пути, и её «проезд насквозь» надо было запомнить.

Разбор: сетки 40 × 40 и 80 × 80

Программа строит сетку улиц со случайной длиной от 1 до 10, сжимает её с ленивой важностью и выполняет по 300 запросов, сверяя каждый ответ с обычной Дейкстрой и считая просмотренные вершины.

import heapq
import random

INF = float("inf")


def poisk(g, s, t=None, bez=None, predel=INF, schet=None):
    dist, h = {s: 0}, [(0, s)]        # Дейкстра: до t, без вершины bez,
    while h:                          # не дальше предела
        d, x = heapq.heappop(h)
        if d > predel:
            break
        if d > dist[x]:
            continue
        if schet is not None:
            schet[0] += 1
        if x == t:
            break
        for y, c in g[x].items():
            if y != bez and d + c < dist.get(y, INF):
                dist[y] = d + c
                heapq.heappush(h, (d + c, y))
    return dist


def opyt(N, rng):
    V = N * N                         # сетка улиц N × N
    graf = [dict() for _ in range(V)]
    for v in range(V):
        for u in (v + 1 if v % N + 1 < N else -1, v + N):
            if 0 <= u < V:
                graf[v][u] = graf[u][v] = rng.randint(1, 10)

    def nuzhnye(g, v):                # сокращения, нужные при сжатии v
        s, otvet = list(g[v].items()), []
        for i, (u, cu) in enumerate(s):
            for w, cw in s[i + 1:]:   # есть ли обход u→w без v?
                obhod = poisk(g, u, w, v, cu + cw).get(w, INF)
                if obhod > cu + cw:
                    otvet.append((u, w, cu + cw))
        return otvet

    g = [dict(x) for x in graf]       # рабочий граф: ещё не сжатые
    vse = [dict(x) for x in graf]     # исходные рёбра + сокращения
    rang, sjato = [0] * V, [0] * V
    vazhn = lambda v: len(nuzhnye(g, v)) - len(g[v]) + sjato[v]
    och = [(vazhn(v), v) for v in range(V)]
    heapq.heapify(och)
    k = sokr = 0
    while och:                        # сжимаем «неважные» первыми
        p, v = heapq.heappop(och)
        tek = vazhn(v)
        if och and tek > och[0][0]:   # важность устарела - назад
            heapq.heappush(och, (tek, v))
            continue
        for u, w, c in nuzhnye(g, v):
            if c < g[u].get(w, INF):
                g[u][w] = g[w][u] = vse[u][w] = vse[w][u] = c
                sokr += 1
        for u in g[v]:
            del g[u][v]
            sjato[u] += 1
        g[v], rang[v], k = {}, k, k + 1
    vverh = [{u: c for u, c in vse[v].items() if rang[u] > rang[v]}
             for v in range(V)]

    osm_d, osm_ch = [0], [0]
    for _ in range(300):
        s, t = rng.randrange(V), rng.randrange(V)
        d = poisk(graf, s, t, schet=osm_d)[t]
        a = poisk(vverh, s, schet=osm_ch)     # оба поиска - только
        b = poisk(vverh, t, schet=osm_ch)     # к более важным
        assert min(a[x] + b[x] for x in a if x in b) == d
    print(f"{V:5} перекрёстков: сокращений {sokr:5}")
    print(f"      вершин на запрос: Дейкстра {osm_d[0] // 300:4},"
          f" иерархия {osm_ch[0] // 300}")


rng = random.Random(16)
for N in (40, 80):
    opyt(N, rng)

 1600 перекрёстков: сокращений  3389
      вершин на запрос: Дейкстра  785, иерархия 120
 6400 перекрёстков: сокращений 15425
      вершин на запрос: Дейкстра 3323, иерархия 272

Все 600 ответов совпали с Дейкстрой. При росте карты вчетверо Дейкстра стала просматривать вчетверо больше вершин, иерархия - чуть больше чем вдвое; на 1600 перекрёстках выигрыш шестикратный, на 6400 - двенадцатикратный. На настоящих дорогах выигрыш ещё больше: там иерархия естественная (дворы, улицы, трассы), и запрос по стране просматривает сотни вершин из десятков миллионов. Цена - сокращения: рёбер стало больше чем вдвое против исходной сетки, и предобработка, которая в Python заняла секунды.

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

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

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

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

Теперь сам

Путь a – b – c с весами 4 и 5, других путей между a и c нет. Ранги: b - 1, a - 2, c - 3. Какие сокращения появятся? Как пройдёт запрос a → c?