Физика плазмы, 2019, T. 45, № 1, стр. 3-13

Осцилляции разряда в стационарном двигателе А. И. Морозова как проявление крупномасштабных мод градиентно-дрейфовой неустойчивости
Е. А. Сорокина, Н. А. Марусов, В. П. Лахин, В. И. Ильгисонис

Е. А. Сорокина ab*, Н. А. Марусов abc, В. П. Лахин a b, В. И. Ильгисонис bd

a НИЦ “Курчатовский институт”
Москва, Россия

b Российский университет дружбы народов
Москва, Россия

c Московский физико-технический институт
Долгопрудный, Россия

d Государственная корпорация по атомной энергии “Росатом”
Москва, Россия

* E-mail: Sorokina_EA@nrcki.ru

Поступила в редакцию 30.05.2018
После доработки 02.08.2018
Принята к публикации 02.08.2018

Полный текст (PDF)

Аннотация

Предложена физическая трактовка явления крупномасштабных осцилляций разряда в стационарном плазменном двигателе А.И. Морозова (СПД), основанная на исследовании глобальных мод градиентно-дрейфовой неустойчивости. Использована идеальная двужидкостная гидродинамическая модель плазмы, включающая эффекты стационарного течения электронов, электронной инерции, а также неоднородностей магнитного поля и плотности плазмы вдоль ускоряющего канала. Для типичных параметров системы типа СПД-100 рассчитаны частоты собственных неустойчивых колебаний и аксиальная структура собственных функций. Полученный спектр характеризуется конечным набором длинноволновых азимутальных мод в нижнегибридном диапазоне частот, локализованных преимущественно в прианодной части двигателя. Показано, что собственные моды могут образовывать волновые пакеты, основные характеристики которых на линейной стадии развития неустойчивости совпадают с параметрами крупномасштабных азимутальных структур типа спиц, наблюдаемых в экспериментах. Исследовано влияние масштабирования геометрии двигателя (длины и ширины канала) на частотные характеристики колебаний и формирование биений.

1. ВВЕДЕНИЕ

Известно, что в процессе работы стационарные плазменные двигатели А.И. Морозова (СПД) “страшно шумят”, испытывая колебания в диапазоне от сверхнизких ~1 кГц до высоких ~10 МГц частот [1]. Для работы плазменных двигателей наиболее значимы низкочастотные (~10 кГц) азимутальные колебания плазмы, наблюдаемые в прианодной области СПД и в других установках холловского типа в виде долгоживущих вращающихся структур типа спиц (spoke в англоязычной литературе). Такие структуры являются потенциальным источником аномальной проводимости электронов поперек внешнего магнитного поля [2] и, как показывают прямые экспериментальные измерения [3], способны переносить до 50% разрядного тока. Помимо малой частоты, эти структуры характеризуются низкими значениями азимутальных волновых чисел $m\sim 1{\kern 1pt} - {\kern 1pt} 5$ и малой скоростью вращения по сравнению со скоростью стационарного вращения электронов в скрещенных полях ${{V}_{\theta }} \ll c[{\mathbf{E}} \times {\mathbf{B}}]{\text{/}}{{B}^{2}}$.

Несмотря на многолетнюю историю исследований, физика генерации крупномасштабных структур в холловской плазме до конца не выяснена. В ранних работах в качестве основного механизма их формирования рассматривались неустойчивости, связанные с ионизационными процессами [4, 5]. Альтернативная концепция связывает появление макроскопических структур с раскачкой длинноволновой части спектра градиентно-дрейфовых колебаний [69].

Градиентно-дрейфовой называют электростатическую неустойчивость, возникающую под действием электрического тока, протекающего поперек магнитного поля при наличии градиента температуры или концентрации компонент плазмы [10, 11]. В англоязычной литературе наряду с термином “градиентно-дрейфовая неустойчивость” (“gradient drift instability”) используется понятие “${\mathbf{E}} \times {\mathbf{B}}$-неустойчивость” (“cross-field instability”), не употребляемое в отечественных публикациях. По своему физическому механизму эта неустойчивость аналогична неустойчивости Рэлея–Тейлора, наиболее известным проявлением которой в физике плазмы служит желобковая неустойчивость. При градиентно-дрейфовой неустойчивости возмущение концентрации электронов в неоднородной плазме и соответствующее возмущение электрического поля нарастают во времени, привнося дополнительный (дрейфовый) вклад в электрический ток поперек магнитного поля. Понятно, что когда этот процесс носит нелокальный характер, в плазме могут формироваться какие-то макроскопические структурные образования (например, неоднородности плотности) или наоборот, могут разрушаться имеющиеся. Эта неустойчивость активно исследуется в физике ионосферы и, в частности, является общепринятым механизмом измельчения нерегулярных структур в F-слое (см., например, [12]).

Градиентно-дрейфовая неустойчивость присуща системе типа СПД [13, 14], и именно из соображений ее подавления профиль магнитного поля в установке выбирают нарастающим от анода к выходу из ускоряющего канала. Локальная теория градиентно-дрейфовой неустойчивости довольно глубоко проработана в литературе, причем как в ставших уже классическими [1517], так и в более современных работах [6, 9, 18] – и включает такие эффекты как конечная скорость ионного потока, инерция электронов, электромагнетизм, конечный ларморовский радиус электронов, неквазинейтральность плазмы. В работах [9, 19] показано, что для описания неустойчивости в холловских двигателях принципиальное значение имеет учет электронной инерции. Инерция электронов приводит к стабилизации коротковолновой части спектра; в противном случае неустойчивым оказывается весь спектр колебаний, а инкремент неустойчивости неограниченно растет с уменьшением длины волны возмущений. В работе [20] в рамках локальной электростатической модели градиентно-дрейфовой неустойчивости рассчитан спектр наиболее неустойчивых мод вдоль оси ускоряющего канала установки СПД-100. Показано, что в холодной плазме неустойчивой является лишь прианодная область канала, характеризуемая раскачкой длинноволновых возмущений. Хотя локальная модель [20] и предсказывает существование длинноволновых азимутальных волн в спектре колебаний, однако частота наиболее неустойчивых мод спектра ~1 MГц значительно превышает частоты экспериментально наблюдаемых структур, что свидетельствует о недостаточном соответствии локальной модели физике реальных процессов, происходящих в плазменных двигателях.

