CC BY
15Работа выполнена при финансовой поддержке РФФИ (проект № 08-01-00863), программы "Ведущие научные школы РФ" (проект НШ-4470.2008.1) и программы фундаментальных исследований Отделения математических наук РАН "Алгебраические и комбинаторные методы математической кибернетики" (проект "Синтез и сложность управляющих систем").
СПИСОК ЛИТЕРАТУРЫ
1. Кузнецов А. В. О средствах для обнаружения невыводимости или невыразимости // Логический вывод. М.: Наука, 1979. 5-33.
2. Касим-Заде О.М. Об одной метрической характеристике неявных и параметрических представлений булевых функций // Матем. вопросы кибернетики. Вып. 6. М.: Наука, Физматлит, 1996. 133-188.
3. Михайлец Е.В. О ранге неявных представлений над классами монотонных функций ^-значной логики // Мат-лы VI Молодежной научной школы по дискр. матем. и ее прил. Ч. II. М., 2007. 26-29.
4. Орехова Е.А. Об одном критерии неявной полноты в трехзначной логике // Матем. вопросы кибернетики. Вып. 12. М.: Наука, Физматлит, 2003. 27-74.
Поступила в редакцию 28.03.2008
УДК 533.6.011.5:532.526:541
ИССЛЕДОВАНИЕ ТЕЧЕНИЯ И ТЕПЛООБМЕНА В МИКРО- И НАНОКАНАЛАХ МЕТОДАМИ МОЛЕКУЛЯРНОЙ ДИНАМИКИ
В. Л. Ковалев, А. Н. Якунчиков
1. Микроэлектронные компоненты уменьшаются в размерах, а количество энергии, рассеиваемое их системами охлаждения, неуклонно увеличивается [1]. В связи с этим вопрос об охлаждении электронных компонентов стоит достаточно остро. Предполагается, что в будущем системы охлаждения будут представлять собой систему микро- и наноканалов, пронизывающую электронный компонент. По каналам будет осуществляться циркуляция охлаждающей жидкости или газа. Прототипы таких систем охлаждения уже появились в исследовательских институтах США. Охлаждающее устройство содержит два базовых элемента — эмиттер и коллектор. Эмиттер — это система заряженных игл с остриями малого диаметра, заряжающая частицы воздуха. Коллектор за счет действия электрического поля на ионизированную среду создает поток, охлаждающий чип. Также поток может создаваться "микронасосом" — осциллирующей стенкой канала.
Течение газов и жидкостей в большинстве случаев исследуется на основе макроскопического подхода, когда среда рассматривается как континуум [2, 3]. Такое описание справедливо в случаях, когда изучаемый объем содержит достаточно большое число молекул и среду можно считать непрерывной. Однако при исследовании течения в микро- и наноканалах, при моделировании физико-химических процессов в газе и на поверхности в ряде случаев необходимо использовать микроскопический подход, основанный на молекулярной динамике и прямом численном моделировании [4-8]. При таком подходе учитывается корпускулярная структура газа, определяются положения и скорости молекул в каждый момент времени, а макроскопические величины отождествляются со средними значениями соответствующих молекулярных величин.
2. Исследовалось течение теплопроводного совершенного газа между двумя пластинами, расположенными на расстоянии Ьу. Течение считалось двумерным, а область течения — симметричной относительно плоскости, равноудаленной от обеих пластин. Задача решалась методом прямого численного моделирования. При этом реальное течение газа описывалось при помощи большого количества моделирующих частиц, изменение координат, скоростей и свойств которых со временем вызвано межмолекулярным взаимодействием и взаимодействием с границами физического пространства. Уравнения движения частиц записывались в виде
(Рщ (М
-1Г =
где Хг, = —--координата и скорость г-й частицы, а правые части этих уравнений ^ задавались в
аЬ
соответствии с выбранной моделью взаимодействия частиц.
Физическое пространство разбивалось на ячейки — малые объемы У0 = (¿о)3 с фиксированными границами. Здесь Ьо — линейный размер ячейки. Изменение времени происходило дискретно с шагом АЬ3. На каждом временном шаге рассчитывались перемещения частиц, проверялось выполнение граничных условий (взаимодействие с поверхностью пластин) и производилось перераспределение частиц по ячейкам. Считалось, что первоначально частицы распределены в физическом пространстве равномерно, а их скорости — в соответствии с равновесной максвелловской функцией распределения.
При описании взаимодействия частиц совершенного газа молекулы представлялись в виде упругих сфер диаметра й\ считалось, что взаимодействие происходит только за счет столкновений (модель твердых сфер). При этом скорости молекул после столкновения вычислялись по скоростям до столкновения с помощью закона сохранения импульса и закона сохранения энергии:
Г ^ + Ш,^ = шХ + Ш, ,
I 2 2 /2 /2
[ Ш1У2 + Ш, V2 = ШЩ + Ш, V, ,
где Шг, VI — масса г-й частицы и проекция ее скорости на линию центров.
В целях экономии вычислительных ресурсов при описании столкновений использовался вероятностный подход [5]. Моделирование столкновений проводилось отдельно в каждой ячейке. Сталкивающиеся частицы выбирались случайным образом. Учитывались только парные столкновения. Относительная ве-
роятность столкновения такой пары pr =
vr
Vm
где vr — скорость одной частицы относительно другой,
а vmax — скорость, заведомо превосходящая любую относительную скорость в этой ячейке. Вводился счетчик времени Ato, значение которого увеличивалось при столкновении на величину ("время одного столкновения")
1 2
h ~ ТЩ ~ Vo(arí2cr)'
Расчеты велись до момента совпадения значения счетчика Ato и шага по времени Ats. При этом число столкновений N, происходящих в единице объема за единицу времени при заданных параметрах потока,
средняя относительная скорость двух частиц, а = nd2 — сечение
1
согласно [5], равно N = — п2асг, где сг -столкновения, n — числовая плотность.
Действительно, пусть и\ — скорость пробной частицы. Рассмотрим частицы с плотностью Ап и скоростью и (частицы класса и). Пробная частица движется относительно этих частиц со скоростью сг = и\ — и. Она столкнется с частицей класса и за интервал времени АЬ, если центр частицы класса и окажется внутри цилиндра с основанием а и высотой сгАЬ (рис. 1). Так как объем цилиндра Ус = асгАЬ = пОРвгАЬ, то количество частиц класса и в объеме Ус равно УсАп = па2сгАЬАп. Поэтому частота столкновений пробной частицы с молекулами класса и равна ис = асг Ап. Средняя частота столкновений вычисляется суммированием по всем скоростям:
An
Рис. 1. Схема сечения столкновения частиц
v = V Апсга = па У" сг-
n
= nacr
Общее число столкновений в единице объема за единицу времени составит N = п^/2 = п2асг/2. Коэффициент 1/2 вводится, чтобы не считать дважды каждое столкновение, так как столкновение включает две молекулы.
Отметим, что при используемом подходе возмущения могут распространяться со скоростью а =
д^-. Поэтому в целях исключения влияния вычислительных эффектов на результаты расчетов Ьо и
выбирались так, чтобы а < а, где а — скорость звука в данной области.
Для описания взаимодействия газа с поверхностью использовалась диффузная модель. Считалось, что скорости каждой из молекул после отражения не зависят от их индивидуальных скоростей падения
и распределены согласно равновесной максвелловской функции распределения в том полупространстве скоростей, где вектор скорости молекул направлен от поверхности. Распределение соответствовало температуре поверхности.
Зеркальное отражение
дш^ийм
Плоскость симметрии
Рис. 2. Схема моделируемой области
3. Схема моделируемой области представлена на рис. 2. Течение считалось симметричным, поэтому
моделирование проводилось по одну сторону от плоскости симметрии. Для этого на ней ставилось условие
зеркального отражения частиц. В области организовывалось течение: газ охлаждался до температуры
Тд = 0,9То и приобретал среднюю скорость г°. Пластина в этой области имела температуру Тт = Тд =
0,9То. Область Б — изучаемая, здесь пластина имела температуру Тт = То.
Исследовалось течение в замкнутом канале: после прохождения изучаемой области газ снова попадал
в область А, охлаждался и приобретал требуемую начальную скорость. Течение считалось стационарным.
Во избежание флуктуаций в результатах проводилось осреднение по времени.
В качестве газа использовался азот N2 с молярной массой Шщ2 = 28 г. Рассматривались условия,
при которых длина свободного пробега Л имеет порядок размера молекулы й = 3,78 • 10_1° м. Средний
4 кг
по сечению расход в канале поддерживался постоянным и равным ро^'о = 1,3 • 10 -^5—. Температура
м2с
торможения Т° = 77 К.
На рис. 3 показано распределение безразмерной скорости - в поперечном сечении канала (Рср —
ср
Л
средняя по сечению скорость) в зависимости от числа Кнудсена Кп = —. При малых числах Кп, когда
Ьу
справедливо описание течения в рамках механики сплошной среды, получен параболический профиль скорости и расчеты совпали с решением уравнений Навье-Стокса для выбранных условий обтекания.
Рис. 3. Профили скорости при различных числах Кнудсена: 1 — аналитическое решение; 2 — Кп = 0,0003; 3 — 0,012; 4 — 0,03; 5 — 0,06
Рис. 4. Профили температуры в нескольких сечениях канала: 1 — х = 0, Кп = 0,006; 2 — х = 0, Кп = 0,03; 3 — х = Ьу, Кп = 0,006; 4 — 0,03; 5 — х = 2Ьу, Кп = 0,006; 6 — х = 2Ьу, Кп = 0,03
х = Ьу, Кп
70
вестн. моск. ун-та. сер. 1, математика. механика. 2008. №5
Для достаточно больших чисел Кнудсена Кп, когда длины свободного пробега сравнимы с линейными размерами задачи, профиль скорости практически выравнивается, а влияние стенки сказывается только в непосредственной близости к ней.
На рис. 4 показано распределение безразмерной тем-Т
пературы
Tw
в нескольких поперечных сечениях канала
Рис. 5. Распределение плотности (с скорости (б) и температуры (в)
при различных числах Кнудсена. Расчеты проводились при Тш = Т°, Тд = 0,9Т°, Кп = 0,006 и Кп = 0,03. Видно, что в канале меньшей ширины газ прогревается быстрее. Вблизи поверхности на начальном участке нагрева безразмерный градиент температуры для Кп = 0,006 и Кп = 0,03 равен соответственно 0,222 и 0,133. Но в размерных координатах градиент температуры (а следовательно, и тепловой поток от поверхности) при Кп = 0,03 больше почти в 3 раза, чем при Кп = 0,006.
На рис. 5 для Кп = 0,006 приведены характерные распределения плотности, скорости и температуры в изучаемой области. О численных значениях параметров следует судить по насыщенности цвета и шкале соответствующих ей значений.
Заключение. В работе развит метод прямого численного моделирования течения и теплообмена в микро-и наноканалах, основанный на подходах молекулярной динамики. Получены распределения плотности, скорости и температуры в канале при различных числах Кнудсена.
СПИСОК ЛИТЕРАТУРЫ
1. Schmidt R.R., Notohardjono B.D. High and server low temperature cooling // IBM J. Res. and Develop. 2002. 46. 739-751.
2. Седов Л.И. Механика сплошной среды. Т. 1, 2. М.: Наука, 1976.
3. Ковалев В.Л. Гетерогенные каталитические процессы в аэротермодинамике. М.: Наука, Физматлит, 2002.
4. Слезкин Н.А. Лекции по молекулярной гидродинамике. М.: Изд-во МГУ, 1981.
5. Bird G.A. Molecular gas dynamics. Oxford: Clarendon, 1976.
6. Ivanov M.S. Computational hypersonic rarefied flows // Annu. Rev. Fluid Mech. 1998. 30. 469-505.
7. Ковалев В.Л, Сазонова В.Ю., Якунчиков А.Н. Динамический метод Монте-Карло моделирования поверхностной рекомбинации // Вестн. Моск. ун-та. Матем. Механ. 2007. № 2. 67-72.
8. Nedea S.V., Markvoort A.J., Frijns A.J., Steenhoven A.A., Hilbers P.A. Hybrid Molecular Dynamics — Monte Carlo method for heat and flow analysis in micro/nanochannels // Proc. 25th Int. Symp. on Rarefied Gas Dynamics. St-Petersburg, July 21-28, 2006.
Поступила в редакцию 04.06.2007
