WWW.DISSERS.RU

БЕСПЛАТНАЯ ЭЛЕКТРОННАЯ БИБЛИОТЕКА

   Добро пожаловать!


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

Давыдов Александр Александрович

Численное моделирование задач газовой динамики на гибридных вычислительных системах

05.13.18 – Математическое моделирование, численные методы и комплексы программ

АВТОРЕФЕРАТ

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

Москва – 2012

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

доктор физико-математических наук, Луцкий Александр Евгеньевич Официальные оппоненты:

доктор физико-математических наук, профессор Елизарова Татьяна Геннадьевна кандидат физико-математических наук, Семенов Илья Витальевич Ведущая организация: Институт теоретической и прикладной мехнаики им. С.А. Христиановича Сибирского отделения РАН

Защита состоится « » 2012 г. в часов на заседании Диссертационного совета Д 002.024.03 при Институте прикладной математи­ ки им. М.В. Келдыша РАН по адресу: 125047, г. Москва, Миусская пл., д. 4.

С диссертацией можно ознакомиться в библиотеке ИПМ им. М.В. Келдыша РАН.

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

Ученый секретарь Диссертационного совета, доктор физико-математических наук Н.В. Змитренко

Общая характеристика работы

Актуальность темы. Диссертационная работа посвящена актуальной создания прикладного программного обеспечения для гибридных в вычисли­ тельных систем с графическими процессорами.

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

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

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

[1]). Простые по внутренней структуре процессоры 20-30-летней давности мед­ ленно выполняли арифметические операции над вещественными числами, на их фоне вспомогательные действия по вычислению адресов этих чисел и вы­ борке их из памяти были мало заметны. Усложнение внутренней структуры процессоров с целью ускорения операций над вещественными числами было магистральным путем развития аппаратной базы высокопроизводительных вычислений.

Современные процессоры в результате многолетнего развития по этому пути выполняют «полезные» операции с плавающими числами так же быст­ ро, как «вспомогательные» операции по вычислению адресов этих чисел, и в десятки раз быстрее, чем эти числа выбираются из памяти. В итоге вспо­ могательные действия доминируют во времени исполнения программы. Рост скорости выполнения операций над вещественными числами перестает при­ водить к значительному ускорению счета. Кроме того, тактовые частоты про­ цессоров практически достигли своего предела 5-7 лет назад.

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

Такие устройства и называются вычислителями с нетрадиционной архитек­ турой.

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

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

Исследование применимости графических процессоров (GPU) для ре­ шения задач газовой динамики.

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

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

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

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

1. Исследована возможность ускорения расчета задач газовой динамики при помощи графических процессоров NVIDIA.

2. Показана эффективность параллельного расчета задач газовой динами­ ки на большом числе графических процессоров.

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

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

Апробация работы. Основные результаты, представленные в диссер­ тационной работе, докладывались и обсуждались на многих российских и международных конференциях:

1. Исследование возможностей ускорения расчета задач аэро-газодинами­ ки с помощью векторных сопроцессоров. V Всероссийская конференция молодых специалистов, СпГУ ИТМО, Санкт-Петербург, 04.2008.

2. Аcceleration of a CFD solver using commodity graphics hardware. XIV International conference of the methods of aerophysical research. Russia, Novosibirsc, July 2008.

3. Применение графических процессоров для расчета задач аэродинами­ ки. XVII Всероссийская конференция «Теоретические основы и констру­ ирование численных алгоритмов для решения задач математической физики с приложением к многопроцессорным системам», посвященная памяти К.И. Бабенко. Дюрсо, 2008 год.

4. Макет гибридного суперкомпьютера «МВС-экспресс». XVII Всероссий­ ская конференция «Теоретические основы и конструирование числен­ ных алгоритмов для решения задач математической физики с приложе­ нием к многопроцессорным системам», посвященная памяти К.И. Ба­ бенко. Дюрсо, 2008 год. Соавторы: С. С. Андреев, С. А. Дбар, А. О.

