Создание и практическое использование алгоритмов сильно зависит от возможности эффективно их реализовать. В лаборатории языковых инструментов JetBrains разрабатывают алгоритмы поиска путей в помеченных графах с дополнительными ограничениями. Эти алгоритмы достаточно естественно выражаются в терминах операций над булевыми матрицами, но в современных высокопроизводительных библиотеках линейной алгебры пока нет полного набора необходимых операций над булевым полукольцом. Поэтому мы решили их реализовать.
Всем привет! Меня зовут Мария Карпенко и в этом году я окончила корпоративную магистерскую программу JetBrains в Университете ИТМО по направлению Software Engineering. Я начала интересоваться программированием, обучаясь в СПбГЭУ на направлении «Математические методы в экономике», и на последнем курсе поступила в Computer Science Center (CSC). Именно в CSC я узнала о существовании корпоративной магистратуры и однозначно определилась с выбором программы. Сейчас я стажируюсь в «Яндексе» в команде Трекера.
В первый год магистратуры студенты активно приобретают новые знания, а на втором курсе начинается работа над дипломом. Я старалась выбрать для себя тему, которая позволит погрузиться в новые технологии и алгоритмы, с акцентом на практическую часть. В итоге я остановилась на теме «Реализация операций с разреженными булевыми матрицами на OpenCL». Работу выполняла под руководством Семёна Григорьева, исследователя лаборатории языковых инструментов JetBrains.
Операции над матрицами и графы
Многие алгоритмы для анализа графов можно выразить через операции с их матрицами смежности. Например, транзитивное замыкание, поиск в ширину, прямое произведение графов. Такие базовые матричные операции входят в стандарт GraphBLAS, разработанный для описания алгоритмов обработки графов в терминах линейной алгебры.
В большинстве практических задач графы оказываются сильно разреженными, что приводит к необходимости хранить матрицы и осуществлять операции с ними в сжатых форматах. В лаборатории языковых инструментов JetBrains разрабатывают алгоритмы анализа графов, сводимые к операциям над разреженными булевыми матрицами. Для подтверждения их практической применимости возникла необходимость в эффективной реализации следующих операций с разреженными булевыми матрицами:
матричное умножение;
матричное сложение;
транспонирование матрицы;
извлечение подматрицы;
редуцирование строк матрицы в вектор;
произведение Кронекера.
Специальные реализации этих операций для матриц, определенных над булевым полукольцом, позволяют, во-первых, экономить память: все форматы для хранения разреженных матриц содержат массив со значениями, который в случае булевой матрицы не нужен. Во-вторых, операции с элементами разреженных булевых матриц сводятся к определению отсутствия или наличия значения на конкретной позиции, что значительно быстрее операций с типами float или double.
Современные высокопроизводительные библиотеки линейной алгебры содержат лишь некоторые из необходимых операций, чаще всего умножение и сложение без специализации для булевых матриц. На момент начала работы библиотека cuBool, разрабатываемая в той же лаборатории языковых инструментов, содержала полный набор необходимых операций с реализацией на CUDA. На CPU отмечу библиотеку SuiteSparse, в которой полностью реализован стандарт GraphBLAS.
Моей задачей было создать полное решение на OpenCL, которое позволит задействовать видеокарты от AMD и другие устройства, поддерживающие стандарт. От меня не требовалось изобретать новые алгоритмы для операций с матрицами, необходимо было выбрать лучшие и обосновать свой выбор.
Знакомство с технологией и сжатыми форматами
Работа началась с погружения в высокопроизводительные вычисления, в чем мне помог курс от CSC «Вычисления на видеокартах». Помимо доступных и понятных объяснений принципов программирования на видеокартах, курс содержит ряд практических советов и приемов, связанных как с организацией кода, так и с поиском ошибок параллельного программирования. В рамках курса мы реализовывали многие необходимые промежуточные операции, которые пригодились мне в дальнейшем.
Также в начале работы мы определились с форматом хранения матриц. Так как нам заранее неизвестна структура матрицы, были выбраны наиболее общие и популярные форматы: COO и CSR (CSC). Именно с этими форматами работают большинство библиотек линейной алгебры. Так как в наших задачах матрицы могут быть сверхразреженными, мы сосредоточились на форматах, в которых матрица занимает O(nnz) памяти, где nnz — количество ненулевых элементов. Таким свойством обладают форматы COO и DCSR, полученный из CSR удалением информации о пустых строках.
Я реализовала операции в разных форматах, выбирая тот или иной формат исходя из особенностей работы с ним. Также были реализованы конвертации между форматами, которые осуществляются достаточно быстро. Некоторые операции по-разному можно реализовать как в COO, так и в DCSR, поэтому основную реализацию я выбирала по результатам замеров.
Выбранные форматы для операций:
COO:
транспонирование;
сложение;
произведение Кронекера.
DCSR:
умножение;
извлечение подматрицы;
редуцирование строк матрицы;
произведение Кронекера.
О реализациях операций
Первыми я реализовала сложение, транспонирование и произведение Кронекера в формате COO. Транспонирование в этом формате получается очень естественно — достаточно переставить массивы строковых и столбцовых индексов и отсортировать результат. Сложение в формате COO — это слияние двух отсортированных массивов с последующим удалением дубликатов.
Произведение Кронекера в COO также достаточно простое, но, как выяснилось позже, реализация через DCSR работает на порядок быстрее.
Для остальных операций однозначно был выбран формат DCSR, так как в них часто требуется индексация строк, которая в COO осуществляется медленнее.
Самая сложная из реализуемых операций — умножение. Алгоритмы умножения матриц развиваются уже несколько десятилетий, и я ознакомилась со многими работами в этой области. Среди всех выделяются работы W. Liu и Y. Nagasaka. Они подтвердили свою практическую применимость в библиотеках, а алгоритм Y. Nagasaka хорошо показал себя в cuBool. Я реализовала оба варианта, и в моей реализации алгоритм Y. Nagasaka оказался лучше и по времени, и по памяти.
Отдельно хочется отметить, что алгоритм Y. Nagasaka значительно проще в реализации: он состоит из однообразных ядер (программ для видеокарты), которые отличаются только распределением работы потоков, размером и размещением структур данных. Также алгоритм Y. Nagasaka активнее использует локальную память видеокарты, скорость обращения к которой в сотни раз выше, чем к глобальной памяти.
Извлечение подматрицы и редуцирование строк матрицы в вектор были реализованы последними в формате DCSR.
Эксперименты
В экспериментальной части мы сосредоточились на операциях, которые есть в большинстве библиотек линейной алгебры: сложение и умножение. Сложение мы тестировали на операции M + M2, а умножение на M2. Для замеров были выбраны десять сильно разреженных матриц из коллекции SuiteSparse, но далее я приведу результаты по трем матрицам, которые можно считать наиболее показательными: две матрицы с равномерно распределенными ненулевыми элементами разных размеров (строки 1 и 3) и матрица с некоторыми плотными строками (строка 2).
Операции выполнялись на трех различных устройствах: две видеокарты (AMD, NVIDIA) и одна FPGA.
Скомпилировать решение на FPGA позволяет технология Intel FPGA SDK for OpenCL. За предоставленный сервер с устройством Arria10 благодарим компанию Selectel.
На FPGA мои реализации сгенерировались в неэффективную архитектуру, и этому есть несколько объяснений. Сами алгоритмы сильно заточены под устройство видеокарт: нагрузка распределяется с учетом организации работы потоков по ворпам и мультипроцессорам. Так как видеокарте не нужно перестраивать свою архитектуру под каждый бинарный файл, мы можем быстро скомпилировать множество ядер, суммарный объем оперативной памяти в которых ничем не ограничен, но в каждом отдельном ядре нужно учитывать ограничения мультипроцессора. С FPGA все иначе: ограничение в 16 КБ затрагивает всю программу, и мне не удалось, например, скомпилировать полный набор ядер для умножения. Пришлось уменьшать число ядер, что привело к неэффективному решению. Время выполнения обеих операций оказалась в сотни раз медленнее CPU-библиотеки SuiteSparse.
На операции сложения в моей реализации разница между результатами AMD и NVIDIA объясняется характеристиками видеокарт — AMD почти в 2 раза мощнее.
В умножении преимущества AMD становится заметно только на матрицах с более плотными строками:
(1) — алгоритм W. Liu
(2) — алгоритм Y. Nagasaka
В таблицах ниже представлено сравнение с некоторыми библиотеками линейной алгебры на видеокарте NVIDIA. Как видно из замеров, сложение через перевод матриц в координатный формат реализуют все библиотеки, кроме cuBool. Отмечу, что по потреблению памяти cuBool серьезно выигрывает, так как в используемом алгоритме не выделяется дополнительная память на некоторые промежуточные вычисления.
Несмотря на один алгоритм и похожую реализацию, на плотных матрицах cuBool работает заметно быстрее моей реализации clBool. Также cuBool и clBool требуют в среднем на треть меньше памяти, чем другие библиотеки.
Заключение
Результат моей работы используется в качестве OpenCL-бэкенда проекта spbla, который был создан для реализации, тестирования и изучения алгоритмов анализа графов, выразимых в терминах операций над булевыми матрицами.
Набор решаемых задач все еще серьезно ограничен по памяти: в ходе итераций алгоритма матрицы становятся плотнее и перестают помещаться в видеопамять. Естественное направление развития проекта — это распределенные алгоритмы для матричных операций и использование виртуальной памяти. Так как библиотеки с подобными решениями появятся нескоро, имеет смысл создавать решения в рамках лаборатории.
Хочется выразить благодарность лаборатории языковых инструментов за предоставление серверов с видеокартами для экспериментов. Отдельно выражаю благодарность научному руководителю Семёну Григорьеву за непрерывное и чуткое руководство, помощь в работе с текстом.
До 9 августа продолжается прием документов на корпоративную магистерскую программу JetBrains «Разработка программного обеспечения». В этом году набор ведется на 30 бюджетных и 5 платных мест. Подробнее о программе можно узнать на сайте или в Telegram-чате для абитуриентов.
