Ускорение кода

Что влияет на скорость работы программ на C++ и как её добиться при высоком уровне кода? Ведущий разработчик библиотеки CatBoost Евгений Петров ответил на эти вопросы на примерах и иллюстрациях из опыта работы над CatBoost для x86_64.


Всем привет. Я занимаюсь оптимизацией для CPU библиотеки машинного обучения CatBoost. Основная часть нашей библиотеки написана на C++. Сегодня расскажу, какими простыми способами мы добиваемся скорости.



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



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

Чтобы эту разницу сгладить, в архитектуре есть несколько уровней кэширования. Самый быстрый и самый маленький L1-кэш. Потом есть кэш побольше и помедленнее, второго уровня. И есть совсем большой кэш, который может быть на десятки мегабайт, кэш третьего уровня, но он самый медленный.



Из-за разной скорости обмена данными вычислений, вычислительный код делится на два класса. Один класс лимитируется пропускной способностью, то есть скоростью обмена данными. Второй класс лимитируется скоростью работы процессора. Граница между ними выставляется в зависимости от числа операций, которые выполняются с одним байтом данных. Это обычно константа для конкретного кода.

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



Из оставшегося компиляторы умеют не всё, так как на их разработку тратится очень ограниченный процент ресурсов. Какие из них сегодня развиваются более-менее активно, то есть поддерживают стандарты, пытаются за ними следить? Это frontend EDG, который используется в различных деривативах, например, компилятор Intel; LLVM; GNU и frontend Microsoft.

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

Что же остается разработчикам? Это можно условно поделить на четыре части. Первая архитектура приложения, компиляторы просто не в состоянии придумать ее за нас.



Параллелизация тоже сложная вещь для компиляторов. Работа с памятью потому что это реально сложно: нужно учитывать и архитектуру, и параллелизацию, и все вместе. Кроме этого, компиляторы не умеют нормально оценивать качество оптимизации, то, насколько быстрым получается код. Это тоже приходится делать нам, разработчикам, принимать решение оптимизировать дальше или остановиться.

В части архитектуры мы посмотрим на амортизацию накладных расходов, виртуальных вызовов, на которых в большом числе случаев архитектура основывается.

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

Использование интерфейсов и абстрактных типов данных один из основных методов проектирования. Рассмотрим подобный вычислительный код из машинного обучения. Это условный код, который обновляет прогноз градиентным методом.



Если заглянуть немножко внутрь и попытаться понять, что внутри происходит, то у нас есть интерфейс IDerCalcer для вычисления производных функции потерь, и функция, которая сдвигает прогноз (наше предсказание) в соответствии с градиентом функции потерь.

Справа на слайде вы видите, что это значит для двумерного случая. А в машинном обучении размер прогноза не два и не три, а миллионы, десятки миллионов элементов. Посмотрим, насколько хорош этот код для вектора порядка 10 млн элементов.



Возьмем в качестве целевой функции среднеквадратичное отклонение и померяем, как она работает, за сколько она сдвинет этот прогноз. Производная этой целевой функции показана на слайде. Время работы на условной машине, которая дальше остается фиксированной, 40 мс.



Попробуем понять, что же тут не так. Первое, что привлекает внимание, виртуальные вызовы. Если смотреть в профилировщик, видно, что в зависимости от числа параметров это примерно пять-десять инструкций. И если, как в нашем случае, вычисление самой производной это всего два арифметических действия, то это запросто может оказаться существенным накладным расходом. Для большого тела при вычислении производных это ок. Для короткого тела, которое вычисляет производную скажем, даже не 500 инструкций, а 20, 50 или еще меньше, это уже будет существенный процент по времени. Что же делать? Попробуем самортизировать вызов виртуальной функции, изменив интерфейс.



Первоначально мы вычисляли производные поточечно, для каждого элемента вектора отдельно. Перейдем от поэлементной обработки к обработке векторами. Возьмем стандартный шаблон C++, который позволяет работать с view на вектор. Или, если ваш компилятор не поддерживает последний стандарт, можно использовать простой самодельный класс, где хранится указатель на данные и размер. Как поменяется код? У нас останется один вызов, который вычисляет производные, и потом нам придется добавить цикл, который будет, собственно, обновлять прогноз.



