Вихреразрешающее моделирование

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

А.С. Монин определил турбулентноcть как "…завихренное течение с очень большим числом возбужденных степеней свободы и с хаотическим распределением дисперсионных соотношений и фазовых сдвигов" (Монин, Яглом 1991). Точный поведения прогноз поведения турбулентного течения невозможен. Однако можно попытаться построить численную модель, воспроизводящую некоторые статистические характеристики турбулентного движения.

Существует три основных подхода к численному моделированию турбулентности:

1. Прямое численное моделирование (в англоязычной литературе - direct numerical, DNS). В рамках этого подхода численно решаются системы алгебраических уравнений, с высокой точностью аппроксимирующие исходную систему дифференциальных уравнений Навье-Стокса (систему уравнений для вязкой жидкости или газа, задающую балансные соотношения для импульса и массы в каждой точке расчетной области). DNS (безусловно самый обоснованный и самый точный метод) требует детального пространственно-временного разрешения и, следовательно, ведет к большим вычислительным затратам. Поэтому DNS редко применяется в практических задачах и чаще всего служит инструментом для получения больших массивов данных о турбулентных потоках простой структуры.

2. Решение систем уравнений, осредненных по Рейнольдсу (Reynolds averaged Navier-Stokes, RANS). В данном случае модель воспроизводит только средние значения скорости (и, при необходимости, скалярных величин, например температуры и влажности воздуха), а влияние всех флуктуаций учитывается при помощи турбулентных замыканий. Под осреднением в RANS понимается осреднение по ансамблю состояний. В силу предполагаемой эргодичности и в случае наличия статистически однородных направлений модели RANS могут быть не только трехмерными, но и двумерными и одномерными. Последние, как правило, используются в качестве блоков, параметризующих турбулентность в геофизических пограничных слоях в глобальных, региональных и мезомасштабных моделях атмосферы и океана. Если турбулентность происходит на фоне какого либо медленного процесса, то возможно построение нестационарной модели RANS (unsteady RANS – URANS).

3. Вихреразрешающее моделирование или моделирование методом крупных вихрей (large eddy simulation, LES). Впервые этот подход был использован в 1970 году Дирдорфом (Deardorff, 1970) для моделирования турбулентного течения в канале, а затем для моделирования атмосферного пограничного слоя (Deardorff, 1973). Основой вихреразрешающего моделирования турбулентных потоков с очень большими числами Рейнольдса (Re 1) является подтверждаемая экспериментально гипотеза о независимости статистических характеристик крупномасштабных турбулентных движений от молекулярной вязкости. Согласно этому предположению, возможно построение численной модели, описывающей нестационарную динамику только относительно крупных вихрей (вихрей, вносящих основной вклад в энергию турбулентного потока и определяющих взаимодействие турбулентности со средним движением). При этом вычислительные затраты необходимые для реализации такой модели не должны зависеть от числа Рейнольдса, поскольку нет необходимости явно и точно рассчитывать все мелкие вихри.

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

Необходимым (но, в ряде случаев, недостаточным) условием корректной работы турбулентных замыканий является способность “подсеточной” модели правильно описывать диссипацию кинетической энергии сглаженных пульсаций скорости – то есть способность замыкания имитировать прямой энергетический каскад от крупных вихрей к мелким. Такой каскад является основным механизмом перераспределения энергии в инерционном интервале трехмерной однородной изотропной турбулентности (в области волновых чисел, в которой энергия не продуцируется и не диссипирует, а только перераспределяется между гармониками). В силу стремления мелкомасштабных флуктуаций скорости к изотропии можно надеяться, что среднестатистический эффект подсеточных взаимодействий можно описать, вводя дополнительный коэффициент вихревой вязкости по аналогии с молекулярной вязкостью. Модели, основанные на введении вихревой вязкости (EVM), составляют широкий класс турбулентных замыканий для LES (обзор возможных подходов и математическое обоснование метода можно найти, например в (Berselly et. al, 2006)).

