Численное моделирование сложных режимов конвекции Рэлея-Бенара тема автореферата и диссертации по механике, 01.02.05 ВАК РФ

005009836

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

Палымский Игорь Борисович

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ СЛОЖНЫХ РЕЖИМОВ КОНВЕКЦИИ РЭЛЕЯ-БЕНАРА

01.02.05 - Механика жидкости, газа и плазмы

Автореферат

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

І С 0ЕЗ «

Пермь - 2012

005009836

Работа выполнена на кафедре общетехнических дисциплин Новосибирского филиала Военного учебно-научного центра Сухопутных войск общевойсковой военной академии ВС РФ (ВУНЦ СВ ОВА ВС РФ, г. Новосибирск)

Официальные оппоненты:

доктор физико-математических наук, профессор Т арунин Евгений Леонидович

доктор технических наук, профессор Цаплин Алексей Иванович

доктор физико-математических наук, профессор Бердников Владимир Степанович

Ведущая организация - Институт вычислительных технологий СО РАН, г. Новосибирск

Защита состоится 20 марта 2012 г. в 1515 на заседании диссертационного совета Д212.189.06 в Пермском государственном национальном исследовательском университете по адресу: 614990, г. Пермь, ул. Букирева, 15, зал заседаний Ученого совета ПГНИУ.

С диссертацией можно ознакомиться в библиотеке Пермского государственного национального исследовательского университета.

Автореферат разослан “¿Г"" января 2012 г.

Ученый секретарь

диссертационного совета Д212.189.06 кандидат физико-математических наук, доцент

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

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

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

В классической задаче о конвекции Рэлея-Бенара ситуация существенно упрощается, так как отсутствие среднего течения обуславливает сравнительно невысокие числа Рейнольдса.

Целью настоящего исследования является изучение различных режимов двумерной и трехмерной конвекции Рэлея-Бенара, от стационарных до сложных, стохастических (переходных). При этом конвекция вязкой и несжимаемой жидкости, возникающая при подогреве снизу, рассматривается в областях простой геометрии и большом диапазоне надкритичности г = Ка/Яасг, где Ка и Яаа - число Рэлея и его критическое значение. Трехмерные расчеты проводятся со свободными от касательных напряжений горизонтальными границами, а двумерные - свободными или жесткими.

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

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

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

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

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

динамики спектров температуры и скорости.

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

В работе вперпыс:

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

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

• Показано, что, несмотря на то, что численный метод для расчета двумерной конвекции строится как метод первого порядка аппроксимации по времени и второго по пространству, при его линеаризации порядок аппроксимации по времени повышается до второго.

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

• Предложено новое представление решения в горизонтальном направ-

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

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

• Показано, что, несмотря на наблюдаемые количественные расхождения между результатами расчетов и экспериментом, трехмерные расчеты со свободными горизонтальными границами дают правильные показатели степенных законов изменения среднеквадратичных пульсаций температуры, вертикальной скорости, числа Рейнольдса, среднеквадратичной скорости и кинетической энергии от надкритичности при / > 150. В двумерных расчетах аналогичное соответствии наблюдается при сравнительно невысокой надкритичности г <250.

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

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

При умеренной горизонтальной протяженности области /х, /у = я:

6

• Установлено, что в двумерных и трехмерных расчетах со свободными границами действие силы плавучести обуславливает стратификационные спектры для скорости, а доминирование конвективного переноса для температуры - спектр пассивной примеси к~5/3. Эти спектры наблюдаются в двумерных и трехмерных расчетах практически во всем рассмотренном диапазоне надкритичности.

При увеличенной до /х = 471 горизонтальной протяженности области в двумерных расчетах со свободными границами:

• При 500<г<4-103 длинноволновый участок спектра скорости состоит из двух ветвей, где идентифицируются степенные законы, отражающие конкурирующие механизмы, связанные с силой плавучести и каскадным процессом переноса энергии, а при 4-103 < г < 104 в спектре скорости видны степенные законы к'3 и к'5/3, соответствующие прямому и обратному каскадным процессам переноса энстрофии и энергии, соответственно. При г> 10* в спектре скорости видна тенденция к установлению единого степенного закона к'2'6 - искаженного закона каскада энстрофии.

• При умеренной надкритичности 500 < г < 2-103 в длинноволновом участке спектра температуры виден стратификационный спектр Болджиано-Обухова к7 5, который при 2-103 < г <6-103 сменяется спектром Бэтчелора к'1. А при г > 6-103 показатель длинноволнового участка спектра температуры совершает колебания в интервале [-1,-0.7], с наиболее вероятным значением -0.8.

• Установлено, что отмеченные выше перестройки спектров температуры и скорости в двумерных расчетах со свободными границами связаны с формированием обратного каскада энергии при г ~ 2-103 (для скорости при г~ 4-103) и уменьшением его энергетической роли при г~ 6-103 (г-104).

Автор защищает:

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

• Вывод о том, что в случае двумерной и трехмерной конвекции со свободными горизонтальными границами при Pr= 1, коэффициенты при шаге

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

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

• Результаты нелинейного анализа, проведенного на модельной нелинейной системе уравнений, которые показывают, что понижение порядка аппроксимации численного метода расчета 2d,free конвекции по времени до первого приводит лишь к незначительному понижению точности вычислений, при этом с первым порядком вычисляется фаза гармонического решения, в то время как его амплитуда - со вторым. И как следствие этого, для практических 2d,free расчетов можно использовать схему первого порядка аппроксимации по времени с вычислением скоростей по функции тока, полученной на первом этапе расщепления.

• Результаты двумерных расчетов с жесткими горизонтальными границами стационарной, валиковой конвекции, которые показывают при невысокой надкритичности 1 < г < 4 удовлетворительное согласование зависимости доминирующей длины волны от надкритичности и хорошее - по числу Нуссельта, профилю средней температуры и изотерме полной температуры.

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

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

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

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

• Вывод о том, что, несмотря на наблюдаемые количественные расхождения между результатами расчетов и экспериментом, трехмерные расчеты со свободными горизонтальными границами дают правильные показатели степенных законов изменения среднеквадратичных пульсаций температуры, вертикальной скорости, числа Рейнольдса, среднеквадратичной скорости и кинетической энергии от надкритичности при г > 150. А в двумерных расчетах аналогичное соответствии наблюдается при сравнительно невысокой надкритичности г <250.

• Вывод о том, что закон роста числа Нуссельта, поведение пульсаций

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

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

При умеренной горизонтальной протяженности области /х, !у = п в двумерных и трехмерных расчетах со свободными границами:

• Результаты исследования пространственных спектров для скорости и температуры, где установлено, что действие силы плавучести обуславливает стратификационные спектры для скорости, а доминирование конвективного переноса для температуры - спектр пассивной примеси к~5/3. Причем указанные спектры наблюдаются в двумерных и трехмерных расчетах практически во всем исследованном диапазоне надкритичности.

11ри увеличенной до /х = 4я горизонтальной протяженности области в

двумерных расчетах со свободными горизонтальными границами:

9

• Результаты исследования спектра скорости, где показано, что при 500 < г < 4-103 длинноволновый участок спектра скорости состоит из двух ветвей, со степенными законами, отражающими конкурирующие механизмы, связанные с силой плавучести и каскадным процессом переноса энергии, а при 4-103</-< 104 в спектре скорости идентифицируются степенные законы к'3 и к'5/3, соответствующие прямому и обратному каскадным процессам переноса энстрофии и энергии, соответственно. А при г> 104 в спектре скорости видна тенденция к установлению единого степенного закона к'2'6 - искаженного закона каскада энстрофии.

• Результаты исследования спектра температуры, показывающие, что при умеренной надкритичности 500 < >■ < 2-103 в спектре температуры виден стратификационный спектр Болджиано-Обухова к'?/5, который при 2-103 < г < 6-103 сменяется спектром Бэтчелора к'1. А при г > 6-103 показатель длинноволнового участка температурного спектра совершает колебания в интервале [-1,-0.7], с наиболее вероятным значением - 0.8.

• Вывод о том, что перестройки спектров температуры (скорости) в двумерных расчетах со свободными границами связаны с формированием обратного каскада энергии при г ~ 2-103 (для скорости при г а 4-103) и уменьшением его энергетической роли при г ~ 6-103 (г =104).

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

Научная и практическая значимость результатов диссертационной работы заключается в том, что в ней разработан и реализован в виде программы метод расчета конвекции Рэлея-Бенара, в том числе и сложных ее режимов при высокой надкритичности. Реализован также и параллельный вариант компьютерной программы для многопроцессорных компьютеров с

общей памятью.

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

С точки зрения возможности сравнения с данными по конвекции с более реальными с экспериментальной точки зрения жесткими горизонтальными границами представляется важным то, что, несмотря на различие в граничных условиях в эксперименте с жесткими границами и трехмерных расчетах со свободными и обусловленные этим количественные расхождения, трехмерные расчеты дают правильные показатели степенных законов изменения среднеквадратичных пульсаций температуры, вертикальной скорости, числа Рейнольдса, среднеквадратичной скорости и кинетической энергии от надкритичности при г> 150. А в двумерных расчетах аналогичное соответствии наблюдается при сравнительно невысокой надкритичности г < 250.

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

Апробация работы. Результаты исследований были представлены на Четвертом сибирском конгрессе по прикладной и индустриальной математике (ИНПРИМ - 2000), Пятой международной конференции “Лаврентьевские чтения по математике, механике и физике” (Новосибирск - 2000), International Conference “Advanced Problems in Thermal Convection” (Perm - 2003),

11

Международной конференции “Нелинейные задачи теории гидродинамической устойчивости и турбулентности” (НеЗаТеГиУс - 2006, 2008), Fourth International Conference on Computational Heat and Mass Transfer (ICCHMT’05, Paris, France - 2005), Seventh International Symposium on Hazards, Prevention and Mitigation of Industrial Explosions (St. Petersburg, Russia - 2008), Всероссийской конференции “Математика в приложениях” (Новосибирск - 2009), International Conference “Dynamics Days of Europe” (DDE09, Gettingen, Germany - 2009), 8* International Symposium on Engineering Turbulence Modelling and Measurements (ETMM8, Marseille, France - 2010), Международной научной конференции “Забабахинские научные чтения” (Снежинск - 2010), 7* International Conference on Computational Heat and Mass Transfer (ICCHMT’ll, Istanbul, Turkey - 2011), Пермском гидродинамическом семинаре им. Г.З. Гершуни и Е.М. Жуховицкого (Пермь 2008, 2010), а также неоднократно на научных семинарах в институтах Теоретической и прикладной механики, Гидродинамики, Вычислительной математики и геофизики, Вычислительных технологий и Вычислительного моделирования СО РАН, Проблем механики РАН и Механики МГУ.

Публикации. Основные результаты диссертации опубликованы в работах [1-21]. Из них статьи [1-14] (в научных журналах, входящих в Перечень ведущих рецензируемых научных журналов и изданий, в которых должны быть опубликованы основные научные результаты диссертаций на соискание ученой степени доктора и кандидата наук), монография [15] и труды международных конференций [16-21].

Личный вклад автора. В работах [5,7,18] автору принадлежит разработка и программная реализация численного метода, проведение расчетов и анализ полученных данных, а в работе [8] автором выполнены расчеты конвекции Рэлея-Бенара и анализ их результатов, при этом расчеты двухдиффузионной конвекции были выполнены соавторами.

Структура и объем работы. Диссертация состоит из введения, трех глав, заключения, двух приложений и списка цитированной литературы, включающего 110 наименований. Общий объем диссертации 208 страниц, включая 111 рисунков.

СОДЕРЖАНИЕ ДИССЕРТАЦИИ

Во введении показана актуальность проблемы, дана общая характеристика работы и сформулирована цель исследования.

В первой главе “Постановка задачи о конвекции Рэлея-Бенара и численные алгоритмы“ описана постановка задачи о конвекции и спектральноразностный метод ее решения.

В п. 1.1 приведены постановки трех рассматриваемых в работе задач о конвекции - двумерной со свободными от касательных напряжений (2d,free), двумерной с жесткими (2d,rigid) и трехмерной со свободными (3d,free) горизонтальными границами. Исходная система уравнений записана в отклонениях от равновесного решения, который характеризуется равными нулю скоростями и линейным профилем температуры. Число Прандгля Рг, за исключением нескольких тестовых расчетов, везде равно 10.

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

В 3d,free расчетах рассматриваются решения вида:

К N М

