Что делают сканом: упаковка, разбиение по биту и ранг в списке

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

Как думаешь, кто куда пишет

Восемь исполнителей держат по числу: 13, 2, 7, 10, 4, 15, 0, 9. Нужно собрать нечётные числа в новый плотный массив, сохранив порядок: 13, 7, 15, 9. Все пишут одновременно. Откуда исполнитель с числом 15 знает, что ему писать в ячейку 2, если он не видит, сколько нечётных было слева?

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

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

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

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

Механизм: позиция из скана

Упаковка (stream compaction): флаг f[i] = 1, если элемент остаётся. Исключающий скан флагов s - сколько единиц левее i. Элемент с флагом пишется в y[s[i]]. Размер ответа - сумма всех флагов, её скан Блеллоха получает даром после подъёма.

Разбиение по биту: нули идут в начало в своём порядке, единицы - за ними. Для нулей скан флагов «бит равен 0» даёт позицию, для единиц - число всех нулей плюс скан флагов «бит равен 1». Порядок внутри групп сохраняется, то есть разбиение устойчиво. Разбиение по битам 0, 1, 2, … от младшего к старшему - и есть поразрядная сортировка: каждый проход устойчив, поэтому порядок по младшим битам не разрушается.

Ранг в списке. Связный список лежит в массиве вперемешку: sled[i] - номер следующего элемента, у последнего он указывает сам на себя. Ранг - расстояние до конца. Прыжки по указателям (Уилли, 1979): у каждого элемента rang сначала 1 (у последнего 0), и в каждом раунде все одновременно делают rang[i] += rang[sled[i]], sled[i] = sled[sled[i]]. После раунда k указатель смотрит на 2ᵏ элементов вперёд, а ранг накопил их число - охват удваивается, как в Хиллисе - Стиле.

Правило

Чтобы много исполнителей записали результаты без очереди, сначала сканом вычисляют позицию каждого, потом пишут все разом. Упаковка - скан флагов; устойчивое разбиение по биту - два скана, и это шаг поразрядной сортировки; ранги в связном списке - прыжками по указателям за ⌈log₂ n⌉ раундов при работе n · log₂ n.

Разбор: упаковка и поразрядная сортировка

Скан здесь написан циклом, чтобы не загромождать программу; в параллельной версии на его месте стоит скан Блеллоха из прошлого урока, а записи элементов в заранее вычисленные ячейки независимы и идут одновременно. Вспомогательные списки программа строит последовательно, а проходы по битам идут один за другим: каждый опирается на результат предыдущего. Движок курса не знает побитовых операций, поэтому бит b числа v берётся арифметикой: v // 2 ** b % 2 - деление на 2ᵇ отбрасывает b младших битов, а остаток от деления на 2 оставляет следующий.

def skan_isklyuch(f):
    # исключающий скан; в параллельной версии - Блеллох из прошлого урока
    s = []
    acc = 0
    for v in f:
        s.append(acc)
        acc = acc + v
    return [s, acc]

def upakovat(x, flag):
    poz, vsego = skan_isklyuch(flag)
    y = [0] * vsego
    for i in range(len(x)):
        if flag[i] == 1:
            y[poz[i]] = x[i]
    return y

