Численное моделирование свободной конвекции вязкой жидкости в канале с нагретым цилиндром Текст научной статьи по специальности «Математика»

УЧЕНЫЕ ЗАПИСКИ КАЗАНСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА

Том 14 7, кн. 3

Физико-математические пауки

2005

УДК 532-681.3

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ СВОБОДНОЙ КОНВЕКЦИИ ВЯЗКОЙ ЖИДКОСТИ В КАНАЛЕ С НАГРЕТЫМ ЦИЛИНДРОМ

А.Б. Мазо

Аннотация

Предложен новый вариант метода конечных элементов для расчета естественной конвекции в канале с нагревателем. Метод основан па оригинальной формулировке граничных условий для нестационарной задачи Навье Стокса в переменных «функция тока вихрь температура».

Введение

Точения со свободной и смешанной конвекцией широко распространены в природе и тех областях техники, где температурные градиенты настолько значительны. что архимедова сила играет существенную роль в ускорении или замедлении потока жидкости. Классические результаты в математическом описании процессов естественной конвекции получены при исследовании теплового и динамического пограничных слоев около вертикальной нагретой стопки и основаны на автомодельных решениях уравнений вязкого течения и теплопроводности. С развитием вычислительных методов появилась возможность теоретического исследования более сложных процессов конвективного теплообмена на основе совместного решения уравнений Навье Стокса и теилопереноса [1 3]. Большинство результатов численного моделирования в данной области относится к смешанной конвекции в вертикальных каналах либо к естественной конвекции в подогреваемых снизу резервуарах. Вместе с тем. некоторые практически важные процессы свободной конвекции в жидкости около горизонтальных нагретых труб изучены не достаточно: результаты экспериментальных исследований в этой области, так же как и конструирование эффективных численных методов решения соответствующих математических задач, остается актуальной проблемой [4 6].

В настоящей статье предлагается метод расчета естественной конвекции несжимаемого газа (жидкости) в открытом вертикальном канале с помещенным в него горячим цилиндром радиуса Я (см. рис. 1). Температура нагревателя постоянна и равна Т, а газ вне канала и в канале в начальном состоянии имеет температуру-То <ТХ.

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

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

^ = + (1) ОЬ р дх

Рис. 1. К постановке задачи

Оу 1 др

БТ

— =аЛТ, Т(х, у, 0) = То! (3)

, ,г Б д д д

= т = т + ид-х + 1%' (4)

Здесь и, V - компоненты вектора скорости V, р - давление; д - ускорение силы тяжести; в ~ коэффициент температурного расширения, для одноатомных газов в = 1/Т; То - начальная температура и температура за пределами капала; а - температуропроводность, V - кинематическая вязкость жидкости. На стенках канала р, Г^ и та поверхности нагревателя 7 ставится уеловие V = 0 (прилипание), а вход Г;п и выход Г01й расположены достаточно далеко от нагревателя, чтобы положить там V = 0. Граничные условия для температуры следующие:

Гт : Т = То;

7 : Т = Т\; (5)

ОТ _

Гош, Г)0р, Гьо1 • —

Формулировка задачи в преобразованных переменных

Так же, как и в задачах без естественной конвекции [7], от уравнений в естественных переменных и,у,р (1)— (4) можно перейти к уравнениям в преобразованных переменных завихренность ^ — функция тока ф:

ду ди дф дф

-Аф = и] ш=---—, и=—~, V =- — ; (6)

дх ду ду дх

Ои г,дТ

— = иАи; + д/З—. (7)

Если ввести безразмерные переменные по формулам

х, У = * = и, V = и0 = ^/д!ЗЩТг - То)] (8)

К К ио

ф Т — То иК 7 ф

1 = ^—^г; ^ = —; Ф

Т\ — То' ио ' Ки

о

то задача будет определяться двумя критериями подобия: числами Рейнольдса Ие и Пекле Ре, а число Ричардсона Ш в рассматриваемом случае свободной конвекции будет равно единице:

Ре=^, т=9'ЩТ\-То) =1. (9)

V а и0

