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

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

Математическое' моделирование реальных газодинамических течений сводится к решению многомерных уравнений газовой динамики [1, 2], что в еще большей мере усложняет теоретические исследования. Поэтому в данной области ведущая роль отводится прямому численному моделированию (DNS — Direct Numeric Simulation) газодинамических течений.

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

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

Основными методами численного решения нелинейных уравнений газовой динамики являются разностные методы [3, 4]. Большинство используемых на сегодняшний день разностных схем для решения систем уравнений газовой динамики [5] конструируются так называемым методом конечных объемов (finite volume method). Суть этого метода в том, что область интегрирования разбивается на множество ячеек и интегралы, выражающие законы сохранения в ячейках, заменяются конечными суммами. При этом предполагается, что сеточные значения в ячейках есть интегральные средние соответствующей величины по объему и потоки на границах между ячейками (дискретные потоки) определяются тем или иным способом. От интегро-интерполяционного метода [3] данный метод отличается тем, что здесь для всех величин используется одна и та же разностная сетка.

Наиболее полное описание существующих на сегодняшний день методов можно найти в работах [5, б].

Множество существующих методов вычисления дискретных потоков основываются на идеях метода С. К. Годунова, предложенного в работе [7]. Здесь описываются разностные схемы для решения нестационарных задач, и предложен способ вычисления потоков на гранях между ячейками дискретной сетки как решений задачи о распаде произвольного разрыва. Решения задачи о распаде разрыва могут находятся или из линейных приближений — при дозвуковом течении, или итерационными методами — при сверхзвуковых течениях. Многие методы вычисления дискретных потоков такие, как методы Куранта-Изаксона-Риса [9], Лакса-Фридрихса [10], Ошера [13, 13], Роу [12] и Хартена-Лакса-Ван Лира [15], по сути являются упрощенными вариантами метода С. К. Годунова. В литературе эти методы носят название методов типа Годунова.

Но метод Годунова имеет только первый порядок точности. В работе [7] доказано, что схемы порядка выше первого не являются монотонными, что приводит к появлению нефизических осцилляций вблизи разрывов (в частности на сильных ударных волнах).

Но сложность решаемых задач требует высокой точности получаемых решений. И первой [5] работой, в которой был предложен метод получения схем с порядком точности большим первого является работа [16]. Здесь и в последующих работах описываются способы перехода от разностных схем высокого порядка аппроксимации к монотонным схемам первого порядка около особенностей решения. Среди таких работ можно выделить работы Ван Лира [17] - описываются «монотонизированные» схемы, Бориса и Бука [18, 19, 20] - описывается алгоритм расчета переноса с коррекцией потоков (РСТ-метод).

В работах Хартена [21, 22] был предложен метод, получивший название ТУБ-метод или метод невозрастания полной вариации решения. Вместо монотонности этот метод обеспечивает невозрастание полной вариации и более точно передает характер поведения разрывных решений. Суть данных методов сводится к использованию разнообразных «монотонизирующих» ограничителей потоков (Нт^егз) с переключателями, зависящими от свойств решений [23, 24].

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

В работе [25] предложен новый способ реконструкции сеточных значений с автоматическим анализом гладкости решения, получивший название ENO схемы (essentially non-oscillatory scheme). В данном методе вместо фиксированного шаблона для интерполяции значений в ячейках сетки используется несколько шаблонов и выбирается тот, на котором решение является наиболее гладким. Этот метод получил дальнейшее развитие в работах [26, 27].

Далее вместо использования при реконструкции значения, полученного на одном шаблоне, в работе [28] было предложено использовать выпуклую линейную комбинацию значений полученных на всех возможных шаблонах, где весовые коэффициенты в линейной комбинации подбираются в зависимости от гладкости решения на каждом шаблоне, используя так называемые индикаторы гладкости. Данный метод получил название WENO (weighted ENO). В работах [29, 30, 31] были предложены новые способы вычисления весовых коэффициентов, а точнее новые способы анализа гладкости решения на шаблоне.

В работе [32] показано, что вблизи гладкого экстремума в ШЕ1Ч0 схемах понижается порядок точности и предложен новый масштабированный метод ^/УЕМО. Где найденные весовые коэффициенты дополнительно масшабируются с помощью специально выбранных функций, что позволяет добиться высокого порядка точности вблизи экстремума.

