МОДЕЛИРОВАНИЕ КОНВЕКТИВНОГО ТЕПЛОПЕРЕНО-СА ВО ВРАЩАЮЩЕЙСЯ ЗАМКНУТОЙ ПОЛОСТИ С ЛОКАЛЬНЫМ ИСТОЧНИКОМ ЭНЕРГИИ Текст научной статьи по специальности «Физика»

ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА

2017 • ФИЗИКА • Вып. 1 (35)

УДК 536.2

PACS 02.60.Cb; 44.25.+f

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

С. А. Михайленкоa, М. А. Шереметb

a Томский государственный университет 634050, Томск, пр. Ленина, 36 email: stepanmihaylenko@gmail.com b Томский государственный университет 634050, Томск, пр. Ленина, 36 email: sheremet@math.tsu.ru

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

Ключевые слова: смешанная конвекция; замкнутая полость; вращение; локальный источник энергии

Поступила в редакцию 19.03.2017; принята к опубликованию 25.04.2017

Simulation of convective heat transfer inside a rotating enclosure with a local heat source

S. A. Mikhailenkoa, M. A. Sheremetb

a Tomsk State University, Lenin Avenue 36, 634050, Tomsk email: stepanmihaylenko@gmail.com b Tomsk State University, Lenin Avenue 36, 634050, Tomsk email: sheremet@math.tsu.ru

Numerical simulation of transient convective heat transfer inside a rotating enclosure with a local heater of constant temperature has been carried out. The boundary-value problem formulated using dimensionless stream function, vorticity and temperature has been solved by finite difference method of the second order accuracy on the basis of uniform mesh. Analysis has been conducted in a wide range of Rayleigh and Taylor numbers as well as time parameter reflecting a formation of periodic flow and heat transfer patterns. The flow and temperature fields illustrating a development of considered process with governing parameters have been shown.

Keywords: mixed convection; enclosure; rotation; local heater

Received 19.03.2017; accepted 25.04.2017 doi: 10.17072/1994-3598-2017-1-19-25

© Михайленко С.А., Шеремет М.А., 2017

1. Введение

Исследования гидродинамики и теплопереноса во вращающихся системах имеют очевидное фундаментальное и прикладное значение в связи с широким спектром технических систем, в которых реализуется такой процесс [1-3].

На сегодняшний день проведено небольшое количество теоретических и экспериментальных исследований конвективного теплопереноса в замкнутых областях в условиях вращения. Так, например, влияние силы Кориолиса на теплопередачу при естественной конвекции в холодной чистой воде проанализировано в [4] при температурах, соответствующих ее максимальной плотности. Показано, что механизм конвекции полностью подавляется при высоких значениях параметра Ко-риолиса, а средняя скорость передачи тепла является нелинейной функцией температурного напора. Конвективный теплоперенос во вращающейся дифференциально обогреваемой пористой полости численно изучен в [5] на основе обобщенной модели Форхгеймера. Установлено, что интенсивность теплообмена увеличивается с ростом пористости внутренней среды и числа Дарси и уменьшается при увеличении числа Тейлора. В [6] проведено численное моделирование естественной конвекции в кубической полости при подогреве снизу и вращении вокруг вертикальной оси симметрии. Показано, что возникновение конвекции всегда происходит колебательным образом с образованием пространственных циркуляций, распространяющихся вдоль боковых стенок. При высоких значениях числа Рэлея наблюдается переход к установившемуся режиму до тех пор, пока параметр Кориолиса не слишком велик. Анализ конвективной теплопередачи внутри радиально вращающегося наклонного изотермического канала квадратного сечения проведен в [7]. Результаты показали, что структуры вихревого потока и, следовательно, явления переноса тепла совершенно отличаются от предыдущих исследований с нулевым углом наклона. При этом канал с углом наклона в 45° дает наилучшую общую эффективность теплопередачи.

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

ности пассивной системы охлаждения элементов электронной аппаратуры в условиях вращения.

2. Математическая модель и метод решения