u(t,x,y,z) = Y,Y,Hui™(t')PtP„P„™s(akx)cos(/}ny)cos(7rmz),

0 й=0 т=0 К-1 N-1 М

v(t,x,y,z) = ZZ Z v*™(Оря sin(akx)sin(0ny)cos(srmz),

Jfc=l ч=1 m-d К-1 N M~I

w (t,x,y,z) = Z Z Z wto,(OP„ sin(akx)cos(finy)sin(?rmz),

k=l n=Q m=l A'-1 N M

p(t,x,y,z) = Z Z Z Ptm(0PKP* sin(akx)cos{finy}cos(^mz),

¿=1 n=0 m=0 K-\ N M-1

Q(t,x,y,z) = Z Z Z QЬт(0ря sin(akx)cos(pny)s\n{7tmz),

л=о jn—l

где и0до=0, a Pk = 0.5 при k = 0, иК = ! при l<k<K-l.

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

G = {0 <х < я/а, 0 <у < Tt/р, 0 <z< 1},

где а и Р - минимальные волновые числа, а /х = я/а и /у = тьф - относительные горизонтальные протяженности области в х и у направлениях, соответственно.

В случае свободных горизонтальных границ все граничные условия следуют из вида решения, например, при х = О, я/а и 0 <у< тс/(3, 0 < г < 1 получаем: их = V = = <3 = 0, что соответствует проницаемой и идеально

теплопроводной стенке. А при у = 0, я/(3 и 0 <х < я/а, 0 < г < 1 находим: иу = V = луу = Ру = 0 - условия на вертикальной границе конвективной ячейки. Некоторая искусственность такой постановки граничных условий обусловлена желанием обеспечить преемственность с двумерной постановкой. А на горизонтальных границах при г = 0,1 имеем: иг = уг = \у = = 0.

В двумерных постановках задачи о конвекции горизонтальное направление у отсутствует, а в случае конвекции приведенное выше гра-

ничное условие для скорости при 2 - 0,1 заменяется на условие прилипания.

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

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

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

На втором этапе расщепления в физическом пространстве учитываются все нелинейные члены или, другими словами, учитывается нелинейное взаи-

14

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

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

В случае жестких границ, метод решения строится, в принципе, аналогично. 1 лавпое отличие заключается в необходимости рассчитывать значение вихря на горизонтальных границах. Используемая здесь формула, основанная на применении алгоритма БПФ, на тестовом примере показала сходимость с первым порядком. Но, по сравнению с известной формулой Тома, имеющей также первый порядок, точность предлагаемой формулы на этом тесте в четыре раза выше. Кроме того, как показывают приведенные здесь результаты методических расчетов, решение задачи о конвекции относительно слабо зависит от числа учитываемых членов в формуле для вычисления вихря.

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

В п. 1.3 описан трехмерный вариант предлагаемого спектральноразностного метода. Исходная система уравнений записана в переменных скорость, давление, температура и переход на следующий слой по времени состоит из трех основных шагов.

Алгоритм расчета на первых двух шагах расщепления, в принципе, аналогичен уже описанному в двумерном случае, но на втором этапе расщепления здесь применена явная схема с направленными разностями первого порядка и поправкой A.A. Самарского, повышающей порядок аппроксимации до второго. На третьем этапе расщепления восстанавливается выполнение уравнения неразрывности, нарушенное на втором этапе расщепления.

Предложенный спектрально-разностный метод имеет первый порядок по времени и, кроме того, в трехмерном его варианте на втором этапе расщепления использована явная схема интегрирования по времени.

15

Но указанные обстоятельства практически не понижают эффективность метода расчета, так как условие корректного моделирования при наличии неустойчивости течения в линейном приближении требует проведения расчетов с достаточно малым шагом по времени. И, кроме того, проведенный во второй главе на модельной системе уравнений нелинейный анализ 2с1,Ггее алгоритма показал, что повышение порядка аппроксимации численного метода до второго не приводит к сколько-нибудь заметному увеличению точности вычислений.

Вторая глава ’’Анализ и тестирование численных алгоритмов” посвящена линейному и нелинейному анализу предложенного спектральноразностного численного метода и описанию результатов тестовых расчетов.

В п.2.1 изложены результаты линейного анализа двумерных алгоритмов.

Для конвекции со свободными горизонтальными границами при Рг = 1 выписаны явные выражения для инкрементов нарастания гармоник в линейных аналогах дифференциальной задачи и предлагаемого численного метода в двумерном и трехмерном случаях. Показано, что коэффициенты при т и Н2 в членах, описывающих схемный эффект, не зависят от числа Рэлея, а определяются только волновыми числами.

А в случае 2с1,п^с1 конвекции сравнивается спектральная кривая, полученная методом ортогонализации из линеаризованной исходной системы уравнений и приближенная, полученная численно из линейного аналога предлагаемого численного метода.

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

Например, для метода расчета 2<1,&ее конвекции формула для инкремента Я при Рг = 1 принимает вид:

К =4 +—(<х6к6 +л6т6)-^-а4к

ту 2 гг:

п1 л4г4 2

96 24 24'

здесь Л5Г и - инкременты нарастания гармоники для предлагаемого численного метода и дифференциальной задачи, а и л - минимальные волновые числа в горизонтальном х и вертикальном направлениях, к и т -номера гармоник, г и Я - шаги по времени и пространству, а X < 0 соответствует нарастанию. Из приведенной формулы видно, что схемный эффект

(все члены правой части, кроме первого) есть величина порядка 0(т? + Н2) и что коэффициенты при т и II2 определяются только волновыми числами.

В п.2.2 изложены результаты нелинейного анализа (на модельной системе уравнений) 2сЦгес алгоритма. Поскольку главной целью нелинейного анализа бьшо исследование аппроксимаций нелинейных членов по времени, то анализ ограничился случаем длинных волн, с дискретизацией только по времени. Изменение вида нелинейных членов в исходной системе уравнений (нелинейные члены вида и их заменены на \и[их, что обуславливает возможное различие в знаке) сделало возможным рассмотрение частных, комплексных решений в виде монохроматической волны в случае дифференциальной системы уравнений и численного метода, с последующим сравнением полученных в виде степенных рядов по г выражений.

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

В п.2.3 изложены результаты линейного анализа трехмерного алгоритма.

Показано, что, как и в двумерном случае, коэффициенты при г и Я2 в выражении для инкремента нарастания гармоники X не зависят от числа Рэлея, а определяются только волновыми числами.

На рис. 1 для трехмерной задачи с учетом 653 гармоник приведен инкремент нарастания Я (первые три старшие моды) при г = Ка/Яасг - 950 и Рг= 1. Здесь ак - волновое число в горизонтальном направлении х (волновые числа в обоих горизонтальных направлениях на рис. 1 положены равными), а область \<0 соответствует нарастанию.

Из рис. 1 видно, что спектральные характеристики численного метода и исходной дифференциальной задачи близки с заметным количественным отличием в области волновых чисел, соответствующей затухающим (А > 0) гармоникам.

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

В п.2.4.1 приведены результаты тестирования метода расчета 2с1,1гее конвекции.

Рис. 1. Инкремент X в дифферещи- Результаты расчета стационар-алъной задаче (кривая 1) и численном ных течений сравниваются с резуль-методе (2) татами других авторов при

надкритичности г < 103 и по более простой модели бесконечного числа Прандтля - при г < 2-103. Отмечено хорошее совпадение по числу Нуссельта.

Проведением при г = 6-103 на разных ПК двух расчетов, с полностью совпадающими компьютерными кодами, трансляторами, числом учитываемых гармоник, начальными и граничными условиями показано, что при сопоставлении параметров сложных конвективных течений имеет смысл, в основном, сопоставление средних (интегральных) характеристик, которые вычисляются достаточно устойчиво. А анализ мгновенных величин можно проводить только при достаточно высокой точности вычислений и на малом отрезке времени.

Устойчивость вычисления интегральных характеристик подтверждает также серия методических расчетов при том же значении надкритичности и различном числе учитываемых гармоник - [129x33], [257x65], [513x129], [1025x257] и [2049x513]. Анализ результатов этой серии расчетов показал, что наблюдается сходимость:

• всех анализируемых интегральных величин,

• профилей средней температуры, среднеквадратичных пульсаций температуры, вертикальной и горизонтальной скорости,

• одномерных пространственных спектров температуры и скорости.

И дополнительно, сравнение временных спектров числа Нуссельта, полученных с учетом [513x129] и [1025x257] гармоник, показало их качественное подобие с совпадением положений максимумов и законов затухания на высоких частотах.

Два 2с1,(гее расчета, проведенные с учетом [129x33] гармоник при г = 6* 103 с двукратной и четырехкратной точностью вычислений также показали практически полное совпадение интегральных характеристик.

Интегральные характеристики, полученные другим методом - модифицированным методом Галеркина (псевдоспектральным), основанным на вычислении нелинейных членов в физическом пространстве на разностной сетке, имеющей в каждом направлении в два раза больше точек, чем гармоник в том же направлении и схеме Рунге-Кутты четвертого порядка для интегрирования по времени, хорошо согласуются с результатами 2сЦгее расчетов стохастических течений двухдиффузионной конвекции и конвекции при надкритичности до г = 2-103.

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

Методические расчеты показали, что при а ~ 1 в 2сЦгее расчетах достаточно учитывать [129x33] гармоник при г < 6-103 и [257x65] - при 6-103 < г < 3.4-104. При увеличении горизонтальной протяженности области 4 = к/а (уменьшении а) количество гармоник в направлении х увеличивается обратно пропорционально а, таким образом, при а = 0.125 расчеты проводились с учетом до [2049x65] гармоник.

Отмечено, что при численном моделировании сложных режимов конвекции стало стандартным требование о примерном равенстве между масштабом диссипации X и пространственным разрешением численного метода Нге (полусумма пространственных шагов в горизонтальном и вертикальном направлениях):

р И

Л =----------гг- = 0.531 гош

" (Ка(Ми-1))ш

в 2сЦгсе расчетах при а = 1. Показано, что Нге < X при всех значениях надкритичности и Яге меньше X на 30% при г = 3.4- Ю4 и числе гармоник [257x65].

В п.2.4.2 приведены результаты тестирования метода расчета 2(1^1<1

конвекции.

Результаты расчетов предлагаемым методом сравниваются с данными других авторов при небольшой и высокой (до г = 6-103) надкритичности. В последнем случае учитывалось [129x129] гармоник - точек, при этом различие по числу Нуссельта не превосходило 3.9% и, более того, проведенная серия расчетов подтвердила при г ~ З-Ю3 смену степенного закона для числа Нуссельта от надкритичности с Ми -1 ~ гш на Л'м -1 ~ /п.

Анализ результатов серии методических расчетов при г = 500 и а = 1 с учетом [33x17], [65x33], [129x65] и [257x129] гармоник-точек показал сходимость всех рассматриваемых интегральных величин.

Анализ результатов методических тестов показал, что при а = 1 в 2с1,пд1(1 расчетах достаточно учитывать [129x65] гармоник-точек при 10<г<2-103, [257x65] - 2-103<г<4-103 и [513x65] - 4-103 < г < 7-103. Показано, что при этом всегда Нге < А.

Сравнение результатов 2<1,г^с1 расчетов с теорией Р. Крайчнана показало, что при 500 < г < 7-103 в профиле пульсаций вертикальной скорости идентифицируются два участка, на которых величина пульсаций приближен-

1/3

но следует предсказанным теоретически степенным законам и , при этом рассчитанное число Нуссельта Л'м = 1.064-г1/3 близко к предсказанному теоретически при 250 < г < 7-103, с практическим совпадением значений при 103 <г <4103.

В п.2.4.3 приведены результаты тестирования метода расчета Зс1,Ггее конвекции.

В качестве первого теста рассчитывается двумерная, стационарная, ва-ликовая конвекция со сравнением результатов ЗсЦгсе и 2(1,й*ее расчетов при г = 2.2, Рг - 6.8, а = 0.09827 и ¡5 = к. В 2(1,Ггее расчете учитывалось [513x17] гармоник, а в Зс1Дгее - [513x5x17] и по координате у проводилось усреднение. Рассматриваемые характеристики течения оказались близки, в частности, относительное среднеквадратичное отклонение полей температуры составило 1.5%, кинетической энергии - 0.36%, при этом профили средней температуры отличались на 0.4% от максимальной температуры, а число Нуссельта - на 1.5%.

Затем расчетом с учетом различного числа гармоник [33x33x33], [65x65x65] и [129x129x129] при а = /3 = ж и г = 950 проверяется сходимость

интегральных характеристик.

Анализ результатов 3d,free расчетов показал, что при надкритичности г < 950, а = р = п и учете [65x65x65] гармоник всегда IIгг < А.

В третьей главе ’’Численное моделирование конвекции Рэлея-Бенара” представлены результаты численного исследования различных режимов конвекции Рэлея-Бенара, от стационарных до сложных, переходных.

В п.3.1 представлены результаты исследования двумерной, стационарной, валиковой конвекции, проведено сравнение экспериментальных и численных (2d,free, 2d,rigid и 3d,free) результатов при г < 4, Рг = 6.8 и а = 0.09827. В двумерных расчетах учитывалось [513x17] гармоник, а в трехмерном- [513x5x17].

При г = 2.2 показано, что результаты всех численных расчетов с графической точностью соответствуют экспериментальным данным по профилю средней температуры и изотерме полной температуры. Полученное в 2d,rigid расчете значение числа Nu на 4% отличается от полученного в эксперименте, а в 2d,free и 3d,free расчетах соответствующие отклонения больше и составляют примерно 21%.

В диапазоне надкритичности 1 < г < 4 результаты 2d,rigid расчета также хорошо соответствуют результатам эксперимента по числу Нуссельта и с удовлетворительной точностью отражают рост доминирующей длины волны стационарной, валиковой конвекции при увеличении надкритичности.

В п.3.2 рассматривается одно- и двухвихревая стационарная конвекция в квадратной области. Эти решения представляют собой интересный пример стационарного и крупномасштабного 2d,free решения, которое удается рассчитать при г < 104 для одновихревого и г < 5-104 - двухвихревого режима конвекции.

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

В п.3.3 рассмотрено поведение вихревого масштаба /те = к/Кт как функции надкритичности (рис. 2). Вычисляемый по кинетической энергии,

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

В трехмерных расчетах при г > 150 вихревой масштаб уменьшается, приближенно следуя степенному закону /те = 1.69т'0084, 1те~ \.1-гшг.

А данные 2с1,Ггее расчетов показывают, что при достаточно высокой надкритичности характерный горизонтальный размер вихря принимает значение порядка размера области Кте ~ а = 1. Учитывая также, что стационарное решение становится периодическим при г = 31.8, заключаем, что наблюдаемая при г >36 в 2сЦгее расчетах тенденция к увеличению вихрево-

В 2d,rigid расчетах (рис. 3) тенденция увеличения вихревого масштаба менее выражена и не так однозначна, из-за наличия конкурирующего механизма генерации завихренности на горизонтальных грани-Рис. 2. Среднее волновое число в го- цах На рис 3 представлеиь1 дашые

ргаонтальном направлении, кривая 2d ng]d расчетов (кривая h

1 - 3d,free и 2 - 2d,free расчеты, ц = ^ З06.г<шз _ да 9„/о уменьшенное

а Р I значение волнового числа, соответст-

вующего наиболее быстрому росту по линейной теории (2), 2.98 г0084 - на 60% увеличенное значение среднего волнового числа, полученного в 3d,free расчетах (3) и асимптотическое значение Кте = 5 (4).

Из рис. 3 видно, что при г ~ 250 происходит перестройка течения, связанная с резким, при-К)’ u>- in5 г мерно на 27% уменьшением сред-

Рис. 3. Среднее волновое число в него волнового числа с дальней-2d, rigid расчетах шим монотонным выходом на

асимптотическое значение Ктс = 5 при г ~ 1.5-103. Отметим также, что среднее волновое число в 2d,rigid расчетах при 400 < г <1.5-103 растет медленнее, чем в трехмерном решении.

В п.3.4 изучаются спектры скорости и температуры, полученные в 2d,free и 3d,free расчетах. Важность такого исследования обусловлена тем,

го масштаба связана с развитием нес

{Ж;

что наличие четко идентифицируемых пространственных спектров указывает на развитый характер течения и на то, какие физические механизмы доминируют. Отдельно рассмотрены случаи умеренной относительной горизонтальной протяженности области к (а = 1) в 2d,free и 3d,free и большой 4л (а = 0.25) - в 2d,free расчетах.

Полученные в 2d,free и 3d,free расчетах временные энергетические спектры температуры в центре конвективной ячейки сопоставляются с экспериментальными данными в газообразном Не при криогенной температуре. В эксперименте и расчетах совпадало значение надкритичности (410 и 6400 в 3d,free и 2d,free расчетах, соответственно) и число Прандгля (Рг = 0.8). Большая консервативность спектров и их слабая зависимость от параметров расчета и размерности позволяет сопоставлять данные эксперимента с результатами двумерного расчета.

I На рис. 4 приведен временной

| спектр температуры, полученный в

I 3d,free расчете. Здесь в качестве харак-

5 7

\ терной частоты / выбрано v/D , а

j fd = VRe3M/D, где D есть расстояние меж. О i

................. ..j ду горизонтальными границами. Видно

хорошее согласование результатов рас-

Рис. 4. Временной спектр тем-

четов с данными эксперимента, причем

пературы, кривая 1 - 3d, free рас-

существенное различие наблюдается чет, 2 - эксперимент по . ,

только на диссипативных частотах f ~fd-конвекции в газообразном Не „

Результат сопоставления качественно не

при 5°К, X.-Zh. etal. Wu, 1990

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

Показано, что при умеренной горизонтальной протяженности области к в 2d,free и 3d,free расчетах наблюдаются качественно одинаковые спектры -стратификационные, связанные с действием силы плавучести для скорости и к'14 и к'5В для температуры. Физическая природа первого из температурных спектров пока неясна, а второй указывает на поведение температуры как пассивной примеси.

На рис. 5 приведен спектр скорости в горизонтальном х направлении при /• =950 и стратификационный спектр Болджиано-Обухова (БО) k'ws.

На рис. 6 приведен температурный спектр в горизонтальном у направлении при г = 950 и спектр пассивной примеси и'5/3. В другом горизонтальном направлении х, в трехмерных расчетах в спектре температуры идентифицируется участок со степенным законом к'2' . Этот спектр отмечен также и в двумерных расчетах при г> 103.

В вертикальном направлении г, в спектре температуры в двумерных и трехмерных расчетах видны участки со степенными законами к'2А, к'513 и к~5ГЗ, соответственно.

Отмеченное качественное соответствие между спектрами в двумерных и трехмерных расчетах наблюдается практически во всем исследуемом диапазоне надкритичности, а именно, при г < 3.4-104 в 2с1, (гее и при г < 950 - в Зс1,1гее расчетах.

Некоторое отличие спектров скорости в 2ё,Ггсе и Зс1,1Уес расчетах состоит в том, что в дополнение к наблюдающемуся в двумерных и трехмерных расчетах в горизонтальном направлении х стратификационному спектру БО к'Ш5, в трехмерных расчетах, в другом горизонтальном направлении у, виден стратификационный спектр Ламли-Шура и'3.

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

Другая картина перестройки спектров наблюдалась при большой гори-

і -

----1.....■-------1 —і----і-----і-----и

ft.fi Ц- *:>*; М

Рис. 5. Спектр скорости, кривая 1 - Зсі,/гее расчет, знак ■ соответствует наиболее быстрому росту в линейном приближении, а вертикальная прямая - границе интервала аппроксимации

■ Г_____1_____________1_____I____1—___I—

Й* • ЧйГ.»: м

Рис. 6. Спектр температуры в горизонтальном у направлении, обозначения как на рис. 5

зонтальной протяженности области 4%. На рис. 7-10 приведены спектры скорости и температуры в горизонтальном направлении х при а = 0.25, знаком

• показаны результаты 2сЦгее расчетов. Диссипация кинетической энергии, ее генерация и потоки к малым и большим масштабам растут при увеличении горизонтальной протяженности области и надкритичности, это обуславливает наблюдаемые перестройки спектров скорости и температуры.

При 500 <г< 4 • 103 длинноволновый участок спектра скорости состоит из двух ветвей, со степенными законами, отражающими конкурирующие механизмы, связанные с силой плавучести и каскадным процессом переноса энергии. Конечно, сосуществование двух спектров обуславливает некоторую неточность определения степенных законов и границы этого сосуществования.

При 4-103 < г < 104 в спектре скорости видны участки с четко идентифицируемыми степенными законами, соответствующими каскадным процессам переноса кинетической энергии из масштаба генерации в малые масштабы с последующей диссипацией (к'3 - прямой каскад энстро-фии) и в большие, с формированием крупномасштабной структуры течения (к5'3 - обратный или красный каскад энергии) (рис. 7). Масштаб генерации здесь соответствует наиболее быстрому в линейном приближении росту.Однако, после формирования обратного каскада энергии его энергетическая роль падает, что обусловлено отсутствием диссипации на больших масштабах. В свою очередь, рост диссипации энергии при увеличении надкритичности и связанное с этим усиление энергетической роли прямого каскада энстрофии обуславливает при г > 104 тенденцию установления единого степенного закона в спектре скорости к'26, близкого к закону энстро-фийного каскада к'3 (рис. 8).

Рис. 7. Спектр скорости при г = 6111

I'

Рис. 8. Спектр скорости при г = 3.4-10'

і:-< В';. ^

Рис. 9. Спектр температуры, г =1.25-11?

А увеличение роли силы плавучести, неизбежное при увеличении горизонтальной протяженности области, обуславливает при 750 < г < 2-103 стратификационный спектр БО к'115 для температуры (рис. 9).

Дальнейшее развитие обратного каскада энергии и связанное с этим формирование крупномасштабной структуры течения сопровождается перекачкой энергии пульсаций скорости в большие масштабы, в результате чего поле скорости становится крупномасштабным при относительно низком уровне пульсаций, и возникает своеобразный вязкоконвективный интервал, где пульсации температуры управляются крупномасштабным полем скорости, что обуславливает появление в спектре температуры участка со спектром Бэтчелора к'\ наблюдаемого в 2d,free расчетах при 2-103<г<6-103(рис. 10).

Уменьшение энергетической роли обратного каскада энергии отражается и на спектре температуры - при г > 6-103 показатель длинноволнового участка спектра температуры совершает колебания в интервале [-1 ,-0.7], с наиболее вероятным значением - 0.8.

В п.3.5 рассмотрено формирование крупномасштабной структуры течения в двумерных расчетах. Проведением серии 2d,free расчетов в областях различной горизонтальной протяженности (ж, 2ж, Аж и 8л) исследовано последовательное включение больших масштабов, а именно, показано, что при г = 103 горизонтальной протяженности области к недостаточно, при г = 6-103 недостаточно длины области 2%, а при г = 1.7-104 - 4тг. Здесь выражение недостаточная длина области обозначает тот факт, что горизонтальная протяженность области становится определяющим решение параметром.

Показано также, что в 2d,free расчетах формирование крупномасштаб-

Ы <.у

Рис. 10.

г = 3-10*

Спектр температуры,

!1||:| //V

Рис. 11. Функция тока при

t = 7.65-1 а4

...

4 '

нои структуры течения является доминирующим и определяющим.

На рис. 11 и 12 при г = 3-104 и 4 = 4л: приведена функция тока в последовательные моменты времени, значения среднего волнового числа равны примерно 1 и 4.26л = 4л на рис. 11 и 12, соответственно. Рис. 12 показывает, что при />0.1 режим конвекции крупномасштабный и практически одновихревой. ________________________

Доминирование процесса формирования крупномасштабной структуры проявляется, в частности, в выходе на степенной закон числа Нуссельта Nu = 0.639-г0 371 при г > 6-103, с показателем заметно выше 2/7, 0.31 и 1/3 и в нефизичном росте среднеквадратичной температуры, пульсаций температуры и вертикальной скорости; а также в близости показателей степенных законов для числа Нуссельта, кинетической энергии, средне-Рис. 12. Функция тока при квадратичной скорости и энстрофии в 2d,free t = 2.29KI1 расчетах и стационарном, одновихревом те-

чении (п.3.2).

Отмечено, что этот процесс менее выражен и не так однозначен (см. рис.

2 и 3 при 4 = л) в 2d,rigid расчетах из-за наличия конкурирующего механизма генерации завихренности на горизонтальных границах.

Например, среднее волновое число при умеренной горизонтальной протяженности области (4 = л), после скачкообразного уменьшения (при г - 250) с монотонным увеличением выходит (г~ 1.5-Ю3) на асимптотическое значение Кте = 5 (рис. 3), а число Нуссельта следует степенному закону Nu = 0.817-г0,367 при г > 2-103, причем показатель этого степенного закона близок к полученному в 2d,free расчетах.

Однако, при конвекции в области малой горизонтальной протяженности (4 = 1) развитие крупномасштабной структуры в 2d,rigid расчетах в значительной степени подавлено малым размером области и это делает возможным рост среднеквадратичного волнового числа при увеличении надкритичности примерно как Ъ.1-гш при г<6-103.

В п.3.6 рассматриваются интегральные характеристики трехмерной и

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

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

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

В п.3.3 отмечено, что основное различие меду трехмерной и двумерной конвекцией связано с разным поведением вихревого масштаба, который рассматривается как функция надкритичности. Результаты 2d,free расчетов показывают, что с ростом надкритичности характерный горизонтальный размер вихря растет до значения порядка размера области, а в трехмерных расчетах характерный горизонтальный размер вихря уменьшается примерно обратно пропорционально корню двенадцатой степени из надкритичности.

Такое поведение вихревого масштаба показывает, что при двумерном моделировании можно ожидать хорошего соответствия вычисленных интегральных характеристик с данными эксперимента только при невысокой (до г ~ 100) надкритичности. Причем, практически полное качественное совпадение пространственных спектров пульсаций температуры и скорости в 2d,free и 3d,free расчетах при конвекции в области относительной горизонтальной протяженности % во всем рассматриваемом диапазоне надкритичности показывает, что для консервативных величин, мало изменяющихся при увеличении надкритичности и пространственного разрешения, интервал соответствия может быть много больше.

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

- '

»»«о

личении надкритичности.

На рис. 13 приведены значения среднеквадратичных пульсаций температуры как функций надкритичности, где кривая 1 - данные 3d,free расчетов,

2 - 2d,free, а 3 - 2d,rigid.

Из рис. 13 видно, что величина среднеквадратичных пульсаций температуры следует степенным законам: q’ = 0.163-г 2/15 при г > 150 в 3d,free и q ’ = 0.316-f02 при 40 < г < 250 - 2d,free расчетах. Оба показателя степенных за-

.-'Я,1 '»'’г,* О' л; ! ; - д*,_.

конов совпадают с полученными Рис. 13. Пульсации температуры в лабораторных экспериментах D.E. как функция надкритичности Fitzjzrrald, 1976 и X.-Zh. Wu,1992, соответственно. Пульсации температуры при г > 250 в 2d,free и во всем рассматриваемом диапазоне надкритичности в 2d,rigid расчетах не следуют какому-либо степенному закону и, более того, пульсации температуры в двумерных расчетах при достаточно большой надкритичности показывают нефизичную тенденцию увеличения.

Деление величины пульсаций скорости на Ргп позволяет сравнивать данные с различным числом Рг. На рис. 14 приведены значения среднеквадратичных пульсаций вертикальной скорости как функций надкритичности, где кривые 1-3 обозначают, как и на рис. 13 данные расчетов, кривые 4-6 -экспериментальные данные J.W. Deardorff,1967; W.V.R. Malkus,1954; F.М. Garon, 1973, соответственно и 7 - результат численного расчета R.M. Кегг,1996.

\

1Л>

50

‘100

1200 г

Рис. 14. Пульсащи вертикальной скорости

Из рис. 14 видно, что величина среднеквадратичных пульсаций вертикальной скорости в трехмерном расчете следует степенному закону:

IV1 = 15.8т0 393 при г > 150, причем показатель этого степенного закона близок к полученным в экспериментах и численном расчете с жесткими горизонтальными границами. В двумерных

расчетах аналогичное соответствие наблюдается при г < 250, а при достаточно высокой надкритичности пульсации вертикальной скорости показывают нефизично быстрый рост.

Нефизично высокие значения пульсаций температуры и вертикальной скорости при достаточно высокой надкритичности в двумерных расчетах обусловлены формированием крупномасштабной структуры течения.

На рис. 15 знаком (•) представлено полученное в ЗсЦгее расчетах число Ке как функция надкритичности.

; , ; чем показатель этого степенного закона

.. .. .....■..Т:м.... ;........ практически совпадает со значением

Рис. 15. Число Ке в Зс1/гее рас- 0.49, полученным в экспериментах по

Простые формулы для среднеквадратичной скорости и кинетической энергии Г’ппз = Ке Рг и Ек= п7-]>г2-Ре2/2 показывают, что соответствие показателей степенных законов наблюдается также и для среднеквадратичной скорости и кинетической энергии.

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

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

3. Показано, что, несмотря на то, что численный метод строится как ме-

Видно, число Рейнольдса в трехмерных расчетах следует степенному закон}' Ке = 1.46-г0497 при г > 150, при-

четах

турбулентной конвекции (Х.-2. 1992; X. СЬауаппе,2001).

ОСНОВНЫЕ РЕЗУЛЬТАТЫ И ВЫВОДЫ

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

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

5. Средствами нелинейного анализа, проведенного на модельной нелинейной задаче, показано, что понижение порядка аппроксимации численного метода расчета 2d,free конвекции до первого по времени приводит лишь к незначительному понижению точности вычислений, при этом с первым порядком вычисляется фаза гармонического решения, в то время как его амплитуда - со вторым. Показано, что для практических вычислений можно использовать схему первого порядка аппроксимации по времени с вычислением скоростей по функции тока, полученной на первом этапе расщепления.

6. Показано, что результаты 2d,rigid расчета стационарной, валиковой конвекции при невысокой надкритичности 1 < г <4 с удовлетворительной точностью отражают зависимость доминирующей длины волны от надкритичности, хорошо согласуясь с данными эксперимента по числу Нуссельта, профилю средней температуры и изотерме полной температуры.

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

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

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

10. Несмотря на наблюдаемые количественные различия между резуль-

31

татами ЗсЦгее расчетов и экспериментом, трехмерные расчеты дают правильные показатели степенных законов изменения среднеквадратичных пульсаций температуры, вертикальной скорости, числа Рейнольдса, среднеквадратичной скорости и кинетической энергии от надкритичности при г > 150. В двумерных расчетах аналогичное соответствие наблюдается при сравнительно невысокой надкритичности (до г ~ 250) по среднеквадратичным пульсациям температуры (в 2с1,Ггее) и вертикальной скорости (в 2сЦгее и 2d,rigid расчетах).

11. Закон нарастания числа Нуссельта, нефизичное поведение пульсаций вертикальной скорости и температуры, а также близость показателей степенных законов для числа Нуссельта, кинетической энергии, среднеквадратичной скорости и энстрофии в 2d,free расчетах и стационарном, одновихревом течении показывают, что в двумерных расчетах со свободными горизонтальными границами при высокой надкритичности формирование крупномасштабной структуры является доминирующим фактором, определяющим характеристики течения.

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

В области умеренной (я) горизонтальной протяженности в 2d,free и Зd,free расчетах:

13. Действие силы плавучести обуславливает стратификационные спектры для скорости, а доминирование конвективного переноса для температуры - к'24 и спектр пассивной примеси к'51ъ. Причем, эти спектры наблюдались в двумерных и трехмерных расчетах практически во всем исследованном диапазоне надкритичности.

В области большой (4л) горизонтальной протяженности в 2d,free расчетах:

14. Показано, что при 500 < г < 4-103 длинноволновый участок спектра скорости состоит из двух ветвей, со степенными законами, отражающими конкурирующие механизмы, связанные с силой плавучести и каскадным

32

процессом переноса энергии, а при 4-103 < г < 104 в спектре скорости идентифицируются степенные законы к'3 и к'513, соответствующие прямому и обратному каскадным процессам переноса энстрофии и энергии. При дальнейшем увеличении надкритичности (r> 104) в спектре скорости видна тенденция к установлению единого степенного закона к'26 - искаженного закона каскада энстрофии.

15. Показано, что при умеренной надкритичности 500 < г < 2-103 в спектре температуры виден стратификационный спектр БО к'7/5, сменяющийся при 2-103 < г < 6-Ю3 спектром Бэтчелора кл. А при г > 6-103 показатель длинноволнового участка температурного спектра совершает колебания в интервале [-1,-0.7], с наиболее вероятным значением - 0.8.

16. Перестройки спектров температуры (скорости) в двумерных расчетах связаны с формированием обратного каскада энергии при г ~ 2-103 (г ~ 4 • 103) и уменьшением его энергетической роли при г = 6Т03 (г~104).

Основные результаты диссертации опубликованы в работах:

1. Пачьшский И.Б. Метод численного моделирования конвективных течений //Выч. техн. 2000. Т.5. N.6. С.53-61.

2. Пачьшский И.Б. Линейный и нелинейный анализ численного метода расчета конвективных течений // Сиб. ж. вычисл. матем. 2004. Т.7. N2. С.143-163.

3. Пачьшский И.Б. Численное моделирование двумерной конвекции при высокой надкритичности //Успехи механики. 2006. N.4. С.3-28.

4. Папьшский И.Б. Численное моделирование двумерной конвекции, роль граничных условий//Изв. РАН. МЖГ. 2007. N.4. С.61-71.

5. Palymskiy I.B., Fomin Р.А. Hieronymus Н. The Rayleigh-Benard convection in gas with chemical reactions // Сиб. ж. вычисл. матем. 2007. T. 10. N.4. С.371-383.

6. Пальшский И.Б. Численное исследование спектров турбулентной конвекции Рэлея-Бенара//Нелинейная динамика. 2008. Т.4. N.2. С.145-156.

7. Palymskiy I.B., Fomin Р.А. Hieroiiymus H. Rayleigh-Benard convection in chemical equilibrium gas II (simulation of surface detonation wave initiation).

33

App. Math. Model. 2008. V.32. N.5. P.660-676.

8. Палымский И.Б., Герценштейн С.Я., Сибгатуллин И.Н. Об интенсивной турбулентной конвекции в горизонтальном плоском слое жидкости//Изв. РАН. ФАО. 2008. Т.44. N.I. С.75-85.

9. Палымский И.Б. Численное исследование спектров трехмерной конвекции Рэлея-Бенара//Изв. РАН. ФАО. 2009. Т.45. N.5. С.691-699.

10. Палымский И.Б. О численном моделировании трехмерной конвекции // Вест. Удмурт, универ. серия 1: матем., мех., компьют. науки. 2009. В.4. С. 118-132.

11. Палымский И.Б. О качественном различии решений двумерной и трехмерной конвекции И Нелинейная динамика. 2009. T.5.N.2. С.183-203.

12. Палымский И.Б. Численное исследование спектров трехмерной турбулентной конвекции // Изв. Сарат. универ. новая серия: матем., мех., информ. 2010. Т. 10. В.1. С.62-71.

13. Палымский И.Б. Численный метод расчета трехмерной конвекции // Сиб. ж. индустр. матем. 2010. Т.13. N.1. С.95-108.

14. Палымский И.Б. О моделировании сложных режимов конвекции Рэлея-Бенара//Сиб. ж. вычисл. матем. 2011. Т. 14. N.2. С. 179-204.

15. Пачымский И.Б. Турбулентная конвекция Рэлея-Бенара. Численный метод и результаты расчетов. Germany: LAP, 2011. 232 с.

16. Palymskiy I. В. Direct numerical simulation of turbulent convection // Fourth International Symposium on Finite Volumes for Complex Applications - Problems and Perspectives. Proceedings in book: Finite Volumes for Complex Applications IV - Problems and Perspectives. Wiley-ISTE, 2005. P.643-654.

17. Palymskiy l.B. Direct numerical simulation of turbulent convection // Fourth International Conference on Computational Heat and Mass Transfer (ICCHMT’05). Proceedings in book: Progress in Computational Heat and Mass Transfer, V.l. Paris, France, 2005. P. 101-107.

18. Palymskiy 1.В., Fomin P.A., Hieronymus H. Rayleigh-Benard convection in chemical equilibrium gas // Fourth International Conference on Computational

34

Heat and Mass Transfer (ICCHMT’05). Proceedings in book: Progress in Computational Heat and Mass Transfer, V. 1. Paris, France, 2005. P.l 16-122.

19. Palymskiy I.B. Numerical investigation of turbulent convection spectrums // International Conference Dynamics Days Europe 2009 (DDE09). Proceedings. Gettingen, Germany, 2009. P.232-233.

20. Palymskiy I.B. Numerical simulation of inverse energy cascade in twodimensional convection // International Conference Dynamics Days Europe 2009 (DDE09). Proceedings. Gettingen, Germany, 2009. P.229-232.

21. Palymskiy I.B. Numerical simulation of turbulent convection // 8th International Symposium on Engineering Turbulence Modelling and Measurements (ETMM8). Proceedings. Marseille, France, 2010. P.894-897.

Подписано в печать 10.01.2012 г. Формат 60x84 х!\6. Уел. печ. л. 2. Тираж 120 экз. Заказ №

630064, г. Новосибирск, ул. Ватутина, 71. ООО “ВАРГУС”

Введение

Глава 1. Постановка задачи о конвекции Рэлея-Бенара и численные алгоритмы

1.1. Постановка задачи

1.2. Метод численного расчета двумерной конвекции

1.2.1. Свободные горизонтальные границы

1.2.2. Жесткие горизонтальные границы 32 1.2.2.1. О вычислении завихренности на границе

1.3. Метод численного расчета трехмерной конвекции 3g

Глава 2. Анализ и тестирование численных алгоритмов

2.1. Линейный анализ двумерных алгоритмов

2.1.1. Свободные границы

2.1.2. Жесткие границы

2.2. Нелинейный анализ 58 2.2.1. Численные эксперименты

2.2.1.1. Модельная система уравнений

2.2.1.2. Исходная система уравнений

2.3. Линейный анализ метода расчета трехмерной конвекции

2.4. Тестирование и методические расчеты

2.4.1. Метод расчета 2d,free конвекции

2.4.1.1. Сравнение с результатами других авторов

2.4.1.2. О сравнении мгновенных и средних величин

2.4.1.3. Проверка сходимости

2.4.1.4. Учет большего числа цифр

2.4.1.5. Сравнение с результатом расчета псевдоспектральным методом

2.4.1.6. Определение необходимой пространственной разрешимости

2.4.2. Метод расчета 2d,rigid конвекции

2.4.2.1. Сравнение с результатами других авторов

2.4.2.2. Проверка сходимости

2.4.2.3. Сопоставление с теорией Р. Крайчнана 102 2.4.3. Метод расчета 3d,free конвекции

2.4.3.1. Сравнение на двумерном решении

2.4.3.2. Проверка сходимости при увеличении пространственной разрешимости

Глава 3. Численное моделирование конвекции Рэлея-Бенара j j

3.1. Стационарная валиковая конвекция

3.2. Стационарная конвекция в квадратной области

3.3. Вихревой масштаб

3.4. Спектры скорости и температуры

3.4.1. Некоторые качественные соображения о динамике спектров

3.4.2. Временной спектр температуры

3.4.3. Пространственные спектры при конвекции в области умеренной горизонтальной протяженности

3.4.3.1. Спектры скорости

3.4.3.2. Спектры температуры

3.4.4. Пространственные спектры при конвекции в области большой горизонтальной протяженности и каскадные процессы

3.4.4.1. Спектры скорости

3.4.4.2. Спектры температуры

3.4.4.3. Диаграмма спектров

3.5. Формирование крупномасштабной структуры течения

3.6. Интегральные характеристики

Классическая задача о конвекции Рэлея-Бенара в различных постановках исследовалась многими авторами численно [1-38] и экспериментально [39-57]. Из-за очевидной связи с прямым численным моделированием турбулентности большой интерес вызывают исследования при высокой надкритичности г = Ra/Raer, где Ra и Ra^ - число Рэлея и его критическое значение, а Рг - число Прандтля.

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

Рассмотрение конвекции со свободными от касательных напряжений горизонтальными границами имеет также и самостоятельный интерес, например, при изучении конвекции в мантии Земли [11,34] или приповерхностном слое океана [60]. В самом деле, в последнем случае формируется своеобразная трехслойная структура - воздух, полностью турбулизированный приповерхностный слой и относительно слабо турбулизированный внутренний слой океана, где изолировано расположенные очаги турбулентности имеют форму горизонтальных блинов [61].

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

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

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

Существенное отличие конвекции со свободными границами от конвекции с жесткими - возможное неравенство нулю значений вертикальной компоненты завихренности [59], что проявляется в различном вращении вокруг вертикальной оси, а направленная горизонтально скорость такого вращения, как следствие, обуславливает возможное расхождение в характеристиках, связанных с горизонтальной скоростью. Заметим, что при двумерном моделировании вертикальная компонента завихренности равна нулю тождественно, что делает возможным согласование среднеквадратичных пульсаций горизонтальной скорости при двумерной конвекции со свободными и жесткими граничными условиями [63]. Интересно, что при этом результаты 2сЦгее расчетов согласуются с экспериментальными данными по величине среднеквадратичных пульсаций горизонтальной скорости и по числу Нуссельта даже при высокой надкритич-ности [63].

Ширина области влияния граничных условий пропорциональна вихревому масштабу, который в трехмерной конвекции уменьшается с ростом надкритич-1 /1 ности как г" [64] и, следовательно, роль граничных условий уменьшается при увеличении надкритичности и сказанное о правомерности сравнения характеристик трехмерной конвекции, связанных с вертикальной скоростью, числом Рейнольдса, среднеквадратичной скоростью, кинетической энергией и температурой, тем более верно при достаточно высокой надкритичности. Утверждение об уменьшении вихревого масштаба с ростом надкритичности в трехмерной конвекции со свободными границами [64], остается справедливым и для трехмерной конвекции с жесткими [3].

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

Между тем, число Рейнольдса (вычисленное по среднеквадратичной скорости и высоте слоя) является относительно медленно растущей функцией надкритичности в конвекции Рэлея-Бенара и Яе ~ 375 при г = 3.4-104 и Рг = 10 [65]. В трехмерных расчетах настоящей работы надкритичность меньше, чем в двумерных и, как следствие, при трехмерном моделировании число Рейнольдса не превышает 44 с максимальным ростом линейных возмущений как е1981 при г = 950 и Рг = 10.

В задаче о конвекции Рэлея-Бенара со свободными границами собственные функции задачи линейной устойчивости выражаются через синусы и косинусы [66] и это обуславливает высокую эффективность применения спектральных (псевдоспектральных) методов. В трехмерной задаче спектральные методы примерно на два порядка эффективнее конечно-разностных, причем, эта порядковая оценка отношения количества точек дискретизации к числу гармоник относительно слабо зависит от конкретных реализаций (в том числе и порядка аппроксимации) конечно-разностного и спектрального методов и определяется тем, что для любого конечно-разностного оператора интервал аппроксимации составляет лишь треть всего спектра и, таким образом, длина волны каждой учитываемой гармоники должна содержать не менее шести интервалов разностной сетки [67].

При расчете конвекции спектральным методом нужно решить две проблеI мы: интегрирование по времени жесткой системы уравнений и вычисление нелинейных членов. Во всех известных автору работах, кроме [15], для интегрирования по времени использовались разностные схемы, что ограничивало шаг по времени и число Рэлея [68]. На наш взгляд, для расчета конвекции со свободными границами наиболее оптимально использование аналитических формул [68,69] или матричной экспоненты [15].

Для проведения расчетов турбулентной конвекции широко используется псевдоспектральный метод в различных его модификациях [1,2,9,13,14,17,2022], недостатками которого являются невыполнение уравнения неразрывности на последнем дробном шаге [17,20] и использование не очень удачной схемы расщепления по физическим процессам, с учетом на одном дробном шаге нелинейных членов и плавучести [1,2,17,20-22], что ограничивает шаг по времени. На наш взгляд, более правильно рассматривать на отдельных дробных шагах линейное и нелинейное развитие возмущений.

Проблема вычисления нелинейных членов рассматривалась в работах [27,62]. В работе [62] рассмотрены три различных способа вычисления нелинейных членов, основанных на применении алгоритма быстрого преобразования Фурье, а в [27] - нелинейные члены вычисляются в смешанном спектрально-физическом пространстве.

Численные алгоритмы с вычислением нелинейных членов в физическом пространстве на разностной сетке, число узлов которой в каждом направлении в два раза превышает число гармоник в том же направлении, обладают высокой точностью, но при этом возникает проблема их эффективности. На наш взгляд, наиболее оптимальный компромисс между эффективностью численного алгоритма и точностью вычислений дает учет нелинейного переноса в физическом пространстве по разностным схемам на отдельном этапе расщепления [68,69].

Для расчета двумерной конвекции со свободными границами возможно применение варианта псевдоспектрального метода (модифицированного метода Галеркина), основанного на вычислении нелинейных членов по аналитическим формулам в физическом пространстве в точках разностной сетки, число узлов которой в каждом направлении в два раза превышает число гармоник в том же направлении и схеме Рунге-Кутта четвертого порядка точности для интегрировании по времени. Подобный численный метод применялся для расчета стохастической двухдиффузионной конвекции [30,31]. Такой метод обладает высокой точностью и результаты расчета по такой методике хорошо согласуются с результатами расчетов предлагаемым спектрально-разностным методом. Но, к сожалению, уже отмеченная выше низкая эффективность подобного численного метода, обусловленная избыточно точным вычислением нелинейных членов и использованием явной схемы для интегрирования по времени, сделала возможным проведение только нескольких тестовых расчетов.

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

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

Применение конечно-разностных методов может быть оправдано при проведении расчетов с небольшими надкритичностями [8,10], в областях со сложной геометрией [4,6] или жесткими стенками [26]. Использование конечно-разностных аппроксимаций возможно также по переменным, по которым применение спектральных разложений затруднительно [70,71].

Для расчета турбулентной конвекции со свободными границами при высокой надкритичности в [12] используется специальный численный метод, основанный на методе характеристик с использованием сплайнов. На наш взгляд, для решения этой задачи более эффективны спектральные методы с разложением по синусам и косинусам, совпадающими с собственными функциями линейной теории устойчивости [66].

Расчеты трехмерной конвекции проводились в [1-9,18] с жесткими и в [1019] со свободными граничными условиями. Использование суперкомпьютеров сделало реальным прямое численное моделирование турбулентной конвекции [1-7,12-14], но, к сожалению, большая сложность такого моделирования, его огромная ресурсоемкость и, как следствие, значительная стоимость обусловили небольшое число подобных работ и проведение надежных и достаточно подробных расчетов пока еще в перспективе. Заметим, что в [1,2,13,14,17] использован псевдоспектральный метод, причем в [2] использована программа [1], а в [3-7] - метод конечных разностей, хотя, как уже отмечалось выше, в таких задачах спектральные (псевдоспектральные) методы более эффективны

72]. В работах [13,14] рассчитывалась трехмерная конвекция в воздухе со сво

3 4 бодными горизонтальными границами при надкритичности г = 9.8-10 и 3.3-10 , соответственно, причем учитываемое число гармоник [96x96x97] кажется явно недостаточным для такой высокой надкритичности.

Очень интересен совместный анализ результатов расчетов, выполненных со свободными и жесткими горизонтальными границами при И.а = 106 и 10 < Рг < 10 [18] (с использованием данных [19]). Показано, что в обоих случаях качественно похожие зависимости чисел Нуссельта и Рейнольдса от числа Прандтля претерпевают перестройки при Рг = 1, причем при Рг < 1 различие чисел Нуссельта и Рейнольдса выражено слабее и более того, числа Рейнольдса (и, следовательно, кинетическая энергия и среднеквадратичная скорость) при Рг —> 0 практически совпадают.

С целью многократно уменьшить используемые вычислительные ресурсы, рассматривают конвекцию жидкости в двумерной постановке. При конвекции в узком канале оси валов параллельны короткой стенке канала и перпендикулярны длинным стенкам. При конвекции в лабораторной модели бесконечно длинного прямого канала - кольцевом канале между коаксиальными цилиндрическими стенками, валы располагаются в радиальном направлении, перпендикулярно стенкам, если канал не слишком широк [59].

Конечно, при надкритичности порядка 10 (Рг = 10) стационарные валы становятся неустойчивыми к трехмерным возмущениям, и течение становится трехмерным [47]. Но, логично ожидать, что в начале своего возникновения трехмерность еще не является доминирующим фактором и законы изменения основных интегральных величин еще мало отличаются от двумерных. Грубую оценку величины надкритичности, выше которой трехмерность является доминирующей можно получить из рассмотрения динамики вихревого масштаба [64].

В самом деле, в работе [64] показано, что кардинальное отличие двумерных течений от трехмерных - различное поведение вихревого масштаба, который в двумерной конвекции растет, а в трехмерной - уменьшается с ростом надкритичности. Указанная тенденция роста вихревого масштаба в двумерной конвекции проявляется при г > гт, гт = 36 (при а = 1 и Рг = 10), а именно, при г < гт масштаб вихрей в двумерной, также как и в трехмерной конвекции уменьшается, а при г > гт - увеличивается до величины порядка горизонтального размера области. Естественно ожидать, что и это различие в поведении вихревого масштаба становится доминирующим не сразу после его появления. Таким образом, можно ожидать примерного соответствия степенных законов изменения средних величин при увеличении надкритичности в двумерной и трехмерной конвекции примерно до значений порядка 100. Причем для консервативных величин, например, пространственных и временных спектров, слабо зависящих от надкритичности, пространственного разрешения и размерности [58,60], интервал соответствия может быть много больше.

Сказанное подтверждает то, что в двумерной и трехмерной конвекции со свободными границами в области умеренной горизонтальной протяженности 71, пространственные спектры пульсационных полей температуры и скорости совпадают практически полностью [58,60].

Отметим также качественное соответствие временных спектров числа Нуссельта при г = 950, полученных в двумерных и трехмерном расчетах и согласование результатов двумерного расчета с экспериментом по энергетическому спектру среднеквадратичных пульсаций температуры в центре ячейки при надкритичности г = 6.4-10 [58].

В работе [73] рассмотрена двумерная конвекция, возникающая при подогреве сбоку, при этом по некоторым характеристикам было получено хорошее согласие расчетных и экспериментальных данных.

В работах [17,27-38] рассматривалась двумерная конвекция со свободными, а в работах [20-26] - жесткими граничными условиями. В работе [26] исследовался начальный этап развития конвекции при внезапном нагреве нижней границы в области большой горизонтальной протяженности при Ид < 106 и Рг = 0.7. Сложные режимы двумерной конвекции при высокой надкритичности исследовались в работах [21-23,35-38]. В [30,31] рассматривались стохастические режимы двухдиффузиоиной конвекции при сравнительно невысокой над-критичности и исследовалась структура стохастического аттрактора. В работах [34,35,38] конвекция при высокой надкритичности рассматривается в приближении бесконечно большого числа Прандтля. В работе [34] рассматриваются стационарные решения, в [35] на основании расчетов двумерной конвекции со свободными граничными условиями получен закон для числа Нуссельта № ~ Иа0'301, а в [38] показано, что интегральные характеристики конвекции сильно зависят от аспектного отношения (длина области отнесенная к высоте), например, число Нуссельта изменяется при этом в несколько раз.

В работах [21,22] были рассчитаны конвективные течения в квадратной области при огромной надкритичности (до г ~ 9.6-104), но разрешимость в горизонтальном направлении при такой надкритичности кажется явно недостаточной (129 и 257 гармоник в работах [21] и [22], соответственно).

В работе [37] рассматриваются периодические в горизонтальном направлении двумерные конвективные течения в области горизонтальной протяженности 3 со свободными горизонтальными границами. Расчеты выполнены псевдоспектральным методом с учетом до [1024x3076] гармоник, причем число

8 14

Рэлея изменяется от 10 до 10 . Некоторые результаты этой работы требуют уточнения или дополнительных пояснений, а именно, странным кажется постоянное при г > 1012 число Нуссельта и явно заниженный показатель степенного закона для числа Рейнольдса (Яе ~ г0'25). Возможно, что указанное несоответствие известным экспериментальным данным и здравому смыслу обусловлено плохой сходимостью тригонометрических рядов, часто наблюдающейся при решении нелинейных задач с учетом чрезмерно большого числа гармоник.

В работе автора [63] описаны результаты численного моделирования двумерной конвекции со свободными и жесткими граничными условиями при высокой надкритичности (до г = 3.4-104 для свободных и г = 7-10 - жестких граничных условий), но результаты расчетов при сравнительно небольшой надкритичности г < 103, а также степенные законы для числа Нуссельта и других интегральных величин потребовали уточнений.

Большое значение имеет исследование динамики пространственных спектров температуры и скорости в двумерных и трехмерных расчетах, так как наличие определенных спектральных законов указывает, во первых, на развитый характер течения и, во вторых, на то, какие конкретные физические механизмы доминируют.

Диссипация и генерация энергии турбулентности растут при увеличении надкритичности примерно как г1'4. При достаточно большой надкритичности большие потоки переносимой энергии обуславливают формирование инерционных интервалов и спектров.

Известно два основных сценария развития конвективной турбулентности [74].

Сценарий Колмогорова, при котором предполагается, что температура ведет себя как пассивная примесь, предполагает наличие двух инерционных интервалов переноса энергий пульсаций температуры и скорости, с формированием одинаковых спектров к"5/3, где к - волновое число в случае зависимости от пространственных переменных либо частота - от времени. Сила плавучести здесь существенной роли не играет.

Напротив, Р. Болджиано и А. Обухов (БО) предположили существование инерционного интервала для переноса энергии пульсаций температуры и в области больших масштабов равенство по порядку величины членов плавучести и нелинейного переноса. Это привело к спектрам к"7/5 и к"п/5 для температуры и для скорости, соответственно. А при большом числе Прандтля возможен баланс между силой плавучести и силой вязкости, что приводит к спектру к"5 для скорости.

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

1 э инерционного интервала переноса энстрофии (к' для температуры и к" для скорости), а в области больших - обратного (красного) каскада переноса энергии пульсаций скорости, направленного от масштаба генерации в область больших масштабов.

В экспериментах по турбулентной конвекции для пульсаций температуры наблюдался спектр пассивной примеси к"5/3, БО к"7/5 и к"2'4 [42-44,49]. Для пульсаций скорости наблюдались спектр БО к"11/5 и к"1'35, но спектр Колмогорова к"5/3 не обнаружен [44,46]. Физические механизмы появления спектров к"2'4 и к"1'35 для пульсаций температуры и скорости, отмеченных в экспериментальных работах [42] и [46], соответственно, пока еще не получили достаточного теоретического обоснования.

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

7/5 1 5аз лучены спектры БО к" и к" [1,4], но спектр пассивной примеси к не обна

СИ 1 <у ружен. А для пульсаций скорости - спектры к" , к" и к [1,4,12], но спектр БО к"11/5 не наблюдался. Отметим также не очевидную идентификацию спектров к"1 и к"5/3 для температуры и скорости, соответственно, в работе [1].

Участки со степенными законами к"1 для температуры и скорости в работе [4] непосредственно предшествуют диссипативному интервалу и при их формировании могли быть существенны численные эффекты, неизбежно возникающие на краю спектров разностных операторов [67].

В [35] проведено моделирование турбулентной конвекции по двумерной

1 ■) модели бесконечного числа Прандтля, получены спектры к", к" для пульсаций температуры и к"2 - скорости.

В работе [36] проведено моделирование двумерной конвекции со свободными от касательных напряжений горизонтальными границами и условием о периодичности в горизонтальном направлении при Иа = 10 и Рг = 1, для температуры получен спектр пассивной примеси к"5/3. Однако, сделанный на основании анализа температурного спектра вывод о существовании красного каскада энергии представляется не достаточно обоснованным, так как анализ спектра скорости и потока энергии по спектру в этой работе не проводился.

В трехмерной турбулентности кинетическая энергия переносится из области генерации в мелкие масштабы, где она диссипируется. А в двумерных течениях возможно появление двух инерционных интервалов, по которым реализуются прямой каскад переноса энстрофии (к" для энергетического спектра скорости), обеспечивающий диссипацию и обратный (красный) каскад кинетической энергии со степенным законом к"5/3, перекачивающий кинетическую энергию из масштаба генерации в область больших масштабов [74].

Обратный каскад энергии можно рассматривать как процесс самоорганизации турбулентности (рождение порядка из хаоса), в результате чего из поля мелкомасштабных пульсаций формируются крупномасштабные когерентные структуры. Такой каскадный процесс наблюдается в двумерных, вращающихся течениях, в плазме и волнах на поверхности жидкости [75]. Очевидна важная роль красного каскада энергии для течений в океане и атмосфере. В самом деле, наряду с вихрями сравнительно небольшого масштаба (порядка 1 км), в океане существуют вихри огромного размера (до 103 км), движения которых квазид-вумерны [61]. Существование таких огромных вихревых образований не представляется реальным без подпитки их энергией из более мелких масштабов.

Каскадные процессы в несжимаемой вязкой жидкости на основе двумерных уравнений Навье-Стокса исследовались в [76-83]. Во всех известных автору работах рассматривалась стационарная однородная двумерная турбулентность, расчеты проводились в квадратной области с периодическими граничными условиями и, как правило, с введением дополнительных членов, обеспечивающих стоки энергии на малых и больших масштабах. Чтобы получить стационарный в среднем процесс, в правую часть вводилась внешняя сила, осуществляющая подкачку энергии в виде белого шума.

Красный каскад энергии в электропроводящей жидкости (водный раствор А1аСГ) исследовался экспериментально в [75,82,84-86]. При этом течение создавалось путем пропускания через жидкость электрического тока, а роль внешней силы играли расположенные под слоем жидкости постоянные магниты.

В работах [76-78,81,82] наблюдался в установившемся (в среднем) течении красный каскад энергии со степенным законом к"5/3. Однако, в [79,80,83] этот степенной закон наблюдался только на начальной стадии расчета, а при выходе на стационарный режим устанавливался спектр к" . Подобную перестройку спектра автор [80] объясняет рождением самоподобных когерентных структур, 5 появление которых и обуславливает спектр к" на больших масштабах.

Эксперименты [75], проведенные в лотке размером 0.18м-0.18м также показали наличие четко идентифицируемого красного каскада энергии со степенси ным законом к" , а на приводимых в работе трасерных фотографиях видно формирование крупномасштабной вихревой структуры.

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

В работе [82] делается попытка совместного численного и экспериментального исследования двумерной турбулентности, при этом в численном моделировании красный каскад виден более определенно.

В эксперименте [86] подкачка энергии создавалась на крупных масштабах большим числом (около 400) постоянных магнитов, расположенных так, чтобы течению не навязывалась какая-либо определенная пространственная структура и электрическим током, меняющим сложным образом свое направление во времени. Такая крупномасштабная генерация турбулентности обусловила наличие продолжительного и четко идентифицируемого спектра энстрофийного каскада к".

В численных исследованиях каскад энстрофии наблюдается менее устойчиво, чем красный [79]. Степенной закон энстрофийного каскада был исследован в [77,79], но рассчитанные в этих работах показатели степенного закона сильно отличаются от -3 и показывают более крутой наклон, а именно: от -5 до -3.5 в [77] и круче -5 - в [79]. Отмеченное расхождение объясняется наличием когерентных вихревых структур, очень устойчивых и имеющих примерно одинаковый размер [74].

Интересное численное исследование проведено также в трехмерной вращающейся турбулентности [87]. Установлено, что при достаточно сильном вращении начинается отток подкачиваемой с помощью внешней силы энергии к большим масштабам, и течение становится квазидвумерным. В зависимости от величины аспектного отношения области (диаметр, отнесенный к высоте)

С/'З -1 наблюдался спектр к" при больших и к" - малых ее значениях. Спектр к' получен также по размерности из предположения, что скорость вращения является определяющим параметром.

В работе [88] описаны результаты совместного численного и экспериментального исследования турбулентности в вертикально стекающей под действием силы тяжести мыльной пленке, где турбулентность генерируется горизонтально и вертикально расположенными цилиндрами. Несмотря на ощутимый разброс в полученном численном и экспериментальном спектрах, автор выделяет инерционные интервалы, соответствующие прямому и обратному каскадам энстрофии и энергии.

Роль и проявление каскадных процессов в конвективной турбулентности пока не исследованы. А между тем, именно различие в каскадных процессах обуславливает качественное различие между двумерной и трехмерной конвективной турбулентностью при высокой надкритичности.

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

7 3

В самом деле, в ряде экспериментов при 11а ~ 10 (г ~ 7-10 ) наблюдалось изменение степенного закона от N11 ~ г1/3 до № ~ г277, что трактуется как переход от режима "мягкой" к режиму 'жесткой' турбулентности [89]. Подобная перестройка наблюдалась и в расчете [21].

С другой стороны, в других экспериментальных [40,41] и численном [5] исследованиях при увеличении надкритичности зафиксирован выход на асим

1 л птотический степенной закон Nu - г , означающий независимость теплообмена (числа Нуссельта) от толщины слоя [89].

И более того, в трехмерных расчетах при высокой надкритичности [4] отмечена обратная последовательность степенных законов - при г ~ 3-Ю5 "жесткая" турбулентность сменяется "мягкой". Отметим также, что другие экспериментальные работы, обзор которых можно найти в [41], содержат и

П 11 другие степенные законы, близкие к Nu ~ г ' .

Интересный и поучительный пример неоднозначности степенных законов конвекции дает совместный анализ трех экспериментальных работ [43,55,57]. Все эти исследования проведены с использованием в качестве рабочей жидкости газообразного гелия при криогенной температуре, в области цилиндрической формы с аспектным отношением 0.5. Однако, полученные в экспериментах [43,55] законы изменения числа Нуссельта как функции надкритичности отличаются качественно от результатов [57].

В самом деле, измеренные в работе [55] числа Нуссельта вплоть до Ra = 1014 следуют закону Nu = 0.17-Ra0'29, близкому к закону "жесткой" турбулентности. В работе [43] при 10б < Ra < 1017 получен закон Nu = 0.124-Ra0'309. А согласно данным [57], степенной закон "жесткой" турбулентности трансформируется в Nu = 0.0225-Ra0375 при 2-Ю11 < Ra < 2-Ю14, что вступает в противоречие с результатами [43,55]. Указанное несоответствие результатов измерений еще не получило достаточно убедительного объяснения и должно стать объектом дальнейшего исследования.

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

Заключение

В двумерной и трехмерной постановках рассмотрена задача о конвекции несжимаемой жидкости между двумя горизонтальными изотермическими плоскостями при подогреве снизу. При трехмерном моделировании горизонтальные границы предполагаются свободными от касательных напряжений, а в двумерном - свободными либо жесткими. За исключением нескольких тестовых расчетов, в которых проводятся сравнения с данными других авторов или экспериментом, число Прандтля равно 10.

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

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

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

В случае конвекции со свободными горизонтальными границами при Рг = 1 выписаны явные выражения для инкрементов нарастания гармоник в линейных аналогах дифференциальной задачи и предлагаемого численного метода. Показано, что коэффициенты при т и Н2, в членах, описывающих схемный эффект, не зависят от числа Рэлея, а определяются волновыми числами.

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

177

Показано, что несмотря на то, что численный метод строится как метод первого порядка по времени и второго по пространству, схемный эффект при линеаризации исходной системы уравнений есть величина порядка ОСт2+Н2) и, следовательно, спектральные характеристики численного метода и исходной дифференциальной системы уравнений близки, что гарантирует правильную динамику бесконечномалых возмущений.

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

В качестве теста, результаты расчета предлагаемым спектрально-разностным методом стационарной 2d,free конвекции сравниваются с результатами других авторов [15,32,33] при надкритичности г < 10 и по более простой модели бесконечного числа Прандтля [34] - при г < 2*10 . Отмечено хорошее совпадение по числу Нуссельта.

Результаты тестового расчета стохастической двухдиффузионной конвекции при RT = Rs = 1.5-104 хорошо согласуются с данными расчетов [30,31], в частности, соответствующие проекции решений на плоскость амплитуд первых гармоник функции тока и температуры представляются совпадающими с графической точностью. По методике, близкой к предложенной в [30,31], в этой серии тестовых расчетов контролировалась невязка, ее относительное значение не превосходило 1%.

Результаты экспериментов по стационарной, валиковой конвекции [45] при надкритичности (1 < г < 4) сопоставляются с данными 3d,free, 2d,free и 2d,rigid расчетов. Во всех расчетах отмечено хорошее согласование профиля средней температуры и изотермы полной температуры и, дополнительно, результаты 2d,rigid расчета хорошо согласуются с данными эксперимента по числу Нуссельта и удовлетворительно - по зависимости доминирующей длины волны от надкритичности.

Проведением при г = 6-10 на разных ПК двух 2d,free расчетов, с полностью совпадающими компьютерными кодами, трансляторами, числом учитываемых гармоник, начальными и граничными условиями, показано, что при сопоставлении параметров сложных конвективных течений имеет смысл, в основном, сопоставление средних (интегральных) характеристик, которые вычисляются достаточно устойчиво. А анализ мгновенных величин можно проводить только при достаточно высокой точности вычислений и на малом отрезке времени. В рассмотренном примере, при проведении вычислений на разных ПК с двойной (Real*8) точностью все анализируемые интегральные величины и профили практически совпадали между собой, а вычисленные мгновенные значения числа Нуссельта существенно отличались при t > 0.18. Качественно аналогичная ситуация наблюдалась и при проведении вычислений с точностью Real* 16. л

Два тестовых расчета, проведенные при г = 6-10 с учетом [129x33] гармоник и точностью вычислений Real* 8 и Real* 16, также показали практическое совпадение интегральных характеристик. л

Была проведена серия методических расчетов при г = 6-10 с учетом [129x33], [257x65], [513x129], [1025x257] и [2049x513] гармоник, анализ результатов которой показал, что наблюдается сходимость:

• всех анализируемых интегральных величин,

• профилей средней температуры, среднеквадратичных пульсаций температуры, вертикальной и горизонтальной скорости,

• одномерных пространственных спектров температуры и скорости.

И дополнительно, сравнение временных спектров числа Нуссельта, полученных с учетом [513x129] и [1025x257] гармоник, показало их качественное подобие с совпадением положений максимумов и законов затухания на высоких частотах.

Интегральные характеристики, полученные другим методом - модифицированным методом Галеркина (псевдоспектральным), основанным на вычислении нелинейных членов в физическом пространстве на разностной сетке, имеющей в каждом направлении в два раза больше точек, чем гармоник в том же направлении и схеме Рунге-Кутта четвертого порядка для интегрирования по времени, хорошо согласуются с результатами расчетов предлагаемым спектрально-разностным методом стохастических течений двухдиффузионной конвекции и конвекции при надкритичности до г = 2-10 .

Для иллюстрации возможности проводить 2d,rigid расчеты стохастических течений предлагаемым методом при высокой надкритичности, сравнивались л числа Нуссельта, полученные в настоящей работе и [21] при г < 6-10 . Среднее отклонение значений было равно 2.8, а максимальное - 3.9% при г = 6-10 . И более того, наши расчеты подтвердили отмеченную в [21] смену степенного закона с 1/3 на 2/7 при г ~ 3-10 в зависимости Nu от г.

Проведено также сравнение временного энергетического спектра (квадрата модуля Фурье преобразования функции времени) пульсаций температуры в центре конвективной ячейки, полученного в 3d,free (при г = 410) и 2d,free (г = 6.4-10 ) расчетах с данными эксперимента по турбулентной конвекции газообразного гелия при криогенной температуре [42]. Надкритичность и число Прандтля в расчетах и эксперименте совпадали. Спектры характеризуются большой консервативностью со слабой зависимостью от надкритичности, временного и пространственного разрешения, а также размерности. Размещение датчика температуры в центре конвективной ячейки в эксперименте обуславливает сравнительно слабую зависимость результатов измерений от геометрии области и граничных условий.

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

Показано, что кардинальное различие 2с1,&ее и ЗсЦгее решений - принципиально разное поведение вихревого масштаба при увеличении надкритично-сти, растущего до значения порядка размера области в двумерной конвекции и уменьшающегося примерно обратно пропорционально корню двенадцатой степени из надкритичности - в трехмерной. Указанная тенденция роста вихревого масштаба в двумерной конвекции проявляется с появлением нестационарности, вихревой масштаб достигает своего минимального значения при г = гт, гт ~ 36 (при а = 1), а затем увеличивается с ростом надкритичности до значения порядка размера области. Незначительное увеличение вихревого масштаба (до ~ 0.9%) наблюдалось и в расчетах стационарных одно- и двухвихревых течений.

В работе исследовались характеристики двумерной 2сЦгее конвекции в областях относительно большой С = 4тг и умеренной С = л горизонтальной протяженности и трехмерной ЗсЦгее конвекции при горизонтальной протяженности области равной я в обеих горизонтальных направлениях.

Различие в характеристиках двумерных 2(1,&ее течений при конвекции в областях умеренной и большой горизонтальной протяженности обусловлено:

• различным вихревым масштабом, по порядку величины равным размеру области,

• увеличением роли силы плавучести, наиболее существенной на больших масштабах,

• ростом отношения масштабов плавучести и генерации и

• увеличением генерации, переноса и диссипации кинетической энергии.

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

В двумерных 2сЦгее и трехмерных ЗсЦгее расчетах вычислялись и анализировались временные и пространственные спектры температуры и скорости. Наблюдаемая в расчетах динамика пространственных спектров укладывается в рамки предложенного в работе качественного физического сценария.

Согласно этому сценарию, при умеренной горизонтальной протяженности области в двумерных и трехмерных расчетах, действие силы плавучести, наиболее существенное на больших масштабах, должно обусловить стратификационный спектр для скорости. В то же время, действие силы плавучести на температуру еще незначительно, что с учетом доминирования конвективного переноса для температуры должно приводить к спектру пассивной примеси к"5/3.

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

С другой стороны, формирование красного каскада энергии при дальнейшем увеличении надкритичности приводит к перекачке энергии пульсаций скорости в большие масштабы, поле скорости становится крупномасштабным при относительно низком уровне пульсаций, и возникает своеобразный вязко-конвективный интервал, где пульсации температуры управляются крупномасштабным полем скорости, в результате чего в спектре температуры вместо стратификационного спектра должен появиться спектр Бэтчелора к"1 [74].

Заметим так же, что после формирования обратного каскада энергии его энергетическая роль неизбежно должна падать, что обусловлено отсутствием диссипации на больших масштабах. В свою очередь, диссипация энергии растет с увеличением надкритичности и горизонтальной протяженности области, обуславливая этим усиление прямого каскада энстрофии и приводя к установлению единого степенного закона в спектре скорости, близкого к закону энст-рофийного каскада к'3.

Описанный физический сценарий подтверждается тем, что при двумерной 2сЦгее и трехмерной ЗсЦгее конвекции в области умеренной горизонтальной протяженности для пульсаций температуры получены спектры к"5/3 и к"2'4. Эти спектры наблюдались в экспериментах по турбулентной конвекции. Ясного физического и теоретического обоснования для второго из перечисленных спектров до сих пор не получено, первый спектр указывает на поведение температуры как пассивной примеси, без доминирующей роли плавучести.

А для пульсаций скорости в расчетах получены стратификационные спектры Болджиано-Обухова (БО) к"11/5 и к"5 и, дополнительно, в трехмерной постановке стратификационный спектр Ламли-Шура (ЛШ) к"3. Все перечисленные здесь спектры имеют стратификационную природу, а спектры БО и ЛШ наблюдались в экспериментах по лабораторной конвекции и атмосфере.

С/-5

Однако, спектр Колмогорова к" , наблюдавшийся для пульсаций скорости в некоторых численных исследованиях с различной степенью убедительности, в трехмерных расчетах настоящей работы получен не был, также как и в немногих известных экспериментальных исследованиях турбулентной конвекции, в которых исследовались спектры пульсаций скорости. Отсутствие спектра Колмогорова в наших трехмерных расчетах объясняется недостаточно высоким значением числа Рейнольдса. Его значение, вычисленное по характерной скорости и высоте слоя не превосходит 44, а значение внутреннего числа Рейнольдса, определенного по характерной скорости и вихревому масштабу - 21.4, в то время как начало развития турбулентности связывается со значением 15.

Подчеркнем, что практически полное совпадение спектров температуры и скорости в 2с1,£гее и ЗсЦгее расчетах в области умеренной горизонтальной протяженности £ = я наблюдается во всем рассмотренном диапазоне изменения надкритичности, при г <3.4-104 в двумерных и г < 950 - трехмерных расчетах.

С описанным выше физическим сценарием согласуется и перестройка спектров температуры и скорости при 2с1,1тее конвекции в области относительно большой горизонтальной протяженности £ = 4п.

При 500 < г < 4-10 на длинноволновом участке спектра скорости наблюдаются одновременно два разных спектра, соответствующие двум конкурирующим механизмам - плавучести и обратному каскаду энергии. Это сосуществование двух различных спектров (к"5/3 и к"11/5) обуславливает их нечеткую идентификацию и в некотором смысле условное проведение верхней границы г = 4-10 этого сосуществования по надкритичности.

При значении надкритичности 4-103 < г < 104 в спектре скорости одновременно наблюдаются два инерционных интервала, соответствующие обратному (красному) каскаду энергии со степенным законом к"5/3, перекачивающему кинетическую энергию из масштаба генерации в крупные масштабы и прямому каскаду энстрофии со степенным законом к", обеспечивающему диссипацию. А при г > 104 в спектре скорости видна тенденция к установлению единого степенного закона к"2'6, близкого к закону энстрофийного каскада.

Перестройка спектра скорости при г ~ 104 подобна отмеченной в численных [79,80,83] и экспериментальном [84] исследованиях двумерной турбулентности и обусловлена формированием крупномасштабной структуры течения. В самом деле, относительное ослабление энергетической роли обратного каскада переноса энергии к большим масштабам (крупномасштабная структура поля скорости уже сформирована, а диссипация энергии пренебрежимо мала на больших масштабах) и усиление прямого каскада энстрофии из-за роста диссипации обуславливает тенденцию установления единого степенного закона к"2'6, близкого к закону энстрофийного каскада. Последнее находится в некотором противоречии с результатами численного [77] и экспериментального [84] исследований двумерной турбулентности, где спектр в интервале переноса энст-рофии заметно круче. Это различие, предположительно, связано с действием присутствующей в конвекции силы плавучести, в результате чего степенной закон энстрофийного каскада заменяется более пологим, но простирающимся в область волновых чисел, меньших частоты генерации.

В длинноволновом спектре температуры при умеренной надкритичности

3 7/5

750 < г < 2-10 наблюдается стратификационный спектр БО к" , который при

Л -3 1 о

2-10 < г <6-10 заменяется на спектр Бэтчелора k" . А при г > 6-10 показатель степени (3 степенного закона к"р принимает значения из интервала [0.7,1], с наиболее вероятным значением 0.8.

Подобная перестройка наблюдалась экспериментально в физически похожей задаче о конвекции в вертикальной мыльной пленке при подогреве снизу [109]. В спектре температуры при dQ < 48° К там наблюдался стратификаци

7/С 1 онный спектр БО k" , а при dQ > 48° К - Бэтчелора к". Авторы [109] связывают перестройку температурного спектра с формированием крупномасштабной структуры течения.

Конечно, перестройка температурного спектра при г ~ 2-103 (спектра скорости при г ~

4-100 связана с формированием крупномасштабной структуры, течения (конденсацией), что обусловленно действием обратного каскада энергии. При меньшей надкритичности красный каскад энергии еще не сформирован и силы плавучести, существенные на больших масштабах, обуславливают стратификационный спектр БО для температуры к"7/5 (к"11/5 для скорости). Но, с повышением надкритичности энергетическая роль красного каскада усиливается и благодаря его действию происходит перекачка энергии пульсаций скорости в большие масштабы, поле скорости становится крупномасштабным при относительно низком уровне пульсаций и возникает своеобразный вязкокон-вективный интервал, где пульсации температуры управляются крупномасштабным полем скорости, что и обуславливает появление спектра Бэтчелора к"1. Понижение энергетической роли обратного каскада энергии при г > 104 обуславливает тенденцию к установлению единого степенного закона в спектре скорости и колебания показателя длинноволнового спектра пульсаций температуры (при г < 6-Ю3).

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

Существенное отличие конвекции со свободными границами от конвекции с жесткими - возможность неравенства нулю значений вертикальной компоненты завихренности, что обуславливает возможное вращение вокруг вертикальной оси, а направленная горизонтально скорость такого вращения и отличие в граничных условиях для нее, обуславливает возможное расхождение в характеристиках, связанных с горизонтальной скоростью. Заметим, что при двумерном моделировании вертикальная компонента завихренности равна нулю тождественно и согласование среднеквадратичных пульсаций горизонтальной скорости при двумерной конвекции со свободными и жесткими граничными условиями лучше, чем вертикальной. Результаты 2с1,&ее расчетов, при этом, хорошо согласуются с экспериментальными данными по числу Нуссельта даже при высокой надкритичности и удовлетворительно - по величине среднеквадратичных пульсаций горизонтальной скорости [63].

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

Уже отмечалось выше, что кардинальное отличие двумерных течений от трехмерных заключается в различном поведении вихревого масштаба при увеличении надкритичности, который в двумерной конвекции растет, а в трехмерной - уменьшается. Указанная тенденция роста вихревого масштаба в двумерной конвекции проявляется при г > гт, где гт ~ 36 (при а = 1). Естественно ожидать, что и это различие в поведении вихревого масштаба становится доминирующим не сразу после его появления. Таким образом, можно ожидать примерного соответствия степенных законов изменения средних величин при увеличении надкритичности в двумерной и трехмерной конвекции примерно до значений порядка 100. Причем результаты вычисления и сравнительного анализа пространственных и временных спектров в двумерных и трехмерных расчетах при а = 1 показывают, что для консервативных величин, слабо зависящих от надкритичности и пространственного разрешения, интервал соответствия может быть много больше.

Сказанное подтверждается тем, что величина среднеквадратичных пульсаций температуры как функция надкритичности в Зё/гее расчетах с хорошей точностью соответствует степенному закону

Я* = 0.163 т2/15 при г> 150, показатель которого совпадает с полученным в эксперименте [53] и близок к расчетному [1].

А результаты 2сЦгее расчетов с хорошей точностью соответствуют степенному закону (при 170 < г < 305 рассматривается верхняя ветвь решения) = 0.316т02 при 40 < г < 305, показатель которого совпадает с полученным в экспериментах по конвекции газообразного Не при криогенной температуре при более высокой надкритичности [55].

Отметим, что в двумерных расчетах (при г ~ 103 в 2d,free и г ~ 250 -2d,rigid) прослеживается качественное изменение зависимости величины среднеквадратичных температурных пульсаций от надкритичности - тенденция убывания сменяется возрастанием. Появление нефизичной тенденции возрастания температурных пульсаций обусловлено формированием крупномасштабной структуры течения.

А величина среднеквадратичных пульсаций вертикальной скорости в 3d,free расчетах с хорошей точностью соответствуют степенному закону

W' = 15.8т0-395 при г >150, с показателем, близким к полученному в экспериментах [48,53] и расчете [1].

В двумерных расчетах соответствие наблюдается при сравнительно невысокой надкритичности (г < 250 в 2d,rigid и г < 305 - в 2d,free расчетах).

А именно, результаты 2d,rigid расчетов и данные эксперимента [48] показывают при надкритичности г < 250 практически совпадающие значения и показатели степенных законов. А степенные законы, аппроксимирующие результаты 3d,free и 2d,free расчетов имеют практически совпадающие при надкритичности г < 137 показатели, несколько завышенные по сравнению с полученными в эксперименте [48] (0.44) и численном расчете [1] (0.46), но значения среднеквадратичных пульсаций вертикальной скорости, полученные в 2d,free расчете при г < 305 согласуются с данными [1].

Отметим так же, что результаты 3d,free и 2d,free расчетов практически совпадают в диапазоне надкритичности 280 < г < 850 (при 170 < г < 305 рассматривается нижняя ветвь решения).

При достаточно большой надкритичности, двумерные расчеты показывают нефизично высокий рост среднеквадратичных пульсаций вертикальной скорости: W ~ г0'57 и W' ~ г0'715 - в 2d,rigid и 2d,free расчетах, соответственно. Это, как и в случае пульсаций температуры, обусловлено формированием крупномасштабной структуры течения.

Различие в пульсациях горизонтальной скорости при двумерном моделировании выражено слабее, чем в пульсациях вертикальной и данные 2d,free расчета удовлетворительно согласуются с экспериментом [63]. Результаты трехмерных и двумерных расчетов близки при г < 250.

Значение числа Nu в 2d,free расчетах хорошо согласуется с экспериментом при большой надкритичности. Причем, при г < 800 значения числа Нуссельта в 2d,rigid расчетах меньше и ближе к экспериментальной кривой, чем в 2d,free, а при г > 800 - наоборот.

В 2d,free и 2d,rigid расчетах настоящей работы при достаточно большой надкритичности наблюдается выход на асимптотический степенной закон, П близкий к Nu ~ г ' . Показатель этого степенного закона близок к показателю степенного закона Nu ~ г0'356, полученного для стационарной, одновихревой конвекции.

Расчет интегральных характеристик стационарного, одновихревого течения в постановке [33], но в более широком диапазоне надкритичности г < 104, показал близость показателей степенных законов для основных интегральных величин в 2d,free расчетах и стационарном, крупномасштабном одно-вихревом течении, при практически полном их совпадении для кинетической энергии и среднеквадратичной скорости.

А в зависимости числа Нуссельта от надкритичности в 3d,free расчетах можно выделить два участка со степенными зависимостями:

Nu = 1.85 т0'302 при 20 < г <137 и Nu = 2.79 т0'230 при 150 < г <950.

Отметим, что при г < 137 показатель закона нарастания числа Нуссельта хорошо согласуется со значением 0.305 (эмпирический закон O'Toole и Silveston, 1961) и с полученным в 2d,free расчетах при г < 305.

В общем, последовательность степенных законов для числа Нуссельта в 3d,free расчетах соответствует полученной в 2d,free, но смена степенных законов в трехмерной конвекции происходит при меньшей надкритичности, чем в двумерной. Но числа Нуссельта, полученные в 2d,free расчетах, хорошо согласуются с экспериментом даже при высокой надкритичности [63]. Это и логика

189 сопоставления степенных законов для числа Нуссельта в 2<3,1тее и ЗсЦгее расчетах указывает на необходимость проведения дальнейших расчетов трехмерной конвекции при более высокой надкритичности.

Суммируя все вышесказанное, сформулируем основные выводы работы:

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

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

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

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

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

6. Показано, что результаты 2d,rigid расчета стационарной, валиковой конвекции при невысокой надкритичности 1 < г < 4 с удовлетворительной точностью отражают зависимость доминирующей длины волны от надкритичности, хорошо согласуясь с данными эксперимента по числу Нуссельта, профилю средней температуры и изотерме полной температуры.

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

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

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

10. Несмотря на наблюдаемые количественные расхождения между результатами 3d,free расчетов и экспериментом, трехмерные расчеты дают правильные показатели степенных законов изменения среднеквадратичных пульсаций температуры, вертикальной скорости, числа Рейнольдса, среднеквадратичной скорости и кинетической энергии от надкритичности при г > 150. В двумерных расчетах аналогичное соответствии наблюдается при сравнительно невысокой надкритичности (до г ~ 250) по среднеквадратичным пульсациям температуры (в 2d,free) и вертикальной скорости (в 2d,free и 2d,rigid расчетах).

11. Закон нарастания при высокой надкритичности числа Нуссельта, нефи-зичное поведение пульсаций вертикальной скорости и температуры, а также

191 близость показателей степенных законов для числа Нуссельта, кинетической энергии, среднеквадратичной скорости и энстрофии в 2сЦгее расчетах и стационарном, одновихревом течении показывает, что в 2с1,£гее расчетах при высокой надкритичности формирование крупномасштабной структуры является доминирующим фактором, определяющим характеристики течения.

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

В области умеренной горизонтальной протяженности л в 2сЦгее и ЗсЦгее расчетах:

13. Действие силы плавучести обуславливает стратификационные спектры для скорости, а доминирование конвективного переноса для температуры -к"24 и спектр пассивной примеси к"5/3. Причем, эти спектры наблюдались в двумерных и трехмерных расчетах во всем рассматриваемом диапазоне надкритичности.

В области большой 4л горизонтальной протяженности в 2сЦгее расчетах:

14. Показано, что при 500 < г < 4-103 длинноволновый участок спектра скорости состоит из двух ветвей, со степенными законами, отражающими конкурирующие механизмы, связанные с силой плавучести и каскадным процессом переноса энергии, а при 4-103 < г < 104 в спектре скорости

4 СИ идентифицируются степенные законы к" и к" , соответствующие прямому и обратному каскадным процессам переноса энстрофии и энергии. При дальнейшем увеличении надкритичности (г > 104) в спектре скорости видна тенденция к установлению единого степенного закона к"2'6 - искаженного закона каскада энстрофии.

15. Показано, что при умеренной надкритичности 500 < г < 2-10 в спектре П температуры виден стратификационный спектр БО к" , сменяющийся при

О 1 о

2-10 < г <6-10 спектромБэтчелора k" . А при г > 6-10 показатель длинноволнового участка температурного спектра совершает колебания в интервале [-1,-0.7], с наиболее вероятным значением -0.8.

16. Перестройки спектров температуры (скорости) в двумерных расчетах связаны с формированием обратного каскада энергии при г ~ 2-103 (для скорости при г ~ 4-10 ) и уменьшением его энергетической роли при г~6-103(г~104).

1. Kerr R.M. Rayleigh number scaling in numerical convection // J. Fluid Mech. 1996. V.310. P.139-179.

2. Hartlep T., Tilgner A., Busse F.H. Large scale structures in Rayleigh-Benard convection at high Rayleigh numbers // Phys. Rev. Lett. 2003. V.91. N.6. P.064501-1 4.

3. Yang H., Zhu Z. Numerical simulation of turbulent Rayleigh-Benard convection // Int. Com. Heat Mass Transfer. 2006. V.33. P.184-190.

4. Verzicco R., Camussi R. Numerical experiments on strongly turbulent thermal convection in a slender cylindrical cell // J. Fluid Mech. 2003. V.477. P. 19-49.

5. Amati G., Koal K., Massaioli F., Sreenivasan K.R., Verzicco R. Turbulent thermal convection at Rayleigh numbers for a Boussinesq fluid of constant Prandtl number // Phys. Fluids. 2005. V.17. P.121701-1 4.

6. Shishkina O., Wagner C. Analysis of thermal dissipation rates in turbulent Rayleigh-Benard convection // J. Fluid Mech. 2006. V.546. P.51-60.

7. Van Reeuwijk M., Jonker H.J.J., Hanjalic K. Identification of the wind in Rayleigh-Benard convection // Phys. Fluids. 2005. V.17. N.4. P.051704-1 4.

8. Grotzbach G. Direct numerical simulation of laminar and turbulent Benard convection // J. Fluid Mech. 1982. V.l 19. P.27-53.

9. Kerr R.M., Herring J.R., Brandenburg A. Large-scale structure in Rayleigh-Benard convection with impenetrable sidewalls // Chaos, Solitons & Fractals. 1995. V.5. N.10. P.2047-2053.

10. Arter W. Nonlinear Rayleigh-Benard convection with square planform // J. Fluid Mech. 1985. V.152. P.391-418.

11. Travis B., Olson P., Schubert G. The transition from two-dimensional to three-dimensional planforms in infinite-Prandtl-number thermal convection // J. Fluid Mech. 1990. V.216. P.71-91.r *

12. Malevsky A.V. Spline-characteristic method for simulation of convective turbulence // J. Comput. Phys. 1996. V.123. N.2. P.466-475.

13. Balachandar S., Maxey M.R., Sirovich L. Numerical simulation of high Rayleigh number convection // J. Sci. Сотр. 1989. V.4. N.2. P.219-236.

14. Cortese Т., Balachandar S. Vortical nature of thermal plumes in turbulent convection // Phys. Fluids. A. 1993. V.5. N.12. P.3226-3232.

15. Thual O. Zero-Prandtl-number convection // J. Fluid Mech. 1992. V.240. P.229-258.

16. Герценштейн С.Я., Родичев Е.Б., Шмидт B.M. Взаимодействие трехмерных волн во вращающемся горизонтальном слое жидкости, подогреваемом снизу // Докл. АН СССР. 1978. Т.238. N.3. С.545-548.

17. Curry J.H., Herring J.R., Loncaric J., Orszag S.A. Order and disorder in two-and three-dimensional Benard convection // J. Fluid Mech. 1984. V.147. P.l-38.

18. Breuer M., Wessling S., Schmalzl J., Hansen U. Effect of inertia in Rayleigh-Benard convection // Phys. Rev. E 69. 2004. P.026302-1 10.

19. Schmalzl J., Breuer M., Hansen U. The influence of the Prandtl number on the style of vigorous thermal convection // Geophys. Astrophys. Fluid Dynamics. 2002. V.96. N.5. P.381-403.

20. Goldhirsch I., Pelz R.B., Orszag S.A. Numerical simulation of thermal convection in a two-dimensional finite box // J. Fluid Mech. 1989. V.199. P. 1-28.

21. DeLuca E.E., Werne J., Rosner R., Cattaneo F. Numerical simulation of soft and hard turbulence: preliminary results for two-dimensional convection // Phys. Rev. Letters. 1990. V.64. N.20. P.2370-2373.

22. Werne J. Structure of hard-turbulent convection in two dimensions: Numerical evidence // Phys. Rev. E. 1993. V.48. N.2. P.1020-1035.

23. Liu J.-G., Wang C., Johnston H. A fourth order scheme for incompressible Bous-sinesq equations //J. of Sci. Сотр. 2003. V.18. N.2. P.253-285.

24. Schneck P., Veronis G. Comparison of some recent experimental and numerical results in Benard convection // Phys. Fluids. 1967. V.10. N.5. P.927-930.

25. Clever R.M., Busse F.H. Transition to time-dependent convection // J. Fluid Mech. 1974. V.65. Pt.4. P.625-645.

26. Полежаев В.И., Яремчук В.П. Численное моделирование двумерной нестационарной конвекции в горизонтальном слое, подогреваемом снизу // Изв. АН СССР. МЖГ. 2001. N.4. С.34-45.

27. Бабенко К.И., Рахманов А.И. Численное исследование двумерной конвекции. Препринт / ИПМ АН СССР им. М.В. Келдыша, N.118.1988. 26 с.

28. Родичева О.В., Родичев Е.Б. О двумерной турбулентности в задаче Рэлея-Бенара // Докл. АН СССР. 1998. Т.359. N.4. С.486-489.

29. Герценштейн С.Я., Шмидт В.М. Нелинейное взаимодействие конвективных волновых движений и возникновение турбулентности во вращающемся горизонтальном слое // Изв. АН СССР. МЖГ. 1977. N.2. С.9-15.

30. Gertsenstein S., Sibgatullin I. Bifurcations, Transition to turbulence and development of chaotic regimes for double-diffusive convection // Wseas Transactions on Applied and Theoretical Mechanics. 2006. V. 1. Is. 1. P. 110-114.

31. Gertsenstein S., Sibgatullin I., Sibgatullin N. Some properties of two-dimensional stochastic regimes of double-diffusive convection in plane layer // Chaos. 2003. V.13. N.4. P.1231-1241.

32. Veronis G. Large-amplitude Benard convection // J. Fluid Mechanics. 1966. V.26. Pt.l. P.49-68.

33. Moore D.R., Weiss N.O. Two-dimensional Rayleigh-Benard convection // J. Fluid Mech. 1973. V.58. Part 2. P.289-312.

34. Schubert G., Anderson C.A. Finite element calculations of very high Rayleigh number thermal convection // Geophys. J. R. Astr. Soc. 1985. V.80. P.576-601.

35. Malevsky A.V., Yuen D.A. Characteristics-based methods applied to infinite Prandtl number thermal convection in the hard turbulent regime // Phys. Fluids. A. 1991. V.3. N.9. P.2105-2115.

36. Vincent A.P., Yuen D.A. Plumes and waves in two-dimensional turbulent convection // Phys. Rev. E. 1999. V.60. N.3. P.2957-2963.

37. Vincent A.P., Yuen D.A. Transition to turbulent thermal convection beyond Ra = 1010 detected in numerical simulations // Phys. Rev. E. 2000. V.61. N.5. P.5241-5246.

38. Hansen U., Yuen D.A., Malevsky A.V. Comparison of steady-state and strongly chaotic thermal convection at high Rayleigh number // Phys. Rev. A. 1992. V.46. N.8. P.4742-4754.

39. Goldstein R.J., Graham D.J. Stability of a horizontal fluid with zero shear boundaries//Phys. Fluids. 1969. V.12. N.6. P.1133-1137.

40. Niemela J.J., Sreenivasan K.R. Turbulent convection at high Rayleigh numbers and aspect ratio 4 // J. Fluid Mech. 2006. V.557. P.411-422.

41. Fleischer A.S., Goldstein R.J. High-Rayleigh-number convection of pressurized gases in a horizontal enclosure // J. Fluid Mech. 2002. V.469. P. 1-12.

42. Wu X.-Z., Kadanoff L., Libchaber A., Sano M. Frequency power spectrum of temperature fluctuations in free convection // Phys. Rev. Lett. 1990. V.64. N.18. P.2140-2143.

43. Niemela J.J., Skrbek L., Sreenivasan K.R., Donnelly R.J. Turbulent convection at very high Rayleigh numbers // Nature. 2000. V.404. N.20. P.837-840.

44. Ashkenazi S., Steinberg V. Spectra and statistics of velocity and temperature fluctuations in turbulent convection // Phys. Rev. Lett. 1999. V.83. N.23. P.4760-4763.

45. Farhadieh R., Tankin R.S. Interferometric study of two-dimensional Benard convection cells // J. Fluid Mech. 1974. V.66. Pt.4. P.739-752.

46. Shang X.-D., Xia K.-Q. Scaling of the velocity power spectra in turbulent thermal convection// Phys. Rev. E. 2001. V.64. P.065301-1 4.

47. Krishnamurti R., Howard L.N. Large-scale flow generation in turbulent convection // Proc. Natl. Acad. Sci. USA. (Applied physical and mathematical sciences). 1981. V.78. N.4. P.1981-1985.

48. Garon A.M., Goldstein R.J. Velocity and heat transfer measurements in thermal convection// Phys. Fluids. 1973. V.16. N.ll. P.1818-1825.

49. Cioni S., Ciliberto S., Sommeria J. Temperature structure functions in turbulent convection at low Prandtl number // Europhys. Lett. 1995. V.32. N.5. P.413-418.

50. Chu T.Y., Goldstein R.J. Turbulent convection in a horizontal layer of water // J. Fluid Mech. 1973. V.60. Pt.l. P.141-159.

51. Deardorff J.W., Willis G.E. Investigation of turbulent thermal convection between horizontal plates // J. Fluid Mech. 1967. V.28. Pt.4. P.675-704.

52. Thomas D.B., Townsend A.A. Turbulent convection over a heated horizontal surface // J. Fluid Mech. 1957. V.2. P.473-492.

53. Fitzjarrald D.E. An experimental study of turbulent convection in air // J. Fluid Mech. 1976. V.73. Pt.4. P.693-719.

54. Denton R.A., Wood I.R. Turbulent convection between two horizontal plates // Int. J. Heat Mass Transfer. 1979. V.22. N.10. P.1339-1346.

55. Wu X.-Z., Libchaber A. Scaling relations in thermal turbulence: the aspect-ratio dependence // Phys. Rev. A. 1992. V.45. N.2. P.842-845.

56. Malkus W.V.R. Discrete transitions in turbulent convection // Proc. Roy. Soc. London. A. 1954. V.225.N.1161. P.185-195.

57. Chavanne X., Chilla F., Chabaud B, Castaing В., Hebral B. Turbulent Rayleigh-Benard convection in gaseous and liquid He II Phys. Fluids. 2001. V.13. N.5. P. 1300-1320.

58. Палымский И.Б. Численное исследование спектров турбулентной конвекции Рэлея-Бенара // Нелинейная динамика. 2008. Т.4. N.2. С.145-156.

59. Гетлинг А.В. Конвекция Рэлея-Бенара. Структуры и динамика. М.: Эдиториал УРСС, 1999. 247 с.

60. Палымский И.Б. Численное исследование спектров трехмерной конвекции Рэлея-Бенара// Известия РАН. ФАО. 2009. Т.45. N.5. С.691-699.

61. Физика океана (ред. Доронина Ю. П.). JL: Гидрометеоиздат, 1978. 294 с.

62. Пономарев С. Г., Стойнов М.И. Решение уравнений Навье-Стокса проекционными методами; вычисление нелинейных членов. Препринт / ИПМ АН СССР им. М.В. Келдыша, N.58. 1987.18 с.

63. Палымский И.Б. Численное моделирование двумерной конвекции, роль граничных условий // Известия РАН. МЖГ. 2007. N.4. С.61-71.

64. Палымский И.Б. О качественном различии решений двумерной и трехмерной конвекции // Нелинейная динамика. 2009. Т.5. N.2. С.183-203.

65. Палымский И.Б. Численное моделирование двумерной конвекции при высокой надкритичности // Успехи механики. 2006. N.4. С.3-28.

66. Гершуни Г.З., Жуховицкий Е.М. Конвективная устойчивость несжимаемой жидкости. М.: Наука, 1972. 392 с.

67. Рождественский Б.Л., Яненко H.H. Системы квазилинейных уравнений и их приложения к газовой динамике. М.: Наука, 1978. 687 с.

68. Палымский И.Б. Линейный и нелинейный анализ численного метода расчета конвективных течений // Сиб. ж. вычисл. матем. 2004. Т.7. N.2. С.143-163.

69. Палымский И.Б. Метод численного моделирования конвективных течений // Вычисл. технологии. 2000. Т.5. N.6. С.53-61.

70. Никитин Н.В., Полежаев В.И. Трехмерные эффекты переходных и турбулентных режимов тепловой гравитационной конвекции в методе Чохральского //ИзвестияРАН. МЖГ. 1999. N.6. С.81-90.

71. Никитин Н.В. Спектрально-конечноразностный метод расчета турбулентных течений несжимаемых жидкостей в трубах и каналах // ЖВМ и МФ. 1994. Т.34. N.6. С.909-925.

72. Рождественский Б.Л., Стойнов М.И. Алгоритмы интегрирования уравнений Навье-Стокса, имеющие аналоги законам сохранения массы, импульса и энергии. Препринт / ИПМ АН СССР им. М.В. Келдыша, N.119.1987. 28 с.

73. Пасконов В.М., Полежаев В.И., Чудов Л.А. Численное моделирование процессов тепло- и массообмена. М.: Наука, 1984. 285 с.

74. Фрик П.Г. Турбулентность: подходы и модели. Москва-Ижевск: Институт компьютерных исследований, 2003. 292 с.

75. Shats M.G., Xia H., Punzmann H. and Falkovich G. Supression of turbulence by self-generated and imposed flows // Phys. Rev. Lett. 2007. V.99. N.16. P.164502-1 4.

76. Frish U., Sulem P.L. Numerical simulation of the inverse cascade in two-dimensional turbulence // Phys. Fluids. 1984. V.27. N.8. P.1921-1923.

77. Babiano A., Dubrulle B. Frick P. Scaling properties of numerical two-dimensional turbulence // Phys. Rev. E. 1995. V.52. N.4. P.3719-3729.

78. Boffetta G., Celani A. Vergassola M. Inverse energy cascade in two-dimensional turbulence: deviations from Gaussian behavior // Phys. Rev. E. 2000. V.61. N.l. P.R29-R32.

79. Tran C.V., Bowman J.C. Robustness of the inverse cascade in two-dimensional turbulence // Phys. Rev. E. 2004. V.69. P.036303-1-7.

80. Borue V. Inverse energy cascade in stationary two-dimensional homogeneous turbulence // Phys. Rev. Lett. 1994. V.72. N.10. P.1475-1478.

81. Smith L.M., Yakhot V. Bose condensation and small-scale structure generation in a random force driven 2D turbulence // Phys. Rev. Lett. 1993. V.71. N.3. P.352-355.

82. Chen S., Ecke R.E., Eyink G.L., Rivera M., Wan M., Xiao Z. Physical mechanism of the two-dimensional inverse energy cascade // Phys. Rev. Lett. 2006. V.96. N.8. P.084502-1 4.

83. Chertkov M., Connaughton C., Kolokolov I., Lebedev V. Dynamics of energy condensation in two-dimensional turbulence // Phys. Rev. Lett. 2007. V.99. N.8. P.084501-1 4.

84. Shats M.G., Xia H., Punzmann H. Spectral condensation of turbulence in plasmas and fluids and its role in low-to-high phase transitions in toroidal plasma // Phys. Rev. E. 2005. V.71. P.046409-1 9.

85. Paret J., Tabeling P. Experimental observation of the two-dimensional inverse cascade // Phys. Rev. Lett. 1997. V.79. N.21. P.4162-4165.

86. Paret J., Jullien M.-C., Tabeling P. Vorticity statistics in the two-dimensional enstrophy cascade//Phys. Rev. Lett. 1999. V.83. N.17. P.3418-3421.c

87. Bruneau С. H., Kellay H. Experiments and direct numerical simulations of two-dimensional turbulence // Phys. Rev. E. 2005. V.71. P.046305-1 5.

88. Фабер Т.Е. Гидроаэродинамика. M.: Постмаркет, 2001. 560 с.

89. Palymskiy I.B., Fomin P.A. Hieronymus H. Rayleigh-Benard convection in chemical equilibrium gas (simulation of surface detonation wave initiation) // App. Math. Model. 2008. V.32. N.5. P.660-676.

90. Palymskiy I.B., Fomin P.A. Hieronymus H. The Rayleigh-Benard convection in gas with chemical reactions // Сиб. ж. вычисл. матем. 2007. Т.10. N.4. С.371-383.

91. Палымский И.Б. О моделировании сложных режимов конвекции Рэлея-Бенара // Сиб. ж. вычисл. матем. 2011. Т.14. N.2. С.179-204.

92. Бахвалов Н. С. Численные методы, T.l. М.: Наука, 1975. 632 с.

93. Роуч П. Вычислительная гидродинамика. М.: Наука, 1980. 618 с.

94. Палымский И.Б. Численное исследование спектров трехмерной турбулентной конвекции // Известия Сарат. универ. Новая серия: матем., мех., информ. 2010. Т.10. В.1. С.62-71.

95. Палымский И.Б. Численный метод расчета трехмерной конвекции // Сиб. ж. индустр. матем. 2010. Т.13. N.l. С.95-108.

96. Палымский И.Б. О численном моделировании трехмерной конвекции // Вест. Удмурт, универ. серия 1: матем., мех., компьют. науки. 2009. В.4. С.118-132.

97. Гольдштик М.А., Штерн В.Н. Гидродинамическая устойчивость и турбулентность. Новосибирск: Наука СО, 1977. 366 с.

98. Палымский И.Б. Качественный анализ разностных схем для модельного нелинейного уравнения со знакопеременной вязкостью // Дифф. уравн. 1992. Т.28. N.12. С.2148-2158.

99. Палымский И.Б. Численное моделирование конвективных течений при высоких числах Вейссенберга. Препринт / ИГ11М СО АН СССР им. С.А. Христиановича, N.15. 1988. 48 с.

100. Senoner J.-M., Garcia M., Mendez S., Staffelbach G., Vermorel O. Growth of rounding errors and repetitively of large-eddy simulations // AIAA J. 2008. V.46. N.7. P.1773-1781.

101. Палымский И.Б., Герценштейн С.Я., Сибгатуллин И.Н. Об интенсивной турбулентной конвекции в горизонтальном плоском слое жидкости // Известия РАН. ФАО. 2008. Т.44. N.l. С.75-85.

102. Заславский Г.М., Сагдеев Р.З. Введение в нелинейную физику. От маятника до турбулентности и хаоса. М.: Наука, 1988. 378 с.

103. Белоцерковский О.М., Опарин A.M. Численный эксперимент в турбулентности. От порядка к хаосу. М.: Наука, 2000, 224 с.

104. Kraichnan R.H. Turbulent thermal convection at arbitrary Prandtl number // Phys. Fluids. 1962. V.5. N.l 1. P.1374-1389.

105. Ландау Л.Д., Лифшиц E.M. Гидродинамика. M.: Наука, 1988. 730 с.

106. Филлипс О.М. Динамика верхнего слоя океана. Л.: Гидрометеоиздат, 1980.319 с.

107. Kraichnan R.H. Inertial-range transfer in two- and three-dimensional turbulence // J. Fluid Mech. 1971. V.47. Pt.3. P.525-535.

108. Zhang J., Wu X.L. Density fluctuations in strongly stratified two-dimensional turbulence // Phys. Rev. Lett. 2005. V.94. N.17. P. 174503-1 4.

109. Xu X., Bajaj K.M.S., Ahlers G. Heat transport in turbulent Rayleigh-Benard convection // Phys. Rev. Lett. 2000. V.84. N.19. P.4357-4360.

110. Вычисление одномерных энергетических спектров

111. Одномерные пространственные энергетические спектры температуры и скорости вычислялись и анализировались в 2сЦгее и ЗсЦгее расчетах в работах 58,60,95. с использованием формул [74].

112. Опишем методику вычисления только спектров температуры, так как способ вычисления спектров скорости отличается лишь несущественными деталями.

113. Рассмотрим сначала вычисление температурных спектров в двумерном случае:

114. По известному полю температуры, осреднением по времени и горизонтальной координате вычисляется профиль средней температуры:1. Кг)=<<2>,здесь угловые скобки означают осреднение по времени и горизонтальной координате х.

115. Вычисляется поле пульсаций температуры в каждый момент времени:

116. Е(д.ы> & = 0,1,.,А72 и га = 0,1,.,М.

117. Проведя суммирование по одному из индексов полученного двумерного энергетического спектра, находим одномерные энергетические спектры:м к/гт=0 к=0

118. Вычисленные одномерные спектры записываются в виде: (2ак,ЕС±к) при к = 0,\,.,К/2 и {7гт,Е0^т) при т = 0,1,.,М,где в круглых скобках первым записано волновое число, а вторым соответствующее ему значение одномерного энергетического спектра.

119. Как и в двумерном случае, полученный одномерный спектр записывается в виде:2ак,Е(±к) при к = 0,1,.,К /2.

120. Аналогичным образом вычисляются одномерные спектры ЕОп и ЕОт в у и г направлениях, соответственно.

121. Для выполнения преобразования Фурье использовались стандартные программы быстрого преобразования Фурье.

122. Определение степенных законов

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

124. Временные спектры при сложных режимах конвекции осцилляторные и это обстоятельство существенно осложняет использование точных методов определения присущих им степенных законов.

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

126. Для определения степенных законов, присущих пространственным спектрам и зависимостям интегральных характеристик от надкритичности, используются две методики.

127. Первая состоит в нахождении параметров степенного закона методом наименьших квадратов с помощью стандартной программы СигуеЕхрег! 1.3 или любой другой аналогичной.

128. Второй метод требует построения компенсационного спектра 43.

129. В качестве примера определения степенного закона рассмотрим спектр:

130. Е(х) = х~2(1 + 0.2Бт(2х)), при Ос* <70.

131. На рис. 1,2 и 3 приведены компенсационные спектры Е(х)-х202, Е(х)-х198 и Е(х)-х2, соответственно. Кривой 1 показаны компенсационные спектры, а 2 -горизонтальные линии.lg(E(x)x2.<>2)005 О-0.05-0.10 -21. Рисlg(E(x)-xl 98) 0.050-0.05 -0.10-2

132. Рис.2. Компенсационный спектр Е(х)-х198.lg(E(x)x2) 0.05О-0.05-0.10 -2

133. По отклонению компенсационных спектров от горизонтали на рис. 1 и 2 видно, что погрешность при определении показателя степенного закона в 1%0 -1.5 -1.0 -0.5 0 0.5 1.0 1б(Х)2 02

134. Определяя теперь недостающий множитель в степенном законе сх методом наименьших квадратов, находим:0.99161 х-2.

135. Использование же метода наименьших квадратов, согласно первой методике, приводит к степенному закону1.0122д;-2-0063. (2)2 0063

136. На рис. 4 приведен компенсационный спектр Е(х)-х (кривая 1), по отклонению от горизонтальной прямой видно, что метод наименьших квадратов показывает в этом тестовом примере меньшую точность.1. Е(х)х20063)0.05О-0.05-0.10 -22 0063

137. Рис.4. Компенсационный спектр Е(х)-х '

138. Рис.5. Исследуемый спектр Е(х) и степенные законы.

139. В настоящей работе использовались обе описанные методики определения степенных законов.