The choice of technology for parallelization of numerical the solution of the convective-diffusion equation on hybrid mu.pdf В настоящее время при исследовании окружающей среды наряду с экспериментальными методами большое внимание уделяется применению методов математического моделирования. Именно математическое моделирование позволяет изучать и прогнозировать развитие атмосферных, гидрологических и геологических процессов на Земле без активного воздействия на них. Основу многих математических моделей механики сплошных сред, гидрометеорологии, экологии и других наук составляют неоднородные нестационарные трехмерные уравнения «конвекции-диффузии», которые описывают изменение таких параметров изучаемой среды, как плотность, компоненты скорости, температура, концентрация, турбулентность. Современные математические модели представляют собой систему из нескольких конвективно-диффузионных уравнений, дополненную замыкающими алгебраическими соотношениями, а также начальными и граничными условиями. Для решения рассматриваемых систем дифференциальных уравнений применяются методы конечных разностей, конечного объема или конечных элементов. Приближенное решение ищется на весьма подробной сетке, состоящей из нескольких миллионов вычислительных узлов, или ячеек. Причем часто при оперативном моделировании атмосферных или гидрологических процессов необходимо в кратчайшие сроки вычислений на компьютере представить результаты прогностического моделирования. Возможность ускорить процесс получения численного решения дает применение многопроцессорной и многоядерной вычислительной техники, снабженной графическими ускорителями [1, 2]. Разнообразие в архитектуре и компоновке высокопроизводительных вычислительных систем порождает разнообразие инструментов для построения параллельных программ. На сегодняшний день в области научных вычислений самыми популярными и узнаваемыми для написания параллельных программ являются технологии MPI, OpenMP и OpenACC. Это подталкивает исследователей к оценке перспективности использования разных подходов и их комбинаций для решения своих задач. В работе [3] на примере модели атмосферы UM (Unified Model UK Met Office) рассматривается эффективность использования гибридного подхода MPI + OpenMP в сравнении с MPI. Также исследуется эффективность выполнения параллельной программы при фиксированном количестве процессов, но при изменяющемся количестве задействованных ядер на каждом процессоре. Это обусловлено тем, что в вычислительных задачах особое внимание уделяется работе с памятью, и увеличение количества активных ядер в процессоре приводит к увеличению их конкуренции в борьбе за ресурсы процессора. Показано, что использование только технологии распараллеливания MPI оказывается эффективнее на всем диапазоне рассматриваемого количества процессов (до 1 536). В то же время использование только 12 из 16 ядер на каждом процессоре также позволяет повысить эффективность вычислений на 10-15%. Это происходит за счет увеличения пропускной способности памяти в расчете на одно активное ядро, сокращения издержек на коммуникации при уменьшении количества процессов MPI и активных ядер, увеличения частоты работы ядер в режиме turbo boost [3]. ВЕСТНИК ТОМСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА 2020 № 51 Управление, вычислительная техника и информатика В работе [4] рассматривается реализация блока переноса примесей в моделях атмосферы и океана на основе гибридного программирования MPI + OpenMP. В работе показано, что при небольшом количестве задействованных MPI-процессов более целесообразно дальнейшее наращивание их количества. В этом случае дополнительное ускорение достигается за счет размещения данных локальных массивов в более быстрой кэш-памяти процессора. При достаточно большом числе задействованных MPI-процессов, наоборот, более эффективным методом ускорения расчетов становится увеличение числа нитей OpenMP. В работе делается акцент на то, что гибридный подход MPI + OpenMP является наиболее универсальным с точки зрения использования ресурсов современных суперкомпьютеров, позволяя оптимизировать число MPI-процессов и число нитей OpenMP с учетом задачи и архитектуры вычислительных элементов. В работе [5] анализируются возможности повышения производительности вычислительных моделей горения на многоядерных компьютерах. Показано, что использование явного распределения данных по потокам оказывается более эффективным по сравнению с непосредственным применением директив OpenMP для автоматического распараллеливания циклов. В работе [6] представлены результаты исследования параллельных реализаций схем расщепления на основе технологий MPI и OpenMP для решения трехмерной задачи теплопереноса. Проводится сравнение следующих параллельных реализаций: MPI, MPI + OpenMP и MPI + OpenMP с выделением потоков-«почтальонов». Под MPI + OpenMP реализацией понимается непосредственно применение автоматического распараллеливания циклов внутри MPI-процессов. Основная идея гибридной MPI + OpenMP реализации с «почтальонами» заключается в выделении на каждом из узлов с общей памятью одного OpenMP потока, отвечающего за выполнение обменов данными между процессами на фоне вычислений. Результаты показали, что использование подхода с выделением потоков-«почтальонов» значительно повышает эффективность простой гибридной реализации, но уступает простой MPI реализации. Таким образом, проведенный анализ использования технологий распараллеливания MPI и OpenMP как отдельно, так и совместно показал, что наибольшей эффективности можно достичь при применении технологии MPI на большом количестве используемых процессов. При небольшом количестве многоядерных центральных процессоров на общей памяти разумнее использовать OpenMP технологию как более простую с точки зрения применения и не менее эффективную при небольшом числе задействованных вычислительных ядер. Постоянное развитие вычислительной техники (процессоров, сопроцессоров, графических процессоров), компиляторов и технологий параллельного программирования требует систематического проведения исследований по выбору наиболее эффективных технологий распараллеливания алгоритмов решения прикладных задач вычислительной математики. Целью данной работы является проведение вычислительных экспериментов на гибридной кластерной вычислительной системе по определению наиболее перспективного использования параллельных вычислительных технологий (MPI, OpenMP, OpenACC) для получения эффективных параллельных программ для решения нестационарного трехмерного неоднородного уравнения «конвекции-диффузии». Полученные результаты будут использованы при создании параллельной версии развиваемой в Томском государственном университете мезомасштабной метеорологической модели высокого разрешения TSUNM3 для прогнозирования опасных погодных явлений и качества атмосферного воздуха над городом [7]. 1. Математическая формулировка и численный метод 1.1. Математическая постановка задачи Рассмотрим обобщенное дифференциальное уравнение конвективно-диффузионного переноса: фФ 8pU Ф 8pV Ф 8pW Ф 8 С 8Ф^ 8 --'---'---'--= -I Kxy-I н- Л 8t 8x 8y 8z 8x l 8x J 8y l Лу 8y С 8Ф^ +ikf|+v (1) K Здесь t, x, y, z - время и координаты (ось Ох направлена на восток, Оу - на север, Oz - вертикально вверх); U, V, W- компоненты вектора скорости; р - плотность; Kxy, Kz - коэффициенты турбулентного обмена; Sф - источниковый член. В (1) могут рассматриваются в качестве Ф = 1 ((1) принимает вид уравнения неразрывности), Ф = U, V, W, ©, qv, qrain (компоненты скорости, потенциальная температура воздуха, абсолютная влажность воздуха и концентрация дождевых капель). Для учета рельефа подстилающей поверхности h(x, y) к уравнению (1) применяется преобразование вида: x' = x; y' = y; z'= IZ^^A.h . (-) H - h(x, y) В качестве граничных условий для (1) на подстилающей поверхности используются граничные условия третьего рода: z = h(x,y): Kz ^8Ф = а(Ф-Фw). (3) 8z На верхней границе - граничные условия второго рода: 8Ф z = H : K - = Уф . (4) 8z На боковых границах - граничные условия «радиационного» типа [8]: x = -L;x = L: ^ + CФ ^^ + C^ ^; (5) 8t 8x 8t 8x y = -L■A = L:™+Ci™=8ФL + cф,8ф:L. 8t 8y 8t 8y Здесь индекс S соответствует параметрам, определяемым из модели глобального масштаба [9]; C(x, C^y - фазовая скорость [8]; -L - размер области по горизонтали. Начальные условия для уравнения (1) определяются в результате интерполяции Фя, вычисляемой на более грубой сетке по модели глобального масштаба [9], на более подробную сетку мезомас-штабной метеорологической модели. 1.2. Численный метод решения задачи Для решения задачи (1)-(6) используются структурированные сетки с равномерным шагом по горизонтальным направлениям Ох и Оу и сгущающиеся к поверхности Земли по вертикальному. При аппроксимации дифференциальной постановки задачи (1)-(6) используется метод конечного объема со вторым порядком аппроксимации по пространственным переменным с применением оригинальной монотонизированной схемы для конвективных членов, построенной с использованием интерполяционных весовых сплайнов [10]. Рассматриваемые явно-неявные аппроксимации во времени (Адамса-Бэшфорда и Кранка-Никольсон) обеспечивали второй порядок точности при интегрировании уравнения (1) по времени. В итоге полученная разностная схема - аналог (1) - имела следующий вид: 8pUФ + 8pVФ + 8pWФ 8 Ik 8Ф^| 8 k ~k+1^isk+1 £ >isk о pP ФР -pP ФР VolP + 3 ж т - Volp С 8Ф^ - S„ K z - 8x 8y 8z 8x \\ 8x J 8y l 8y dpUo dpVo dpwф df do^ d --'---'----1 Kxy - dx dy dz dx I dx (7) k Щ K, do dz I dz - 2 jjj 2 VolP z - 2 f 2 VolP dxdydz = 0. k+1 4 K do dz I dz k-1 dO - 2 ж 2 VolP - SO dxdydz - K dy I dy Здесь k - номер слоя по времени; (xp, yp, zp) - координаты середины конечного объема; VoIp - значение конечного объема; здесь и далее нижний индекс p (point) указывает, что переменная определяется в центре конечного объема. Величина шага интегрирования по времени т выбиралась из условия устойчивости разностной схемы. Рассматриваемая разностная схема является экономичной и для вычисления значений сеточной функции ^Pk+1} позволяет использовать метод прогонки. Следует отметить, что неявная аппроксимация для вертикального диффузионного переноса, который имеет важное значение в пограничном слое атмосферы, позволяет уйти от более жесткого ограничения на шаг интегрирования по времени. Вычислительное ядро программы расчета значений сеточной функции ^Pk+1} на новом (k + 1)-м слое по времени представлено на рис. 1. Пунктирными линиями выделен блок, который может быть преобразован в подпрограмму TDMA(I, J) для вычисления {Ф^1} на отдельной вертикальной сеточной линии с помощью метода прогонки. Для оценки времени работы программной реализации (7) рассматривалась задача о переносе примеси от постоянно действующего источника в области с размерами L = 25 000 м и H = 600 м. Источник примеси расположен в точке с координатами (0,014L, 0,014L, 0,0101H). Компоненты скорости считались постоянными: U = 1 м/с, V= 1 м/с, W = 0 м/с. Коэффициенты диффузии также принимались постоянными Kxy = Kz = 100 м2/с. Шаг по времени выбирался равным 6 с, а время окончания моделирования TimeFin = 500т. Расчет проводился на сетке 256 х 256 х 32 на узлах вычислительного кластера ТГУ Cyberia со следующими характеристиками: 48Gb RAM, 2xIntel® Xeon® CPU 5670, 2,93 ГГц. Использовался компилятор Intel Fortran. Время счета последовательной программы (см. рис. 1) составило 107 с. Точность проведения расчета контролировалась по известному аналитическому решению рассматриваемой задачи. Time = 0 Do While (Time < TimeFin) Do I = 1,Nx Do J = 1,Ny !-----------------------TDMA(I,J)----------------------- P(0) = !из разностного аналога формулы (3) Q(0) = !из разностного аналога формулы (3) Do K = 1,Nz !расчет коэффициентов разностной схемы и !прогоночных коэффициентов P и Q из разностной схемы (7) End Do !цикл K F(Nz+1,I,J) = Q(Nz)/(1.0 - P(Nz)) Do K = Nz,0,-1 F(K,I,J) = P(K) * F(K+1,I,J) + Q(K) End Do !цикл K !задание боковых граничных условий из (5) и (6) !-----------------------TDMA(I,J)----------------------- End Do !цикл J End Do !цикл I F0(:,:,:) = F(:,:,:) !обмен значениями Time = Time + tau End Do !While Рис. 1. Вычислительное ядро численного метода (7) решения уравнения переноса Fig. 1. The computational kernel of the numerical method (7) of solving the transport equation 2. Параллельные алгоритмы для численного метода 2.1. Технология Message Passing Interface При применении этой технологии в качестве основного подхода распараллеливания для вычислительной системы с распределенной памятью выбрана двумерная (по осям Ох и Оу) декомпозиция сеточной области на подобласти (рис. 2). При такой декомпозиции каждому процессорному элементу (ПЭ), задействованному в вычислениях, распределяются значения сеточной функции {Фр} принадлежащие локальной подобласти, для вычисления значений на следующем слое по времени {Фр^1} с использованием метода прогонки. Из-за используемого шаблона явно-неявной разностной схемы (7) (см. рис. 2, выделено в нижнем правом углу) для вычисления очередного приближения в приграничных узлах каждой локальной подобласти необходимо знать значения сеточной функции {Фр4} соседней граничащей подобласти. Для этого создаются фиктивные ячейки (на рис. 2 светлые значки) для хранения данных с процессорного элемента, рассчитывающего {Фр^1} в соседней подобласти декомпозиции, и организуются пересылки этих граничных значений, необходимых для обеспечения однородности вычислений в каждой локальной подобласти [11]. Для передачи данных другим процессорным элементам и получения необходимых для продолжения вычислений данных от них в данной работе используется интерфейс передачи сообщений MPI. Рис. 2. Двумерная декомпозиция сеточной области вдоль направлений Ох и Оу и организация передачи данных между подобластями Fig. 2. Two-dimensional decomposition of the grid area along the Ox and Oy directions and the organization of data transfer between the subdomains При разработке параллельной программы дополнительно используется ряд функций MPI для создания и тестирования двумерной декартовой логической топологии, создания новых типов данных для обеспечения эффективной пересылки между процессорами, а также функция SENDRECV (в процедуре MPE_SEND_F) для проведения парных блокирующих пересылок значений сеточной функции в соседних подобластях декомпозиции (рис. 3). Time = 0 Do While (Time < TimeFin) Do I = 2,Px-1 !- число узлов в подобласти по х (Рис.2) Do J = 2,Py-1 !- число узлов в подобласти по у (Рис.2) CALL TDMA(I,J) End Do !цикл J End Do !цикл I Call MPE_SEND_F !обмен значений между ПЭ F0(:,:,:) = F(:,:,:) !обмен значениями между слоями k и k+1 Time = Time + tau End Do !While Рис. 3. Фрагмент параллельной части MPI-программы для решения уравнения переноса Fig. 3. Fragment of the parallel part of the MPI program for solving the transport equation Проведенные на вычислительном кластере ТГУ Cyberia (использовались узлы с указанными выше характеристиками) расчеты показали достаточно высокую эффективность (не ниже 70% на 128 процессах) выбранного способа распараллеливания явно-неявного вычислительного алгоритма (рис. 4, светлые значки). Для сравнения на этом же рисунке приведены значения ускорения и эффективности, полученные с использованием технологии OpenMP (см. п. 2.2). Значения ускорения рассчитывались по формуле Sp = TJTp, а эффективность параллельной программы вычислялась как Ep = Sp/p. Здесь T1, Tp - время счета параллельной программы на 1 и p процессах. Из рис. 4 также видно, что время счета MPI-программы снижается практически линейно. При количестве процессов параллельной программы p = 128 время счета снижается в 90 раз по сравнению с временем счета последовательной программы. « 10 Число процессов р 10 100 Число процессов р 1000 100 10 Число процессов р 1000 1 -г 100 Рис. 4 Время счета, ускорение и эффективность MPI (o) и OpenMP (•) программ. Штриховая линия представляет график идеального ускорения параллельных программ Fig. 4 Time of calculation, acceleration and efficiency of MPI (o) and OpenMP (•) programs. The dashed line represents the graph of the ideal acceleration of parallel programs Этот результат чрезвычайно важен при разработке программ для проведения прогностических расчетов, когда требуется за минимальное время выполнить прогноз опасного погодного явления (ливень, интенсивный снегопад, заморозки, град и т.п.), тем самым повышая уровень заблаговременной предупрежденности и оставляя больше времени для подготовки и снижения последствий. 2.2. Технология Open Multi-Processing (OpenMP) Данная технология распараллеливания применяется для многопроцессорных (многоядерных) вычислительных систем с общей оперативной памятью. Основная задача, которая стоит перед программистом, выбравшим эту технологию распараллеливания, заключается в распределении вычислений между параллельно выполняемыми потоками обработки общих разделяемых данных. Это реализуется с помощью специальных директив препроцессора используемого языка программирования, вызова библиотечных процедур и задания переменных окружения [12]. В рассматриваемом случае при создании параллельного многопоточного приложения особое внимание нужно уделить основному вычислительному блоку (см. рис. 1), в котором производится перевычисление значений сеточной функции {Фр^1} на каждом следующем слое по времени, а также распределению вычислений, выполняемых для каждого фиксированного значения переменных циклов (I, J), между задействованными в расчетах многоядерными процессорами. Для этого в последовательную версию программы в соответствующих местах нужно разместить соответствующие директивы OpenMP (рис. 5). Заметим, что трудоемкость преобразования последовательной программы в параллельную с помощью технологии OpenMP значительно ниже, чем при разработке MPI параллельной программы. Более того, для того чтобы размещенные в OpenMP-программе директивы заработали, нужно при компиляции программы использовать некоторые дополнительные опции. Без них директивы OpenMP при компиляции будут игнорироваться. Time = 0 Do While (Time < TimeFin) $OMP PARALLEL PRIVATE(I,J) $OMP DO Do I = 1,Nx Do J = 1,Ny CALL TDMA(I,J) End Do !цикл J End Do !цикл I $OMP END DO !обмен значениями $OMP DO Do I = 0,Nx+1 Do J = 0,Ny+1 Do K = 0,Nz+1 F0(K,I,J) = F(K,I,J) End Do !цикл K End Do !цикл J End Do !цикл I $OMP END DO $OMP END PARALLEL Time = Time + tau End Do !While Рис. 5. Фрагмент параллельной части OpenMP-программы для решения уравнения переноса Fig. 5. Fragment of the parallel part of the OpenMP program for solving the transport equation При проведении расчетов на одном вычислительном узле кластера ТГУ Cyberia с двумя процессорами Intel Xeon 5670, каждый из которых имеет 6 ядер, были получены следующие результаты ускорения и эффективности разработанной OpenMP-программы (см. рис. 4). Из рисунка видно, что в относительных величинах ускорения Sp и эффективности Ep, а также по возможности привлечения других вычислительных узлов кластера OpenMP-программа уступает MPI-программе, однако в абсолютных значениях время работы OpenMP-программы для каждого параллельного запуска несколько меньше (различие уменьшается с ростом количества процессов). Тем не менее с использованием относительно простой технологии параллельного программирования OpenMP на одном вычислительном узле кластера Cyberia можно получить 10-кратное ускорение параллельной программы по сравнению с последовательной. 2.3. Технология Open Accelerators (OpenACC) Данная технология параллельного программирования используется для создания параллельных программ на гетерогенных вычислительных системах с центральными и графическими процессорами. Эта технология параллельного программирования, как и технология OpenMP, с помощью директив компилятора позволяет отмечать участки программы, которые следует выполнять параллельно на графическом процессоре, обозначает общие и локальные переменные [13]. Здесь также при параллельном программировании большое внимание уделяется распределению циклически выполняющихся вычислений графическому процессору (рис. 6). Time = 0 !$ACC DATA COPYIN(FO) COPYOUT (F) Do While (Time

Скачать электронную версию публикации