Параллельная реализация метода естественных соседей

Существует широкий круг задач механики жидкости, требующих для своего решения вычислительных ресурсов, значительно больших, чем может предоставить персональный компьютер. Многим из этих задач необходимо не только высокое быстродействие, но также обработка и хранение большого объема информации, что предъявляет повышенные требования к оперативной памяти. Использование для решения таких задач высокопроизводительных систем и систем с параллельной архитектурой значительно расширяет возможности исследователей, занимающихся компьютерным моделированием сложных физических процессов. Решение задачи в параллельном режиме требует правильного выбора модели распараллеливания программного кода, тесно связанной с аппаратными и системными средствами параллельной обработки данных. В настоящей работе рассматривается ряд вопросов по моделированию задач вязкой и идеальной жидкости методами естественных соседей (Natural element method) на многопроцессорных вычислительных системах с распределенной памятью. Обсуждаются методы блочной декомпозиции матрицы системы, способы распределения данных по процессорам, а также особенности параллельной реализации численных методов, применяемых для решения подзадач. На примере задачи о колебаниях жидкости в бассейне исследуется производительность программного кода в зависимости от размера задачи и числа процессоров. 1. Метод естественных соседей и его распараллеливание Метод естественных соседей (Natural element method, NEM) предложен Л. Траверсони в 1994 г. для решения задач теории упругости. Данный метод представляет собой разновидность бессеточного метода Галеркина. При формировании дискретной системы © Институт вычислительных технологий Сибирского отделения Российской академии наук, 2008. 56 Параллельная реализация метода естественных соседей… 57 уравнений используется метод взвешенных невязок. Для интерполяции значений неизвестных используются функции формы Сибсона и Лапласа, основанные на диаграммах Вороного первого и второго порядков. В работе метод естественных соседей адаптирован для решения задач идеальной (NEM) и вязкой (General NEM) несжимаемой жидкостей с большими деформациями расчетной области. Течения вязкой и идельной несжимаемой жидкости описываются системами уравнений Навье — Стокса (1) и Эйлера (2) соответственно, а также уравнением неразрывности (3): Du 2„- — Vp + pX/zu + pF; (1) (2) (3) Dt p- = — Vp + pF-Vw = 0. РВ случае идеальной жидкости на твердых границах задается условие непротекания йп = 0, на свободных — р = patm; для вязкой жидкости условие непротекания заменяется условием прилипания и = 0, а на свободных границах задаются плотности поверхностных сил рп = t^Uj — рщ = ani. Здесь n^, rij — компоненты вектора внешней нормали к твердой границе; т^- — тензор вязких напряжений; / = (/ь/г) — вектор внешних сил; р — плотность жидкости. На рисунке приведена блок-схема метода GNEM. Для решения системы уравнений (1), (3) и (2), (3) используется метод дробных шагов. После ввода начальных данных строится диаграмма Вороного и определяются границы расчетной области (шаг 2). Затем вычисляется предиктор скорости u* = (u*,v*) (шаг 3), из решения уравнения Пуассона вычисляется давление, которое используется для определения гидродинамических нагрузок (шаг 4), и определяются компоненты вектора скорости и, v на текущем временном слое (шаг 5). Шаг 6 — стабилизация функции давления методом конечных приращений — выполняется в методе GNEM для моделирования вязких течений. При переходе к дискретной форме уравнений (шаги 3-6) формируются системы линейных алгебраических уравнений (СЛАУ), которые решаются методом сопряженных градиентов. Для исследования производительности приведенного алгоритма методов NEM и GNEM рассматривалась задача о колебаниях жидкости в прямоугольном бассейне. В табл. 1 приведены временные затраты на выполнение основных операций. На одном шаге по времени формируется несколько однотипных блоков, которые включают численное интегрирование, сбор матрицы, внедрение граничных условий, решение полученной СЛАУ. В рамках модели вязкой жидкости число таких блоков равно семи, в случае идеальной — трем. В табл. 1 приводится суммарное время реализации блоков Input x(0), u(0),v(0) щCalculate Voronoi, a-shape re- 5»Calculate predictor u*, v* ft Calculatepressure /1 — P [4 |t = t + dt| ¦ calculate load . PsP7«Ј Stabilize pressure \6«Ј Calculate velocity stop с 5 Блок-схема расчета одного шага по времени методом естественных соседей 58 , Таблица 1. Временные затраты на реализацию блоков последовательной программы, с. Число узлов, NДискретизация областиСбор матрицы (вязк./идеал.)Внедрение граничных условий (вязк./идеал.)Решение СЛАУ (вязк./идеал.) 12000.051.21 / 0.510.15 / 0.061.08 / 0.40 29300.174.74 / 1.860.56 / 0.235.77 / 2.43 41200.317.39 / 3.122.24 / 0.9910.25 / 4.29 62080.7620.01 / 8.824.55 / 2.0146.24 / 19.96 одного временного шага. Из данных табл. 1 видно, что при увеличении числа расчетных узлов в области процентное соотношение времени сбора матрицы к общему времени одного временного шага уменьшается с 48 до 28 %, а соотношение времени решения СЛАУ увеличивается с 43 до 65 %. Внедрение граничных условий занимает 5-12 %, а дискретизация области (построение диаграммы Вороного и определение границ) — 1-2 %. Параллельная программа представляет собой систему процессов, взаимодействующих посредством передачи сообщений (Message Passing Interface, MPI) . В рамках данной модели для создания эффективного приложения необходимо снижать долю последовательных вычислений и минимизировать пересылки данных между процессорами. Подход к распараллеливанию NEM и GNEM основан на декомпозиции матрицы СЛАУ на блоки, число которых равняется числу процессоров, расчете каждым процессором своего блока, обмене данными между процессорами на каждом шаге по времени. В результате достигается равномерная загрузка данными и вычислительными операциями каждого процессора. Дополнительным преимуществом подобной декомпозиции является значительная экономия оперативной памяти: при построении СЛАУ выделяется матрица размером N х N/np, где пр — число узлов кластера. Вычисление нового положения узлов осуществляется каждым процессором в пределах своей подобласти. Затем производится обмен векторов координат расчетных узлов для синхронизации данных. Одним из трудоемких шагов в NEM и GNEM является вычисление функций формы Сибсона и Лапласа и построение СЛАУ. Матрицы и векторы правой части формируются в результате численного интегрирования уравнений (1), (3) для задач вязкой и (1), (2) для задач идеальной жидкостей. После декомпозиции области операции по вычислению интегралов осуществляются только над теми узлами, которые хранятся на данном процессоре. Каждая результирующая строка матрицы строится независимо от других, при этом полностью отсутствуют коммуникационные операции. Таким образом, степень параллелизма данного участка алгоритма является максимальной. 2. Решение СЛАУ Важное значение имеют распределение данных по процессорам, синхронизация вычислений и маршрутизация данных, неудачная организация которых приводит к неверному счету или большим задержкам, связанным с ожиданием и пересылкой данных и бездействием процессоров. Дополнительно организация данных на каждом процессоре должна соответствовать схеме алгоритма численного метода, применяемого для решения СЛАУ Полученная система линейных алгебраических уравнений решается методом сопряженных градиентов (МСГ) . Распараллеливание МСГ осуществляется следующим Параллельная реализация метода естественных соседей… 59 образом. Используется одномерная декомпозиция по столбцам, т. е. каждому процессору выделяется определенное число столбцов, для которых он будет производить вычисления. Алгоритм метода состоит из операций умножения матрицы на вектор и вычисления скалярных произведений, которые достаточно легко распараллеливаются. Структура массивов на каждом процессоре без труда позволяет организовать вычисления так, что обмен требуется лишь для пересылки вектора неизвестных на каждой итерации. 3. Тестирование Для количественной оценки эффективности распараллеливания используются следующие характеристики. 1. . Ускорение Sp = Т\/Тр, где Т\ — время выполнения программы на одном процессоре, Тр — время выполнения программы в системе из р процессоров. 2. . Эффективность Ер = Sp/p . Для определения эффективности и ускорения реализованного параллельного алгоритма метода естественных соседей проведена серия расчетов на учебном кластере кафедры ЮНЕСКО по НИТ Кемеровского государственного университета и СКИФ Таблица 2. Характеристики используемого кластера ХарактеристикаКластер кафедры ЮНЕСКО по НИТСКИФ Cyberia, ТГУ Процессорная архитектурах86×86 с поддержкой 64 Кол-во узлов / проц.8/8286 / 566 (1132 ядра) ПроцессорPentium III, 1 ГГцIntel Xeon 5150, 2.66 ГГц Суммарный объем ОЗУ2 Гб1136 Гб Суммарный объем ПЗУ160 Гб22 Тб Коммуникационная сетьFast EthernetInfiniBand Операционная системаRed Hat Linux 9.0SUSE Enterprise Server 10.0 КомпиляторFortran Intel 9.0 + MPI 2Fortran Intel 9.1 + MPI 2 Таблица З. Ускорение и эффективность на кластере КемГУ УскорениеЭффективность пр248248 Nвязкидвязкидвязкидвязкидвязкидвязкид 12000.530.750.440.630.380.540.260.370.110.150.040.06 41200.861.090.861.090.851.070.430.540.210.270.10.13 62081.651.961.722.051.72.070.820.980.430.510.210.25 Таблица 4. Ускорение на кластере СКИФ Cyberia пр2481632128 Nвязкидвязкидвязкидвязкидвязкидвязкид 12001.151.141.181.171.171.161.161.151.141.130.980.97 41201.361.341.431.41.611.571.691.631.651.611.661.61 62081.371.351.481.451.71.651.871.831.761.71.821.7 10 1601.441.421.551.511.991.922.212.122.312.212.292.2 60 , Таблица 5. Эффективность на кластере СКИФ Cyberia пр2481632128 Nвязкидвязкидвязкидвязкидвязкидвязкид 12000.570.570.290.290.140.140.070.070.030.030.0070.007 41200.680.670.350.350.020.190.10.10.050.050.0130.012 62080.680.670.370.360.210.20.120.110.060.050.0140.013 10 1600.720.710.380.370.250.240.140.130.070.070.0170.017 Cyberia Томского государственного униве ситета. В табл. 2 приведены основные характеристики указанных кластеров. В табл. 3 приведены значения ускорения и эффективности параллельного алгоритма, полученные на кластере КемГУ. В табл. 4 и 5 приведены значения ускорения и эффективности параллельного алгоритма для различного числа процессоров на кластере СКИФ Cyberia при решении одного шага по времени в задачах о колебании вязкой и идеальной несжимаемой жидкостей в бассейне. Заключение Проведено тестирование программного кода и его составных частей на многопроцессорных вычислительных системах. Вследствие разделения ресурсов переход от однопроцессорной конфигурации к многопроцессорной системе из пр процессоров не приводит к сокращению времени счета в пр раз. С ростом числа процессоров эффективность уменьшается. В результате тестирования параллельного алгоритма установлено, что наилучшее значение ускорения достигается в тех случаях, когда на каждом узле кластера хранится не менее 2500-3000 расчетных узлов. Данная зависимость объясняется особенностью МСГ, который занимает до 65 % процессорного времени: при малом числе расчетных узлов на каждом процессоре пересылка данных занимает гораздо большее время, чем непосредственные вычисления с ними. Распараллеливание алгоритма методов NEM и GNEM позволяет проводить расчеты для большого числа узлов за счет разделяемой памяти. Sukumar N., Moran В., Belytschko T. The natural element method in solid mechanics // Intern. J. Num. Methods Eng. 1998. Vol. 43, N 5. P. 839-887. , , , Метод естественных соседей на основе интерполяции Сибсона // Вест. Том. гос. ун-та. 2006. № 19. С. 210-219. Многопроцессорные вычислительные системы и параллельное программирование / , СВ. Стуколов, , . Кемерово: Кузбассвуз-издат, 2003. Saad Y. A basic tool kit for sparse matrix computation // Tech. Rep. CSRD TR 1029, CSRD, Univ. of Illinois, Urbana, IL., 1990.