Физика плазмы, 2020, T. 46, № 8, стр. 721-729

1. ВВЕДЕНИЕ

Исследования процессов энергетического обмена в неоднородных системах взаимодействующих частиц вызывает значительный интерес в различных областях науки и техники (физике плазмы, биологии, физике полимеров и т.д.) [1–5]. Ряд актуальных вопросов касается особенностей физических характеристик в ансамблях неидентичных частиц, имеющих различный характер парного взаимодействия, заряды, размеры, диэлектрическую проницаемость и т.д. [1–8].

Пылевая (комплексная) плазма представляет собой ионизованный газ, содержащий заряженные частицы вещества микронных размеров (пыль). Такая плазма широко распространена в природе и образуется в ряде технологических процессов [1–3]. Большинство теоретических и численных работ, посвященных исследованию свойств пылевой плазмы, имеют дело с идентичными пылевыми частицами, поскольку такие системы легче поддаются математическому описанию и более просты для понимания. Тем не менее, в реальных условиях пылевые структуры редко содержат идентичные частицы. Даже в случае лабораторных исследований монодисперсные пылевые частицы могут иметь различную величину заряда и/или стохастической кинетической энергии в зависимости от их пространственного положения [1, 2].

Большинство лабораторных исследований пылевой плазмы проводится в газовых разрядах различных типов [9–12]. Обычно, в центре газоразрядных камер наблюдается некоторое превышение концентрации ионов плазмы над концентрацией ее электронной компоненты [13]. Данное обстоятельство приводит к формированию эффективных ловушек для отрицательно заряженных частиц пыли [1, 2]. Стохастическая кинетическая энергия пылевых частиц (их “кинетическая температура”) в таких условиях может достигать ~0.2–5 эВ, что значительно выше температуры окружающего их газа. Механизмы такого “аномального разогрева” пылевых частиц обычно связывают с временными и/или пространственными изменениями их зарядов, или положения в объеме неоднородной плазмы [14–19]. Так как заряд пылевой частицы определяется локальными параметрами плазмы в ее окрестности, мощность источников подкачки энергии, а, соответственно, и “кинетическая температура” частицы могут существенно изменяться в пространстве [1, 2]. Источниками неравномерного нагрева системы частиц также могут являться неоднородное распределение температуры окружающего газа, лазерное излучение, используемое для диагностики и т.д.

Флуктуации зарядов пылевых частиц, вызванные случайной природой ионных и электронных токов, заряжающих эти частицы, присущи любым типам плазмы [1, 2]. В условиях лабораторной газоразрядной плазмы дополнительная стохастическая кинетическая энергия, ΔТ_f, для отдельной пылевой частицы, связанная с этими флуктуациями может быть записана в виде [1, 14‒16]

Здесь Q, М и ν – заряд, масса и коэффициент трения частицы, α ≈ 0.5 – параметр, отвечающий за амплитуду флуктуаций их зарядов, υ ∝ а_d – характерная частота этих флуктуаций для частицы радиусом а_d, а E – напряженность электрического поля в анализируемой системе, необходимая для равновесного положения пылевой частицы в поле действующих сил. В условиях микрогравитации величина ΔТ_f определяется флуктуациями зарядов окружающих частиц пылевого облака, E ≅ E_int ~ Q/d², где d – среднее расстояние между пылевыми частицами. Для наземных лабораторных экспериментов, где основной внешней неэлектрической силой является сила тяжести: E ≅ E_int + E_ext, где E_ext ≅ gM/Q.

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

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

