Моделирование микротечений методом решеток Больцмана
Рассмотрен подход к анализу микро- и макротечений с использованием метода решеток Больцмана. Получены профили скорости на начальном участке для микроканала без проскальзывания и с условием проскальзывания на стенке. Проанализировано влияние числа Кнудсена на длину начального гидродинамического участ...
Saved in:
| Published in: | Промышленная теплотехника |
|---|---|
| Date: | 2011 |
| Main Authors: | , , , |
| Format: | Article |
| Language: | Russian |
| Published: |
Інститут технічної теплофізики НАН України
2011
|
| Subjects: | |
| Online Access: | https://nasplib.isofts.kiev.ua/handle/123456789/60312 |
| Tags: |
Add Tag
No Tags, Be the first to tag this record!
|
| Journal Title: | Digital Library of Periodicals of National Academy of Sciences of Ukraine |
| Cite this: | Моделирование микротечений методом решеток Больцмана / А.И. Тыринов, А.А. Авраменко, Б.И. Басок, Б.В. Давыденко // Промышленная теплотехника. — 2011. — Т. 33, № 2. — С. 11-18. — Бібліогр.: 7 назв. — рос. |
Institution
Digital Library of Periodicals of National Academy of Sciences of Ukraine| id |
nasplib_isofts_kiev_ua-123456789-60312 |
|---|---|
| record_format |
dspace |
| spelling |
Тыринов, А.И. Авраменко, А.А. Басок, Б.И. Давыденко Б.В. 2014-04-14T06:29:32Z 2014-04-14T06:29:32Z 2011 Моделирование микротечений методом решеток Больцмана / А.И. Тыринов, А.А. Авраменко, Б.И. Басок, Б.В. Давыденко // Промышленная теплотехника. — 2011. — Т. 33, № 2. — С. 11-18. — Бібліогр.: 7 назв. — рос. 0204-3602 https://nasplib.isofts.kiev.ua/handle/123456789/60312 532.5: 536.24 Рассмотрен подход к анализу микро- и макротечений с использованием метода решеток Больцмана. Получены профили скорости на начальном участке для микроканала без проскальзывания и с условием проскальзывания на стенке. Проанализировано влияние числа Кнудсена на длину начального гидродинамического участка и гидравлическое сопротивление канала. Розглянуто підхід до аналізу мікро- й макротечій з використанням методу решіток Больцмана. Отримано профілі розгінної течії для мікроканалу без проковзування та з умовою проковзування на стінці. Проаналізовано вплив числа Кнудсена на довжину початкової гідродинамічної ділянки й гідравлічний опір каналу. The analysis for micro- and macroflows using the lattice Boltzmann method is considered. The velocity profiles in the inlet length of a microchannel with slip and nonslip boundary conditions are obtained. The Knudsen number effect on the inlet length and channel hydraulic resistance is analysed. ru Інститут технічної теплофізики НАН України Промышленная теплотехника Тепло- и массообменные процессы Моделирование микротечений методом решеток Больцмана Modelling of microflows by the method of Boltzman lattices Article published earlier |
| institution |
Digital Library of Periodicals of National Academy of Sciences of Ukraine |
| collection |
DSpace DC |
| title |
Моделирование микротечений методом решеток Больцмана |
| spellingShingle |
Моделирование микротечений методом решеток Больцмана Тыринов, А.И. Авраменко, А.А. Басок, Б.И. Давыденко Б.В. Тепло- и массообменные процессы |
| title_short |
Моделирование микротечений методом решеток Больцмана |
| title_full |
Моделирование микротечений методом решеток Больцмана |
| title_fullStr |
Моделирование микротечений методом решеток Больцмана |
| title_full_unstemmed |
Моделирование микротечений методом решеток Больцмана |
| title_sort |
моделирование микротечений методом решеток больцмана |
| author |
Тыринов, А.И. Авраменко, А.А. Басок, Б.И. Давыденко Б.В. |
| author_facet |
Тыринов, А.И. Авраменко, А.А. Басок, Б.И. Давыденко Б.В. |
| topic |
Тепло- и массообменные процессы |
| topic_facet |
Тепло- и массообменные процессы |
| publishDate |
2011 |
| language |
Russian |
| container_title |
Промышленная теплотехника |
| publisher |
Інститут технічної теплофізики НАН України |
| format |
Article |
| title_alt |
Modelling of microflows by the method of Boltzman lattices |
| description |
Рассмотрен подход к анализу микро- и макротечений с использованием метода решеток Больцмана. Получены профили скорости на начальном участке для микроканала без проскальзывания и с условием проскальзывания на стенке. Проанализировано влияние числа Кнудсена на длину начального гидродинамического участка и гидравлическое сопротивление канала.
Розглянуто підхід до аналізу мікро- й макротечій з використанням методу решіток Больцмана. Отримано профілі розгінної течії для мікроканалу без проковзування та з умовою проковзування на стінці. Проаналізовано вплив числа Кнудсена на довжину початкової гідродинамічної ділянки й гідравлічний опір каналу.
The analysis for micro- and macroflows using the lattice Boltzmann method is considered. The velocity profiles in the inlet length of a microchannel with slip and nonslip boundary conditions are obtained. The Knudsen number effect on the inlet length and channel hydraulic resistance is analysed.
|
| issn |
0204-3602 |
| url |
https://nasplib.isofts.kiev.ua/handle/123456789/60312 |
| citation_txt |
Моделирование микротечений методом решеток Больцмана / А.И. Тыринов, А.А. Авраменко, Б.И. Басок, Б.В. Давыденко // Промышленная теплотехника. — 2011. — Т. 33, № 2. — С. 11-18. — Бібліогр.: 7 назв. — рос. |
| work_keys_str_mv |
AT tyrinovai modelirovaniemikrotečeniimetodomrešetokbolʹcmana AT avramenkoaa modelirovaniemikrotečeniimetodomrešetokbolʹcmana AT basokbi modelirovaniemikrotečeniimetodomrešetokbolʹcmana AT davydenkobv modelirovaniemikrotečeniimetodomrešetokbolʹcmana AT tyrinovai modellingofmicroflowsbythemethodofboltzmanlattices AT avramenkoaa modellingofmicroflowsbythemethodofboltzmanlattices AT basokbi modellingofmicroflowsbythemethodofboltzmanlattices AT davydenkobv modellingofmicroflowsbythemethodofboltzmanlattices |
| first_indexed |
2025-11-24T17:56:53Z |
| last_indexed |
2025-11-24T17:56:53Z |
| _version_ |
1850490422161309696 |
| fulltext |
ISSN 0204-3602. Пром. теплотехника, 2011, т. 33, №2 11
ТЕПЛО- И МАССООБМЕННЫЕ ПРОЦЕССЫ
Розглянуто підхід до аналізу
мікро- й макротечій з використан-
ням методу решіток Больцмана.
Отримано профілі розгінної течії
для мікроканалу без проковзу-
вання та з умовою проковзуван-
ня на стінці. Проаналізовано
вплив числа Кнудсена на довжи-
ну початкової гідродинамічної
ділянки й гідравлічний опір кана-
лу.
Рассмотрен подход к анализу
микро- и макротечений с использо-
ванием метода решеток Больцмана.
Получены профили скорости на на-
чальном участке для микроканала
без проскальзывания и с условием
проскальзывания на стенке. Проана-
лизировано влияние числа Кнудсена
на длину начального гидродинами-
ческого участка и гидравлическое
сопротивление канала.
The analysis for micro- and
macroflows using the lattice Boltzmann
method is considered. The velocity
profiles in the inlet length of a
microchannel with slip and nonslip
boundary conditions are obtained. The
Knudsen number effect on the inlet
length and channel hydraulic resistance
is analysed.
УДК 532.5: 536.24
Тыринов А.И., Авраменко А.А., Басок Б.И., Давыденко Б.В.
Институт технической теплофизики НАН Украины
МОДЕЛИРОВАНИЕ МИКРОТЕЧЕНИЙ МЕТОДОМ РЕШЕТОК БОЛЬЦМАНА
e – координаты положений решетки;
f – функция распределения;
f e – функция равновесного распределения;
h – полуширина канала;
k – постоянная Больцмана;
L – ширина канала;
m – масса молекулы;
p – давление;
R – универсальная газовая постоянная;
t – время;
T – температура;
u – скорость течения;
w – весовые коэффициенты;
x,y – координаты;
Λ – длина свободного пробега;
ν – кинематическая вязкость;
ρ – плотность;
τ – время релаксации;
υ – скорость молекулы;
Kn=Λ/L – число Кнудсена.
Типы переменных обозначаются
диакритическими знаками:
u~ – размерная переменная;
u – безразмерная переменная;
u – «решетчатая» переменная.
Введение
В настоящее время, в связи с развитием на-
нотехнологий, резко возрос интерес к модели-
рованию микро- и нанопотоков. Разработан ряд
методов, которые соответствуют различным
линейным и временным масштабам потока.
Можно отметить следующую классификацию.
В микродиапазоне (microscale), который услов-
но имеет верхнюю границу 10 nm и 1 ns исполь-
зуются методы молекулярной динамики (MD),
в которых учитываются только консервативные
(потенциальные) межмолекулярные валентные
и невалентные силы [1,2]. В мезоскопическом
диапазоне (mesoscale – 10 nm …1 μm и 1 ns
…10 ms) используется метод динамики дисси-
пативных частиц (DPD), который учитывает не
только невалентные потенциальные силы (по-
тенциал Ленарда-Джонса), но и диссипатив-
ные и стохастические взаимодействия [2]. Для
этого же диапазона масштабов используется и
метод решеток Больцмана, основанный на ре-
шетчатом уравнении Больцмана [3,4]. В настоя-
щем исследовании данный метод используется
для расчета гидродинамики микропотока в ка-
нале со скользящими граничными условиями
на стенке. Такой режим течения соответствует
диапазону чисел Кнудсена Kn = 0,1…0,001.
Метод решеток Больцмана (Lattice-
Boltzmann Method, LBM), являющийся развити-
ем идеи клеточных автоматов, появился в конце
20-го века. Одним из его достоинств является
возможность выбора оператора столкновений,
отражающего особенности микродинамики ис-
следуемой среды. Таким образом, удается мо-
ISSN 0204-3602. Пром. теплотехника, 2011, т. 33, №212
ТЕПЛО- И МАССООБМЕННЫЕ ПРОЦЕССЫ
делировать явления, с которыми плохо справ-
ляются традиционные методы, основанные
непосредственно на имитации макроскопиче-
ских свойств жидкости. С другой стороны, этот
метод можно рассматривать как дискретизацию
уравнения Больцмана [5], что позволяет решать
его численно.
Решетчатое уравнение Больцмана
Метод решеток Больцмана основан на дис-
кретизации уравнения Больцмана в приближе-
нии BGK [6]
( )1 e
i
i
f f f f
t x
∂ ∂
+ υ ⋅ = −
∂ ∂ τ
% %
% %%
% % %
, (1)
названном так по фамилиям авторов
(Bhatnagar-Gross-Krook). В уравнении (1) τ% –
время релаксации, а ef% – равновесное распре-
деление молекул, соответствующее распреде-
лению Максвелла
( )2
3/2
exp
2 2
e
m umf
kT kT
υ− = ρ − π
r r% %
% % . (2)
Уравнение (1) в виде конечных разностей будет выглядеть следующим образом
( ) ( ) ( ) ( ), , , , , , , ,
.
e
i i i i i i i i i i i
i
i
f x x t t f x x t f x x t f x t f f
t x
+ ∆ υ + ∆ − + ∆ υ + ∆ υ − υ −
+ υ =
∆ ∆ τ
% % % % % %% % % % %% % % %% % % % % % %
%
% % %
(3)
Для i-мерного пространства число векторов скоростей iυ% является конечным. Принимая, что
положение i ix t+ υ ∆%%% должно занимать фиксированное положение в пространственной решетке,
i ix t∆ = υ ∆%%% . Кроме того, так как вектор скорости фиксирован, то ( ) ( ), , ,i i if x t f x tυ =% %% %%% % . С учетом сказан-
ного, формула (3) примет вид
( ) ( ), , e
i i if x x t t f x t f f
t
+ ∆ + ∆ − −
=
∆ τ
% % % %% % %% % %
% %
. (4)
Раскладываем левую часть выражения (4) в ряд Тейлора до второго порядка
( )31
2
e
i i j i
i i i
f f f f f f f ft
t x t t x t t x
∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ −
+ υ + ∆ + υ + υ + υ +Ο ∂ = ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ τ
% % % % % % % %
%% % % %
% % % % %% % % %
. (5)
После подстановки из (1) значения для f% получим
( )3
2
e e e e
i i j i
i i i
f f t f f f f f f
t x t t x t t x
∂ ∂ ∆ ∂ ∂ ∂ ∂ ∂ ∂ − + υ − τ− + υ + υ + υ +Ο ∂ = ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ τ
% % % % % % % %%
% % % %%
% % % % %% % % %
. (6)
Умножив уравнение (6) на iυ% , проинтегрируем полученное уравнение по υ%
2
e e e e
j i i j i
i i i
f f t f f f f f fd d
t x t t x t t x
∂ ∂ ∆ ∂ ∂ ∂ ∂ ∂ ∂ − υ + υ − τ− + υ + υ + υ υ = υ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ τ
∫ ∫
% % % % % % % %%
% % % % % % %%
% % % % %% % % %
. (7)
Используя для ef% выражение из (2) и проинтегрировав (7) получим
( ) 1
2
ji i i
j
j i j j i
uu u uRTu tRT
t x x x x x
∂∂ ∂ ∂∂ ρ ∂ ρ +ρ = − + ρ τ− ∆ + ∂ ∂ ∂ ∂ ∂ ∂
%%% % %
%% % %%
% % % % % %
, (8)
где / tτ = τ ∆%% , R = k/m. Правая часть уравнения (7) будет равна нулю, так как согласно [5] интеграл
от интеграла столкновений, умноженного на инвариант столкновения υ% равен нулю. При выводе
уравнения (8) учитывались соотношения
ISSN 0204-3602. Пром. теплотехника, 2011, т. 33, №2 13
ТЕПЛО- И МАССООБМЕННЫЕ ПРОЦЕССЫ
, .fd u fdρ = υ ρ = υ υ∫ ∫% %% % % % %% (9)
Сравнивая уравнение Навье-Стокса с (8), мож-
но увидеть, что кинематическая вязкость опре-
деляется как
1
2
tRT ν = τ − ∆
%% . (10)
Обезразмеривая (8) по микромасштабам решет-
ки с учетом того, что средняя скорость частиц
3 3A
kT RT
m
υ = =% , (11)
получим безразмерную форму в «решетчатых»
переменных для несжимаемого потока
1 1 ,
3 2
i i
j
j i
ji
j j i
u u pu
t x x
uu
x x x
∂ ∂ ∂
+ = − +
∂ ∂ ∂
∂∂∂ + τ − + ∂ ∂ ∂
(12)
где
2, , ,i i
i i
A A
u xp tu p x t
x t
= = = =
υ ρυ ∆ ∆
%% %%
%% % % %
. (13)
Из (12) получаем, что вязкость для дискретной
«решетчатой» модели
1 1 2 1
3 2 6
τ − ν = τ − =
. (14)
Связь между безразмерными и «решетчатыми» переменными
Для определения времени релаксации τ необходимо определить связь между ним и макромас-
штабами потока. Для этого обезразмерим уравнение Навье – Стокса
22 2
0 0 0 0 00 0 0 0
2
0
ji i i
j
j i j j i
uu u up
uU U U U UU U U U
x x x x xt L U L L
L LL L L
∂ ∂ ∂ ∂ ∂ ρ ν ∂ + = − + + Ω ∂ ∂∂ ∂ ∂ ∂ Ω
%% % %%
% % % % %% % % %% % %
% % % % %% % % % % %
% %% % % %
, (15)
где 0 , ,U LΩ% % % – макромасштабы скорости, времени и длины соответственно. Макромасштабы вре-
мени и длины связаны с «решетчатыми» масштабами
/ , /iter t xt N x L N L∆ =Ω = Ωδ ∆ = = δ% % % % . (16)
После преобразований выражения (15) получим
0
1
Re
ji i i
j
j i j j i
uu u uL pu
U t x x x x x
∂∂ ∂ ∂∂ ∂
+ = − + + Ω ∂ ∂ ∂ ∂ ∂ ∂
%
% %
, (17)
где
0
2
0 0
Re , , , ,i i
i i
U L u xp tp u t x
U U L
= = = = =
ν ρ Ω
% % %% %%
% % % %% %
. (18)
Примем, что 0U LΩ = и перейдем в (18) от безразмерных переменных к «решетчатым».
2 2 2
1 1
Re
jx i x i x i
j
t t j t i x t j j i
uu u upu
t x x x x x
∂δ ∂ δ ∂ δ ∂∂ ∂
+ = − + + δ ∂ δ ∂ δ ∂ δ δ ∂ ∂ ∂
. (19)
Сравниваем (12) и (19) с учетом (14), получим
}0
0
2
1 1 1 1 1 1 /
Re Re Re Re / Re
A
U
t t
x x x x
U Nt L L N
x x t
υ
δ δ ∆ Ω
ν = = = = =
δ δ δ Ω ∆ δ ∆ ∆
&&&
%
% % %%
% %% %123
, (20)
где 0 0 / AU U= υ% % .
ISSN 0204-3602. Пром. теплотехника, 2011, т. 33, №214
ТЕПЛО- И МАССООБМЕННЫЕ ПРОЦЕССЫ
Таким образом, проведенный анализ пока-
зывает эквивалентность «решетчатого» уравне-
ния Больцмана (1) безразмерному уравнению
Навье-Стокса.
Метод решеток Больцмана
Для решения уравнения (1) был предложен
метод решеток Больцмана. Одни из первых ра-
бот, где был рассмотрен данный метод, это ра-
боты [7, 8].
Основная идея метода решеток Больцмана
состоит в том, что положение частиц предпо-
лагается дискретным, т.е. частица может за-
нимать положение только в узлах заданной
решетки. Форма ячеек решетки может быть
выбрана произвольно, но для декартовых коор-
динат наиболее приемлемой является прямоу-
гольная решетка.
При описании решеток используются обо-
значения в виде DnQi, где n – размерность ре-
шетки, а i – количество скоростных каналов.
Для двумерных задач обычно используется дву-
мерная решетка D2Q9 с девятью скоростными
каналами (рис. 1), а для трехмерных – D3Q15 с
пятнадцатью.
Для дискретной решетки с количеством
Рис. 1. Скоростные каналы решетки D2Q9.
скоростных каналов равным α выражения для
расчета макропараметров имеют следующий
вид
2
, ,
1 ( ) .
3
f u f
p u f
α α α
α α
α α
α
ρ = ρ = υ
= υ −
∑ ∑
∑
(21)
Для выполнения расчетов с помощью моде-
ли решеток Больцмана необходимо определить
выражение для равновесного распределения f e,
зависящее от выбранной пространственной ре-
шетки.
Распределение Максвелла (2) раскладыва-
ем в ряд Тейлора по текущей скорости потока
u, так как она намного меньше скорости моле-
кул (υ >> u). Полученное выражение имеет вид
3/2 2 2 2 2
exp 1
2 2 2 2
e m m u u uf
kT kT RT RT RT
υ υ⋅ υ = ρ − − + + π
. (22)
Согласно [9] вычисление таких макровеличин, как плотность, скорость и энергия в трехмерном
пространстве эквивалентно вычислению интеграла
3/2 2 2 2 2
( ) ( ) exp 1
2 2 2 2
e m m u u uI f d d
kT kT RT RT RT
υ υ⋅ υ = ψ υ υ = ρ ψ υ − − + + υ π
∫ ∫ , (23)
где ψ(υ) является полиномиальной функцией от υ. Для решетки D2Q9, показанной на рис. 1, при-
мем функцию ψ(υ) в виде
n
y
m
xnm υυυψ =)(, , (24)
где υx и υy – пространственные x и y компоненты скорости молекул υ соответственно. После под-
становки (24) в (23) интеграл (23) примет вид
( ) ( ) 2 22
1 1 2 1 1 2 22 2
2 1
2 2
m n x m n y m n x m n x y m n m n
m n
u I I u I I u I I u u I I u I IuI RT I I
RT RTRT
+ + + + + + +
+ + + ρ = − + + π
, (25)
где
2
,
2
m
mI e d
RT
+∞ −ζ
−∞
υ
= ζ ζ ζ =∫ . (26)
ISSN 0204-3602. Пром. теплотехника, 2011, т. 33, №2 15
ТЕПЛО- И МАССООБМЕННЫЕ ПРОЦЕССЫ
Используя формулу Эрмита для интегрирова-
ния, получим
3
1
m
m j j
j
I
=
= ω ζ∑ . (27)
Для трех значений переменной
1 2 33 / 2, 0, 3 / 2 ζ = − ζ = ζ = (28)
весовые коэффициенты будут равны
1 2 3, 2 ,
6 3 6
π π π
ω = ω = ω = . (29)
Таким образом, интеграл (25) примет следую-
щий вид
2 23
, ,
, 2
, 1
( ) ( )
( ) 1
2( ) 2
i j i j
i j i j
i j
u u uI
RT RT RT=
υ ⋅ υ ⋅ρ
= ωω ψ υ + + −
π
∑ , (30)
где , ( , ) 2 ( , )i j i j i jRTυ = υ υ = ζ ζ . Следовательно, сравнивая выражения (21) и (30), функцию равно-
весного распределения можно записать в виде выражения
2 2
, ,
, 2
( ) ( )
1
2( ) 2
i j i j i je
i j
u u uf
RT RT RT
ω ω υ ⋅ υ ⋅
= ρ + + −
π
. (31)
Определим координаты положений решетки
(0,0), 0
(cos ,sin ) , ( 1) / 2, 1, 2,3, 4
( 5) / 2 / 4, 5,6,7,82(cos ,sin ) ,
e c
c
α α α α
αα α
α =
= Θ Θ Θ = α − π α =
Θ = α − π + π α =Θ Θ
, (32)
и весовые коэффициенты для каждого направ-
ления
4 / 9, 0
1/ 9, 1,2,3,4
1/ 36, 5,6,7,8
i j
w
α
α =
ωω = = α =π α =
. (33)
Тогда, учитывая выражение для скорости звука
в решетке
2
3
2
s
cRT c = = , (34)
получим выражение для равновесного распре-
деления частиц при использовании решетки
D2Q9
2 2
2 4 2
3( ) 9( ) 31
2 2
e e u e u uf w
c c c
α α
α α
⋅ ⋅
= ρ + + −
. (35)
Граничные условия
Выбор граничных условий при численном
моделировании имеет первостепенное значе-
ние, поскольку они влияют на его точность и
сходимость. Граничные условия являются не-
простой проблемой в методе решеток Больцма-
на. Трудности возникают из того, что известна
только макроскопическая информация на гра-
ницах расчетной области (например, отсут-
ствие скольжения на стенках). Для реализации
этих граничных условий необходимо перевести
эту макроскопическую информацию в соответ-
ствующую микроскопическую функцию рас-
пределения.
В настоящее время существуют два типа
граничных условий для метода решеток Боль-
цмана. В первом подходе граничные узлы яв-
ляются «мокрыми», то есть являются частью
жидкости. Следовательно, совокупности ча-
стиц такого узла удовлетворяют результатам
разложения Чепмена-Энскога. Они могут быть
разложены на равновесные и неравновесные
и связаны с макроскопическими переменны-
ми потока. В подходе "обратного отскока" гра-
ничные узлы расположены вне жидкости. Они
осуществляют динамику обратного отскока,
т.е. значение вероятности для каждого направ-
ISSN 0204-3602. Пром. теплотехника, 2011, т. 33, №216
ТЕПЛО- И МАССООБМЕННЫЕ ПРОЦЕССЫ
ления частицы копируется в противоположную
вероятность соседней частицы. Поскольку эти
узлы не являются частью жидкости, для них
невозможно вычислить значения макроскопи-
ческих переменных.
В 1997 году был предложен способ задания
значений скорости и давления на «мокрых»
границах [9,10]. Рассмотрим подробно прин-
цип задания значений скорости с компонента-
ми (ux, uy) для нижней горизонтальной границы
прямоугольной расчетной области. После рас-
чета течения значения распределений для f0, f1,
f3, f4, f7, f8 известны. Необходимо определить
выражения для плотности ρ и f2, f5, f6. Из (21)
для решетки D2Q9 можно записать следующие
зависимости
0 1 2 3 4 5 6 7 8f f f f f f f f fρ = + + + + + + + + , (36)
1 3 5 6 7 8xu f f f f f fρ = − + − − + , (37)
2 4 5 6 7 8yu f f f f f fρ = − + + − − . (38)
Формулу для расчета плотности можно най-
ти из (36) и (38)
[ ]0 1 3 4 7 8
1 2( )
1 y
f f f f f f
u
ρ = + + + + +
− , (39)
но для решения системы относительно f2, f5, f6
необходимо еще одно уравнение. Его можно
получить из предположения, что метод обрат-
ного отскока применим для неравновесной ча-
сти распределения частиц в направлении нор-
мали к поверхности. В данном случае
2 2 4 4
e ef f f f− = − . (40)
Учитывая (32) – (38) , получим выражения
для f2, f5, f6
2 4
5 7 1 3
6 8 1 3
2 ,
3
1 1 1( ) ,
2 2 6
1 1 1( ) .
2 2 6
y
x y
x y
f f u
f f f f u u
f f f f u u
= + ρ
= − − + ρ + ρ
= + − − ρ + ρ
(41)
Аналогичный подход позволяет задать значе-
ние давления на границе.
Граничное условие, известное как «откры-
тый выход», реализует нулевой градиент ско-
рости на выходе канала. Допустим, что выход
расположен на вертикальной границе с нор-
мальным компонентом скорости ux и номером
ячейки n. Для первого порядка точности необ-
ходимо установить 1n n
x xu u −= и затем пересчи-
тать функции распределения для f3, f6, f7 , ис-
пользуя (42). Если необходим второй порядок
точности, то значение скорости определяется
из выражения 1 2(4 ) / 3n n n
x x xu u u− −= − .
3 1
6 8 2 4
7 5 2 4
2 ,
3
1 1( ) ,
2 6
1 1( ) .
2 6
n
x
n
x
f f
f f f f u
f f f f u
= + ρ
= − − − ρ
= + − − ρ
(42)
Несжимаемые течения
Описанная выше модель решеток Больцма-
на, основанная на равновесном распределении
(35), является сжимаемой и при моделировании
несжимаемых сред дает погрешность резуль-
тата, вызванную сжимаемостью. Для расчета
несжимаемых течений (таких, например, как
течение Пуазейля) необходимо модифициро-
ванное уравнение равновесного распределе-
ния, имеющее вид
2 2
2 4 2
3( ) 9( ) 3
2 2
e e u e u uf w
c c c
α α
α α
⋅ ⋅
= ρ+ + −
. (43)
Тогда значения макропараметров будут вычис-
ляться следующим образом
, ,f u fα α α
α α
ρ = = υ∑ ∑ (44)
а значение плотности при задании скорости на
границе вместо выражения (39) будет опреде-
ляться как
[ ]0 1 3 4 7 82( ) .yu f f f f f fρ = + + + + + + (45)
Моделирование микротечений
При моделировании микротечений значе-
ние времени релаксации определяется не по
вязкости среды, а по числу Кнудсена. Исполь-
зуя уравнение Максвелла для вязкости
ISSN 0204-3602. Пром. теплотехника, 2011, т. 33, №2 17
ТЕПЛО- И МАССООБМЕННЫЕ ПРОЦЕССЫ
3
AΛυ
ν =% (46)
и выражение (10), находим связь между числом
Кнудсена и временем релаксации
1
2yKnNτ = + . (47)
где Ny – количество элементов решетки попе-
рек канала.
Кроме того, в диапазоне значений числа
Кнудсена Kn = 0,001÷0,1 условие прилипания
на стенках не выполняется. Скорость проскаль-
зывания на стенке определяется выражением
x
xw
uu
y
∂
= Λ
∂
. (48)
В «решетчатых» переменных значение скорос-
ти проскальзывания вычисляется следующим
образом
(1
y
xw x
y
Kn N
u u
Kn N )
′=
+
, (49)
где u'x – продольная компонента скорости в
соседней со стенкой ячейке.
Результаты моделирования
Было проведено моделирование гидродина-
мики разгонного течения в микроканале с соот-
ношением сторон, равном 15 (решетка 600х40
элементов). Исследовался диапазон чисел
Кнудсена Kn = 0,001…0,1. Входная скорость в
«решетчатых» единицах – 0,01.
Метод решеток Больцмана был реализо-
ван в среде GNU Octave, совместимой со сре-
дой MATLAB. Время счета одного варианта
на 4-ядерном компьютере с частотой 2,4 ГГц и
объемом оперативной памяти 4 Гб составило
16 минут.
На рис. 2 представлены профили скорости
разгонного течения без проскальзывания (а) и
для значения числа Кнудсена Kn = 0,1 (б).
Из рисунков видно, что стабилизация по-
тока при увеличении числа Кнудсена насту-
пает раньше. Это обусловлено ослаблением
взаимодействия потока со стенками канала.
Исследования показали, что для Kn = 0,1 уча-
Рис. 2. Профили скорости разгонного
течения без проскальзывания (а)
и при значении Kn = 0,1 (б) для разных
расстояний от входа, нормированных
шириной канала:
1 – x/L = 0; 2 – 0,25; 3 – 1,25; 4 – 15.
а)
б)
сток гидродинамической стабилизации потока
уменьшается примерно на 30 % по сравнению
со случаем течения без проскальзывания.
Следует также отметить, что по мере раз-
вития потока скорость на стенке канала умень-
шается, что вызвано уменьшением локального
градиента скорости по мере удаления от входа
канала.
ISSN 0204-3602. Пром. теплотехника, 2011, т. 33, №218
ТЕПЛО- И МАССООБМЕННЫЕ ПРОЦЕССЫ
На рис. 3 показаны установившиеся профи-
ли скорости на выходе канала для различных
значений числа Кнудсена. Из рисунка видно,
что при увеличении числа Кнудсена значения
максимальной скорости на выходе уменьша-
ются, а скорости проскальзывания на стенках
растут (скорость на стенке пропорциональна
значению числа Кнудсена). При этом среднеин-
тегральная скорость в канале сохраняется.
дравлического сопротивления была получена
на основе приближенного аналитического ре-
шения в работе [11].
Выводы
Для расчета микротечений предложен ме-
тод решеток Больцмана. Данный метод явля-
ется альтернативой подходу, использующему
уравнение Навье-Стокса. Он позволяет моде-
лировать как макротечения (Kn = 0), так и ми-
кротечения (Kn > 0).
Используя метод решеток Больцмана, про-
веден расчет разгонного микротечения в канале
с соотношением сторон 1:15. На основе полу-
ченных результатов построены профили ско-
ростей в различных сечениях канала. Установ-
лено, что значение числа Кнудсена уменьшает
длину начального гидродинамического участ-
ка. Показано, что увеличение степени разре-
женности потока ведет к уменьшению гидрав-
лического сопротивления.
ЛИТЕРАТУРА
1. J.M. Haile. Molecular dynamics simulation,
NY.: Wiley, 1992. – 512 p.
2. G. Karniadakis, A. Beskok, Aluru N.
Microflows and Nanoflows Fundamentals and
Simulation, NY:. Springer, – 2005. – 817 p.
3. Brian J.N. Wylie Application of two-
dimensional cellular automaton lattice-gas models
to the simulation of hydrodynamics. – University
of Edinburgh, 1990.
4. Maxwell J.B. Lattice Boltzmann methods
for interfacial wave modeling. – University of
Edinburgh, 1997.
5. Черчиньяни К. Теория и приложения
уравнения Больцмана. – М.: Мир. – 1978. –
492 c.
6. P.L. Bhatnagar, E.P. Gross, M. Krook
(1954). "A Model for Collision Processes in
Gases. I. Small Amplitude Processes in Charged
and Neutral One-Component Systems". Physical
Review 94. – Р. 511–525.
7. Shan X, Chen H. (1993) Lattice Boltzmann
model for simulating flows with multiple phases
and components. Phys Rev E 47. –Р. 1815–1819.
Рис. 3. Установившиеся профили скорости
для различных значений чисел Кнудсена:
1 – Kn = 0,1; 2 – 0,05; 3 – 0,01.
На основе полученных результатов была
проведена оценка значений коэффициента ги-
дравлического сопротивления для значений
числа Кнудсена Kn = 0,0…0,1. При этом разни-
ца давлений рассчитывалась по формуле (21).
Анализ результатов показал, что для относи-
тельного коэффициента гидравлического соп-
ротивления справедливо соотношение
0
1
1 6Kn
λ
=
λ +
, (50)
где λ0 – коэффициент гидравлического сопро-
тивления для течения без проскальзывания на
стенке. Из формулы (50) следует, что рост ско-
рости проскальзывания ведет к уменьшению
коэффициента гидравлического сопротивления
по гиперболической зависимости. Подобная
формула для относительного коэффициента ги-
|