CC BY
42
УДК 517.958:532.5:536.25
И.А. Ермолаев, А.С. Шаповалов, В.Б. Байбурин
МОДЕЛИРОВАНИЕ СМЕШАННОЙ ТЕРМОГРАВИТАЦИОННОЙ КОНВЕКЦИИ В ОБЛАСТИ С НЕРЕГУЛЯРНОЙ ГЕОМЕТРИЕЙ И НЕОДНОРОДНЫМИ УСЛОВИЯМИ
НА ГРАНИЦАХ
Методом конечных элементов Галёркина решается задача смешанной термогравитационной конвекции в области сложной формы с неоднородными условиями на
границе. Задача усложняется наличием выходной границы, условия на которой не определены в явном виде и вычисляются на основе эмпирических моделей.
Математическое моделирование, метод конечных элементов, смешанная термогравитационная конвекция
I.A. Ermolaev, A.S. Shapovalov, V.B. Baiburin
MODELING OF MIXED THERMOGRAVITATIONAL CONVECTION IN THE AREA WITH THE IRREGULARLY GEOMETRY AND NON-UNIFORM BOUNDARIES CONDITIONS
The results of the finite elements received using Galerkin numerical simulation method of mixed thermogravitational convection in the area with the irregularly geometry and non-uniform boundaries conditions are represented.
Mathematical simulation, finite elements method, mixed thermogravitational convection
Введение. Известно значительное количество работ, посвященных численному моделированию естественной термоконвекциии в замкнутых областях [1, 2]. Чаще всего объектом исследования является конвекция в областях регулярной формы - квадратной или прямоугольной. Нерегулярность геометрии и неоднородность условий на границах существенно влияют на свойства конвективных течений и осложняет их моделирование [3-6].
Вместе с тем наиболее сложными являются задачи, связанные с анализом динамики комбинированных конвективных течений в сочетании с нерегулярной геометрией и неоднородными условиями на границах. Одним из видов комбинированных течений является смешанная термоконвекция, обладающая двумя взаимодействующими механизмами неустойчивости: конвективным и гидродинамическим. Численные исследования таких течений проводились в основном для полостей с относительно простой геометрией [7, 8].
В настоящей работе методом конечных элементов моделируется смешанное термоконвективное течение в двумерной области сложной формы с неоднородными условиями на границе. Расчетная область представляет собой схематически (рис. 1) полуоткрытый корпус некоторого типа радиоэлектронных устройств с внешним охлаждением и локальным источником тепла в виде нагретого выступа. Часть верхней стенки является выходной границей. На выступе задана плотность теплового потока, температура остальных границ фиксирована и равна температуре окружающей среды. Вынужденное охлаждение осуществляется через подводящие каналы внизу корпуса. Подобная геометрия представляется также типичной для ряда задач обогрева и вентиляции помещений специального назначения и др.
Постановка задачи, модель, метод решения. Рассматривалось двухмерное смешан- ______________________t...t...,t..t....t ____________
ное конвективное течение вязкой, термически сжимаемой жидкости (газа), для которой справедливо приближение Буссинеска. В начальный момент времени поле температур считалось однородным, температура равна температуре окружающей среды.
При решении задачи использовались нестационарные двухмерные уравнения Буссинеска [1]. В качестве масштабов расстояния, времени, скорости и температуры были выбраны
Н, H2N, v/H, qoH/X. Безразмерные переменные были равны соответственно X = x/H, Y = y/H, Рис. 1. Схема расчетной области
т = Vt/H2, U = uH/v, V = UH/V, 0 = ХЪ/qoH. Здесь
x, y - координаты, t - время, V - коэффициент кинематической вязкости, u, U - составляющие скорости в проекции на оси x, y соответственно, Ф = T - T0, T0 - начальная температура; X - коэффициент теплопроводности, q0 - плотность теплового потока. Безразмерные уравнения Буссинеска в переменных «вихрь скорости - функция тока - температура» записывались как
х ft! Ttt
дю дш Эю дш Эю ^ Э0
— + —--------------—— = Лю-вг —,
дт дУ дХ дХ дУ у дХ
Лш = ю,
д0 дш д0 дш д0 1 л
— + —----------------------------— = —Л0.
дт дУ дХ дХ дУ Рг
(1)
(2)
(3)
где ю, у - вихрь скорости, функция тока соответственно; вгу = gУвqoH4/Xv2 - число Грасгофа (вгх = 0); Рг = у/% - число Прандтля; gy - составляющая ускорения силы тяжести в проекции на ось у = 0); в - температурный коэффициент объемного расширения; % - коэффициент температуропроводности. Граничные условия для системы (1)-(3) имели вид (рис. 1): на твердых непроницаемых границах
а д0
■ 0,---= 1; на
дЫ
Ш = = 0, 0 = 0ср; на входных границах ш = Ре, ю = 0, 0 = 0ср; на нагревателе ш =
дЫ дЫ
выходной границе ш ■
дш дю п д0 п. та. „ „
—— = 0, ------= 0, -----= Бі, где Ві - критерий Био, вычисляемый по эмпириче-
дЫ дЫ дЫ
ским формулам [9] для условий, соответствующих внешней, по отношению к задаче, термоконвекции.
В начальный момент времени ю(Х, У, 0) = у(^, У, 0) = 0(Х, У, 0) = 0. Вихрь скорости на стенках определялся по формуле Вудса [10]. Задача решалась методом конечных элементов Галеркина [11]. Температура, вихрь скорости и функция тока аппроксимировались линейной комбинацией не зависящих от времени функций формы на линейных треугольных конечных элементах. Для временной аппроксимации использовалась неявная двухслойная схема.
Уравнения (1)-(3) решались последовательно, каждый временной шаг начинался с вычисления поля температуры (3), затем определялись граничные условия для вихря скорости, и решалось уравнение (1), далее поле вихря скорости корректировалось и определялось поле функции тока (2). Расчеты проводились по конечно-элементной программе ЕРЕМА1_ТБРМО, реализующей данный алгоритм [12]. Стационарные решения были получены методом установления, критерием установления являлось неравенство
І0к+1 -0*1 + |юк+1 -юЧ + |ш
т т т т т
-шт < ерє.
— і
где 0т, ют, шт - экстремальные значения температуры, вихря скорости и функции тока. Индекс к -номер шага по времени, величина ерз изменялась в интервале 10-5^10-6. Шаг по времени 10-3. Расчеты проводились на неравномерной сетке. Проверка на более подробной сетке показала, что относительное изменение максимума температуры не превосходит 1%.
Применение процедуры метода конечных элементов Галеркина для решения системы (1)-(3) (ортогонализация невязки относительно базисных функций, определенных локально - на конечных элементах) приводит к системе из трех матричных уравнений вида
д [с ]{Ф} + [к ]{ф}+М=0,
дґ
где [С] - матрица демпфирования, [.К] - матрица жесткости, {Ф} - вектор узловых значений, {^} -вектор нагрузки. Элементы матриц [С], [К] и вектора {^}: для вихря скорости
4 Лт I 1 7
дш лг дЫ7 ,дш лг дЫ7 „дЫ; дЫ7 „дЫ1 дЫ7
*7 = Л^) еЫ^ - (^) N
дУ
дх дх ' ду
д0
дх
-V-
дх дх
-V-
дУ дУ
іІБ ,
ґі =-[Мі (+ЛТ {ыю-1^ ;
для функции тока
*7 =- {
Сі] = 0,
( ЩдЫ± +ЭЫ1Щ ^
дх дх дУ дУ
/і = { Ы,юе йБ ;
йБ ,
Б
Б
Б
для температуры
Рис. 2. Мгновенные значения изотерм (а) Рис. 3. Изменения максимума безразмерной
и линий тока (Ь) температуры во времени
Здесь Су - элемент матрицы демпфирования, Ат - шаг по времени, N N - функции формы, Бе - площадь конечного элемента, е - индекс конечного элемента, ку - элемент матрицы жесткости, / - элемент вектора нагрузки; п - индекс шага по времени. Величины с нижним индексом е усредняются по каждому конечному элементу, величины с верхним индексом п-1 относятся к предыдущему временному шагу. Матрицы [С] и [й], а также вектор формируются суммированием по всем конечным элементам. При решении системы алгебраических уравнений применялся метод Гаусса для ленточных матриц.
Результаты. Численные расчеты проводились для среды, характеризуемой критерием Пранд-тля Рг = 0,7, что соответствует воздуху, в диапазоне параметров Грасгофа 5106 < вг < 5109 и Рейнольдса 10 < Ие < 410 . Данные диапазоны соответствуют рабочим режимам прибора.
В таких условиях реализуется заведомо турбулентный режим конвективного течения воздуха, характеризуемый колебаниями температуры, вихря скорости и функции тока вблизи некоторых средних величин. Симметричное стационарное решение неустойчиво при заданной совокупности параметров в данной геометрии, поэтому устойчивыми являются колебания интенсивности циркуляции и температуры, сказывающиеся на структуре потока.
На рис. 1 представлены мгновенные значения поля температур (а) и поля течения (Ь) для случая вг = 210 и Ие = 4-10 . Видно, что нагретый вблизи поверхности выступа воздух не циркулирует внутри полости, а поднимается к выходной границе. Кроме того, видна также циркуляция холодного воздуха, поступающего в полость извне через верхнюю границу.
На рис. 2 показаны колебания максимума безразмерной температуры. Эти колебания устанавливаются примерно через 6.5 единиц безразмерного времени. Отклонения температуры невелики -менее одной безразмерной единицы. Тем не менее эти колебания связаны с заметными колебаниями потока воздуха (рис. 1 Ь).
Анализ поля температур на рис. 1а показывает, что значительная часть температурного поля расположена в пределах одного шага изотерм, т.е. значительные массы воздуха вблизи боковых стенок практически изотермические и имеют температуру, близкую к температуре стенок.
Градиент температуры сосредоточен в непосредственной близи поверхности источника тепла. Максимум температуры находится в верхней части источника. Таким образом, область наибольшей тепловой нагрузки расположена только в непосредственной близи от источника тепла, на расстоянии, меньшем половины его характерного размера сбоку, либо порядка его характерного размера сверху. Расположение термочувствительных элементов и узлов возможно вдоль боковых и нижних стенок, учитывая указанные расстояния.
Из анализа поля функции тока на рис. 1Ь также следует, что значительные области вблизи боковых стенок представляют собой зоны, где конвективное течение практически отсутствует. Таким образом, любое заполнение этих областей конструктивными элементами слабо скажется на конвективном течении и, как следствие, на поле температур.
Заключение. Таким образом, полученные результаты моделирования структур комбинированного конвективного течения и теплообмена в области со сложной геометрией, имитирующей корпус РЭА, позволяют сформулировать рекомендаций по его конструктивному заполнению и расположению термочувствительных элементов. Кроме того, результаты моделирования позволяют дать сравнительные характеристики различных режимов теплообмена, связанных с вариацией параметров вг и Ие.
ЛИТЕРАТУРА
1. Гершуни Г.З. Конвективная устойчивость несжимаемой жидкости / Г.З. Гершуни, Е.М. Жу-ховицкий. М.: Наука, 1972. 392 с.
2. Математическое моделирование конвективного тепломассообмена на основе уравнений Навье-Стокса / В.И. Полежаев, А.В. Бунэ, Н.А. Верезуб и др. М.: Наука, 1987. 271 с.
3. Karyakin Yu.E. Transient natural convection in prismatic enclosures of arbitrary cross-section / Yu.E. Karyakin // Heat Mass Transfer. 1988, Vol. 31. No. 6. P. 1177-1187.
4. Ramos E. Natural convection in a two-dimensional square loop / E. Ramos, A. Castrejon // Heat Mass Transfer. 1990. Vol. 33. No. 5. Р. 917-930.
5. Ермолаев И. А. Моделирование естественной термогравитационной конвекции в горизонтальном канале с сечением нерегулярной формы / И.А. Ермолаев, А.И. Жбанов, В.С. Кошелев // Инженерно-физический журнал. 2003. Т. 76. № 4. С. 134-137.
6. Acharya S. Heat transfer due to buoyancy in a partially divided square / S. Acharya, R. Jetli // Heat Mass Transfer. 1990. Vol. 33. No. 5. P. 931-942.
7. Yerkes K.L. An experimental and numerical simulation of mixed convection in large baffled rectangular chambers / K.L. Yerkes, A. Faghri // Heat Mass Transfer. 1991. Vol. 34. № 6. Р. 1525-1542.
8. Бернер С. Обтекание перегородок в теплообменнике / С. Бернер, Ф. Дарст, Д. Макэлигот // Теплопередача. 1984. Т. 106. №4. С. 48-55.
9. Мартыненко О.Г. Свободно-конвективный тепло- и массообмен: библиогр. указатель (17971981) / О.Г. Мартыненко, Ю.А. Соковишин. Минск: АН БССР, ИТМО им. А.В. Лыкова, 1982. Ч. 1. 390 с., 1983. Ч. 2. 418 с.
10. Тарунин Е.Л. Вычислительный эксперимент в задачах свободной конвекции / Е.Л. Тару-нин. Иркутск: Изд-во Иркут. ун-та, 1990. 225 с.
11. Зенкевич О. Конечные элементы и аппроксимация / О. Зенкевич, К. Морган. М.: Мир, 1986. 288 с.
12. Ермолаев И.А. Решение двумерной нестационарной задачи тепло- и массопереноса методом конечных элементов / И.А. Ермолаев, А.И. Жбанов, В.С. Кошелев // Вопросы прикладной физики. Саратов: Изд-во СГУ, 2002. Вып. 8. С. 60.
Ермолаев Игорь Анатольевич -
кандидат физико-математических наук, ведущий инженер кафедры прикладной физики Саратовского государственного университета имени Н.Г. Чернышевского
Шаповалов Александр Степанович -
доктор физико-математических наук, профессор, заведующий кафедрой прикладной физики Саратовского государственного университета имени Н.Г. Чернышевского
Igor A. Ermolaev -
PhD, Associate Professor Department of Applied Physics Chernyshevsky Saratov State University
Alexander S. Shapovalov -
Dr. Sc., Professor
Head: Department of Applied Physics Chernyshevsky Saratov State University
Байбурин Вил Бариевич -
доктор физико-математических наук, профессор заведующий кафедрой «Программное обеспечение вычислительной техники и автоматизированных систем»
Саратовского государственного технического университета имени Гагарина Ю.А.
Vil B. Baiburin -
Dr. Sc., Professor
Head: Department of Software of the Computing
Machinery and Automated Systems
Yu. Gagarin Saratov State Technical University
Статья поступила в редакцию 21.10.11, принята к опубликованию 15.11.11