Лацис, Е. А. Плоткина 5. Численное моделирование задач аэро-газодинамики на гибридном су­ перкомпьютере «МВС-Экспресс». Международная научная конферен­ ция «Параллельные вычислительные технологии 2009», Нижний Нов­ город, 30 марта - 1 апреля 2009 года. Сборник тезисов.

6. Программный комплекс для расчета задач газовой динамики на гибрид­ ном суперкомпьютере «МВС-Экспресс». XVIII Всероссийская конфе­ ренция «Теоретические основы и конструирование численных алгорит­ мов для решения задач математической физики с приложением к много­ процессорным системам», посвященная памяти К. И. Бабенко. Дюрсо, 2010 год.

7. Numerical Simulation of the Continuous Media Mechanics Problems on the Hybrid Supercomputer «MVS-Express». XV International conference of the methods of aerophysical research. Russia, Novosibirsc, November 2010.

8. О моделях и технологиях программирования суперкомпьютеров с нетра­ диционной архитектурой. Научный сервис в сети Интернет: суперком­ пьютерные центры и задачи: Труды Международной суперкомпьютер­ ной конференции (20-25 сентября 2010 г., г. Новороссийск). - М.: Изд-во МГУ, 2010. - 694 с. Соавторы: С. С. Андреев, С. А. Дбар, А. О. Лацис, Е. А. Плоткина Структура работы. Диссертация состоит из введения, 5 глав, заклю­ чения и списка цитируемой литературы.

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

В первой главе приводится описание физического устройства графиче­ ского процессора NVIDIA и аппаратно программных средств CUDA (Compute Unified Device Architecture). Представлена реализация численного метода С.

К. Годунова для системы трехмерных уравнений Эйлера для архитектуры CUDA.

Несмотря на то, что в целом архитектура CUDA хорошо подходит для реализации явных сеточных методов, применимость ее для таких методов как метод С. К. Годунова и алгоритмов на его основе не вполне очевидна.

В действительности, в методе С. К. Годунова на каждой грани расчетной сетки решается задача о распаде произвольного разрыва. В классическом подходе эта задача решается итерационным способом, и, скорость сходимо­ сти итерационного процесса меняется от грани к грани. Архитектура CUDA устроена таким образом, что параллельные потоки разделяются на так на­ зываемые варпы (warps). Внутри варпа все потоки выполняют идентичный набор инструкций. Это означает, что если итерационный процесс сходится с разной скоростью в разных потоках одного варпа, то время выполнения вар­ па будет определяться временем выполнения самого медленно сходящегося потока. Тем не менее, на практике оказывается, что даже в случае расчетов течений с большими градиентами физических величин, ускорение расчетов на GPU относительно CPU падает незначительно. Это связано, по видимому, тем, что увеличение числа итераций не требует дополнительного чтения из медленной глобальной памяти, а все вычисления производятся с данными из сверхбыстрой регистровой памяти.

На основании проведенных в главе 1 исследований можно сделать вывод, что аппаратно–программные средства CUDA могут быть успешно использо­ ваны для решения задач газовой динамики. Так, при расчете трехмерных уравнений Эйлера методом С. К. Годунова был получены ускорения расчета вплоть до 70 раз по сравнению с одним ядром универсального процессора (Рис. 1).

Вторая глава посвящена описанию программного комплекса «Express 3D». Программный комплекс «Express 3D» представляет собой структуру па­ мяти для хранения многоблочных индексных сеток, набор методов для об­ мена данными между узлами вычислительной системы, набор методов для обмена данными между графическими ускорителями и универсальными про­ цессорами. В комплекс «Express 3D» также входят утилиты для подготовки начальных данных, оптимального разбиения блоков между процессорами.

