CC BY
10
DOI: 10.15514/ISPRAS-2022-34(5)-14
Моделирование сопряженного теплообмена в микроканалах в среде OpenFOAM
Е.С. Байметова, ORCID: 0000-0002-4534-0936 <baimetova.e.s@gmail.com>
М.Е. Хвалько, ORCID: 0000-0001-8524-4749 <izhgazproekt@mail.ru> А.Ю. Армянин, ORCID: 0000-0001-7541-7239 <armynin-izhevsk@mail.ru> Ижевский государственный технический университет имени М.Т. Калашникова, 426069, Россия, г. Ижевск, ул. Студенческая, д. 7
Аннотация. В работе с использованием инструментов пакета OpenFOAM проведено численное моделирование вынужденного конвективного теплообмена в кремниевом микроканальном радиаторе. В качестве теплоносителя использовалась однофазная жидкость - вода. Модель микроканального радиатора представлена в виде кремниевой подложки длиной 10 мм, с прямоугольными микроканалами шириной 57 мкм и глубиной 180 мкм. расположенными по всей длине радиатора. Выполнен сравнительный анализ в виде кроссплатформенной верификации полученных численных результатов с данными сторонних авторов. Анализ полученных данных показал хорошую сходимость результатов исследования и возможность применения пакета OpenFOAM в качестве расчетной среды для численного моделирования физических процессов, протекающих в канальных радиаторах.
Ключевые слова: численное моделирование; микроканал; теплообмен; сравнительный анализ
Для цитирования: Байметова Е.С., Хвалько М.Е., Армянин А.Ю. Моделирование сопряженного теплообмена в микроканалах в среде OpenFOAM. Труды ИСП РАН, том 34, вып. 5, 2022 г., стр. 205214. DOI: 10.15514/ISPRAS-2022-34(5)-14
Modeling of coupled heat transfer in microchannels in OpenFOAM
E.S. Baymetova, ORCID: 0000-0002-4534-0936 <baimetova.e.s@gmail.com>
M.E. Hval'ko, ORCID: 0000-0001-8524-4749 <izhgazproekt@mail.ru> A.Yu. Armyanin, ORCID: 0000-0001-7541-7239 <armynin-izhevsk@mail.ru> Kalashnikov Izhevsk State Technical University 7, Studencheskaya st., Izhevsk, 426069, Russia
Abstract. In this paper, a numerical simulation of forced convective heat transfer in a silicon microchannel heat sink has been performed using the OpenFOAM tools. A single-phase fluid - water - was used as a heat transfer medium. The model of microchannel heat sink is represented as silicon substrate with length 10 mm, with rectangular microchannels 57 microns wide and 180 microns deep located along the full length of the heat sink. A comparative analysis in the form of cross-platform verification of the numerical results obtained with data from third-party authors was performed. The analysis of the obtained data has shown a good convergence of the study results and the possibility of using the OpenFOAM package as a computational environment for the numerical simulation of the physical processes occurring in channel radiators.
Keywords: numerical simulation; microchannel; heat exchange; comparative analysis.
For citation: Baymetova E.S., Hval'ko M.E., Armyanin A.Yu. Modeling of coupled heat transfer in microchannels in OpenFOAM. Trudy ISP RAN/Proc. ISP RAS, vol. 34, issue 5, 2022. pp. 205-214 (in Russian). DOI: 10.15514/ISPRAS-2022-34(5)-14
1. Введение
В настоящее время применение микромасштабных охлаждающих устройств является очень актуальной задачей теплоотвода и поддержания оптимальной температуры в микромеханических системах, для охлаждения электроники, в электрических цепях и других устройствах, которые подвергаются высоким тепловым нагрузкам. Микроканальные радиаторы позволяют обеспечивать высокие показатели коэффициента теплопередачи при течении жидкости в небольших объемах. Отличительная особенность микроканальных радиаторов заключается в том, что хладагент протекает через большое количество микроканалов, которые зачастую бывают связаны между собой сложным образом, что не позволяет разбивать теплообменник на отдельные составляющие и расчет приходится проводить целиком, что ведет за собой большую трату вычислительных ресурсов. В работе [1] для решения задачи создания эффективного численного алгоритма для моделирования сопряженного теплообмена в микротеплообменниках предлагается использовать гибридный алгоритм, основанный на сочетании методов ТГЦ и пространственных методов вычислительной гидродинамики (CFD). В большинстве случаев решение гибридных задач осуществляется путем разделения пространственной и сетевой частей на разные расчетные модули и последующей их сшивки. Такой подход может быть неустойчивым и привести к значительному увеличению времени расчета. В случае с сопряженным теплообменом (твердое тело-жидкость) реализация гибридного подхода тем более затруднительна, так как сетевые элементы непосредственно включены в пространственную область. Авторами статьи [1] предложен метод, в котором подобная гибридная задача решается на основе общего уравнения на поправку давления. В данной статье представлены результаты тестирования методики на задачах гидродинамики и теплообмена в микроканалах.
В работе [2] рассматривается процесс прогрева стенки микроканала, для решения управляющих уравнений авторами был разработан конечно-разностный численный код с трехдиагональным матричным алгоритмом (TDMA), на основе которого проводился верификационный анализ при помощи решения задачи в среде OpenFOAM с использованием оригинальной схемы и применением исходных теплофизических характеристик исследуемых сред.
Численным расчетам в области микроканального моделирования посвящен ряд работ. В [3] рассматриваются процессы, которые необходимо учитывать при уменьшении диаметра каналов, включая проблемы с распределением потока в одиночном, параллельном и разделенном потоках, нестабильность потока в параллельных каналах, процессы зарождения и эффекты проводимости стенки.
В работе [4] коллектив авторов решал сопряженную задачу теплопередачи для микроканальных радиаторов, с целью получить подробные пространственные распределения температуры поперечного сечения радиатора по длине каналов в твердой подложке и жидкости, приведен расчетный алгоритм для выбора размеров теплообменников. В [5] показан эксперимент по измерению коэффициента трения ламинарного потока деионизированной воды в гладких кремниевых микроканалах трапециевидного сечения с гидравлическими диаметрами в диапазоне 25,9-291,0 мкм. Доказано, что на постоянную трения этих микроканалов большое влияние оказывает соотношение сторон поперечного сечения микроканала. Установлено, что экспериментальные данные хорошо согласуются с существующим аналитическим решением для несжимаемого, полностью развитого ламинарного течения при граничном условии без скольжения. Подтверждено, что уравнения Навье-Стокса по-прежнему справедливы для ламинарного течения деионизированной воды в гладких кремниевых микроканалах с гидравлическим диаметром 25,9 мкм. Применение методов численного моделирования для исследования и верификации физических процессов и результатов вычислений в технических устройствах является
актуальной задачей. Существует ряд пакетов прикладных программ для решения задач МСС, такие как SolidWorks, Ansys, FlowVision, OpenFOAM и др. Из представленного списка пакет OpenFOAM является полноценным, завершенным, свободно распространяемым программным продуктом с открытым исходным кодом и набором библиотек, на основе которых можно разработать полноценные численные алгоритмы для решения задач МСС. Зачастую проведение натурных экспериментов невозможно или экономически не целесообразно, в связи с этим применяются численные эксперименты, в результате которых можно получить достаточно физичную картину процессов, протекающих в исследуемой задаче. Численный эксперимент должен быть подкреплен верификацией с натурным экспериментом или кроссверифицирован.
Целью работы является проведение кроссплатформенной верификации результатов численного эксперимента, проведенного в OpenFOAM, с результатами исследований представленных в работе [2].
2. Метод решения
Для решения задачи сопряженного теплообмена используется решатель chtMultiRegionSimpleFoam, который является объединением программ расчета тепломассопереноса в текучих средах buoyantSimpleFoam и процессов теплопроводности в твердом теле на основе уравнения Лапласа laplacianFoam. Основные уравнения сохранения, решаемые с помощью метода контрольных объемов можно записать в следующей форме [68]:
V • (ри) = 0, (1)
V • (рии) = -Vp + pg + V • (2^effD(u)) - V (^VeffQV • и)), (2)
где р - плотность среды, и( и, v, w) - вектор скорости, р - статическое давление, g -ускорение свободного падения, ¡xeff - эффективная вязкость, представляющая собой сумму молекулярной вязкости (динамического коэффициента вязкости жидкости) ц и турбулентной вязкости pvt, D(u) - тензор скорости деформации, определяющийся как D(u) = -(Vu + (Vu)T). 2
Уравнение энергии может быть выражено с использованием внутренней энергии е следующим образом:
V • (рие) + V • (риК) + V • (ри) = V • (aeffVe) + ри^ g, (3)
где К = \и\2/2 - кинетическая энергия на единицу массы, aeff - эффективная температуропроводность: aeff = pvt/Prt + ц/Рг = pvt/Prt + к/ср, где к -теплопроводность, ср - удельная теплоемкость, Рг- число Прандтля, Prt- турбулентное число Прандтля.
Кондуктивный теплообмен через твердое тело может быть определен с помощью решателя laplacianFoam, для которого основное уравнение записывается следующим образом:
дТ
= V • (aVT), (4)
где а- коэффициент термической диффузии а = к/(рср), р - плотность тела, ср - удельная теплоемкость твердого тела. Когда решается стационарная задача, слагаемое в левой части уравнения (4), зависящее от времени не учитывается.
Алгоритм решения задачи в виде блок-схемы приведен на рис. 1. С помощью утилиты построения сетки blockMesh создается сетка для общего домена. Затем, проводится разделение областей, занятых жидкостью и твердым телом с использованием утилит topoSet и splitMeshRegions. На этом этапе проводится определение границ контакта жидкости и
твердого тела, создаются отдельные директории в проекте, позволяющие задавать разные физические свойства веществ, начальные и граничные условия.
Рис. 1. Ход решения задачи Fig. 1. Work flow of task
3. Постановка задачи
Объектом исследования являются процессы сопряженного теплообмена в системе микроканалов кремниевого радиатора. По системе параллельных микроканалов, вытравленных в кремниевой пластине, которая примыкает к нагреваемому устройству, пропускается охлажденная жидкость (рис. 2а). Задача решается для единичного канала с условиями симметрии по боковым стенкам. Расчетная область с размерами приведена на рис. 2б. Геометрия радиатора, расположение и размеры каналов соответствуют работе [2].
Рис. 2. Геометрическая постановка задачи: (а) схема микроканального радиатора, (б) схема расчетной области Fig. 2. Geometric problem description: (а) microchannel radiator scheme, (b) calculation domain scheme
Геометрические размеры радиатора, представленного на рис. 2: длина Lx, высота Ly, ширина Lz равны 10 мм, 0.9 мм и 0.1 мм соответственно. Высота каналов охлаждения H равна 0.18 мм, ширина W = 0.057 мм, расстояние между каналами составляет I = 0.043 мм. Корпус радиатора выполнен из кремния с плотностью 2330 кг/м3, теплоемкостью 714 Дж/(кг*К), теплопроводностью 149 Вт/(м*К). В качестве охлаждающей среды принята вода с плотностью 1000 кг/м3, теплоемкостью 4183 Дж/(кг*К), молекулярной вязкостью 1.004*10" 3 м2/с, число Прандтля в расчетах бралось равным 7.02.
Расчетная область (рис.3) включает границу подвода тепла (1), границы симметрии (2), адиабатические непроницаемые стенки (3), границы контакта жидкости и твердого тела (4), границы подвода (5) и отвода (6) жидкости.
Для твердого тела:
1: Y = 0, qw = 90 Вт/см2;
2: Z = 0, Z = Lz, 0;
z' дп
дТ
3: Y= Ly, X = 0, X = Lx (за исключением границ 5 и 6), — = 0.
Граница контакта жидкости с твердым телом:
4: 0<Х <LX, Y = Ly/2 - H, Y = Ly/2, Z = 1/2Z = l/2 + W, и = v = w = 0, -Àw ^ =
—Af-^,TW = Tf, где индексы «/» и «w» относятся к жидкости и твердому телу, соответственно.
Для жидкости:
5: X = 0, Ly/2 -H<Y < Ly/2, р = 50 кПа, Т = 293 К, 6: X = Lx, Ly/2 -H <Y < Ly/2, p = 0 Па.
Рис. 3. Граничные условия Fig. 3. Boundary conditions Жидкость двигается по каналу под действием перепада давления 50 кПа. При таком перепаде согласно [2] поток остается ламинарным и характеризуется числом Рейнольдса равным 144, которое определяется по эквивалентному диаметру D и средней скорости жидкости. Для
данного канала значение D, определяемое по формуле D =
- равно 0.8658 мкм.
2(№+Н) -
Ввиду простоты геометрии расчетной области и регионов, выделенных под жидкость и твердое тело, сетка состояла из множества параллелепипедов. На объем жидкости пришлось 90 тыс. ячеек, на твердое тело - 720 тыс. ячеек. Внутри жидкого региона ячейки были одинакового размера. Сетка твердого домена строилась так, чтобы на границах сопряжения сетки были полностью согласованными.
4. Анализ результатов численного моделирования
В результате численного моделирования получены распределения полей физических величин в рабочих объемах микроканалов радиатора. Так, на рис. 4 приведены линии равных скоростей в выходном поперечном сечении микроканала. В целом, наблюдается установившееся равномерное течение жидкости, подчиняющееся классическим представлениям о течении жидкости в прямолинейных каналах прямоугольной формы [9]. Видно хорошее качественное соответствие изотах, полученных с использованием матмодели (1)-(4) и средств ОреиРОАМ (черные линии на рис. 4) с результатами численного моделирования, представленными в работе [2] (цветные линии). Средняя скорость воды на выходе равняется 1.1026 м/с, а максимальная скорость составляет 1.97985 м/с, что также соответствует результатам работы [2]. Рассчитанное по данным значениям число Рейнольдса равно 95.548 и 171.415 для средней и максимальной скорости. Расхождение с экспериментальными и теоретическими значениями чисел Рейнольдса из работы [2] не превышает 0.5%.
Распределение линий равных температур в поперечном выходном сечении микроканала представлено на рис. 5. Видна некоторая асимметрия изотерм со смещением локального минимума вверх по поперечной координате, что обусловлено разностью температур между нижней и верхней стенками радиатора. При этом полученные с использованием сформированной математической модели (1)-(3) и средств ОреиРОАМ (черные линии на рис. 5) хорошо согласуются с результатами [2] в центральной зоне и на периферии вблизи верхней стенки канала. В работе [2] отмечено локальное повышение температуры жидкости в углах
канала, что также показано и на рис. 5. Однако, характер, а именно угол наклона, изотерм, полученных в результате используемого в настоящей работе подхода отличен от изотрем, представленных в работе [2]. В частности, наблюдается более интенсивный отвод тепла в стенки радиатора, что может быть связано с расхождением в задании физических свойств кремния. Максимальный прогрев воды в опорной работе составляет 32°С. В текущих расчетах максимальная средняя температура воды равна 31.08°С. В этом случае относительная погрешность составляет 2.8%.
Рис. 4. Сопоставление распределения изотах на выходе из канала: черные линии - численный расчет в
OpenFOAM; цветные - расчет [2] Fig. 4. Comparison of the distribution of isotopes at the outlet of the channel: black lines - numerical calculation in OpenFOAM; color lines - results from work [2]
Рис. 5. Сопоставление распределения изотерм в канале (черные линии - численный расчет в
OpenFOAM; цветные - расчет [2]) Fig. 5. Comparison of the distribution of isotherms in the channel (black lines - numerical calculation in
OpenFOAM; color lines - calculation [2])
На рис. 6 приведено распределение температуры по продольным сечениям расчетной области. Рис. 6а показывает температурное поле на грани симметрии в плоскости Z = 0 или Z = Lz. Рис. 6б соответствует средней плоскости Z = Lz/2 и иллюстрирует распределение не только в твердой стенке, но и в жидкости. Охлаждение кремниевого радиатора происходит неравномерно. Наиболее интенсивно охлаждается первая треть объема, затем интенсивность снижается ввиду прогрева жидкости. Движение жидкости влияет на температурное поле радиатора, это видно по наклону линий изотерм в нижней и верхней части радиатора. Наиболее сильный градиент соответствует первой трети канала, что говорит о сильном отводе тепла от нижней нагреваемой части радиатора.
Рис. 6. Распределение температуры в продольных сечениях радиатора: (а) плоскость симметрии, Z = 0, Z = Lz; (б) центральная плоскость, Z = Lz/2 Fig. 6. Temperature field in the radiator cross section: (a) symmetry plane at Z = 0, Z = Lz; (b) central plane at Z = Lz/2 Изменение температуры жидкости и в стенках радиатора по каналу показано на графиках рис. 7. Три линии на рис. 7а показывают температуру воды вдоль нижней стенки, на середине канала и на верхней стенке канала. Вначале канала температура воды выше на верхней стенке, по всей видимости, сказывается нагрев радиатора между каналами. В середине канала в точке Z = 0.005 м верхняя и нижняя стенки имеют одинаковую температуру, далее нижняя стенка начинает прогреваться быстрее. Однако разница температур к концу канала достигает всего 1 К. Температура в ядре потока всегда ниже, чем температура жидкости на стенках, максимальная разница в центре канала составляет 3 К. На рис. 7б видно, что верхняя, относительно канала с водой, часть радиатора имеет более низкую температуру по сравнению с нижней частью. По длине канала эта разница несущественно, но увеличивается. В начале канала разница составляет 2 К, в конце - 3 К.
(а) (б)
Рис. 7. График изменения температуры вдоль радиатора для: (а) воды:1 - нижняя стенка, 2 - середина канала, 3 - верхняя стенка, (б) твердого тела: 1 - верхняя часть радиатора (над каналом), 2 - нижняя часть (под каналом)
Fig. 7. Temperature plots along radiator: (a) for water in channel: 1 - bottom wall, 2 - middle channel, 3 - top wall, (b) for solid: 1 - top radiator part, 2 - bottom radiator part Температурное поле во входном и выходном сечениях радиатора приведено на рис.8. Температура стенки, к которой осуществляется подвод тепла равна 300 К во входном и 308 К в выходном сечении. Температура верхней стенки равна 297.5 К и 304.8 К соответственно. Таким образом во входном сечении разница температур на верхней и нижней стенках радиатора составляет 2.5 К, а в выходном сечении - 3.2 К. Средняя температура
прогреваемой стенки составляет 305.4 К, а верхней стенки радиатора - 302 К, т.е. температура в среднем понижается на 3.4 К. Температура воды в ядре потока меняется с 293 К до 302 К. Средний прогрев жидкости в микроканале составляет 11 К. Количество тепла, отводимое водой, составило 14 Вт/см2.
У
2?3. 294.9 Ж,? 296.6 300.
У
I I I ^^^ III I
(б)
Рис. 8. Распределение температуры во входном (а) и выходном (б) сечении радиатора Fig. 8. Temperature field in the inlet (a) and the outlet (b) section of the radiator
5. Заключение
В работе исследован процесс сопряженного теплообмена в микроканальном радиаторе. Выявлено, что реализуемые в микроканалах течения подчиняются законам континуальной физики сплошных сред и могут быть корректно смоделированы с применением предложенной математической модели с учетом вязких эффектов рабочей жидкости. Проведенная кроссплатформенная верификация показала хорошую сходимость результатов численного моделирования. Это подтверждает корректность предложенных подходов в целом и возможность применения пакета OpenFOAM, в частности, для математического моделирования процесса сопряженного теплообмена в микроканальных радиаторах. В результате проведенных расчетов, была получена общая картина температурного режима в объеме кремниевого радиатора с микроканалами, по которым прокачивается охлаждающая жидкость. Получены локальные и интегральные характеристики работы радиатора.
Список литературы / References
[1] Филимонов С.А., Дектерев А.А. и др. Моделирование сопряженного теплообмена в системе микроканалов при помощи гибридного алгоритма. Сибирский журнал индустриальной математики, том 18, вып. 3, 2015 г., стр. 86-97 / Filimonov S.A., Dekterev A.A. et al. Simulation of conjugate heat transfer in a microchannel system by a hybrid algorithm. Journal of Applied and Industrial Mathematics, vol. 9, issue 4, 2015, pp. 469-479.
[2] Li J., Peterson G.P., Cheng P. Three-dimensional analysis of heat transfer in a micro-heat sink with single phase flow. International Journal of Heat and Mass Transfer, vol. 47, issues 19-20, 2004, pp. 4215-4231.
[3] Kendall G.E., Griffith P. et al. Small diameter effects on internal flow boiling. In Proc. of the 2001 ASME International Mechanical Engineering Congress and Exposition, 2001, pp. 1-17.
[4] Weisberg A., Bau H.H., Zemel J.N. Analysis of microchannels for integrated cooling. International Journal of Heat and Mass Transfer, vol. 35, issue 10, 1992, pp. 2465-2474.
[5] Wu H.Y., Cheng P. Friction factors in smooth trapezoidal silicon microchannels with different aspect ratios. International Journal of Heat and Mass Transfer, vol. 46, issue 14, 2003, pp. 2519-2525.
[6] Li Yu. Implementation of multiple time steps for the multi-physics solver based on chtMultiRegionFoam. In Proceedings of CFD with OpenSource Software, 2016, 50 p.
[7] Koroleva M.R., Mishchenkova O.V. et al. A Theoretical research of the internal gas dynamics processes of measurements of hot air curtain with cross-flow fan. MM Science Journal, June 2020, pp. 3966-3972.
[8] Байметова Е.С., Гиззатуллина А.Ф., Пушкарев Ф.Н. Решение задачи сопряженного теплообмена в оребренной трубке с использованием OpenF oam. Химическая физика и мезоскопия, том 23, вып. 2, 2021 г., стр. 154-164 / Baimetova E.S., Gizzatullina A.F., Pushkarev F.N. Solving the conjugate heat transfer problem in the ribbed tube with OpenFoam. Chemical physics and mesoscopy, vol. 23, issue 2, 2021, pp. 154-164 (in Russian).
[9] Лойцянский Л.Г. Механика жидкости и газа. М., Наука, 1970 г., 904 стр. / Loytsyansky L.G. Mechanics of liquid and gas. Moscow, Nauka, 1970, 904 p. (in Russian).
Информация об авторах/ Information about authors
Елена Сергеевна БАЙМЕТОВА - старший преподаватель кафедры «Тепловые двигатели и установки». Сфера научных интересов: вычислительная гидрогазодинамика и интенсификация процессов теплообмена.
Elena Sergeevna BAIMETOVA is a senior lecturer at the Department of Heat Engines and Plants. Research interests: computational fluid dynamics and intensification of heat exchange processes. Михаил Евгеньевич ХВАЛЬКО - аспирант кафедры «Тепловые двигатели и установки». Сфера научных интересов: газодинамика, гидродинамика, конвективный теплообмен.
Mikhail Evgenevich HVAL'KO is a Postgraduate Student at the Department of Heat Engines and Plants. Research interests: gas dynamics, hydrodynamics, convective heat transfer. Алексей Юрьевич АРМЯНИН - аспирант кафедры «Тепловые двигатели и установки». Сфера научных интересов: гидродинамика, конвективный теплообмен.
Aleksej Yurevich ARMYANIN is a Postgraduate Student at the Department of Heat Engines and Plants. Research interests: hydrodynamics, convective heat transfer.