Бурное развитие высокопроизводительных вычислительных систем и потребности современной науки в моделировании сложных газодинамических течений способствуют развитию параллельных алгоритмов численного решения уравнений газовой динамики [33, 34]. В настоящее время технологии параллельных вычислений посвящено достаточно много работ [35, 36, 37, 38, 39], также много и электронных ресурсов освещающихактуальные вопросы данного направления [40, 41, 42].

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

• построение существенно неосциллирующих разностных схем высокого порядка точности для решения многомерной системы уравнений газовой динамики в Эйлеровых переменных;

• тестирование построенных схем на известных одномерных и двумерных модельных задачах;

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

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

Научная новизна работы заключается в следующем:

• в работе предложены новые существенно неосциллирующие схемы высокого порядка точности с масштабированием весовых коэффициентов;

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

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

• проведено прямое численное моделирование развитой турбулентной стадии неустойчивости Рихтмайера-Мешкова с использованием разработанного программного комплекса на кластере, составленном из персональных ЭВМ;

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

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

• построенные вычислительные алгоритмы апробированы на тестовых задачах и показывают работоспособность для различных типов задач;

• полученные численные результаты хорошо согласуются с известными экспериментальными и расчетными данными;

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

Автор данной работы выносит на защиту:

• построенную существенно неосциллирующую разностную схему высокого порядка точности для решения многомерной системы уравнений газовой динамики в переменных Эйлера;

• параллельный вычислительный алгоритм для построенных схем;

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

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

По главам содержание работы распределено следующим образом. В первой главе описан класс существенно неосциллирующих разностных схем высокого порядка точности для численного решения Эйлеровой системы уравнений газовой динамики. Представлены результаты численного решения с помощью построенных схем одномерной задачи о распаде разрыва с начальными данными Сода [49] и с начальными данными Лакса [50].

Вторая глава включает в себя описание параллельного вычислительного алгоритма для построенных разностных схем. Описан разработанный программный комплекс ОВ" У^Е1ЮРА11ЗВ для вычислений па многопроцессорных системах, реализующий построенные алгоритмы.

Проведен анализ эффективности распараллеливания для кластера, построенного на базе персональных ЭВМ.

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

Результаты исследования докладывались и обсуждались на:

• семинарах Средневолжского математического общества под руководством профессора Е. В. Воскресенского;

• XXXIV, XXXV, XXXVI Огаревских чтениях в Мордовском государственном университете им. Н. П. Огарева (г. Саранск, 2005, 2006, 2007 гг.);

• XIV, XV, XVI научных конференциях молодых ученых в Мордовском государственном университете им. Н. П. Огарева (г. Саранск, 2005, 2006, 2007, 2008 гг.);

• Второй Международной научной школе «Математическое моделирование, численные методы и комплексы программ» (г. Саранск, 2005 г.);

• VII Международной конференции «Дифференциальные уравнения и их приложения» (г. Саранск, 2006 г.);

• на семинаре отдела № 4 Института математического моделирования РАН (г. Москва, 2007 г.);

• Третьей Международной научной школе «Математическое моделирование, численные методы и комплексы программ» (г. Саранск, 2007 г.);

• Восьмой Международной конференции «Диффернциальные уравнения и их приложения» (г. Саранск, 2008 г.).

Основные результаты опубликованы в работах [52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63].

Заключение

.

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

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

2. Разработан параллельный программный комплекс GDWEN0PAR3D для прямого численного моделирования газодинамических течений с использованием построенных разностных схем. Построен параллельный вычислительный алгоритм на основе технологии MPI, для выполнения расчетов по построенным разностным схемам, оптимизированный под вычисления на кластере, составленном из объединенных сетью Ethernet персональных ЭВМ. Рассмотрены различные способы оценки эффективности параллельного алгоритма. Проведен анализ эффективности алгоритма при различных способах декомпозиции расчетной области. Выполнен анализ объемов межпроцессорных обменов при различных способах декомпозиции.

3. Выполнено прямое численное моделирование развития неустойчивости Рихтмайера-Мешкова с использованием программного комплекса СО? ЕГЮРА11ЗВ. Выполнено сравнение результатов расчетов с экспериментальными данными на основе значений пульсаций продольной компоненты скорости. Показано, что предлагаемый вычислительный алгоритм дает более приближенные к экспериментальным результаты по сравнению с вычислительными экспериментами других авторов.