Проведем элементарные оценки для воды, когда То ~ 20°С Т ~ 100°С, в = = 2 • 10-4К-1 [1], которая окружает трубу радиуса К = 10-2 м. По формулам (8), (9) при V = 10-6 м2 /с, а = 4 • 10-7 м2 /с получаем и0 = 4 • 10- 2 м/с, Ие « 4 • 10 2, Ре « 2 • 103. Данные оценки показывают, что изучаемый класс течений вполне может быть описан в рамках моделей ламинарного течения и теплообмена.

Система определяющих уравнений принимает вид (черта над безразмерными переменными опущена)

л / дг' ди д^ д^ ит

-Аф = и\ Ш=---—, и=—~, V =- — ; (10)

дх ду ду дх

ди ди ди 1 ОТ

т+иЖ:+г% = ^ + ^ (11)

дТ дТ дТ 1

+ = (12)

Граничные условия (5) для температуры, а также для вихря и функции тока теперь можно записать так:

Ь = 0 : Т = 0; ф =0; и = 0;

Г„ : Т = 0; дф/дп = ут = 0; ди/дп = 0;

Гоиг : дТ/дп =0; дф/дп = ут = 0; ди/дп = 0; (13)

ГЬог : дТ/дп = 0; ф = 0; дф/дп = ут = 0;

Г0р : дТ/дп = 0; ф = (; дф/дп = Ут = 0;

7 : Т =1; ф = ф7; дф/дп = ут = 0.

Здесь и далее индексы пит указывают на нормаль и касательную к границе. Начальные условия свидетельствуют, что при Ь = 0 газ в канале — холодный и неподвижный. «Мягкие» условия па входе и выходе моделируют установившееся течение, параллельное оси канала: однако предполагаем, что течение всегда направлено от Г;п к Г0^, и во входное сечение поступает холодный газ (Т = 0), тогда как на выходе из капала температура практически ие меняется. Стопки теплоизолированы и так же, как поверхность нагревателя 7, являются линиями тока.

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

(

ф7 в каждый момент времени будем использовать подход [8], согласно которому для произвольной гармонической функции п справедливо равенство

J^cT-J ф^-сТ+ =0, = У ил]сЮ. (14)

г

г

и

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

Здесь Г - объединение всех границ, Б - область расчета. В силу граничных условий (13) первый интеграл в (14) равен нулю, а второй интеграл можно представить в виде

ф^с1Г + <2 [ ^с/Г + ф^ [ (15)

Г;„ + Гоц1

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

Ащ = 0; Г,ор : = 1; Г;п, ГоиЬ 7 : ^ = 0; ГЬо| : т = (1С)

Тогда из (14) (16) получаем

д = п(т,ш)/\ти>р (17)

| Г| Г

дп2 дп2

Д?у2 =0; 7 : = 1; Г;п, ГоиЬ Г,ор : = 0; ГЬо| : т = 0, (18)

то по аналогии получаем искомую формулу для значения функции тока на границе цилиндра

ф.7 = П(п2,ш)/ \7\. (19)

Вспомогательная функция ц2(х,у) показана на рис. 2.

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

Ащ = 0; 7 : ^ = 1; Г;п, Гоц, : ^ = 0; ГЬо|, Г,ор : щ = 0, дп дп

решение которой удовлетворяет условию

/ / /^ = 2 Г =

.] дп ] дп ] дп ] дп ] дп ] дп

