Слияние с галопом: когда прыгать выгоднее, чем шагать

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

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

Две отсортированные серии по сто тысяч чисел. В первой - блоки 0-99, 200-299, 400-499…, во второй - 100-199, 300-399… Обычное слияние сравнивает голову с головой и делает около 2n сравнений: здесь около двухсот тысяч. Сколько сделает sorted, получив a + b?

Ответь себе до чтения. Если «двести тысяч» - разбор покажет 27 992 сравнения на слияние: в семь раз меньше.

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

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

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

Механизм: обрезка краёв, по одному, галоп

Слияние соседних серий A (левая) и B (правая) в merge_at:

1. Обрезка краёв. Элементы A, не большие B[0], уже стоят на месте - их ищут галопом и пропускают. Так же отбрасывают хвост B, не меньший последнего элемента A. Если вся A не больше B[0], слияния нет вовсе; если вся A больше всей B, слияние сводится к двум галопам и копированию кусков. Меньшая из оставшихся частей копируется во временный буфер: памяти нужно min(|A|, |B|), а не n. Этого хватает: если в буфере A, результат пишется слева направо на её прежнее место, и свободных ячеек впереди всегда столько, сколько элементов ещё лежит в буфере, - запись не догонит непрочитанную часть B (если в буфере B - то же самое справа налево).
2. По одному. Сравнить B[j] < A[i]: при равенстве берётся элемент A - это и даёт устойчивость. Считать, сколько раз подряд выиграла одна серия.
3. Галоп - когда счётчик дошёл до min_gallop (вначале 7). По очереди: найти в A все элементы, не большие B[j], и забрать их куском; забрать B[j]; найти в B все элементы, меньшие A[i], - забрать куском; забрать A[i]. Пока хоть один кусок длиной от 7 - продолжать. Каждый круг галопа уменьшает min_gallop на 1, выход из галопа увеличивает.

Экспоненциальный поиск позиции x начиная с L[s]: сравнить с L[s], L[s+1], L[s+3], L[s+7]… до первого «перескочили», потом двоичный поиск в последнем промежутке.

Правило

Слияние в Timsort идёт по одному, пока одна серия не выиграет min_gallop раз подряд (вначале 7); тогда - галоп: кусок длины k ≥ 1 находится экспоненциальным поиском за 2⌊log₂ k⌋ + 2 сравнения вместо k + 1 (если серия кончилась раньше, сравнений меньше). Окупившийся галоп понижает порог, неокупившийся - повышает. Перед слиянием края, уже стоящие на месте, обрезаются, и буфер нужен на меньшую из частей.

Разбор: модель слияния на четырёх входах

Повторим merge_lo из CPython на модели и посчитаем сравнения: по одному и с галопом. Серии по 400 элементов.

import random

schet = [0]

def menshe(x, y):
    schet[0] = schet[0] + 1
    return x < y

def galop(L, s, x, pravyj):
    # сколько элементов L[s:] идут раньше x: правый - те, что <= x,
    # левый - те, что < x. Прыжки 1, 3, 7, 15, ..., потом двоичный поиск
    if pravyj:
        ok = not menshe(x, L[s])
    else:
        ok = menshe(L[s], x)
    if not ok:
        return 0
    prosh = 0
    ofs = 1
    while s + ofs < len(L):
        if pravyj:
            ok = not menshe(x, L[s + ofs])
        else:
            ok = menshe(L[s + ofs], x)
        if not ok:
            break
        prosh = ofs
        ofs = 2 * ofs + 1
    lo = prosh + 1
    hi = min(ofs, len(L) - s)
    while lo < hi:
        mid = (lo + hi) // 2
        if pravyj:
            ok = not menshe(x, L[s + mid])
        else:
            ok = menshe(L[s + mid], x)
        if ok:
            lo = mid + 1
        else:
            hi = mid
    return lo

def slit(A, B, s_galopom):
    schet[0] = 0
    i = 0
    j = 0
    if s_galopom:            # обрезка краёв, как в merge_at
        i = galop(A, 0, B[0], True)
        if i == len(A):
            return schet[0]  # A целиком не больше B[0]: сливать нечего
        konec_b = galop(B, 0, A[len(A) - 1], False)
    else:
        konec_b = len(B)
    min_galop = 7
    while i < len(A) and j < konec_b:
        pobed_a = 0
        pobed_b = 0
        while i < len(A) and j < konec_b:          # по одному
            if menshe(B[j], A[i]):
                j = j + 1
                pobed_b = pobed_b + 1
                pobed_a = 0
            else:
                i = i + 1
                pobed_a = pobed_a + 1
                pobed_b = 0
            if s_galopom and max(pobed_a, pobed_b) >= min_galop:
                break
        while i < len(A) and j < konec_b:          # галоп
            if min_galop > 1:
                min_galop = min_galop - 1
            ka = galop(A, i, B[j], True)
            i = i + ka
            if i == len(A):
                break
            j = j + 1                              # сам B[j]
            if j == konec_b:
                break
            kb = galop(B, j, A[i], False)
            j = j + kb
            if j == konec_b:
                break
            i = i + 1                              # сам A[i]
            if ka < 7 and kb < 7:
                min_galop = min_galop + 2          # штраф: не окупился
                break
    return schet[0]