Общим недостатком всех EVM является очень низкая корреляция (~0.2) моделируемых и реальных турбулентных напряжений и неспособность EVM описывать взаимодействия, приводящие к локальной генерации энергии крупных вихрей за счет мелкомасштабной подсеточной турбулентности. Это ограничивает область применения таких моделей условием “попадания” пространственного разрешения в инерционный интервал. Это условие не всегда достижимо, например, при моделировании турбулентности в областях, ограниченных твердыми границами, где генерация турбулентной энергии происходит на масштабах, сравнимых с расстоянием от стенки. Похожие проблемы возникают при моделировании турбулентных пограничных слоев при наличии устойчивой плотностной стратификации в поле силы тяжести.

Название "подсеточные” утвердилось за замыканиями благодаря тому, что на начальном этапе развития вихреразрешающих моделей считалось, что роль фильтра выполняет пространственная дискретизация, связывающая минимальный размер вихря с размером сеточной ячейки. Однако, вскоре было обнаружено (см. (Ghosal, 1996)), что численные ошибки аппроксимации нелинейных членов в конечно-разностных моделях могут существенно превосходить величину "подсеточных" турбулентных напряжений. Для конечно-разностных или конечно-элементных схем низкого порядка точности эти ошибки заведомо превышают вклад "подсеточных сил" практически во всем спектральном диапазоне, воспроизводимом моделью независимо от пространственного разрешения. Вполне естественно, что наилучшим выбором для LES моделей являются спектральные численные методы, наиболее точно аппроксимирующие весь воспроизводимый диапазон невозможностью их использования для расчетов турбулентных потоков в областях сложной конфигурации. Поэтому задача построения конечно-разностных или конечно-элементных LES моделей, а также разработки турбулентных замыканий, учитывающих особенности дискретизации в таких моделях, остается актуальной. Следует отметить, что в силу малости “подсеточных/подфильтровых” сил по сравнению с численными ошибками использование неконсервативных по энергии численных методов в LES нежелательно, поскольку систематическая ошибка, вносимая в энергию системы, может существенно исказить коротковолновый участок спектра. Известно, что линейные противопотоковые схемы, даже высокого порядка точности, неприменимы в LES.

Уменьшить роль численных ошибок можно, моделируя только относительно крупные вихри на достаточно мелкой сетке. При этом считается, что пространственный фильтр задан заранее, а его свойства известны. Такая концепция носит название “явной фильтрации” (Ghosal, 1996, Lund and Kaltenbach, 1995). Общепринятого подхода к построению LES, основанных на явной фильтрации, и замыканий к ним на данный момент нет. Известно, что модели, основанные на явной фильтрации более чувствительны к выбору замыкания. Отметим, что в данном случае, более уместно говорить не о “подсеточном”, а о “подфильтровом” турбулентном замыкании.

Как правило, турбулентные замыкания содержат набор констант. Эти константы могут быть либо определены из физических соображений, либо эмпирически подобраны по результатам численных расчетов. Однако, такие оценки не могут быть универсальными и даже в тех случаях, когда само замыкание оказывается применимым, параметрические константы приходится подбирать исходя из особенностей моделируемого потока. Кроме того, в силу большой чувствительности LES к численным схемам, физические соображения не всегда оправданы. Альтернативой могут быть динамические замыкания, впервые предложенные в работе (Germano et al, 1991). В динамических подсеточных и подсеточных/подфильтровых моделях турбулентности используется предположение об инвариантности констант, входящих в замыкание, относительно эффективной ширины фильтрации. Иными словами, предполагается, что одна и та же подсеточная модель без какого либо изменения может применяться как при более грубом, так и при более мелком пространственном разрешении.

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

Идея подобия масштабов, лежащая в основе динамического подхода, часто используется при конструировании самого вида турбулентных замыканий. Простейшая модель, построенная из соображений подобия (scale similarity model SSM), была впервые предложена и проверена на априорных тестах в работе Bardina et al, 1980). В этой модели тензор турбулентных напряжений аппроксимируется путем формальной замены неизвестной (нефильтрованной скорости) на известную (фильтрованную). Оказалось, что такая простая модель обладает рядом преимуществ по сравнению с моделями вихревой вязкости. В частности, она способна обеспечивать локальную генерацию турбулентной энергии, имитируя обратный каскад энергии от мелких вихрей к крупным. Это свойство SSM оказывается важным для моделирования анизотропных турбулентных потоков при реализации процесса самоорганизации мелкомасштабной турбулентности в квазикогерентные крупномасштабные структуры.

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