Рассмотрены условия перераспределения стохастической кинетической энергии в двумерных ансамблях с двумя разделенными фракциями неидентичных частиц разных размеров и температур. Причины формирования раздельных фракций для частиц разных размеров обычно связывают с наличием термофоретических сил, ${{F}_{{\text{T}}}} = {{C}_{{\text{T}}}}a_{d}^{2}\nabla {{T}_{n}}$, или сил ионного увлечения, ${{F}_{{\text{I}}}} = {{C}_{{\text{I}}}}a_{d}^{2}{{u}_{i}}$. Здесь a_d – радиус частицы, ∇T_n – градиент температуры нейтрального газа (который направлен в сторону уменьшения T_n), u_i – дрейфовая скорость ионов (которая направлена по внешнему электрическому полю), а коэффициенты C_T и C_I зависят от внешних условий и параметров окружающей плазмы [1, 20–22]. Указанные силы (F_I, F_Т) зависят от размера частицы иным образом, чем сила электрического поля, F_Е ≈ QE ∝ a_d, действующая на заряженную пылевую частицу. За счет направления электрических полей и градиентов температур газа в обычных условиях газоразрядных камер обе упомянутых силы (F_I, F_Т) будут давать дополнительный вклад в баланс сил, ослабляющий воздействие внешних электрических сил F_Е на пылевую частицу [13]. И этот вклад будет тем больше, чем больше размер пылевых частиц. Таким образом, в случае системы, состоящей из двух фракций частиц, более крупные частицы будут располагаться на оболочке пылевого облака.

В настоящей работе внешнее электрическое поле полагалось линейным. Для представленных здесь теоретических и численных расчетов мы полагали, что плотность материала, ρ, одинакова для различных фракций частиц, ρ₁= ρ₂, т.е. отношение их масс M₁/M₂ ∝ (a_d1/a_d2)³. Заряды частиц задавались согласно приближению ограниченных орбит (Orbit Motion Limited): Q_i ∝ a_di [1, 2]; а их коэффициенты трения согласно свободномолекулярному приближению: ${{\nu }_{i}} \propto a_{{di}}^{2}{\text{/}}{{M}_{i}}$ [23]. С учетом принятых приближений отношение величин дополнительной энергии (1) за счет флуктуаций зарядов пылевых частиц разного размера (a_d1, a_d2) можно оценить как: $\Delta Т_{f}^{1}{\text{/}}\Delta Т_{f}^{2}$ ~ (d₂ /d₁)⁴ для E_int ≫ E_ext; или $\Delta Т_{f}^{1}{\text{/}}\Delta Т_{f}^{2}$ ~ (a_d1/a_d2)² для E_int ≪ E_ext, если основной действующей неэлектрической силой в системе является сила тяжести.

2. ДВЕ ЧАСТИЦЫ (ОСНОВНЫЕ СООТНОШЕНИЯ)

Рассмотрим систему линеаризованных уравнений движения, описывающих отклонения двух частиц (с зарядами Q₁₍₂₎ , массой M₁₍₂₎ и попарным взаимодействием) от их положения равновесия (ξ₁, ξ₂) в поле внешних сил под действием случайной силы F_b1(2), которая является источником стохастической (тепловой) энергии частиц:

где ν₁₍₂₎ – коэффициенты трения частиц за счет их столкновений с нейтралами окружающего газа, а коэффициенты a₁₍₂₎, b зависят от физики решаемой задачи и рассматриваемой степени свободы.