random.seed(2)
n = 400
vse = list(range(2 * n))
random.shuffle(vse)
chet = list(range(0, 2 * n, 2))
nechet = list(range(1, 2 * n, 2))
nabory = [["вперемешку", sorted(vse[:n]), sorted(vse[n:])],
          ["через один", chet, nechet]]
ryad = list(range(2 * n))
A = []
B = []
for k in range(2 * n):
    if (k // 50) % 2 == 0:
        A.append(k)
    else:
        B.append(k)
nabory.append(["блоками по 50", A, B])
nabory.append(["B целиком больше", ryad[:n], ryad[n:]])
for imya, A, B in nabory:
    print(imya, "- по одному:", slit(A, B, False),
          "с галопом:", slit(A, B, True))

вперемешку - по одному: 799 с галопом: 814
через один - по одному: 799 с галопом: 815
блоками по 50 - по одному: 750 с галопом: 192
B целиком больше - по одному: 400 с галопом: 16

На перемешанных и чередующихся сериях галоп почти никогда не включается, и проигрыш - 2%: лишние сравнения на обрезку краёв и редкие неудачные заходы. На блоках по 50 галоп забирает каждый блок за десяток сравнений вместо пятидесяти: 192 против 750. На непересекающихся сериях хватает обрезки краёв: 16 сравнений вместо 400, а слияния нет вовсе.

Разбор: где галоп проигрывает

Почему не галопировать всегда? Посчитаем сравнения, если кусок, который можно забрать, имеет длину k. Обычный путь - k + 1 сравнение: k «да» и одно «нет». Галоп: сравнения с позициями 0, 1, 3, 7, … до первого «нет», потом двоичный поиск в промежутке.

длина куска k:      0  1  2  3  4  5  6  7  8  15  100
по одному:          1  2  3  4  5  6  7  8  9  16  101
галоп:              1  2  4  4  6  6  6  6  8   8   14

До k = 5 галоп нигде не выигрывает, а при k = 2 и k = 4 проигрывает: 4 сравнения вместо 3 - на треть больше, 6 вместо 5 - на пятую. Дальше выигрыш неограничен: кусок из ста - 14 сравнений вместо 101. Угадать длину куска заранее нельзя, поэтому Timsort ждёт улики - семь побед подряд, - а потом подстраивает порог: на случайных данных он растёт, и галоп почти не включается; на блочных падает до единицы.

Тот же эффект на настоящем sorted: две серии по сто тысяч, разбитые на блоки разной длины. Обёртка считает сравнения; 2n − 1 из них уходит на поиск двух серий, остальное - на слияние.

class Klyuch:
    schet = 0

    def __init__(self, x):
        self.x = x

    def __lt__(self, drugoj):
        Klyuch.schet += 1
        return self.x < drugoj.x


n = 100000
for blok in [1, 10, 100, 1000, n]:
    # блоки с чётным номером - в первую серию, с нечётным - во вторую
    a = [k for k in range(2 * n) if (k // blok) % 2 == 0]
    b = [k for k in range(2 * n) if (k // blok) % 2 == 1]
    Klyuch.schet = 0
    sorted(a + b, key=Klyuch)
    print("блоки по", blok, "- сравнений", Klyuch.schet,
          "из них на слияние", Klyuch.schet - (2 * n - 1))

блоки по 1 - сравнений 399998 из них на слияние 199999
блоки по 10 - сравнений 359992 из них на слияние 159993
блоки по 100 - сравнений 227991 из них на слияние 27992
блоки по 1000 - сравнений 203985 из них на слияние 3986
блоки по 100000 - сравнений 199999 из них на слияние 0

Через один - 2n сравнений на слияние, галоп не помогает и почти не мешает. Блоки по 100 - в семь раз меньше, блоки по 1000 - в пятьдесят. Блоки во всю длину - это уже одна серия, слияния нет.

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

Слияние «как в учебнике» на блочных данных. Склейка двух журналов, где события идут пачками по секундам, или объединение отсортированных выгрузок с разных серверов: голое слияние делает 2n сравнений, Timsort - в разы меньше. Самописное слияние на Python здесь проигрывает sorted(a + b) и по сравнениям, и по скорости интерпретатора.

Буфер на весь массив. Учебное слияние копирует обе половины. Timsort копирует только меньшую часть после обрезки краёв - на почти упорядоченных данных это килобайты, а не мегабайты.

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

Экспоненциальный поиск - самостоятельный инструмент: поиск в отсортированном потоке неизвестной длины, пересечение отсортированных списков в поисковых индексах (Тим Петерс взял галоп из работы Демейна, Лопес-Ортиса и Манро о пересечении множеств, а сама идея восходит к Бентли и Яо, 1976), поиск в журнале по времени, когда искомое скорее всего близко к началу. Везде, где ответ «обычно рядом, но иногда далеко», прыжки 1, 2, 4, 8 дают O(log k) вместо O(log n) - цена зависит от расстояния до ответа, а не от размера данных.

Теперь сам

1. Сколько сравнений сделает галоп, если забрать можно ровно 31 элемент, а за ними в серии есть ещё?
2. Серии A = [1, 2, 3, 10, 11] и B = [4, 5, 6, 7]. Какие элементы останутся после обрезки краёв?
3. Почему при равенстве B[j] и A[i] слияние берёт элемент A?