69.5x 55.3x 49.3x 35.2x 12.9x 253 503 1003 1503 20Размер расчетной сетки Рис. 1. Зависимость ускорения расчета на GPU по сравнению с одним ядром CPU зависи­ мости от размерности сетки Модульная структура комплекса позволяет в короткие сроки добавлять к нему новые математические модели и алгоритмы. Реализованные в пер­ вой главе алгоритмы для решения уравнений Эйлера, а так же алгоритмы решения квазигазодинамических и модифицированных квазигазодинамиче­ ских уравнений, о которых речь пойдет в последующих главах естественным образом были включены в состав комплекса в виде отдельных модулей.

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

Время обработки блока на универсальном процессоре практически ли­ нейно зависит от его размерности (общего количества элементов). Это поз­ воляет балансировать загрузку процессоров по объему обрабатываемых дан­ ных. Фактически это можно делать до начала расчета. Напротив, ускорение расчета на графическом ускорителе по сравнению с универсальным процессо­ ром существенно зависит от объема обрабатываемых ядром данных и о неко­ торых других параметров, таких как размерности блока по направлениям. В таблице 1 приведено ускорение расчета блоков характерного размера при ис­ Ускорение Размерность блока сетки 4050 11070 12690 19215 37515 43005 52521 101199 117547 1975Время GPU, с. 11 18 18 17 23 23 23 44 46 Время CPU, с. 29 79 96 142 283 323 407 810 978 17Ускорение 2,6 4,4 5,3 8,4 12,3 14,0 17,7 18,4 21,3 24,Таблица 1 : Ускорение расчета блока на GPU по сравнению с одним ядром CPU в зависимости от размерности блока сетки пользовании графического процессора. Таким образом, время обработки бло­ ка на GPU зависит от размерности блока существенно нелинейно. Более того, заранее трудно спрогнозировать с какой скоростью будет обрабатываться тот или иной блок. Все это приводит к необходимости балансировать загрузку не по объему данных, а по временам обработки каждого блока.

В комплексе «Экспресс-3D» применяется следующий алгоритм распреде­ ления блоков про процессорам (или по графическим процессорам) обеспечи­ вающий равномерность загрузки: сначала производится разбиение блоков по объему данных на некоторое количество процессоров (ускорителей). Это мо­ жет быть и один процессор (ускоритель), главное, чтобы объема оперативной памяти хватило для размещения данных задачи. Далее проводится расчет некоторого числа шагов с целью определения среднего времени обработки каждого блока. После чего производится разбиение блоков на процессоры (графические процессоры) уже про временам обработки блоков с использова­ нием алгоритма оптимального распределения камней по ящикам предложен­ ного А. Шараховым [11], который является переработанным и дополненным алгоритмом решения "Задачи о куче камней"[9].

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

Комплекс «Экспресс-3D» позволяет использовать для расчета не только графические процессоры, но и ядра универсального процессора. Так, если узел вычислительной системы состоит из 8 ядер универсального процессора и 2 графических ускорителей (например, СК «Ломоносов», МГУ), то рабо­ ту следует разбивать на 8 частей - на 6 универсальных ядер и на 2 пары CPU+GPU. Конечно, сложность балансировки загрузки при этом существен­ но возрастает.

В работе предложен эффективный алгоритм для распределения блоков между процессорами и ускорителями. Основная идея алгоритма состоит в том, чтобы мелкие блоки, которые не очень хорошо ускоряются на GPU, но число которых, как правило, велико, считать на CPU. Применение такого подхода позволяет на 30-40% повысить эффективность узла вычислительной системы (по сравнению с вычислениями только на ускорителях). В таблице 2 приведены времена расчета и ускорение расчета на одном и двух узлах установки К-100, работающей в ИПМ им. М. В. Келдыша РАН [12], при ис­ пользовании только GPU и при расчете на GPU и CPU одновременно.