Остановимся на случае вертикальной конфигурации частиц, см. рис. 1. Тогда для их вертикальных смещений ξ₁₍₂₎= z₁₍₂₎: коэффициент a₁₍₂₎ = –(Q₁₍₂₎β_z – F') и b = –F'. Для радиальных смещений ξ₁₍₂₎= r₁₍₂₎: a₁₍₂₎ = –(Q₁₍₂₎β_r- F/d) и b = = F/d. Здесь β_r и β_z – величины градиентов электрического поля в радиальном и вертикальном направлении, соответственно, d – расстояние между частицами, F – сила взаимодействия между частицами, а F' ее производная в направлении z. Для кулоновского взаимодействия: F = Q₁Q₂/d²; F  ' = –2Q₁Q₂/d³.

Задача об устойчивой конфигурации двух идентичных частиц, взаимодействующих с различными типами потенциалов, рассматривалась в работах [24–27]. Было показано, что устойчивость их вертикальной конфигурации определяется соотношением β_r > β_z, в обратном случае (β_r < β_z) – формируется горизонтальная конфигурация частиц. Для двух частиц с разной массой для условий наземных экспериментов (где нельзя пренебречь силой тяжести) условия баланса сил дают [8]

При этом условие формирования неустойчивости в такой системе можно записать как [8]:

Для M₁₍₂₎ ≡ M и Q₁₍₂₎ ≡ Q данное условие приобретает вид β_r < β_z , в соответствие с критерием, изложенным выше. Исследование условий развития неустойчивостей (4) для различных параметров задачи представлены в работе [8].

Для поиска корреляторов скоростей и смещений частиц в системе (2а)–(2б) отметим, что корреляторы случайной силы F_b1(2) подчиняются уравнениям: 〈F_b1〉 = 〈F_b2〉 ≡ 0, 〈F_b1F_b2〉 = 0, 〈F_b1V₂〉 = = 〈F_b2V₁〉 ≡ 0, 〈F_b1ξ₂〉 = 〈F_b2ξ₁〉 ≡ 0, 〈F_b1ξ₁〉 = 〈F_b2ξ₂〉 ≡ 0, 〈F_b1V₂〉 = 〈F_b2V₁〉 ≡0. Здесь и далее угловые скобки 〈 〉 обозначают усреднение по времени при t → ∞. При движении частиц по ограниченным траекториям: 〈ξ₁V₁〉 = 〈ξ₂V₂〉 ≡ 0, а 〈V₁₍₂₎F_b1(2)〉 = $\nu T_{{1\left( 2 \right)}}^{0}$ , где $T_{{1\left( 2 \right)}}^{0}$ – температура тепловых источников [28–31]. Тогда уравнения для корреляторов можно представить в виде [30, 31]

Здесь ${{T}_{1}} = T_{1}^{0} + \delta {{T}_{1}}$, ${{T}_{2}} = T_{2}^{0} + \delta {{T}_{2}}$, где T₁, T₂– температура частиц для равновесного состояния системы, $T_{1}^{0}$, $T_{2}^{0}$ – энергия источников (которая при численном моделировании задачи соответствует их заданной/начальной температуре), а δT₁, δT₂– приращение температуры в процессе установления равновесия.

Обозначим $\Delta T = T_{2}^{0} - T_{1}^{0}$, тогда

где

Для двух идентичных частиц (Q₁₍₂₎ = Q, M₁₍₂₎= = M, ν₁₍₂₎ = ν, a₁₍₂₎ = а) решение системы уравнений (5а)–(5д) для их вертикальной конфигурации можно записать как [8]

В этом случае δT₁ = –δT₂, а при b² ≫ –ν²Мa величина |δT₁₍₂₎| → |△T|/2, т.е. энергия равномерно распределяется между частицами системы, а при ν → ∞ величина |δT₁₍₂₎| → 0.

Иллюстрация зависимостей δT₁/ΔT и δT₂/ΔT от ν₁/ω, где $\omega = {{(Q_{1}^{2}{\text{/}}{{d}^{3}}{{М}_{1}})}^{{1/2}}}$, описывающая перераспределение энергии в вертикальном и радиальном направлениях, для двух частиц с кулоновским взаимодействием при g = 0, β_r/β_z = 4, $T_{1}^{0} < T_{2}^{0}$ показана на рис. 2.

В заключение данного раздела рассмотрим два частных случая для частиц равной массы M₁₍₂₎ = = M при $T_{2}^{0} \ne T_{1}^{0}$: (i) ν₁₍₂₎ = ν, a₁ ≠ a₂; (ii) ν₁ ≠ ν₂, a₁₍₂₎ = а.

В первом случае (ν₁₍₂₎ = ν, a₁ ≠ a₂) получим следующее уравнение баланса:

При этом общая энергия системы будет сохраняться, δT₁ + δT₂ = 0, а при b² ≫ –ν²М(a₁ + a₂), величина |δT₁₍₂₎/△T| → (2 + (a₁ – a₂)²/2b²)^–1.

Отметим, что случай a₁ ≠ a₂ может реализоваться для пылевых частиц в плазме даже в случае их попарного взаимодействия, например, за счет разницы в зарядах частиц, или отличия градиентов внешнего поля в точке их равновесного положения. Данные обстоятельства определяются возможными пространственными изменениями параметров окружающей их плазмы (например, концентраций и температур ионов/электронов).

Учтем уравнение баланса сил (3) для двух частиц при Q₁ ≠ Q₂ и M₁₍₂₎ = M (для g = 0). Тогда величина |δT₁₍₂₎/ΔT| → 0 при (Q₁/Q₂– Q₂/Q₁)²→ ∞. Это связано с уменьшением сил взаимодействия в рассматриваемой системе с увеличением разницы между Q₁ и Q₂, что является следствием условия баланса сил (3), необходимого для равновесного состояния системы. На рис. 3 показаны отношения δT_z/ΔT и δT_r/ΔT в зависимости от Q₁/Q₂, описывающие перераспределение энергии в вертикальном и радиальном направлениях, соответственно. Расчеты были выполнены при β_r/β_z = 4 для частиц с кулоновским взаимодействием, расположенными на расстоянии d ≈ 0.1 см.

Рассмотрим второй случай: M₁₍₂₎ = M, a₁₍₂₎ = а, $T_{2}^{0} \ne T_{1}^{0}$, ν₁ ≠ ν₂. Неравная величина коэффициентов трения, ν₁ ≠ ν₂, в этом случае может возникать, например, за счет градиентов температур и/или давлений окружающего буферного газа. Тогда для частиц с равными зарядами Q₁₍₂₎ ≡ Q получим уравнения баланса в виде

Общая энергия системы не сохраняется: δT₁ + + δT₂ ≠ 0. При b² ≫ –ν₁ν₂ Мa: δT₁ + δT₂ → ΔT (ν₂ – – ν₁)/(ν₁+ ν₂). Поскольку $\Delta T = T_{2}^{0} - T_{1}^{0}$, это означает, что при $T_{2}^{0} > T_{1}^{0}$ и ν₁> ν₂ – система будет терять энергию, в случае $T_{2}^{0} < T_{1}^{0}$ и ν₁> ν₂ – будет приобретать ее. Это связано с тем, что при меньших коэффициентах трения одной из частиц в рассматриваемой системе, она менее активно обменивается энергией с другой частицей. Легко увидеть, что диапазон изменений кинетической температуры всей системы при любых ν₁ и ν₂ будет варьироваться в пределах от $T_{2}^{0}$ до $T_{1}^{0}$.

Зависимость |δT₁/ΔT| и |δT₂/ΔT| от ν₁/ν₂ при β_r/β_z = 4 для частиц с кулоновским взаимодействием, расположенными на расстоянии d ≈ ≈ 0.1 см, показана на рис. 4.

Предельные случаи для двух последних из рассмотренных задач (Q₁/Q₂ ≫ 1 и ν₁/ν₂ ≫ 1) не могут реализоваться в лабораторной практике. Тем не менее, иллюстрации, представленные на рис. 3, 4 качественно демонстрируют изменения величины перераспределяемой энергии в зависимости от параметров задачи.

3. РЕЗУЛЬТАТЫ ЧИСЛЕННОГО МОДЕЛИРОВАНИЯ И ИХ ОБСУЖДЕНИЕ

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

Моделирование проводилось для двух отдельных частиц и для двумерных ансамблей с разделенными фракциями частиц разных размеров и температур (при g = 0). Отношения масс, M₂/M₁, и отношения температур тепловых источников, $T_{2}^{0}{\text{/}}T_{1}^{0}$ (или $T_{2}^{0}{\text{/}}T_{1}^{0}$), для частиц различных фракций изменялись от одного до десяти. Заряды и коэффициенты трения задавались как Q_i ∝ a_di и ${{\nu }_{i}} \propto a_{{di}}^{{ - 1}}$, соответственно; i = 1, 2. Отношение ω₁₍₂₎/ν₁₍₂₎ варьировалось от ~0.7 до ~7, где ${{\omega }_{{1\left( 2 \right)}}} = {\text{ }}$ $ = {\text{ }}{{(Q_{{1\left( 2 \right)}}^{2}{\text{/}}{{d}^{3}}{{М}_{{1\left( 2 \right)}}})}^{{1/2}}}$. Здесь для двумерных бинарных систем: d – среднее расстояние между частицами отдельной фракции. Число частиц в каждой из фракций (N₁, N₂) менялось от 50 до 400.

Температура тепловых источников $T_{{1\left( 2 \right)}}^{0}$ варьировалась в пределах от ~0.2 эВ до ~2 эВ и задавалась одинаковой по степеням свободы: $T_{z}^{0} = T_{x}^{0} = T_{y}^{0}$. В процессе моделирования начальная стохастическая кинетическая энергия (энергия источников) перераспределялась от более “горячих” частиц к менее “горячим”. Во всех случаях наблюдаемые распределения скоростей частиц были близки к Максвелловским функциям. При этом для случая двух частиц фиксировалось неравномерное перераспределение энергий по степеням свободы T_z ≠ T_x = T_y.

Иллюстрации отдельных результатов численных исследований для двух частиц при β_r/β_z = 4, ω₁/ν₁ ≈ 1.4 и d ≈ 0.1 см представлены на рис. 5–7 совместно с аналитическими решениями задачи. Так на рис. 5 показаны функции распределения скоростей для двух идентичных частиц f₁(V₁) и f₂(V₂) в вертикальном и в радиальном направлениях для при $T_{1}^{0}{\text{/}}T_{2}^{0}$ = 10. Согласно результатам, представленным в разд. 2, полная энергия моделируемой системы не менялась (T₂ + T₁) = = ($T_{2}^{0} + T_{1}^{0}$), (δT₁ + δT₂) = 0.

На рис. 6, 7 оказаны функции f₁(V₁) и f₂(V₂) для М₂ = 2М₁, Q₂/Q₁ = 2^1/3, ν₁/ν₂ = 2^1/3 при $T_{2}^{0}{\text{/}}T_{1}^{0}$ = 10 и $T_{1}^{0}{\text{/}}T_{2}^{0}$ = 10, соответственно. В обоих случаях (δT₁ + δT₂) ≠ 0. При этом для $T_{2}^{0}{\text{/}}T_{1}^{0}$ = 10, см. рис. 6, приращение энергии было отрицательным (δT₁ + δT₂) < 0, а (T₂ + T₁)/($T_{2}^{0} + T_{1}^{0}$) ≈ 0.97. Для случая $T_{1}^{0}{\text{/}}T_{2}^{0}$ = 10, см. рис. 7, приращение энергии было положительным (δT₁ + δT₂) > 0, а (T₂ + T₁)/($T_{2}^{0} + T_{1}^{0}$) ≈ 1.02. Величина |δT₁/δT₂| ≅ ≅ ν_2/ν₁, см. формулы (10а), (10б).

Таким образом, результаты численного моделирования для двух частиц согласуются выводами, представленными в разд. 2. Функции распределения скоростей f₁(V₁) и f₂(V₂) полностью соответствовали их температурам, найденным из аналитических соотношений.

Перейдем к двумерным ансамблям, состоящим из двух разделенных фракций частиц разного размера. Частицы в моделируемых системах находились под воздействием электрических сил ловушки, F_Е ∝ a_d, и сторонних сил ${{F}_{{{\text{T(I)}}}}} \propto a_{d}^{2}$, направленных противоположно силе F_Е. Здесь мы рассматриваем упомянутые силы только в первом (линейном) приближении. Поскольку основным предметом исследований является баланс энергии, а подробный анализ природы сил, механизмов и условий разделения частиц по размерам выходит за рамки настоящей работы.

Расчеты проводились для различных отношений |F_T(I)/F_Е|, которые варьировались от 5% до 25% для мелких частиц, и от 10% до ~30% для крупных частиц, соответственно, в зависимости от отношения размеров крупных и мелких частиц. При этом отношение градиентов ∇F₁/∇F₂ суммарных сил, F₁₍₂₎, действующих на частицы изменялось от ~1.05 до ~1.2. Во всех случаях наблюдалось заметное разделение системы по фракциям и более крупные частицы находились на оболочке двумерного облака, см. рис. 8.

На рис. 9, 10 показаны результаты численного моделирования для частиц с параметрами, близкими к параметрам экспериментов в условиях микрогравитации [21, 22]: М₂/М₁ = 8, Q₂/Q₁ = 2, ν₁/ν₂ = 2, $T_{1}^{0} + T_{2}^{0}$ = 4 при ∇F₁/∇F₂ ~ 1.1. Средние расстояния между частицами крупной фракции составляли d ≈ 0.07 см, между легкими частицами – d ≈ 0.035 см. Отношение ω₁/ν₁ ~ 6.8, ω₂/ν₂ ~ 3.4, где ω₁₍₂₎ = ${{(Q_{{1\left( 2 \right)}}^{2}{\text{/}}{{d}^{3}}{{М}_{{1\left( 2 \right)}}})}^{{1/2}}}$. На рис. 9 показана вероятность, f(r), нахождения частиц различной массы М₁₍₂₎ на расстояниях, r, от центра ловушки, а на рис. 10 представлена радиальная зависимость их кинетических температур, Т(r).

Результаты исследования баланса энергии в системе частиц, состоящей из двух фракций, представлены на рис. 11–14 для случая М₂/М₁ = 2, Q₂/Q₁ = 2^1/3, ∇F₁/∇F₂ ~ 1.1, ν₁/ν₂ = 2^1/3 при ν₁ = 5 с^–1 (ω₁/ν₁ ~ 5.6, ω₂/ν₂ ~ 4.5) и ν₁ = 10 с^–1 (ω₁/ν₁ ~ 2.8, ω₂/ν₂ ~ 2.25).

Так на рис. 11, 12 показана вероятность, f(r), нахождения частиц различной массы на разных расстояниях, r, от центра ловушки и радиальная зависимость их кинетических температур, Т(r), при $T_{1}^{0} - T_{2}^{0}$ = 10. А на рис. 13, 14 представлены функции  f(r) и Т(r) для случая $T_{2}^{0} - T_{1}^{0}$ = 10.

Анализ результатов численного моделирования для двумерных ансамблей частиц показал, что вдоль линий раздела двух фракций, когда |δT₁₍₂₎/T₁₍₂₎| > 3% (т.е. ошибок численного эксперимента), величина отношения |δT₁(r)/δT₂(r)| = = ${\text{|}}({{T}_{1}}(r) - T_{1}^{0}){\text{/}}({{T}_{2}}(r) - T_{2}^{0}){\text{|}}$ ≈ ν_2/ν₁, см. рис. 10, 12, 14; где r – расстояние до центра ловушки. При этом дополнительная энергия δT₁₍₂₎ ∝ ΔT = $T_{2}^{0} - T_{1}^{0}$, и уменьшалась с ростом ν₁₍₂₎, см. рис. 12, 14.

4. ЗАКЛЮЧЕНИЕ

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

Выполнен численный анализ процессов перераспределения энергии в двумерных ансамблях неидентичных частиц, состоящих из двух раздельных фракций, формирующихся в электрических полях ловушки (F_Е ∝ a_d) под воздействием дополнительных сил пропорциональных квадрату радиуса частиц (${{F}_{{{\text{T(I)}}}}} \propto a_{d}^{2}$), например, термофоретических сил, или сил ионного увлечения.

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

Данная работа частично поддержана Российским фондом фундаментальных исследований (грант № 18-38-20175), а также Программой Президиума РАН.