Г ГЬо(; Г4ор 7 Гор 7

Отсюда следует, что

Г

Подставив это выражение в (14). (15). получим выражение для искомой константы

С3

^=2+ПтР" ( 5

При стационарном или симметричном течении завихренность антисимметрична относительно оси канала (т. е. ш(х,у) = —ш(х, —у) ), поэтому интеграл в (20) будет равен нулю; отклонение «центральной» линии тока ф7 от Q/2 можно ожидать лишь при периодическом режиме обтекания (вихревая дорожка Кармана).

Метод численного решения

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

• (10) для функции тока ф с граничными условиями (13), по которой дифференцированием определялись скорости и,у, а также вычислялось распределение завихренности на стенках с помощью уравнения вида ш = —△ф [8];

• (12) для температуры Т, то которой вычислялся «источник» дТ/ду в уравнении для ш;

• (11) для завихренности ш;

• с помощью соотношений (17) и (20) уточнялись значения полного расхода Q и функции тока на цилиндре ф7.

При записи конечноэлементных аппроксимаций конвективные слагаемые брались как полусумма значений на нижнем временном слое и на предыдущей итерации, для остальных слагаемых в уравнениях использовалась неявная схема. Решение соответствующих систем алгебраических уравнений находилось прямым методом, причем разложения Холецкого строились один раз, поскольку матрицы систем уравнений но зависят от времени. Также одни раз, до временного цикла, решались задачи (16), (18) для вспомогательных гармонических функций щ, 'ц2-

В качестве иллюстрации работы предложенного численного метода приведем некоторые результаты моделирования естественной конвекции в канале с размерами —20 < х < 60, —10 < у < 10 при И.е = 100, Ре = 100. Как показывают расчеты, газ, нагреваясь у поверхности цилиндра, поднимается вверх в виде тонкой струн непосредственно над ним. Эта струя вызывает вторичное течение газа в канале, которое имеет вид пары крупных вихрей с замкнутыми линиями тока, которые медленно перемещаются вверх (рис. 3). Поскольку единственным источником движения является поперечный градиент температуры дТ/ду в уравнении (11), течение локализовано в окрестности восходящей струи, а в верхней и особенно в нижней частях канала среда остается неподвижной до тех пор, пока туда не распространится тепловое возмущение.

-10 0 10 X

Рис. 3. Липли тока при свободной конвекции

10

-10

_1_

_1_

0 10

Рис. 4. Температурное поле над цилиндрическим нагревателем

Температурное поле имеет характерный вид «гриба», поднимающегося над нагревателем. как это показано на рис. 4. Полученная картина свободной конвекции вполне согласуется с численными результатами [4]. полученными при моделировании аналогичных процессов около горизонтальной трубы эллиптического сечения.

5

0

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

Summary

А. В. Mazu, Numerical simulation of natural convection of viscous fluid in a channel with a hot cylinder.

A new variation of finite element method for calculation of natural convection in a channel containing a heater is offered. The method is based 011 the original formulation of the boundary conditions for unsteady Navier Stokes problems in variables "stream function vort.icit.y temperature".

Литература

1. Себиси Т., Бредшоу П. Конвективный теплообмен. М.: Мир. 1987. 592 с.

2. Математическое моделирование конвективного теплообмена па основе уравнений На-вье Стокса / Под ред. B.C. Авдуевского. М.: Наука, 1987. 272 с.

3. Ши Д. Числеппые методы в задачах теплообмена. М.: Мир. 1988. 544 с.

4. Mahfouz F.M., Kocabiyik S. Transient numerical simulation of buoyancy driven flow adjacent to an elliptic tube // Int. J. of Heat and Fluid Flow. 2003. No 24. P. 864-873.

5. Srivastava A., Dutta S., Panigrahi P.K., Muralidhar K. Laser sclilieren measurement of vertical flow past a heated horizontal circular cylinder // Int. Communications in Heat and Mass Transfer. 2005. No 32. P. 520-528.

6. Corcione M. Correlating equations for free convection heat transfer from horizontal isothermal cylinders set in a vertical array // Int. J. of Heat and Mass Transfer. 2005. No 48. P. 3660-3673.

7. Флапче.р К. Вычислительная гидродинамика. Ч. 2. М.: Мир. 1991. 552 с.

8. Мазо А.Б., Даутоо Р.З. О граничных условиях для уравнений Навье Стокса в переменных функция тока - завихренность при моделировании обтекания системы тел // ИФЖ. 2005. Т. 78, Л» 2. С. 75 79.

Поступила в редакцию 14.10.05

Мазо Александр Венцианович доктор физико-математических паук, профессор кафедры аэрогидромехапики Казанского государственного университета. Е-шаП: amazoQksu.ru