Вязкая теплопроводная жидкость находится в замкнутой полости, имеющей высоту Н и длину Ь. Сила тяжести направлена вертикально вниз по оси у (рис. 1). Полость вращается относительно оси, проходящей перпендикулярно плоскости рисунка через центр области решения, с постоянной угловой скоростью ю0. На нижней границе полости находится локальный источник постоянной температуры Тк. Вертикальные стенки (х = 0 и х = Ь) поддерживаются при постоянной температуре Тс (Тк > Тс), остальные стенки являются адиабатическими. При численном моделировании рассматривается случай квадратной полости (Ь = Н).

ччч\ччч<

дТ/ду = 0 Рис. 1. Область решения задачи

Дифференциальные уравнения в частных производных, описывающие процесс переноса массы, импульса и энергии в рассматриваемой области, имеют вид нестационарных уравнений Обербека-Буссинеска [8]:

VV = 0,

р^ ^^ + юх(юх г ) + 2 (юхУ =

= -Ур + у^2У + р£, ят__

— + (¥ ут ) = аУ2Т.

(2.1)

(2.2)

(2.3)

Здесь V - вектор линейной скорости; р - плотность; t - время; ю - вектор угловой скорости вращения; г - радиус вектор положения; р - статическое давление; ц - коэффициент динамической вязкости; ^ - вектор ускорения свободного

падения; T - температура; a - коэффициент температуропроводности.

В уравнении движения (2.2) статическое давление можно представить в виде суммы гидростатического давления и гидродинамического давления следующим образом:

Р = Ре + Ра ,

ре (юх(юх Г )) = ^Ре +РеЯ, УРе = Ре [ 8 -Юх(Юх Г)] • Из соотношения (2.5) получаем: Фр • / \ 2

= -РеЯ • вШ(Ю0/) + РеХЮо ,

др» / \ 2 —е = -Ре8 • СОв (Юо^) + РеУЮ •

ду

или

Тогда гидростатическое давление ре примет следующий вид:

Ре = РеЯ^Ш (<Ю/) Ре8УСОв (Ю0?) + + 0^5рю( х2 + 0^5РЮ(у2 •

Подставляя соотношения (2.4) и (2.5) в уравнение движения (2.2), получим:

РI + юх(юх г ) + 2 (ю х V ) 1 =

Ж

или

■-УРа + цУ2— + ря --Ре [ 8 -Юх(Юх Г)]

Р|ддг+(—]=-ура+

+ цУ— -(Ре -р) 8 +

+ (Ре -р)^(юх(юхГ))-2р-(юх V)•

вращательная

выталкивающая

сила:

переменных «функция тока - завихренность» можно переписать следующим образом:

д2 Т

д2 Т

(2.4)

дХ2 дУ дО дТ дО

дх дУ дХ

= -О,

(2.10)

где гидростатическое давление Ре определяется из уравнения Эйлера с учетом вращения полости:

1

^д2 О

л/га

дТ дО _ дУ дУ "

д2 ол

чдХ2 дУ' ,

(2.11)

Яа

(2.5)

Рг • Та

д© . , д© . , ' — сов (х)--в1П (х)

дХ w дУ w

д© дТ д© дТ д©

дх дУ дХ дХ дУ

1

^д2 ©

(2.6)

Ргл/Та

д2 ©^

дХ2 дУ2

(2.12)

(2.7)

(2.8)

(2.9)

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

(ре-р)-(юх(юх г )) и сила Кориолиса: ю х V). Использование приближения Бус-синеска дает следующее соотношение: Ре =Р[1 -Р(Те - Т)] иЛи Ре -Р = РР(Т - Т) . Тогда уравнения (2.1), (2.3) и (2.9) в декартовой системе координат в безразмерных преобразованных

Здесь X, Y - безразмерные координаты, соответствующие координатам х, у, х - безразмерное время; © - безразмерная температура; Т - безразмерная функция тока (и = дТ/дУ, V = -дТ/дХ) ; Ц

V - безразмерные составляющие скорости в проекции на оси X, Y соответственно; О - безразмерная завихренность скорости (О = д—/дХ -ди/дУ) ; Рг = ц/(ра) - число

Прандтля; Та = р2ю((X4/ц2 - число Тейлора; Яа = ряР(ГА - Тс) X3/(ца) - число Рэлея; в - температурный коэффициент объемного расширения.

Начальные и граничные условия для сформулированной задачи математической физики (2.10)-(2.12) рассматривались в следующем виде.

В начальный момент времени предполагалось, что жидкость, заполняющая полость, неподвижна, поэтому Т (Х, У ,0) = О (Х, У ,0) = 0. Начальная

температура - ©(Х ,У ,0) = 0.

Граничные условия:

• на границах Y = 0 и Y = 1:

Т = 0, дТ/дУ = 0, д©/дУ = 0;

• на границах Х = 0 и Х = 1: Т = 0, дТ/дХ = 0, © = 0;

• на поверхности источника энергии: Т = 0, дТ/Ш = 0, © = 1.

Сформулированная краевая задача для дифференциальных уравнений в частных производных (2.10)-(2.12) с соответствующими начальными и граничными условиями решена методом конечных разностей [9-11] на равномерной сетке. Подробное описание используемого численного метода представлено в [11].

Разработанный метод решения был протестирован на модельной задаче конвективного тепло-переноса в дифференциально-обогреваемой вра-

+

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

Рис. 2. Изотермы при углах поворота 0°, 270° и 180°в сравнении с данными: б - [12], в - [13]

п-1-1-1-1-г

О 10 20 30 40 50 60 I

Рис. 3. Зависимость среднего числа Нус-сельта от времени и размерности разностной сетки

На рис. 3 представлены временные зависимости среднего числа Нуссельта на поверхности источника энергии от размерности разностной сетки при Ra = 105, Pr = 0.7, Ta = 105. Из рис. 3 видно, что процесс становится периодическим после пятого оборота полости. Использование сеток размерности 50x50, 100x100 и 150x150 не приводит к

значительным расхождениям. Поэтому основные исследования были проведены на равномерной сетке 100x100.

3. Результаты численного моделирования

Численный анализ проведен при следующих значениях безразмерных комплексов: Ra = 104-106; Ta = 104-106; Pr = 0.7; 0 < т < 70. Исследовано влияние этих параметров на распределения изолиний функции тока и изотерм, а также на среднее число Нуссельта по поверхности источника энергии.

На рис. 4 представлены линии тока и изотермы для Ra = 105, Ta = 105 при различных углах поворота полости после 10 полных оборотов. Как было отмечено выше на основании рис. 3, начиная с шестого полного оборота наблюдается формирование периодической картины как в теплопереносе, так и в гидродинамике. Распределения, представленные на рис. 4, отражают главным образом влияние инерционной составляющей движения, в связи с чем данные поля не совпадают с аналогичными распределениями при отсутствии вращения, но при наличии соответствующего наклона полости. Поворот полости сопровождается слиянием или разделением вихрей, что в динамике можно описать как постоянная эволюция вихря, зарождающегося в зоне (0 < X <0.5). Причем в пределах одного оборота полости такая эволюция вихря происходит дважды, а именно при изменении угла поворота в диапазоне от 0 до 3л/4 наблюдается зарождение и развитие вихря в области (0 < X <0.5), что проявляется в деформации и ослаблении конвективной ячейки в зоне (0.5 < X <1.0).

Дальнейший поворот полости от 3л/4 до 2л вновь характеризуется развитием вихревой структуры со стороны (0 < X <0.5) и подавлением ранее образованного вихря.

При этом поле температуры позволяет проследить эволюцию теплового факела над правой (при нулевом угле поворота) угловой точкой источника энергии с координатами (0.6,0.2). Один такт развития вихря при изменении угла поворота от 0 до 3л/4 сопровождается интенсификацией отмеченного выше теплового факела, при этом тепловой факел, сформированный при нулевом угле поворота полости, в области слева от источника энергии ослабевает и исчезает. Следующий такт при изменении угла поворота полости от 3л/4 до 2л характеризует обратный процесс - усиливающийся ранее тепловой факел начинает ослабевать, а ослабленный факел развивается. Отмеченные особенности можно детально проследить по зависимости среднего числа Нуссельта на источнике энергии от угла поворота (рис. 5).

в

Рис. 4. Изолинии функции тока Т и температуры © при Ra = 105, Ta = 105 и различных углах поворота: а - 0, б - л/4, в - л/2, г - 3л/4, д - л, е - 5л/4, ж - 3л/2, з - 7л/4, и - 2л

Рис. 5. Зависимость среднего числа Нус-сельта от угла поворота полости после десятого полного оборота

Описанная выше граница такта (ф = 3л/4) отражает максимальный теплоотвод с поверхности источника энергии. При этом в пределах первого такта происходит интенсификация теплоотвода, а на втором такте интенсивность теплообмена уменьшается.

На рис. 6 представлено влияние числа Рэлея на среднее число Нуссельта при Та = 105. Рост величины выталкивающей силы проявляется в интенсификации теплообмена. Следует также отметить, что с увеличением Яа изменяется характер колебаний вследствие более интенсивного течения.

Рис. 6. Зависимость среднего числа Нус-сельта от времени и числа Рэлея

Влияние числа Тейлора на среднее число Нус-сельта при Ra = 105 показано на рис. 7. Увеличение безразмерной угловой скорости вращения полости приводит к снижению интенсивности теплоотвода вследствие большой инерционности среды, что проявляется также в снижении линейной скорости движения внутри полости.

Рис. 7. Зависимость среднего числа Нус-сельта от времени и числа Тейлора

4. Заключение

Численно проанализировано влияние вращения полости на конвективный теплоперенос при нали-

чии локального источника энергии. Исследования проведены в широком диапазоне изменения определяющих параметров: Ra = 104-106; Ta = 104-106; Pr = 0.7; 0 < х < 70. Установлено, что при Ra = Ta = 105 после пятого полного оборота полости эволюцию гидродинамики и теплопереноса можно разбить на два такта: при 0 < ф < 3л/4 происходит полное развитие вихря, формирующегося в зоне (0 < X <0.5), от момента зарождения до диссипации, что сопровождается ростом интенсивности теплоотвода от источника энергии; а при 3л/4 < ф < 2л наблюдается полная эволюция новой вихревой структуры, зарождающейся также в зоне (0 < X <0.5), при этом интенсивность теплоотвода снижается. Показано, что с ростом числа Рэлея не только происходит увеличение среднего числа Нус-сельта, но и изменяется характер зависимости Nuavg = f (х) . С ростом числа Тейлора интенсивность теплоотвода снижается.

Работа выполнена при финансовой поддержке Совета по грантам Президента РФ для молодых российских ученых (грант МД-2819.2017.8).

Список литературы

1.

Гореликов А. В., Ряховский А. В., Фокин А. С. Численное исследование некоторых нестационарных режимов естественной конвекции во вращающемся сферическом слое // Вычислительная механика сплошных сред. 2012. Т. 5. № 2. С. 184-192.

2. Вяткин А. А., Иванова А. А., Козлов В. Г., Сабиров Р. Р. Конвекция тепловыделяющей жидкости во вращающемся горизонтальном цилиндре // Известия РАН. Механика жидкости, газа и плазмы. 2014. № 1. С. 21-31.

3. Baig M. F., MasoodA. Natural convection in a two-dimensional differentially heated square enclosure undergoing rotation // Numerical Heat Transfer, Part A. 2001. Vol. 40. N. 2. P. 181-202. Sandaravadivelu K., Kandaswamy P. Nonlinear convection in a rotating square cavity // Acta Me-chanica. 2000. Vol. 144. N. 1. P. 119- 125. Saleh H., Alhashash A. Y. N., Hashim I. Rotation effects on non-Darcy convection in an enclosure filled with porous medium // International Communications in Heat and Mass Transfer. 2013. Vol. 43. P. 105 - 111.

Sedelnikov G. A., Busse F. H., Lyubimov D. V. Convection in rotating cubical cavity // European Journal of Mechanics B/Fluids. 2012. Vol. 31. P. 149-157.

7. Fann S., Yang W. J. Convective heat transfer in rotating square channel with oblique cross section // Computational Mechanics. 1994. Vol. 14. N. 5. P. 513-527.

4.

5.

6.

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

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

10. Гибанов Н. С., Шеремет М. А. Влияние формы и размеров локального источника энергии на режимы конвективного теплопереноса в квадратной полости // Компьютерные исследования и моделирование. 2015. Т. 7. № 2. С. 271-280.

11. Астанина М. С., Шеремет М. А. Моделирование термогравитационной конвекции с переменной вязкостью в замкнутой полости с локальным источником энергии // Вестник Пермского Университета. Серия: Физика. 2015. Вып. 3 (31). С. 52-58.

12. Hamady F. J., Lloyd J. R., Yang K. T., Yang H. Q. A study of natural convection in a rotating enclosure // J. Heat Transfer. 1994. Vol. 116. P. 136143.

13. Tso C. P., Jin L. F, Tou S. K. W. Numerical segregation of the effects of body forces in a rotating, differentially heated enclosure // Numerical Heat Transfer. 2013. Vol. 51. P. 85-107.

References

1. Gorelikov A. V., R'akhoskii A. V., Fokin A. S. Numerical investigation of some transient natural convection modes in a rotating spherical layer. Computational Continuum Mechanics, 2012, vol. 5, no. 2, pp. 184-192.

2. Vyatkin A. A., Ivanova A. A., Kozlov V. G., Sabi-rov R. R. Convection of a heat-generating fluid in a rotating horizontal cylinder. Fluid Dynamics, 2014, vol. 49, no. 1, pp. 17-25. Baig M. F., Masood A. Natural convection in a two-dimensional differentially heated square enclosure undergoing rotation. Numerical Heat Transfer, Part A, 2001, vol. 40, no. 2, pp. 181-202.

3. Sandaravadivelu K., Kandaswamy P. Nonlinear convection in a rotating square cavity. Acta Me-chanica, 2000, vol. 144, no. 1, pp. 119- 125.

4. Saleh H., Alhashash A. Y. N., Hashim I. Rotation effects on non-Darcy convection in an enclosure filled with porous medium. International Communications in Heat and Mass Transfer, 2013, vol. 43, pp. 105 - 111.

5. Sedelnikov G. A., Busse F. H., Lyubimov D. V. Convection in rotating cubical cavity. European Journal of Mechanics B/Fluids, 2012, vol. 31, pp. 149-157.

6. Fann S., Yang W. J. Convective heat transfer in rotating square channel with oblique cross section.

Computational Mechanics, 1994, vol. 14, no. 5, pp. 513-527.

7. Gershuni G. Z., Zhukhovitskii E. M. Convective stability of incompressible fluids. Jerusalem, Israel: Keter Publishing House, 1976, 330 p.

8. Paskonov V. M., Polezhaev V. I., Chudov L. A. Chislennoe modelirovanie processov teplo- i mas-soobmena (Numerical simulation of heat and mass transfer processes). Moscow: Nauka, 1984, 288 p. (In Russian).

9. Gibanov N. S., Sheremet M. A. Vlijanie formy i razmerov lokal'nogo istochnika jenergii na rezhimy konvektivnogo teploperenosa v kvadratnoj polosti (Effect of shape and sizes of a local heat source on convective heat transfer in a square cavity). Computer Research and Modeling, 2015, vol. 7, no 2, pp. 271-280 (In Russian).

10. Astanina M. S., Sheremet M. A. Simulation of natural convection with variable viscosity in an enclosure with a local heat source. Bulletin of Perm University. Series: Physics, 2015, no. 3 (31), pp. 52-58 (In Russian).

11. Hamady F. J., Lloyd J. R., Yang K. T., Yang H. Q. A study of natural convection in a rotating enclosure. Journal of Heat Transfer, 1994, vol. 116, pp. 136-143.

12. Tso C. P., Jin L. F., Tou S. K. W. Numerical segregation of the effects of body forces in a rotating, differentially heated enclosure. Numerical Heat Transfer, 2013, vol. 51, pp. 85-107.

Просьба ссылаться на эту статью в русскоязычных источниках следующим образом:

Михайленко С. А., Шеремет М. А. Моделирование конвективного теплопереноса во вращающейся замкнутой полости с локальным источником энергии // Вестник Пермского университета. Физика. 2017. № 1 (35). С. 19-25. doi: 10.17072/1994-3598-2017-1-19-25

Please cite this article in English as:

Mikhailenko S. A., Sheremet M. A. Simulation of convective heat transfer inside a rotating enclosure with a local heat source // Bulletin of Perm University. Physics, 2017, no. 1 (35), pp. 19-25. doi: 10.17072/1994-3598-20171-19-25