Интерполяционный поиск

Когда в словаре ищут «яблоко», его не открывают посередине — открывают у самого конца, потому что «я» — последняя буква. Двоичный поиск так не умеет: он делит отрезок пополам, не глядя на то, что ищет. Интерполяционный поиск смотрит: прикидывает место по значениям на краях отрезка и пробует там. Ниже пять отсортированных рядов из 45 227 городов мира: хеши названий, сами названия как в телефонной книге, широты, население и придуманный ряд с «великаном». Окно сцены следит за отрезком, который остался, и хорда между его краями всегда идёт по диагонали окна.

На хешах цель находится за три-четыре пробы, а двоичному поиску на тех же 42 311 числах нужно четырнадцать-шестнадцать. Хеш перемешивает биты, значения ложатся ровно, и кривая в любом окне — почти диагональ. На населении поиск ползёт: город в 22 026 жителей он находит за 2139 проб. Кривая населения — длинная полка из десятков тысяч городков и обрыв вверх, где миллионники. Хорда идёт высоко над полкой и показывает почти на левый край, а при приближении картина повторяется.

Идея

Интерполяционный поиск устроен как двоичный: проба, сравнение, отрезок сужается. Отличается выбор места. Двоичный поиск берёт середину, а интерполяционный проводит прямую через края отрезка — точки (lo, a[lo]) и (hi, a[hi]) — и пробует там, где она достигает высоты цели x.

Цель стоит выше нижнего края на долю (x − a[lo]) / (a[hi] − a[lo]) всего подъёма. Если значения растут равномерно, на такую же долю отрезка она и отстоит от его начала:

доля = (x − a[lo]) / (a[hi] − a[lo]),
pos = lo + ⌊доля · (hi − lo)⌋.

Дальше всё как у двоичного поиска. Если на месте pos цель — готово; если там меньше — ищем правее pos, если больше — левее. Поиск идёт, пока цель лежит между значениями на краях: при x < a[lo] или x > a[hi] её в отрезке нет. Это условие нужно не только для того, чтобы вовремя остановиться.

Без проверки краёв доля выходит за [0, 1], и проба улетает за пределы массива. Когда a[lo] = a[hi], делить нельзя, но тогда все значения отрезка равны и пробовать можно любое место — листинг берёт lo.

Почему так мало проб

Возьмём модель, в которой метод хорош. Пусть n ключей — независимые случайные числа, равномерно разбросанные по отрезку [A, B], и края A и B известны. Ищем число x; пусть K — сколько ключей не больше x, то есть место, куда x встанет. Первая проба угадывает его как np, где p = (x − A)/(B − A).

На миллионе ключей первая проба в среднем промахивается не больше чем на 500 мест — а двоичный поиск после первой пробы знает лишь, в какой половине цель. Дальше та же картина в меньшем масштабе, и тут уже прикидка, а не доказательство. После пробы новая хорда начинается от значения, которое близко к x, и доля p в новой задаче маленькая — порядка ошибки, делённой на длину отрезка m. Дисперсия mp(1 − p) тогда порядка прошлой ошибки, и новая ошибка — порядка корня из прошлой. Отрезок, где ещё может быть цель, сжимается как n, √n, n^1/4, … и после k проб он порядка n^1/2k. До нескольких мест это доходит, когда 2^k ≈ log₂ n, то есть при k ≈ log₂ log₂ n. Для миллиона ключей это 4,3 пробы, для миллиарда — 4,9.

Строго это доказали в 1978 году Иехошуа Перл, Алон Итаи и Хаим Авни: в среднем хватает около log₂ log₂ n проб, и сделать заметно больше почти невероятно. На хешах из сцены log₂ log₂ 42 311 ≈ 3,94, а листинг ниже намерил в среднем 4,0.

Где он ломается

Всё держится на том, что кривая значений близка к прямой. Если это не так, хорда врёт, и можно построить ряд, где она врёт каждый раз.