В настоящей работе на примере классического представителя линейки СПД – двигателя СПД-100 – исследуется устойчивость глобальных градиентно-дрейфовых волн в холловской плазме. В рамках идеальной двужидкостной гидродинамической модели, включающей эффекты стационарного течения электронов, электронной инерции, а также неоднородности магнитного поля и плотности плазмы вдоль ускоряющего канала, рассматриваются электростатические колебания в диапазоне частот между ионной ωBi и электронной ωBe циклотронными частотами, распространяющиеся строго поперек внешнего магнитного поля.

Анализ, проведенный в настоящей работе, основан на численном решении граничной задачи для потенциала возмущений вдоль оси ускоряющего канала с нулевыми граничными условиями на электродах. Решениями такой задачи на собственные значения являются собственные функции, описывающие глобальную структуру собственных мод, и собственные числа, формирующие спектр собственных мод. Основной особенностью этого спектра служит наличие нескольких близких собственных значений (частот), имеющих наибольшие мнимые части (инкременты). Совместная раскачка таких мод сопровождается биениями, имеющими существенно меньшую (разностную) частоту по сравнению с частотой исходных колебаний – собственных мод. Тем самым показано, что крупномасштабные низкочастотные азимутальные структуры могут формироваться на линейной стадии развития градиентно-дрейфовой неустойчивости из высокочастотных мод в виде волновых пакетов с частотой огибающих ~10 кГц.

Работа построена следующим образом. Раз-дел 2 посвящен выводу дифференциального уравнения на собственные частоты глобальных градиентно-дрейфовых мод. В разд. 3 представлены характерные профили параметров плазмы в системе типа СПД-100, используемые в ходе расчетов. Раздел 4 содержит результаты численного решения уравнения глобальных мод. Представлены спектр неустойчивых колебаний и структура собственных функций, определена область пространственной локализации мод внутри ускоряющего канала. Также исследована временная эволюция волнового пакета, формируемого глобальными модами, на линейной стадии развития неустойчивости. Представлена картина биений. Раздел 5 посвящен исследованию зависимости характеристик спектра от геометрии ускоряющего канала с целью экстраполяции полученных результатов на другие двигатели семейства СПД. Краткие выводы работы представлены в заключительном разделе.

2. УРАВНЕНИЕ ГЛОБАЛЬНЫХ МОД

Конструкция СПД и основные принципы его работы подробно описаны в [13]. Для упрощения рассмотрения будем считать ускоряющий канал плоским и пренебрежем радиальной неоднородностью физических величин. Введем декартову систему координат $\{ x,y,z\} $ с осью x, направленной вдоль оси ускоряющего канала. Координата y будет обозначать азимутальное, а координата z – радиальное направления соответственно – см. рис. 1. Будем рассматривать базовую модель СПД с чисто аксиальным электрическим полем – ${\mathbf{E}} = E(x){{{\mathbf{e}}}_{x}},$ и радиальным магнитным полем в виде ${\mathbf{B}} = {{B}_{x}}(x,z){{{\mathbf{e}}}_{x}} + {{B}_{z}}(x,z){{{\mathbf{e}}}_{z}} \approx {{B}_{z}}(x){{{\mathbf{e}}}_{z}}$, причем ${\text{rot}}{\mathbf{B}} \approx 0$. Остальные невозмущенные величины также будем считать зависящими только от координаты x.

Рис. 1.

Сечение ускоряющего канала системы типа СПД в плоскости $y = {\text{const}}$. Здесь ${{R}_{1}}$ и ${{R}_{2}}$ – внутренний и внешний радиусы ускоряющего канала, d – длина ускоряющего канала.

Исследование устойчивости проведем в рамках идеальной двужидкостной гидродинамической модели холловской плазмы. Напомним, что магнитное поле в СПД подбирается таким образом, чтобы размер ларморовского радиуса ионов ${{\rho }_{i}}$ значительно превышал длину ускоряющего канала d: ${{\rho }_{i}} \gg d$, в то время как для электронов должно выполняться обратное условие: ${{\rho }_{e}} \ll d$. В этом случае ионы свободно ускоряются электрическим полем вдоль оси канала (в x-направлении) и создают тягу, электроны же дрейфуют в азимутальном направлении y со скоростью электрического дрейфа ${{{\mathbf{V}}}_{E}} = c[{\mathbf{E}} \times {\mathbf{B}}]{\text{/}}{{B}^{2}}$. Замагниченность электронов в такой системе позволяет избежать ограничений, связанных с образованием объемного заряда (классический закон Чайлда–Ленгмюра).

Динамика незамагниченной ионной компоненты плазмы описывается уравнением движения

(1)
$\frac{{\partial {{{\mathbf{v}}}_{i}}}}{{\partial t}} + ({{{\mathbf{v}}}_{i}} \cdot \nabla ){{{\mathbf{v}}}_{i}} = \frac{e}{{{{m}_{i}}}}{\mathbf{E}},$
и уравнением непрерывности

(2)
$\frac{{\partial {{n}_{i}}}}{{\partial t}} + {\text{div}}\left( {{{n}_{i}}{{{\mathbf{v}}}_{i}}} \right) = 0.$

Здесь ${{{\mathbf{v}}}_{i}}$ и ${{n}_{i}}$ – скорость и концентрация ионов, E – напряженность электрического поля, e – элементарный заряд и ${{m}_{i}}$ – масса иона. В уравне-нии (1) помимо слагаемого с силой Лоренца опущены силы ионного давления, что оправдано в случае холодной плазмы, когда фазовая скорость рассматриваемых возмущений ${{v}_{{ph}}}$ существенно превосходит тепловую скорость ионов ${{v}_{{Ti}}}$, ${{v}_{{ph}}} \gg {{v}_{{Ti}}}$.

