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

УДК 532.685

МОДЕЛИРОВАНИЕ ПРОЦЕССОВ ГИДРОДИНАМИКИ ТЕЧЕНИЯ ОХЛАДИТЕЛЯ В НАНОСТРУКТУРАХ НА ОСНОВЕ НИТЕВИДНЫХ КРИСТАЛЛОВ КРЕМНИЯ

Д.А. Коновалов, И.Г. Дроздов, И.Н. Лазаренко, Д.П. Шматов

Данная работа посвящена моделированию гидродинамики течения охладителя в наноструктурах на основе нитевидных кристаллов кремния. Предложено два варианта матрицы из нитевидных кристаллов кремния, которые образуют развитую поверхность теплообмена в микроканальных теплообменниках. Произведен расчет с применением четырех моделей турбулентности и выбран ее наиболее предпочтительный вариант. Моделирование течения охладителя проводилось с помощью программного комплекса АЫБУБ СЕХ

Ключевые слова: гидродинамика, микроканальные теплообменники, наноструктуры, моделирование

Современное развитие энергетики и

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

Одним из способов интенсификации

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

Исследованию этого вопроса посвящено

множество работ отечественных и зарубежных авторов, например [1].

Одним из наиболее эффективных способов защиты считается пористое охлаждение [2], которое предусматривает изготовление энергетических конструкций или деталей из пористых материалов, через поры которых во время работы подается охладитель (газ или жидкость) на поверхность, соприкасающуюся с нагретой рабочей средой.

Широкий диапазон структурных,

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

малоэффективными, в условиях действия высоких тепловых нагрузок и температур.

Коновалов Дмитрий Альбертович - ВГТУ, канд. техн. наук, ст. преп., тел. (473) 277-27-55, e-mail:

[email protected]

Дроздов Игорь Геннадьевич - ВГТУ, д-р техн. наук, доцент, тел. (473) 277-27-55, e-mail: [email protected] Лазаренко Игорь Николаевич - ВГТУ, аспирант, тел. (473) 277-27-55, e-mail: [email protected] Шматов Дмитрий Павлович - ВГТУ, канд. техн. наук, ст. преп., тел. (473) 277-27-55, e-mail: [email protected]

При чрезвычайно высоких тепловых потоках, больших температурных напорах, сложной конструкции или малой доступности поверхности пористое и микроканальное охлаждение - наиболее надежный метод тепловой защиты. Так, например, в системах тепловой защиты жидкостных ракетных двигателей, сверхзвуковых огнеструйных резаков, плазмотронов, микропроцессоров и микросхем тепловой поток составляет 100 Вт/см2 и выше, а применение иных систем охлаждения ограничено размерами тепловыделяющего элемента [2 - 9].

Следует отметить, что на процесс теплопереноса оказывает непосредственное влияние не только конвективная составляющая, но и теплофизические свойства непосредственно теплоотводящего элемента, а также термическое сопротивление между «горячей» поверхностью и охладителем. Так, например, при величине коэффициента внутрипорового теплообмена равным ку = 109...1011 Вт/(м3 К), эффективный коэффициент теплопроводности медной матрицы составит 160 Вт/(мК), а удельное термическое

сопротивление между нагретой поверхностью и теплообменником составляет порядка

10-3, 10-2 (м К)/Вт, что в целом сводит на «нет» эффективность пористого охлаждения.

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

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

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

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

Фотографии макета базового варианта теплоотводящих элементов из нитевидных кристаллов кремния выращенных на подложке полупроводника представлены на рис. 1.

Рис. 1. Фото матрицы нитевидных кристаллов (высота -200 мкм), увеличение 65 крат

Был разработан и создан теплообменный аппарат с теплоотводящими элементами (рис. 2) в виде матрицы нитевидных кристаллов кремния, которые выращены на подложке полупроводника и образуют развитую поверхность теплообмена с системой коллекторов и проницаемых перегородок [9]. В таком варианте элементы охлаждения образуют монолитную конструкцию вместе с тепловыделяющим элементом, при этом исключается термическое сопротивление, которое присутствует в варианте при раздельном исполнении теплообменника и тепловыделяющего элемента.