Кроме того, что добавляется цикл, нам еще придется второй раз посмотреть на данные, то есть второй раз прочитать сам вектор прогнозов и вектор градиента, который мы только что вычислили.



Померяем снова на той же машине и увидим, что получилось хуже, что-то пошло не так. Давайте разбираться, что случилось с кодом.



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



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



Когда посчитан вектор производных, дело доходит до цикла, который сдвигает прогноз. Перед этим циклом в быстром кэше первого уровня, который работает на частоте процессора, останется только очень маленькая часть данных. На слайде это зеленый цвет на светофоре. Остальные данные будут вытеснены из кэша в память, и когда цикл пойдет обновлять прогнозы, данные придется второй раз читать из памяти. А она у нас работает, в общем, весьма небыстро, со скоростью пешехода.



Когда мы обновляем прогноз, нам не обязательно читать все производные сразу. Достаточно считать их большими пачками, чтобы самортизировать виртуальные вызовы. Поэтому есть смысл раздробить вычисление производных и обновление прогноза на небольшие блоки и перемешать эти два действия. К чему это приведет, если посмотреть на то, откуда будут читаться данные?



Это приведет к тому, что мы будем все время брать данные, и к тому, что данные будут оставаться в L1-кэше и не успеют уйти в медленную память. А дальше нам нужно понять, кто же нам скажет этот размер блока.



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



Вот он, внешний по блокам.



А вот внутренний по элементам блока.



Мы учитываем, что последний блок может быть неполным.



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



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

Чтобы это исправить, нужно вынести аллокацию из цикла. Для этого нам придется еще раз поменять интерфейс, перестать возвращать векторы и начать записывать производные в память, которую мы получаем от вызывающего кода.



Это стандартный прием вынос всех манипуляций с ресурсами из узких мест в вычислительном коде. Добавим к методу CalcDer еще один параметр, view на вектор, куда должны попасть производные.



Код тоже изменится очевидным образом. Вектор производных будет один, снаружи от всех циклов, и к методу просто добавится новый параметр.



Смотрим. Получается, что мы выиграли еще где-то восемь процентов по сравнению с предыдущей версией, а по сравнению с базовой уже 15%.

Понятно, что оптимизации не исчерпываются амортизацией накладных расходов, что узкие места бывают и других видов.



Для иллюстрации того, как искать узкие места, понадобится еще один простой подопытный код. Например, я взял транспонирование матрицы. У нас есть матрица approx и матрица approxByCol, куда мы должны положить транспонированные данные. И простое гнездо из двух циклов. Здесь нет никаких виртуальных вызовов, создания векторов. Это просто перекладывание данных. Цикл относительно удобен для компилятора.

Померяем, как работает этот код на достаточно большой матрице и на конкретной машине.



Для примера я взял число строк 1000, число колонок 100 000. Машина интеловский сервер, одно ядро. Память именно такая, нам это важно, потому что вся работа с памятью и скорость будет зависеть от скорости работы памяти. Замерили и получили 1,4 с. Много это или мало? Что мы успеваем сделать за это время?



Мы успеваем прочитать 800 мегабайт, это не транспонированная матрица, а исходная. А также прочитать и записать 1,6 ГБ, это уже транспонированная матрица. Процессор выполняет вспомогательное чтение перед записью, чтобы проинициализировать данные в кэше.



Посчитаем, сколько пропускной способности мы использовали с пользой. Получается, пропускная способность нашего кода составила 1,7 ГБ/с.



Это был теоретический расчет. Дальше можно взять профилировщик, который умеет мерить скорость работы с памятью. Я взял VTune. Посмотрим, что он покажет. Показывает похожую цифру 1,8 ГБ. В принципе хорошо согласуется, потому что в нашем расчете мы не учитывали, что приходится читать адреса строк и адреса колонок. Плюс к этому VTune регистрирует фоновую деятельность в операционной системе. Поэтому наша модель согласуется с реальностью.

Чтобы понять, 1,7 ГБ это много или мало, нужно выяснить, какая максимальная пропускная способность нам доступна.