Динамику электронной компоненты удобно описывать при помощи уравнения непрерывности, редуцированного на случай достаточно медленных процессов с характерным временем τ, таким что ${{\tau }^{{ - 1}}} \ll {{\omega }_{{Be}}}$ (см. например, [9]):

(3)
$\begin{gathered} \left( {\frac{\partial }{{\partial t}} + {{{\mathbf{V}}}_{E}} \cdot \nabla } \right){{n}_{e}} - 2{{n}_{e}}{{{\mathbf{V}}}_{E}} \cdot \nabla lnB + \\ + \;{\text{div}}\left( {\frac{{{{n}_{e}}}}{{{{\omega }_{{Be}}}}}\left[ {\left( {\frac{\partial }{{\partial t}} + ({{{\mathbf{V}}}_{E}} \cdot \nabla )} \right){{{\mathbf{V}}}_{E}}} \right] \times {\mathbf{b}}} \right) = 0. \\ \end{gathered} $

Здесь ${{n}_{e}}$ – концентрация электронов, ${{\omega }_{{Be}}} = $ $ = eB{\text{/}}{{m}_{e}}c$, ${{m}_{e}}$ – масса электрона, c – скорость света, ${\mathbf{b}} = {\mathbf{B}}{\text{/}}B$ – единичный вектор вдоль направления магнитного поля. Электронная компонента плазмы также считается холодной. В работе [20] показано, что для типичных параметров холловских ускорителей эффекты электронного давления и бесстолкновительной вязкости качественно не влияют на характер неустойчивости азимутальных возмущений внутри ускоряющего канала, и их учет приводит лишь к незначительному (~10%) снижению частоты и волнового числа неустойчивых колебаний. Это обусловлено низкой температурой разряда, ${{T}_{e}}\sim 10$ эВ, и большой величиной электронного тока в плазме СПД. Величина этого тока, выступающего в роли “драйва” рассматриваемой неустойчивости, значительно превышает пороговое значение, необходимое для начала развития неустойчивости, тогда как эффекты бесстолкновительной вязкости существенны лишь непосредственно вблизи самого порога [9, 18].

Далее будем рассматривать только потенциальные возмущения, поэтому для замыкания системы уравнений (1)–(3) используем уравнение Пуассона

(4)
${{\nabla }^{2}}\Phi = 4\pi e({{n}_{e}} - {{n}_{i}}),$
где Φ – электростатический потенциал плазмы (${\mathbf{E}} = - \nabla \Phi $).

Интересуясь линейной стадией развития неустойчивости, представим возмущенные величины в виде стандартного фурье-разложения по времени и азимутальной координате. В частности, для возмущения потенциала