Рис. 2. Микроканальные теплоотводящие элементы на основе матрицы нитевидных монокристаллов кремния

расположением теплоотводящих элементов (шипов) из нитевидных кристаллов кремния. На их основе были созданы 2 модели области течения охладителя через теплоотводящие элементы (рис. 3, 4). Данные модели были созданы в графическом пакете SolidWorks.

Рис. 3. Модель области течения со сплошным расположением теплоотводящих элементов

Рис.4. Модель области течения с зигзагообразным расположением теплоотводящих элементов

Расход охладителя (воды) Q через обе матрицы был выбран исходя из параметров помп, применяемых в системах охлаждения радиоэлектронной техники и экспериментальном стенде, на котором будут проведены дальнейшие исследования.

Геометрические размеры области течения матрицы со сплошным расположением теплоотводящих элементов указаны на рис. 3. Теплоотводящие элементы расположены в шахматном порядке, расстояние между центрами по оси х равно 600 мкм, шаг по оси у равен 600 мкм (рис. 5).

При проектировании сложных технических систем математическое моделирование является простым средством, позволяющим на этапе проектирования рассмотреть поведение

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

специализированных расчетных комплексов, предназначенных для численного интегрирования уравнений движения жидкости в заданной пользователем расчетной области, таких как ANSYS СБХ.

В материалах статьи рассмотрены два варианта матриц со сплошным и зигзагообразным

§

о

о

со

600 мкм

ТтлТт

ТI ТI т IТI

гггтгт-

Рис. 5. Подложка со сплошным расположением теплоотводящих элементов

Геометрические размеры области течения матрицы с зигзагообразным расположением теплоотводящих элементов указаны на рис. 4. Теплоотводящие элементы расположены в шахматном порядке, расстояние между центрами по оси х равно 400 мкм, шаг по оси у равен 400 мкм

(рис. 6). В отличие от области течения матрицы со сплошным расположением шипов, в матрице с зигзагообразным расположением теплоотводящих элементов сплошные области равные четырем заполненным рядам чередуются с пустыми двумя рядами.

400мкм

Вход охладителя

Рис. 6. Подложка с зигзагообразным расположением теплоотводящих элементов

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

Наиболее широкое применение получили модели на основе метода сеток.

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

Схематично процесс решения задачи можно представить следующим образом.

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

Рис. 7. Расчетная область со сплошным расположением теплоотводящих элементов

Геометрические размеры элемента расчетной области течения с зигзагообразным расположением теплоотводящих элементов показаны на рис. 8.

Вход охладителя

Выход охладителя

Рис. 8. Расчетная область с зигзагообразным расположением теплоотводящих элементов

Для расчетных областей была сгенерирована сетка. Для области течения со сплошным расположением теплоотводящих элементов сетка имеет параметры:

- метод построения сетки: patch conforming (построение сетки снаружи модели);

- тип сетки: tetrahedrons (тетрагональная);

- размеры ячеек: min - 3,0131 10-6м, max -

2.3327-10-5м;

- параметр достаточности сетки (aspect ratio): min - 1,1578, max - 10,94.

При построении сетки этот параметр (aspect ratio) должен составлять значение в диапазоне от 9 до 11 для получения сходимости и точного результата.

Для расчетной области с зигзагообразным расположением теплоотводящих элементов сетка имеет параметры:

- метод построения сетки: patch conforming (построение сетки снаружи модели);

- тип сетки: tetrahedrons (тетрагональная);

- размеры ячеек: min - 2,4529 10-6м, max -

2.3327-10-5м;

- параметр достаточности сетки (aspect ratio): min - 1,1582, max - 10,412.

При построении математических моделей был принят ряд допущений:

- рабочий агент считается вязкой несжимаемой средой (вода);

- течение потока трехмерное стационарное;

- теплофизические свойства потока

принимаются постоянными и равными средним значениям в исследуемом интервале температур;

- на входе в расчетную область имеется полностью развитое течение с изотропной турбулентностью;

- теплообмен с окружающей средой отсутствует (на внешней стороне выполняется условие адиабатности) [10].

Система уравнений Навье-Стокса для

описания турбулентного движения вязкой

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

V-V = 0,

pD = ~VP + V'( mVV) •

(1)

(2)