Для этого нужно почитать спеки по процессору. Есть специальный сайт ark.intel.com, где про любой процессор можно все узнать. Если мы смотрим конкретно на наш сервер, то видим, что у него восемь ядер и для самой быстрой памяти DDR3, которую он поддерживает, обеспечивается передача со скоростью примерно 60 ГБ/с в одну сторону.



Но тут надо учесть, что мы используем только одно ядро и память у нас помедленнее, то есть нужно масштабировать эти 60 ГБ на наши условия пропорционально числу ядер и частоте памяти.

Получается, что наш код мог бы использовать 5,3 ГБ в одну сторону. А поскольку параллельно можно читать и записывать, то в идеале, если бы мы просто копировали данные с места на место, достигли бы 10,6. Поскольку у нас два чтения и одна запись, то должно быть примерно 8 ГБ/с. Мы помним, что у нас получилось 1,7. То есть мы использовали где-то 20%.

Почему так получается? Снова нужно разбираться с архитектурой. Дело в том, что данные между памятью и кэшем передаются не произвольными пакетами, а по 64 байта ровно, ни больше ни меньше. Это первое соображение.



Второе соображение: мы записываем транспонированные данные не последовательно, а произвольно, поскольку строки матрицы расположены в памяти непредсказуемым образом.

Получается, что перед тем, как записать одно вещественное число, нам приходится вычитывать 64 байта данных. Если обозначить размер матрицы N, то вместо оптимального времени работы (N/5,3 + N/10,6) у нас получается (8*N/5,3 + N/10,6). Где-то в четыре-пять раз больше, что и объясняет эту эффективность в 20%.



Что с этим делать? Нужно перестать записывать данные по одному столбцу и начать записывать столько колонок, сколько укладывается в одну кэш-линию (64 байта). Для этого разобьем цикл по столбцам на цикл по кэш-линиям и вложенный цикл по элементам кэш-линии.



Вот они, итерации по кэш-линиям.



И вот они, итерации внутри кэш-линии. Тут мы для простоты считаем, что данные выровнены на границу кэш-линии. Теперь проверим с помощью VTune, что получится.



Видим, что получилось близко к расчетным восьми гигабайтам в секунду 7,6. Но еще не факт, что все эти 7,6 полезная работа. Может быть, сколько-то из них накладные расходы.

Чтобы понять, сколько пользы мы получили, измерим время работы после оптимизации. Получается 0,5 с на той же самой машине. Пропускная способность, которая приходится на само транспонирование, стала 4,8 ГБ/с. Видно, что еще остался запас, который мы не выбрали, но все равно, мы из 20-процентной эффективности получили 60-процентную.

С помощью профилировщика можно разобраться, почему мы не получили 80% или 95%.



Дело в том, что мы храним матрицы в виде вектора векторов, то есть используем доступ в память с двойной косвенностью.



С помощью VTune можно увидеть, какие инструкции сгенерированы для доступа к элементам массива. Желтым цветом слева выделены инструкции, которые вычитывают адреса колонок транспонированной матрицы. И это, во-первых, дополнительные инструкции, а во-вторых, дополнительные передачи данных. Но дальше оптимизировать уже не будем, остановимся и подведем итоги.



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

Вынос аллокаций из узких мест это тоже их амортизация. И тоже повышение скорости доступа за счет фиксации временных буферов в памяти.

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

На этом у меня все. Если вы пользуетесь CatBoost или первый раз про него услышали и хотите узнать, что это такое, приходите к нам на GitHub, пишите в телеграм. Большое спасибо за внимание.

Множество (Set) структура данных, которая позволяет достаточно быстро (в зависимости от реализации) применить операции add, erase и is_in_set. Но иногда этого не достаточно: например, невозможно перебрать все элементы в порядке возрастания, получить следующий / предыдущий по величине или быстро узнать, сколько элементов меньше данного есть в множестве. В таких случаях приходится использовать Упорядоченное множество (ordered_set). О том, как оно работает, и какие реализации есть для питона далее.

Стандартный Set