В сцене это ряд «великан»: 1 … 4999 и 10¹⁸. Окно не сужается — слева отщипывается по одному месту, и великан в углу всё время тянет хорду вниз. Настоящие данные ломают поиск мягче, но по той же причине.

• Население. Медианный город из 45 227 — 25 247 жителей, у 73 % городов меньше 50 000, а самый большой — почти 25 миллионов. Кривая — полка с обрывом, и обрыв играет роль великана: в среднем 1587 проб, в худшем 2211.
• Телефонная книга. Название — число из первых четырёх букв по основанию 27, как если бы книгу открывали по алфавиту. Но буквы неравны: на S начинается 4898 названий, на X — 211. Интерполяция делает в среднем 17,6 пробы, в худшем 60, а двоичный поиск — 13,1.
• Широты. Почти две трети городов лежат между 20° и 60° северной широты, к полюсам их почти нет. Хорда, проведённая через крайние широты, промахивается — 16,0 пробы против 14,0.

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

Страховка

Лечат это просто — чередованием. Нечётная проба — по хорде, чётная — в середину отрезка. Половинные шаги гарантируют, что отрезок сжимается, а хорды сохраняют удачу там, где кривая ровная.

Худший случай — вдвое хуже двоичного поиска, зато в среднем на наших рядах выходит так:

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

Кривая из прямых кусков

Разрежем ряд на куски так, чтобы в каждом прямая через первый ключ угадывала место любого ключа с ошибкой не больше ε. Кусок берёт ключи по одному, пока такая прямая есть: каждый ключ разрешает свой отрезок наклонов, и их пересечение сужается, как конус с вершиной в первом ключе. Пересечение опустело — начинается новый кусок. Индекс хранит начала кусков и наклоны. Поиск находит кусок двоичным поиском по началам, прикидывает место по его прямой и ищет цель среди 2ε + 1 соседних мест.

Население городов при ε = 64 укладывается в 18 кусков. Закладок через каждые 129 ключей, которые дали бы то же окно, понадобилось бы 238. Хешам хватает шести кусков вместо 328, широтам — 25 вместо 252. Проб выходит примерно столько же, сколько у двоичного поиска, — выигрыш в памяти: у населения индекс меньше в 13 раз, у хешей — в 55, а последние пробы идут внутри 129 соседних ключей — для 64-битных ключей это один килобайт памяти.

Так устроены выученные индексы баз данных. FITing-tree в 2019 году резал кривую ключей жадным «сужающимся конусом», PGM-index в 2020-м — так, чтобы кусков было как можно меньше, и строил над началами кусков такие же куски. Это интерполяционный поиск, у которого вместо одной хорды на всю кривую — по хорде на каждый почти прямой кусок.

Код

Листинг строит из 45 227 городов те же пять рядов, что и сцена, и ищет в каждом каждое 50-е значение тремя способами: по хорде, по очереди и двоичным поиском. Для каждого способа он печатает среднее и наибольшее число проб.

sorted row                 size  interpolation       by turns         binary
                                    mean   max     mean   max     mean   max
hashes of the names       42311      4.0     9      6.7    15     14.4    16
names, a phone book       17205     17.6    60     12.9    23     13.1    15
latitudes                 32407     16.0    43     10.2    17     14.0    15
populations               30576   1587.2  2211     12.3    23     13.9    15
1 .. 4999 and one giant    5000   2476.0  4951      4.9    15     11.3    13

Наибольшее число проб у великана, 4951, — это утверждение о худшем случае в действии: цель на месте 4950 ищется пробами по местам 0, 1, …, 4950. Поиск по очереди на тех же целях делает не больше 15 проб, а его граница — 2⌊log₂ 5000⌋ + 2 = 26.

Пробы и время

Проба не единица времени. В пробе двоичного поиска сравнения и сдвиг, а в пробе по хорде — ещё вычитания, умножение и деление, в Go к тому же 128-битные. Замер ниже ищет случайные ключи среди n случайных 63-битных — на них кривая ровная, как на хешах.