1 узел К-100 2 узла К-13GPU 3GPU+8CPU 6GPU 6GPU+16CPU Время 65.98 47.72 38.8 28.ускорение 1.38 1.Таблица 2 : Сравнение времени расчета на GPU и на CPU+GPU.

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

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

В комплексе «Express 3D» в качестве коммуникационной среды и моде­ ли параллельного программирования была использована технология Shmem [13].

Модель программирования shmem [15] подразумевает взаимодействие независимых процессов, каждый – со своей локальной памятью, занумерован­ ных по порядку от нуля, и в этом она похожа на модель программирования MPI [6]. Отличие от MPI состоит в том, что обмены данными между процес­ сами являются односторонними. Чтобы данные были переданы из процесса А в процесс Б, требуется не согласованная активность обоих процессов, как в MPI, а лишь «желание» одного из участников. Например, процесс А может «насильно» послать данные процессу Б, без какой-либо ответной активности с его стороны. Процесс Б, в свою очередь, может «насильно» прочитать дан­ ные из процесса А. В англоязычных описаниях этой модели ее принято назы­ вать «модель put/get», в отличие от принятой в MPI «модели send/receive».

Выбор модели программирования Shmem был обусловлен следующими факторами. Во-первых, первый гибридный кластер, появившийся у нас в стране - «МВС-Экспресс» [3], долгое время не имел другой модели парал­ лельного программирования. Во-вторых, модель односторонних обменов, по мнению автора, существенно упрощает понимание программы и сокращает время разработки. В-третьих, существуют программные реализации shmem основанные на MPI [8, 16], которые хотя и несколько снижают производитель­ ность относительно реализаций shmem привязанных к оборудованию, позво­ ляют упростить программирование, сохранив при этом переносимость кода на многие вычислительные системы.

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

Это также вносит вклад в снижение параллельной эффективности прило­ жений. Системы программирования которые позволили бы избежать такого двухуровневого обмена данными в настоящий момент только разрабатыва­ ются.

В таблице 3 приведены времена расчетов уравнений Эйлера (обтекание крылатого тела, см. Главу 4) для различного числа универсальных и гра­ фических процессоров. Количество блоков сетки – 132, общее число узлов 8 · 106.

1x CPU 11xCPU 3x 6x 9x 12x 24x 30x core cores GPU GPU GPU GPU GPU GPU Время, с 13510 1424 360 211 136 108 60 Ускорение 1,0 9,5 37,5 64,0 99,3 125,6 226,7 291,Эффективность CPU 100% 86% Эффективность GPU 100% 85% 88% 84% 77% 78% Таблица 3 : Время расчета уравнений Эйлера на различном числе CPU и GPU.

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

В таблице 4 приведены времена расчета системы квазигазодинамических уавнений (задача о каверне с движущейся крышкой, см главу 3). Количество блоков сетки – 216, размерность каждого блока 50 50 50.

Число GPU 4 6 9 18 27 36 54 72 1Время, с 11090 7580 5080 2670 1800 1360 940 730 5Эффективность, % 100 97,5 97,0 92,3 91,3 90,6 87,4 84,4 77,Таблица 4 : Время расчета квазигазодинамических уравнений на различном числе GPU.

В Третьей главе описано построение численного метода для решения квазигазодинамической системы (1)-(2). Описание и физический смысл кото­ рой приведены в [5] + jmi = 0, (1) t xi ui + jmjui + p = ji, (2) t xj xi xj E + p E + jmj + qj = jiui. (3) t xj xj xj Здесь и далее по повторяющимся индексам подразумевается суммирова­ ние, E = u2 2 + – полная энергия. Поток массы задается как jmi = (ui - wi), wi = ujui + p, (4) xj xi Также рассматривается модифицированная квазигазодинамическая си­ стема (5)-(7):