В языке Python есть стандартная стукрура set, реализованная с помощью хэш-таблиц. Такую структуру обычно называют unordered_set. Данный метод работает так: каждый элемент присваивается какому-то классу элементов (например, класс элементов, имеющих одинаковый остаток от деления на модуль). Все элементы каждого класса хранятся в одтельном списке. В таком случае мы заранее знаем, в каком списке должен находиться элемент, и можем за короткое время выполнить необходимые операции. Равновероятность каждого остатка от деления случайного числа на модуль позволяет сказать, что к каждому классу элементов будет относиться в среднем size / modulo элементов.

Но хэш-таблица не позволяет выполнить операцию count_lower или подобные, поэтому придётся использовать другие структуры данных.

Что есть в других языках

В языке c++ есть структура std::set, которая поддерживает операции изменения, проверку на наличие, следующий / предыдущий по величине элемент, а также for по всем элементам. Но тут нет операций получения элемента по индексу и индекса по значению, так что надо искать дальше (индекс элемента количество элементов, строго меньших данного)

И решение находится достаточно быстро: tree из pb_ds. Эта структура в дополнение к возможностям std::set имеет быстрые операции find_by_order и order_of_key, так что эта структура именно то, что мы ищем.

Она реализована с помощью красно-чёрных деревьев. Смысл этой струкруты в том, что все элементы образуют собой двоичное дерево поиска, которое балансируется так, чтобы высота не превышала логарифм. Нам это даёт возможность с помощью одного спуска по дереву выполнить необходимые операции. Также с этой задачей может справиться Декартово дерево (Дерамида) по неявному ключу или AVL дерево.

Таким образом, целью этой статьи станет поиск аналога этой структуры в Python.

Как будем тестировать скорость работы структур данных

Для оценки времени работы я написал программу, которая будет выполнять последовательно несколько типов операций:

