Авторефераты по всем темам >>
Авторефераты по разное
На правах рукописи
Палымский Игорь Борисович
ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ СЛОЖНЫХ РЕЖИМОВ КОНВЕКЦИИ РЭЛЕЯ-БЕНАРА
01.02.05 - Механика жидкости, газа и плазмы
Автореферат диссертации на соискание ученой степени доктора физико-математических наук
Пермь - 2012
Работа выполнена на кафедре общетехнических дисциплин Новосибирского филиала Военного учебно-научного центра Сухопутных войск общевойсковой военной академии ВС РФ (ВУН - СВ ОВА ВС РФ, г. Новосибирск)
Официальные оппоненты:
доктор физико-математических наук, профессор Тарунин Евгений Леонидович доктор технических наук, профессор Цаплин Алексей Иванович доктор физико-математических наук, профессор Бердников Владимир Степанович Ведущая организация - Институт вычислительных технологий СО РАН, г. Новосибирск
Защита состоится _________________2012 г. в 1515 на заседании диссертационного совета Д212.189.06 в Пермском государственном национальном исследовательском университете по адресу: 614990, г. Пермь, ул. Букирева, 15, зал заседаний Ученого совета ПГНИУ.
С диссертацией можно ознакомиться в библиотеке Пермского государственного национального исследовательского университета.
Автореферат разослан У_____Ф ____________________ 2012 г.
Ученый секретарь диссертационного совета Д212.189.кандидат физико-математических наук, доцент В.Г. Гилев
ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ
Актуальность проблемы. Конвекция, связанная с неоднородным нагревом является одним из самых распространенных видов течений газа и жидкости в природе, поэтому исследование свойств конвективных течений имеет большую важность в виду многочисленных приложений в технических устройствах, астрофизике, геофизике, атмосфере и океане.
Но режим течения жидкости во многих практически важных задачах сложный, турбулентный. Хотя в настоящее время уже не вызывает сомнения тот факт, что турбулентные режимы могут быть получены как решения уравнений гидродинамики, при практическом воплощении идеологии прямого численного моделирования турбулентности возникают огромные технические и методические трудности, связанные, в основном, с большими числами Рейнольдса и, как следствие, огромным числом степеней свободы, которые необходимо учитывать при корректном моделировании.
В классической задаче о конвекции Рэлея-Бенара ситуация существенно упрощается, так как отсутствие среднего течения обуславливает сравнительно невысокие числа Рейнольдса.
Целью настоящего исследования является изучение различных режимов двумерной и трехмерной конвекции Рэлея-Бенара, от стационарных до сложных, стохастических (переходных). При этом конвекция вязкой и несжимаемой жидкости, возникающая при подогреве снизу, рассматривается в областях простой геометрии и большом диапазоне надкритичности r = Ra/Racr, где Ra и Racr - число Рэлея и его критическое значение. Трехмерные расчеты проводятся со свободными от касательных напряжений горизонтальными границами, а двумерные - свободными или жесткими.
Основные трудности при численном моделировании конвекции при высокой надкритичности связаны с наличием растущих в линейном приближении возмущений с огромными инкрементами, например, при r = 3.4104 (в двумерных расчетах со свободными границами) существуют возмущения, растущие как e1367t. Последнее обстоятельство накладывает серьезные ограничения на численные методы и обуславливает актуальность методической работы, направленной на разработку, тестирование, линейный и нелинейный анализ методов расчета конвективных течений.
Решение задачи о конвекции существенно упрощается, делая возможным расчет при существенно большей надкритичности, если конвекция рассматривается со свободными (от касательных напряжений) горизонтальными границами, так как в этом случае возможно применение Фурье декомпозиции во всех направлениях. Однако при этом возникает проблема интерпретации полученных численных результатов, так как подавляющее большинство экспериментальных исследований проведено с жесткими границами и экспериментальные данные по конвекции со свободными горизонтальными границами практически отсутствуют. Таким образом, необходимость проведения корректного сравнения полученных численных результатов с данными эксперимента делает актуальной задачу о сравнительном исследовании характеристик конвекции со свободными и жесткими горизонтальными границами.
Представляется важным исследование динамики вихревого масштаба с ростом надкритичности в двумерных и трехмерных течениях, так как его различное поведение (убывание в трехмерных и возрастание в двумерных) обуславливает качественное различие характеристик двумерных и трехмерных решений при достаточно высокой надкритичности.
Большое значение имеет также исследование динамики пространственных спектров температуры и скорости в двумерных и трехмерных расчетах в областях умеренной и большой горизонтальной протяженности, так как наличие определенных степенных законов в спектрах температуры и скорости указывает, во-первых, на развитый характер течения и, во-вторых, на то, какие физические механизмы доминируют. Причем, при двумерной конвекции в области большой горизонтальной протяженности должны проявиться черты, присущие только двумерным течениям - например, обратный каскад энергии, обусловленный перекачкой энергии из масштаба генерации в большие масштабы.
Цель работы - разработка и программная реализация специального спектрально-разностного метода, его линейный и нелинейный (на модельной нелинейной системе уравнений) анализ, проведение тестовых расчетов и расчетов различных (от стационарных до стохастических, переходных) режимов трехмерной и двумерной конвекции со сравнительным анализом законов изменения интегральных величин, а также исследование каскадных процессов, динамики спектров температуры и скорости.
Научная новизна работы состоит в разработке и программной реализации специального спектрально-разностного метода расчета двумерной и трехмерной конвекции, его линейного и нелинейного анализа (на модельной нелинейной системе уравнений) и проведении расчетов различных (в том числе и стохастических, переходных) режимов трехмерной и двумерной конвекции со сравнительным анализом законов изменения интегральных величин, а также исследование каскадных процессов, динамики спектров температуры и скорости.
В работе впервые:
Х Предложена и реализована в виде компьютерной программы новая реализация спектрально-разностного численного метода расчета течений трехмерной и двумерной конвекции, первого порядка по времени и второго по пространству.
Х В случае двумерной и трехмерной конвекции со свободными горизонтальными границами при Pr = 1 выписаны явные выражения для инкрементов нарастания гармоник в линейных аналогах дифференциальной задачи и предлагаемого численного метода. Показано, что коэффициенты при шагах по времени и пространству, в членах, описывающих схемный эффект, определяются только волновыми числами и не зависят от числа Рэлея.
Х Показано, что, несмотря на то, что численный метод для расчета двумерной конвекции строится как метод первого порядка аппроксимации по времени и второго по пространству, при его линеаризации порядок аппроксимации по времени повышается до второго.
Х Средствами нелинейного анализа, проведенного на модельной нелинейной задаче, показано, что понижение порядка аппроксимации численного метода по времени до первого приводит лишь к незначительному понижению точности вычислений, при этом с первым порядком вычисляется фаза гармонического решения, в то время как его амплитуда - со вторым. Показано, что для практических вычислений можно использовать схему первого порядка аппроксимации по времени с вычислением скоростей по функции тока, полученной на первом этапе расщепления.
Х Предложено новое представление решения в горизонтальном направлении x и связанная с ним численная реализация граничных условий на проницаемой стенке. Это делает возможным неравенство завихренности нулю на боковых границах и существенно облегчает получение стохастических решений в задаче о двумерной конвекции со свободными горизонтальными границами.
Х Показано, что результаты двумерных расчетов с жесткими горизонтальными границами стационарной, валиковой конвекции при невысокой надкритичности 1 r < 4 с удовлетворительной точностью отражают зависимость доминирующей длины волны от надкритичности, хорошо согласуясь при этом с данными эксперимента по числу Нуссельта, профилю средней температуры и изотерме полной температуры.
Х Показано, что, несмотря на наблюдаемые количественные расхождения между результатами расчетов и экспериментом, трехмерные расчеты со свободными горизонтальными границами дают правильные показатели степенных законов изменения среднеквадратичных пульсаций температуры, вертикальной скорости, числа Рейнольдса, среднеквадратичной скорости и кинетической энергии от надкритичности при r 150. В двумерных расчетах аналогичное соответствии наблюдается при сравнительно невысокой надкритичности r 250.
Х Показано, что закон нарастания при числа Нуссельта, поведение пульсаций вертикальной скорости и температуры, а также близость показателей степенных законов для числа Нуссельта, кинетической энергии, среднеквадратичной скорости и энстрофии в двумерных расчетах со свободными границами и стационарном, одновихревом решении свидетельствуют о том, что при высокой надкритичности формирование крупномасштабной структуры является доминирующим фактором, определяющим характеристики течения.
Х Отмечено, что процесс формирования крупномасштабной структуры течения в двумерных расчетах с жесткими горизонтальными границами менее выражен и не так однозначен из-за наличия конкурирующего механизма образования завихренности на горизонтальных границах. Показано, что в областях малой пространственной протяженности процесс формирования крупномасштабной структуры блокируется малым размером области, что делает возможным рост среднеквадратичного волнового числа.
При умеренной горизонтальной протяженности области lx, ly = :
Х Установлено, что в двумерных и трехмерных расчетах со свободными границами действие силы плавучести обуславливает стратификационные спектры для скорости, а доминирование конвективного переноса для температуры - спектр пассивной примеси k-5/3. Эти спектры наблюдаются в двумерных и трехмерных расчетах практически во всем рассмотренном диапазоне надкритичности.
При увеличенной до lx = 4 горизонтальной протяженности области в двумерных расчетах со свободными границами:
Х При 500 r < 4103 длинноволновый участок спектра скорости состоит из двух ветвей, где идентифицируются степенные законы, отражающие конкурирующие механизмы, связанные с силой плавучести и каскадным процессом переноса энергии, а при 4103 r 104 в спектре скорости видны степенные законы k-3 и k-5/3, соответствующие прямому и обратному каскадным процессам переноса энстрофии и энергии, соответственно. При r > 1в спектре скорости видна тенденция к установлению единого степенного закона k-2.6 - искаженного закона каскада энстрофии.
Х При умеренной надкритичности 500 r 2103 в длинноволновом участке спектра температуры виден стратификационный спектр Болджиано-Обухова k-7/5, который при 2103 < r 6103 сменяется спектром Бэтчелора k-1. А при r > 6103 показатель длинноволнового участка спектра температуры совершает колебания в интервале [-1,-0.7], с наиболее вероятным значением -0.8.
Х Установлено, что отмеченные выше перестройки спектров температуры и скорости в двумерных расчетах со свободными границами связаны с формированием обратного каскада энергии при r 2103 (для скорости при r 4103) и уменьшением его энергетической роли при r 6103 (r 104).
Автор защищает:
Х Предложенный и реализованный в виде компьютерной программы специальный спектрально-разностный численный метод расчета сложных, переходных течений трехмерной и двумерной конвекции, первого порядка по времени и второго по пространству.
Х Вывод о том, что в случае двумерной и трехмерной конвекции со свободными горизонтальными границами при Pr = 1, коэффициенты при шаге по времени и пространству в выражении для инкремента нарастания гармоник в линейном аналоге предлагаемого численного метода не зависят от числа Рэлея, а определяются только волновыми числами.
Х Вывод о том, что, несмотря на то, что численный метод для расчета двумерной конвекции строится как метод первого порядка аппроксимации по времени, при его линеаризации порядок аппроксимации повышается до второго.
Х Результаты нелинейного анализа, проведенного на модельной нелинейной системе уравнений, которые показывают, что понижение порядка аппроксимации численного метода расчета 2d,free конвекции по времени до первого приводит лишь к незначительному понижению точности вычислений, при этом с первым порядком вычисляется фаза гармонического решения, в то время как его амплитуда - со вторым. И как следствие этого, для практических 2d,free расчетов можно использовать схему первого порядка аппроксимации по времени с вычислением скоростей по функции тока, полученной на первом этапе расщепления.
Х Результаты двумерных расчетов с жесткими горизонтальными границами стационарной, валиковой конвекции, которые показывают при невысокой надкритичности 1 r < 4 удовлетворительное согласование зависимости доминирующей длины волны от надкритичности и хорошее - по числу Нуссельта, профилю средней температуры и изотерме полной температуры.
Х Результаты тестовых расчетов с учетом различного числа гармоник, с разным числом значащих цифр, а также результат сравнения с данными расчета другим (псевдоспектральным) методом. Перечисленные тесты свидетельствуют об устойчивости вычисления исследуемых интегральных величин и о сходимости их значений при увеличении числа учитываемых гармоник.
Х Сравнение оцененного по результатам расчетов диссипативного и разрешаемого масштабов, где показано, что при расчете двумерной и трехмерной конвекции разрешаемый масштаб всегда меньше диссипативного.
Х Новое представление решения в горизонтальном направлении x и связанную с ним численную реализацию граничных условий на проницаемой стенке. Это существенно облегчает получение стохастических решений в задаче о двумерной конвекции со свободными горизонтальными границами.
Х Вывод о том, что, несмотря на наблюдаемые количественные расхождения между результатами расчетов и экспериментом, трехмерные расчеты со свободными горизонтальными границами дают правильные показатели степенных законов изменения среднеквадратичных пульсаций температуры, вертикальной скорости, числа Рейнольдса, среднеквадратичной скорости и кинетической энергии от надкритичности при r 150. А в двумерных расчетах аналогичное соответствии наблюдается при сравнительно невысокой надкритичности r 250.
Х Вывод о том, что закон роста числа Нуссельта, поведение пульсаций вертикальной скорости и температуры, а также близость показателей степенных законов для кинетической энергии, среднеквадратичной скорости и энстрофии в двумерных расчетах со свободными границами и стационарном, одновихревом решении свидетельствуют о том, что при высокой надкритичности формирование крупномасштабной структуры является доминирующим фактором, определяющим характеристики течения.
Х Замечание о том, что процесс формирования крупномасштабной структуры течения в двумерных расчетах с жесткими горизонтальными границами менее выражен и не так однозначен из-за наличия конкурирующего механизма образования завихренности на горизонтальных границах. В частности, в областях малой пространственной протяженности процесс формирования крупномасштабной структуры блокируется малым размером области, что делает возможным рост среднеквадратичного волнового числа.
При умеренной горизонтальной протяженности области lx, ly = в двумерных и трехмерных расчетах со свободными границами:
Х Результаты исследования пространственных спектров для скорости и температуры, где установлено, что действие силы плавучести обуславливает стратификационные спектры для скорости, а доминирование конвективного переноса для температуры - спектр пассивной примеси k-5/3. Причем указанные спектры наблюдаются в двумерных и трехмерных расчетах практически во всем исследованном диапазоне надкритичности.
При увеличенной до lx = 4 горизонтальной протяженности области в двумерных расчетах со свободными горизонтальными границами:
Х Результаты исследования спектра скорости, где показано, что при 500 r < 4103 длинноволновый участок спектра скорости состоит из двух ветвей, со степенными законами, отражающими конкурирующие механизмы, связанные с силой плавучести и каскадным процессом переноса энергии, а при 4103 r 104 в спектре скорости идентифицируются степенные законы k-3 и k-5/3, соответствующие прямому и обратному каскадным процессам переноса энстрофии и энергии, соответственно. А при r > 104 в спектре скорости видна тенденция к установлению единого степенного закона k-2.6 - искаженного закона каскада энстрофии.
Х Результаты исследования спектра температуры, показывающие, что при умеренной надкритичности 500 r 2103 в спектре температуры виден стратификационный спектр Болджиано-Обухова k-7/5, который при 2103 < r 6103 сменяется спектром Бэтчелора k-1. А при r > 6103 показатель длинноволнового участка температурного спектра совершает колебания в интервале [-1,-0.7], с наиболее вероятным значением - 0.8.
Х Вывод о том, что перестройки спектров температуры (скорости) в двумерных расчетах со свободными границами связаны с формированием обратного каскада энергии при r 2103 (для скорости при r 4103) и уменьшением его энергетической роли при r 6103 (r 104).
Достоверность результатов диссертационной работы подтверждается количественным соответствием полученных в работе численных данных с результатами расчетов других авторов, в том числе и при высокой надкритичности, с результатами расчетов другим, псевдоспектральным методом, тестовыми расчетами с различным числом учитываемых гармоник и разным числом значащих цифр, сравнением с данными эксперимента по стационарной, валиковой конвекции и другими экспериментальными данными при более высокой надкритичности, а также сравнением данных, полученных в двумерных расчетах с жесткими горизонтальными границами и предсказанными теорией Р. Крайчнана.
Научная и практическая значимость результатов диссертационной работы заключается в том, что в ней разработан и реализован в виде программы метод расчета конвекции Рэлея-Бенара, в том числе и сложных ее режимов при высокой надкритичности. Реализован также и параллельный вариант компьютерной программы для многопроцессорных компьютеров с общей памятью.
С научной точки зрения представляется важным вывод о том, что в двумерных расчетах со свободными границами при высокой надкритичности формирование крупномасштабной структуры, сопровождающееся ростом вихревого масштаба до горизонтального размера области, является доминирующим фактором, определяющим характеристики течения, в отличие от трехмерных расчетов, где вихревой масштаб уменьшается с ростом надкритичности. Отмеченное различие в поведении вихревого масштаба в двумерных и трехмерных течениях обуславливает и существенно большую роль граничных условий при двумерном моделировании.
С точки зрения возможности сравнения с данными по конвекции с более реальными с экспериментальной точки зрения жесткими горизонтальными границами представляется важным то, что, несмотря на различие в граничных условиях в эксперименте с жесткими границами и трехмерных расчетах со свободными и обусловленные этим количественные расхождения, трехмерные расчеты дают правильные показатели степенных законов изменения среднеквадратичных пульсаций температуры, вертикальной скорости, числа Рейнольдса, среднеквадратичной скорости и кинетической энергии от надкритичности при r 150. А в двумерных расчетах аналогичное соответствии наблюдается при сравнительно невысокой надкритичности r 250.
Большое значение для понимания физических механизмов развития течения имеет исследование пространственных спектров скорости и температуры, в частности, формирование крупномасштабной структуры в двумерных расчетах со свободными границами в области большой горизонтальной протяженности связано со спектром обратного каскада энергии для скорости. А в области умеренной горизонтальной протяженности, в двумерных и трехмерных расчетах типичны стратификационные спектры для скорости, обусловленные действием силы плавучести и спектр пассивной примеси - для температуры.
Апробация работы. Результаты исследований были представлены на Четвертом сибирском конгрессе по прикладной и индустриальной математике (ИНПРИМ - 2000), Пятой международной конференции УЛаврентьевские чтения по математике, механике и физикеФ (Новосибирск - 2000), International Conference УAdvanced Problems in Thermal ConvectionФ (Perm - 2003), Международной конференции УНелинейные задачи теории гидродинамической устойчивости и турбулентностиФ (НеЗаТеГиУс - 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), 8th International Symposium on Engineering Turbulence Modelling and Measurements (ETMM8, Marseille, France - 2010), Международной научной конференции УЗабабахинские научные чтенияФ (Снежинск - 2010), 7th International Conference on Computational Heat and Mass Transfer (ICCHMTТ11, 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) горизонтальными границами. Исходная система уравнений записана в отклонениях от равновесного решения, который характеризуется равными нулю скоростями и линейным профилем температуры. Число Прандтля Pr, за исключением нескольких тестовых расчетов, везде равно 10.
Приближенное решение задач о конвекции со свободными горизонтальными границами разыскивается в виде конечной суммы собственных функций линейной теории устойчивости, а в случае жестких горизонтальных границ - в смешанном спектрально-физическом пространстве, с использованием метода конечных разностей в вертикальном направлении.
В 3d,free расчетах рассматриваются решения вида:
K N M u(t, x, y, z) = uknm (t) cos( kx) cos( ny) cos( mz), k n m k =0 n=0 m=K -1 N -1 M v(t, x, y, z) = v (t) sin( kx) sin( ny) cos( mz), knm m k =1 n=1 m=K -1 N M -w(t, x, y, z) = w (t) sin( kx) cos( ny) sin( mz), knm n k =1 n=0 m=K -1 N M p(t, x, y, z) = pknm (t) sin( kx) cos( ny) cos( mz), n m k =1 n=0 m=K -1 N M -Q(t, x, y, z) = Q (t) sin( kx) cos( ny) sin( mz), knm n k =1 n=0 m=где u0,0,0 = 0, а k = 0.5 при k = 0, и K = 1 при 1 k K - 1.
Трехмерная конвекция рассматривается в прямоугольном параллелепипеде единичной высоты G = {0 x /, 0 y /, 0 z 1}, где и - минимальные волновые числа, а lx = / и ly = / - относительные горизонтальные протяженности области в x и y направлениях, соответственно.
В случае свободных горизонтальных границ все граничные условия следуют из вида решения, например, при x = 0, / и 0 y /, 0 z 1 получаем: ux = v = w = Q = 0, что соответствует проницаемой и идеально теплопроводной стенке. А при y = 0, / и 0 x /, 0 z 1 находим:
uy = v = wy = Qy = 0 - условия на вертикальной границе конвективной ячейки.
Некоторая искусственность такой постановки граничных условий обусловлена желанием обеспечить преемственность с двумерной постановкой. А на горизонтальных границах при z = 0,1 имеем: uz = vz = w = Q = 0.
В двумерных постановках задачи о конвекции горизонтальное направление y отсутствует, а в случае 2d,rigid конвекции приведенное выше граничное условие для скорости при z = 0,1 заменяется на условие прилипания.
Заметим, что вид решения и граничные условия в горизонтальном направлении х отличаются от обычно рассматриваемых и являются новыми (рассматриваемое в большинстве работ представление решения получается заменой синусов на косинусы и наоборот в направлении х), что делает возможным неравенство нулю завихренности на боковых границах при x = 0, / и существенно облегчает получение стохастических решений в задаче о двумерной конвекции со свободными границами.
В п.1.2 описан двумерный вариант предлагаемого спектральноразностного метода. Исходная система уравнений записана в переменных функция тока, вихрь, температура и в случае свободных границ переход на следующий слой по времени состоит из двух основных шагов.
На первом дробном шаге расщепления в спектральном кЦпространстве рассчитывается рост гармоник без учета их взаимодействия, здесь учитываются все линейные члены исходной системы уравнений. Из соображений устойчивости вычислений вязкость, разбитая на две половинки, учитывается на первом и втором этапах расщепления. Для амплитуд гармоник выписана система из двух обыкновенных уравнений, которая решается по аналитическим формулам, подобным, полученным в линейной теории устойчивости.
На втором этапе расщепления в физическом пространстве учитываются все нелинейные члены или, другими словами, учитывается нелинейное взаимодействие гармоник. Здесь используется неявная схема продольнопоперечной прогонки в реализации В.И. Полежаева, ранее использованная при расчете сложных режимов двумерной конвекции при подогреве сбоку.
Пересчет искомых полей из спектрального в физическое пространство и обратно производится по стандартным программам быстрого преобразования Фурье по косинусам и синусам.
В случае жестких границ, метод решения строится, в принципе, аналогично. Главное отличие заключается в необходимости рассчитывать значение вихря на горизонтальных границах. Используемая здесь формула, основанная на применении алгоритма БПФ, на тестовом примере показала сходимость с первым порядком. Но, по сравнению с известной формулой Тома, имеющей также первый порядок, точность предлагаемой формулы на этом тесте в четыре раза выше. Кроме того, как показывают приведенные здесь результаты методических расчетов, решение задачи о конвекции относительно слабо зависит от числа учитываемых членов в формуле для вычисления вихря.
На первом дробном шаге расщепления, где в спектральном пространстве рассчитывается линейный рост гармоник без учета их взаимодействия, система уравнений записана в переменных функция тока, температура. Выписанная система уравнений решается численно по неявной схеме КранкаНиколсона. В остальном методы расчета конвекции со свободными и жесткими границами не отличаются.
В п.1.3 описан трехмерный вариант предлагаемого спектральноразностного метода. Исходная система уравнений записана в переменных скорость, давление, температура и переход на следующий слой по времени состоит из трех основных шагов.
Алгоритм расчета на первых двух шагах расщепления, в принципе, аналогичен уже описанному в двумерном случае, но на втором этапе расщепления здесь применена явная схема с направленными разностями первого порядка и поправкой А.А. Самарского, повышающей порядок аппроксимации до второго. На третьем этапе расщепления восстанавливается выполнение уравнения неразрывности, нарушенное на втором этапе расщепления.
Предложенный спектрально-разностный метод имеет первый порядок по времени и, кроме того, в трехмерном его варианте на втором этапе расщепления использована явная схема интегрирования по времени.
Но указанные обстоятельства практически не понижают эффективность метода расчета, так как условие корректного моделирования при наличии неустойчивости течения в линейном приближении требует проведения расчетов с достаточно малым шагом по времени. И, кроме того, проведенный во второй главе на модельной системе уравнений нелинейный анализ 2d,free алгоритма показал, что повышение порядка аппроксимации численного метода до второго не приводит к сколько-нибудь заметному увеличению точности вычислений.
Вторая глава ФАнализ и тестирование численных алгоритмовФ посвящена линейному и нелинейному анализу предложенного спектральноразностного численного метода и описанию результатов тестовых расчетов.
В п.2.1 изложены результаты линейного анализа двумерных алгоритмов.
Для конвекции со свободными горизонтальными границами при Pr = выписаны явные выражения для инкрементов нарастания гармоник в линейных аналогах дифференциальной задачи и предлагаемого численного метода в двумерном и трехмерном случаях. Показано, что коэффициенты при и H2 в членах, описывающих схемный эффект, не зависят от числа Рэлея, а определяются только волновыми числами.
А в случае 2d,rigid конвекции сравнивается спектральная кривая, полученная методом ортогонализации из линеаризованной исходной системы уравнений и приближенная, полученная численно из линейного аналога предлагаемого численного метода.
В двумерном случае показано, что в линейном приближении порядок аппроксимации по времени предлагаемого численного метода повышается до второго.
Например, для метода расчета 2d,free конвекции формула для инкремента при Pr = 1 принимает вид:
H1 H2 sr = d + (6k6 + m6) - 4k4 - m4, 96 24 здесь sr и d - инкременты нарастания гармоники для предлагаемого численного метода и дифференциальной задачи, и - минимальные волновые числа в горизонтальном x и вертикальном направлениях, k и m - номера гармоник, и H - шаги по времени и пространству, а < 0 соответствует нарастанию. Из приведенной формулы видно, что схемный эффект (все члены правой части, кроме первого) есть величина порядка О(2 + H2) и что коэффициенты при и H2 определяются только волновыми числами.
В п.2.2 изложены результаты нелинейного анализа (на модельной системе уравнений) 2d,free алгоритма. Поскольку главной целью нелинейного анализа было исследование аппроксимаций нелинейных членов по времени, то анализ ограничился случаем длинных волн, с дискретизацией только по времени. Изменение вида нелинейных членов в исходной системе уравнений (нелинейные члены вида uux заменены на |u|ux, что обуславливает возможное различие в знаке) сделало возможным рассмотрение частных, комплексных решений в виде монохроматической волны в случае дифференциальной системы уравнений и численного метода, с последующим сравнением полученных в виде степенных рядов по выражений.
Показано, что понижение порядка аппроксимации численного метода по времени до первого приводит лишь к незначительному понижению точности вычислений, а именно, с первым порядком при этом вычисляется фаза гармонического решения, в то время как его амплитуда - со вторым. Показано, что для практических вычислений можно использовать схему первого порядка аппроксимации по времени с вычислением скоростей по функции тока, полученной на первом этапе расщепления.
В п.2.3 изложены результаты линейного анализа трехмерного алгоритма.
Показано, что, как и в двумерном случае, коэффициенты при и H2 в выражении для инкремента нарастания гармоники не зависят от числа Рэлея, а определяются только волновыми числами.
На рис. 1 для трехмерной задачи с учетом 653 гармоник приведен инкремент нарастания (первые три старшие моды) при r = Ra/Racr = 950 и Pr = 1. Здесь k - волновое число в горизонтальном направлении х (волновые числа в обоих горизонтальных направлениях на рис. 1 положены равными), а область < 0 соответствует нарастанию.
Из рис. 1 видно, что спектральные характеристики численного метода и исходной дифференциальной задачи близки с заметным количественным отличием в области волновых чисел, соответствующей затухающим ( > 0) гармоникам.
В п.2.4 изложены результаты тестирования предлагаемых алгоритмов и результаты методических расчетов.
В п.2.4.1 приведены результаты тестирования метода расчета 2d,free конвекции.
Результаты расчета стационарРис. 1. Инкремент в дифференциных течений сравниваются с резульальной задаче (кривая 1) и численном татами других авторов при методе (2) надкритичности r 103 и по более простой модели бесконечного числа Прандтля - при r 2103. Отмечено хорошее совпадение по числу Нуссельта.
Проведением при r = 6103 на разных ПК двух расчетов, с полностью совпадающими компьютерными кодами, трансляторами, числом учитываемых гармоник, начальными и граничными условиями показано, что при сопоставлении параметров сложных конвективных течений имеет смысл, в основном, сопоставление средних (интегральных) характеристик, которые вычисляются достаточно устойчиво. А анализ мгновенных величин можно проводить только при достаточно высокой точности вычислений и на малом отрезке времени.
Устойчивость вычисления интегральных характеристик подтверждает также серия методических расчетов при том же значении надкритичности и различном числе учитываемых гармоник - [129x33], [257x65], [513x129], [1025x257] и [2049x513]. Анализ результатов этой серии расчетов показал, что наблюдается сходимость:
Х всех анализируемых интегральных величин, Х профилей средней температуры, среднеквадратичных пульсаций температуры, вертикальной и горизонтальной скорости, Х одномерных пространственных спектров температуры и скорости.
И дополнительно, сравнение временных спектров числа Нуссельта, полученных с учетом [513x129] и [1025x257] гармоник, показало их качественное подобие с совпадением положений максимумов и законов затухания на высоких частотах.
Два 2d,free расчета, проведенные с учетом [129x33] гармоник при r = 6103 с двукратной и четырехкратной точностью вычислений также показали практически полное совпадение интегральных характеристик.
Интегральные характеристики, полученные другим методом - модифицированным методом Галеркина (псевдоспектральным), основанным на вычислении нелинейных членов в физическом пространстве на разностной сетке, имеющей в каждом направлении в два раза больше точек, чем гармоник в том же направлении и схеме Рунге-Кутты четвертого порядка для интегрирования по времени, хорошо согласуются с результатами 2d,free расчетов стохастических течений двухдиффузионной конвекции и конвекции при надкритичности до r = 2103.
В частности, проекции решения на плоскость амплитуд первых гармоник функции тока и температуры, полученные в расчетах стохастической двухдиффузионной конвекции, представляются совпадающими с графической точностью. В этом расчете контролировалась невязка выполнения всех уравнений исходной системы, ее относительное значение не превосходило 1%.
Методические расчеты показали, что при = 1 в 2d,free расчетах достаточно учитывать [129x33] гармоник при r 6103 и [257x65] - при 6103 < r 3.4104. При увеличении горизонтальной протяженности области lx = / (уменьшении ) количество гармоник в направлении x увеличивается обратно пропорционально , таким образом, при = 0.125 расчеты проводились с учетом до [2049x65] гармоник.
Отмечено, что при численном моделировании сложных режимов конвекции стало стандартным требование о примерном равенстве между масштабом диссипации и пространственным разрешением численного метода Hre (полусумма пространственных шагов в горизонтальном и вертикальном направлениях):
Pr1/ Hre, = = 0.531 r-0.3(Ra(Nu -1))1/в 2d,free расчетах при = 1. Показано, что Hre < при всех значениях надкритичности и Hre меньше на 30% при r = 3.4104 и числе гармоник [257x65].
В п.2.4.2 приведены результаты тестирования метода расчета 2d,rigid конвекции.
Результаты расчетов предлагаемым методом сравниваются с данными других авторов при небольшой и высокой (до r = 6103) надкритичности. В последнем случае учитывалось [129x129] гармоник - точек, при этом различие по числу Нуссельта не превосходило 3.9% и, более того, проведенная серия расчетов подтвердила при r 3103 смену степенного закона для числа Нуссельта от надкритичности с Nu -1 ~ r1/3 на Nu -1 ~ r2/7.
Анализ результатов серии методических расчетов при r = 500 и = 1 с учетом [33x17], [65x33], [129x65] и [257x129] гармоник-точек показал сходимость всех рассматриваемых интегральных величин.
Анализ результатов методических тестов показал, что при = 1 в 2d,rigid расчетах достаточно учитывать [129x65] гармоник-точек при 10 r < 2103, [257x65] - 2103 r < 4103 и [513x65] - 4103 r 7103. Показано, что при этом всегда Hre < .
Сравнение результатов 2d,rigid расчетов с теорией Р. Крайчнана показало, что при 500 r 7103 в профиле пульсаций вертикальной скорости идентифицируются два участка, на которых величина пульсаций приближенно следует предсказанным теоретически степенным законам ~z и ~z1/3, при этом рассчитанное число Нуссельта Nu = 1.064r1/3 близко к предсказанному теоретически при 250 r 7103, с практическим совпадением значений при 103 r 4103.
В п.2.4.3 приведены результаты тестирования метода расчета 3d,free конвекции.
В качестве первого теста рассчитывается двумерная, стационарная, валиковая конвекция со сравнением результатов 3d,free и 2d,free расчетов при r = 2.2, Pr = 6.8, = 0.09827 и = . В 2d,free расчете учитывалось [513x17] гармоник, а в 3d,free - [513x5x17] и по координате y проводилось усреднение.
Рассматриваемые характеристики течения оказались близки, в частности, относительное среднеквадратичное отклонение полей температуры составило 1.5%, кинетической энергии - 0.36%, при этом профили средней температуры отличались на 0.4% от максимальной температуры, а число Нуссельта - на 1.5%.
Затем расчетом с учетом различного числа гармоник [33x33x33], [65x65x65] и [129x129x129] при = = и r = 950 проверяется сходимость интегральных характеристик.
Анализ результатов 3d,free расчетов показал, что при надкритичности r 950, = = и учете [65x65x65] гармоник всегда Hre < .
В третьей главе ФЧисленное моделирование конвекции Рэлея-БенараФ представлены результаты численного исследования различных режимов конвекции Рэлея-Бенара, от стационарных до сложных, переходных.
В п.3.1 представлены результаты исследования двумерной, стационарной, валиковой конвекции, проведено сравнение экспериментальных и численных (2d,free, 2d,rigid и 3d,free) результатов при r 4, Pr = 6.8 и = 0.09827. В двумерных расчетах учитывалось [513x17] гармоник, а в трехмерном - [513x5x17].
При r = 2.2 показано, что результаты всех численных расчетов с графической точностью соответствуют экспериментальным данным по профилю средней температуры и изотерме полной температуры. Полученное в 2d,rigid расчете значение числа Nu на 4% отличается от полученного в эксперименте, а в 2d,free и 3d,free расчетах соответствующие отклонения больше и составляют примерно 21%.
В диапазоне надкритичности 1 r < 4 результаты 2d,rigid расчета также хорошо соответствуют результатам эксперимента по числу Нуссельта и с удовлетворительной точностью отражают рост доминирующей длины волны стационарной, валиковой конвекции при увеличении надкритичности.
В п.3.2 рассматривается одно- и двухвихревая стационарная конвекция в квадратной области. Эти решения представляют собой интересный пример стационарного и крупномасштабного 2d,free решения, которое удается рассчитать при r 104 для одновихревого и r 5104 - двухвихревого режима конвекции.
В обоих случаях показано, что учет [513x513] гармоник обеспечивает достаточную точность. В расчетах были определены степенные законы зависимостей основных интегральных величин от надкритичности и показано, что наблюдаемое в расчетах увеличение среднего волнового числа с ростом надкритичности сопровождается увеличением среднеквадратичной температуры.
В п.3.3 рассмотрено поведение вихревого масштаба lme = /Kme как функции надкритичности (рис. 2). Вычисляемый по кинетической энергии, вихревой масштаб характеризует горизонтальный масштаб, на котором сосредоточена основная часть кинетической энергии.
В трехмерных расчетах при r > 150 вихревой масштаб уменьшается, приближенно следуя степенному закону lme = 1.69r-0.084, lme 1.7r-1/12.
А данные 2d,free расчетов показывают, что при достаточно высокой надкритичности характерный горизонтальный размер вихря принимает значение порядка размера области Kme ~ = 1. Учитывая также, что стационарное решение становится периодическим при r = 31.8, заключаем, что наблюдаемая при r > 36 в 2d,free расчетах тенденция к увеличению вихревого масштаба связана с развитием нестационарности.
В 2d,rigid расчетах (рис. 3) тенденция увеличения вихревого масштаба менее выражена и не так однозначна, из-за наличия конкурирующего механизма генерации завихренности на горизонтальных граниРис. 2. Среднее волновое число в гоцах. На рис. 3 представлены данные ризонтальном направлении, кривая 2d,rigid расчетов (кривая 1, 1 - 3d,free и 2 - 2d,free расчеты, = 1), 3.06r0.123 - на 9% уменьшенное = = значение волнового числа, соответствующего наиболее быстрому росту по линейной теории (2), 2.98r0.084 - на 60% увеличенное значение среднего волнового числа, полученного в 3d,free расчетах (3) и асимптотическое значение Kme = 5 (4).
Из рис. 3 видно, что при r 250 происходит перестройка течения, связанная с резким, примерно на 27% уменьшением средРис. 3. Среднее волновое число в него волнового числа с дальней2d,rigid расчетах шим монотонным выходом на асимптотическое значение Kme = 5 при r ~ 1.5103. Отметим также, что среднее волновое число в 2d,rigid расчетах при 400 r 1.5103 растет медленнее, чем в трехмерном решении.
В п.3.4 изучаются спектры скорости и температуры, полученные в 2d,free и 3d,free расчетах. Важность такого исследования обусловлена тем, что наличие четко идентифицируемых пространственных спектров указывает на развитый характер течения и на то, какие физические механизмы доминируют. Отдельно рассмотрены случаи умеренной относительной горизонтальной протяженности области ( = 1) в 2d,free и 3d,free и большой 4 ( = 0.25) - в 2d,free расчетах.
Полученные в 2d,free и 3d,free расчетах временные энергетические спектры температуры в центре конвективной ячейки сопоставляются с экспериментальными данными в газообразном He при криогенной температуре. В эксперименте и расчетах совпадало значение надкритичности (410 и 6400 в 3d,free и 2d,free расчетах, соответственно) и число Прандтля (Pr = 0.8). Большая консервативность спектров и их слабая зависимость от параметров расчета и размерности позволяет сопоставлять данные эксперимента с результатами двумерного расчета.
На рис. 4 приведен временной спектр температуры, полученный в 3d,free расчете. Здесь в качестве характерной частоты f выбрано /D2, а fd = VRe3/4/D, где D есть расстояние между горизонтальными границами. Видно хорошее согласование результатов расРис. 4. Временной спектр темчетов с данными эксперимента, причем пературы, кривая 1 - 3d,free рассущественное различие наблюдается чет, 2 - эксперимент по только на диссипативных частотах f ~ fd.
конвекции в газообразном He Результат сопоставления качественно не при 5K, X.-Zh. et al. Wu, 19изменяется, если данные эксперимента сравниваются с результатом двумерного расчета.
Показано, что при умеренной горизонтальной протяженности области в 2d,free и 3d,free расчетах наблюдаются качественно одинаковые спектры - стратификационные, связанные с действием силы плавучести для скорости и k-2.4 и k-5/3 для температуры. Физическая природа первого из температурных спектров пока неясна, а второй указывает на поведение температуры как пассивной примеси.
На рис. 5 приведен спектр скорости в горизонтальном x направлении при r = 950 и стратификационный спектр Болджиано-Обухова (БО) k-11/5.
На рис. 6 приведен температурный спектр в горизонтальном y направлении при r = 950 и спектр пассивной примеси n-5/3. В другом горизонтальном направлении x, в трехмерных расчетах в спектре температуры идентифицируется участок со степенным законом k-2.4. Этот спектр отмечен также и в двумерных расчетах Рис. 5. Спектр скорости, кривая при r 103.
1 - 3d,free расчет, знак соотВ вертикальном направлении z, в ветствует наиболее быстрому спектре температуры в двумерных и росту в линейном приближении, трехмерных расчетах видны участки со а вертикальная прямая - границе степенными законами k-2.4, k-5/3 и k-5/3, интервала аппроксимации соответственно.
Отмеченное качественное соответствие между спектрами в двумерных и трехмерных расчетах наблюдается практически во всем исследуемом диапазоне надкритичности, а именно, при r 3.41в 2d,free и при r 950 - в 3d,free расчеРис. 6. Спектр температуры в тах.
горизонтальном y направлении, Некоторое отличие спектров скорообозначения как на рис. сти в 2d,free и 3d,free расчетах состоит в том, что в дополнение к наблюдающемуся в двумерных и трехмерных расчетах в горизонтальном направлении x стратификационному спектру БО k-11/5, в трехмерных расчетах, в другом горизонтальном направлении y, виден стратификационный спектр Ламли-Шура n-3.
Наблюдаемое здесь качественное различие физической природы спектров скорости и температуры обусловлено тем, что сила плавучести входит только в уравнение для скорости, действуя на температуру только опосредованно, через нелинейные члены. Отсюда следует, что при последовательном увеличении надкритичности стратификация должна первоначально проявится в спектре скорости.
Другая картина перестройки спектров наблюдалась при большой горизонтальной протяженности области 4. На рис. 7Ц10 приведены спектры скорости и температуры в горизонтальном направлении x при = 0.25, знаком показаны результаты 2d,free расчетов. Диссипация кинетической энергии, ее генерация и потоки к малым и большим масштабам растут при увеличении горизонтальной протяженности области и надкритичности, это обуславливает наблюдаемые перестройки спектров скорости и температуры.
При 500 r < 4103 длинноволновый участок спектра скорости состоит из двух ветвей, со степенными законами, отражающими конкурирующие механизмы, связанные с силой плавучести и каскадным процессом переноса энергии. Конечно, сосуществование двух спектров обуславливает некоторую неточность определения степенных законов и границы этого сосуществования.
При 4103 r 104 в спектре скорости видны участки с четко идентифицируемыми степенными законами, соответствующими каскадным процессам переноса кинетической энергии из масштаба генерации в малые масштабы с последующей диссипацией (k-3 - прямой каскад энстроРис. 7. Спектр скорости при фии) и в большие, с формированием r = 61крупномасштабной структуры течения (k-5/3 - обратный или красный каскад энергии) (рис. 7). Масштаб генерации здесь соответствует наиболее быстрому в линейном приближении росту.Однако, после формирования обратного каскада энергии его энергетическая роль падает, что обусловлено отсутствием диссипации Рис. 8. Спектр скорости при на больших масштабах. В свою очередь, r = 3.41рост диссипации энергии при увеличении надкритичности и связанное с этим усиление энергетической роли прямого каскада энстрофии обуславливает при r > 104 тенденцию установления единого степенного закона в спектре скорости k-2.6, близкого к закону энстрофийного каскада k-3 (рис. 8).
А увеличение роли силы плавучести, неизбежное при увеличении горизонтальной протяженности области, обуславливает при 750 r < 2103 стратификационный спектр БО k-7/5 для температуры (рис. 9).
Дальнейшее развитие обратного каскада энергии и связанное с этим формирование крупномасштабной структуры течения сопровождается перекачкой энергии пульсаций скорости в большие масштабы, в результате чего Рис. 9. Спектр температуры, поле скорости становится крупномасr = 1.25103 штабным при относительно низком уровне пульсаций, и возникает своеобразный вязкоконвективный интервал, где пульсации температуры управляются крупномасштабным полем скорости, что обуславливает появление в спектре температуры участка со спектром Бэтчелора k-1, наблюдаемого в 2d,free расчетах при 2103 r 6103 (рис. 10).
Уменьшение энергетической роли обратного каскада энергии отражается и на спектре температуры - при r > 6103 показатель длинноволнового участка спектра температуры совершает колебания в интервале [-1,-0.7], с наиболее вероятным значением - 0.8.
Рис. 10. Спектр температуры, В п.3.5 рассмотрено формирование r = 31крупномасштабной структуры течения в двумерных расчетах. Проведением серии 2d,free расчетов в областях различной горизонтальной протяженности (, 2, 4 и 8) исследовано последовательное включение больших масштабов, а именно, показано, что при r = 103 горизонтальной протяженности области недостаточно, при r = 61недостаточно длины области 2, а при r = 1.7104 - 4. Здесь выражение недостаточная длина области обозначает тот факт, что горизонтальная протяженность области становится определяющим решение параметром.
Показано также, что в 2d,free расчетах формирование крупномасштабной структуры течения является доминирующим и определяющим.
На рис. 11 и 12 при r = 3104 и lx = 4 приведена функция тока в последовательные моменты времени, значения среднего волнового числа равны примерно 1 и 4.26 4 на рис. 11 и 12, соответственно. Рис. 12 показывает, что при t > 0.1 режим конвекции крупномасштабный и практически одновихревой.
Доминирование процесса формирования крупномасштабной структуры проявляется, в частности, в выходе на степенной закон числа Нуссельта Nu = 0.639r0.371 при r 6103, с поРис. 11. Функция тока при казателем заметно выше 2/7, 0.31 и 1/3 и в t = 7.6510-нефизичном росте среднеквадратичной температуры, пульсаций температуры и вертикальной скорости; а также в близости показателей степенных законов для числа Нуссельта, кинетической энергии, среднеквадратичной скорости и энстрофии в 2d,free Рис. 12. Функция тока при расчетах и стационарном, одновихревом теt = 2.2910-чении (п.3.2).
Отмечено, что этот процесс менее выражен и не так однозначен (см. рис.
2 и 3 при lx = ) в 2d,rigid расчетах из-за наличия конкурирующего механизма генерации завихренности на горизонтальных границах.
Например, среднее волновое число при умеренной горизонтальной протяженности области (lx = ), после скачкообразного уменьшения (при r 250) с монотонным увеличением выходит (r ~ 1.5103) на асимптотическое значение Kme = 5 (рис. 3), а число Нуссельта следует степенному закону Nu = 0.817r0.367 при r 2103, причем показатель этого степенного закона близок к полученному в 2d,free расчетах.
Однако, при конвекции в области малой горизонтальной протяженности (lx = 1) развитие крупномасштабной структуры в 2d,rigid расчетах в значительной степени подавлено малым размером области и это делает возможным рост среднеквадратичного волнового числа при увеличении надкритичности примерно как 3.7r1/8 при r 6103.
В п.3.6 рассматриваются интегральные характеристики трехмерной и двумерной конвекции как функции надкритичности в области умеренной горизонтальной протяженности .
К сожалению, в литературе практически отсутствуют экспериментальные данные по конвекции со свободными от касательных напряжений горизонтальными границами, и это обстоятельство заставляет нас обращаться к данным по конвекции с жесткими.
При сопоставлении данных по конвекции со свободными и жесткими граничными условиями следует ожидать хорошего соответствия по пульсационным характеристикам температуры и вертикальной скорости, так как для них граничные условия на горизонтальных границах правильные и их динамика определяется, в основном, подогревом и гравитацией. Соответствие должно наблюдаться также для числа Рейнольдса, среднеквадратичной скорости и кинетической энергии, так как их динамика определяется, главным образом, вертикальной скоростью.
В п.3.3 отмечено, что основное различие меду трехмерной и двумерной конвекцией связано с разным поведением вихревого масштаба, который рассматривается как функция надкритичности. Результаты 2d,free расчетов показывают, что с ростом надкритичности характерный горизонтальный размер вихря растет до значения порядка размера области, а в трехмерных расчетах характерный горизонтальный размер вихря уменьшается примерно обратно пропорционально корню двенадцатой степени из надкритичности.
Такое поведение вихревого масштаба показывает, что при двумерном моделировании можно ожидать хорошего соответствия вычисленных интегральных характеристик с данными эксперимента только при невысокой (до r ~ 100) надкритичности. Причем, практически полное качественное совпадение пространственных спектров пульсаций температуры и скорости в 2d,free и 3d,free расчетах при конвекции в области относительной горизонтальной протяженности во всем рассматриваемом диапазоне надкритичности показывает, что для консервативных величин, мало изменяющихся при увеличении надкритичности и пространственного разрешения, интервал соответствия может быть много больше.
И, наоборот, при трехмерном моделировании из-за уменьшения вихревого масштаба различие в граничных условиях в расчете и эксперименте постепенно нивелируется и можно ожидать улучшения соответствия вычисленных интегральных характеристик с данными эксперимента при увеличении надкритичности.
На рис. 13 приведены значения среднеквадратичных пульсаций температуры как функций надкритичности, где кривая 1 - данные 3d,free расчетов, 2 - 2d,free, а 3 - 2d,rigid.
Из рис. 13 видно, что величина среднеквадратичных пульсаций температуры следует степенным законам:
qТ = 0.163r-2/15 при r 150 в 3d,free и qТ = 0.316r-0.2 при 40 r 250 - 2d,free расчетах. Оба показателя степенных законов совпадают с полученными Рис. 13. Пульсации температуры в лабораторных экспериментах D.E.
как функция надкритичности Fitzjzrrald,1976 и X.-Zh. Wu,1992, соответственно. Пульсации температуры при r > 250 в 2d,free и во всем рассматриваемом диапазоне надкритичности в 2d,rigid расчетах не следуют какомулибо степенному закону и, более того, пульсации температуры в двумерных расчетах при достаточно большой надкритичности показывают нефизичную тенденцию увеличения.
Деление величины пульсаций скорости на Pr1/3 позволяет сравнивать данные с различным числом Pr. На рис. 14 приведены значения среднеквадратичных пульсаций вертикальной скорости как функций надкритичности, где кривые 1Ц3 обозначают, как и на рис. 13 данные расчетов, кривые 4Ц6 - экспериментальные данные J.W. Deardorff,1967; W.V.R. Malkus,1954;
F.M. Garon, 1973, соответственно и 7 - результат численного расчета R.M.
Kerr,1996.
Из рис. 14 видно, что величина среднеквадратичных пульсаций вертикальной скорости в трехмерном расчете следует степенному закону:
WТ = 15.8r0.395 при r 150, причем показатель этого степенного закона близок к полученным в экспериментах и численном расчете с жесткими гориРис. 14. Пульсации вертикальной зонтальными границами. В двумерных скорости расчетах аналогичное соответствие наблюдается при r < 250, а при достаточно высокой надкритичности пульсации вертикальной скорости показывают нефизично быстрый рост.
Нефизично высокие значения пульсаций температуры и вертикальной скорости при достаточно высокой надкритичности в двумерных расчетах обусловлены формированием крупномасштабной структуры течения.
На рис. 15 знаком () представлено полученное в 3d,free расчетах число Re как функция надкритичности.
Видно, число Рейнольдса в трехмерных расчетах следует степенному закону Re = 1.46r0.497 при r 150, причем показатель этого степенного закона практически совпадает со значением Рис. 15. Число Re в 3d,free рас- 0.49, полученным в экспериментах по четах турбулентной конвекции (X.-Z. Wu, 1992; X. Chavanne,2001).
Простые формулы для среднеквадратичной скорости и кинетической энергии Vrms = RePr и Ek = 2Pr2Re2/2 показывают, что соответствие показателей степенных законов наблюдается также и для среднеквадратичной скорости и кинетической энергии.
ОСНОВНЫЕ РЕЗУЛЬТАТЫ И ВЫВОДЫ 1. Предложен и реализован в виде комплекса программ специальный спектрально-разностный численный метод расчета сложных, переходных течений трехмерной и двумерной конвекции, первого порядка по времени и второго по пространству.
2. В случае двумерной и трехмерной конвекции со свободными горизонтальными границами при Pr = 1 выписаны явные выражения для инкрементов нарастания гармоник в линейных аналогах дифференциальной задачи и предлагаемого численного метода. Показано, что коэффициенты при шагах по времени и пространству не зависят от числа Рэлея, а определяются только волновыми числами.
3. Показано, что, несмотря на то, что численный метод строится как метод первого порядка по времени, в линейном приближении в случае двумерной конвекции порядок аппроксимации по времени повышается до второго.
4. Показано, что спектральные характеристики предлагаемого численного метода в двумерном и трехмерном случаях и исходной дифференциальной задачи хорошо согласуются, заметные количественные отклонения наблюдаются только для инкрементов затухающих гармоник.
5. Средствами нелинейного анализа, проведенного на модельной нелинейной задаче, показано, что понижение порядка аппроксимации численного метода расчета 2d,free конвекции до первого по времени приводит лишь к незначительному понижению точности вычислений, при этом с первым порядком вычисляется фаза гармонического решения, в то время как его амплитуда - со вторым. Показано, что для практических вычислений можно использовать схему первого порядка аппроксимации по времени с вычислением скоростей по функции тока, полученной на первом этапе расщепления.
6. Показано, что результаты 2d,rigid расчета стационарной, валиковой конвекции при невысокой надкритичности 1 r < 4 с удовлетворительной точностью отражают зависимость доминирующей длины волны от надкритичности, хорошо согласуясь с данными эксперимента по числу Нуссельта, профилю средней температуры и изотерме полной температуры.
7. Тестовым расчетом с учетом различного числа гармоник, с учетом различного числа значащих цифр, а также сравнением с результатом расчета другим (псевдоспектральным) методом, показано, что все исследуемые интегральные величины вычисляются устойчиво и при увеличении числа учитываемых гармоник наблюдается сходимость их значений.
8. Сравнением величин разрешаемого и оцененного по результатам расчетов диссипативного масштабов показано, что во всех двумерных и трехмерных расчетах разрешаемый масштаб меньше диссипативного.
9. Предложено новое представление решения в горизонтальном направлении x и связанная с ним численная реализация постановки граничных условий на проницаемой стенке. Это делает возможным неравенство завихренности нулю на боковых границах при x = 0, / и существенно облегчает получение стохастических решений в задаче о двумерной конвекции со свободными горизонтальными границами.
10. Несмотря на наблюдаемые количественные различия между результатами 3d,free расчетов и экспериментом, трехмерные расчеты дают правильные показатели степенных законов изменения среднеквадратичных пульсаций температуры, вертикальной скорости, числа Рейнольдса, среднеквадратичной скорости и кинетической энергии от надкритичности при r 150. В двумерных расчетах аналогичное соответствие наблюдается при сравнительно невысокой надкритичности (до r ~ 250) по среднеквадратичным пульсациям температуры (в 2d,free) и вертикальной скорости (в 2d,free и 2d,rigid расчетах).
11. Закон нарастания числа Нуссельта, нефизичное поведение пульсаций вертикальной скорости и температуры, а также близость показателей степенных законов для числа Нуссельта, кинетической энергии, среднеквадратичной скорости и энстрофии в 2d,free расчетах и стационарном, одновихревом течении показывают, что в двумерных расчетах со свободными горизонтальными границами при высокой надкритичности формирование крупномасштабной структуры является доминирующим фактором, определяющим характеристики течения.
12. Процесс формирования крупномасштабной структуры течения в 2d,rigid расчетах менее выражен из-за наличия конкурирующего механизма образования завихренности на горизонтальных границах. Показано, что в областях малой горизонтальной протяженности процесс формирования крупномасштабной структуры блокируется малым размером области и это делает возможным рост среднеквадратичного волнового числа.
В области умеренной () горизонтальной протяженности в 2d,free и 3d,free расчетах:
13. Действие силы плавучести обуславливает стратификационные спектры для скорости, а доминирование конвективного переноса для температуры - k-2.4 и спектр пассивной примеси k-5/3. Причем, эти спектры наблюдались в двумерных и трехмерных расчетах практически во всем исследованном диапазоне надкритичности.
В области большой (4) горизонтальной протяженности в 2d,free расчетах:
14. Показано, что при 500 r < 4103 длинноволновый участок спектра скорости состоит из двух ветвей, со степенными законами, отражающими конкурирующие механизмы, связанные с силой плавучести и каскадным процессом переноса энергии, а при 4103 r 104 в спектре скорости идентифицируются степенные законы k-3 и k-5/3, соответствующие прямому и обратному каскадным процессам переноса энстрофии и энергии. При дальнейшем увеличении надкритичности (r > 104) в спектре скорости видна тенденция к установлению единого степенного закона k-2.6 - искаженного закона каскада энстрофии.
15. Показано, что при умеренной надкритичности 500 r 2103 в спектре температуры виден стратификационный спектр БО k-7/5, сменяющийся при 2103 < r 6103 спектром Бэтчелора k-1. А при r > 6103 показатель длинноволнового участка температурного спектра совершает колебания в интервале [-1,-0.7], с наиболее вероятным значением - 0.8.
16. Перестройки спектров температуры (скорости) в двумерных расчетах связаны с формированием обратного каскада энергии при r 21(r 4103) и уменьшением его энергетической роли при r 6103 (r 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 P.A. Hieronymus H. The Rayleigh-Benard convection in gas with chemical reactions // Сиб. ж. вычисл. матем. 2007. Т.10.
N.4. С.371Ц383.
6. Палымский И.Б. Численное исследование спектров турбулентной конвекции Рэлея-Бенара // Нелинейная динамика. 2008. Т.4. N.2. С.145Ц156.
7. 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.
8. Палымский И.Б., Герценштейн С.Я., Сибгатуллин И.Н. Об интенсивной турбулентной конвекции в горизонтальном плоском слое жидкости // Изв. РАН. ФАО. 2008. Т.44. N.1. С.75Ц85.
9. Палымский И.Б. Численное исследование спектров трехмерной конвекции Рэлея-Бенара // Изв. РАН. ФАО. 2009. Т.45. N.5. С.691Ц699.
10. Палымский И.Б. О численном моделировании трехмерной конвекции // Вест. Удмурт. универ. серия 1: матем., мех., компьют.
науки. 2009. В.4. С.118Ц132.
11. Палымский И.Б. О качественном различии решений двумерной и трехмерной конвекции // Нелинейная динамика. 2009. Т.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.B. 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 I.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.1. Paris, France, 2005. P.101Ц107.
18. Palymskiy I.B., Fomin P.A., Hieronymus H. Rayleigh-Benard convection in chemical equilibrium gas // Fourth International Conference on Computational Heat and Mass Transfer (ICCHMTТ05). Proceedings in book: Progress in Computational Heat and Mass Transfer, V.1. Paris, France, 2005.
P.116Ц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 20(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.