CC BY
89
УДК 519.612
ЭФФЕКТИВНЫЙ ПАРАЛЛЕЛЬНЫЙ РАСЧЕТ ЗАДАЧ МЕХАНИКИ СПЛОШНОЙ СРЕДЫ ДЛЯ СИСТЕМ ГИБРИДНОЙ АРХИТЕКТУРЫ
М.И. КЛИМОВ, И.С. МЕНЬШОВ
В статье предлагается общий подход к эффективному параллельному решению задач механики сплошных сред для систем гибридной архитектуры. Ввиду своих архитектурных особенностей эффективность расчета на GPU сильно зависит от многих факторов. Предложенный подход к построению параллельных программ позволит достичь корректного двухуровневого параллелизма: на уровне декомпозиции области и поэлементого параллелизма, при этом обеспечив должную производительность на обоих уровнях.
Ключевые слова: параллельные алгоритмы, высокопроизводительные вычисления, GPGPU, численные методы.
Введение
В настоящее время требования для расчета задач аэродинамики к вычислительным системам очень высоки: помимо сложной геометрии, сложных физических моделей, современные требования включают в себя расчет на больших сетках и ресурсоемкую вычислительную часть. Всё это заставляет строить огромные вычислительные комплексы с уникальной архитектурой. Обеспечить большую производительность таких систем до недавних пор позволяло увеличение тактовой частоты процессора, но на сегодняшний день учеными практически достигнут "кремниевый тупик". Рост производительности таких систем в последнее время происходит благодаря увеличению числа процессоров в системе или числа ядер на процессоре. Таким образом, ядра центральных процессоров вычислительных узлов параллельно решают необходимую задачу. Параллельное взаимодействие в таком случае, в основном, происходит на основе технологий MPI[2] и/или OpenMP.
Графические карты появились еще в 1981 году и были призваны освободить часть нагрузки на CPU по работе с графикой. В связи с особенностью работы с текстурной графикой GPU быстро превратился в простой SIMD-процессор, грубо говоря, процессор, одинаково обрабатывающий несколько потоков различных данных (рис. 1).
Входящий поток Выходящий поток
SIMD
Рис. 1. SIMD-модель
Из-за своей архитектуры, максимально нацеленной на увеличение скорости расчёта текстур и сложных графических объектов, а также ограниченного набора команд графический процессор долгое время занимался лишь ускорением трехмерной графики и не мог использоваться в качестве устройства для математических расчетов. Это стало возможным совсем недавно благодаря добавлению программируемых шейдерных блоков и более высокой арифметической точности растровых конвейеров. Ряд крупных производителей графических карт выпустили свои программные инструменты для возможности проведения общих вычислений на GPU. Таким образом, появилась техника GPGPU (англ. General-purpose graphics processing units - «GPU общего назначения») - техника использования графического процессора видеокарты для общих вычислений, которые обычно проводит центральный процессор. Программные реализации GPGPU доступные сегодня:
1. CUDA - технология кампании NVidia;
2. OpenCL - открытая разработка;
3. FireStream — продукт кампании AMD;
4. DirectCompute — реализация от корпорации Microsoft.
CUDA и OpenCL на данный момент доминируют в секторе GPGPU с большим преимуществом технологии CUDA.
Рис. 2. Теоретическая вычислительная производительность и пропускная способность памяти
За счет наличия сотен скалярных ядер, одновременно обрабатывающих данные, пиковая производительность GPU оказывается на порядки выше, чем CPU. Память графического ускорителя обладает многоуровневой системой. Основной памятью служит так называемая глобальная память видеокарты, которая представлена наибольшим объемом. Но глобальная память латентна, и работа с ней часто является «узким местом» программы. Для уменьшения количества обращений к глобальной памяти применяются шаблоны запросов, при которых несколько обращений могут объединяться в одно (coalescing). Некоторые данные по производительности можно наблюдать на рис. 2.
Программирование с помощью техники GPGPU сложно из-за большого количества специфичных для архитектуры ограничений и особенностей и может существенно отличаться от программирования для CPU. Ввиду таких архитектурных особенностей для успешной адаптации кода под видеокарты необходимо учесть множество факторов, прямо или косвенно влияющих на производительность. Среди таких факторов можно выделить следующие:
• наличие ветвлений снижает производительность;
• специальным образом организованная память может существенно повысить эффективность;
• увеличение количества нитей может компенсировать бесполезное ожидание данных, но уменьшает количество регистровой памяти, выделенной потоку.
Учитывая подобные факторы при разработке, можно получить ускорение на порядки по сравнению с аналогичным расчетом на одном ядре CPU.
Исходная задача
с Ч 1 с 1 I
Подзадача Подзадача Подзадача
Рис. 3. Модель двухуровневого параллелизма
Безусловно, эффективность GPGPU высока, но ограничена. Для дальнейшего увеличения производительности используют метод разбиения общей задачи на подзадачи, решаемые параллельно. В таком случае обычно решение каждой подзадачи происходит на отдельном вычислительном объекте. Взаимодействие между процессами происходит с помощью MPI. Данная техника применима как для систем с общей памятью, так и для систем с раздельной памятью, что делает ее самой популярной на сегодняшней день. Каждый параллельный процесс при этом может вести свой параллельный расчет на GPU, формируя второй уровень параллелизма (рис. 3). Связка MPI и OpenMP является еще одним примером двухуровневого распараллеливания.
1. Численные модели механики сплошных сред
Несмотря на разнообразие численных и физических моделей механики сплошных сред, среди них можно выделить несколько общих концептуальных положений:
• законы сохранения (массы, импульса, энергии);
• происходит пространственная дискретизация области (строится структурированная или неструктурированная сетка);
• определяющие уравнения можно представить в дивергентной форме —Q + divF = H ;
dt
• вводится понятие вектора состояния ячейки Q = Q ( X, t) = (р, pUk,...) ;
• связи дискретной модели локальны.
При численном моделировании задач механики сплошных сред требуются знания не только физики и математики, но и программирования. Разработчику приходится решать проблемы организации памяти, хранения данных, параллельных вычислений. Подобного рода программные задачи, вообще говоря, не зависят от конкретной задачи механики сплошных сред и могут быть реализованы отдельно, составив так называемый «базис». В свою очередь, выбранные для решения физическая и математическая модели, опирающиеся на приведенные выше общие положения, могут быть также программно выделены в независимый модуль - так называемое «вычислительное ядро». Таким образом, становится возможным разработать обобщенную программную модель, в которой возможно разделение реализации вычислительного ядра от остальной части. Это обеспечит возможность более легкой разработки вычислительного ядра, а также его простую замену, не касаясь чисто программных проблем, позволив конечному разработчику вычислительного ядра сосредоточиться на физической и численной моделях. При этом «базис» должен обеспечить в первую очередь корректность и высокую эффективность при параллельном расчете.
Исходя из общих принципов дискретных моделей, стала возможна разработка обобщенной программной модели для современных параллельных вычислений, а также реализация модели в виде программного пакета базовых функций.
2. Обобщенная программная модель
Разработка обобщенной программной модели для эффективного решения задач механики сплошных сред сводится к формированию ряда правил и рекомендаций по созданию системы, обеспечивающего максимальную универсальность, расширяемость и, главное, корректность и высокую производительность параллельных вычислений. Следуя данной модели, становится возможным разделение разработки вычислительного ядра от остального кода, снабдив разработчика ядра необходимыми инструментами для проведения успешного и параллельного расчета. Таким инструментом может стать программная библиотека, содержащая широкий набор функций. Далее следует описание подхода по созданию представленной выше системы.
На расчетной области строится структурированная или неструктурированная сетка с некоторым количеством элементов. Каждый элемент окружен конечным количеством
соседей. Каждой ячейке приписан вектор состояния / некоторой размерности, содержащий в качестве компонент примитивные переменные расчета и другие переменные, необходимые для расчета (рис. 4).
Рис. 4. Пример двумерной дискретной модели на регулярной сетке f = (р,Vx,V ,p)
В процессе расчета выбранной ячейки производятся математические операции с вектором состояния ячейки и ее соседей. Обновление вектора состояния может происходить лишь в текущей ячейке. Доступ к данным расчетной сетки не зависит от вычислительного ядра и может быть реализован независимо, на основании общих положений.
Использование параллельных технологий влечет за собой значительное усложнение системы и порождает целый ряд особенностей и специфичных для параллельного программирования проблем. Их решение нетривиально и может отнять у разработчика много сил. Далее, в следующих разделах будут рассмотрены примеры возникнувших проблем и методы их решения. В рамках данной работы необходимо обеспечить 2 уровня параллелизма: основанное на межузловом взаимодействие и поэлементное распараллеливание.
При этом согласно обобщенной программной модели разработчик вычислительного ядра не должен самостоятельно заниматься вопросами эффективного распараллеливания, представления данных в памяти компьютера. Его задача состоит в том, чтобы предоставить полный алгоритм расчета одной ячейки.
Конечно же, данные рассуждения вносят некоторые ограничения в математические модели. Во-первых, в шаблон включены лишь соседи первого уровня. Во-вторых, не рассматриваются адаптивные сетки. В-третьих, численные методы, применимые в данном подходе, не должны требовать обращения глобальной матрицы - это итерационные методы, такие как явные методы, LU-SGS, GMRES, метод Якоби и другие.
3. Раскраска сетки
Данный раздел описывает общий подход к проблеме поэлементного распараллеливания. На данном уровне параллелизма может быть применена не только техника GPGPU, но и другие различные технологии, основанные на многопоточности, допустим OpenMP или Intel Cilk Plus. В дальнейшем на втором уровне параллелизма будет подразумеваться использование GPGPU как наиболее специфичной технологии параллельного программирования.
Подход GPGPU, хотя и похож на работу с SIMD-моделью, имеет некоторые отличия. Нити разбиваются на группы по 32 нити, называемые warp'ами. Только нити в пределах одного warp'а выполняются физически одновременно, при этом быстрое управление нитями прозрачно осуществляет мультипроцессор GPU. В рамках такой многопоточной модели возможно возникновение состояния гонки, когда результат зависит от того, в каком порядке исполняются участки кода. Такие случаи приводят к коллизиям данных, называемые коллизиями типа «read-after-write».
Рис. 5. Пример возникновения коллизий данных К примеру, пусть ячейки А и В имеют векторы состояния и соответственно, как
показано на рис. 5. В процессе расчета векторы состояния обновляются / —/ / —/и. При
'в ся
расчете ячейки А запрашивается вектор состояния ячейки В. Какой вектор получим, ]В или/В,?
Результат неоднозначен, а значит, такой расчет некорректен. Для решения проблемы гонок данных предлагается применить «раскраску» сетки, т.е. каждой ячейке назначается цвет таким образом, чтобы никакие ячейки одного цвета не были соседними друг с другом. Простой схематичный пример раскраски сетки изображен на рис. 6. Так как каждая расчетная сетка представляется графом, вершины которого взаимно-однозначно отображаются в элементы сетки, а инцидентные ребру вершины соответствуют соседним ячейкам, то раскраска сетки является аналогичной для раскраски графа процедурой.
Рис. 6. Раскраска простой сетки
Очевидно, что при параллельном расчете ячеек выбранного цвета коллизий данных происходить не будет. Запустив расчет последовательно по ячейкам каждого цвета и обработав тем самым все ячейки, получим корректный результат.
Понятно, что ввиду последовательного запуска процедуры расчета по элементам каждого цвета для улучшения производительности необходимо уменьшить количество цветов в раскраске. Минимальное количество цветов среди всевозможных раскрасок называется хроматическим числом. Поиск одной из «минимальных раскрасок», т.е. раскрасок с хроматическим числом цветов, относится к ряду сверхсложных задач с экспоненциальным временем решения (зависимость от числа элементов). Поэтому раскраску будем искать приближенно с числом цветов, близким к хроматическому. К. Аппель и В. Хакен в 1976 году доказали известную проблему четырех красок, гласящую, что хроматическое число любого плоского графа не превосходит четырех. Так как любая двумерная сетка может быть отображена плоским графом, то данная теорема говорит о том, что минимальное число цветов двумерной сетки также не превосходит четырех. Существует ряд алгоритмов [7], обеспечивающих относительно быструю раскраску произвольного графа в заданное количество цветов, при этом согласно теореме о четырех красках, это количество разумно выбрать равным четырем. Самым простым и почти самым эффективным методом раскраски является жадный алгоритм с возвратами [6]. Примененный к структурированной сетке жадный алгоритм, к сожалению, с большой вероятностью даст сетку, раскрашенную в четыре цвета; такой результат не устраивает нас, поскольку хотелось бы в таком случае наблюдать очевидную минимальную «шахматную» раскраску с хроматическим числом, равным двум. Также данный алгоритм (как и другие известные алгоритмы) не заботится о количестве элементов каждого цвета, что может привести к результату, когда количество ячеек определенного цвета будет много меньше
количества ячеек других цветов. Такой расклад не сможет обеспечить расчет должной эффективностью. Поэтому в данной работе был разработан метод раскраски, базирующийся на жадном алгоритме, адаптированный к структурированным сеткам и производящий коррекцию раскраски, уравновешивающую количество ячеек каждого цвета. Результат применения разработанного алгоритма на некоторых сетках можно видеть на рис. 7.
а
б
Рис. 7. Пример раскраски с помощью разработанного алгоритма: а - неструктурированная сетка; б - структурированная сетка
Алгоритм реализуется в два этапа. На первом этапе все элементы обходятся по соседству с запоминанием «маршрута». При возникновении ситуации, когда у текущего элемента все соседи были пройдены, необходимо продолжить обход из последней ячейки маршрута, для которой это возможно. На втором этапе можно действовать обычным жадным алгоритмом раскраски, при этом выбор ячейки производя неслучайно, а по порядку следования в маршруте. Очевидно, что для структурированных сеток такой проход по первоначальному маршруту даст требуюмую раскраску в два цвета.
4. Разбиение на подобласти
Параллелизм первого уровня, основанный на межузловых взаимодействиях, происходит за счет разбиения расчетной области на подобласти (рис. 8). Каждый вычислительный объект (в данном случае GPU или ядро CPU) производит расчет своей подобласти.
Рис. 8. Пример разбиения области на подобласти
Такой расчет должен сопровождаться обменами данных между объектами, так как в общем случае подобласти не являются независимыми. Локальность связей и ограничение на соседство первого уровня между элементами сетки, однако, позволяет организовать синхронизацию данных за счет их обмена только из граничных ячеек, т.е. ячеек, хотя бы один сосед которых находится в другой подобласти. При практической реализации удобно достроить к подобласти дополнительный ряд ячеек, называемых ghost-ячейками, копирующими граничные ячейки соседних подобластей (рис. 9). Тогда синхронизация сводится к передаче информации от граничных ячеек подобластей и ее приему в ghost-ячейки соседних подобластей.
I
Рис. 9. Дополнение подобласти ghost-ячейками
5. Алгоритм расчета
Итерация расчетного цикла будет состоять из следующих шагов:
1. Выбираем цвет, который еще не рассчитан.
2. Рассчитываем граничные ячейки выбранного цвета во всех подобластях.
3. Синхронизация: передаем необходимые данные из посчитанных граничных ячеек в ghost-ячейки соседних блоков (можно выполнять асинхронно, т.е. в фоновом режиме).
4. Расчет внутренних ячеек выбранного цвета.
5. Выбранный цвет считается рассчитанным. Если все цвета рассчитаны, то итерация завершается, иначе переходим к пункту 1.
6. Распределение памяти
Память графических карт имеет многоуровневую систему. К сожалению, основная память GPU, называемая глобальной, латентна, что заставляет разработчиков предпринимать попытки уменьшить количество обращений к ней. Одним из способов, о котором рассказывалось ранее, является coalescing, т.е. объединение нескольких запросов к памяти в один. GPU аппаратно поддерживает несколько шаблонов оптимального доступа к глобальной памяти. Более подробно об этом написано в документации CUDA [8]. В данном пункте рассказывается о рекомендациях по быстрой и удобной организации памяти для осуществления эффекта coalescing и возможности быстрой межузловой синхронизации.
Как известно, для представления массивов структур в CUDA наилучшим образом подходит хранение в памяти в виде структуры массивов. Поэтому данные всех ячеек желательно хранить блочно покомпонентно. В таком случае хранение векторов состояния элементов сетки можно визуально представить рис. 10.
Компонента Расположение ячеек в памяти
1 2 3 4 N-1 N
1 2 3 4 N-1 N
Рис. 10. Покомпонентная организация памяти
Ввиду последовательного расчета по цветам ячеек для возникновения максимального количества coalesced-запросов к памяти, необходимо переупорядочить элементы в памяти по цветам. Так как порядок обхода ячеек внутри цвета неважен, то в таком случае, отдав соседним нитям соседние ячейки, можно добиться того, что все запросы к данным текущих обрабатываемых ячеек будут объединены. Схематично упорядоченные ячейки в блоке изображены на рис. 11.
Рис. 11. Упорядочивание ячеек по цвету
Согласно алгоритму расчета, граничные ячейки каждого цвета обрабатываются отдельно. Поэтому внутри подблока каждого цвета желательно выделить граничные ячейки и расположить их также линейным подблоком. Так как на каждой итерации будет происходить их пересылка между узлами, то такое линейное расположение также ускорит процесс синхронизации. Выделение ghost-ячеек в отдельный подблок позволит быстро принимать ячейки во время обменов. Такое распределение ячеек внутри подблока одного цвета одной компоненты можно представить (рис. 12).
Рис. 12. Упорядочивание ячеек по типу
Тогда рекомендованное расположение ячеек в памяти можно проиллюстрировать
схемой (рис. 13).
подблок "белых" ячеек
подблок "черных" ячеек
подблок "красных" ячеек
1 .....!......и 1 1 L ! i
V V V подблок подблок подблок внутренних приграничных дополнительных ячеек ячеек ячеек
ММ : :
: О г M ОГНЯ
1 M ! I M : |
Рис. 13. Рекомендованное схематичное представление данных в памяти
7. Заключение
В рамках данной работы разработана обобщенная программная модель для эффективного решения задач механики сплошных сред. Предложенная схема ведения параллельного расчета помогает добиться высокой производительности, используя особенности как программирования для систем с распределенной памятью, так и многопоточного программирования. В том числе уделено особое внимание актуальному направлению: расчету с использованием графических плат на основе технологии CUDA. На основе данной работы возможно создать высокоэффективную программную библиотеку, полностью охватывающую сложные параллельные процессы синхронизации, вопросы распределения памяти и хранения данных, оставляя своему пользователю лишь проблему программирования физической модели.
ЛИТЕРАТУРА
1. Антонов А.С. Параллельное программирование с использованием технологии MPI. - М.: МГУ, 2004.
2. Воеводин В.В. Параллельные вычисления. - СПб.: БХВ-Петербург, 2002.
3. Информационно-аналитический центр по параллельным вычислениям [Электронный ресурс]. URL: http://parallel.ru.
4. Men'shov I., Nakamura Y. Numerical Simulations and Experimental Comparisons for High Speed Nonequili-brium Air Flows -Fluid Dynamics Research Journal, v. 27, No 5, 2000, pp. 305-334.
5. Бьерн Страуструп. Язык программирования C++. Специальное издание / пер. с англ. - М.: Бином, 2011.
6. Кристофидес Н. Теория графов. Алгоритмический подход. - М.: Мир, 1978.
7. Родионов В.В. Методы четырехцветной раскраски вершин плоских графов. - М.: КомКнига, 2005.
8. NVIDIA CUDA C Programming Guide [Электронный ресурс]. URL: http://developer.nvidia.com.
EFFECTIVE PARALLEL COMPUTING ON HYBRID ARCHITECTURE SYSTEMS FOR PROBLEMS OF CONTINUUM MECHANICS
Klimov M.I., Menshov I.S.
A common core approach for effective parallel computing problems of continuum mechanics on hybrid architecture systems is proposed. The efficiency of GPGPU computations highly depends on several factors due to specifics of the architecture features. The proposed approach provides general means to develop parallel programs based on a correct and highperformance two-level parallel programming model that employs domain decomposition with message passing between computing nodes on the first level and element-by-element parallelism with multithread technologies on the second one.
Key words: parallel methods, high-performance computing, GPGPU, numerical methods.
Сведения об авторах
Климов Максим Игоревич, 1989 г.р., окончил МГУ им. М.В. Ломоносова (2012), аспирант МГУ им. М.В. Ломоносова, область научных интересов - численные методы, высокопроизводительные вычисления, параллельное программирование, механика сплошных сред.
Меньшов Игорь Станиславович, 1957 г.р., окончил МГУ им. М.В. Ломоносова (1979), доктор физико-математических наук, ведущий научный сотрудник ИПМ им. М.В. Келдыша РАН, профессор кафедры высшей математики МГТУГА и вычислительной механики МГУ им. М.В. Ломоносова, автор 47 научных работ, область научных интересов - численные методы, математическое моделирование, механика сплошных сред.
