2009
ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА
Математика. Механика. Информатика Вып.з(29)
УДК 532.5
Устойчивость осесимметричных течений в цилиндрической области
Е. Л. Тарунин, А. М. Шарапова
Пермский государственный университет, 614990, Пермь, ул. Букирева, 15
Численно исследовано ветвление осесимметричных решений тепловой конвекции в цилиндрической геометрии при подогреве снизу. Характер ветвления был вызван конечным возмущением температуры на нижней границе. Вычислительные эксперименты показали, что рассмотренное возмущение существенно повышает амплитуду предпочтительного решения.
Введение. Явление тепловой конвекции широко распространено в природе и технике [1]. Из метеорологии известно, что при нагревании земной поверхности обычно возникает течение двух типов.
Первый тип характеризуется восходящей тепловой струей воздуха вблизи центра нагретой области.
Второй тип течения отличается опускным движением холодного воздуха вблизи центра нагретой области. Такой характер тепловой конвекции обнаруживают летчики при пролете над нагретым участком.
При вхождении в такую зону, охваченную тепловой конвекцией, на самолет действует вертикальная компонента скорости воздуха. Вначале он испытывает восходящую струю, потом опускную, а затем снова восходящую. Такой тип течения называют "Костром Каина" (для плоской области соответствующие течения иногда называли течениями в виде "кошачьи глаза"). Средствами вычислительных экспериментов выясняются условия, при которых реализуется и устойчиво течение, похожее на Костер Каина.
Постановка задачи. Исследуется естественная конвекция в цилиндрической полос-
ти при подогреве снизу (см. рис. 1). При этом основное внимание уделено вопросам устойчивости двух форм осесимметричного течения при конечных изменениях температуры на нижней границе. Для плоской области подобное исследование было выполнено в работе [5].
Используемые уравнения естественной конвекции в приближении Буссинеска [1] в переменных двухполевого метода [2, 3] имеют вид
дф _ 1
дt г
дщ дф дщ дф 1 дщ
дг дг дг дг г дг
+
+
д 2ф д ф \дф_ф_
дг2 дг2 г дг г2
+С дт,
дг
(1)
д 2щ д 2щ 1 дщ
г ф =--------------------— л-----------
дг дг г дг
(2)
Т- 1
дХ г
1
дщ дТ дщ дТ дг дг дг дг
+
+ -
Рг
д 2Т д 2Т 1 дТ —т л----------------т л---------
дг дг г дг
(3)
Исследования были поддержаны грантами РФФИ (Урал) 07-01-96040 и CRDF-PE-009
© Е. Л. Тарунин, А. М. Шарапова, 2009
Уравнения записаны в цилиндрических координатах в предположении осевой симметрии решения. Функция тока связана с компонентами скорости соотношениями
1 дщ
°г = ' л г д2
1 дщ г дг
(4)
£(0) = 0.
Т(0) = 1 - г .
(6)
Уравнения написаны в безразмерных переменных. В качестве единиц расстояния, времени, скорости, температуры выбраны соответственно: Н - высота цилиндра, Н2 /V (V - коэффициент кинематической вязкости ), V Н, 0= тах Т , г,0).
Г
Значения температуры на нижней границе соответствуют подогреву снизу. В том случае, когда температура на нижней границе постоянная, возможно равновесное решение [1]. Случай с нарушением условия равновесия для плоской полости описывался в работе [3].
Г еометрия рассматриваемой области изображена на рис. 1.
Устойчивость такого решения рассматривалась в [1]. При 0 равновесное решение невозможно.
Уравнения (1)—(3) содержат два безразмерных параметра - числа Грасгофа и Пран-дтля:
gP®H3
О = -
Рг = -.
Рис. 1. Геометрия области
Нижняя граница полости предполагается твердой, неподвижной и непроницаемой.
Боковая и верхняя границы области считаются свободными от напряжения (штриховые линии на рис. 1 отражают этот факт), значение вихря скорости на этих участках границы и на оси цилиндра полагается равным нулю. Боковая граница области предполагается теплоизолированной. Температура на верхней границе области принята за начало отсчета, а температура на нижней границе определяется следующим образом:
— • г
Т(ї, г,0) = 1 + єєоб(——) . (5)
Вариант при є = 0 соответствует постоянному ("невозмущенному") значению температуры на нижней границе, при котором возможно равновесное решение:
В данной работе исследуются нелинейные конвективные движения при различных значениях е и О. При этом выясняются границы перестроек течений: перестройка течения с восходящей струей по центру в течение с нисходящей струей и наоборот.
В качестве начального условия для получения первого решения при е = 0 использовалось равновесное решение и к нему добавлялось возмущение для вихря скорости с амплитудой А(А=±100):
Ф(0)= Аг (Я - г)г(1 - £). (7)
В случае е ^ 0 использовались различные начальные условия, которые указаны при обсуждении результатов вычислительных экспериментов. При выборе метода продолжения по параметру в качестве начального условия для новых параметров применялось стационарное решение, полученное при предыдущих значениях параметров (е ,О).
Метод решения. Поставленная задача решалась методом сеток на квадратной сетке с шагом Аг = Аг = Н = 1/ N. Эволюционные уравнения для вихря и температуры аппроксимировались двухслойной схемой с использованием центральных разностей для первых производных по пространственным координатам. Уравнение (2), связывающее функцию тока с вихрем скорости, решалось методом последовательной верхней релаксации (ПВР) с параметром релаксации, соответствующим оптимальному значению для прямоугольной области [3]. Настроечные вычислительные эксперименты подтвердили близость этого значения параметра релаксации к оптимальному. Шаг по времени вычислялся из соображений устойчивости с небольшим запасом. Вычисление вихря скорости на нижней твердой границе области осуществлялось по формулам Тома. Температура на вертикальных границах (г = 0, г = Я) вычислялась по формулам, обес-
2
V
печивающим второй порядок аппроксимации по пространственным координатам.
Результаты. В расчетах фиксированными были число Прандтля Р=0.7(воздух) и безразмерное значение радиуса цилиндра R=1. Число Грасгофа изменялось до 8000, а величина Є от -0.9 до +0.9.
Стационарное решение находилось методом установления. Точность стационарности решения определялась по установлению суммы S модулей относительных изменений трех интегральных характеристик течения -числа Нуссельта, кинетической энергии, максимального значения функции тока. Решение считалось стационарным при выполнении неравенства Б <Аї ■ 10-4.
0 5 0 6 0 7 0 8 0 9
Кроме того, для запаса надежности требовалось, чтобы безразмерное время счета (0 было не менее 0.5.
Вычислительный эксперимент показал, что при 0=8000 неравенство для суммы относительных изменений выполняется уже при г « 0.25. Соответствующие линии тока и изотермы изображены на рис. 2. В этом случае теплая жидкость поднимается по краям полости, а холодная опускается в центре полости. Это течение похоже на Костер Каина. Заметим, что для подобных течений в плоской области употреблялось название - течение типа "кошачьи глаза".
При смене знака возмущения (А=-100) стационарное решение описывает другое направление вращения жидкости (рис. 3).
Рис. 2. Изолинии функции тока и изотермы
Рис. 3. Изолинии функции тока и изотермы
Оценка погрешности полученного реше- грешность Е при Н=И20 составляет около ния по кинетической энергии Е по идее Рунге- 2%.Отметим, что значение 0=8000 является Ромберга показывает, что относительная по- максимальным из использованных в расчетах.
Обработка результатов решений при значениях числа Грасгофа
О = 5000 + к • 500 (k=0,1,...,6) методом наименьших квадратов позволила выяснить, что для максимума функции тока выполняется корневой закон Ландау [4]. Следствием этого закона являются линейные зависимости для кинетической энергии:
Ек = 0,0150(0 - 2980), (А=100),
(8)
Ек = 0,0204(0 - 3114). (А=-100).
Критическое число О* в этих формулах отличается от среднего 3050 примерно на 2.3% при И=0.1 (при И =1/30 это отличие менее 0.5%).
Перейдем к обсуждению результатов вычислительных экспериментов, в которых исследуется влияние параметра е на интенсивность течения. Выяснено, что положительные е увеличивают интенсивность течения с восходящей струей воздуха по центру цилиндрической области, а отрицательные -уменьшают. Для течений с нисходящей струей воздуха зависимости интенсивности течения от е противоположные (увеличение интенсивности течения происходит при отрицательных значениях параметра, а уменьшение при положительных).
Соответствующее увеличение модулей значений функции тока показано на рис. 4 (штриховая линия соответствует е = 0).
Рис. 4. Зависимость экстремальных значений функции тока от числа Грасгофа
Перейдем к обсуждению вариантов, в Итоговые результаты для этих вычис-
которых изменение величины е приводит лительных экспериментов изображены на к"перекрутке" течения. рис. 5.
Рис. 5. Зависимость экстремальных значений функции тока от числа Грасгофа при є = 0.01 ■ т(т=0Д,2,...)
Тесно расположенными звездочками на рис. 5 указана величина экстремального значения функции тока при различных е . Эти значения соответствуют стационарным решениям. При критических значениях е происходила перестройка течения (стрелками показаны соответствующие переходы).
Дадим пояснение происходящим перестройкам для фиксированного значения числа 0=8000. Вначале было получено решение с нисходящей струей при е = 0 . Затем стационарное решение рассчитывалось методом продолжения по параметру при значениях ет+\ =ет + 0.01 до т=63. При т=64
(е = 0.63 — 0.64) произошла резкая перестройка течения: при е = 0.63 экстремальное
значение функции тока равно +1.171, а при е = 0.64 равно -5.755.
Критические значения еч в интервале значений от 3000 до 6000 описывается линейной зависимостью
ем« 0.13 -10—3 О - 0.28. (9)
Аналогичные расчеты были выполнены для выяснения критических значений е*2, при которых происходит переход с восходящей струи в нисходящую. Выяснено, что при О < 6000 можно считать, что е*2 также описывается зависимостью (9), но при О > 6000 значения |е.2| меньше еп. Одна из причин отсутствия симметрии критических величин (| еп | ^ | е*21), на наш взгляд, неравенство
площадей с положительной и отрицательной добавками температуры на нижней границе при е Ф 0.
Обсудим карты изолиний функции тока
Рис. 6. Линии тока и
Рис. 7. Линии тока и
Выводы
Показано, что возникшее течение типа костра Каина довольно устойчиво по отношению к изменению температуры подстилающей поверхности. Однако для его возникновения требуется особое (е > 0) распределение температуры подстилающей поверхности. Поэтому естественно, что такое течение наблюдается не так часто. Найдены критические значения еч и е*2, при которых происходит смена типа течения, и выяснено, что переход с одной ветки на другую имеет широкий гис-терезисный интервал по параметру е .
и изотерм при 0=6000 до перестройки и после перестройки. Соответствующие карты изолиний при е = 0.49 и е = 0.5 показаны на рис.6, 7.
изотерм при е = 0.49
изотермы при е = 0.5
Список литературы
1. Гершуни Г.З. Конвективная устойчивость несжимаемой жидкости / Г.З.Гершуни, Е.М.Жуховицкий. М.: Наука, 1972. 392 с.
2. Гидродинамика: учен. зап. № 327. Вып. VI. Пермь, 1975. 187 с.
3. Тарунин Е.Л. Вычислительный эксперимент в задачах свободной конвекции / Е.Л.Тарунин. Иркутск: Изд-во Иркут. ун-та, 1990. 228 с.
4.Ландау Л.Д. Гидродинамика / Л.Д.Ландау, Е.М. Лифшиц. М.: Наука, 1966.
5. Тарунин Е.Л. Конвекция в замкнутой полости, подогреваемой снизу, при нарушении условий равновесия / Е.Л.Тарунин // Изв. АН СССР. Механика жидкости и газа. 1977. № 2. С.203-207.
Е. fl. TapyHUH, A. M. ^apanoea
Stability of axis symmetrical flows in a cylindrical domain
E. L. Tarunin, A. M. Sharapova
Perm State University, 614990, Perm, Bukirev st., 15
Numerical methods were used for investigation of branching solutions of equations of heat or free convection in a cylindrical domain. A character of the branching was created by final perturbation of temperature on the bottom the boundary. Numerical experiments showed that the considered perturbation has an essential influence for stability of two flows - with up and down velocity in the central part of the domain.