2 + + jmi = 0, (5) 2 t2 t xi 2ui + ui + jmjui + p = ji, (6) 2 t2 t xj xi xj 2E E + p + E + jmj + qj = jiui. (7) 2 t2 t xj xj xj Тензор вязких напряжений П, тепловой поток q и соотношения, замы­ кающие систему, определяются следующими соотношениями:

1 ij = NS + ui uk uj + p + ij uk p + p uk, (8) ij xk xj xk xk 2 NS = +µ uj + ui - ij uk, (9) ij xi xj 3 xk 1 NS NS qi = qi - uiuj + p, qi = - T, (10) xj xj xi p p = ( - 1), T =. (11) R Динамическая вязкость µ, коэффициент теплопроводности и релакса­ ционный параметр , имеющий размерность времени, имеют следующий вид T R µ µ = µ, = µ, =, (12) T ( - 1) Pr pSc КГД система (1) – (3) получена из балансного уравнения в предполо­ жении малости времени по сравнению с характерным газодинамическим временем tгаз. Учитывая следующий член разложения функции распределе­ ния в ряд Тейлора, придем к варианту системы (5) – (7).

Эта система имеет гиперболический тип [7]. Важным свойством КГД системы уравнений (5) – (7) является ее близость к классическим уравнениям Навье – Стокса (НСУ). Связь между этими двумя системами может быть представлена в виде КГД = НСУ + O Kn2, (13) где Kn – число Кнудсена. Следует отметить, что сами уравнения Навье – Стокса получены из уравнения Больцмана с точностью до членов второго порядка малости по числу Кнудсена. Дополнительные диссипативные члены в КГД системе, содержащие множитель , обращаются в нуль в тех обла­ стях течения, где решение удовлетворяет стационарным уравнениям Эйлера.

Заметим также, что наличие вторых производных по времени в уравнениях (5) – (7) окажется важным фактором при формировании вычислительного алгоритма.

Конечно разностная аппроксимация систем (1) – (3) и (5) – (7), строит­ ся на основе метода контрольных объемов. Все газодинамические параметры относятся к центрам ячеек разностной сетки, а в качестве контрольного объ­ ема берется сама ячейка. Интегрируя уравнения (1) – (3) по объему ячейки, мы получаем законы сохранения в интегральной форме. При этом изменение газодинамических параметров в ячейке определяется суммой потоков консер­ вативных переменных (плотности , импульса u и полной энергии E) через все ее грани. Для аппроксимации пространственных производных, входящих в выражения для потоков, используются центральные разности, а значения газодинамических переменных в центрах граней вычисляются с помощью ли­ нейной интерполяции. Для обеспечения устойчивости счета системы (1) – (3) по явной схеме к релаксационному параметру в (12) добавляется слагаемое, пропорциональное шагу пространственной сетки.

Для построения численного алгоритма для гиперболического варианта КГД–системы уравнений (5) – (7) перепишем ее в компактной форме:

2 f U + U = FQGD. (14) t2 t Здесь U = (, u1, u2, u3, E)T – вектор консервативных переменных, FQGD – потоки, соответствующие системе (1) – (3), f – новый релаксаци­ онный параметр. Физический смысл дополнительных слагаемых со второй производной по времени, а с ними и параметра f становится понятным, если переписать (14) в форме релаксации потоков:

U = F, f F = FQGD - F. (15) t t В соответствии с (15) потоки консервативных переменных не могут в от­ личие от системы (1) – (3) мгновенно достичь значений FQGD, они стремятся к ним, стартуя со значений в предыдущий момент времени, и f – характер­ ное время релаксации потоков. При малых значениях f можно пренебречь левой частью второго из уравнений (15), и мы получим систему (1) – (3), поскольку в этом случае F = FQGD.