(5)
$\Phi {\text{'}}(x,y,t) = \sum {\Psi (x)exp( - i\omega t + i{{k}_{y}}y)} .$

Здесь ${{k}_{y}}$ – азимутальная компонента волнового вектора. Выражая в (1)–(4) все возмущенные величины через $\Phi '$, приходим к искомому дифференциальному уравнению на собственные частоты глобальных градиентно-дрейфовых мод

(6)
$\begin{gathered} \frac{{{{d}^{2}}\Psi }}{{d{{x}^{2}}}} + {{\kappa }_{n}}\frac{{d\Psi }}{{dx}} - k_{y}^{2}\Psi - \\ - \;\frac{{{{\omega }^{2}}}}{{{{\omega }^{2}} - \omega _{{lh}}^{2}}}\left[ {2{{\kappa }_{B}}\frac{{d\Psi }}{{dx}} - \frac{{{{k}_{y}}\Lambda (x)}}{{\omega - {{\omega }_{E}}}}\Psi } \right] = 0. \\ \end{gathered} $

Здесь ${{\omega }_{{lh}}} = \sqrt {{{\omega }_{{Be}}}{{\omega }_{{Bi}}}} $ – нижнегибридная частота, ${{\kappa }_{{(n,B)}}} = dln(n(x),B(x)){\text{/}}dx$ – параметры, характеризующие аксиальные градиенты равновесной концентрации плазмы $n = {{n}_{i}} = {{n}_{e}}$ (плазма считается квазинейтральной) и магнитного поля, соответственно, ${{\omega }_{E}} = {{k}_{y}}u$, $u = - cE{\text{/}}B$ – скорость равновесного вращения электронов, в рассматриваемом случае равная ${{V}_{E}}$, и

(7)
$\Lambda (x) = {{\omega }_{{Be}}}({{\kappa }_{n}} - 2{{\kappa }_{B}}) + ({{\kappa }_{n}} - 2{{\kappa }_{B}})\left( {\frac{{du}}{{dx}}} \right) + \frac{{{{d}^{2}}u}}{{d{{x}^{2}}}}.$

Уравнение (6) описывает аксиальную структуру глобальных электростатических градиентно-дрейфовых мод, вращающихся в направлении стационарного дрейфа электронов в скрещенных полях в неоднородной плазме с замагниченными электронами и незамагниченными ионами, и учитывает эффекты электронной инерции и неоднородности магнитного поля, а также шир равновесной скорости азимутального вращения электронов. При выводе уравнения (6) явно предполагалось, что плазма является достаточно плотной, ${{\omega }_{{pe}}} \gg {{\omega }_{{Be}}}$, где ${{\omega }_{{pe}}} = \sqrt {4\pi n{{e}^{2}}{\text{/}}{{m}_{e}}} $ – плазменная электронная частота (для типичных максимальных значений $n \simeq 0.5 \times {{10}^{{12}}}$ см–3 и $B \simeq 200$ Гс внутри канала СПД ${{\omega }_{{pe}}} \simeq 50$ ГГц и ${{\omega }_{{Be}}} \simeq 3.5$ ГГц), а расстояние пролета ионов за время $1{\text{/}}\omega $ мало по сравнению с длиной волны возмущений. Последнее предположение позволяет пренебречь влиянием равновесной ионной скорости ${{v}_{i}}$ на неустойчивость. Также не учитывались столкновения электронов с нейтральными атомами. Можно показать (см., например, [21]), что для азимутальных волн столкновения важны при ${{k}_{y}} \gtrsim ({{\kappa }_{n}} - 2{{\kappa }_{B}}){{\omega }_{{Be}}}{\text{/}}{{\nu }_{{ne}}}$, где ${{\nu }_{{ne}}}$ – частота столкновений электронов с нейтралами, т.е. ввиду ${{\omega }_{{Be}}} \gg {{\nu }_{{ne}}} \sim 1$ МГц только для очень коротковолновых возмущений. Отдельно стоит отметить, что второй и третий члены в выражении (7), вообще говоря, пренебрежимо малы в силу условия ${{\tau }^{{ - 1}}} \ll {{\omega }_{{Be}}}$, однако в случае ${{\kappa }_{n}} - 2{{\kappa }_{B}} \approx 0$ слагаемое $\sim {\kern 1pt} {{d}^{2}}u{\text{/}}d{{x}^{2}}$ играет ключевую роль при описании неустойчивостей типа Кельвина–Гельмгольца, возникающих в гидродинамических и плазменных системах с существенно неоднородным профилем вращения среды.

Уравнения, схожие по виду с (6), хорошо известны в теории плазменных и гидродинамических неустойчивостей [15, 22, 23]. Например, в случае ${{\kappa }_{n}} = {{\kappa }_{B}} = 0$ и $\omega \gg {{\omega }_{{lh}}}$ уравнение (6) в точности принимает вид классического уравнения Рэлея [23], описывающего структуру собственных колебаний плоскопараллельных течений несжимаемой жидкости:

(8)
$\frac{{{{d}^{2}}\Psi }}{{d{{x}^{2}}}} - k_{y}^{2}\Psi + \frac{\Psi }{{\omega {\text{/}}{{k}_{y}} - u}}\frac{{{{d}^{2}}u}}{{d{{x}^{2}}}} = 0.$

В случае обычной жидкости величина Ψ определяется как возмущение функции тока $\zeta = {{\zeta }_{0}} + \Psi $, где ${{\zeta }_{0}}$ – функция равновесного тока жидкости: $u = - d{{\zeta }_{0}}{\text{/}}dx$. Плазменный аналог (8) описывает случай желобковых колебаний газа заряженных частиц [22]. Другим примером высокочастотной неустойчивости ($\omega \gg {{\omega }_{{lh}}}$), описываемой уравнением (6), является диокотронная неустойчивость пространственно неоднородного электронного облака.

3. СТАЦИОНАРНЫЕ ПРОФИЛИ ПАРАМЕТРОВ СИСТЕМЫ

Уравнение динамики ионной компоненты (1) в случае одномерного движения ${{{\mathbf{v}}}_{i}} = {{v}_{i}}(x){{{\mathbf{e}}}_{x}}$ дает явную связь между равновесной ионной скоростью и электрическим полем в виде

${{v}_{i}}\frac{{d{{v}_{i}}}}{{dx}} = \frac{e}{{{{m}_{i}}}}E.$

Такое равновесие принято называть баллистическим. Рассматривая стационар уравнения (2) $n(x){{v}_{i}}(x) = {\text{const}}$ ($n = {{n}_{i}} = {{n}_{e}}$), легко получить связь между величиной электрического поля и градиентом концентрации плазмы

(9)
$\frac{1}{n}\frac{{dn}}{{dx}} = - \frac{1}{{{{v}_{i}}}}\frac{{d{{v}_{i}}}}{{dx}} \equiv \frac{{{{V}_{E}}{{\omega }_{{Bi}}}}}{{v_{i}^{2}}}.$

Устойчивость равновесия (9) была исследована в работе [17]. В действительности, экспериментальные данные, полученные на различных ускорителях холловского типа, демонстрируют существенную неоднородность величины ионного потока $n(x){{v}_{i}}(x)$ вдоль ускоряющего канала. Это отклонение связывают с эффектами радиального расхождения ионного потока, процессами ионизации, столкновениями, взаимодействием со стенками ускоряющего канала и т.п. (см. обсуждение в работе [6]). Поэтому далее для описания волн в системе СПД вместо баллистического равновесия мы будем использовать наиболее типичные профили параметров плазмы ($n(x)$, $B(x)$ и $\Phi (x)$), полученные в ходе различных экспериментальных и численных исследований [2426], аппроксимируя их простыми аналитическими зависимостями.

Хорошо известно, что в типовых конструкциях СПД конфигурация магнитного поля подобрана таким образом, чтобы модуль поля был максимален вблизи среза ускоряющего канала. Следуя результатам работы [24], определим профиль магнитного поля в виде стандартного гауссова профиля

(10)
$B(x) = {{B}_{m}}exp\left[ { - {{\nu }_{1}}\mathop {\left( {\frac{x}{d} - 1} \right)}\nolimits^2 } \right].$

Здесь ${{B}_{m}}$ – максимальное значение напряженности магнитного поля, достигаемое в сечении $x = d$, и ${{\nu }_{1}}$ – коэффициент, определяющий пространственный масштаб функции $B(x)$.

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

(11)
$n(x) = {{n}_{m}}exp\left[ { - {{\nu }_{2}}\mathop {\left( {\frac{x}{{{{l}_{1}}}} - 1} \right)}\nolimits^2 } \right].$

Здесь ${{n}_{m}}$ – максимальное значение концентрации в точке $x = {{l}_{1}}$ (${{l}_{1}} < d$) и ${{\nu }_{2}}$ – коэффициент, определяющий пространственный масштаб функции $n(x)$.

Аппроксимация профилей магнитного поля и концентрации плазмы в виде гауссовых функций крайне удобна, т.к. параметры неоднородности ${{\kappa }_{n}}$ и ${{\kappa }_{B}}$ в этом случае имеют простую линейную зависимость от $x$

(12)
${{\kappa }_{B}}(x) = - \frac{{2{{\nu }_{1}}}}{d}\left( {\frac{x}{d} - 1} \right),\quad {{\kappa }_{n}}(x) = - \frac{{2{{\nu }_{2}}}}{{{{l}_{1}}}}\left( {\frac{x}{{{{l}_{1}}}} - 1} \right).$

Для дальнейшего анализа выберем следующие значения коэффициентов: ${{\nu }_{1}} = 2.5$, ${{\nu }_{2}} = 4$, ${{l}_{1}} = $ $ = 0.62d$. При таких значениях выражения (10) и (11) с хорошей точностью описывают сглаженные профили магнитного поля и концентрации плазмы, рассчитанные в работе [25] для двигателя СПД-100. Рассматриваемые профили $B(x)$, $n(x)$, ${{\kappa }_{B}}(x)$ и ${{\kappa }_{n}}(x)$ представлены на рис. 2

Рис. 2.

Зависимости: а) – $B(x){\text{/}}{{B}_{m}}$ и $n(x){\text{/}}{{n}_{m}}$; б) – ${{\kappa }_{B}}(x)$ и ${{\kappa }_{n}}(x)$ от координаты $x{\text{/}}d$.

