Моделирование течений жидкости и газа с поверхностью раздела сред, турбулентностью и стратификацией тема автореферата и диссертации по механике, 01.02.05 ВАК РФ
|
Яковенко, Сергей Николаевич
АВТОР
|
||||
|
доктора физико-математических наук
УЧЕНАЯ СТЕПЕНЬ
|
||||
|
Новосибирск
МЕСТО ЗАЩИТЫ
|
||||
|
2014
ГОД ЗАЩИТЫ
|
|
01.02.05
КОД ВАК РФ
|
||
|
|
||||
На правах рукописи
Яковенко Сергей Николаевич
МОДЕЛИРОВАНИЕ ТЕЧЕНИЙ ЖИДКОСТИ И ГАЗА С ПОВЕРХНОСТЬЮ РАЗДЕЛА СРЕД, ТУРБУЛЕНТНОСТЬЮ И СТРАТИФИКАЦИЕЙ
01.02.05. — Механика жидкости, газа и плазмы
АВТОРЕФЕРАТ диссертации на соискание ученой степени доктора физико-математических наук
15 ЯКЗ 2015
005557441
Новосибирск - 2014
005557441
Работа выполнена в ФГБУН Институте теоретической и прикладной механики им. С.А. Христиановича СО РАН
Официальные оппоненты:
Исаев Сергей Александрович - д.ф.-м.н., профессор,
ФГБОУ ВПО Санкт-Петербургский государственный университет
гражданской авиации, профессор кафедры механики
Терехов Виктор Иванович - д.т.н., профессор,
ФГБУН Институт теплофизики им. С.С. Кутателадзе СО РАН,
заведующий отделом термогазодинамики
Платов Геннадий Алексеевич - д.ф.-м.н.,
ФГБУН Институт вычислительной математики и математической геофизики СО РАН, ведущий научный сотрудник лаборатории математического моделирования процессов в атмосфере и гидросфере
Ведущая организация:
Федеральное государственное автономное образовательное учреждение высшего профессионального образования
Московский физико-технический институт (государственный университет)
Защита состоится 27 марта 2015 г. в 1400 на заседании диссертационного совета 003.035.02 в ФГБУН Институте теоретической и прикладной механики им. С. А. Христиановича СО РАН по адресу: 630090, Новосибирск, ул. Институтская, 4/1
С диссертацией можно ознакомиться в библиотеке ФГБУН Института теоретической и прикладной механики им. С. А. Христиановича СО РАН и на сайте http://www.itam.nsc.ru/ru/thesis/
Автореферат разослан «_» декабря 2014 г.
Ученый секретарь диссертационного совета доктор технических наук
ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ
Актуальность темы. В большинстве течений жидкости и газа присутствуют поверхности раздела (твердые стенки, поверхности раздела текучих сред), имеющие сложную геометрию, за исключением самых простых случаев. Кроме того, окружающая среда, как правило, стратифицирована по плотности за счет изменений температуры и концентрации примеси. При некоторых условиях, например, при больших числах Рейнольдса (Re), Прандтля, Шмидта (Se), градиенты плотности оказываются достаточно резкими, приводя к образованию «поверхностей раздела» слоев разной плотности в однофазной среде. В течениях почти всегда есть не только стратификация, но и сдвиги скорости, например, в пограничном слое или в слое смешения потоков разной скорости.
Эти особенности течений жидкости и газа, безусловно, затрудняют их теоретическое, экспериментальное и численное исследование. Однако наиболее сложный вызов — свойство течений к переходу от спокойного ламинарного состояния к нерегулярному турбулентному, реализуемое в большинстве окружающих нас явлений, в частности, в упомянутых выше сдвиговых слоях. Аналитическое решение законов сохранения (для частиц среды) получено лишь в некоторых простых ситуациях, в основном при малых Re. Физические эксперименты по изучению турбулентности в окружающей среде — чрезвычайно сложная задача из-за очевидных ограничений возможностей измерительной аппаратуры в натурных условиях и упрощений в лабораторных условиях.
В частности, турбулентность в области обрушения внутренних волн исследовалась в основном в гидродинамических каналах с затопленным препятствием, двигаемым вдоль свободной поверхности [Castro I.P., Snyder W.H. // JFM. 1993. Vol. 255. P. 195-211]. Продолжительность опыта, однако, недостаточна для достижения квазистационарной структуры развитой турбулентности: через некоторое время после его начала эффекты отражения волн от входной и выходной границы искажают реальную картину. То есть в каналах конечной длины не воспроизводится ситуация обрушения внутренних волн в атмосфере и океанах (где горизонтальные размеры много больше высоты). С другой стороны, численное моделирование этого явления и анализ полученных новых данных могут дать ключ к пониманию слабо изученных механизмов возникновения и поддержания геофизической турбулентности, наблюдаемой на больших высотах от подстилающей поверхности и элементов топографии.
Развитие компьютеров в последние полвека сделало возможным численное решение уравнений механики жидкости и газа, на основе различных типов моделирования, для получения недостающих данных о природе физических процессов. Подход DNS разрешает явно все турбулентные масштабы, но ограничен небольшими Re из-за повышенных требований к ресурсам компьютера. Модели RANS, наоборот, вводят при осреднении по всем турбулентным масштабам, но при этом теряется информация о нестационарных свойствах турбулентных пульсаций, хотя знания осредненных величин обычно достаточно в практических задачах. Полученные RANS-методом распределения часто отличаются от данных опытов и DNS, стимулируя полуэмпирические коррекции ги-
потез замыкания. Промежуточный между DNS и RANS, подход LES явно разрешает крупные масштабы и включает подсеточные модели для аппроксимации мелкомасштабных движений, позволяя воспроизвести течения при больших Re, чем DNS. Однако в LES нет адекватного обратного контроля эффектов мелких масштабов на крупные, и возникают сильные требования к измельчению сетки в слоях у стенки. Для преодоления трудностей описания отрывных течений при высоких Re, где ни LES, ни RANS не способны в одиночку справиться, предложены гибридные RANS-LES подходы, находящиеся на ранней стадии развития [Spalart P.R. // Ann. Rev. Fluid Mech. 2009. Vol. 41. P. 181-202].
Обзор существующих методов и сопутствующих проблем численной реализации и адекватного представления данных показывает, что в настоящее время развитие универсального, эффективного и точного инструмента для воспроизведения течений жидкости и газа по-прежнему является нерешенной актуальной проблемой. Тем не менее «классические» методы (RANS, LES, DNS) способны помочь в решении многих задач механики сплошной среды, которые могут усложняться эффектами поверхностей раздела, стратификации, закрутки, фазовых переходов, химических реакций, радиационного переноса и т.д.
Цели работы:
- развитие иерархии адекватных RANS-моделей уравнений переноса для моментов второго порядка, их верификация в канонических тестовых течениях с устойчивой стратификацией, свободной поверхностью жидкости, твердой стенкой, отрывом и рециркуляцией из-за наличия препятствия;
- исследование путем DNS и LES процесса обрушения внутренних волн в стратифицированном течении с препятствием, развития неустойчивости и турбулентности; анализ величин, найденных при осреднении мгновенных полей скорости и плотности, их сравнение с результатами RANS-моделей;
- развитие адекватных способов, позволяющих получить разрешение поверхности раздела несмешивающихся текучих сред, и апробация этих способов в задачах обрушения плотины и неустойчивости Рэлея-Тейлора;
- эффективная организация серий предварительных и параметрических расчетов, усовершенствование алгоритмов решения и способов обработки полученных данных для повышения качества (точности) моделирования.
Научная новизна работы:
1. Развита иерархия аппроксимаций RANS-метода, проведена их верификация в канонических течениях с устойчивой стратификацией, стенкой, свободной поверхностью, отрывом и рециркуляцией, вводимой препятствием. Показано, в частности, что модифицированная модель второго порядка описывает противоградиентный характер турбулентных потоков тепла и импульса, наблюдаемый в устойчиво стратифицированном течении в открытом канале.
2. Путем решения уравнений Навье-Стокса и уравнения для функции объемной фракции, с помощью континуальной модели поверхностного натяжения, адекватно воспроизведена эволюция неустойчивости Рэлея-Тейлора (НРТ) на поверхности раздела несмешивающихся текучих сред, в том числе для реальных случаев «вода - воздух», «вода - бензол». Показано, что при небольшом перепаде плотности на нелинейном этапе наблюдаются эффекты вторичной не-
устойчивости Кельвина-Гельмгольца (НКГ) с грибовидными конвективными структурами, а при большом перепаде плотности более плотная среда тонкими струями проникает глубоко вниз между толстыми колоннами менее плотной среды, поднимающимися вверх. Полученные новые результаты моделирования НРТ подтверждают и уточняют данные предшествующих работ.
3. Выполнено прямое численное моделирование (DNS) и моделирование крупных вихрей (LES) обрушения внутренних волн, генерируемых препятствием в устойчиво стратифицированном течении при Re = 4000, намного большем, чем в предыдущих работах, и 1 < Se < 700. Впервые показано, что в начальной стадии развиваются мелкомасштабные возмущения НРТ с длиной волны, зависящей от Se и согласующейся с теоретическим значением для наиболее неустойчивой моды на поверхности раздела слоев разной плотности.
4. Изучена сложная структура течения в области опрокидывания внутренних волн, с изогнутыми квазидвумерными элементами, генерирующими ряд вихревых трубок на поверхности цилиндра, который опоясывает эту область и содержит резкие градиенты плотности. Взаимодействие матрицы структур НРТ и НКГ в области обрушения волн порождает развитую турбулентность с протяженным инерционным интервалом «—5/3 » на спектрах полей скорости и скаляра. На поздних стадиях развития идентифицированы долгоживущие тороидальные вихри, которые соответствуют крупномасштабной моде, представляющей собой наиболее неустойчивое возмущение исходной системы двумерной ламинарной пары вихрей в области опрокидывания волн.
5. Показано, что возникающая при обрушении волн турбулентная область является квазистационарной и имеет баланс порождения сдвигом, адвекции и диссипации энергии турбулентности. В рассматриваемой области выполнены анализ статистических моментов, полученных из данных DNS, и априорная проверка RANS-моделей уравнений напряжений Рейнольдса, а также сделаны оценки геофизически важных величин, в частности, эффективности перемешивания, равной 0.2, как и в океанологических приложениях.
Достоверность полученных результатов подтверждена при сравнении с аналитическими решениями, данными измерений и расчетов других авторов.
Практическая значимость работы состоит в применимости развитых моделей и методов к самым разнообразным приложениям техники и окружающей среды, в частности, в задачах гидрометеорологии на микромасштабах (с обтеканием сооружений и рельефа местности) и мезомасштабах (с учетом физико-химических процессов в атмосфере и водоемах), прогнозирования изменений в биосфере под воздействием антропогенных нагрузок, проектирования различных технических устройств с проточными трактами, нефте- и газопроводов, летательных и плавательных аппаратов, автомобилей, ж/д транспорта.
Основные положения, выносимые на защиту:
— результаты развития аппроксимаций RANS-метода и их верификации в канонических течениях со стратификацией, свободной поверхностью, стенкой, в частности, воспроизведение модифицированной моделью второго порядка противоградиентных турбулентных потоков тепла и импульса в устойчиво стратифицированном течении в открытом канале;
— результаты эволюции неустойчивости Рэлея-Тейлора на поверхности раздела несмешивающихся текучих сред, в том числе и для реальных случаев «вода - воздух», «вода - бензол», показывающие эффекты влияния перепада плотности сред и поверхностного натяжения;
— результаты DNS/LES обрушения внутренних волн, инициированных препятствием в стратифицированном потоке при 1 < Sc < 700 и Re = 4000;
— сценарий развития неустойчивости при обрушении волн, который включает рост мелкомасштабных возмущений неустойчивости Рэлея-Тейлора в начальной стадии, появление вытянутых и изогнутых квазидвумерных элементов на поверхности цилиндра (опоясывающего область опрокидывания волн) и ряда соответствующих вихревых трубок в промежуточной стадии, цепочки долгоживущих тороидальных структур на поздних стадиях;
— получение состояния развитой турбулентности с протяженным инерционным интервалом в возникающей при обрушении волн квазистационарной области, где имеется глобальный баланс скорости порождения сдвигом, адвекции и диссипации энергии турбулентности, а эффективность перемешивания согласуется с наблюдаемой в геофизических приложениях.
Личный вклад автора. Результаты, представленные в работе, получены лично автором. Исследования, результаты которых приведены в диссертации, выполнены в соавторстве с А.Ф. Курбацким, В.В. Ларичкиным, К.С. Chang, I.P. Castro, T.G. Thomas. Автор написал программы расчета, провел численные исследования, участвовал в постановке задач, анализе литературных ссылок и полученных данных, подготовке публикаций к печати.
Апробация работы. Результаты диссертационной работы докладывались и обсуждались на семинарах ИТПМ СО РАН, ИВТ СО РАН, ИГиЛ СО РАН, ИВЭП СО РАН, ИПМ РАН, НИИМ МГУ, МФТИ, Университета Саутгемптона, видеосеминаре по аэромеханике ЦАГИ - ИТПМ СО РАН - СПбГПУ - НИИМ МГУ, а также на всероссийских съездах по теоретической и прикладной механике (1991, 2001, 2011), международных конференциях по методам аэрофизических исследований (1992-2014), вычислительной гидродинамике (1999, 2002), по турбулентности (2000, 2009, 2011, 2012, 2013), гидромеханике (2008, 2010), международных симпозиумах по тепломассопереносу (2000, 2009), вычислительному теплопереносу (2012), 1-м Международном симпозиуме по явлениям турбулентности и сдвиговых течений (1999) и других научных мероприятиях. Работа по теме диссертации поддержана интеграционными проектами СО РАН, грантами РФФИ, Инженерно-физического научно-исследовательского совета Великобритании и других научных фондов. Соискатель получил премию СО РАН им. H.H. Яненко (1995) за цикл работ в области прикладной и вычислительной математики, стипендию Московского центра фундаментальной физики в Москве (1995-1996), государственную научную стипендию (1997-2000).
Публикации. По теме диссертации имеется 108 печатных работ, в том числе 29 статей в рецензируемых журналах из списка ВАК.
Структура и объем работы. Диссертация состоит из введения, 4 глав, заключения и списка литературы из 375 наименований. Общий объем работы составляет 340 страниц, включая 126 рисунков и 4 таблицы.
СОДЕРЖАНИЕ РАБОТЫ
Во введении обсуждается актуальность и практическая значимость численных исследований течений с поверхностью раздела, турбулентностью и стратификацией. Отмечено, что до сих пор нет универсальной модели, способной с необходимой точностью воспроизвести детали интересующих течений и применимой к большому кругу задач. На основе обзора существующих методов расчета турбулентных течений обоснована необходимость применения «классических» подходов (RANS, LES, DNS) и их развития для получения приемлемых результатов. В частности, для задачи обрушения внутренних волн в стратифицированном течении анализ полученных новых данных моделирования мог бы способствовать пониманию процессов возникновения и поддержания геофизической турбулентности. В таком анализе автор впервые обнаружил и исследовал явления НРТ при опрокидывании волн и полагает полезным сравнить их с эволюцией неустойчивости на поверхности раздела несмешивающих-ся сред, изучение которой также является одной из целей работы.
Явление НРТ при различных условиях активно исследовалось, начиная с работ Тейлора, Льюиса, Юнгса, в том числе и отечественными авторами. В частности, Руев, Федоров, Фомин и Демьянов, Долуденко, Иногамов, Сон изучали некоторые режимы НРТ в аномальных условиях, Герценштейн и др. моделировали НРТ в ячейке Хеле-Шоу. Вопросы численного моделирования обрушения внутренних волн в отечественной литературе не представлены (насколько известно автору), что указывает на сложность проблемы для численной реализации. С другой стороны, имеется множество работ по общей тематике моделирования стратифицированных течений (в том числе и с внутренними волнами) различными методами. В частности, Чашечкин, Байдулов, Гущин, Матюшин и др. провели численные, аналитические и экспериментальные исследования в течениях с внутренними волнами вокруг сферы, цилиндра и в других конфигурациях. Черных, Воропаева, Мошкин, Фомина, Дружинин и др. моделировали турбулентные стратифицированные следы и другие потоки путем RANS и DNS. Ткаченко, Гурьев и др. выполнили расчеты стратифицированных течений, в том числе в морской среде с подводной возвышенностью и внутренними приливными волнами, при использовании URANS и LES.
Помимо целей работы и краткого обзора исследований других авторов дана информация о новизне результатов, структуре и апробации работы.
В первой главе приведены уравнения для описания течений несжимаемой среды и иерархия методов их решения. В 1.1.2 кратко изложены основные типы подсеточных моделей LES, с акцентом на их использование для стратифицированных течений, а в § 1.2-1.4 более подробно рассмотрены аппроксимации RANS-моделей, вопросы их развития и применения.
Уравнения для осредненных по Рейнольдсу величин (давления, скорости Uh скаляра F) в приближении Буссинеска записываются как
ртг птг агг 8U.U. я я/7
~^-ccgiF, (1)
р Эх,
PF Dt '
8F 8U.F
!--1--'—■■
5t 5X;
5_ dx,
. 8F -тт.
Л--uj
dx,, J
где для напряжений Рейнольдса и турбулентных потоков скаляра (тепла, вещества) необходима модель турбулентности, замыкающая (1)-{2).
Простейшие гипотезы замыкания основаны на градиентной гипотезе Бус-синеска с изотропными коэффициентами турбулентного обмена V, и X,:
du, dUj —- + —'-
дх, дх,
3 дх, ст, дх,
(3)
где к = и\и\12 - турбулентная кинетическая энергия (ТКЭ), <х, = V, /Я, - турбулентное число Прандтля (Шмидта), а v, и X, определяются из моделей турбулентности различного уровня замыкания. В практических приложениях получила распространение, в основном, двухпараметрическая к-е-модель:
Dk = Dt ~
De Dt
дх,
_д_ дх,
v + -
К
v + -
'kj
дк_ дх,
де дх,
+ P + G-S, Р = -и\и'^, G = —agfU'jf,
dXj
+ {Ctl{P + Ct3G)-Ct2e}j, У,=С^к2/е,
(4)
(5)
где а, полагают константой или эмпирической функцией, а стандартные значения коэффициентов [Методы расчета турбулентных течений. М.: Мир, 1984] получены эмпирически при расчетах свободных и пристеночных-течений.
Однако двухпараметрические модели, в частности ¿-е-модель (1)-(5), не описывают нелокальность турбулентного переноса в сложных течениях (напр., при наличии массовых сил, стенки, свободной поверхности), которая воспроизводится с привлечением уравнений переноса для моментов второго порядка:
Duft _ д Dt ~ дх„
du^u'j дх1
Du, f Dt
dT1 _
Dt
дх,
О»
1, ч du\f -Iv + Â)—'— 2K ' дх,
Qijf
+ Pif+Gif+7tif-£if
■agiU'jf -agjU'J', Gif=-agY\
df'2 dUjf*
dt
dx,
_d_
' dx,
дх,
+pf
pf=
-2
dx,
(6)
(7)
(8)
P^-uy^dUjdx^-uy^dUjdx,), Pif =-u'it'l{dFldxl)-U'f{dUildx]),
где производную по времени и конвекцию (D/Dt) балансируют молекулярная и турбулентная диффузия (в квадратных скобках), порождения сдвигом скорости и скаляра (Ру, Pif, Pj) и силами плавучести (G,y, G,/), корреляции с пульсациями давления (ху, л-,/) и диссипативные члены (%, е,у, £/). Уравнения (1)-(2), (6)-(8) образуют полную модель турбулентного переноса второго порядка после введения замыкающих аппроксимаций для неизвестных статей баланса в (6)-(8), в частности, корреляций щ и жиграющих ключевую роль в модели. В 1.3.1-1.3.5 приведены стандартные выражения для е,у, е,у, щ, щ/ (модель 2) и их модифи-
кации, удовлетворяющие предельному переходу турбулентности в двумерное состояние (модель 1). В предположении локальной изотропии при больших турбулентных числах Рейнольдса обычно принимают г,у = (2/3)е<5,у , 'н/ = 0. Решение уравнения Пуассона для пульсаций давления дает интегралы для щ и гг,/, с вкладами различных физических процессов. В однородной турбулентности 7г,у имеет только слагаемое 4', обусловленное взаимодействием пульсаций скорости. Аппроксимация Ротта, линейная по тензору анизотропии а,у, имеет вид
4" = (Ф)Щ - 0/ЗН*} - , (9)
где С\ - эмпирическая константа. При наличии значительных сдвигов средней скорости и массовых сил (стратификации) выражение (9) «стремления к изотропии» обычно дополняется аппроксимациями «изотропизации порождения»:
42) + 43' = -С2 [ру - (2/3)4^] - Сз \рц ~ (2/3)5^]. (10)
По аналогии с л:) = 4" + + л]р моделируется корреляция ж,/ в (7):
л, = 4" + 42) + л? = -Сч {е!к)ЦГ' + С2Г^Т'(ди,1дх]) - С3/Су. (11) Дефекты выражения (9) — неуниверсальность С] и неучет нелинейных свойств 7г,у. Для описания демпфирования турбулентности и ее перехода к сильно анизотропному, приближенно двумерному состоянию под действием стенки, свободной поверхности, устойчивой стратификации и других эффектов доминирующий вклад я,'1' должен обращаться в нуль, т.е. С\ —> 0 (при подавлении
компоненты и'п члены и\и'п -уравнения или их комбинация должны исчезать).
Модифицированная модель 1 отличается от стандартной модели 2 более корректными соотношениями для л^ и еу. Недостатки линейного выражения (9) преодолеваются в «квазиквадратичной» аппроксимации вида
4° =-С114АА2Е[а0+ С1в(а,ка^-{1/3)5^)], Л2=ам1Г Аг=а^а]каи, (12) где параметр анизотропии /1 = 1 — (9/8)(А2 - А3) и инварианты А2 и Аъ нелинейно зависят от а,у. В модели 1 также применяется анизотропная аппроксимация:
г- 2 / и'м', +3-11 .и' +8,„и'.и' +дИи'м'
*;=" <13>
Нелинейные выражения (12)—(13) сводятся к стандартным из модели 2 при слабой анизотропии (А —<► 1) и удовлетворяют предельному поведению при и'п—> 0 (А —> 0). То есть, модифицированная модель 1 способна без дополнительной эмпирической коррекции воспроизвести анизотропию напряжений Рейнольдса при подавлении пульсации и'п свободной (твердой) поверхностью или массовыми силами. В стандартной модели 2 влияние свободной поверхности жидкости (стенки) учитывается «приповерхностными» поправками в 7г,у (и к,]) с демпфирующей функцией /т{хп), нейтрализующими дефекты линейной связи (9) и изотропной модели е,у= (2/3)е<5,у. Такой способ оказывается неподходящим для течений, ограниченных поверхностями сложной топографии, где проблематично выделить расстояние до поверхности х„, имеющее физический смысл.
Члены турбулентной диффузии, включающие тройные корреляции, выражаются в 1.3.6 через моменты второго порядка при помощи градиентных моделей. Если процессы конвекции и диффузии в (6)-(8) несущественны, то ими пренебрегают или выражают алгебраически через известные величины в предположении локального равновесия для /'2, пропорциональности переноса г/'н'
переносу ТКЭ и переноса вектора it'J' - переносу величины \]к f12 :
Pf-Ef=0, Ъ+^е* P + G-e P+G-e
7 f u]u) ,J к u\f f 2 к
где a{j = 1 в неравновесном приближении и atj = 0 в локально-равновесном. Из
(14) моменты второго порядка определяются с учетом точного вида порождений и аппроксимаций Жу, жу, Еу , с,у, е/. Стандартная модель 2 дает (без учета стенки и свободной поверхности, т.е. при f^ =0) искомые величины в виде
Jз. I ('-^fri-p/^j-i'-^fe-pm6] П5)
к 3 iJ C,£ + a..[P + G-f]
к (i-c2/)(-^7)(a^,)+(i-c3/)G,/ к 8F
UJ =----7-—77i7v\—\~Б~Гп-т-' J =~2K~UJ (1Ь>
е C{fe + (\l2)ajfyP + G-e\ s dxi
где R = rfjт = (f'2/ef)-(s/к) = const, согласно аппроксимации для r.f, обседаемой в 1.3.7. Выражения (15)-(16) замыкают (1)-(2) с помощью уравнений вида (4)-(5). Такая модель с алгебраическими уравнениями для моментов второго порядка (модель 3), где обычно а,у = S:j и atf = 0, является более общей, чем
стандартная ¿-е-модель (1}-(5) с изотропными коэффициентами турбулентного обмена. В плоских (осесимметричных) сдвиговых течениях модель 3 приводит к сложным анизотропным функциям v,= СДР,е,.. .)i2/c и Л, = vt/at{P,E,...), которые учитывают эффекты массовых сил, поверхностей и других факторов, отражаемых источниками в исходных уравнениях переноса (6)-(8).
Во второй главе проведена верификация развиваемых RANS-моделей, сформулированных в § 1.2-1.3, при сравнении с данными опытов для течений со стратификацией (в слое смешения, поверхностной струе, открытом канале) и с препятствием квадратного сечения (в пограничном слое, закрытом канале).
В § 2.1 модели 1-3 упрощаются для плоских установившихся (с точки зрения осредненных величин) тонких сдвиговых слоев. Коэффициенты моделей соответствуют найденным предшественниками в сравнении с данными опыта в канонических тестах и путем численной оптимизации. Например, основные параметры С] и С2 моделей 2 и 3, влияющие на перераспределение турбулентных напряжений, удовлетворяют связи 1 - С2 ~ 0.2 Ci для согласия с величинами ау, измеренными в свободном локально-равновесном течении со слабым сдвигом скорости. Значения Сц = 2.5 и C\q = 1.2 определены при калибровке модифицированной модели 1 в нейтральном течении в канале. Для искомых функций заданы краевые условия в соответствии с типом границ: (1) свободная граница, отделяющая турбулизированную от невозмущенной среды; (2) свободная no-
верхность, т.е. поверхность раздела жидкости и газа; (3) стенка (дно канала). Система определяющих уравнений и граничных условий представлена в безразмерном виде при помощи входных параметров. Дифференциальные уравнения переноса, имеющие параболизованный вид в переменных «координата -функция тока», решены численно методом контрольного объема.
В § 2.2 приведены результаты расчета слоя смешения жидкостей разной плотности (солености). Модель, промежуточная по сложности между моделями 2 и 3, включает уравнения для горизонтальной скорости С/, функции плавучести В(= F) = —g(p-р0)/р0, к, е (как в ¿-е-модели 3), уравнения (6), (9)—(10) для напряжений Рейнольдса (как в модели 2 второго порядка) и выражения (16) для корреляций и'р, Ь2 (с пульсацией плавучести /' = Ь) при а^ = 1. Возможность
использования (16) вместо (7)-(8) подтвердили тестовые расчеты.
Сначала моделируется слой смешения постоянной плотности - канонический тест для развития моделей турбулентности. Получено, что применение выражения е где масштаб Ь пропорционален характерной ширине слоя,
дает лучшее согласие с экспериментом на высокоскоростной стороне слоя, чем уравнение переноса для е. В условиях устойчиво стратификации вычисленные профили удовлетворительно описывают измеренные (напр., см. рис. 1, где профили В(у) и Ь'(у) = (Ь2)хп подтверждают обнаруженную в опытах асимметрию смешения, вызываемого крупномасштабными турбулентными движениями).
Рис. 1. Профили средней плавучести и интенсивности пульсаций плавучести: точки - измерения, линии - настоящий расчет в различных сечениях слоя смешения.
Развитие слоя смешения с ростом координаты х характеризуют толщины 6и=и0/\ди1ду\п^, ^я=50/|55/ау|тах, где В0 =-g(pí-р0)/р0, 170 - скорость потока плотности р, над слоем смешения. На начальном этапе при слабом влиянии стратификации в измерениях и расчете наблюдается, как и при постоянной плотности, линейный рост этих величин. При х>Ь5 под действием сил плавучести скорость роста сдвигового слоя замедляется. Вычисленное поведение ди и 8В не согласуется с наблюдаемым в опыте уменьшением (коллапсом) толщины слоя. Последний эффект связан, вероятно, с неразвитостью мелкомасштабной структуры турбулентного течения и выбросом из слоя смешения
части вовлеченной жидкости, не перемешанной до молекулярных масштабов. Рассматриваемое течение носит перемежающийся характер, учет которого в ЯА^-моделях весьма нетривиален и сдерживается из-за отсутствия необходимых опытных данных. Вычисленные зависимости максимальных по поперечному сечению уровней пульсаций качественно согласуются с измеренными, в частности, воспроизводя степенные законы при х > 5Ьх, подобные законам затухания однородной турбулентности за решеткой. Медленное падение с ростом х и задержка выхода вычисленных кривых на степенные функции связаны с расхождением толщин слоя в расчете и опыте. Другая возможная причина медленного уменьшения и'т и Ъ'т - неуниверсальность алгебраических выражений (уравнений переноса) для Е/ и е, «не чувствующих» стратификацию.
В § 2.3 представлены результаты расчета устойчиво стратифицированного течения в канале со свободной поверхностью, нагретой сверху. Моделирование проведено в два этапа в соответствии с условиями опытов: а) получение развитого изотермического течения (динамическая задача), б) наложение стратификации на развитое течение (сопряженная динамическая и тепловая задача). Модель 3 с демпфирующей функцией, учитывающей воздействие свободной поверхности и стенки, описывает развитое течение в канапе при нейтральных условиях. Для воспроизведения противоградиентных потоков в условиях сильной устойчивости применены модели 1 и 2, которые в изотермическом течении дают, как &-е-модель 3, приемлемое согласие с опытом, в частности, демпфирование поверхностью вертикальных турбулентных пульсаций скорости (рис. 2).
Рис. 2. Характеристики развитого течения в канале. Сплошная, штриховая линии, пунктир, штрихпунктир - настоящие расчеты по моделям 1, Г, 2, 3 соответственно (модель Г включает, вместо аппроксимации (10), кубическую модификацию я^2').
При решении сопряженной задачи начальные профили средней скорости и температуры взяты в виде функциональных зависимостей, соответствующих измеренным распределениям в приближенно развитом устойчиво стратифицированном течении в канале. Это позволяет адекватно учесть влияние стратификации, проследив эволюцию турбулентных потоков импульса и тепла с ростом х, их изменения по поперечной координате у и локальному градиентному числу
Ричардсона Ri = ag(dT/8y)/(dU/dy)2. Полные модели второго порядка (1 и 2) воспроизводят противоградиентность (рис. 3) потока тепла vO, а модифицированная модель 1-й потока импульса uv, зафиксированную в опыте и DNS. Двухпараметрическая модель 3 не отражает эти эффекты сильно устойчивой стратификации даже качественно. Модель 1 с нелинейной модификацией (12) в основном правильно описывает динамику корреляций с пульсациями давления: уменьшение их величин в условиях сильно устойчивой стратификации по сравнению с нейтрально и слабо устойчивой, изменение знака луг при Ri > 0.5. Поведение пульсационных характеристик качественно согласуется с результатами DNS, обнаруживая подавление стратификацией с ростом Ri, немонотонность по х при Ri > 0.5 и достижение «стабилизированного» состояния при х/д ~ 20.
эе/<У
ю
Z Rc 4
Рис. 3. Корреляционные коэффициенты в нагретом сверху открытом канале, полученные по модели 1 при разных Ш в зависимости от х/б, где 6 - глубина канала (а); по моделям 1, 2, 3 при 0О2 = 0.002(- дТ/ду) (сплошная линия, штрихпунктир, длинный штрих) и модели 1 при = 0 (короткий штрих) в сечении х = 55 в зависимости от Ш (б).
В § 2.4 модели 1-3 верифицированы для струи, распространяющейся по поверхности более плотной жидкости. Вычисленные по моделям 1 и 2 поперечные профили характеристик струи соответствуют измеренным. Модифицированная модель 1 адекватно описывает данные опыта для зависимости от х максимумов средних величин ит,Вт, толщины 8и = у1/2 (где и(уу2) = ит/2) и стабилизацию среднего течения под влиянием устойчивой стратификации (рис. 4). Расширение струи происходит так же, как и для слоя смешения, за исключением стадии коллапса, отсутствующей в струе из-за более развитой мелкомасштабной структуры турбулентности за счет наличия свободной поверхности. Модели 2' (с /т,(у) = 0), 2 и 3 занижают асимптотические уровни максимумов
средних величин, не описывая выход 8и и на постоянные зна-
чения. Все модели 1-3 дают более медленный, как и для слоя смешения, закон затухания и'т(х), который, однако, близок к степенной зависимости, зафиксированной в однородной турбулентности за решеткой. Одно из соображений по расхождению расчета и опыта состоит в том, что при достижении течением стабилизированного состояния энергия турбулентности может расходоваться как на вязкую диссипацию, так и на излучение внутренних волн. Последний механизм не учитывается полуэмпирическими моделями 11А№-подхода.
Рис. 4. Эволюция интегральных характеристик плавучей струи при Ш0 к 0.007: □ - опыт;
Рис. 5. Схема течения: 1-3 - области рециркуляции, 4 - слой смешения, 5 — линия '// = 0, б - точка присоединения потока к стенке, 7 - восстановление пограничного слоя.
В § 2.5-2.7 рассматриваются вопросы моделирования течений в пограничном слое или канале с двумерным выступом квадратного сечения (рис. 5). Сформулированы четыре подхода: ¿-¿'-модель (1)-(5) с выражением Буссинеска для напряжений Рейнольдса (МБ); алгебраическая модель нормальных напряжений (АМН), где сдвиговое напряжение Рейнольдса находится из (3), а нормальные напряжения - из алгебраических соотношений, полученных из уравнений (6) с учетом (14); дифференциальная модель нормальных напряжений (ДМН), где сдвиговое напряжение определяется из (3), а нормальные напряже-
ния - из (6); полная дифференциальная модель второго порядка (ДМ), где все напряжения Рейнольдса находятся из уравнений (6). Последняя модель эквивалентна модели 2, примененной для течений со стратификацией. Вязкий подслой у стенки явно не разрешался, а краевые условия заданы в некоторой ее окрестности, исходя из степенных или логарифмических «законов стенки». Дискретизация уравнений выполнена на смещенной сетке для предотвращения рассогласования полей скорости и давления. Поле давления в эллиптических задачах для течений с препятствием (в отличие от рассматриваемых в § 2.1-2.4 параболических задач для тонких сдвиговых слоев) находится методом одновременных итераций давления и скорости [Viecelli J.A. // JCP. 1971. Vol. 8. P. 119-143], при помощи сформулированных в работе модификации для неравномерной сетки, критериев итерационной и аппроксимационной сходимости.
В § 2.6 проведена отладка алгоритма на основе модели МБ для течения с выступом в пограничном слое. Изучены некоторые эффекты алгоритма, влияющие на точность получаемых решений. В частности, установлены оптимальные положения границ области расчета, практически не влияющие на структуру течения с препятствием и позволяющие ставить краевые условия как условия на бесконечности. Разработана процедура растяжения сетки при удалении от выступа, экономящая ресурсы компьютера без потери точности счета. Помимо сходимости итераций расчета поля среднего давления и сходимости по времени, продемонстрирована аппроксимационная сходимость с измельчением сетки. На мелкой сетке с минимальным шагом Axmin = HI48 у препятствия выявлены малые области рециркуляции размерами в согласии с данными опытов.
В 2.6.2 в расчетах по модели МБ изучено влияние глубины погружения выступа квадратного сечения (высотой Н) в пограничный слой (толщиной 3) на структуру отрывного течения. Показано, что с уменьшением глубины погружения длина зоны рециркуляции в следе за препятствием растет в согласии с данными сопутствующих опытов (рис. 6). При этом применен оригинальный метод экстраполяции точек, соответствующих вычислениям с различным разрешением, при Лхт[п/Н->0, т.е. на случай бесконечно мелкой сетки.
В § 2.7 в расчетах течения вокруг препятствия на дне канала с большим коэффициентом загромождения выполнена верификация моделей напряжений Рейнольдса различного уровня описания. Модель МБ с линейным выражением (3) для нормальных напряжений дает сильные количественные расхождения и
даже отрицательные нефизические значения для и2 (рис. 7, а). Нелинейная по градиентам средней скорости, модель АМН исправляет этот недостаток и позволяет точнее описать характеристики течения с препятствием. Полная дифференциальная модель (ДМ) наиболее точно воспроизводит поле средней скорости перед выступом, над ним и в зоне восстановления за ним (рис. 7). ДМ также описывает эффект слияния зон рециркуляции над и за препятствием (рис. 8), наблюдаемый в опытах, что не удается получить по более простым моделям.
Третья глава посвящена моделированию поверхности раздела несмеши-вающихся текучих сред при решении уравнений Навье-Стокса и функции объемной фракции/ с континуальной моделью поверхностного натяжения (CSF):
Рис. 6. Зависимость длины зоны рециркуляции от относительной толщины пограничного слоя (слева): 1-4 - настоящие расчеты, 5-15 - данные опытов с различными параметрами; экстраполяция длин областей рециркуляции с измельчением сетки (справа).
а б
Рис. 7. Поперечные профили 1.5и2Д/а2„ (о), и/1/т (б) в потоке с выступом в канале: о - опыт; длинный, короткий штрих, пунктир и сплошные линии - модели МБ, АМН, ДМН и ДМ.
Рис. 8. Линии тока в канале с высотой, в два раза большей высоты препятствия: о - опыт, сплошные линии - модель ДМ с противопоточной схемой для членов конвекции в (1)—(6), пунктир - модель ДМ с ОШСК-схемой в (1) и схемой «против потока» в (4)-(6).
du дх,
du, ди.и.
1 Эр+I
д_ рдхк
p=pj+a(i-/).
CJK
р[г] 1J л
<Э? р дх)
H = t\f + Я2О-/).
дп,
du.
I—L
9/
dutf.
-A = l k = -
dt дхf
dfs
(17)
(18)
ЙГ;
= (19)
и, dx,
где и /¿i — плотность и динамическая вязкость более плотной среды (/= 1), а р2 и Ц2 - то же для менее плотной среды (f= 0), Fsi - объемная сила поверхностного натяжения [Brackbill J.U. et al. //JCP. 1992. Vol. 100. P. 335-354], a - коэффициент поверхностного натяжения, к - кривизна поверхности, й,. - единичный вектор нормали к поверхности, fs - сглаженная функция объемной фракции. Последнюю можно взять в первом приближении равной объемной фракции, найденной численно (т.е. fs = f в «упрощенной» CSF-модели), в связи с размыванием из-за численных схем. Корректное представление fs предполагает плавное изменение поперек поверхности конечной толщины за счет свертки характеристической функции / с интерполяционной функцией сглаживающего ядра К в «модифицированной» CSF-модели. Сформулированы выражения для fs (с функцией К в виде полинома восьмого порядка) применительно к плоским двумерным течениям с границами в виде стенок и плоскостей симметрии.
Для решения (17)-(19) взята та же, что и в § 2.5-2.7, процедура одновременных итераций поля скорости и давления, основанная на методе искусственной сжимаемости и дискретизации на смещенной сетке. В работе формулируется коррекция алгоритма расчета давления (с учетом переменной плотности) и критериев для шага по времени (с учетом поверхностного натяжения).
tit, = 1.2
т. = 2.4
ut. = г.4
xfW
Рис. 9. Изолинии объемной фракции с различными схемами адвекции: слева, сверху вниз -схема «против потока», QUICK, c4mm в (17) и «донор-акцептор» в (18); справа - c4mm.
В § 3.2 выполнена проверка аппроксимаций членов адвекции для задачи разрушения плотины (без учета поверхностного натяжения). Показано (рис. 9), что схема «против потока» дает нереалистично толстую поверхность раздела из-за численной диффузии, а QUICK-схема существенно искажает поверхность из-за ложных осцилляций. Процедура «донор-акцептор» [Hirt C.W., Nichols B.D. // JCP. 1981. Vol. 39. P. 335-354], обеспечивая резкую поверхность, демон-
стрирует ложные образования в виде «плавающих осколков». Наиболее приемлемый результат с весьма тонкой и гладкой поверхностью раздела обеспечила MUSCL-схема с QUICK-интерполянтами и TVD-ограничителями "compressive minmod" (c4mm) [Kelecy F.J., Pletcher R.H. // JCP. 1997. Vol. 138. P. 201-225].
С помощью нормализованных погрешностей основных характеристик течения исследована независимость решения от сетки и шага по времени. Показано, что реалистичные результаты достигаются как с более адекватными схемами, так и с измельчением сетки и модификацией пристеночных условий, особенно вблизи движущихся контактных линий. Интегральные характеристики, найденные в расчетах, хорошо согласуются с измеренными (рис. 10).
В § 3.3 исследована неустойчивость Рэлея-Тейлора (НРТ) на поверхности раздела (с наложенными малыми возмущениями) между двумя слоями тяжелой (сверху) и легкой сред. Показано (рис. 11, 12), что средняя амплитуда поверхности yf имеет начальный экспоненциальный рост, соответствующий стадии линейной устойчивости со скоростью роста « = const, а вязкость и поверхностное натяжение демпфируют развитие неустойчивости в согласии с данными теории и опытов и на этапе линейной устойчивости, и на нелинейном этапе.
Для близких по плотности сред (вода и бензол) поведение поверхности раздела служит аналогом развития НРТ в стратифицированной среде. Результаты даны на рис. 13 для условий опыта Льюиса [Lewis D.J. // Proc. R. Soc. Lond. A. 1950. Vol. 202. P. 81-96] в сравнении сданными линейной теории Тейлора
Г
■ ' ПО О
Рис. 10. Эволюция переднего фронта воды х/ (при у = 0) и высоты столба воды yf (х = 0) при обрушении плотины: 7и2-расчет, 3 и 4 - измерения при 11е = 4.3-104 и 12.2-104.
0.2
0.1
250
500 Re_ 750
Рис. 11. Эволюция НРТ без учета поверхностного натяжения при р\1р2 = 2. Слева - изолинии
объёмной фракции при Rem = gmXmlv = 400 и tltr = 3.2, 4.8, 6.4, 8.0, где tr = (Llg)m, X = 21. Справа - скорость п* — n(v/g2)] , где п = d(\r\ yj)/dt при -3 < \n(j>//L) < -2, в зависимости от числа Рейнольдса: 1 - теория, 2 - расчет (Kelecy & Pletcher 1997), 3 - настоящие вычисления.
Рис. 12. HPT с учетом поверхностного натяжения: эволюция In(y//L), где а ■ 10s = 0, 1, 2, 3, 3.5, 3.8, 4.0, 4.13, 4.2 - кривые 1-9 (слева); скорость п в зависимости от где щ = п{Ф= 0): теория (1), измерения (2), другие расчеты (3), настоящие вычисления по схеме c4mm с упрощенной (4) и модифицированной (5) версиями CSF-модели (справа).
поверхности раздела^ (справа): 1 - теория, 2 - измерения, 3 - невязкие расчеты Пуллина при Ф = 0.333 ays0/L = 0.08; 4-6 - настоящие расчеты при Rem= 1.03Т05, Ф = 0 и yJL = 0.01, 0.03,0.08; 7 - при Ф = 0.0422 и ysJL = 0.03 (для условий опыта), 8 - Ф = 0.333 и ys0/L= 0.08.
[Taylor G.I. // Proc. R. Soc Lond. A. 1950. Vol. 201. P. 192-196] и расчета Пуллина [Pullin D.I. // JFM. 1982. Vol. 119. P. 507-532]. Кривая 3 расчета Пуллина рано обрывается в отличие от настоящих результатов (кривые 4-8), показывающих согласие с данными Тейлора (7), Льюиса (2) и демпфирующее влияние вязкости и поверхностного натяжения. Для случая «вода - воздух», также изученного в опытах Льюиса, эволюция (рис. 14) отличается от предыдущего случая, однако, веер кривых 5-7, вычисленных для различных условий опыта, попадает в центр разброса данных измерений (2-3), а невязкий расчет Пуллина (4) рано обрывается из-за проблем численной реализации и занижает амплитуду поверхности.
Эволюция неустойчивости на начальном этапе не зависит от перепада плотности на поверхности раздела и согласуется с линейной теорией Тейлора. Если перепад плотности не слишком велик (р\/р2 < 5), на нелинейном этапе из-за усиления сдвига скорости между противонаправленными потоками сред наблюдаются эффекты неустойчивости Кельвина-Гельмгольца, приводящие к
появлению характерных грибовидных конвективных структур со спиралевидной закруткой потока (рис. 11, 13, 15). При большом перепаде плотности (как в случае «вода — воздух») тяжелая среда тонкими струями проникает глубоко вниз в легкую среду, приобретающую вид толстых колонн, поднимающихся вверх (рис. 14, 15). При pjpz < 2 картина развития неустойчивости сохраняется симметричной до значительных моментов, а с ростом р\/рг усиливается асимметрия НРТ. Поверхностное натяжение не только демпфирует развитие возмущений, но и препятствует искажению (фрагментации) поверхности раздела (рис. 14), способствуя утолщению кончика струи тяжелой среды и даже возврату к грибовидным структурам при p\!pi ~ 104 (рис. 15). Новые результаты моделирования НРТ подтверждают и уточняют данные предшествующих работ.
Рис. 14. НРТ в системе «вода - воздух» (слева направо): изолинии объёмной фракции if = 0.5)
при Ф = 0.04, t* = 0.5, 1.0, 1.5, 1.9; Ф = 0, t*= 1.5,1.9; амплитуда поверхности: / -теория, 2 и 5 — пределы разброса опытных данных, 4 - расчет Пуллина, 5-7 - настоящие вычисления.
О x/L О дЛ 0 x/L x/L 0 x/L 0 x/L x/L 0 x/L 0 x/L 0 x/L 1
Q
0 x/L 1 0 x/L 1
Рис. 15. Изолинии объёмной фракции воды (/= 0.5) на поверхности раздела воды и менее плотной среды (слева направо): при Ф = 0.0422, г* =1.9 и р\!рг = 1.01, 1.1,2,5, 10,20,50, 102, 103, 104; при Ф = 0.422, /* = 2.1,2.6 и />,/р2 = Ю4.
В четвертой главе представлены результаты DNS- и LES-моделирования внутренних волн с конвективным опрокидыванием, генерируемых двумерным препятствием высоты h в устойчиво стратифицированном потоке с постоянными значениями скорости U и градиента плотности (dp/dz)0 на входе. Определяющие «отфильтрованные» уравнения в приближении Буссинеска имеют вид
Зй „ Зй Вййк 32П Зтк 1 ар -т.
—L = 0, —'- + —'-¡- = v—f--¡¡----—-afgt, (20)
cxj dt дхк дхк дхк р ôxj
+ = J = p = gPnBL. (21)
dt дх, dxf dxt pJiN2
В DNS подсеточные напряжения (г,*) и потоки скаляра (qï) исчезают, а в LES выражаются моделями вихревой вязкости и диффузии типа Смагоринского.
Для сохранения эффективности алгоритма, сформулированного на равномерной декартовой сетке, топография задана неявно: в правую часть (20) введен дополнительный член объемной силы «сопротивления» степенного вида
Fl=-A(U,-,"/h)ïïl(ntnky"/2 при z<hs(x,y), F=0 при z>hs{x,y), где А - интенсивность силы, hs(x,y) - поверхность препятствия; m = 0 и А = 1000 выбраны после оптимизации, включавшей многочисленные расчеты.
Для устранения влияния отражения волн от входной и выходной границы применены поглощающие слои: в правые части (20)-(21) у границ вводятся дополнительные члены вида -о(х)(Ф - Фг) с переменной интенсивностью поглощения а(х). Эти члены вынуждают численное решение для Ф = (г7,,/) следовать предписанному решению Фг (значениям Ф на входе), эффективно абсорбируют отраженные от границ волны, препятствуя искажению картины при их возврате к препятствию, и экономят ресурсы компьютера, позволяя выбрать более короткую область расчета. В работе тестировались зоны поглощения с различными интенсивностью, толщиной, формой и положением и найдены оптимальные параметры поглощающих слоев косинусоидальной формы.
В § 4.1 изложены основные результаты многочисленных тестов для течений в каналах с пассивной/активной примесью, выступом и стратификацией, где верифицирован алгоритм и определены численные, геометрические и физические параметры, позволяющие получить эффекты обрушения волн наиболее отчетливо и с приемлемой точностью. Базовые расчеты реализованы при числе Фруда F;, = U/(Nh) = 0.6 (где N = [-ig/po)(dp/dz)0]m - частота Брента-Вяйсяля, Ро - характерная плотность), числе Прандтля-Шмидта Se = v/Â = 1 и 700, числе Рейнольдса Re = Uhlv = 4000 - большем, чем в предшествующих работах. Для интенсивного обрушения взят двумерный «холм» косинусоидальной формы длиной L = 3.56/г (на полувысоте), обеспечивающий «резонансную настройку» на теоретическую длину подветренных волн. Достаточное разрешение (шаг сетки Ах = 0.02Л < 3//, где Колмогоровский микромасштаб ц оценивался в точке максимальной диссипации) позволяет воспроизвести мелкомасштабные детали развития неустойчивости и возникающую в итоге турбулентную область.
В § 4.2 обсуждаются данные DNS обтекания холма при Sc = 1 (применительно к атмосфере). Развитие обрушения волн иллюстрируется траекториями
частиц в плоскости^ = 0 (рис. 16). После начального этапа происходит конвективное опрокидывание, приводящее к паре двумерных ламинарных вихрей, быстро разрушающихся при ? > 25 и переходящих в турбулентность. В глобальной картине течения (рис. 16, д) видна вторичная область обрушения вниз по потоку от первичной. Обе эти области генерируются одинаково и разделены примерно одной длиной подветренной волны. Цепочка приподнятых зон обрушения наблюдалась и в опытах, соответствуя «захваченной» подветренной волне, но вторая и последующие области обрушения были более нестабильными. Турбулентный ротор в следе холма (отрывной пузырь, обусловленный отрывом пограничного слоя, происходящим из-за большого положительного продольного градиента давления на поверхности холма) появляется задолго до приподнятой турбулентной области, возникающей при обрушении волн.
Рис. 16. Траектории частиц вокруг зоны обрушения при г = 17.5 (а), 22.5 (б), 27.5 (в), 40.0 (г); в большей области при г = 37.5 (д), где прямоугольник показывает границу графиков (а-г).
Рис. 17. Эволюция (безразмерных) осредненных по пространству корреляций: А-(и'2), □-(v'2), 0-(W2), <~(f'2), 0-e = ((Vu'f)/Re, V-sf = 2((V/')2)/(ReSc) .
Переход к турбулентности инициирован введением «белого шума» малой амплитуды к полю скаляра при t = 7.5. Возмущение постепенно затухает, а течение остается ламинарным до t ~ 20. Затем энергия извлекается из оставшегося шума и возбуждает трехмерную неустойчивость определенной длины волны в трансверсальном направлении (Ху = 0.5h при t ~ 25), приводящую к характерным для эффектов НРТ (изучаемых в § 3.3, 4.4) конвективным структурам. Характеристики турбулентных пульсаций, полученные осреднением данных DNS в покрывающей зону обрушения области (см. прямоугольник на рис. 16, г), показывают экспоненциальный рост на переходном этапе и установление квазистационарного состояния при 35 < t < 55 (рис. 17). На этом периоде отношение между диссипациями полей скорости и скаляра оказывается примерно тем же, что и отношение между самими пульсациями скорости и скаляра, что прямо свидетельствует о гипотезе R = const в RANS-моделях, обсуждаемой в 1.3.7.
Трансверсальные спектры, построенные в зависимости от волнового числа ку = h/Xy (см. ниже) также указывают, что при t ~ 35 стартует стационарный период. Спектры при t > 35 в различных местах турбулентной области содержат не менее декады с наклоном «—5/3», указывающим на наличие протяженного инерционного интервала и на отсутствие существенных эффектов плавучести на турбулентность в области перемешивания, возникающей при обрушении. Осредненные по у интенсивности турбулентности демонстрируют квазистационарное поведение при 35 < t < 55 и достигают, как и в опытах, примерно 20% от входной скорости U в некоторых точках турбулентной области, а напряжения Рейнольдса, осредненные по квазистационарному периоду, показывают подобие с наблюдаемыми в дальнем поле классического плоского следа (рис. 18).
Эволюция статей баланса ТКЭ, осредненных «глобально», как и величины на рис. 17, показывает (рис. 19) относительно малый вклад порождения плавучести во время квазистационарного периода. То есть существует приблизительный баланс между порождением ТКЭ сдвигом средней скорости и расходованием за счет диссипации и адвекции. Турбулентная диффузия (за счет тройных корреляций и пульсаций давления) очень мала и остаточный член баланса dk/dt также незначителен на квазистационарном интервале 35 < / < 55.
Рис. 19. Глобальный баланс ТКЭ, осредненный по турбулентной области (см. рис. 19, г): >-к„ V-(D,+Dv), <-Dp, □ -Р, A-G, О-(-*), 0-(-Ак).
Геофизически важные величины «эффективности смешения» Г = G/s и Гр = Ер/s имеют глобальные (осредненные по области обрушения и интервалу 35 < t < 55) значения 0.2 и 0.13 соответственно. Во многих океанологических приложениях также получено Г — 0.2. С другой стороны, величина /"меняется во времени на протяжении квазистационарного периода и локально в турбулентной области, а Гр со временем изменяется несущественно вокруг осреднен-ного значения Гр = 0.13. Это меньше, чем величина 0.2, найденная в некоторых океанологических приложениях, но скалярная эффективность смешения Г„ , по данным других работ, заметно меняется для различных случаев и параметров.
В рассматриваемой турбулентной области проведен анализ статистических моментов, полученных из данных DNS. В частности, выполнена априорная оценка общепринятых гипотез RANS-замыкания (которые применялись, например, в расчете стратифицированного слоя смешения в главе 2) для основных слагаемых в уравнениях напряжений Рейнольдса. В этих уравнениях, член Л7, обеспечивает перераспределение нормальных напряжений согласно аппроксимации стремления к изотропии (9)-(10), а тензор диссипации приближенно
Рис. 18. Осредненные по размаху холма профили (при х = 2.5) продольной интенсивности турбулентности в разные моменты квазистационарного периода (слева); осредненные при 35 < ?< 55 моменты: А-(м'2}, о-(^2}, 0-(и/2}, V-(«V), <-{р,2), 0-(и'р') (справа).
1
-0.01
следует локальной изотропии, 8у = (2/3)е0у . Поведение ей к похоже (рис. 20), т.е. приемлемыми аппроксимациями будут s = kl г* или е = kV2/L* (аналогичное выражение применено в § 2.2 вместо е-уравнения при расчете слоя смешения со стратификацией), где г* = 6.4 и L* = 1.45 - масштабы времени и длины крупных турбулентных вихрей при х ~ 2.2 и z ~ 3.0, где и к, и е имеют максимумы.
В § 4.3 обсуждаются результаты LES обтекания холма при Se = 700 (применительно к диффузии соленой воды в водоемах, напр., эстуариях рек). Расчет без подсеточных моделей показал: взятая из § 4.2 сетка 5120x512x512 недостаточна для разрешения турбулентных масштабов поля плотности, а скалярный диссипативный микромасштаб оказывается гораздо меньше шага сетки и кол-могоровского микромасштаба, поэтому некорректно воспроизводится диссипа-тивная подобласть спектров поля плотности, генерирующего численный шум. С другой стороны, привлечение подсеточных моделей типа Смагоринского со стандартным значением С„ = 0.1 и турбулентным числом Прандтля-Шмидта SCsgs 0.3 позволяет устранить ложный численный шум, обеспечивая корректные инерционный и диссипативный интервалы спектра. LES при Se = 700 с подсеточной моделью, действующей на малых масштабах, способен также воспроизвести детали перехода и квазистационарной области турбулентности. Вихревая вязкость достаточно мала, а энергосодержащие и инерционные интервалы в расчетах DNS при Se = 1 и LES при Se = 700 получаются идентичными, поэтому вклад подсеточных пульсаций в статистические моменты в первом приближении можно считать малым и результаты близки к аналогам из § 4.2, как показывает сравнение мгновенных и осредненных характеристик. Лишь на ранних стадиях перехода к турбулентности (когда подсеточные модели еще не включаются, а эффекты влияния различной молекулярной диффузии при Se = 1 и 700 существенны) проявляются различия, рассмотренные в 4.4.4.
В § 4.4 отдельно рассматриваются вопросы развития неустойчивости при опрокидывании внутренней волны (рис. 21). В момент t ~ 10 в поле плотности появляется неустойчивый градиент при х ~ 4 и z ~ 2, а из изолиний плотности в трансверсальной плоскости при х = 4 видны периодические структуры с длиной волны порядка Лу= 0.3-Ю.5/?, едва заметные на фоне «белого шума». Длину волны наиболее неустойчивого возмущения можно оценить как = 47rv2/3(ga)~1/3, где а = (р\ — р2)/(р 1 + Pi), согласно классической теории возмущений НРТ с различными длинами волн, налагаемых на плоскую поверхность раздела несжимаемых сред с плотностями рх и р2. Для рассматриваемого течения при 10 < t < 15 оценка дает Хт = 0.45А, что довольно близко к данным визуализации структуры с периодичностью по у в настоящей работе. При этом, амплитуды неустойчивости все еще очень малы, с перепадом скаляра 0.01-Ю.05% от его изменения по вертикали в области расчета. При t ~ 22 (см. рис. 21) наблюдается взаимное проникновение «языков» среды с разной плотностью, и формируется довольно резкая «поверхность раздела» с неустойчивым градиентом плотности, т.е. условия для появления неустойчивости Рэлея-Тейлора, рассматриваемой в главе 3. Действительно, при t > 21 происходит экспоненциальный рост вторичных возмущений, как видно из распределений плотности в трансверсальной плоскости y-z (рис. 22), показывающих формирование грибовидных конвективных эле-
ментов HPT при t > 24, и из эволюции «высоты перемешивания» S(t) (рис. 23). Определение последней величины предложено в настоящей работе по аналогии со средней амплитудой возмущений НРТ в § 3.3.
2 3x4 2 3*4
Рис. 20. Изолинии ТКЭ (слева) и ее диссипации (справа), осредненных по^ и при 35 < t < 55.
Рис. 21. Изолинии скаляра в сечении у = 0, показывающие развитие опрокидывания волны.
Рис. 22. Изолинии скаляра в DNS с Sc = 1 в сечении х = 2.5, показывающие развитие НРТ.
1
0.6
0.001
0.01
<5
0.1
0.5
0.4
0.1
0.3
0.2
0
лр
14 16 18 20
22 24 26
t
Рис. 23. Эволюция высоты <5 (о,«) и перепада безразмерной плотности Ар (Д, А) в сечении х = 2.5: о,Д - DNS при Se = 1; А - LES при Se = 700.
Из-за продолжающегося вторжения слоя тяжелой среды при t> 24 помимо закручивающегося вверх языка появляется (в точке торможения) и второй язык, закручивающийся вниз (рис. 24), следуя вращению вихрей, структура которых усложняется по сравнению с вихревой парой при t = 22.5 (рис. 16, б). Изолинии скаляра в плоскости x-z содержат отпечаток структур НРТ, развивающихся на поверхности раздела тяжелой и легкой сред (рис. 22). Изолинии скаляра при t = 25 свидетельствуют о наличии в области обрушения ряда слоев с более и менее плотной средой. Некоторые локальные максимумы плотности появляются из-за взаимного проникновения по вертикали двух элементов НРТ. Однако довольно резкие поверхности раздела при 1.8 < z < 2.1 (для х ~ 2.5) свидетельствуют именно о наличии двух слоев, простирающихся по размаху препятствия и являющихся следствием крупномасштабного опрокидывания волны, сопровождаемого вторжением тяжелого языка вверх по течению (рис. 21) и его расщеплением на два языка с разворотом на 180° (рис. 24). Поверхность раздела с резким положительным градиентом плотности изогнута, наклонена в плоскости (x,z) и содержит грибовидные элементы НРТ, видимые в различных сечениях (рис. 25). То есть НРТ на нелинейной стадии развития приводит к ряду противоположно вращающихся вихревых трубок («жгутов»), с завихренностью, направленной тангенциально поверхности вихревого цилиндра (рис. 26).
Природу структур (квазидвумерных на мелких масштабах порядка длины волны возмущений НРТ) иллюстрирует трехмерная визуализация инварианта Q = —ui jUj il2 (рис. 26). Вихревые структуры получаются сильно вытянутыми, а доминирование завихренности и характерные величины Q в них быстро растут со временем. Трансверсальная длина волны, соответствующая элементам НРТ при 22 < t < 25, согласно рис. 22, 25-26 составляет в среднем около 0.5h, в согласии с информацией из спектров (рис. 27). Из изолиний плотности и векторов скорости для отдельных элементов НРТ (рис. 28) видны зоны интенсивной завихренности (показанные окружностями) с закруткой в различных направлениях, расположенные на кромках шляпок грибов, как для НРТ на поверхности раздела несмешивающихся сред с небольшим перепадом плотности в § 3.3. Эти
зоны соответствуют вторичной неустойчивости Кельвина-Гельмгольца (НКГ). Такие же НКГ-вихри визуализируются на «шляпках грибов» и в других сечениях (рис. 25) вдоль цилиндрической «поверхности раздела», так что в итоге составляют вытянутые вихревые структуры, показанные на рис. 26. Рассмотренные особенности развития НРТ не наблюдались в предшествующих работах.
Рис. 24. Изолинии скаляра (слева) и линии тока (справа) при / = 25 в сечении у - 0. зГ
-5 -4 -3 -2 -1 0 у 1 2 3 4 5
Рис. 25. Изолинии плотности при t = 25 в трех сечениях (2.35 <х< 3.04).
Рис. 26. Инвариант градиента скорости при ? = 24.5, 2 = 5 (слева) и I = 25, (? = 20 (справа).
Рис. 27. Безразмерные трансверсальные спектры дисперсии пульсаций скаляра Еа = Е(//'),
осредненные при 1.25 <х < 5.00, 1.25 <г< 3.75. Вертикальные линии обозначают Лу/И = 0.3, 0.5, 2.5; числа у кривых - моменты времени; штриховая линия - наклон «-5/3».
Из трансверсальных спектров пульсаций плотности в различные моменты видно (рис. 27), что возмущение «белого шума» доминирует при t = 10 для большинства ку (= h/Xy), кроме наибольших, где влияет молекулярная диффузия. Позже, при 12 < t < 24 происходит экспоненциальное усиление первой моды, соответствующей ширине области расчета (Ау = 10К). При t > 21 появляется доминирующая мелкомасштабная мода с длиной волны Ху = 0.5h, обсуждаемая выше. При t ~ 25 эта мода генерирует наиболее мощный пик, соответствующий активизации периодических структур НРТ, но полностью исчезает при t = 27.5. К моменту t — 29 спектр приобретает классический вид, соответствующий турбулентности с инерционным интервалом «—5/3». При этом ТКЭ и ее компоненты при t ~ 29 достигают максимума, затем немного уменьшаются при релаксации к квазистационарному состоянию при / > 35 (рис. 17). Однако при t > 22 появляется пик, соответствующий крупномасштабной моде Ху = 2.5И, доминирующей при 30 < t < 45 (при осреднении спектра по всей области обрушения), особенно в момент t ~ 32.5. Эту моду можно полагать результатом преимущественного усиления произвольного малого поперечного возмущения исходного двумерного ламинарного движения (вихревой пары на рис. 16, б), возникающего при опрокидывании волн. Картины плотности и завихренности при t > 24 показывают также вероятность слияния грибовидных элементов НРТ с Ху~ 0.5h в более крупные образования. Мода 1У = 2.5h наблюдалась в опытах и расчетах других авторов, где идентифицированы связанные с ней крупномасштабные тороидальные структуры, т.е. сценарий укрупнения мелкомасштабных элементов НРТ не упоминался, да и сами структуры НРТ не выявлены (по-видимому, из-за недостаточного разрешения), хотя предполагалось, что исходная неустойчивость является конвективной и подобна НРТ-возмущениям.
Тороидальные вихри, обнаруженные на поздних этапах развития, когда область обрушения оказывается почти полностью турбулентной, и имеющие наклон около 40° к горизонтали в плоскости (х, z), выявляются и в настоящей работе (рис. 29). Хотя поток становится турбулентным, эти вихри остаются заметными и при t > 32.5, т.е. являются относительно долгоживущими. Область обрушения в плоскости (х, z) при t = 32.5 похожа на исходную вихревую пару, обусловленную крупномасштабным конвективным опрокидыванием и наблюдаемую на более раннем этапе (рис. 16, б), когда течение все еще ламинарное и практически двумерное. Во время последующего (t > 22) развития неустойчивости и перехода к турбулентности, внешняя глобальная система внутренних волн остается неизменной. Следовательно, эта система порождает путем усиления моды (Ху = 2.5h) глобальной трансверсальной неустойчивости исходной двумерной вихревой пары и одновременного процесса сложного взаимодействия и укрупнения мелкомасштабных структур НРТ и НКГ, фактически, единственные структуры в виде тороидальных вихрей диаметром в пределах 2.5/7.
3 У 3.5
Рис. 28. Изолинии скаляра и векторы скорости при ? = 25 в сечении х = 2.5.
Рис. 29. Траектории частиц и векторы скорости при I = 32.5 в сечениях^ = -0.2 и г = 2.26. Свидетельства наличия тороидальной структуры обозначены эллипсом и линией, соединяющей центры противоположно вращающихся вихрей.
В 4.4.4 обсуждаются результаты LES при Se = 700, где выявлены те же процессы развития неустойчивости, что и рассматриваемые в 4.4.1-4.4.3. Однако это развитие оказывается более быстрым и характеризуется меньшими длинами волн периодических структур НРТ (из-за более слабых эффектов диффузии и более резких поверхностей раздела между слоями разной плотности), чем в DNS при Se = 1. Тем не менее, можно полагать, что переход к турбулентности в атмосферных и океанических внутренних волнах, включающих одно и то же явление конвективного опрокидывания, происходит по близким сценариям.
В заключении представлены основные результаты диссертации:
1. С целью развития иерархии адекватных аппроксимаций RANS-метода проведена верификация стандартных и модифицированных моделей уравнений переноса для моментов второго порядка в канонических течениях с устойчивой стратификацией, свободной поверхностью, стенкой, отрывом и рециркуляцией, вводимой препятствием. Показано, в частности, что модифицированная модель второго порядка воспроизводит (в отличие от /г-е-модели с алгебраическими выражениями для корреляций второго порядка) противоградиентный характер турбулентных потоков тепла и импульса, наблюдаемый в устойчиво стратифицированном течении в канале со свободной поверхностью.
2. Исследована эволюция поверхности раздела двух несмешивающихся сред при помощи решения уравнений Навье-Стокса и уравнения для функции объемной фракции, а также континуальной модели поверхностного натяжения. Сформулированный метод адекватно воспроизводит эволюцию неустойчивости Рэлея-Тейлора (НРТ) на этапах и линейной, и нелинейной устойчивости, в том числе и для реальных случаев «вода - воздух», «вода - бензол». В последнем случае поведение поверхности раздела может служить аналогом развития НРТ в стратифицированной среде. Показано, что при небольшом перепаде плотности на нелинейном этапе наблюдаются эффекты вторичной неустойчивости Кельвина-Гельмгольца (НЮ") с грибовидными конвективными структурами, а при большом перепаде плотности более плотная среда тонкими струями проникает глубоко вниз между толстыми колоннами менее плотной среды, поднимающимися вверх. Поверхностное натяжение приводит к демпфированию развития неустойчивости и препятствует искажению (фрагментации) поверхности раздела. Новые результаты, полученные при моделировании НРТ, подтверждают и уточняют результаты работ других авторов.
3. Выполнено DNS/LES-моделирование обрушения внутренних волн, инициированных препятствием в устойчиво стратифицированном течении при числе Шмидта 1 < Se < 700 и числе Рейнольдса Re = 4000, намного большем, чем в предыдущих работах. Впервые показано, что в начальной стадии происходит рост мелкомасштабных возмущений НРТ, длина волны которых зависит от числа Шмидта и согласуется с теоретическим значением для наиболее неустойчивой моды на поверхности раздела слоев разной плотности. В отличие от классических двухслойных систем, сложная структура течения в области опрокидывания волн приводит к вытянутым и изогнутым квазидвумерным элементам НРТ, расположенным на поверхности цилиндра, опоясывающего эту область и содержащего резкие градиенты плотности. Эффекты НКГ приводят к
появлению ряда вихревых трубок, следующих поверхности упомянутого цилиндра. Сложное взаимодействие матрицы структур НРТ и НКГ в области обрушения приводит к генерации энергии во всех масштабах движения и возникновению турбулентности с протяженным инерционным интервалом «—5/3» на спектрах. На поздних стадиях развития доминирует крупномасштабная мода, представляющая наиболее неустойчивое возмущение исходной системы двумерной пары вихрей в области опрокидывания и соответствующая цепочке долгоживущих тороидальных вихрей. Подобные структуры обнаружены предшественниками, но механизмы их возникновения не изучены.
4. Показано, что возникающая при обрушении волн турбулентная область является квазистационарной, так что наблюдается приблизительный глобальный баланс порождения сдвигом, адвекции и диссипации энергии турбулентности. В рассматриваемой области выполнен анализ статистических моментов, полученных из данных DNS, и априорная оценка RANS-аппроксимаций в уравнениях переноса напряжений Рейнольдса. Этот анализ позволил сделать оценки геофизически важных величин: в частности, эффективность перемешивания получилась равной 0.2, как и в океанологических приложениях.
Полученные сценарии развития неустойчивости и характеристики квазистационарной перемешанной области при обрушении внутренних волн проливают свет на природу возникновения геофизической турбулентности в стратифицированной среде вдали от подстилающей поверхности.
Естественным продолжением исследований обрушения волн будет изучение поведения ротора в следе за препятствием и системы вторичных приподнятых зон обрушения, решение задач дальнейшей эволюции течения и перехода к реальным условиям в геофизических приложениях. Представляется важным создание карты режимов развития неустойчивости и турбулентности при обрушении волн в зависимости от чисел Рейнольдса, Прандтля-Шмидта, Фру-да. Увеличение Re и Se неизбежно приведет к введению моделей подсеточных напряжений и потоков скаляра (и связанных с этим погрешностей). Решение этих задач является весьма трудоемким и потребует гораздо больших компьютерных ресурсов, чем использованные в настоящей работе.
Автор выражает благодарность проф. А.Ф. Курбацкому за поддержку и внимание к вопросам RANS-моделирования, проф. К.С. Чану за мотивацию к развитию способов разрешения поверхности раздела сред, проф. И.П. Кастро и д-ру Т.Г. Томасу за помощь в получении данных DNS/LES и их обсуждение.
ОСНОВНЫЕ ПУБЛИКАЦИИ ПО ТЕМЕ ДИССЕРТАЦИИ
1. Курбацкий А.Ф., Яковенко С.Н. К моделированию турбулентного слоя смешения в устойчиво стратифицированной жидкости У/ Известия СО АН СССР, серия технических наук. 1989. Вып. 4. С. 59-65.
2. Курбацкий А.Ф., Яковенко С.Н. Моделирование слоя смешения переменной плотности (Ч. 1 : без учета перемежаемости) // Сиб. физ.-техн. журнал (Известия СО АН СССР). 1991. Вып. 4. С. 69-76.
3. Курбацкий А.Ф., Яковенко С.Н. Моделирование слоя смешения переменной плотности (4.2: с учетом перемежаемости) // Сиб. физ.-техн. журнал (Известия СО АН СССР). 1991. Вып. 4. С. 77-84.
4. Курбацкий А.Ф., Яковенко С.Н. Моделирование турбулентного переноса импульса и тепла в устойчиво стратифицированном течении со свободной поверхностью // Сиб. физ.-техн. журнал (Известия СО РАН). 1992. Вып. 4. С. 55-63.
5. Курбацкий А.Ф., Яковенко С.Н. Моделирование турбулентной струи, распространяющейся по поверхности более плотной жидкости // Сиб. физ.-техн. журнал (Известия СО РАН). 1993. Вып. 1. С. 50-62.
6. Яковенко С.Н. Верификация моделей второго порядка в развитом турбулентном течении в открытом канале // Теплофизика и аэромеханика. 1994. Т. 1, № 2. С. 151-158.
7. Курбацкий А.Ф., Поросева C.B., Яковенко С.Н. Вычисление статистических характеристик турбулентного потока во вращающейся круглой трубе // Теплофизика высоких температур. 1995. Т. 33, № 5. С. 738-748.
8. Курбацкий А.Ф., Поросева C.B., Яковенко С.Н. О моделировании поведения одноточечных моментов поля скорости четвертого порядка в развитом турбулентном потоке в круглой трубе // Теплофизика и аэромеханика. 1996. Т. 3, № 1.С. 41-52.
9. Курбацкий А.Ф., Яковенко С.Н. Численное исследование турбулентного течения вокруг двумерного препятствия в пограничном слое // Теплофизика и аэромеханика. 1996. Т. 3, № 2. С. 145-163.
Ю.Заец П.Г., Курбацкий А.Ф., Онуфриев А.Т., Поросева C.B., Сафаров H.A., Сафаров P.A., Яковенко С.Н. Экспериментальное изучение и математическое моделирование статистических характеристик турбулентного потока во вращающейся относительно продольной оси прямой круглой трубе // Прикладная механика и техническая физика. 1998. Т. 39, № 2. С. 101-114.
П.Зиновьев А.Т., Яковенко С.Н. Моделирование вертикального турбулентного обмена в пристенном стратифицированном водоеме // Прикладная механика и техническая физика. 1998. Т. 39, № 6. С. 57-64.
12. Курбацкий А.Ф., Яковенко С.Н. Моделирование структуры турбулентного потока вокруг препятствия с острыми кромками в плоском канале. Модели турбулентности // Теплофизика высоких температур. 1998. Т. 36, № 6. С. 927-932.
13.Яковенко С.Н., Курбацкий А.Ф. Моделирование структуры турбулентного потока вокруг препятствия с острыми кромками в плоском канале. Результаты численного моделирования // Теплофизика высоких температур. 1999. Т. 37, № 1.С. 98-105.
14. Курбацкий А.Ф., Яковенко С.Н. Диффузия пассивной примеси от линейного источника в нейтральном атмосферном приземном слое // Физика атмосферы и океана. 1999. Т. 35, № 4. С. 506-515.
15. Kurbatskii A.F., Yakovenko S.N. Diffusion of passive contaminant from a line source in a neutrally stratified turbulent boundary layer // Wind and Structures. 2000. Vol. 3, No. l.P. 11-21.
16. Kurbatskii A.F., Yakovenko S.N. Turbulence closure schemes suitable for air pollution and wind engineering // Wind Engineering and Industrial Aerodynamics. 2000. Vol. 87. P. 231-241.
17. Kurbatskii A.F., Yakovenko S.N. Computational study of complex turbulent flows with sharp-edged boundary and passive scalar sources // Computational Fluid Dynamics Journal. 2001. Vol. 9, Special Number. P. 459-469.
18. Ilyushin B.B., Yakovenko S.N. Testing of models for fourth-order cumulant, third- and second-moment models in stationary and rotating pipe flows // Journal of Engineering Thermophysics. 2002. Vol. 11, No. 1. P. 45-71.
19.Ларичкин В.В., Яковенко С.Н. Влияние толщины пограничного слоя на структуру пристенного течения с двумерным выступом // Прикладная механика и техническая физика. 2003. Т. 44, № 3. С. 76-84.
20. Yakovenko S.N., Chang К.С. Performance examination of geometry-independent near-wall second-moment closures in simple and backstep flows // Numerical Heat Transfer, Part B: Fundamentals. 2007. Vol. 51, Issue 2. P. 179-204.
21. Яковенко C.H., Чан К.С. Аппроксимация потока объемной фракции в течении двух жидкостей // Теплофизика и аэромеханика. 2008. Т. 15, № 2. С. 181-199.
22. Яковенко С.Н., Чан К.С. Применение неразрывной модели для силы поверхностного натяжения к задаче неустойчивости Рэлея-Тейлора // Теплофизика и аэромеханика. 2011. Т. 18, № 3. С. 449-461.
23. Yakovenko S.N., Thomas T.G., Castro I.P. A turbulent patch arising from a breaking internal wave // J. Fluid Mech. 2011. Vol. 677. P. 103-133.
24. Яковенко С.Н. Исследование неустойчивости Рэлея-Тейлора в задачах механики многофазных и стратифицированных сред // Вестник Нижегородского университета им. Н.И. Лобачевского. 2011. № 4, часть 3. С. 1280-1282.
25. Яковенко С.Н. Бюджет уравнений для напряжений Рейнольдса в области турбулентности, возникающей при обрушении внутренних волн // Вестник Новосиб. гос. ун-та. Серия: Физика. 2012. Т. 7, вып. 4. С. 87-95.
26. Yakovenko S.N. Rayleigh-Taylor instability in two-fluid and stratified media // Computational Thermal Sciences. 2012. Vol. 4, Issue 5. P. 399-409.
27. Яковенко С.Н. Моделирование неустойчивости Рэлея-Тейлора в задачах механики многофазных и стратифицированных сред // Известия вузов. Физика. 2013. Т. 56, № 6/3. С. 179-181.
28. Яковенко С.Н. Влияние перепада плотности и поверхностного натяжения на поверхности раздела текучих сред на развитие неустойчивости Рэлея-Тейлора // Известия РАН. Механика жидкости и газа. 2014. № 6. С. 54-69.
29. Yakovenko S.N., Thomas T.G., Castro I.P. Transition through Rayleigh-Taylor instabilities in a breaking internal lee wave // J. Fluid Mech. 2014. Vol. 760. P. 466-493.
Ответственный за выпуск С.Н. Яковенко
Подписано в печать 18.12.2014 Формат бумаги 60x84/16, Усл. печ. л. 2.0 Уч.-изд. л. 2.0, Тираж 150 экз., Заказ № 12
Отпечатано в типографии ООО «Параллель» 630090, Новосибирск, Институтская, 4/1