Обратим внимание, что идея релаксации потоков не нова. Еще в се­ редине прошлого века физики обратили внимание на нефизичное свойство уравнения теплопроводности, основанного на законе Фурье: бесконечная ско­ рость распространения возмущений, свойственная параболическим уравнени­ ям. Этот факт приводит к явному расхождению с экспериментальными дан­ ными при исследовании быстро протекающих процессов с большими градиен­ тами температур, таких, как течения разреженных газов, низкотемператур­ ный теплоперенос в твердых телах, электронная теплопроводность в плазме и др. Чтобы избежать этого и обеспечить конечную скорость распространения возмущений, было предложено гиперболическое уравнение теплопроводности [16].

При построении разностной схемы для гиперболической системы (5) – (7) используем ту же пространственную аппроксимацию потоков FQGD, как и для системы (1) – (3). Второе уравнение в (15) является линейным ОДУ первого порядка, и может быть проинтегрировано аналитически на интервале времени tj, tj+1 :

Fj+1 = FjD + FQGD(1 - D), D = exp (-t/f), t = tj+1 - tj. (16) Относительно выбора параметра f, можно сказать, что с одной стороны, он должен быть достаточно большим, чтобы обеспечить наилучшую устойчи­ вость схемы, а с другой – достаточно малым, чтобы получаемое решение не слишком сильно отличалось от решения исходной системы.

Переход на новый временной слой осуществляется в следующем поряд­ ке. Сначала по значениям газодинамических величин в момент времени tj вычисляются потоки FQGD для уравнений (1) и (3), описывающих изменение скалярных переменных. Затем в соответствии с (16) находим новые значения потоков этих переменных Fj+1. Используя простейшую аппроксимацию по времени первого уравнения (15), Uj+1 = Uj + t · Fj+1, вычисляем новые значения плотности и энергии. После этого та же процедура производится с уравнением импульса (2), причем потоки FQGD для этого уравнения вычис­ ляются по вновь полученным значениям j+1 и j+1. Переход на слой tj+завершается вычислением новых значений скоростей uj+1.

µ h = + , (17) pSc c где c – локальная скорость звука, h = x2 + x2 + x2, x1, x2, x3 – 1 2 шаги прямоугольной пространственной сетки, – числовой параметр поряд­ ка единицы, подбираемый экспериментально.

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

Для различных значений релаксационного параметра f были числен­ но определены максимальные значения , обеспечивающие вычислительную устойчивость. Результаты вычислений для случая Ma = 0.05 изображены на рис. 2 Точки на оси ординат соответствуют f = 0, т.е. решению исходной системы (3) – (5). Видно, что значения max оказываются значительно боль­ ше, чем при решении исходной системы и практически не зависят от шага пространственной сетки. Таким образом, мы имеем практически Курантов­ ское условие устойчивости для существенно дозвукового течения вязкого га­ за. При этом оптимальное значение f, при котором достигается max, почти пропорционально шагу h.

Рис. 2. Зависимость максимально допустимого числа Куранта max от релаксационного параметра f для различных пространственных сеток.

При решении той же задачи для случая Ma = 0.1 были получены анало­ гичные результаты. Максимально допустимые значения max близки к еди­ нице, а оптимальные значения f, при которых эти значения достигаются удвоились по сравнению с вариантом Ma = 0.05. В Таблице 5 представлено сравнение оптимальных значений параметра f и соответствующих им чисел Куранта max для различных вариантов расчета.

Таблица 5 : Оптимальные значения параметра f и соответствующие значениия max Анализируя эти результаты, можно заметить, что для всех вариантов расчета оптимальное значение f близко к величине Ma · h, т.е. релаксаци­ онный параметр должен быть одного порядка с шагом по времени. С фи­ зической точки зрения этот результат представляется вполне естественным.

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

Таким образом, достигнутое увеличение шага по времени приводит к реаль­ ной экономии времени расчета задачи.

Отметим, что установившиеся при max распределения газодинамиче­ ских величин и структура течения, полученные по модифицированной си­ стеме, отличаются от полученных по системе (1) – (3) на доли процента.

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