В скалярно-тензорной форме уравнения неразрывности и изменения количества движения записываются так:

duj

ex.

= 0;

(3)

5uk д ¡ \ iep i e

—- +-------(uu, )=-------------+--------Tk• (4а)

k k

,v pdXk p dXj

дик дик 1др 1 д

—- + и—— =--------— +------т к. (4)

дt } дхрдхк р дх !

С учетом уравнения неразрывности (1.3) уравнение (1.4) может быть представлено в виде

1др 1 д

. -1 и и, I = — '

дt дх

В уравнениях (1.1) - (1.4) используемые индексы определяют направления декартовой системы координат х) (здесь у = 1,2,3; к = 1,2,3;

ик, и у - декартовые составляющие скорости в направлении соответствующих осей; р - давление; t - время; р - плотность жидкости; т/к -составляющие тензора вязких напряжений; т]к = т(диу /дхк + дик /дху); ¡л - коэффициент

динамической (молекулярной) вязкости; V - вектор местной скорости потока; V = в[и1; в1 - единичные векторы; V - оператор Г амильтона; -

производная по времени.

Вопросы замыкания полученной системы уравнений решаются на различном уровне сложности. Простейший путь - использование эмпирической информации о характеристиках турбулентности, наиболее сложный заключается в выводе уравнений относительно рейнольдсовых напряжений [11].

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

эксперимента:

- на входе в расчетную область задается

постоянный расход теплоносителя V (r) = V0,

степень начальной турбулентности 5 %;

- на выходе из расчетной области задаются «мягкие» граничные условия (условие продолжения решения);

- на всех боковых поверхностях условия прилипания для уравнения движения и неразрывности u(Г) = 0;

- для турбулентных величин предполагают равенство нулю турбулентной энергии и нормальной производной скорости диссипации на стенках канала: k = delдп = 0 ;

- задаются граничные условия 1-го рода: постоянная температура теплоносителя t = 25 0C .

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

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

- Модель должна быть протестирована для рассматриваемого класса задач для того, чтобы иметь уверенность в удовлетворительном качественном и количественном описании турбулентных эффектов.

- Модель должна быть относительно простой и

удобной в вычислительном плане, легко алгоритмизироваться и опираться на экспериментальные параметры, которые

определяются точно и однозначно.

- Модель должна характеризоваться высокой вычислительной эффективностью, т.е. она не должна существенно увеличивать время счета.

Для численного моделирования были применены 4 модели:

1) Модель прямого численного моделирования Laminar, (DNS - Direct Numerical Simulation). Подразумевает полное решение уравнений Навье -Стокса и уравнений неразрывности в зависимости от времени. С практической точки зрения оно может применяться, чтобы проверить предложенные аппроксимации замыкания в инженерных моделях. На самом фундаментальном уровне DNS может использоваться для понимания структуры турбулентности. Результаты такого моделирования могут также рассматриваться как дополнительный источник экспериментальных данных, которые недоступны для измерительных методов. Однако для реализации прямого численного моделирования необходима чрезвычайно мелкая сетка, достаточная для разрешения вихрей наименьших масштабов. Для разрешения турбулентного течения необходима сетка в 10 узлов по каждому направлению для

каждого турбулентного вихря. В настоящее время расчет модели из 100 элементов занимает несколько часов работы современных ЭВМ. Кроме того, при прямом численном решении уравнений Навье-Стокса потребуется значительное уменьшение временного шага до величины, сравнимой с частотой турбулентных пульсаций, что значительно увеличит время расчета.

2) Модель к — т . Относится к модели с двумя

дифференциальными уравнениями. Содержит уравнение переноса кинетической энергии турбулизации к и удельной скорости диссипации пульсаций. Данная модель используется преимущественно при моделировании

пристеночных течений и низких числах Рейнольдса.

3) Модель к — е . Также относится к модели с двумя дифференциальными уравнениями. Она

развитых турбулентных течений при больших числах Рейнольдса.

4) Двухслойная SST модель. Относится к модели с двумя дифференциальными уравнениями. Учитывает в себе модели к -w и к - е . Данная модель обеспечивает плавный переход от к -w модели в пристеночных областях к модели к -е вдали от твердых стенок [11, 12].