Электрическое поле зададим посредством электростатического потенциала Φ

(13)
$\Phi (x) = {{\Phi }_{m}}{{a}_{1}}\left( {{\text{arcctg}}\left[ {{{\nu }_{3}}\left( {\frac{x}{{{{l}_{2}}}} - 1} \right)} \right] - {{a}_{2}}} \right),$
где ${{\Phi }_{m}}$ – значение потенциала на аноде; ${{l}_{2}}$ – некоторая точка внутри ускоряющего канала (${{l}_{2}} < d$), в которой электрическое поле достигает своего максимального значения; ${{a}_{1}} = 1{\text{/}}({\text{arcctg}}[ - {{\nu }_{3}}] - {{a}_{2}})$ – нормировочный коэффициент; ${{\nu }_{3}}$ – коэффициент, регулирующий градиент функции $\Phi (x)$. Тогда профиль электрического поля ($E = - d\Phi {\text{/}}dx$) определяется выражением
(14)
${{E(x) = {{E}_{m}}} \mathord{\left/ {\vphantom {{E(x) = {{E}_{m}}} {\left\{ {1 + \nu _{3}^{2}\mathop {\left( {\frac{x}{{{{l}_{2}}}} - 1} \right)}\nolimits^2 } \right\}}}} \right. \kern-0em} {\left\{ {1 + \nu _{3}^{2}\mathop {\left( {\frac{x}{{{{l}_{2}}}} - 1} \right)}\nolimits^2 } \right\}}},$
где ${{E}_{m}} = {{\Phi }_{m}}{{a}_{1}}{{\nu }_{3}}{\text{/}}{{l}_{2}}$. При ${{\nu }_{3}} = 7.5$, ${{l}_{2}} = 0.84d$ и ${{a}_{2}} = 0.1$ выражение (13) описывает сглаженный профиль потенциала, рассчитанный в работе [25]. Профили $\Phi (x)$ и $E(x)$ представлены на рис. 3.

Рис. 3.

Зависимости $\Phi (x){\text{/}}{{\Phi }_{m}}$ и $E(x){\text{/}}{{E}_{m}}$ от координаты $x{\text{/}}d$.

Для анализа неустойчивости далее будем использовать следующие абсолютные значения параметров системы: ${{B}_{m}} = 200$ Гс, ${{n}_{m}} = 0.5 \times {{10}^{{12}}}$ см–3, ${{\Phi }_{m}} = 300$ В. При этом будем считать, что рабочим веществом является ксенон.

4. СПЕКТР ГЛОБАЛЬНЫХ МОД

Решение уравнения (6) для профилей равновесных параметров (10)–(14) и нулевых граничных условий для возмущений потенциала плазмы на электродах $\Psi (0) = \Psi (d) = 0$ находилось численно с применением алгоритма матричной дискретизации задачи на собственные значения, описанного в работе [27]. Используемые граничные условия соответствуют фиксированным значениям потенциала на аноде и виртуальном катоде, т.е. магнитной поверхности, проходящей через катод-нейтрализатор на открытом конце системы. Таким образом, рассматриваются возмущения, локализованные внутри ускоряющего канала двигателя. Для расчета длина и полурадиус ускоряющего канала были выбраны равными $d = 2.5$ см и $R = ({{R}_{1}} + {{R}_{2}}){\text{/}}2 = 4$ см (${{R}_{1}} = 3$ см и ${{R}_{2}} = 5$ см – внутренний и внешний радиусы канала), что примерно соответствует геометрическим размерам канала установки СПД-100. Спектр неустойчивых колебаний (зависимость $\gamma (f)$, где γ – инкремент, f – линейная частота: $\omega = 2\pi f + i\gamma $), полученный для данной геометрии, представляет собой набор из пяти длинноволновых неустойчивых мод с азимутальными волновыми числами $m \equiv {{k}_{y}}R = 1{\kern 1pt} - {\kern 1pt} 5$ в нижнегибридном диапазоне частот $f = - (0.26{\kern 1pt} - {\kern 1pt} 0.67)$ МГц и $\gamma {\text{/}}2\pi = (0.26{\kern 1pt} - {\kern 1pt} $ 0.38) МГц – см. рис. 4а. Моды с $m > 5$ стабилизируются электронной инерцией [9, 19]. Знак минус у частоты означает, что направление распространения волны совпадает с направлением стационарного ${\mathbf{E}} \times {\mathbf{B}}$-вращения электронов.

Рис. 4.

Спектр собственных частот неустойчивых мод с азимутальными волновыми числами $m = {{k}_{y}}R$ (а); аксиальная структура собственных функций (б). Здесь $d = 2.5$ см и $R = 4$ см.