В Четвертой главе с помощью разработанного программного комплек­ са проводиться расчет аэродинамических характеристик летательного аппа­ рата (крылатого тела) рис. 3 обтекаемого сверхзвуковым потоком газа под различными углами атаки при различных числах Маха набегающего потока.

Расчетная сетка состояла из 132 блоков. Общее число ячеек сетки 1.8·106.

Для каждого из углов атаки = 5, 10, 15 моделировались течения с различными числами Маха M = 1.5, 2.0, 2.5 набегающего потока. Получены зависимости сопротивления подъемной силы в зависимости от углов атаки и чисел Маха набегающего потока.

Общая картина течения для M = 2.0, = 10 изображена на рисунке На рисунке 5 представлены зависимости подъемной силы от угла атаки для различных чисел Маха. Для всех числе Маха зависимость подъемной силы от угла атаки практически линейно, что хорошо согласуется с теорией [14].

Давление торможения определяется параметрами набегающего потока, Рис. 3. Модель крылатого тела и расположение блоков индексной сетки.

Рис. 4. Общая картина течения, уровни давления.

и не должно зависеть от угла атаки, при этом положение точуи торможения, конечно, меняется. В таблице 6 приведено давление торможение для числа Маха набегающего потока M = 2.5 и погрешность относительно аналитиче­ ского решения P = 8.52613589 [14].

Рис. 5. Зависимость подъемной силы от угла атаки.

Погрешность Угол Давление относительно атаки торможения точного решения 5 8.6378 1.31 % 10 8.5646 0.45 % 15 8.4986 0.32 % Таблица 6 : Давление торможения для M = 2.В Пятой главе в рамках модели модифицированных квазигазодинами­ ческих уравнений проводится моделирование истечения тяжелой жидкости из резервуара.

Модель резервуара представляет из себя куб с равномерной простран­ ственной сеткой. В начальный момент жидкость покоится рис. (6). В одной из боковых стенок имеется отверстие прямоугольной формы совпадающее с ячейками сетки. Для отверстия задаются условия внешнего давления, плот­ ность и компоненты скорости “сносятся” из ячейки резервуара. Считая что объем вытекающей жидкости много меньше объема резервуара положим на свободной верхней границе входящий поток, равный по массе вытекающему Рис. 6. Модель резервуара.

Рис. 7. Линии тока установившегося течения.

через отверстие.

На рисунке 7 показаны линии тока установившегося течения.

В работе были исследованы течения с различными положениями отвер­ стия при различных формах резервуаров.

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

ОСНОВНЫЕ РЕЗУЛЬТАТЫ 1. Разработан и реализован комплекс программ «Express-3D» для реше­ ния задач газовой динамики на многоблочных индексных сетках. Ком­ плекс позволяет решать задачи, как на универсальных кластерах, так и на гибридных вычислительных системах. Модульная структура ком­ плекса позволяет в короткие сроки добавлять к нему новые математи­ ческие модели и алгоритмы. Показана высокая параллельная эффек­ тивность комплекса при расчетах на большом числе графических про­ цессоров.

2. В рамках разработанного комплекса реализованы алгоритмы для реше­ ния уравнений Эйлера, алгоритмы для расчета квазигазодинамических и модифицированных (гиперболизованных) квазигазодинамических урав­ нений для архитектуры CUDA.

3. При помощи программного комплекса «Express-3D» Проведено экспериментальное исследование устойчивости численного алгоритма для модифицированных (гиперболизованных) квазигазоди­ намических урвнений. Определены зависимости максимально допусти­ мых (при которых сохраняется устойчивость) чисел Куранта от значе­ ний релаксационного параметра.

Проведено численное исследование аэродинамических характеристик летательного аппарата (крылатого тела) при различных углах атаки и числах Маха набегающего потока. Получены зависимости подъемной силы от угла атаки и чисел Маха.

Построена модель и исследованы параметры истечения тяжелой жид­ кости из резервуара через отверстие. Определены особенности течения в зависимости от положения отверстия и формы резервуара.