Далее было произведено решение поставленной задачи. В результате расчета были получены: графики распределения давлений

(рис. 9, 12), величины потери давления в образцах (табл. 1, 3), значения максимальных скоростей (табл. 2, 4) и линии тока двух моделей

турбулентности SST и k-Epsilon обеих расчетных областей (рис. 10, 11, 13. 14).

основана

AP, Па 450 -г 400 350 300

250 -|------

200

150

100

50

0

О laminar ■ k-Epsilon -А- - k-Omega _ — * -SST

N

' *•

H-

4-

H-

4-

L, мм

0 4 8 12 16 20

Рис. 9. Графики распределения давлений в области течения со сплошным расположением теплоотводящих элементов

ШШБ>

Л

_________________0.0015__________________________0.003 (т)

0.00075 0.00225

Рис. 10. Линии тока области течения со сплошным расположением теплоотводящих элементов

(модель турбулентности k-Epsilon)

Velocity Streamline 1

МШ

ч

Рис. 11. Линии тока области течения со сплошным расположением теплоотводящих элементов

(модель турбулентности SST)

Результаты расчета потери давления в области течения со сплошным расположением теплоотводящих элементов приведены в табл. 1.

В табл. 2 показаны значения максимальных скоростей в области течения со сплошным расположением теплоотводящих элементов.

Таблица 1

Таблица 2

Наименование модели Потери давления в образце, Па

Laminar 278,568

K-Epsilon 382,172

K-Omega 260,566

SST (shear stress transport) 258,414

Наименование модели Максимальная скорость в образце, м/с

Laminar 0,282

K-Epsilon 0,216

K-Omega 0,277

SST (shear stress transport) 0,277

ÁP, Па 9 -г ■к

V

_ _ laminar

■ k-Epsilon' А k-Omega _ — * SST

Л.

ч-

н-

н-

Ж L, мм

0 0,48 0,96 1,44 1,92 2,4

Рис. 12. Графики распределения давлений в области течения с зигзагообразным расположением теплоотводящих

элементов

5

4

3

Velocity Streamline 1 0.033

0.025

0.017

0.008

0.000 [m sA-1]

ШВ

Л

Рис. 13. Линии тока области течения с зигзагообразным расположением теплоотводящих элементов

(модель турбулентности k-Epsilon)

Velocity Streamline 1 ГГ 0.044

0.033

0.022

0.011

0.000 [m sA-1]

■A

Рис. 14. Линии тока области течения с зигзагообразным расположением теплоотводящих элементов

(модель турбулентности SST)

Результаты расчета потери давления в области течения с зигзагообразным расположением теплоотводящих элементов приведены в табл. 3.

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

Таблица 3

Таблица 4

Наименование модели Потери давления в образце, Па

Laminar 6,332

K-Epsilon 8,302

K-Omega 5,439

SST (shear stress transport) 5,415

Наименование модели Максимальная скорость в образце, м/с

Laminar 0,044

K-Epsilon 0,033

K-Omega 0,044

SST (shear stress transport) 0,044

По полученным данным расчетов видно, что модели турбулентности K-Omega и SST дают идентичный результат, Laminar очень близкий к ним, а k-Epsilon, описывающая развитые турбулентные течения, значительно отличающийся от всех остальных моделей, что можно судить по значениям перепада давления и значениям максимальных скоростей в образцах (табл. 1 - 4, рис. 9, 12). Необходимо заметить, что модель турбулентности SST учитывает пристеночные течения и открытую турбулентность, k-Omega только пристеночные течения, поэтому можно сделать вывод о том, что развитой турбулентности внутри расчетных областей течения нет. Модель Laminar не учитывает особенности обтекания и пристеночные течения, поэтому расчеты с применением данной модели имеют некую погрешность, но для инженерных расчетов они являются допустимыми.

Следовательно, в дальнейших исследованиях модель SST будет взята за основу и приоритетна для исследования гидродинамики течения охладителя и теплообмена в наноструктурах на основе нитевидных кристаллов кремния.

Работа выполнена в рамках госзадания по НИР № ГЗ 6/12 П.

Лицензия Ansys CFX Customers: Voronezh State Technical University (618632).

Литература