Все моды спектра обладают одинаковой аксиальной структурой собственных функций, представленной на рис. 4б. Как видно, неустойчивые моды локализованы преимущественно в прианодной области ускоряющего канала, что согласуется с результатами локального анализа устойчивости [20]. Отметим, что возмущение электрического поля не равно нулю на аноде, но исчезает к выходу из ускоряющего канала. Это свидетельствует о нечувствительности решения к граничному условию на открытом конце системы.

Спектр неустойчивости имеет характерный “параболический” вид с двумя преобладающими гармониками $m = 3$ и $m = 4$ с частотами ~0.5 МГц, что многократно превышает частоту крупномасштабных структур, наблюдаемых в СПД. При этом все остальные линейные характеристики глобальных мод (область локализации, длина волны, направление вращения) весьма хорошо согласуются с параметрами наблюдаемых структур. Между тем, близость инкрементов и частот неустойчивых гармоник указывает на возможность образования биений (волновых пакетов) с заведомо меньшей частотой огибающих.

Для иллюстрации азимутальной динамики волнового пакета введем функцию

(15)
${{\Phi }_{b}}(y,t) = \sum\limits_{m = 1}^5 exp\left[ { - i(2\pi {{f}_{m}} + i({{\gamma }_{m}} - {{\gamma }_{3}}))t + im\frac{y}{R}} \right],$
входящую в возмущение потенциала (5) в виде

