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

На правах рукописи

ГОРОБЕЦ Андрей Владимирович

ПАРАЛЛЕЛЬНЫЕ ТЕХНОЛОГИИ МАТЕМАТИЧЕСКОГО МОДЕЛИРОВАНИЯ ТУРБУЛЕНТНЫХ ТЕЧЕНИЙ НА СОВРЕМЕННЫХ СУПЕРКОМПЬЮТЕРАХ

Специальность 05.13.18 - математическое моделирование, численные методы и комплексы программ

АВТОРЕФЕРАТ диссертации на соискание ученой степени

доктора физико-математических наук 1 |-|0Я 2015

Москва - 2015

005564519

Работа выполнена в Федеральном государственном бюджетном учреждении науки Институте прикладной математики им. М. В. Келдыша РАН

Научный консультант: доктор физико-математических наук

Татьяна Константиновна Козубская

Официальные оппоненты: Воеводин Владимир Валентинович

доктор физико-математических наук, профессор, член-корр. РАН, НИВЦ МГУ, заместитель директора

Волков Андрей Викторович

доктор физико-математических наук, доцент,

ФГУП «ЦАГИ», начальник отделения

Дерюгин Юрий Николаевич доктор физико-математических наук, старший научный сотрудник, РФЯЦ-ВНИИЭФ, начальник отдела

Ведущая организация: Федеральное государственное бюджетное

учреждение науки Вычислительный центр им. А. А. Дородницына Российской академии наук

Защита диссертации состоится «3» декабря 2015 г. в_на заседании

диссертационного совета № Д 002.024.03 при ИПМ им. М. В. Келдыша РАН по ' адресу 125047, г. Москва, Миусская пл., 4.

С диссертацией можно ознакомиться в библиотеке ИПМ им. М.В.Келдыша РАН, а также на сайте: http://keldysh.ni/council/3/D00202403/defence2.htm

Автореферат разослан «_»_2015 г.

Ученый секретарь Диссертационного совета Д 002.024.03

кандидат физико-математических наук м- А- Корнилина

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

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

Непрерывный рост производительности суперкомпьютеров открывает все более широкие возможности перед вычислительным экспериментом. Уже существуют системы, обладающие пиковой производительностью в несколько десятков PFLOPS (1 PFLOPS = 1015 вычислительных операций в секунду). Но если в недалеком прошлом рост производительности поддерживался в основном за счёт увеличения числа процессоров в системе, то теперь рост сопряжен с существенным усложнением архитектур, программных моделей и их разнообразием. Это делает создание расчетных кодов, в полной мере использующих возможности современных вычислительных систем, сложной и актуальной научной проблемой.

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

С другой стороны, рост производительности достигается за счет использования массивно-параллельных ускорителей. К таким ускорителям относятся, в частности, графические процессоры GPU (Graphics Processing Unit) и ускорители Intel Xeon Phi существенно-многоядерной архитектуры MIC

3

(Many Integrated Core). Среди десяти самых мощных суперкомпьютеров мира уже можно видеть несколько систем такой гибридной архитектуры (Tianhe-2, Tianhe-1, Cray Titan, Intel Stampede). Гетерогенные вычисления, в частности, с использованием графических процессоров, сами по себе являются новой областью, появившейся всего несколько лет назад. Средства программирования и технологии вычислений находятся в настоящее время в активном развитии, постоянно дополняется функциональность и возможности низкоуровневых интерфейсов программирования CUD A, OpenCL, и высокоуровневых директивных средств ОрепАСС и OpenMP 4.0. Несмотря на все возрастающую популярность гетерогенных вычислений, развитие газодинамических алгоритмов, ориентированных на их использование, все еще далеко от стадии зрелости, активно ведутся исследования в этой области.

Таким образом, разработка параллельных алгоритмов, которые могут одновременно удовлетворять всё возрастающим требованиям по степени параллелизма и масштабируемости, и быть адаптированными к еще более сложной параллельной модели, является крайне актуальной научной проблемой. Также важной и актуальной представляется разработка параллельных комплексов программ, позволяющих проводить расчеты широкого круга задач газовой динамики, аэроакустики, тепломассопереноса с использованием десятков тысяч процессорных ядер и различных типов ускорителей, включая архитектуры Intel, AMD и NVIDIA. Только на основе эффективных алгоритмов и программных комплексов, отвечающих современным тенденциям в эволюции суперкомпьютеров, можно достичь прогресса в математическом моделировании турбулентных течений.

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

Цели и задачи диссертационной работы

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

4

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

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

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

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

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

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

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

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

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

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

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

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

Методы исследования. В данной работе вычислительный эксперимент является методом исследования турбулентного течения. Для численного исследования применяются конечно-объемные и конечно-разностные методы повышенного порядка аппроксимации. Для моделирования турбулентности используются нестационарные вихреразрешающие подходы LES (Large eddie simulation - моделирование крупных вихрей), семейство гибридных подходов DES (Detached eddie simulation - моделирование отсоединенных вихрей), и метод прямого численного моделирования DNS (Direct numerical simulation -прямое численное моделирование). Параллельные алгоритмы и программные комплексы, реализующие численные методы, основываются на объектно-ориентированном подходе и многоуровневой параллельной модели, сочетающей различные типы параллелизма.

Положения, выносимые на защиту

• Многоуровневый параллельный алгоритм расчета сжимаемых турбулентных течений с учетом явлений аэроакустики.

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

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

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

• Результаты расчётов по моделированию аэродинамического шума от турбулентного следа и от взаимодействия турбулентности с твердым телом.

• Результаты расчетов ряда фундаментальных задач по моделированию несжимаемых турбулентных течений и тепломассопереноса при естественной конвекции.

• Результаты расчетов ряда фундаментальных задач по моделированию несжимаемых турбулентных течений при вынужденной конвекции.

Научная новизна

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

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

