Порядок слияний: стек серий, сломанный инвариант и powersort

Если ты знаешь, почему Timsort сливает только соседей, что было не так со старым правилом стека и как powersort считает силу границы, листай до «Теперь сам».

Как думаешь, в каком порядке сливать

Серии длин 1, 200, 200, 2, 2. Слить их надо попарно соседей, цена слияния - сумма длин, ведь каждый элемент переписывается (это модель: обрезку краёв и галоп не учитываем). Самое простое - копить слева: 1 + 200, потом результат с 200 и так далее. Сколько это стоит, и насколько лучше можно? (Длины здесь - модель: в настоящем Timsort при n = 405 короткие серии добились бы вставками до minrun = 51. Порядок слияний считаем для уже готовых длин.)

Ответь себе до чтения. Копить слева - 1410 переписанных элементов. Лучший порядок - 814. Старое правило Timsort - 1211, а нынешнее - 814.

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

Сборка большой пачки документов из стопок, лежащих в ряд на длинном столе. Стопки можно объединять только с соседкой - иначе перепутаются документы с одинаковыми номерами. Объединение стоит столько, сколько листов в обеих стопках. Если таскать одну растущую стопку через весь стол, самые первые листы перекладываются десятки раз. Разумно сначала объединять мелкие стопки между собой, а крупные - в самом конце.

Механизм: стек и сила границы

Проход Timsort кладёт серии на стек по одной и ещё не знает, какие будут дальше. Решать «сливать сейчас или подождать» приходится по тому, что уже видно.

Почему только соседей. Серии A, B, C идут подряд. Слить A с C в обход B нельзя: если в A, B и C есть равные элементы, их взаимный порядок перепутается. Допустимы только (A + B) + C и A + (B + C).

Старое правило (до 3.11, merge_collapse): после каждой новой серии сливать верх стека, пока длины не растут вниз быстрее Фибоначчи - len[-3] > len[-2] + len[-1] и len[-2] > len[-1]. Из этого следует, что стек короткий: длины растут как числа Фибоначчи, и серий на нём не больше логарифма. Размер стека посчитали из инварианта и выделили заранее.

Powersort (с 3.11, found_new_run). Представим массив отрезком [0, 1) и идеальное дерево слияний, которое делит его пополам, четвертями, восьмыми… Граница между соседними сериями получает силу - номер первого уровня, на котором между серединами этих серий проходит линия деления: середина массива - сила 1, четверти - 2, восьмые - 3. Правило: перед тем как положить новую серию, слить на стеке всё, у чего сила больше силы новой границы. Силы на стеке растут снизу вверх, и стек тоже логарифмической глубины.

Правило

Сливать можно только соседние серии - ради устойчивости. Цена слияния - сумма длин. Powersort (Python 3.11+): сила границы - глубина первого деления пополам, разделяющего середины двух соседних серий; перед каждой новой серией сливать верхние серии, пока сила под вершиной больше силы новой границы. Стоимость powersort доказанно не больше n·H + 2n, где n·H = Σ ℓᵢ · log₂(n / ℓᵢ) - нижняя граница любого порядка.

Откуда эта граница. Порядок слияний - двоичное дерево с сериями в листьях: элемент серии переписывается столько раз, какова глубина её листа, и цена - сумма длин, умноженных на глубины. Это задача об оптимальном алфавитном дереве - родственнице кода Хаффмана из главы о жадных алгоритмах, только листья обязаны остаться по порядку (сливать можно лишь соседей). Нижняя граница Σ ℓᵢ log₂(n / ℓᵢ) у них общая: та же, что у средней длины кода.

Разбор: четыре порядка на одних длинах

Повторим оба правила CPython и сравним их с «копить слева» и с лучшим порядком, найденным перебором точки последнего слияния (как в задаче о порядке умножения матриц из раздела про динамическое программирование).

import random

def slit_na(st, k, cena):
    # слить серии k и k + 1 стека; цена слияния - сумма длин
    cena[0] = cena[0] + st[k][1] + st[k + 1][1]
    st[k] = [st[k][0], st[k][1] + st[k + 1][1], st[k][2]]
    st.pop(k + 1)

def dokonchit(st, cena):
    # конец сортировки: сливать сверху, выбирая меньшего соседа
    while len(st) > 1:
        k = len(st) - 2
        if k > 0 and st[k - 1][1] < st[k + 1][1]:
            k = k - 1
        slit_na(st, k, cena)

def staroe_pravilo(serii):
    # Timsort до 3.11: len[-3] > len[-2] + len[-1] и len[-2] > len[-1]
    st = []
    cena = [0]
    for r in serii:
        st.append([0, r, 0])
        while len(st) > 1:
            n = len(st) - 2
            d = []
            for z in st:
                d.append(z[1])
            if ((n > 0 and d[n - 1] <= d[n] + d[n + 1]) or
                    (n > 1 and d[n - 2] <= d[n - 1] + d[n])):
                if d[n - 1] < d[n + 1]:
                    n = n - 1
                slit_na(st, n, cena)
            elif d[n] <= d[n + 1]:
                slit_na(st, n, cena)
            else:
                break
    dokonchit(st, cena)
    return cena[0]

def sila(s1, n1, n2, n):
    # глубина узла между сериями в идеальном дереве слияний
    # удвоенные середины серий: так деление пополам остаётся в целых
    a = 2 * s1 + n1
    b = a + n1 + n2
    r = 0
    while True:
        r = r + 1
        # a/n и b/n - середины в долях массива; берём их двоичные
        # разряды по одному: >= n - разряд 1 (вычитаем n), иначе 0
        if a >= n:
            a = a - n
            b = b - n
        elif b >= n:
            return r
        a = a * 2
        b = b * 2