Будучи недостаточно диссипативными SSM используются в практических расчетах в комбинации с моделями вихревой вязкости. Такие модели принято называть смешанными. Смешанные модели возможно реализовать в динамической постановке (dynamic mixed model, DMM). Подробный обзор замыканий, базирующихся на идее подобия масштабов, приводится в работе (Meneveau and Katz , 1999). Вихреразрешающая модель нейтрально стратифицированнго пограничного слоя атмосферы с использованием смешанного динамического замыкания тестировалась в работе (Глазунов, 2006).

Интересно отметить, что в моделях подобия и смешанных моделях в отличие от моделей вихревой вязкости оператор фильтрации по пространству используется в явном виде. В том случае, когда фильтр имеет гладкую в спектральном пространстве функцию пропускания, SSM параметризует не только часть “подсеточных” взаимодействий, но и часть взаимодействий в разрешаемом диапазоне, подавленных фильтрацией. Это позволяет рассматривать LES, использующие смешанные замыкания, как модели с явной фильтрацией (Winckelmans and Jeanmart, 1999). При реализации таких замыканий в конечно-разностных LES снижается роль ошибок аппроксимации.

Моделирование сдвиговой пристеночной турбулентности при очень больших числах Рейнольдса

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

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

Рис.1 Схема численного эксперимента по моделированию турбулентного потока в канале, ораниченном двумя шероховатыми пластинами.

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

Рис.3 Нормированная дисперсия "реконструированных" пульсаций скорости, полученных в LES, близка к дисперсии пульсаций скорости в DNS (при максимально достижимых числах Рейнольдса) и слабо зависит от пространственного разрешения модели .

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

Рис.4 Моделирование обтекания куба, расположенного на шероховатой поверхности, турбулентным потоком воздуха. Слева - одна из мгновенных реализаций поля скорости; справа - осредненный по времени, динамически определенный коэффициент Смагоринского.

Вихреразрешающее моделирование геофизических пограничных слоев (пограничного слоя атмосферы и верхнего слоя океана)

Пограничный слой атмосферы (ПСА) является ключевым звеном климатической системы, поскольку с помощью механизмов турбулентного переноса он регулирует преобразование энергии солнечной радиации, поглощенной подстилающей поверхностью, в энергию движений в атмосфере, океане и других водных системах, контролирует уровень диссипации кинетической энергии и обеспечивает тепловлагоперенос между атмосферой и деятельным слоем суши. Известно, что наряду с мелкомасштабной хаотической трехмерной турбулентностью в пограничном слое атмосферы присутствуют когерентные движения с пространственными масштабами, сравнимыми с толщиной пограничного слоя. Крупномасштабные структуры в ряде случаев определяют основную часть потоков тепла, влаги и импульса в пограничном слое. Мелкомасштабная хаотическая трехмерная турбулентность оказывается в этом случае важной лишь постольку, поскольку она влияет на формирование и статистику крупномасштабных вихрей
В настоящее время для исследования ПСА широко используется подход, основанный на применении методологии вихреразрешающего моделирования. При этом основные энергонесущие масштабы движений в ПСА описываются явным образом, а влияние мелкомасштабной хаотической трехмерной турбулентности параметризуется. Аналогичный подход может применяться и для моделирования верхнего слоя океана (ВСО), а также для совместного моделирования двух взаимодействующих пограничных слоев ПСА и ВСО (Глазунов, Лыкоссов 2003).

Рис.5 Совместная конвективная циркуляция двух взаимодействующих пограничных слоев, рассчитанная при помощи совместной модели ПСА и ВСО. Стрелками указано направление ветра и течений, цветом обозначены аномалии температуры воздуха и воды.

Параллельные вычисления в задачах вихреразрешающего моделирования