def razbit_po_bitu(x, b):
    # нули по порядку в начало, единицы по порядку за ними
    bit = []
    for v in x:
        bit.append(v // 2 ** b % 2)
    nul = []
    for v in bit:
        nul.append(1 - v)
    poz0, nulej = skan_isklyuch(nul)
    poz1, edinic = skan_isklyuch(bit)
    y = [0] * len(x)
    for i in range(len(x)):
        if bit[i] == 0:
            y[poz0[i]] = x[i]
        else:
            y[nulej + poz1[i]] = x[i]
    return y

x = [13, 2, 7, 10, 4, 15, 0, 9]
flag = []
for v in x:
    if v % 2 == 1:
        flag.append(1)
    else:
        flag.append(0)
print("нечётные:", upakovat(x, flag))
for b in range(4):
    x = razbit_po_bitu(x, b)
    print("бит", b, x)

нечётные: [13, 7, 15, 9]
бит 0 [2, 10, 4, 0, 13, 7, 15, 9]
бит 1 [4, 0, 13, 9, 2, 10, 7, 15]
бит 2 [0, 9, 2, 10, 4, 13, 7, 15]
бит 3 [0, 2, 4, 7, 9, 10, 13, 15]

После разбиения по биту 0 все чётные стоят впереди, причём в исходном порядке: 2, 10, 4, 0. После бита 1 впереди те, у кого бит 1 нулевой: 4 (100), 0 (000), 13 (1101), 9 (1001) - и внутри этой группы порядок по биту 0 сохранился. Четыре прохода по четырёхбитным числам - массив отсортирован. Каждый проход - два скана и одна волна записей: глубина O(log n) на бит, работа O(n).

Разбор: ранг в перемешанном списке

Список 5 → 2 → 7 → 0 → 3 → 6 → 1 → 4 лежит в массиве так, что соседи по списку не соседи по памяти. Пройти его по указателям - n зависимых шагов. Прыжки делают это за логарифм раундов.

import random

def rang_pryzhkami(sled):
    # sled[i] - следующий элемент списка, у последнего sled[i] = i
    n = len(sled)
    rang = []
    for i in range(n):
        if sled[i] == i:
            rang.append(0)
        else:
            rang.append(1)
    sled = sled[:]
    raundov = 0
    operacij = 0
    gotovo = False
    while not gotovo:
        gotovo = True
        novyj_rang = rang[:]
        novyj_sled = sled[:]
        for i in range(n):
            # указатель ещё не упёрся в хвост - нужен ещё раунд
            if sled[i] != sled[sled[i]]:
                gotovo = False
            novyj_rang[i] = rang[i] + rang[sled[i]]
            novyj_sled[i] = sled[sled[i]]
            operacij = operacij + 1
        if gotovo:
            break
        rang = novyj_rang
        sled = novyj_sled
        raundov = raundov + 1
    return [rang, raundov, operacij]

# список 5 -> 2 -> 7 -> 0 -> 3 -> 6 -> 1 -> 4, в массиве вперемешку
sled = [3, 4, 7, 6, 4, 2, 1, 0]
print(rang_pryzhkami(sled))
random.seed(3)
poryadok = list(range(1000))
random.shuffle(poryadok)
sled = [0] * 1000
for k in range(999):
    sled[poryadok[k]] = poryadok[k + 1]
sled[poryadok[999]] = poryadok[999]
r = rang_pryzhkami(sled)
print(r[0][poryadok[0]], r[0][poryadok[500]], r[1], r[2])

[[4, 1, 6, 3, 0, 7, 2, 5], 3, 32]
999 499 10 11000

Ранги в выводе стоят по номерам элементов массива, а не по порядку списка: ранг элемента 0 - 4, элемента 5 (головы) - 7, элемента 4 (хвоста) - 0. Раунды кончаются, когда у каждого элемента указатель и указатель соседа совпадают: все смотрят на хвост, и ранги готовы. Три раунда на восьми элементах, десять на тысяче: ⌈log₂ 1000⌉ = 10. Операций 11 000 - десять раундов по тысяче плюс проход, на котором программа видит, что все указатели уже упёрлись в хвост. Работа n · log₂ n против n у обхода одним исполнителем - прыжки неработоэффективны. Работоэффективные алгоритмы ранжирования есть (сначала выбрасывают из списка независимое множество элементов, ранжируют остаток и возвращают выброшенные), но они заметно сложнее.

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

Общий счётчик вместо скана. poz = atomic_fetch_add(counter, 1) - неделимое «вернуть старое значение и прибавить» (подробно - в главе о разделяемой памяти) - тоже раздаёт уникальные адреса, но в порядке прихода, а не в порядке массива: упаковка перестаёт быть устойчивой, и поразрядная сортировка на ней ломается. Плюс счётчик - одна горячая ячейка, куда пишут все.

Прыжки на месте. Если обновлять sled[i] в том же массиве, из которого читают соседи, часть элементов в раунде прыгнет дальше, чем на 2ᵏ, часть - меньше, и ранги перепутаются. Как и в Хиллисе - Стиле, раунд читает старые значения и пишет новые.

Неустойчивое разбиение в поразрядной сортировке. Если внутри группы порядок не сохраняется, проход по биту 1 разрушает порядок по биту 0 - и после четырёх проходов массив не отсортирован.

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

Фильтр на видеокарте (выбрать строки таблицы по условию, частицы внутри области, рёбра для следующего фронта поиска в ширину) - упаковка сканом. Поразрядная сортировка через разбиение - основная сортировка коротких ключей на GPU (cub::DeviceRadixSort), а в базах данных тот же приём раскладывает строки по секциям при параллельном хеш-соединении. Прыжки по указателям нужны, когда данные - связные структуры: обход Эйлера по дереву, связные компоненты, корни леса указателей.

Теперь сам

1. Флаги 0, 1, 1, 0, 1, 0, 0, 1. Куда попадёт элемент с номером 4 при упаковке, и какого размера ответ?
2. Разбей по биту 0 массив 6, 3, 8, 5, 1, 4 и выпиши позиции каждого элемента.
3. Сколько раундов прыжков нужно списку из 5000 элементов? А если список - уже по порядку в массиве?