Проектирование и эффективная реализация средств численного моделирования в задачах вихревой гидроаэродинамики С.В. Гувернюк, Д.А. Григоренко, А.И. Гирча, П.Р. Андронов, Г.Я. Дынникова НИИ механики МГУ

Плоские нестационарные течения вязкой (или идеальной) несжимаемой жидкости Re = 2000 Связанная задача: совместное решение уравнений гидродинамики и динамики при движении тела в жидкости

Используемые математические модели В качестве базовых моделей среды используется две: модель идеальной несжимаемой жидкости (двумерные уравнения Эйлера) модель вязкой несжимаемой жидкости (двумерные уравнения Навье-Стокса) Во всех случаях начальное состояние среды – покой во всем пространстве Движущиеся тела – двумерные, произвольной формы конечной толщины Граничные условия: условие прилипания (в вязкой жидкости) условие непротекания с априорным заданием точек отрыва (в идеальной жидкости)

«Вихревая» формулировка плоской сопряжённой задачи динамики и аэрогидродинамики R В начальный момент среда и тело покоятся в неограниченном пространстве. Тело начинает двигаться под действием внешних сил, моментов и наложенных связей, испытывая также дополнительное влияние аэродинамических сил и моментов со стороны среды, возмущённой движением тела. - гидродинамический импульс среды (Г. Ламб, 1947) - вращательный импульс среды (Г. Ламб, 1947) - P.G.Saffman. Vortex Dynamics (1992) Для безвихревого движения теория развита Томсоном, Тэтом, Кирхгофом (1869).

Лагранжевы методы Метод дискретных вихрей (МДВ) С.М. Белоцерковский и др. Модификации Кратные цепочки С.Я. Герценштейн, С.В. Гувернюк, Х. Исванд Метод вязких вихревых доменов (ВВД) П.Р. Андронов, С.В. Гувернюк, Г.Я. Дынникова Уравнения Навье-Стокса Генерация вихрей Диффузионная скорость Y. Ogami, T. Akamatsu

Схема алгоритма решения задач нестационарной гидродинамики 1. Вычисление циркуляции, обеспечивающей условие непротекания, нахождение скорости тел 2. Вычисление скорости движения жидкости в точках нахождения вихрей 3. Вычисление скорости смещения вихря относительно жидкости, связанной с диффузией (диффузионной скорости Vd) 4. Перемещение вихрей 5. Аннигиляция близко расположенных вихрей противоположного знака

Требования Необходимо обобщение Различные типы тел Различные типы взаимодействия Эффективность Переносимость Расширяемость Отчуждаемость

Структура комплекса Средства подготовки заданий Скрипты Средства расчета Средства обработки результатов Скрипты Дополнительные средства (GUI) Сторонние программы

Обобщенный алгоритм Единые обозначения величин Подробное описание, вспомогательные формулы Замкнутые двумерные тела произвольной формы Случаи идеальной и вязкой среды Число степеней свободы 0, 1, 2, 3 Возможны внешние силы Различные сценарии схода вихрей Несколько не сталкивающихся тел Множество рассчитываемых параметров (Cx,Cy,F,V,P,t…)

Архитектура Счетная программа: ООП С++ Графический интерфейс: Python Qt Концепция: Open Source GPL

Структура программы Вспомогательные классы Application Body Medium Scope Container

Эффективность Факторы: Прозрачность кода Переносимость Основная проблема: расчет поля скоростей Пути повышения: Оптимизация Профилирование Использование возможностей компилятора Распараллеливание Изменение алгоритма

Алгоритмы Точные Приближенные (до ) Отличия от классической задачи n тел: Аннигиляция вихрей Нет симметрии

Использование Т-системы Автоматическая балансировка загрузки Гетерогенность (сериализация) Новые возможности (в сравнении с MPI) Недостатки: Лишние пересылки данных Тонкие ограничения на C++ Модификатор static

Время работы MPI 1248 Без объединения вихрей 30 шагов % % % % С объединением вихрей 30 шагов % % % % С объединением вихрей 60 шагов % % % %

T-система vs MPI Вычислительные машины: AMD Athlon(tm) 64 Processor Intel(R) Celeron(R) CPU 2.53GHz Результаты T-система (Celeron+Athlon): 70 секунд MPI (Celeron+Athlon): 105 секунд Athlon: 73 секунды Celeron: 171 секунда Идеально: 51 секунда

Переносимость Windows (MSVC, Cygwin + g++) Linux (g++) AIX (xlC) MPI T-система Трудности при переносе возможны

Достоверность результатов Теоретическое обоснование Проверка инвариантов Сравнение с расчетами Библиотека GMP Другие методы Сравнение с экспериментом

Нестационарное обтекание цилиндра квадратного сечения при Re = 57. Мгновенные линии тока. Слева – расчёт В.П. Шкадовой и др. (2004) сеточным методом при Re=56,6; справа – расчёт методом ВВД. Интервалы безразмерного времени между последовательными картинками слева и справа одинаковые.

Р езультаты тестовых расчётов на основе метода ВВД при тех же начальных данных для случая Re=1000 при разной степени дискретизации по времени и по пространству: чёрная линия – при разбиении контура пластины на 164 отрезка и при шаге по времени 0,025, красная линия – при разбиении контура пластины на 328 отрезков и при шаге по времени 0,0125. Максимальное наблюдавшееся отличие угловой скорости при разной степени дискретизации составляет около 5%.

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

Расчётная (слева) и экспериментальная (справа) картины течения при авторотации вертушки с двумя лопастями для различных фаз поворота (поток слева направо, вращение против часовой стрелки) Sh=0,064 Sh=0,089 Re = НИИ механики МГУ, аэродинамическая труба А-8; гг.

Сравнение с экспериментом (3 лопасти) Кратные цепочки вихрей Автоколебания Скорость потока U=30 м/с, длина лопасти h= 1,4 R

Сравнение с экспериментом (3 лопасти) Кратные цепочки вихрей Авторотация Скорость потока U=30 м/с, длина лопасти h= 1,4 R

Практическое применение Расчет мощности ветроэнергетических установок Ветровые нагрузки на высокие строения Нестационарные аэродинамические нагрузки на крыло Нестационарная теплопередача при обтекании нагретых цилиндров

Расчёт методом ВВД эволюции поля завихренности при авторотации пластины для случая: Re = 1000; безразмерный момент инерции J = 20; начальный угол поворота α = 10 градусов. Шкала от синего до красного цвета соответствует изменению завихренности от –0,5 до +0,5. а)а)б) в)г) д)е)

Перспективы Дальнейшая профилировка и отладка Пересмотр структуры программы Применение приближенных методов задачи n тел Дополнительные связи на тело Теплообмен и теплоперенос Источники и стоки Реализация вспомогательных средств комплекса

Спасибо за внимание Авторы выражают благодарность Васенину В.А., Водомерову А.Н.