ОСНОВНЫЕ ПУБЛИКАЦИИ ПО ТЕМЕ РАБОТЫ Статьи в рецензируемых журналах, рекомендованных ВАК [1] Давыдов А. А. Применение графических сопроцессоров на суперкомпью­ тере «МВС Экспресс» для расчета задач аэро-газодинамики// Научно­ технический вестник СПГУ ИТМО, № 54, Технологии высокопроизво­ дительных вычислений и компьютерного моделирования, сс. 178-180, 2008 г.

[2] Давыдов А. А. Численное моделирование задач аэро-газодинамики на ги­ бридном суперкомпьютере «МВС-Экспресс»// Журнал Математическое моделирование, том 22, № 4, 2010 г., с. 90-98.

[3] Давыдов А. А., Лацис А. О., Луцкий А. Е., Смольянов Ю. П., Чет­ верушкин Б. Н., Шильников Е. В. Многопроцессорная вычислительная система гибридной архитектуры «МВС-Экспресс»// ДОКЛАДЫ АКА­ ДЕМИИ НАУК, 2010, том 434, № 4, с. 459-4[4] Давыдов А. А., Четверушкин Б. Н., Шильников Е. В. Моделирование течений несжимаемой жидкости и слабосжимаемого газа на многоядер­ ных гибридных вычислительных системах// Журнал вычислительной математики и математической физики, 2010, том 50, № 12, с. 2275-22Список литературы [5] Елизарова Т. Г. Квазигазодинамические уравнения и методы расчета вязких течений. – М: Научный мир, 2007. – 352 с.

[6] Воеводин В. В., Воеводин Вл. В. Параллельные вычисления. - СПб.:

БХВ-Петербург, 2002. - 600 с.

[7] Злотник А.А., Четверушкин Б.Н. Параболичность квазигазодинами­ ческой системы уравнений, гиперболичность одной ее модификации и устойчивость малых возмущений для них // Ж. вычисл. мат. и ма­ тем. физ., 2008, т. 49, №3, сс. 445 - 472.

[8] Корж А. А. Масштабирование Data-Intensive приложений с помощью библиотеки DISLIB на суперкомпьютерах Blue Gene/P и “Ломоносов” // Труды конференции “Научный сервис в сети Интернет-2011” [9] Романовский И. В. Алгоритмы решения экстремальных задач, “Нау­ ка”, М., 1977, с.248.

[10] Четверушкин Б. Н. Кинетические схемы и квазигазодинамическая си­ стема уравнений. М.: Макс Пресс, 2004.

[11] http://delphiworld.narod.ru/base/stones_in_boxes.html [12] http://www.kiam.ru/MVS/resourses/k100.html [13] Климов Ю. А., Лацис А. О. Руководство по использованию сети МВС­ Экспресс на К-100.

http://www.kiam.ru/MVS/documents/k100/mvseuserguide.html [14] Ландау Л. Д., Лифшиц Е. М. Теоретическая физика: Учебное пособие.

В 10 т. Т. VI. Гидродинамика. – 3-е изд. перераб. – М: Наука. Гл. ред физ-мат. лит., 1986. – 736 с.

[15] http://www.shmem.org/ [16] http://www.kiam.ru/MVS/research/ [17] http://www.nvidia.com Подписано к печати 23.04.2012. Формат 60х84/16. Усл. печ. л. 1,5. Тираж 100 экз. Заказ.

ИПМ им.М.В.Келдыша. 127047, Москва, Миусская пл.,







© 2011 www.dissers.ru - «Бесплатная электронная библиотека»

Материалы этого сайта размещены для ознакомления, все права принадлежат их авторам.
Если Вы не согласны с тем, что Ваш материал размещён на этом сайте, пожалуйста, напишите нам, мы в течении 1-2 рабочих дней удалим его.