1. Резник С. В. Математическое обеспечение экспериментальных исследований теплообмена в пористых материалах тепловой защиты многоразовых космических аппаратов / С.В. Резник, П.В. Просунцов, В.П. Тимошенко // Дисперсные потоки и пористые среды. Интенсификация теплообмена: труды IV Рос. нац. конф. по теплообмену. М.: МЭИ, 2006. - Т. 6. С. 109-112.

2. Калинин Э.К. Эффективные поверхности теплообмена. / Э.К. Калинин, Г.А. Дрейцер, И.З. Копп. -М.: Энергоатомиздат, 1998. - 408 с.

3. Дзюбенко Б. В. Интенсификация тепло- и массообмена на макро-, микро- и наномасштабах / Б. В.

Дзюбенко, Ю.А. Кузма-Кичта, А.И. Леонтьев, И.И. Федик, Л.П. Холпанов - Монография. - М.: ФГУП «ЦНИИАТОМИНФОРМ». 2008. - 532 с.

4. Калинин Э.К. Интенсификация теплоотдачи в каналах / Э.К. Калинин, Г.А. Дрейцер, С. Л. Ярхо. - М.: Машиностроение, 1990. - 200 с.

5. Поляев В.М. Гидродинамика и теплообмен в пористых элементах конструкций летательных аппаратов / В.М. Поляев, В.А. Майоров, Л.Л. Васильев. - М.: Машиностроение, 1988. - 168 с.

6. Пат. 51441 Российская Федерация, МПК

Устройство охлаждения электронных компонентов / Дроздов И.Г., Мозговой Н.В., Шматов Д.П., Кожухов Н.Н.; заявитель и патентообладатель Воронежский государственный технический университет.

№2005121936/22; заявл. 11.07.2005; опубл. 10.02.2006, Бюл. № 4.

7. Пат. 51449 Российская Федерация, МПК. Устройство для охлаждения катода плазматрона / Дроздов И.Г.; Кожухов Н.Н.; Габасова Э.Р.; Коновалов Д.А., Шматов Д.П. №2005121935/22; заявл. 11.07.2005; опубл. 10.02.2006, Бюл. № 4.

8. Пат. 58788 Российская Федерация. U1. Устройство охлаждения для электронных компонентов / Дроздов И.Г., Кожухов Н.Н., Мозговой Н.В., Коновалов Д.А., Шматов Д.П. - № 2006113838/22; заявл. 24.04.2006; опубл. 27.11.2006, Бюл. № 33.

9. Пат. 2440641 Российская Федерация. Устройство отвода теплоты от кристалла полупроводниковой микросхемы / Савинков А.Ю., Дроздов И.Г., Шматов Д.П., Дахин С.В., Коновалов Д.А., Кожухов Н.Н., Небольсин В.А. - № 2010146036/28; заявл. 10.11.2010; опубл. 20.01.2012, Бюл. № 2.

10. Menter F.R., Esch T. Advanced Turbulence Modelling in CFX // CFX Update - Spring 2001. - No. 20. -P. 4-5.

11. Белов И. А. Моделирование турбулентных течений: учеб. пособие И. А. Белов, С. А. Исаев // Балтийский государственный технический университет "Военмех", Типография БГТУ, Санкт-Петербург, 2001. -106 с.

12. Морозов И. И. Введение в численные методы вычислительной гидроаэродинамики: учеб. пособие И. И. Морозов, А.С. Ляскин // Издательство СГАУ, 2011. -65 с.

Воронежский государственный технический университет

MODELING HYDRODYNAMICS COOLANT FLOW IN NANOSTRUCTURES BASED ON

THREADLIKE CRYSTALS OF SILICON

D.A. Konovalov, I.G. Drozdov, I.N. Lazarenko, D.P. Shmatov

This work is devoted to the modeling hydrodynamics coolant flow in nanostructures based on threadlike crystals of silicon. Two options of a matrix are offered on threadlike crystals of silicon, which form the developed surface of heat exchange in microchannel heat exchangers. The calculation using the four turbulence models and selected its most preferred option. Modeling system of cooling was spent by means of program complex ANSYS CFX

Key words: hydrodynamics, microchannel heat exchangers, nanostructures, modeling