На миллионе ключей поиск по хорде делает в среднем 4,4 пробы, двоичный — 19,0: в 4,3 раза меньше. А по времени в Go разница почти вдвое, в Python — в 1,4 раза. На тысяче ключей в Python поиск по хорде не выигрывает: интерпретатор тратит на умножение и деление больших чисел столько же, сколько сэкономил на пробах. Выигрыш растёт с n: когда массив не помещается в кеш, дорогой становится каждая проба, а не вычисления в ней.

Похожее видели и другие. Перл, Итаи и Авни в 1978 году писали, что в памяти машины интерполяционный поиск идёт вровень с двоичным и чуть быстрее лишь на файлах больше 5000 записей; выигрыш в пробах важен, когда каждая проба — чтение с диска. В git интерполяцию по хешам пробовали с 2007 года: хеши SHA-1 равномерны, лучшего случая не придумать. В 2017-м её убрали — на git rev-list --objects --all в репозитории ядра Linux она давала 37,4 секунды против 35,4 без неё. А в 2019 году Питер Ван Сандт, Яннис Хронис и Джигнеш Патель показали, что аккуратная реализация, которая переиспользует наклон прошлой хорды, бывает до четырёх раз быстрее двоичного поиска на равномерных данных.

Сколько это стоит

Последняя строка — два двоичных поиска: по m началам кусков и по окну из 2ε + 1 мест. Нарезка на куски жадным конусом — один проход, O(n). Средний случай по хорде — log₂ log₂ n только для равномерных ключей. На других распределениях он может быть каким угодно, вплоть до линейного, и без страховки этот поиск разумно ставить только туда, где форма данных известна заранее.

Где работает

• git. В индексе пакета объектов есть таблица из 256 чисел: сколько хешей начинаются с байта не больше данного. По первому байту хеша таблица сразу даёт отрезок длиной около n/256, и двоичный поиск идёт уже по нему. Это интерполяция, сделанная заранее и таблицей. А функция oid_pos в коде git делает по хорде только первую пробу — по паре байтов хеша, — а дальше делит пополам.
• Выученные индексы. FITing-tree, PGM-index и их потомки хранят вместо страниц B-дерева куски прямых с гарантированной ошибкой. Об этом задача статьи.
• Словари. Вырубки по буквам на обрезе — та же таблица, что в git: доля страниц на каждую букву известна заранее, и интерполировать остаётся внутри буквы.
• Корни уравнений. Метод хорд ищет корень монотонной функции так же: проводит прямую через края отрезка и берёт точку, где она пересекает ноль. И болезнь та же — на выпуклой функции один край застревает, и сходимость становится медленной. Лекарство тоже: метод Брента, например в scipy.optimize.brentq, перемежает хорды и параболы с делением пополам.

История

Поиск по прикидке места описал Уэсли Питерсон в 1957 году, в статье о том, как находить записи в памяти с произвольным доступом. Через четыре года он придумает циклический избыточный код CRC. Питерсон полагал, что в среднем поиску нужно не меньше половины log₂ n обращений. Оценка оказалась слишком пессимистичной. В 1976 году Эндрю Яо и Фрэнсис Яо сложным комбинаторным рассуждением доказали, что на случайных равномерных таблицах хватает log₂ log₂ n и что заметно меньше в среднем не сможет никакой поиск. В 1977-м Перл и Эдвард Рейнгольд дали наглядное объяснение, а в 1978-м Перл, Итаи и Авни — своё доказательство и опыт: на 400 000 равномерных ключах в среднем 4,28 обращения при log₂ log₂ n = 4,21.

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

Задачи

В задаче нужно нарезать на прямые куски моменты землетрясений, население и широты городов — тем самым сужающимся конусом — и сосчитать куски. Ключи до 2·10¹⁸, и наклоны приходится сравнивать точно.