Классической областью применения вихреразрешающих моделей ПСА является статистический анализ структурной турбулентности (квазиупорядоченных вихревых структур) с целью построения параметризаций для мезо-масштабных и глобальных моделей атмосферы. В таких задачах, как правило, применяются относительно небольшие горизонтальные размеры расчетной области (до нескольких километров) в несколько раз превышающие масштаб вихревых структур. При этом число узлов сетки в модели относительно невелико (около 10^5). Поэтому такие задачи не требуют больших вычислительных ресурсов и могут быть успешно решены на последовательных компьютерах. С другой стороны, вихреразрешающее моделирование может рассматриваться как самостоятельная прикладная проблема. Модели ПСА, основанные на этой методологии, применимы для таких задач как перенос примесей в атмосфере на масштабах до нескольких десятков километров, обтекание крупных препятствий (особенностей рельефа, зданий) и т.д. В этом случае возникает необходимость увеличения расчетной области модели до размера, определяемого конкретной прикладной постановкой, а шаг сетки модели не может произвольно увеличиваться в связи с требованием «попадания» в инерционный интервал. Соответственно сказанному, число расчетных узлов сетки модели значительно увеличивается. Например, для решения задачи о переносе загрязняющих примесей в ПСА на площади 30 x 30 км^2 необходимо использовать модель, содержащую порядка 10^7 расчетных точек сетки. Это налагает высокие требования к эффективности численных алгоритмов и вызывает необходимость использования параллельных вычислений в задачах подобного рода. Параллельные реализации вихреразрешающих моделей, разработанных в ИВМ РАН, ориентированы на использование суперкомпьютеров с распределенной памятью. Параллельные вычисления организованы следующим образом (см. рис. 6).

Рис.6 Организация параллельных вычислений в вихреразрешающих моделях ПСА и ВСО.

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

В ИВМ РАН разрабатывается параллельная версия модели ПСА, предназначенная для расчета обтекания городской застройки турбулентными потоками воздуха. Эта модель позволяет за приемлемое время провести серию расчетов с целью определения статистических характеристик ветрового потока между зданиями при различных значениях и направлениях фонового ветра. Модель включает в себя блок расчета лагранжева переноса частиц примеси. Этот блок позволяет проводить расчеты траекторий десятков миллионов частиц одновременно с моделированием динамики турбулентного потока.

Литература.

Монин А.С., Яглом А.М. Статистическая гидромеханика: теория турбулентности. // М.: Наука. Гл. ред. физ.-мат. лит., 1992. 695 с.

Deardorff J.W. A numerical study of three-dimensional turbulent channel flow at large Reynolds numbers. // J. Fluid Mech., 1970, 41, pp. 453-480

Deardorff J.W. The use of subgrid transport equations in a three-dimensional model of atmospheric turbulence // Journal of Fluids Engineering. 1973. V. 9. P. 429-438.

Berselli L.C., Iliescu T., Layton, W.J. Mathematics of Large Eddy Simulation of Turbulent Flows. // Springer. Series: Scientific Computation. 2006, XVIII, 348 p.

Ghosal S. An analysis of numerical errors in large-eddy simulations of turbulence. // J. Comput. Phys., 1996, 125, 187–206.

Lund, T. S., Kaltenbach, H.-J. Experiments with explicit filtering for LES using a finite-difference method. // Center for Turbulence Research, Annual Research Briefs 1995, pp. 91-105

Germano M., Piomelli U., Moin P. and Cabot W. H. A dynamic subgrid-scale eddy viscosity model // Phys. Fluids. A, 1991, 3, pp. 1760-1765.

Bardina J, Ferziger J.H., Reynolds W.C. Improved subgrid scale models for large-eddy simulation. // Am. Inst. Aeronaut. Astronaut., 1980, Paper 80-1357.

Meneveau C, Katz J. Dynamic testing of subgrid models in LES based on the Germano identity.// Phys. Fluids, 1999 11 pp. 245–47

Winckelmans, G. S., and H. Jeanmart Assessment of some models for LES without/with explicit filtering.// Direct and Large-Eddy Simulation IV, B. J. Geurts, R. Friedrich, and O.M'etais, Eds., Kluwer, 2001, pp. 55–66.

Glazunov A.V. , V.N. Lykossov. Large-eddy simulation of interaction of ocean and atmospheric boundary layers. - Russ. J. Numer. Anal. Math. Modelling, 2003, 18, pp. 279-295.

Глазунов А.В. Моделирование нейтрально стратифицированного турбулентного потока воздуха над горизонтальной шероховатой поверхностью., 2006. Изв. РАН, ФаиО, 42, No. 3, 307–325.