def powersort(serii):
    n = sum(serii)
    st = []                          # [начало, длина, сила]
    cena = [0]
    s = 0
    for r in serii:
        if len(st) > 0:
            p = sila(st[-1][0], st[-1][1], r, n)
            while len(st) > 1 and st[-2][2] > p:
                slit_na(st, len(st) - 2, cena)
            st[-1][2] = p
        st.append([s, r, 0])
        s = s + r
    dokonchit(st, cena)
    return cena[0]

def optimum(serii):
    # лучший порядок слияний соседей: перебор точки последнего слияния
    k = len(serii)
    pref = [0]
    for x in serii:
        pref.append(pref[-1] + x)
    D = []
    for i in range(k):
        D.append([0] * k)
    for dl in range(2, k + 1):
        for i in range(k - dl + 1):
            j = i + dl - 1
            luchshe = -1
            for m in range(i, j):
                c = D[i][m] + D[m + 1][j]
                if luchshe == -1 or c < luchshe:
                    luchshe = c
            D[i][j] = luchshe + pref[j + 1] - pref[i]
    return D[0][k - 1]

def podryad(serii):
    cena = 0
    nakopleno = serii[0]
    for x in serii[1:]:
        nakopleno = nakopleno + x
        cena = cena + nakopleno
    return cena

for serii in [[1, 200, 200, 2, 2], [200, 5, 3, 200],
              [64, 64, 64, 64, 64, 64, 64, 64]]:
    print(serii)
    print("  подряд", podryad(serii), "старое", staroe_pravilo(serii),
          "powersort", powersort(serii), "оптимум", optimum(serii))

random.seed(11)
hudshee = [0, 0]
for t in range(150):
    serii = []
    for i in range(random.randint(3, 8)):
        serii.append(random.choice([1, 3, 8, 20, 50, 120, 300]))
    o = optimum(serii)
    hudshee = [max(hudshee[0], staroe_pravilo(serii) / o),
               max(hudshee[1], powersort(serii) / o)]
print("худшее отношение к оптимуму на 150 наборах: старое",
      round(hudshee[0], 2), "powersort", round(hudshee[1], 2))

[1, 200, 200, 2, 2]
  подряд 1410 старое 1211 powersort 814 оптимум 814
[200, 5, 3, 200]
  подряд 821 старое 624 powersort 816 оптимум 624
[64, 64, 64, 64, 64, 64, 64, 64]
  подряд 2240 старое 1536 powersort 1536 оптимум 1536
худшее отношение к оптимуму на 150 наборах: старое 1.33 powersort 1.17

Восемь равных серий - так выглядят случайные данные после добивания до minrun - оба правила сливают идеальным деревом: 1536, каждый элемент переписан ровно три раза. «Копить слева» - 2240.

На [1, 200, 200, 2, 2] старое правило дорого: оно слило единицу с первой двухсоткой (201), вторую двухсотку с двойкой (202), потом их вместе (403) и в конце с последней двойкой (405) - каждая двухсотка переписана трижды. В оптимуме сначала сливается мелочь справа (2 + 2, потом 200 + 4), затем 1 + 200, и большие куски встречаются один раз, в самом конце: каждая двухсотка переписана дважды. Powersort видит, что граница двухсоток - середина массива (сила 1), и откладывает это слияние на самый конец: ровно оптимум 814.

Разбор: и powersort не всегда лучший

Второй набор - [200, 5, 3, 200]. Здесь старое правило попало в оптимум 624, а powersort - 816, на 31% хуже. Середина массива проходит между 5 и 3, поэтому граница 5|3 получила силу 1 и сливается последней, а мелочь прилипает к большим сериям по отдельности. Это граница правила: powersort не оптимален на каждом входе, его гарантия - про худший случай. Оптимальный порядок требует знать все длины заранее (простая динамика из разбора - O(r³), алгоритм Ху - Таккера - O(r log r)), а Timsort видит серии по одной. Зато для старого правила известны входы, где оно примерно в полтора раза хуже n·H, а powersort выше n·H + 2n не поднимается никогда. На 150 случайных наборах худшее отношение к оптимуму: 1,33 у старого правила против 1,17 у powersort.

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

Инвариант, проверенный не везде. До 2015 года старое правило проверяло условие только для трёх верхних серий (модель урока - уже исправленная версия из Python 3.10). В 2015 году де Гау, Рот, де Бур, Бубель и Хенле пытались формально доказать Timsort из Java в системе KeY - и нашли вход из множества коротких серий, где инвариант глубже по стеку нарушается. Стек, размер которого посчитали из инварианта, переполнялся: java.util.Arrays.sort падала с ArrayIndexOutOfBoundsException на массиве из 67 108 864 элементов (JDK-8072909). В Python исправили саму проверку - смотреть и четвёртую сверху серию: в коде 3.10 это второе условие в merge_collapse, n > 1 and len[n-2] <= len[n-1] + len[n]. В OpenJDK вместо этого увеличили запас стека.

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

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

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

Теперь сам

1. Серии длин 4, 4, 4, 4. Какова цена «копить слева» и цена идеального порядка?
2. Массив длины 100, в начале серии длин 30 и 20 (позиции 0-29 и 30-49). Какая сила у границы между ними?
3. Почему слить сначала A и C в обход B - ошибка, даже если это дешевле?