$\Phi {\text{'}}(x,y,t) = \Psi (x){{\Phi }_{b}}(y,t)exp({{\gamma }_{3}}t).$

Здесь ${{\gamma }_{3}}$ инкремент наиболее неустойчивой моды $m = 3$, и мы воспользовались эквивалентностью аксиальной структуры собственных функций для мод с разными m. Подобная нормировка позволяет исключить из рассмотрения экспоненциальное нарастание амплитуды возмущений и исследовать “внутреннюю структуру” биения. Очевидно, что на больших временах возмущение потенциала должно выйти на наиболее неустойчивую моду, т.е. $\mathop {lim}\limits_{t \to \infty } {{\Phi }_{b}}(y,t) = exp\left( { - 2\pi i{{f}_{3}}t + 3iy{\text{/}}R} \right)$, с частотой ~0.5 МГц. Однако на начальной стадии развития неустойчивости при относительно небольших амплитудах возмущений, когда, собственно, и применима линейная теория, совместная раскачка различных мод выглядит как крупномасштабная осцилляция.

Эволюция величины ${{\Phi }_{b}}$ в азимутальном сечении $y = {\text{const}}$ изображена на рис. 5. Видно, что возмущения равной начальной амплитуды релаксируют к биениям, образованным модами $m = 3$ и $m = 4$, с огибающими $\sim {\kern 1pt} cos(2\pi {{f}_{b}}t)$, где ${{f}_{b}} = $ $ = ({{f}_{3}} - {{f}_{4}}){\text{/}}2$. Частота огибающих ${{f}_{b}} \simeq 52$ кГц, что соответствует наблюдаемому диапазону частот крупномасштабных осцилляций.

Рис. 5.

Зависимость функции ${{\Phi }_{b}}$ от времени t в сечении $y = {\text{const}}$. Амплитуда функции нормирована на единицу. Штрихованные кривые описывают огибающие волнового пакета $\Phi _{b}^{*}$.

Для иллюстрации полученных возмущений вернемся к реальной (цилиндрической) геометрии СПД и рассмотрим в этой геометрии динамику волнового пакета в аксиальном сечении $x = {\text{const}}$, т.е. на срезе двигателя: именно в этой плоскости в экспериментах удается визуально наблюдать крупномасштабные структуры при помощи высокоскоростной видеофиксации11 – см., например, [28]. Выберем параболическую зависимость возмущения потенциала $\Phi {\text{'}}$ от радиуса z: $\Phi '(x,y,z,t) = [z_{0}^{2} - {{(z - R)}^{2}}]\Phi {\text{'}}(x,y,t)$, обеспечивающую зануление возмущения на внешнем ${{R}_{2}}$ и внутреннем ${{R}_{1}}$ радиусах ускоряющего канала при ${{z}_{0}} = 1$ см. Следует отметить, что для полученных нами решений радиальная (от z) зависимость потенциала может быть произвольной. Мы не решаем точную трехмерную задачу, поскольку учет продольного (вдоль направления магнитного поля) движения электронов не вносит существенного вклада в дисперсионное уравнение [15, 17].

Пространственно-временная эволюция волнового пакета в сечении $x = {\text{const}}$ представлена на рис. 6 в виде линий уровня огибающей функции $\Phi _{b}^{*}$, построенных в разные моменты времени (координаты y* и z* – декартовы координаты прямоугольной системы в плоскости $x = {\text{const}}$ – приведены на рис. 6 для сопоставления с реальными размерами СПД). В начальный момент времени $t = 0$ волновой пакет представляет из себя локализованную структуру $m = 1$ с прослеживающимся вкладом от более высоких гармоник – см. рис. 6а. Релаксируя к биениям на двух гармониках, наблюдаемая структура растягивается и смещается в направлении ${\mathbf{E}} \times {\mathbf{B}}$-дрейфа электронов – см. рис. 6б, в. Групповая скорость волнового пакета, рассчитанная по движению максимума огибающей, равна ${{V}_{\theta }} \simeq 5.5 \times {{10}^{5}}$ см/с, что по порядку величины составляет $\sim {\kern 1pt} 0.02{{V}_{E}}$ (${{V}_{E}}$ рассчитана в точке $x = 0.15d$, соответствующей максимальному значению амплитуды моды вдоль оси ускоряющего канала – см. рис. 4б).

Рис. 6.

Линии уровня огибающей волнового пакета $\Phi _{b}^{*}$ в сечении $x = {\text{const}}$ для различных моментов времени t: a) – 0.0, б) – 3.0 мкс, в) – 6.0 мкс. На рисунке в плоскости $x = {\text{const}}$ введена декартова система координат $\{ y{\text{*}},z{\text{*}}\} $.

5. ВЛИЯНИЕ МАСШТАБИРОВАНИЯ ГЕОМЕТРИИ СПД НА ХАРАКТЕРИСТИКИ НЕУСТОЙЧИВОСТИ

Рассматриваемая картина биений формируется лишь в случае, когда спектр собственных колебаний обладает хотя бы двумя модами с близкими по величине инкрементами. Частота огибающей волнового пакета и его групповая скорость при этом напрямую зависят от “удаленности” мод (т.е. от разности частот колебаний с примерно равными γ). Форма спектра (на плоскости $f - \gamma $), в свою очередь, оказывается весьма чувствительной к изменению геометрических размеров СПД.

Для начала исследуем изменение формы спектра колебаний в зависимости от радиуса ускоряющего канала при фиксированной длине канала $d = 2.5$ см. Спектры неустойчивых мод для разных значений R представлены на рис. 7. Видно, что с увеличением радиуса число неустойчивых гармоник растет: $m = 1{\kern 1pt} - {\kern 1pt} 5$ при $R = 4$ см, $m = 1{\kern 1pt} - {\kern 1pt} 8$ при $R = 6$ см и $m = 1{\kern 1pt} - {\kern 1pt} 11$ при $R = 8$ см – а форма спектра остается неизменной. Также увеличивается и азимутальный номер моды, обладающей наибольшим инкрементом: $m = 3$ при $R = 4$ см, $m = 5$ при $R = 6$ см и $m = 7$ при $R = 8$ см, – при этом сами значения инкрементов наиболее неустойчивых мод практически не изменяются по величине: ${{\gamma }_{3}}{\text{/}}2\pi \simeq 0.38$ МГц, ${{\gamma }_{5}}{\text{/}}2\pi \simeq 0.39$ МГц и ${{\gamma }_{7}}{\text{/}}2\pi \simeq 0.39$ МГц. Самым существенным следствием увеличения радиуса является сближение частот неустойчивых мод с близкими γ, что приводит к уменьшению частоты биений (n-кратное увеличение радиуса приводит к n-кратному уплотнению спектра, т.е. ${{f}_{b}}\sim 1{\text{/}}R$). Так, например, при радиусе $R = 8$ см частота огибающей волнового пакета, образованного модами m = 6–8, равна ${{f}_{b}} \simeq 23$ кГц. Тенденция уменьшения частоты спиц с ростом радиуса ускоряющего канала явно наблюдается в экспериментах [29].

Рис. 7.

Спектры собственных частот неустойчивых колебаний при $R = 4$ см, $R = 6$ см и $R = 8$ см; $d = $ = 2.5 см. Стрелка указывает направление увеличения азимутального волнового числа $m$. Аксиальная структура собственных функций – рис. 4б.

Изменение длины ускоряющего канала d при фиксированном радиусе ($R = 4$ см) оказывает более сложное влияние на спектр колебаний – см. рис. 8а – поскольку сопровождается изменением градиентов равновесных параметров плазмы и магнитного поля. Как и при увеличении радиуса, с ростом длины канала растет число неустойчивых гармоник: $m = 1{\kern 1pt} - {\kern 1pt} 5$ при $d = 2.5$ см, $m = 1{\kern 1pt} - {\kern 1pt} 10$ при $d = 3.5$ см и $m = 1 - 12$ при $d = 4.5$ см, причем в последнем случае в спектре неустойчивости появляются две дополнительные моды $m = 1$, 2 (на рисунке отмечены как ${{m}^{ + }}$) c более мелкомасштабной аксиальной структурой собственных функций – см. рис. 8б. Увеличение d также приводит к росту азимутального номера наиболее неустойчивой моды: $m = 3$ при $d = 2.5$ см, $m = 9$ при $d = $ = 3.5 см и $m = 11$ при $d = 4.5$ см, что, однако, сопровождается и изменением инкремента: ${{\gamma }_{3}}{\text{/}}2\pi \simeq $ $ \simeq \;\,0.38$ МГц, ${{\gamma }_{9}}{\text{/}}2\pi \;\, \simeq \;\,0.51$ МГц и ${{\gamma }_{{11}}}{\text{/}}2\pi \;\, \simeq $ $ \simeq 0.49$ МГц. Кроме того, с ростом d увеличивается интервал между частотами и инкрементами наиболее неустойчивых гармоник. Это усложняет процесс образования биений и повышает их частоту.

Рис. 8.

Спектры собственных частот неустойчивых колебаний при $d = 2.5$ см, $d = 3.5$ см и $d = 4.5$ см; $R = 4$ см (а); аксиальная структура собственных функций мод ${{m}^{ + }} = 1$, 2 (б). Аксиальная структура собственных функций мод mрис. 4б. На рисунке (а) стрелка указывает направление увеличения волнового числа m.

Дополнительные расчеты показывают, что для формирования спектра, наиболее неустойчивые моды в котором обладают соразмерными инкрементами (разность которых не превышает 10%) и близкими частотами ${{f}_{{m + 1}}} - {{f}_{m}} \lesssim 100$ кГц, соотношение между радиусом и длиной ускоряющего канала должно превосходить некоторое критическое значение $R{\text{/}}d \geqslant 1.7$.

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

В работе исследована устойчивость глобальных градиентно-дрейфовых волн в стационарном плазменном двигателе. Расчеты проведены в рамках идеальной двужидкостной гидродинамической модели холловской плазмы; использованы аксиальные профили параметров в СПД-100, аппроксимированные аналитическими выражениями (10)–(14). Показано, что собственные неустойчивые колебания в двигателе характеризуются крупномасштабными азимутальной $m = 1{\kern 1pt} - {\kern 1pt} 5$ и аксиальной (см. рис. 4б и 8б) структурами и локализованы в прианодной области ускоряющего канала. Ключевую роль в формировании спектра играет инерция электронов, стабилизирующая коротковолновые возмущения [9, 19], а пространственная локализация мод объясняется выполнением локального критерия неустойчивости ${{V}_{E}}({{\kappa }_{n}} - 2{{\kappa }_{B}}) < 0$ лишь в прианодной части ускорителя [20].

Показано, что для характерного “параболического” спектра (см. рис. 4a) глобальные моды на линейной стадии развития неустойчивости образуют волновые пакеты $m = 1$ с характерной частотой огибающей ~50 кГц и групповой скоростью $\sim {\kern 1pt} 0.02{{V}_{E}}$. Варьирование параметров СПД показало, что биения легче формируются в ускоряющих каналах c малой длиной d и большим радиусом R, причем их частота $\sim {\kern 1pt} 1{\text{/}}R$. Соответствие между основными характеристиками волновых пакетов и параметрами крупномасштабных азимутальных структур, наблюдаемых в двигателях семейства СПД, как-то: прианодная локализация, частоты ~10 кГц, низкие m, медленное вращение в направлении $[{\mathbf{E}} \times {\mathbf{B}}]$-дрейфа, обратная пропорциональность частоты радиусу ускоряющего канала – позволяет сделать вывод о том, что именно градиентно-дрейфовая неустойчивость может служить вероятным механизмом возникновения макроскопических структур в виде спиц в СПД.

Работа частично поддержана грантом Российского фонда фундаментальных исследований № 16-02-00640, при финансовой поддержке которого получены результаты, приведенные в разд. 5, и грантом Российского научного фонда № 17-12-01470. Публикация подготовлена при поддержке Программы РУДН “5-100”.

Список литературы

  1. Choueiri E.Y. // Phys. Plasmas. 2001. V. 8. P. 1411.

  2. Плазменные ускорители / Под ред. Арцимови-ча Л.А., Гришина С.Д., Гродзовского Г.Л., Леско-ва Л.В., Морозова А.И., Дороднова А.М., Падал-ки В.Г., Пергамента М.И. М.: Машиностроение, 1973.

  3. Ellison C.L., Raitses Y., Fisch N.J. // Phys. Plasmas. 2012. V. 19. P. 013503.

  4. Janes G.S., Lowder R.S. // Phys. Fluids. 1966. V. 9. P. 1115.

  5. Есипчук Ю.В., Морозов А.И., Тилинин Г.Н., Трофимов А.В. // ЖТФ. 1973. Т. 43. № 7. С. 1466.

  6. Frias W., Smolyakov A.I., Kaganovich I.D., Raitses Y. // Phys. Plasmas. 2012. V. 19. P. 072112.

  7. Lakhin V.P., Ilgisonis V.I., Smolyakov A.I., Soroki-na E.A. // Phys. Plasmas. 2016. V. 23. P. 102304.

  8. Smolyakov A.I., Chapurin O., Frias W., Koshkarov O., Romadanov I., Tang T., Umansky M., Raitses Y., Kaganovich I.D., Lakhin V.P. // Plasma Phys. Control. Fusion. 2017. V. 20. P. 052108.

  9. Lakhin V.P., Ilgisonis V.I., Smolyakov A.I., Soroki-na E.A., Marusov N.A. // Phys. Plasmas. 2018. V. 25. P. 012106.

  10. Simon A. // Phys. Fluids. 1963. V. 6. P. 382.

  11. Hoh F.C. // Phys. Fluids. 1963. V. 6. P. 1184.

  12. Hargreaves J.K. The Solar-Terrestrial Environment, Cambridge Univ. Press, 1992.

  13. Morozov A.I., Savelyev V.V. // Rev. Plasma Phys. / Ed. by B.B. Kadomtsev and V.D. Shafranov (Consultants Bureau, New York, 2000), V. 21. P. 203.

  14. Boeuf J.P. // Front. Phys. 2014. V. 2. P. 74.

  15. Михайловский А.Б. Теория плазменных неустойчивостей. Т. 2. Неустойчивости неоднородной плазмы. М.: Атомиздат, 1971.

  16. Морозов А.И., Есипчук Ю.В., Капулкин А.М., Невровский В.А., Смирнов В.А. // ЖТФ. 1972. Т. 42. № 3. С. 612.

  17. Есипчук Ю.В., Тилинин Г.Н. // ЖТФ. 1976. Т. 46. № 4. С. 718.

  18. Lakhin V.P., Ilgisonis V.I., Smolyakov A.I., Sorokina E.A., Marusov N.A. // Phys. Plasmas. 2018. V. 25. P. 012107.

  19. Гордеев А.В., Гречиха А.В. // Физика плазмы. 1992. Т. 18. С. 3.

  20. Marusov N.A., Sorokina E.A., Lakhin V.P., Ilgisonis V.I., Smolyakov A.I., 2018 Plasma Sources Sci. Technol. in press https://doi.org/10.1088/1361-6595/aae23d.

  21. Escobar D., Ahedo E. In 34th Int. Electric Propulsion Conference, IEPC-2015-371, 2015.

  22. Тимофеев А.В. Резонансные явления в колебаниях плазмы. М.: Физматлит, 2009.

  23. Ландау Л.Д., Лифшиц Е.М. Гидродинамика. М.: Наука, 1986.

  24. Boeuf J., Garrigues L. // J. Appl. Phys. 1998. V. 84. P. 3541.

  25. Hofer R., Mikellides I., Katz I., Goebel D. In 43rd AIAA/ASME/SAE/ASEE Joint Propulsion Conference & Exhibit, AIAA 2007-5267, 2007.

  26. Kronhaus I., Kapulkin A., Balabanov V., Rubanovich M., Guelman M., Natan B. // J. Phys. D: Appl. Phys. 2012. V. 45. P. 175203.

  27. Сорокина Е.А. // Физика плазмы. 2009. Т. 35. С. 472.

  28. Parker J.B., Raitses Y., Fisch N.J. // Appl. Phys. Lett. 2010. V. 97. P. 091501.

  29. McDonald M.S., Gallimore A.D. In 32nd Int. Electric Propulsion Conf., IEPC-2011-242, 2011.

Дополнительные материалы отсутствуют.