• Предложен новый алгоритм с многоуровневым распараллеливанием для расчетов сжимаемых течений на основе экономичных EBR схем повышенного порядка на неструктурированных сетках (Abalakin I.V., Bakhvalov P.A., Kozubskaya Т.К. Edge-based reconstruction schemes for prediction of near field flow region in complex aeroacoustic problems // Int. J. Aeroacoust. 2014. V.13. p. 207-234). Параллельный алгоритм рассчитан на системы с сотнями тысяч ядер и адаптирован к использованию гибридных суперкомпьютеров с ускорителями Intel Xeon Phi.

• Представлен новый комплекс программ NOISEtte для математического моделирования сжимаемых турбулентных течений методами повышенного порядка точности на неструктурированных сетках. В нем впервые создана параллельная реализация экономичных численных схем семейства EBR.

• Предложен оригинальный параллельный метод решения уравнения Пуассона на гибридных системах с массивно-параллельными ускорителями. Новый метод за счет сочетания прямых и итерационных методов решения СЛАУ наиболее эффективен для моделирования несжимаемых турбулентных течений в задачах с одним периодическим направлением. Метод также может применяться для расчетов трехмерных задач без периодического направления, благодаря новому расширению на основе многосеточного подхода. Для задач с пространственной симметрией предложен новый способ расширения границ масштабирования.

• На основе этого метода создан новый параллельный алгоритм и реализующий его программный комплекс для математического моделирования несжимаемых турбулентных течений на традиционных суперкомпьютерах с числом процессоров порядка сотни тысяч и на гибридных суперкомпьютерах с различными типами ускорителей (GPU NVIDIA, GPU AMD, Intel Xeon Phi)

• Выполнено численное исследование механизмов генерации шума от турбулентного следа и взаимодействия турбулентных структур с твердым телом на задаче об обтекании тандема квадратных цилиндров. Получены новые данные о применимости гибридных RANS-LES подходов на неструктурированных сетках к математическому моделированию

аэродинамического шума и о качестве численного решения в сравнении с экспериментальными данными.

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

• Выполнена серия DNS расчетов конвекции Рэлея-Бенара на сетках с числом узлов до 600 миллионов. Получены новые численные результаты и подробные данные о топологии течения.

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

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

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

Достоверность результатов. Разработанные параллельные комплексы программ надежно и тщательно верифицированы на широком круге задач.

Выполнено сравнение с аналитическим решением на модельных задачах, сравнение с экспериментальными данными, сравнение с численными результатами других авторов. Корректность реализации и порядки сходимости дискретных операторов подтверждены, в том числе, на основе широко известного метода MMS (Method of Manufactured Solutions). Эффективность и производительность параллельных вычислений подтверждается серией тестов, выполненных на многопроцессорных системах различных архитектур с использованием до 24000 процессоров.

Основные публикации. По теме диссертации в печати опубликовано 28 работ в журналах, рекомендованных ВАК для опубликования научных результатов докторских диссертаций, включая 21 статью в международных журналах, входящих в реферативную базу Scopus. Полный список публикаций приводится в конце автореферата.

Апробация работы. Результаты, входящие в данную диссертационную работу, были представлены в 70-ти докладах на конференциях (более 50 международных). В том числе в докладах:

• А.В. Горобец, С.А. Суков, П.Б. Богданов, Ф.Х. Триас, Конечно-объемные алгоритмы для моделирования турбулентных течений на гибридных суперкомпьютерах различной архитектуры, XV международная конференция «Супервычисления и математическое моделирование», 13-17 октября 2014, г. Саров.

• Gorobets, F.X. Trias and A. Oliva, Fighting against massively parallel accelerators of various architectures for the efficiency of finite-volume parallel CFD codes, 26 Parallel CFD, 2014, 20-22 of May, Trondheim, Norway.

• Sergey Soukov, Andrey Gorobets and Pavel Bogdanov, OpenCL Implementation of Basic Operations for a High-order Finite-volume Polynomial Scheme on Unstructured Hybrid Meshes, Parallel CFD, 2013, May 20-24, Changsha, China.

• Andrey Gorobets, Francesc Xavier Trias Miquel and Assensi Oliva, An OpenCL-based Parallel CFD Code for Simulations on Hybrid Systems with Massively-parallel Accelerators, Parallel CFD, 2013, May 20-24, Changsha, China.

• И.А. Абапакин, П.А. Бахвалов, A.B. Горобец, А.П. Дубень, Т.К. Козубская, Комплекс программ NOISETTE для моделирования задач аэродинамики и аэроакустики, XXIV Научно-техническая конференция по аэродинамике, п. Володарского, 28 февраля - 1 марта 2013 г.

• А. В. Горобец, С. А. Суков, На пути к освоению гетерогенных супервычислений в газовой динамике, Первый Национальный Суперкомпьютерный Форум (НСКФ-2012), Россия, Переславль-Залесский, ИПС имени А.К. Айламазяна РАН, 29-30 ноября 2012 года.

• Andrey Gorobets, F. Xavier Trias, Assensi Oliva, A parallel OpenCL-based solver for large-scale DNS of incompressible flows on heterogenous systems, Parallel CFD 2012, Atlanta, USA, May 2012.

• Gorobets, F. X. Trias, R. Borrell, M. Soria and A. Oliva, Hybrid MPI+OpenMP parallelization of an FFT-based 3D Poisson solver that can reach 100000 CPU cores, Parallel CFD 2011, Barcelona, Spain, 16-20 May 2011.

• V. Gorobets, S. A. Soukov, P. B. Bogdanov, A. O. Zheleznyakov and B. N. Chetverushkin, Extension with OpenCL of the two-level MPI+OpenMP parallelization for large-scale CFD simulations on heterogeneous systems, Parallel CFD 2011, Barcelona, Spain, 16-20 May 2011.

• A.B.Горобец, C.A.Суков, А.О.Железняков, П.Б.Богданов, Б.Н.Четверушкин, Применение GPU в рамках гибридного двухуровневого распараллеливания MPI+OpenMP на гетерогенных вычислительных системах, Параллельные вычислительные технологии (ПаВТ) 2011, 28 марта - 1 апреля, Москва.

Реализация и внедрение результатов работы. Работа выполнялась в рамках научных планов ИПМ им. М. В. Келдыша РАН и поддерживалась грантами Российского фонда фундаментальных исследований, проектами Министерства образования и науки РФ, проектами совета по грантам при президенте РФ. Методики, методы, алгоритмы, программные средства, результаты расчетов, представляемые в работе к защите, использовались в проектах, совместных научных исследованиях и договорных работах со следующими организациями: ЦАГИ, ЦНИИМаш, ОАО «Авиадвигатель», ОАО «Камов», ОАО «ОКБ Сухого», МФТИ, РФЯЦ-ВНИИЭФ, НИИСИ РАН, ИБРАЭ РАН, Санкт-Петербургский государственный политехнический университет, Технический университет Каталонии (Испания), Университет Гронингена (Нидерланды), Termo Fluids S.L. (Испания), а также в проектах РНФ, Минобрнауки и в проекте 7-й рамочной программы Евросоюза VALIANT.

Структура и объем работы. Диссертация состоит из введения, 6 глав, заключения и списка литературы. Объем составляет 226 машинописных страниц, текст содержит 131 рисунок и 11 таблиц.

СОДЕРЖАНИЕ РАБОТЫ

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

Первая глава посвящена технологии математического моделирования турбулентного течения посредством крупномасштабных вычислительных экспериментов с использованием конечно-объемных и конечно-разностных численных методов повышенной точности. Технология представляет собой свод согласованных между собой подходов для эффективной программной реализации алгоритмов вычислительной газовой динамики, рассчитанной на современные суперкомпьютеры, в том числе гибридные. Результаты данной главы опубликованы в [1] и по вычислениям на ускорителях в [3-5,7,8,11].

В § 1.1 рассматривается многоуровневая модель параллельных вычислений и проблема выбора средств разработки для различных ее уровней. Для соответствия современным архитектурам суперкомпьютеров предлагается использовать гибридную многоуровневую параллельную модель, сочетающую модели с общей и распределенной памятью, параллелизм M1MD (Multiple Instruction, Multiple Data) и SIMD (Single Instruction, Multiple Data). Уровни упорядочены следующим образом: 1) распределенная память, MIMD; 2) общая память, MIMD; 3) потоковая обработка (Stream processing); 4) векторизация, SIMD. В качестве средства разработки на первом уровне, объединяющем узлы суперкомпьютера, выбран интерфейс прикладного программирования MPI (Message passing interface). На втором уровне для распараллеливания с общей памятью выбран ОрепМР - открытый стандарт, поддерживаемый основными компиляторами языка С++. На третьем уровне предлагается использовать открытый вычислительный стандарт OpenCL, поддерживаемый основными производителями процессоров и ускорителей, включая NVIDIA, AMD, Intel.

