CC BY
30
ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ
«наука. инновации. технологии», № 3, 2017
УДК 621.1.016.4 Коновалов Д.А. [Konovalov D.A.]
МОДЕЛИРОВАНИЕ
теплогидравлических характеристик микроканальных
ТЕПЛООБМЕННых ЭЛЕМЕНТОВ НА ОСНОВЕ МАТРИцы МОНОКРИСТАЛЛОВ КРЕМНИЯ
Modeling of heat-hydraulic characteristics of microchannel heat transfer elements based on the matrix of silicon monocrystals
Предложена конструкция теплообменного элемента на основе матрицы нитевидных кристаллов кремния для систем термостабилизации миниатюрных источников тепловыделения с удельной мощностью до 100 Вт/см2, работающих в широком диапазоне температур окружающей среды. На основе разработанной математической модели конвективного переноса теплоты в микроканальном компактном теплообменнике с развитой поверхностью теплообмена проведено численное моделирование процессов гидродинамики и теплообмена для различных конфигураций микроканальных вставок. Получены поля давлений, скоростей течения, температур охладителя и матрицы из монокристаллов кремния в широком диапазоне расходов охладителя, определены критериальные зависимости для числа Нус-сельта, а также потерь давления различных геометрических конфигураций теплообменников. Исследованы критические режимы работы, предложены направления оптимизации. По разработанной технологии изготовлены опытные образцы для проведения испытаний.
The construction of a heat exchange element based on a matrix of silicon whiskers for thermal stabilization systems of miniature heat sources with specific power up to 100 W / cm2 operating over a wide range of ambient temperatures is proposed. Based on the developed mathematical model of convective heat transfer in a microchannel compact heat exchanger with a developed heat exchange surface, numerical simulation of the hydrodynamics and heat transfer processes for various configurations of microchannel insertions was carried out. Fields of pressures, flow velocities, coolant temperatures and matrix from silicon single crystals have been obtained in a wide range of coolant flow rates, criteria dependencies for the Nusselt number and pressure losses of various geometric configurations of heat exchangers have been determined. Critical operation modes are investigated, optimization directions are proposed. According to the developed technology, prototypes for testing have been manufactured.
Ключевые слова: монокристалл кремния, конвективный теплообмен, критериальная зависимость, интенсификация теплообмена, критические режимы работы.
Key words: monocrystalline silicon, convective heat transfer, the criteria dependence, intensification of heat exchange, critical modes of operation.
Введение
Активное развитие таких отраслей как энергетика, электроника ведет к росту энергонапряженности отдельных элементов энергетического и электронного оборудования и росту тепловыделения на фоне повышенных требований к компактности и миниатюризации систем теп-
ловой защиты [1—3]. Для обеспечения устойчивой работы указанного оборудования необходимо поддержание рабочих тепловых режимов, что ведет к необходимости создания систем термостабилизации, работающих в условиях как высоких, так и низких температур окружающей среды. Так, например, для наземных систем управления спутниковыми комплексами требуется поддержание высокой тепловой стабильности (± 0,1°С), что может потребовать применения различного рода высокоэффективных теплообмен-ных элементов. Миниатюризация устройств и повышение нагрузок требует интенсификации теплообмена. Большое распространение получили различного рода пористые и микроканальные рекуперативные теплообменники. Широкий диапазон структурных, теплофизических, гидравлических, химических, оптических и других свойств пористых и микроканальных материалов, простота изготовления из них элементов конструкций, высокая интенсивность теплообмена - все это дает возможность использовать указанные теплообменные элементы в условиях действия высоких тепловых нагрузок и температур, там, где другие виды охлаждения конструкций оказываются малоэффективными. Так, например, в системах тепловой защиты жидкостных ракетных двигателей, сверхзвуковых огнеструйных резаков, плазмотронов, а также блоков питания, микросхем, современных процессоров, станций базовой, спутниковой и космической связи были разработаны и успешно апробированы пористые теплообменники с межканальной транспи-рацией охладителя.
Одним из перспективных способов подвода/отвода тепла с энергонапряженной поверхности является применение микроканальных теплообменников на основе матрицы нитевидных кристаллов кремния. Это позволяет существенно увеличить коэффициент теплоотдачи от тепловыделяющей поверхности к теплоносителю, а удельный теплосъем может составлять до 100 Вт/ см2. Следует отметить, что на процесс теплопереноса оказывает непосредственное влияние не только конвективная составляющая, но и теплофизичес-кие свойства самого теплоотводящего элемента, а также термическое сопротивление между «горячей» поверхностью и охладителем. Так, при величине коэффициента внутрипорового теплообмена равным = 109 ... 1011 Вт/(м3К), эффективный коэффициент теплопроводности медной матрицы составит 160 Вт/(мК), а удельное термическое сопротивление между нагретой поверхностью и теплообменником составляет порядка от 10-3 до 10-2 (мК)/Вт, что в целом сводит на «нет» эффективность пористого охлаждения.
В свою очередь, технология получения монокристаллов кремния позволяет в некоторых случаях выращивать теплообменный элемент с регулярной структурой на тепловыделяющей поверхности, что исключает термическое сопротивление в месте контакта. Следовательно, поддержание заданной рабочей температуры энергонапряженных устройств, точное прогнозирование режимов их работы является актуальной задачей.
Материалы и методы исследований
На основе технологии [4], был разработан и создан макет базового варианта теплоотводящего элемента из нитевидных монокристаллов кремния, выращенных на кремниевой подложке (рис. 1а). Нитевидные кристаллы кремния выращиваются на кремниевых монокристаллических подложках в печи с горизонтальным расположением трубчатого кварцевого реактора в открытой хлоридно-водородной системесреде. После разращива-ния кристаллов, подача тетрахлорида кремния в реакционную зону прекращается, а реактор с выращенными образцами нитевидных кристаллов охлаждается до комнатной температуры. Морфологические исследования выполняются методами сканирующей зондовой микроскопии. Предложенная технология позволяет выращивать нитевидные кристаллы кремния с различ-
1 V-/ / \
ной конфигурацией их расположения (сплошной ряд, змейка и пр.)
Теплообменник состоит из медного основания, в котором располагается теплообменный элемент из матрицы монокристаллов кремния и верхней крышки со штуцерами подвода и отвода охладителя. Основание и крышка стягиваются винтами. Для предотвращения утечки охладителя предусмотрена резиновая прокладка. Охладитель подается в теплообменник (рис. 1б) через штуцер подвода охладителя, проходит через область матрицы нитевидных монокристаллов кремния, нагревается от них, охлаждая теплонапряжен-ный элемент, на котором установлен теплообменник, и отводится через штуцер отвода охладителя.
Кремниевая подложка имеет размеры 20 х 20 мм (рис. 2а), толщина подложки на которой выращены шипы 0,2 мм, высота шипов 1 мм, диаметр шипов 0,1 мм. Теплоотводящие элементы (шипы) расположены в шахматном
а - матрица нитевидных б - общий вид теплообменника.
монокристаллов кремния;
Рис. 1. Теплообменник с матрицей из нитевидных монокристал-
лов кремния.
ш
200 мкм
400 мкм
ГТгТГ+т4_
600 мкм
»
* - гТ-тТт4|4--ьььь-
а - 3D модель; б - схема расположения шипов на подложке.
Рис. 2. Кремниевая подложка.
порядке, расстояние между их центрами по оси х и у равно 200, 400, 600 мкм соответственно (рис. 2б). В таком варианте теплоотводящие элементы образуют монолитную конструкцию вместе с подложкой, что исключает термическое сопротивление.
В настоящее время практически отсутствуют результаты математического и численного моделирования процессов тепломассопереноса в микроканальных структурах подобного рода. Применение математических моделей и классических критериальных уравнений тепломассопереноса для макромоделей в том числе обтекание пучков труб в теплообменниках для оценки процессов в микромасштабах дает неудовлетворительные результаты. Использование модели идеальной пористой среды для прогнозирования работы мик-роканльных элементов также является несовершенной, что заключается в серьезном расхождении результатов численного и экспериментального моделирования.
Это связано с тем, что необходимо исследование особенностей процессов гидродинамики и теплообмена в микроканалах, когда следует учитывать пристеночные течения вблизи стенок и шипов, а также развитые течения в свободном пространстве. Существенный градиент температур между шипом и охладителем также вносит свои коррективы в точность моделирования процессов гидродинамики, в том числе и на нестационарных режимах работы. Некорректность существующих моделей ведет к потере точности в оценке процессов тепломассопереноса в разрабатываемых элементах тепловой защиты, что может привести к ее разрушению при критических режимах работы.
Моделирование конвективного теплообмена в подобного рода средах основывается на феноменологических уравнениях Дарси-Бринкмана-Фор-чхеймера [5-6]. Для получения аналитического решения выделяется элементарный объем для одиночного шипа, производится объемное осреднение, а затем совместно применяется одностороннее интегрального преобразования Лапласа и конечное интегральное преобразования Фурье. Для подтверждения корректности и адекватности предложенной математической модели, уточнения гидродинамической обстановки было проведено численное моделирование в среде ANSYS.
ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ
Моделирование теплогидравлических характеристик.
Важным этапом при построении математической модели для исследования характеристик течения охладителя является выбор системы уравнений для расчета. От этого зависит качество и точность полученных результатов. В ANSYS Fluent эти системы уравнений представлены в моделях турбулентности.
Система уравнений для неизотермического (с теплообменом) течения несжимаемой жидкости в декартовых прямоугольных координатах будет состоять из уравнений неразрывности (1), движения (2) и энергии (3).
4U = 0, (1)
— + (U-V)U = J--VP + vAU, (2) 5т р
- + (U-V)T = aM + -^, (3) 5т р с.
где и - скорость течения охладителя;
р - давление;
т - время;
р - плотность жидкости;
О - результирующий вектор массовых сил,
/л - коэффициент динамической (молекулярной) вязкости;
V - кинематическая вязкость среды (у = /л / р),
а - коэффициент температуропроводности (а = X /рср).
С целью ускорения вычислительного процесса был использован модельный элемент подложки нитевидных монокристаллов кремния, показанной на рисунке 3 а. Геометрические размеры модели составили 2 х 20 мм. Остальные геометрические размеры и схема расположения шипов оставались неизменными. Расчетная область вместе со сгенерированным дискретным разбиением показана на рисунке 3б.
Выбор модели турбулентности при использовании пакета ANSYS следует учитывать, что модель должна быть протестирована для аналогичного класса задач, должна быть инвариантна к использованию геометрических и теплофизических параметров, характеризующих выбранный тип теплооб-менных элементов с регулярной матрицей и использовать имеющиеся в наличии экспериментальные данные по гидротермическим полям, а также отвечать критерию устойчивости и сходимости вычислительного процесса.
В качестве математической модели гидродинамики течения охладителя была выбрана двухслойная модель SST, которая обеспечивает хорошую сходимость с достаточной скоростью с учетом структуры турбулентного потока
а - участок подложки с шипами б - 3D модель расчетной области
Рис. 3. Модельный элемент.
в ядре и пограничном слое вблизи поверхности шипов. Преимущество данной модели состоит в том, что данная модель учитывает одновременно вих-реобразование (модель k — ю) и интенсивность пульсаций (модель k — е). Такой выбор коррелируется с исследованиями модели SST, приведенными в работах [7-8].
При анализе представленного программного образа для пакета ANSYS были использованы следующие допущения:
- теплоноситель являлся ньютоновской несжимаемой жидкостью;
- течение теплоносителя носит стационарный и трехмерный характер;
- теплофизические свойства потока и матрицы нитевидных монокристаллов кремния рассчитывались при средней температуре процесса;
- начальный гидродинамический участок на входе теплооб-менного элемента отсутствовал;
- внутренние тепловыделения отсутствуют;
- тепловой поток подводится через нижнюю стенку, боковые и верхняя стенка теплоизолированы.
Для проведения математического моделирования использовался специализированный расчетный комплекс Ansys Fluent v. 17, предназначенный для численного решения уравнений движения жидкости и теплообмена в интересуемой расчетной области, что позволило снизить трудоемкость и сократить длительность расчетов. Алгоритм использования Ansys Fluent для решения поставленной задачи заключается в следующем. На первом этапе осуществляется импорт построенной Зд-модели исследуемого объекта в решатель Fluent; на втором - декомпозиция расчетной области и построение сетки, оценивается ее качество; на третьем этапе необходимо определить граничные условия, выбрать параметры расчета, свойства материалов,
выбор дополнительных моделей для моделирования турбулентности; на четвертом производится решение поставленной задачи. Не учитывалось влияние подводящих и отводящих коллекторов.
Для проведения расчетов была сгенерирована объемная модель объекта, после чего декомпозировалась расчетная область в виде плоских треугольных элементов для подложки и в виде тетраэдров для области течения с использованием граничного разделения. Сгенерированная сетка для расчетных областей имела следующие параметры:
- для области с шагом расположения шипов 200 х 200 мкм: тип сетки - тетрагональная; размеры ячеек:
min 3,793176 х 10-17 м, max 2,098848 х 10-11 м; количество ячеек - 17736744 шт.;
- для области с шагом расположения шипов 400 х 400 мкм тип сетки - тетрагональная; размеры ячеек:
min 4,980473 х 10-17 м, max 3,267845 х 10-11 м; количество ячеек - 16999847 шт.;
- для области с шагом расположения шипов 600 х 600 мкм: тип сетки - тетрагональная; размеры ячеек:
min 3,569304 х 10-17 м, max 3,999074 х 10-11 м; количество ячеек - 16534874 шт.
В расчетах варьировались геометрические размеры матрицы, а именно межцентровое расстояние шипов 200 х 200, 400 х 400 и 600 х 600 мкм. Диапазон расходов, температура на входе и подводимый тепловой поток были одинаковы для всех конструкций подложек.
Результаты исследований и их обсуждение
Для проведения вычислительного эксперимента были выбраны следующие исходные данные. В качестве теплоносителя была выбрана вода. Массовый расход охладителя изменялся в диапазоне от 0,0006 кг/с до 0,0100 кг/с. На поверхности нагрева задаются граничные условия второго рода, удельный тепловой поток с тепловыделяющей поверхности q = 100 Вт/ см2. Отсутствует тепловой поток через боковые и верхнюю стенки. Коэффициент теплопроводности кремния Xs = 130 Вт/мК; коэффициент теплопроводности воды Xf = 130 Вт/м •К. На входе в теплообменник заданы граничные условия 1 рода, а температура охладителя составила 20 °С.
В качестве примера расчета на рисунках 4-5 представлена визуализация расчетов ANSYS для матрицы нитевидных кристаллов кремния с шагом 600 х 600 мкм.
Скорость, м/с
а - поле скоростей
~~ ^С^Г^Г^ЕГ : Г г г -г г - т- - £ - г
Ч=~- "_ 4:-.' С- Г -¿-Г > С-Г £
__ _й ---о-- -й- -5 й----ег-- -¿Г" ^ ¡Г Е _ Ь _ о
— ~~ _ ^ ¿Г 1 £ -Г _ 'С _ С
п?
1
Скорость, м/с
б - линии тока
Рис. 4. Гидродинамическая картина течения для матрицы 600х600 мкм.
«$>■ V <$>■ & & & >§> У У <$>■
Температура. оС
а - температура теплоносителя
ц
ш
Температура, оС
I I !
б - то же вид сверху
I М I ' I ! I ! I
Температура. оС
в - температура матрицы нитевидных монокристаллов кремния
Рис. 5.
Распределение температуры микроканальном элементе (600x600 мкм):
79
Т, °С
G, кг/с
69
59
49
39
29
Каркас
0,000
0,002
0,004
0,006
0,01
Рис. 6.
Зависимость температуры потока и каркаса от расхода охладителя (600 х 600 мкм).
19
600 х 600
10 А^ °С \ G, кг/с
8
6 \\ 200 к 200
4 ччч
2 400 х
0 0,002 0,004 0,006 0,008 0,01
Рис. 8. Сравнительный график зависимости перепада температу-
ры на входе и выходе из матрицы от расхода охладителя.
180_Ар, кПа _G, кг/с
Рис. 9. Сравнительный график зависимости перепада давления
от расхода охладителя.
На рисунке 6 приведены зависимости температур матрицы и охладителя в исследуемом диапазоне расходов.
На рисунках 7-9 представлены зависимости теплогидравлических характеристик в зависимости от расхода теплоносителя.
Обработка результатов вычислительных экспериментов представлена в классическом критериальном виде для чисел Нуссельта и потерь давления для примененного спектра структуры используемых матриц. Результаты обработки результатов для:
- подложки с межцентровым расстоянием шипов 200 х 200 мкм N4 = 12,035 • Яе0,3876, Ар = 154,295 • С1,4888;
- подложки с межцентровым расстоянием шипов 400 х 400 мкм N4 = 3,173 • Яе0,4515, Ар = 24,841 • С1,5422;
- подложки с межцентровым расстоянием шипов 600 х 600 мкм N4 = 1,9698 • Яе0,4833, Ар = 28,645 • С1 6057.
где N4 = а^й/к; Яе = МЫ;
а; - коэффициент теплоотдачи, Вт/(м2-К);
й - характерный размер, м.
Анализа результатов вычислительного эксперимента показал, что течение теплоносителя симметричный характер относительно стенок и ядра потока. При расходах 0,0006-0,0030 кг/с течение теплоносителя носит ламинарный характер. Отсутствуют завихрения потока, как в целом, так и вблизи локальных элементов (шипов). При расходах теплоносителя свыше 0,0030 кг/с наблюдается ярко выраженный переходный режим течения от ламинарного к турбулентному с локальным возникновением отрывных течений за шипами, что свидетельствует о возникновении застойных зон, величина которых достигает расстояния между шипами только при высоких расходных характеристиках. Наблюдается существенное различие между температурой каркаса матрицы нитевидных монокристаллов кремния и температуры охладителя. Это свидетельствует о корректности применения двухтемпературно-го подхода при построении математической модели. Из полученных данных следует, что с увеличением расхода охладителя не наблюдается возрастания теплосъема, т.к. прогрев теплоносителя происходит только в нижней его части. Уменьшение шага шипов также не способствует увеличению эффективности охлаждения, и к тому же значительно увеличивает потери давления охладителя через матрицу нитевидных монокристаллов кремния.
Найдена рациональная высота шипа ~ 0,6 мм и межцентровое расстояние между шипами 600 х 600 мкм, которые позволяют максимально охлаждать подложку при несущественном увеличении потерь давления на прокачку теплоносителя. На начальном тепловом участке X ~ 0,21 (I - длина теплообменника) высота шипа может быть снижена до величины от 0,2 до 0,3 мм, что положительно скажется на гидравлическом сопротивлении, без снижения
тепловых показателей. Разработанная математическая модель позволяет оценить локализацию по полю температур теплоносителя критическую изотерму и сделать вывод о возможном фазовом переходе с учетом получаемого на данной изотерме абсолютного давления.
Выводы
Вычислительный эксперимент на базе пакета ANSYS позволил уточнить локальную гидродинамическую структуру течения теплоносителя и показать, что влияние частоты шахматного расположения шипов в диапазонах 0,2 х 0,2; 0,4 х 0,4 и 0,6 х 0,6 мм не оказывает существенного влияния на перемешивание теплоносителя при его течении через теплообмен-ный элемент. Подтверждено использование модели идеального вытеснения при описании теплообмена в теплообменном элементе с регулярной матрицей из нитевидных монокристаллов кремния. Результаты аналитического и численного решения показывают существенную разность температур между каркасом и охладителем, а также существенную нелинейность распределения температуры по высоте шипов. При моделировании микроканальных тепло-обменных элементов с матрицей нитевидных монокристаллов кремния необходимо учитывать данные факты для исключения локального перегрева. Результаты моделирования могут рассматриваться как основа для инженерной методики расчета и инструментарий для выбора конструктивных параметров проектируемых теплообменных систем и идентификации рациональных режимов их функционирования.
Библиографический список
1. Интенсификация тепло- и массообмена на макро-, микро- и на-номасштабах: монография / Б. В. Дзюбенко, Кузма-Кичта Ю. А., Леонтьев А. И. и др. // М.: ФГУП «ЦНИИАТОМИНФОРМ», 2008. 532 с.
2. Разработка и моделирование микроканальных систем охлаждения: монография / Д. А. Коновалов, Дроздов И. Г., Шматов Д.П. и др. // Воронеж: Воронеж. гос. техн. ун-т, 2013. 222 с.
3. Коновалов Д.А., Лазаренко И.Н., Дроздов И.Г., Шматов Д.П. Современные подходы к разработке и созданию элементов систем тепловой защиты радиоэлектронных компонентов. Вестник ВГТУ Т. 10. № 1. 2014. С. 97-104.
4. Небольсин В.А., Иевлева Е.В., Шмакова С.С. О взаимосвязи электронного строения и каталитических свойств металлов -катализаторов роста нитевидных кристаллов кремния // Вестник ВГТУ. 2012. Т. 8. № 7. С. 37-42.
5. Bear J, Bachmat Y. Introduction to modeling of transport phenomena in porous media.Netherlands: Kluwer Academic Publishers, 1991 -553 p.
6. Chen G.M., Tso C.P. A two-equation model for thermally developing forced convection in porous medium with viscous dissipation. International Journal of Heat and Mass Transfer. 2011. v. 54, №2526. pp 5406-5414.
7. Белов И. А., Исаев С. А. Моделирование турбулентных течений: учеб. пособие Балтийский государственный технический университет «Военмех». Типография БГТУ. Санкт-Петербург, 2001. 106 с.
8. Гарбарук, А.В., Стрелец М.Х., Шур М.Л. Моделирование турбулентности в расчетах сложных течений: учеб. пособие. СПб.: Изд-во Поли-техн. ун-та, 2012. 88 с.
References
1. Dzyubenko B.V., Kuzma-Kichta Yu.A., Leont'ev A.I., Fedik I.I., Khol-panov L.P. Intensifikatsiya teplo- i massoobmena na makro-, mikro-i nanomasshtabakh. [Intensification of heat and mass transfer at the macro-, micro - and nanoscale]. Moscow, FGUP «TsNIIATOMIN-FORM», 2008. 532 p.
2. Konovalov D.A., Drozdov I.G., Shmatov D.P., Dakhin S.V., Kozhuk-hov N.N. Razrabotka i modelirovanie mikrokanal'nykh sistem okhla-zhdeniya [Development and modeling of microchannel cooling sys-tems].Voronezh, VSTU publ., 2013. 222 p.
3. Konovalov D.A., Lazarenko I.N., Drozdov I.G., Shmatov D.P. Sovre-mennye podkhody k razrabotke i sozdaniyu elementov sistem teplovoy zashchity radioelektronnykh komponentov [Modern approaches to the development and creation of elements of system of thermal protection of electronic components]. Voronezh, VSTU publ., 2014, vol. 10, no. 1, pp. 97-104. (In Russ).
4. Nebol'sin V.A., Ievleva E.V., Shmakova S.S. O vzaimosvyazi elek-tronnogo stroeniya i kataliticheskikh svoystv metallov - katalizato-rov rosta nitevidnykh kristallov kremniya [On the relationship of the electronic structure and catalytic properties of metal catalysts for the growth of silicon whiskers]. Voronezh, VSTU publ., 2012, vol. 8, no. 7.2, pp. 37-42. (In Russ).
5. Bear J., Bachmat Y. Introduction to Modeling of Transport Phenomena in Porous Media. Dordrecht, Kluwer Academic Publishers, 1991. 553 p.
6. Chen G.M., Tso C.P. A two-equation model for thermally developing forced convection in porous medium with viscous dissipation. International Journal of Heat and Mass Transfer, 2011. vol. 54, no. 25-26, pp. 5406-5414.
7. Belov I.A., Isaev S.A. Modelirovanie turbulentnykh techeniy. [Modeling of turbulent flows] Sankt-Peterburg, BGTU pub., 2001. 106 p.
8. Garbaruk, A.V., Strelets M.Kh., Shur M.L. Modelirovanie turbulent-nosti v raschetakh slozhnykh techeniy [Modeling turbulence in complex flows calculations: proc. manual]. Sankt-Peterburg. Politechnic pub. 2012. 88 p.