1. Добавление в множество миллиона случайных чисел (при данном сиде среди них будет 999'936 различных)
2. Проверка миллиона случайных чисел на присутствие в множестве
3. Прохождение циклом по всем элементам в порядке возрастания
4. В случайном порядке для каждого элемента массива узнать его индекс (а, соответственно, и количество элементов, меньше данного)
5. Получение значения i-того по возрастанию элемента для миллиона случайных индексов
6. Удаление всех элементов множества в случайном порядке

from SomePackage import ordered_setimport randomimport timerandom.seed(12345678)numbers = ordered_set()# adding 10 ** 6 random elements - 999936 uniquelast_time = time.time()for _ in range(10 ** 6):    numbers.add(random.randint(1, 10 ** 10))print("Addition time:", round(time.time() - last_time, 3))# checking is element in set for 10 ** 6 random numberslast_time = time.time()for _ in range(10 ** 6):    is_element_in_set = random.randint(1, 10 ** 10) in numbersprint("Checking time:", round(time.time() - last_time, 3))# for all elementslast_time = time.time()for elem in numbers:    now_elem = elemprint("Cycle time:", round(time.time() - last_time, 3))# getting index for all elementslast_time = time.time()requests = list(numbers)random.shuffle(requests)for elem in requests:    answer = numbers.index(elem)print("Getting indexes time:", round(time.time() - last_time, 3))# getting elements by indexes 10 ** 6 timesrequests = list(numbers)random.shuffle(requests)last_time = time.time()for _ in range(10 ** 6):    answer = numbers[random.randint(0, len(numbers) - 1)]print("Getting elements time:", round(time.time() - last_time, 3))# deleting all elements one by onerandom.shuffle(requests)last_time = time.time()for elem in requests:    numbers.discard(elem)print("Deleting time:", round(time.time() - last_time, 3))

SortedSet.sorted_set.SortedSet

Пакет с многообещающим названием. Используем pip install sortedset

К сожалению, автор не приготовил нам функцию add и erase в каком-либо варианте, поэтому будем использовать объединение и вычитание множеств

Использование:

from SortedSet.sorted_set import SortedSet as ordered_setnumbers = ordered_set()numbers |= ordered_set([random.randint(1, 10 ** 10)])  # добавлениеnumbers -= ordered_set([elem])  # удаление

Протестируем пока на множествах размера 10'000:

Как так получилось? Давайте загляем в исходный код:

def __init__(self, items=None):    self._items = sorted(set(items)) if items is not None else []def __contains__(self, item):    index = bisect_left(self._items, item)

Как оказалось, это обычный массив, в котором наличие элемента определяется бинпоиском. Это действительно отсортированное множество, но очень ленивое.

Вывод: почти бесполезно, несколько строчек кода завернули в класс

sortedcontainers.SortedSet

Внеший пакет, для установки можно использовать pip install sortedcontainers. Посмотрим же, что он нам покажет

Но, не смотря на это, кажется мы нашли то, что искали! Все операции выполняются за приличное время. По сравнению с ordered_set некоторые операции выполняются дольше, но за то операция discard выполняется не за o(n), что очень важно для возможности использования этой структуры.

Также пакет нам предлагает SortedList и SortedDict, что тоже может быть полезно.

И как же оно работает?

На странице пакета мы можем прочитать, что реализована структура не так, как мы предполагали в начале статьи.

Из-за особенностей реализации языка Python, в нём быстро работают list, а также bisect.insort (найти бинарным поиском за o(log n) место, куда нужно вставить элемент, а потом вставить его туда за o(n)). Insert работает достаточно быстро на современных процессорах. Но всё-таки в какой-то момент такой оптимизации не хватает, поэтому структуры реализованы как список списков. Создание или удаление списков происходит достаточно редко, а внутри одного списка можно выполнять операции даже за быструю линию.

Если говорить кратко, то принцип действия похож на корневую оптимизацию.

Проблема с ordered_set

Что вообще такое упорядоченное множество? Это множество, в котором мы можем сравнить любые 2 элемента и найти среди них больший / меньший. В течение всей статьи под операцией сравнения воспринималась операция сравнения двух элеметнов по своему значению. Но все пакеты называющиеся ordered_set считают что один элемент больше другого, если он был добавлен раньше в множество. Так что с формулировкой ordered_set нужно быть аккуратнее и уточнять, имеется ввиду ordered set или sorted set.

Bintrees

Так есть же модуль bintrees! Это же то, что нам нужно? И да, и нет. Его разработка была приостановлена в 2020 году со словами Use sortedcontainers instead.

Пакет предлагает нам несколько структур. К сожалению, ни одна из них не поддерживает операции find_by_order и подобные, так что эти струкруты являются аналогами std::set. Посмотрим же, на что они способны:

pip install bintrees

Название AVLTree говорит само за себя, RBTree красно-чёрное дерево, BinaryTree несбалансированное двоичное дерево, префикс Fast означает реализацию на Cython (соответственно, необходимо наличие Visual C++, если используется на Windows).

Результаты тестирования отчётливо показывают нам, почему использовать деревья поиска на Python плохая идея в плане производительности. А вот в интеграции с Cython всё становится намного лучше.

Оказывается, эта структура и SortedSet очень похожи по производительности. Все 3 Fast версии структур bintrees достаточно близки, поэтому будем считать, что оттуда мы используем FastAVLTree.

Как мы видим, AVL в полтора раза быстрее в скорости добавления элементов и почти в 2 раза быстрее в операциях удаления. Но он в те же 2 раза медленнее в проверке на наличие и цикле по всем элементам. К тому же не стоит забывать, что 2 операции он выполнять не умеет, то есть не является тем ordered_set, что мы ищем.

Использование:

import bintreesnumbers = bintrees.FastAVLTree()numbers.insert(value, None)  # второй параметр - значение, как в словаре

Что же выбрать

Мои рекомендации звучат так: если вам нужны операции find_by_order и order_of_key, то ваш единственный вариант sortedcontainers.SortedSet. Если вам нужен только аналог std::map, то выбирайте на своё усмотрение между SortedSet и любым из fast контейнеров из bintrees, опираясь на то, каких операций ожидается больше.

Можно ли сделать что-то быстрее

Скорее нет, чем да. Использование Cython один из самых мощных способов оптимизации, а AVL считается очень быстрым решением исходной задачи. Про остальные операции ordered_set можно сказать, что модификация красно-чёрного дерева так, чтобы оно поддерживало эти операции, вряд ли будет быстрее SortedContainers, так что смысла изобретать велосипед я не вижу.

---

Облачные VPS серверы от Маклауд быстрые и безопасные.

Зарегистрируйтесь по ссылке выше или кликнув на баннер и получите 10% скидку на первый месяц аренды сервера любой конфигурации!