В §1.2 предлагается разделить алгоритм на составляющие его базовые операции, выделить реализации этих операций из основного расчетного кода и снабдить каждую операцию тестовой оболочкой для быстрого многократного запуска, профилирования и проверки корректности вычислений. Каждая базовая операция реализуется в трех вариантах: 1) вариант для CPU в исходной

11

структуре данных; 2) переходный вариант, работающий на CPU, но в структуре данных, адаптированной к ускорителю и потоковой обработке; 3) вариант для ускорителя на OpenCL. В итоге существенная часть отладки остается на CPU, разработка значительно упрощается. Конечная версия для ускорителя на основе переходной реализации создается с существенно меньшими трудозатратами. Корректность вычислений гарантируется сравнением выходных данных с исходной реализацией. Далее, работая с операциями по отдельности в составе тестовой оболочки, реализующей проверку корректности, быстрый запуск с разными наборами входных данных, точное профилирование, выполняется оптимизация вычислений. Такой подход упрощает разработку, повышает надежность кода, ускоряет отладку.

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

Описывается реализация обмена данными одновременно с вычислениями. Для этого алгоритмическая операция над набором элементов некоторого типа, зависящая от обмена данными, разделяется на два этапа: расчет внутренних элементов и расчет интерфейсных элементов. Расчет внутренних элементов может выполняться одновременно с обменом данными, который нужен только интерфейсным элементам. Для упрощения гетерогенной реализации предлагается использовать планировщик задач (Bogdanov Р. В., Efremov А. А. Programming infrastructure of heterogeneous computing based on OpenCL and its applications // GPU Technology Conference GTC-2013), который по заданному графу зависимостей между вычислительными и коммуникационными операциями алгоритма автоматически обеспечивает передачу данных одновременно с вычислениями на ускорителях.

В §1.4 речь идет о MIMD распараллеливании с общей памятью. Проблемы для распараллеливания операций, представимых, как правило, в виде цикла по набору элементов расчетной области, следуют из зависимости по данным между итерациями цикла. При одновременной обработке итераций, имеющих зависимости, нарушается целостность доступа к данным, возникает состояние гонки, что приводит к некорректному результату. Показаны характерные примеры таких операций: расчет потоков через грани контрольных объемов, расчет узловых градиентов, и т.д. Предлагаются подходы к устранению пересечений по данным, в частности, подход на основе многоуровневой декомпозиции для распараллеливания на сотни нитей, подходящий для ускорителей Intel Xeon Phi. В этом подходе подобласть MPI процесса разделяется далее на подобласти 2-го уровня по нитям, затем декомпозиция применяется к интерфейсным элементам 2-го уровня. Операция выполняется в 3 этапа с точками синхронизации между ними: расчет внутренних элементов 2-го уровня; расчет внутренних элементов 3-го уровня; расчет интерфейса 3-го уровня. Рассматривается проблема улучшения доступа к памяти. При работе с неструктурированными сетками большое влияние оказывает нумерация сеточных элементов, ячеек и граней. Чтобы повысить компактность данных и снизить число промахов в кэш, используется алгоритм Cuthill-McKee.

§1.5 посвящен проблеме адаптации алгоритмов к потоковой обработке для вычислений на ускорителях. Приведены способы устранения зависимости по данным для парадигмы потоковой обработки. Первый способ состоит в разделении операции на два этапа: 1) расчет по элементам (или граням) и хранение промежуточного результата в массиве по элементам; 2) сборка результатов с элементов в узлы (или ячейки). Второй способ (менее эффективный, если первый применим) состоит в разделении операции на этапы путем раскраски графа зависимостей между элементами расчетной области. Исходный набор элементов делится на поднаборы, в которых нет элементов, имеющих пересечение по данным. Также рассматривается проблема оптимизации доступа к памяти ускорителя и оптимизация вычислений.

В §1.6-1.7 представлена методика выполнения крупных расчетов, охватывающая различные аспекты выполнения численного эксперимента. Рассматриваются проблемы оптимального выбора конфигурации многоуровневой параллельной модели. Уделено внимание проблеме эффективности работы на системах коллективного доступа, с учетом простоя в системе очередей. Рассматривается проблема надежности вычислений. Для

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

Внешняя файловая система

Резервная копия текущей точки восстановления расчета 1 ВР 1/2

полный интервал осреднения

»ЩогЗЕК

установление режима течения (факт.) интервал

'установление режима течения (прогноз) осреднения (ИО)

Обозначения записей:

ВР - восстановление расчета, ВО - восстановление осреднения. ИО - интервал осреднения

Рис. 1. Схема записи точек восстановления и интервалов осреднения

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

14

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

В §2.1 описывается используемая математическая модель. Течение моделируется системой безразмерных уравнений Навье-Стокса в дивергентной форме относительно вектора консервативных переменных - плотности, трех компонент импульса и полной энергии, <? = {р,ри,ру,р\л/,ЕУ.

(1)

8Q 5F,(Q) cF,(Q) 5F4(Q) _ 5F,"(Q) 3Ff(Q) 5Ff(Q)

--1--r--1--Z.--1 ~--^ -i a

ox oy óz

(2)

dt ' дх ' dy dz

Конвективные потоки имеют вид:

F,(Q) = (pK,piT + p,puv,puw,u(E + р)) , F2(Q) = (pv,p¡<v,pv2 + p,pvw,v{E + p)j , F3 (Q) = (pw,puw, pvw, pw2 + p,w(E + p)) . где и - (и, v, w) - компоненты скорости; р,р - плотность и давление; Е = р(и2 + v2 + w2)/2 + ре - полная энергия, г - внутренняя энергия. Система уравнений (1)-(2) замыкается уравнением состояния идеального газа р = ps(y — 1), где у - показатель адиабаты.

Диссипативные потоки:

F1j®(Q) = (0,T„,TV,X„,IITX[ +vi v + yn„ + qx)T,

Ff (Q) = (О , VW«^ + vr^ + wr^ + qyf, (3)

F3ra(Q) = (0>T„,VT='MT» + vtw +yvx2Z+q2)T.

Компоненты вязкого тензора напряжений -г и вектора теплового потока q. представляются в виде следующих выражений:

И (~ди 2 ^ \i ( dv 2 \

т =— 2---diva , т =— 2-—-diva ..

" Reí дх 3 j " ReU 3 )

u (ди dv\ _ fM =T =iL

Тд= = T= = Re v & дх )' " ^ ~ Reí дТ у дТ y

Re

т = т = —

™ Re

ц I 2 т =— 2---diva

dz 3

f dv 8w — + —

dz dy дТ

Я:

Рг(у -1) дх ' Рг(у -1) ду ' ^ Рг(у -1) & Далее описывается используемый численный метод решения задачи (1). Метод повышенной точности для моделирования течения сжимаемого газа основан на конечно-объемном подходе, который имеет второй порядок аппроксимации. Для повышения порядка до 6-го включительно используются

экономичные EBR (Edge-Based Reconstruction) схемы с квазиодномерной реконструкцией переменных. Вязкие члены в уравнениях Навье-Стокса аппроксимированы, в частности, на основе конечно-элементного метода 2-го порядка точности. Пример пространственного шаблона и контрольных объемов показан на рис. 2. Для интегрирования по времени используется явная 4-х шаговая схема Рунге-Кутты 4-го порядка точности или неявная трехслойная схема второго порядка точности с поточечной линеаризацией по Ньютону пространственной разностной системы уравнений.

Рис. 2. Пример пространственного шаблона (слева) и контрольных объемов

В §2.2 описываются базовые операции над наборами элементов расчетной области (узлы, грани ячеек, сеточные элементы до 6 граней) и операции линейной алгебры (матрично-векторное произведение с разреженной матрицей мелкоблочной структуры, скалярное произведение, линейная комбинация векторов и т.д.). Все операции, составляющие алгоритм, подобраны исходя из условий совместимости со всеми уровнями многоуровневой параллельной модели, включая потоковую обработку, что обеспечивает применимость на современных гибридных архитектурах с массивно-параллельными ускорителями, включая и GPU, и многоядерные сопроцессоры архитектуры MIC (Many integrated core). Для каждой из операций, составляющих алгоритм, предложено распараллеливание в рамках параллельной модели с распределенной памятью (на десятки тысяч процессоров), существенно-многопоточное MIMD распараллеливание с общей памятью на сотни нитей, SIMD распараллеливание в рамках потоковой обработки.

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

Рис. 3. Схема работы алгоритма для моделирования сжимаемых турбулентных течений

В случае неявной схемы, получающаяся в результате линеаризации линейная система уравнений решается стабилизированным методом бисопряженных градиентов (ВЮОЗТАВ) с локальным предобуславливателем. Предложен оригинальный вариант решателя В1СС8ТАВ, адаптированный к структуре матрицы, состоящей из блоков размерности числа переменных в узлах, каждый из которых состоит из двух диагональных подблоков (для основных переменных и для переменных модели турбулентности). По сути, сразу решается две системы - для основной и турбулентной части, при этом все операции сгруппированы, что позволяет объединить обмены данными и избежать потерь на латентность сети.

В третьей главе представлен разработанный программный комплекс NOISEtte [14], реализующий новый параллельный алгоритм для моделирования сжимаемых турбулентных течений.

В §3.1 резюмируются реализованные в NOISEtte математические модели: базовые системы уравнений Навье-Стокса и Эйлера для сжимаемого газа, и их линеаризованные аналоги; модели турбулентности - RANS, LES, DES, DDES, IDDES; погруженные граничные условия; метод расчета акустики в дальнем поле FW/H, и т.д. Кратко описаны реализованные численные методы на неструктурированных гибридных сетках: EBR схемы с квазиодномерной реконструкцией переменных, полиномиальная конечно-объемная схема; римановские солверы Годунова, Роу, Хуанг; предобуславливатель Туркеля для маломаховых течений; поддержка разрывных решений - схемы WENO-EBR, MUSCL-TVD EBR; интегрирование по времени - явная схема Рунге-Кутгы до 4-го порядка, неявная схема до 2-го порядка с линеаризацией по Ньютону разностной системы уравнений; граничные условия - твердые стенки, неотражающие условия, условия характеристического типа на свободных границах, периодические условия.

В §3.2 описаны особенности программной архитектуры и реализации. Исходный код NOISEtte написан на языке С++, выбор подходов к программной реализации был обусловлен следующими требованиями: переносимость -совместимость с различными операционными системами (ОС), включая Windows, Linux; совместимость с различными типами вычислительных систем от рабочих станций до крупных суперкомпьютеров, с различными архитектурами процессоров, включая Intel, AMD, IBM и ускорители MIC; совместимость с различными компиляторами С++, включая Microsoft VS С++, Intel С++, GNU С++, IBM XL С++; автономность - базовая конфигурация может работать без подключения каких либо дополнительных библиотек, отличных от стандартных библиотек С++. Этот список дополняют такие естественные требования, как удобство работы с исходным кодом, надежность и защищенность кода от внесения ошибок, высокая производительность вычислений. Параллельная конфигурация использует MPI для модели с распределенной памятью, версия стандарта 1.2, и ОрепМР для модели с общей памятью, стандарт 3.0. На верхнем уровне программной архитектуры NOISEtte состоит из вычислительного ядра и обширной инфраструктуры, включающей в себя средства для постановки задачи, подготовки сеточных данных, обработки результатов. Вычислительное ядро реализовано в виде библиотеки статической

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

§3.3 посвящен препроцессору, основу которого составляют параллельные средства для работы с неструктурированными сетками большого объема (с числом элементов более миллиарда). Общая схема работы препроцессора показана на рис. 4. Для генерации сеток используются либо коммерческие генераторы, либо собственный сеточный генератор.

Рис. 4. Схема работы средств препроцессора

В §3.4 описывается состав и структура вычислительного ядра, которое условно разделено на функциональные модули, включая: модуль граничных условий; модуль точных решений; геометрический модуль - построение геометрии контрольных объемов; модуль формирования начальных данных и среднего поля; модуль интегрирования по времени; модуль ввода-вывода; модуль численной схемы (пространственная аппроксимация); модуль диссипативных потоков; параллельный модуль; турбулентный модуль, реализующий набор моделей турбулентности RANS, LES, DES, DDES, IDDES. Отдельное внимание уделено новому блочному варианту параллельного решателя СЛАУ на основе метода бисопряженных градиентов с локальным предобуславливателем. Также описана инфраструктура вычислительного ядра, включающая менеджер памяти, систему ручного инструментального профилирования, интерпретатор пользовательского ввода, счетчик FLOPs.

19

§3.5 посвящен постпроцессору, основу которого составляют параллельные средства для обработки результатов расчета сверхбольшого объема на неструктурированных сетках (рис. 5).

Вычислительное ядро

................i................ .................J.......... >

-j Записи моментальных полей ! | Запись восстановления счета ; | Запись контрольных поверхностей j

Л

Записи осредненных полей

Динамические данные

1.............

Суммирование интервалов

X

Сборка интервалов записи

Расчет статистики

]

Расчет дальнего поля

Средства Фурье анализа

Сборкл^запио^и^аспределенного формата |

i

¥

А

Спектральные характеристики, диаграммы направленности, АЧХ, эволюция подъемной силы, сопротивления, и т. д.

; Бинарный файл в формате Tecplot |

Рис. 5. Схема работы средств постпроцессора

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

В §3.7 представлено многоуровневое распараллеливание средствами MPI и ОрепМР. Подробно рассматриваются базовые операции, составляющие численный алгоритм, представлены различные способы реализации параллельных вычислений и исключения конфликтов по доступу к данным. Приводятся показатели производительности вычислений, измеренные встроенным в код счетчиком операций с плавающей точкой и встроенной системой таймирования. Так, чистая фактическая производительность составила 12-15% и 10-12% от пиковой производительности CPU ядра для явной и неявной схемы, соответственно. Это является достаточно хорошим результатом для таких алгоритмов на неструктурированных сетках, которым

20

1°2 Число ядер 10 Число ядер

Рис. 6. Ускорение на суперкомпьютере Ломоносов, явная схема, сетка 16 млн. узлов (слева); на МВС-10П, неявная схема, модель турбулентности, сетка 8 млн. узлов (справа)

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

Представлены показатели производительности и параллельной эффективности для MPI и ОрепМР распараллеливания. Полученное с ОрепМР ускорение на сетке порядка 1 млн на 8-ядерном процессоре Intel Xeon Е5-2690 составило 6.8 для явной схемы, и 4.7 для неявной. Подробно исследована параллельная эффективность MPI распараллеливания в суперкомпьютерных расчетах различных задач. На рис. 6 и 7 показаны примеры ускорения на суперкомпьютерах Ломоносов, МВС-10П, МВС-100К.

— МР1+ОрепМРх8

- Линейное ускорение

Эффективность 80%

—^---Неявная схема. MUSCL

- Линейное ускорение

............ Эффективность 80%

1500 -

1000 8000 Число СР^°°° 24000

Рис. 7. Тесты на исчерпание параллелизма, явная схема, слева - ускорение на сетке 1 млн узлов, МВС-100К; справа - ускорение на суперкомпьютере Ломоносов, сетка 16 млн узлов

В тестах на 60-ядерных ускорителях Intel Xeon Phi среднее внутреннее ускорение относительно 1 ядра ускорителя составило порядка 90 раз (отдельных операций до 125) и до 1.2 раз относительно 8-ядерного процессора.

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

Рис. 8. Примеры расчетов (слева направо): турбулентное течение вокруг клиновидного тела с обратным уступом, расчет интегральных характеристик модели летательного аппарата, расчет характеристик и акустических нагрузок для винта в кольце

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

В §4.1 описывается математическая модель и численный метод. Рассматривается течение несжимаемой жидкости с постоянными физическими свойствами в расчетной области, представленной параллелепипедом высотой шириной Ьу, и глубиной Ьх. Приближение Буссинеска используется для зависимости плотности от температуры. Тепловым излучением пренебрегается. В этом случае для поля скорости, и, и температуры, в, безразмерная несжимаемая система уравнений Навье-Стокса (N8) с уравнением переноса тепла имеет вид:

V-u = 0,

5,u+C(u,u) = £>u-V/?+f, д,в+С(и,9) = ?х-х D6,

(4)

(5)

(6)

где конвективный член определен как С(\\,ф) = (х\-У)ф. В случае моделирования естественной конвекции диффузионный член представляется в виде Du = Рггде Рг - число Прандтля, а вектор массовых сил f =(0,0,Рг<?), Ra - число Рэлея (по высоте Lz). В случае вынужденной конвекции с пренебрежением изменениями плотности температура становится пассивным скаляром, вектор массовых сил f = 0, а диффузионный член представляется в виде Du = Re"1 V2u, где Re - безразмерное число Рейнольдса.

Уравнения (4-6) представлены в дискретном виде на разнесенных декартовых сетках посредством сохраняющей симметрию схемы (Verstappen R.W.C.P., Veldman А. Е. P., Symmetry-Preserving Discretization of Turbulent Flow // J. Comput. Phys. 2003. Vol. 187. p. 343-368). Конечно-объемная пространственная дискретизация для уравнения (5) записывается в виде:

где М - дискретный оператор дивергенции; О - симметрическая, положительно определенная матрица диффузии, представляющая интеграл диффузионного потока через грани контрольных объемов; диагональная матрица О - матрица размеров контрольных объемов; конвективная матрица С -кососимметрическая: С(и;,) = -Сг(",,)■ Дискретизация уравнения (6) выполнена аналогично. Для дискретизации по времени конвективных и диффузионных членов используется явная схема второго порядка. Для связи давления и скорости применяется классический проекционный метод: сначала в явном виде вычисляется предиктор и£ без учета градиента давления, затем, накладывая ограничение несжимаемости, Ми£+1 = 0,,, получаем уравнение Пуассона, = Ми£, которое надо решать на каждом шаге по времени.

Оператор Лапласа представлен симметричной отрицательной полуопределенной матрицей: Ь=-МОМг.

В §4.2 представлен разработанный параллельный метод решения уравнения Пуассона, совместимый с гибридными системами с ускорителями различной архитектуры. Дискретизация уравнения Пуассона приводит к системе линейных уравнений А}°х}° = Ь3°, где А3° = -ЬеЕл''"1', а N = NxxNyxNz-

Q^l+C(uJu, +Du„-MTPh = 0А>

число узлов в сетке (индекс 3D означает связи между неизвестными в трех пространственных направлениях). Решатель для уравнения Пуассона, обозначенный KSFD (Krylov-Schur-Fourier Decomposition) [22], предназначен для расчетных областей с одним периодическим направлением. Наличие однородного направления с постоянным шагом сетки и периодическими граничными условиями позволяет применить быстрое преобразование Фурье (БПФ) для декомпозиции исходной 3D системы на набор независимых 2D систем: Â*Dx;D =b;D, ; = i,--,Nx. Системы, соответствующие низким частотам, хуже обусловлены, а система для самой низкой частоты особенно проблематична для итерационного метода, поскольку связывает всю расчетную область. Идея метода состоит в применении прямого метода на основе дополнения Шура (Soria M. et al. A direct parallel algorithm for the efficient solution of the pressure-correction equation of incompressible flow problems using loosely coupled computers // Numerical Heat Transfer, part B. 2002. Vol. 41 p. 117138) для нескольких наиболее низкочастотных систем, в то время как остальные системы решаются итерационным методом сопряженных градиентов с предобуславливателем, совместимым с потоковой обработкой.

В §4.3 описывается новое MPI+OpenMP распараллеливание решателя уравнения Пуассона [19]. Расчетная область разделена на Рх частей в периодическом направлении и на Pyz = P,,PZ частей в двух других направлениях. MPI группа разделена на подгруппы: одномерную в периодическом направлении, и двухмерную в двух других направлениях. Для выполнения БПФ в подгруппах вдоль периодической оси используется групповой обмен данными для восстановления подвекторов размера Nx вдоль этой оси. Само БПФ реплицировано внутри этих одномерных групп. MPI распараллеливание по периодическому направлению ограничено групповым обменом и заметно теряет эффективность после Рх=8. Распараллеливание по двухмерным подгруппам ограничено прямым методом на основе дополнений Шура, также имеющим групповой обмен, объем которого возрастает с числом процессоров. По опыту применения Pyz ограничено примерно 200-300 для схемы 4-го порядка, и 500-800 для схемы 2-го порядка. Таким образом, эффективный диапазон числа MPI процессов ограничен несколькими тысячами. Второй уровень распараллеливания с применением ОрепМР увеличивает число процессоров в Pt раз, где Pt - число параллельных нитей на каждый MPI процесс. Применение ОрепМР в дополнение к MPI на суперкомпьютере, состоящем из двухпроцессорных узлов с 8-ядерными процессорами,

увеличивает максимальное число процессорных ядер как минимум до -40 тысяч (для схемы 4-го порядка). Затем подробно описывается распараллеливание каждой из вычислительных операций алгоритма. Схема разбиения расчетной области показана на рис. 9.

MPI

ОрепМР

Pt

Physical space, explicit part Fourier space, Poisson solver

Рис. 9. Двухуровневая декомпозиция расчетной области

В §4.4 представлена модификация решателя, позволяющая для расчетных областей, имеющих одну или две пространственных симметрии, расширить диапазон числа процессоров еще в 2 или 4 раза, соответственно. В случае симметрии, исходная система может быть разделена далее на решение независимых подсистем.

В §4.5 представлено разработанное расширение решателя для задач без периодического направления. Расширение работает по принципу многосеточного метода и обозначено МО-К8РО [27]. В качестве приближения на втором уровне используется постановка, позволяющая применить базовый решатель К^Б: игнорируется присутствие твердых тел в потоке, шаг сетки -равномерный по одной оси, сетка с тем же числом узлов в каждом направлении, как и в исходной сетке. Операторы перехода между сетками равны единичной матрице. Показана хорошая сходимость для полностью трехмерных задач, если соблюдать некоторые ограничения по сгущению сетки. При этом сохраняется параллельная эффективность исходного К^Б метода, поскольку добавляются только операции, имеющие более высокую параллельную эффективность.

Наконец, в §4.6 представлен обший алгоритм шага интегрирования по времени и приводятся данные по тестированию параллельной эффективности его программной реализации в различных режимах и на различных

7

с. 6

о 5 ОТ

г» 4

N

В з

вычислительных системах. На рис. 10 показано ускорение для ОрепМР, число MPI процессов фиксировано и равно 64, число нитей варьируется от 1 до 8 на системе МВС-ЮОК, сетка 11 миллионов узлов.

х

Рис. Î0. Ускорение с ОрепМР, 64 MPI процесса, сетка 11 млн узлов, МВС-ЮОК Пример ускорения в режиме MPI+OpenMP на суперкомпьютере Ломоносов показан на рис. 11. Отдельно измерялось ускорение по периодической оси и по двум другим осям. В первом случае Pyz=200, Рх=8 фиксированы, Рх меняется от 1 (1600 ядер) до 8 (12800 ядер - максимально доступное для тестов), во втором случае Рх=1, Р,=8, Pyz меняется от 200 до 800.

2000 4000 6000 8000 10000 12000 ,nnn ,nnn .nnn cnn„

Number of CPU rares 2000 3000 4000 5000 6000

iNumoer oi cores Number of CPUs

Рис. 11. Ускорение на суперкомпьютере Ломоносов для MPI распараллеливания по

периодической оси (слева) и по двум другим осям (справа)

Из результатов видно, что параллелизм по периодической оси исчерпывается на Рх= 8, и на Ру=800 для сетки 330 млн узлов, при этом на сетке с 1 миллиардом узлов эффективность еще достаточно высока. Беря полученные

Poisson solver Overall algorithm Linear Speedup 80% efficiency 60% efficiency

--- I000M mesh

330Mmesh - Linear 100%

оценки 8 и 800 как оценку снизу предела параллелизма по MPI, получим оценку снизу общего числа процессоров как 8*8*800 = 51200 для узлов с двумя 4-ядерными процессорами. Для современных 8-ядерных процессоров оценка оставляет более 100 тысяч ядер, соответственно.

В пятой главе представлен разработанный программный комплекс STG-CFD&HT для крупных DNS расчетов на гибридных суперкомпьютерах.

В §5.1 приводятся реализованные математические модели и численные методы: система уравнений Навье-Стокса для несжимаемого газа, приближение Буссинеска для зависимости плотности от температуры, использованные LES модели турбулентности, модели турбулентности на основе регуляризации, сохраняющая симметрию дискретизация на структурированных сетках.

В §5.2 кратко описаны структура и особенности программной архитектуры, которые во многом совпадают с главой 3. Основное отличие состоит в применении OpenCL на третьем уровне распараллеливания для вычислений на массивно-параллельных ускорителях различной архитектуры. Также описана параллельная инфраструктура для обработки результатов расчета в распределенном формате.

В §5.3 описываются особенности программной реализации, представление данных расчетной области на структурированной разнесенной сетке, базовые структуры данных. Алгоритм составлен из 6 базовых операций, полностью адаптированных к потоковой обработке: линейный оператор (диффузия), нелинейный оператор (конвекция), линейная комбинация векторов, матрично-векторное произведение, скалярное произведение, БПФ.

В §5.4 описана OpenCL реализация каждой из базовых операций. Для ускорения позиционирования по данным используются константы препроцессора OpenCL: подставляются расстояния между соседними элементами в каждом направлении, смещение в каждом направлении с учетом гапо-элементов, смещения сетки и т. д. Данные размещаются в виде 3D или 2D полей. Коэффициенты шаблона численной схемы, а также матриц Л;", размещаются в виде отдельных 2D полей, развернутых в линейные массивы с учетом выравнивания, как показано на рис. 12. Для оператора конвекции используются специальные структуры данных для ускорения вычисления коэффициентов схемы (рис. 13).

§5.5 посвящен исследованию производительности на различных вычислительных архитектурах, включая GPU AMD и NVIDIA. В таблице 1 приводятся результаты по ускорению отдельных базовых операций на GPU

относительно последовательных вычислений на CPU. В таблице 2 представлены данные по фактической производительности для отдельных базовых операций на CPU и GPU. На рис. 14 показано сравнение результирующей производительности относительно 6-ядерного процессора для всего алгоритма в рабочем режиме.

2D поля коэффициентов для каждого из узлов шаблона схемы

.......... Тх......................"....................

аьэ

шаблон схемы

запись в виде одномерного массива

расширенная подобласть; подобласть

фиктивная зона для выравнивания

; граница выравнивания

Рис. 12. Представление 20 поля в расчетная области ускорителя в рамках распараллеливания верхнего уровня

шаблон схемы

Рис. 13. Схема вычисления коэффициентов шаблона для нелинейного оператора (снизу) Таблица 1 Ускорение на GPU относительно 1 ядра Intel Xeon Х5670 на сетке 1.3 миллиона

Операция Время на CPU, Ускорение на Ускорение на

сек. NVIDIA С2050 AMD 7970

Диффузия 0.08 22.5 45.5

Конвекция 0.25 17.6 43.0

Сложение векторов: d = aa + ßh + хс 0.01 18.0 96

Матрично-векторное произведение 0.042 19.2 44.7

FFT 0.035 7.7 12.7

Скалярное произведение 0.002 7.6 11.6

Таблица 2. Достигаемая производительность на различных операциях

Операция 1 ядро Intel Xeon Х5670 NVIDIA С2050 AMD 7970

GFLOPS % от пика GFLOPS % от пика GFLOPS % от пика

Диффузия 2.4 20 42 8.1 85 9

Конвекция 2.0 17 36 7 88 9.3

Сложение Зх векторов 1.1 9 19.3 3.7 143 15.1

М.-в.п. 0.76 6.3 14.6 2.8 34 3.6

FFT 1.5 12 11.5 2.2 19 2.0

Скалярное произв. 1.2 10 9 1.7 14 1.5

.■III

6-core CPU Intel 8-core CPU Intel GPU NVIDIA GPU NVIDIA GPU AMD XeonXS690 Xeon E5-2690 C2090 TITAN HE HREPR0 9100

Рис. 14. Относительная производительность на различных устройствах

В шестой главе представлены крупномасштабные расчеты фундаментальных задач по моделированию турбулентных течений, выполненные по разработанной технологии с помощью созданных программных комплексов. Использовались суперкомпьютеры Ломоносов, МВС-10П, МВС-100К, MareNostrum.

В §6.1 представлено моделирование дозвукового турбулентного течения вокруг тандема цилиндров с помощью программного комплекса NOISEtte, описанного в 3-й главе. Исследуются механизмы генерации аэродинамического шума от шасси самолета при заходе на посадку. Для задачи имеются подробные экспериментальные данные, полученные в NLR (Нидерланды) в рамках работ по проекту 7-й рамочной программы Евросоюза VALIANT. Конфигурация состоит их двух квадратных цилиндров, выставленных по одной линии. Расстояние между центрами цилиндров шириной D=0.04m составляет S=0.16m. Скорость потока составляет 70 м/с (М = 0.2). Число Рейнольдса (по ширине цилиндров) равно 182,000. Вдоль оси цилиндров используются периодические граничные. Ширина по периодике равна трем диаметрам цилиндров. Для пространственной дискретизации использовались сетки различных конфигураций, содержащие до 13 млн узлов, 72 млн тетраэдров. Для моделирования сжимаемого вязкого турбулентного течения использован

гибридный подход, DDES (Delayed Detached Eddy Simulation). Подробно представлены результаты расчетов для углов атаки 0° и 10° в сравнении с экспериментом и численными результатами других авторов: средние поля течения и сравнение с экспериментом (PIV измерения), профили физических величин по центральной линии между цилиндрами, распределение коэффициента давления по поверхности цилиндров, спектры пульсаций давления в контрольных точках на поверхности цилиндров и в дальнем поле. Примеры показаны на рис. 15. Для моделирования акустических возмущений в дальнем поле применялся метод Фокса-Уилльямса-Хокингса, нижняя граница диапазона разрешенной частоты 18.5 Гц.

О 500 ГЦ 1000 1500 о 500 Гц 1000 1500

Рис. 15. Сравнение с экспериментом: средние поля продольной (а) и вертикальной (б) компонент скорости, спектр в дальнем поле (в), спектр пульсаций давления на цилиндре (г)

В §6.2 представлены результаты расчета естественной конвекции на вертикальных стенках в каверне с соотношением сторон 5:1, число Рэлея Ra = 4.5xl010. Такая постановка исследуется с середины 80-х, в том числе экспериментально, и достаточно широко используется для валидации моделей. Здесь и далее использовался программный комплекс STG-CFD&HT (глава 5). Полученные результаты расчета [13] лучше подходят для валидации моделей турбулентного течения, чем экспериментальные данные, поскольку явление естественной конвекции изолировано, а вклад излучения, который существенно влияет на течение в физическом эксперименте, отсутствует в DNS расчете. Схема расчетной области и картина течения показаны на рис. 16.

В §6.3 представлен расчет отбойной струи в канале, число Рейнольдса равно 20000 (по ширине щели), сетка 103 млн. узлов. Постановка задачи показана на рис. 17 слева: на стенках, ортогональных оси Z, -периодические

условия, на нижней стенке - изотермические, на верхней - адиабатические. Результаты расчета опубликованы в [16], где они использовались для исследования и сравнения множества ЯА^ моделей турбулентности. Моментальная картина течения показана на рис. 17 справа.

Рис. 16. Схема расчетной области (слева) и моментальная картина течения в закрытой каверне с соотношением сторон 5:1 (справа), показано поле температуры

Рис. 17. Схема расчетной области для отбойной струи в канале (слева), поле модуля скорости в центральном сечении (справа)

В §6.4 представлен расчет течения вокруг бесконечного цилиндра с квадратным сечением. Число Рейнольдса равно 22000 (по ширине цилиндра), сетка 320 млн узлов. Данный расчет, на момент завершения, судя по открытой печати, являлся самым крупным в мире для данной конфигурации. Моментальная картина течения - рис. 18.

Рис. 18. Моментальная картина течения вокруг квадратного цилиндра, показано поле давления, в зоне цилиндра (справа) видно развитие нестабильности Кельвина-Гельмгольца

В §6.5 представлены расчеты, моделирующие течение в квадратной бесконечной трубе (рис. 19). Число Рейнольдса до Яе,, = 1200 (что соответствует ~21 ООО по ширине канала). Схема 4-го порядка, сетки до ] 70 млн узлов. Результаты расчетов опубликованы в [2], где также представлены расчеты течения в присутствии твердых частиц и в [28].

Рис. 19. Моментальная картина течения в квадратном канале - поле модуля завихренности в поперечном сечении (слева) и изоповерхности модуля завихренности (справа)

В §6.6 представлен расчет течения вокруг вмонтированного в стену куба [27]. Число Рейнольдса (по высоте куба) равно 7235, сетка 16 млн узлов. Это полностью 3D постановка, в которой задействовано расширение MG-KSFD (§4.6). Схема расчетной области и картина течения показаны на рис. 20.

В §6.7 приведены DNS расчёты, выполняемые автором в настоящее время. Для этих расчётов ещё не завершена обработка результатов и анализ полученных данных, поэтому в контексте данной работы расчёты приводятся только в качестве примеров, демонстрирующих применимость, производительность и эффективность программного комплекса STG-CFD&HT.

Рис. 20. Вмонтированный в стену куб в канапе, схема (слева) и моментальная картина течения - линии тока и изоповерхности модуля скорости (справа)

Первый расчёт (рис. 21) моделирует аналогичную §6.2 конфигурацию, ограниченную стенками во всех трёх пространственных направлениях. Конфигурация имеет следующие параметры: соотношение высоты к ширине 3.84, соотношение глубины к ширине 0.86, Ra=1.2xl0H. Цель расчёта дополнить экспериментальные и LES данные результатами DNS расчёта на подробной сетке из 104 млн узлов. Исследование подразумевает сопоставление результатов с новым LES расчётом, который ещё не завершен коллегами из LIMSI CNRS, Франция. Здесь расчёт представляет интерес как наиболее крупный пример эффективного применения MG-KSFD решателя (см. §4.5).

Рис. 21. Моментальная картина течения в закрытой каверне - показаны поля (а) завихренности и (б) температуры, осредненное поле температуры в верхней половине (в)

Второй расчет (рис. 22) моделирует естественную конвекцию Рэлея-Бенара. Расчетная область - параллелепипед с соотношением ширины, высоты и длины 1:1: На коротких вертикальных стенках заданы периодические граничные условия. Число Рэлея Ю10, сетка 600 млн узлов, схема 4-го порядка. На данный момент это самый крупный расчёт с 5ТС-СРВ&ПТ.

Рис. 22. Моментальная картина течения конвекции Рэлея-Бенара (поле температуры)

В заключении резюмируются результаты диссертационной работы.

ОСНОВНЫЕ РЕЗУЛЬТАТЫ РАБОТЫ

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

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

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

расчетах задач с одним периодическим направлением при использовании схем произвольного порядка точности.

• На основе разработанной технологии и алгоритмов созданы программные комплексы для крупномасштабных расчетов с использованием суперкомпьютеров с десятками тысяч процессоров, а также гибридных суперкомпьютеров:

- NOISEtte для математического моделирования турбулентного течения сжимаемого газа и создаваемых им акустических полей с повышенной точностью на неструктурированных сетках;

- STG-CFD&HT для математического моделирования турбулентного течения несжимаемой жидкости и явлений теплопереноса на основе алгоритма повышенной точности на структурированных сетках.

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

ПУБЛИКАЦИИ АВТОРА ПО ТЕМЕ ДИССЕРТАЦИИ

1. Горобец A.B. Параллельная технология численного моделирования задач газовой динамики алгоритмами повышенной точности // Журнал выч. математики и матем. физики. 2015. том. 4.

2. Нао Zhang, F. Xavier Trias, Andrey Gorobets, Dongmin Yang, Assensi Oliva, Yuanqiang Tan, Yong Sheng. Effect of collisions on the particle behavior in a turbulent square duct flow // Powder Technology. 2015. Vol. 269. pp. 320-336.

3. Богданов П.Б., Ефремов A.A., Горобец A.B., Суков С.А. Применение планировщика для эффективного обмена данными на суперкомпьютерах гибридной архитектуры с массивно-параллельными ускорителями // Вычислительные методы и программирование. 2013. т.14. стр. 122-134.

4. Горобец A.B., Суков С.А., Богданов П.Б. На пути к освоению гетерогенных супервычислений в газовой динамике // Информационные технологии и вычислительные системы. 2013. т.4. с. 15-26.

5. G. Oyarzun, R. Borrell, A. Gorobets, A. Oliva. MPI-CUDA sparse matrix-vector multiplication for the conjugate gradient method with an approximate inverse preconditioner// Computers and Fluids. 2014. Vol. 92. pages 244-252.

6. Andrey Gorobets, F. Xavier Trias, Assensi Oliva. An OpenCL-based Parallel CFD Code for Simulations on Hybrid Systems with Massively-parallel Accelerators // Procedía Engineering. 2013. Vol. 61. Pages 81-86.

7. S.A. Soukov, A.V. Gorobets, P.B. Bogdanov. OpenCL Implementation of Basic Operations for a High-order Finite-volume Polynomial Scheme on Unstructured Hybrid Meshes // Procedia Engineering. 2013. Vol. 61. Pages 76-80.

8. G. Oyarzun, R. Borrell, A. Gorobets, O. Lehmkuhl, A. Oliva. Direct Numerical Simulation of Incompressible Flows on Unstructured Meshes Using Hybrid CPU/GPU Supercomputers // Procedia Engineering. 2013. Vol. 61. Pages 87-93.

9. F. Xavier Trias, Andrey Gorobets, Hao Zhang, Assensi Oliva. New Differential Operators and Discretization Methods for Eddy-viscosity Models for LES // Procedia Engineering. 2013. Vol. 61. Pages 179-184.

10. F.X. Trias, A. Gorobets and A. Oliva. A simple approach to discretize the viscous term with spatially varying (eddy-)viscosity // Journal of Computational Physics. 2013. vol. 253. pp 405^17

11. Суков С. А., Горобец А. В., Богданов П. Б. Адаптация и оптимизация базовых операций газодинамического алгоритма на неструктурированных сетках для расчетов на массивно-параллельных ускорителях // Журнал выч. математики и матем. физики. 2010. том 53. № 8. с. 1360-1373.

12. A. Gorobets, F.X. Trias, A. Oliva. A parallel MPI+OpenMP+OpenCL algorithm for hybrid supercomputations of incompressible flows // Computers and Fluids 2013. Vol. 88. Pages 764-772.

13. F.X. Trias, A. Gorobets, C.D. Perez-Segarra, A. Oliva. DNS and regularization modeling of a turbulent differentially heated cavity of aspect ratio 5 // International Journal of Heat and Mass Transfer. 2013. vol. 57(1). pp. 171-182.

14. Абалакин И.В., Бахвалов П.А., Горобец A.B., Дубень А.П., Козубская Т.К. Параллельный программный комплекс NOISETTE для крупномасштабных расчетов задач аэродинамики и аэроакустики // Вычислительные методы и программирование. 2012. т.13. стр. 110-125.

15. F.X. Trias, A. Gorobets, C.D. Perez-Segarra, A. Oliva. Numerical simulation of turbulence at lower costs: Regularization modeling // Computers and Fluids. 2013. vol. 80. pp 251-259.

16. J. E. Jaramillo, F. X. Trias, A. Gorobets, C. D. Perez-Segarra, A. Oliva. DNS and RANS modelling of a Turbulent Plane Impinging Jet // International Journal of Heat and Mass Transfer. 2012. Volume 55. Issue 4. Pages 789-801.

17. F. X. Trias, A. Gorobets, R. W. C. P. Verstappen, A. Oliva. Symmetry-preserving regularization of wall-bounded turbulent flows // Journal of Physics: Conference Series. 2011. Volume 318. Section 4. Issue 4. p.42-60.

18. Горобец A.B., Суков С.А., Железняков A.O., Богданов П.Б., Четверушкин Б.Н. Применение GPU в рамках гибридного двухуровневого распараллеливания MPI+OpenMP на гетерогенных вычислительных

системах // Вестник ЮУрГУ. Серия: Вычислительная математика и информатика. 2011. №25 (242). стр. 76-86.

19. A. Gorobets, F. X. Trias, R. Borrell, О. Lehmkuhl, A. Oliva. Hybrid MPI+OpenMP parallelization of an FFT-based 3D Poisson solver with one periodic direction// Computers and Fluids. 2011. vol. 49. pp. 101-109.

20. F.X. Trias, R.W.C.P. Verstappen, A. Gorobets, M. Soria, A. Oliva. Parameterfree symmetry-preserving regularization modelling of a turbulent differentially heated cavity // Computers & Fluids. 2010. Vol. 39. Issue 10. P.1815-1831.

21. Горобец А. В., Суков С. А., Триас Ф. X. Проблемы использования современных суперкомпьютеров при численном моделировании в гидродинамике и аэроакустике // Ученые записки ЦАГИ. 2010. том XLI, №2, стр. 217-225.

22. A. Gorobets, F. X. Trias, М. Soria and A. Oliva. A scalable parallel Poisson solver for three-dimensional problems with one periodic direction // Computers and Fluids. 2010. Vol.39, p. 525-538.

23. F. X. Trias, A. Gorobets, M. Soria, A. Oliva. Direct numerical simulation of a differentially heated cavity of aspect ratio 4 with Ra-number up to 1011 - Part I: Numerical methods and time-averaged flow // International Journal of Heat and Mass Transfer. 2010. Vol. 53. p. 665-673

24. F. X. Trias, A. Gorobets, M. Soria, A. Oliva. Direct numerical simulation of a differentially heated cavity of aspect ratio 4 with Ra-number up to 10u - Part II: Heat transfer and flow dynamics // International Journal of Heat and Mass Transfer. 2010. Vol. 53. p. 674-683

25. Савин Г.И., Четверушкин Б.Н., Суков C.A., Горобец А.В., Козубская Т.К., Вдовикин О.И., Шабанов Б.М. Моделирование задач газовой динамики и аэроакустики с использованием ресурсов суперкомпьютера МВС-ЮОК // Доклады академии наук. 2008. том 423. №3. с. 312-315.

26. A.V. Gorobets, Т.К. Kozubskaya, S.A. Soukov. On efficiency of supercomputers in CFD simulations // Lecture Notes in Computational Science and Engineering. 2011. Vol.74. Part 9. p. 347-354.

27. A. Gorobets, F. X. Trias, M. Soria, C. D. Perez-Segarra and A. Oliva. From extruded-2D to fully-3D geometries for DNS: a Multigrid-based extension of the Poisson solver // Lecture Notes in Computational Science and Engineering. 2011. Vol. 74. Part 6. p. 219-226.

28. H. Zhang, F.X. Trias, A. Gorobets, and A. Oliva. Direct numerical simulation of a fully developed turbulent square duct flow up to Re_tau=1200 // International Journal of Heat and Fluid Flow, 54:258-267, 2015.

Подписано в печать 02.09.2015. Формат 60x84/16. Усл. печ. л. 1,1. Тираж 70 экз. Заказ А-3. ИПМ им. М. В. Келдыша РАН. 125047, Москва, Миусская пл., 4.