CC BY
25
ТЕОРЕТИЧЕСКАЯ И МАТЕМАТИЧЕСКАЯ ФИЗИКА
Электронная вайбелевская неустойчивость плазмы с температурной
анизотропией
Л. В. Бородачёв0, Д. О. Коломиец6
Московский государственный университет имени М. В. Ломоносова, физический факультет, кафедра математики. Россия, 119991, Москва, Ленинские горы, д. 1, стр. 2.
E-mail: а [email protected], ь [email protected]
Статья поступила 12.10.2009, подписана в печать 19.11.2009
Методом дискретного моделирования самосогласованной безызлучательной плазмы получена общая картина вайбелевской неустойчивости, связанной с анизотропным распределением температуры в электронной компоненте. Исследованы зависимости от величины исходной анизотропии среды основных параметров неустойчивости (характерного периода и максимальной энергии магнитного поля).
Ключевые слова: вайбелевская неустойчивость, метод макрочастиц, модель Власова-Дарвина.
УДК: 517.9:533.9. PACS: 52.35.Qz.
Введение
Как известно, анизотропное распределение электронов по скоростям в однородной разреженной плазме может вызывать развитие так называемой вайбелевской неустойчивости (ВН) [1], суть и механизм которой легко понять из следующей одномерной по пространству и трехмерной по скоростям физической картины.
Пусть в однородной бесстолкновительной плазме с положительным (ионным) фоном и плотностью электронов п0 задано распределение скоростей частиц вида
щ ехр
т =
ТГ3/'2ихиуиг
(1)
где иг — тепловая скорость электронов, движущихся вдоль оси г, много выше, чем тепловые скорости по другим осям (движение вдоль оси у для простоты можно вообще не рассматривать). Обычно указанную анизотропию среды количественно выражают с помощью показателя (коэффициента) Ло = (и2г/и2х - 1).
В этом случае плазма, очевидно, обладает избытком быстрых частиц по оси 2. Однако в силу симметрии скоростного распределения электронов полный ток системы равен нулю.
Теперь допустим, что в области случайно (например, из шумового фона) возникло ненулевое магнитное поле Во = еуВуъ эт^д-л'), поперечное к акцентированной составляющей скорости частиц. Тогда сила Лоренца р1 = х В0/с искривит траекторию частиц, двигавшихся вдоль 2, так, как это показано на рис. 1, что приведет в конечном итоге к образованию пространственно разделенных токовых слоев. Локализация этих слоев будет совпадать с областями нулевых значений магнитного поля, вызывающего отклонения частиц. Следовательно, в нашем случае можно ожидать образования двух токовых слоев (с разными знаками) на каждую длину волны возмущения, усиливающих исходное магнитное поле, в свою очередь увеличивающего плотность токовых слоев. Развитие неустойчивости
Рис. 1. Качественное представление линейной фазы развития вайбелевской неустойчивости. Образование токовых слоев J вокруг нулей Ву (на рисунке изображена в виде поверхности)
(рост амплитуды, а вместе с ней и энергии поля) продолжается до тех пор, пока частицы в основной своей массе не становятся существенно замагниченными.
Характерной чертой описанного процесса является наличие существенно нелинейного режима, при котором резкий рост амплитуды поля, обусловленный образованием и развитием токовых слоев (пучков в многомерном случае) в начальном периоде, сменяется их (пучков) динамической трансформацией, стабилизирующейся вместе с энергией магнитного поля на стадии насыщения.
Таким образом, уже из этого короткого экскурса в механизмы генерации и этапы развития вайбелевской неустойчивости становится очевидно, что ее детальное изучение (особенно нелинейной стадии) предполагает активное применение численного анализа, и в частности его наиболее эффективного инструментария — компьютерного эксперимента на основе кинетического моделирования по методу частиц [2]. В этой связи можно отметить ряд примечательных работ [3-5], позволивших получить достаточно полную картину клас-
сической (электронной) ВН. Вместе с тем влияние ряда факторов, и прежде всего величины исходной анизотропии среды, на развитие неустойчивости представляется не достаточно изученным, что определило основное направление настоящих исследований.
Нужно подчеркнуть, что подавляющее большинство работ по ВН, связанных с дискретным моделированием (особенно многомерным), выполнено в рамках полного электромагнитного описания бесстолкновительной плазмы. Притом что сравнительно низкочастотный характер ВН позволяет использовать для ее изучения физически адекватное безызлучательное (магнитоиндук-ционное или дарвинское) представление, существенно более экономичное в вычислительном плане. Этот аспект магнитоиндукционного моделирования по методу частиц становится особенно актуальным при исследовании неустойчивости на временах, когда активную роль начинают играть ионы, т. е. в масштабах, на порядки превосходящих характерные времена классической (электронной) ВН, где ионная компонента в силу инертности обычно рассматривается в виде нейтрализующего положительного фона.
Цели настоящих численных экспериментов состояли в получении классической картины многомерной ВН в рамках безызлучательного приближения разреженной плазмы и изучении зависимости основных параметров неустойчивости (характерного периода развития и максимальной энергии магнитного поля) от исходной анизотропии среды.
1.Основные уравнения модели
В настоящей работе исследования вайбелевской неустойчивости были выполнены в рамках бесстолкно-вительного кинетического описания плазмы [61:
dfs „г Я („ . v х В
ot т
V„/s = О,
где q, т — соответственно заряд и масса частицы, s = /, е — сорт частиц, f{r,v,t) — одночастичная функция распределения. Самосогласованные поля при этом определяются приближением Дарвина [7]:
„ _ 4тг _ 1 дЕр VxB=
с с at
VB = О,
р 1 дв
с at
VEP = 4 -p. [E = EP+EV, V x £p = 0, VEv = 0,
где Ep и Ev — соответственно потенциальная (продольная) и соленоидальная (поперечная) компоненты электрического поля, а плотности заряда и тока выражаются через моменты функции распределения:
р(г) = £
s
qsfs(r, v) dvxdvydvz, qsvfs(r, v) dvxdvydvz,
Как легко видеть, предложенный полевой формализм отличается от максвелловского лишь опущенной поперечной составляющей тока смещения, что физически равносильно пренебрежению излучением и переходу к мгновенному дальнодействию [8]. При этом часть индукционных эффектов (связанных с законом Фарадея) и тождественное удовлетворение уравнения непрерывности в виде
4тгJ„
Щ,
dt
= 0
в системе сохраняются.
Нужно отметить, что в случае любой корректной разностной аппроксимации характерная особенность дарвинского приближения — совмещение незапаздыва-ющего характера и гиперболической формы представления — требует его предварительной переформулировки, исключающей производные по времени из уравнений поля [9-11].
Используемый в настоящих компьютерных экспериментах безызлучательный алгоритм представляет собой численную интерпретацию по методу крупных частиц приведенной выше интегро-дифференциальной модели Власова-Дарвина с дискретной функцией распреде-
ления
N„
fD(r,v,t) = Y|R(r,rp(t))5(v^Vp(t)),
Р=1
где Ыр — общее число макрочастиц, а формфактор отдельной р-й макрочастицы (с координатой гр и линейными размерами Нх, ку) имеет вид
R{r^rp) = {
1
hxhy'
0,
hx 2
hr
^ \У-Ур\ ^ if5
hy 2 hu
*p| > -J и\У^Ур\ > у-
При его построении была использована эллиптическая редакция на основе подхода, предложенного в работе [12], сводящего исходное полевое описание к следующей системе краевых задач в заданной области О:
' А(р = —\т:р.
(л, У(р)\т = (п,Еь), Ер =
4тг
АВ= —-V х I,
Щ Ш = 0' B\m = Bb,
' АЕи = Gu + (рЕи)и, VEv = 0,
Ev | qq Е^
"Plan'
где
G = t + pEp + Cx В,
s s s
4mu
ГПсС3
Js,
8 ВМУ. Физика. Астрономия. № 2
8 = I, е — сорт частиц, — их масса, и — плотности заряда и тока соответственно, Д, — дивергенция тензора переноса тока, Еь — граничное поле.
Вычисления полей проводились с использованием методики, описанной в работе [13], на равномерной (с шагами Их, Иу) прямоугольной сетке 0(г,-,г/), краевые узлы которой совпадают с границей рассматриваемой области <90.
Разностное интегрирование динамических уравнений частиц (являющихся уравнениями характеристик кинетического уравнения Власова с дискретной функцией распределения)
= урхВ{гр,{)\
аА ' аА тр \ с )
р=1,2.....Мр, (2)
осуществлялись по оптимизированной неявной схеме [14, 15] с временным шагом т, выбираемым из физических соображений и условий численной корректности расчета.
2. Численные эксперименты
Рассмотрим основные результаты численного моделирования вайбелевской неустойчивости в 2.5-мерной (х, у,их,'Оу,иг) постановке, где акцентированная компонента тепловой скорости частиц выбрана перпендикулярно плоскости моделирования х - у. Попутно заметим, что развитие неустойчивости теперь будет сопровождаться образованием уже токовых жгутов, а не слоев, как в описанном выше 1.5-мерном случае.
Итак, зададим равномерное пространственное распределение электронов и однозарядных ионов, причем для последних реализуем фоновый режим. В выражении (1) положим иг = 0.1с, их = и,у = 0.0316с. При этих значениях тепловых скоростей показатель исходной анизотропии среды Л0 = 9. Отметим, что вариация базового значения параметра Ло в серии экспериментов осуществлялась изменением величины (л\ у)-компонент тепловой скорости при фиксированной величине иг.
Выбор линейного размера области в нашей постановке был продиктован тем соображением, чтобы на этой длине могла реализоваться конкуренция возбуждаемых мод магнитного поля, имеющих разные волновые числа и инкременты возрастания [3]. В частности, выбор Ьх = 2Ьс/шре, Ьу = 25с/шре позволил детально проследить развитие и структурную перестройку на нелинейной стадии исходной многомодовой системы возмущений магнитного поля.
Наконец, следуя уже упомянутой работе [3], для полей и частиц были приняты периодические граничные условия по х и у.
Динамика нарастания неустойчивости для различных исходных значений показателя анизотропии хорошо видна на графиках зависимости средней по пространству плотности энергии магнитного поля от времени, приведенных на рис. 2. Начальное значение средней плотности энергии близко к нулю, но уже в момент времени ¡шре = 100, соответствующий появлению областей локализации тока, наблюдается ее заметное увеличение, тем большее, чем больше исходная анизотропия среды. Как было отмечено выше, линейная стадия неустойчивости заканчивается,
х10~4
Рас. 2. Зависимость величины средней плотности энергии магнитного поля от времени для различных значений начальной анизотропии (Ло = 1-10) при постоянном значении акцентированной тепловой скорости. Размер счетной области: Lx = Ly = '25с/шре; иг = 0Лс\ размер пространственной сетки 128 х 128; 500 частиц каждого сорта на ячейку; шаг по времени г = O.'25/Шре
когда в системе образуется ярко выраженная токовая структура, а сами частицы становятся существенно замагниченными (момент, когда t ■ шре = 150). Далее процесс носит нелинейный характер и по сути соответствует стадии насыщения, характеризуемой относительной стабилизацией величины энергии магнитного поля.
Для уточнения полученной картины весьма полезным представляется рассмотрение временной эволюции средней плотности энергии частиц, а также магнитного и электрического полей, приведенной на рис. 3 для базового (Ло = 9 и иг = 0.1с) варианта.
Здесь прежде всего отметим, что, как и следует из теории электромагнитных неустойчивостей [16], исходная анизотропия системы генерирует в ней в основном именно магнитное поле, энергия которого быстро растет (с согласованным падением начальной кинетической
хЮ 4
Рас. 3. Динамика средней плотности энергии магнитного поля ги<в, электрического поля ги<£, кинетической энергии ги% и полной энергии системы (пунктир) для А = 9, иг=0Лс
энергии частиц) на активной (линейной) стадии вайбе-левской неустойчивости. Отметим при этом, что, как легко видеть из графика, доля переходящей кинетической энергии совпадает с оценкой (10%), приведенной в работе [3], где, в отличие от настоящих численных экспериментов, на базе дарвинского (безызлучательно-го) приближения самосогласованных полей проведено дискретное моделирование ВН в рамках полного электромагнитного описания плазмы.
Далее (см. нижнюю часть графика), в соответствии с линейной теорией неустойчивости энергия электрического поля остается практически постоянной и при этом крайне малой величины (по сути на уровне шума). Последнее объясняется тем фактом, что чистые электромагнитные неустойчивости (и вайбелевская, в частности) практически не влияют на исходно однородную плотность заряда, оставляя ее возмущение на уровне ленгмюровских колебаний.
Как было указано выше, одна из задач настоящих экспериментов состояла в выяснении зависимости основных характеристик неустойчивости: ее инкремента (характерного периода) и максимальной плотности энергии магнитного поля от показателя исходной анизотропии среды.
В этой связи вернемся к рис. 2, из которого усматривается ожидаемый эффект: повышение максимальной плотности энергии поля и инкремента неустойчивости при увеличении исходного коэффициента анизотропии.
Менее тривиален, однако, характер этой зависимости, представленной на рис. 4. Здесь, как мы видим, резкий рост кривой плотности энергии на отрезке Ло < 10, сменяется столь же резким выполаживанием графика при Ло > 10. Более того, значение Ло в области 25 можно считать своеобразным порогом неустойчивости, поскольку дальнейшее увеличение степени исходной анизотропии среды (при фиксированной величине акцентированной компоненты тепловой скорости частиц) практически не влияет на базовый параметр ВН — максимальную плотность энергии магнитного поля.
Отметим также, что, если характерное время (7)) неустойчивости определить как период достижения максимальной плотности энергии магнитного поля, то из того же рис. 2 легко оценить и характер зависимости Т[ от Л0. При этом наблюдается хорошее
хЮ"4 4.0
| 3.0
¡2.0
1.0
О 10 20 30 40 50
Ай = {и2, / м^-1)
Рис. 4. Зависимость максимального значения средней плотности энергии магнитного поля от величины начальной анизотропии при фиксированном значении иг = 0.1с (те же расчетные параметры, что и для рис. 2)
совпадение численного значения безразмерного инкремента неустойчивости (0.034) с прогнозом линейной теории (0.037) в рассмотренном выше случае Л0 = 9 и иг = 0.1с.
Далее, если ввести понятие энергии анизотропии
шл = ^гп-е ^и2г - ^ (и2 + и2у + и2) ^ , (3)
то ее можно рассматривать как своего рода свободную энергию системы, идущую на развитие вайбелевской неустойчивости. При том, как следует из расчетов, переходящая доля этой энергии зависит от степени исходной анизотропии среды и варьируется от примерно 8% (для Л0 = 9) до примерно 3.5% (для Л0 = 2).
Эти данные уточняют результаты цитируемой выше работы [3], где получена лишь верхняя оценка переходящей в энергию неустойчивости доли начальной кинетической энергии электронов плазмы.
Заключение
Проведенные компьютерные эксперименты позволили получить многомерную картину классической (электронной) вайбелевской неустойчивости в рамках низкочастотного (дарвинского) приближения самосогласованной плазмы. При этом удалось дополнительно выявить характер зависимостей времени развития ВН и максимального значения плотности энергии магнитного поля от начальной величины анизотропии среды. Установить характер связи доли «свободной» энергии системы, идущей на развитие неустойчивости, и исходного показателя анизотропии; оценить его пороговое значение.
Попутно отметим, что апробированный в настоящих численных исследованиях прикладной пакет численного моделирования разреженной плазмы без излучения DarWin является программной реализацией в технологии параллельных вычислений представленного выше многомерного дискретного дарвинского алгоритма. Анализ полученных результатов как в содержательном (физическом), так и в техническом (вычислительном) аспектах позволяет уверенно говорить о возможности его эффективного использования при изучении на кинетическом уровне широкого круга явлений плазмофи-зики, имеющих сравнительно низкочастотный и магни-тоактивный характер.
Список литературы
1. Weibel E.S. // Phys. Rev. Lett. 1959. 2. P. 83.
2. Хокни P., Иствуд Дж. Численное моделирование методом частиц. М., 1987.
3. Morse R.L., Nielsen C.W. // Phys. Fluids. 1971. 14, N 11.
P. 830.
4. Davidson R.S., Hammer DA. // Phys. Fluids. 1972. 15,
N 2. P. 317.
5. Lemons D.S., Winske D., Gary S.P. // J. Plasma Phys.
1979. 21, pt 2. P. 287.
6. Власов A.A. Теория многих частиц. M.; Jl., 1950.
7. Darmin С.G. 11 Phil. Mag. 1920. 39. P. 537.
8. Kaufman A.N., Rostler P.S. // Phys. Fluids. 1971. 14.
P. 446.
9 ВМУ. Физика. Астрономия. 2
9. Нильсон К., Льюис Г. Модели укрупненных частиц в безызлучательном пределе // Управляемый термоядерный синтез. М., 1980. С. 395.
10. Бородачёв Л.В. // Вестн. Моск. ун-та. 1993. № 3. С. 87.
11. Бородачёв Л.В. // Матем. моделирование. 2005. 17, № 9. С. 53.
12. Бородачёв Л.В., Мингалёв И.В., Мшгалёв О.В. // Матем. моделирование. 2006. 18, № 11. С. 117.
13. Бородачёв Л.В., Коломиец Д.О., Литвинюк В.В. // Вестн. Моск. ун-та. 2006. № 6. С. 14.
14. Бородачёв Л.В. // ЖВММФ. 1991. 30, № 6. С. 934.
15. Бородачёв Л.В., Коломиец Д.О. // Тез. докл. междунар. конф. «Тихонов и современная математика». М., 2006. С. 37.
16. Михайловский А.Б. Теория плазменных неустойчиво-стей. М., 1975.
Electron Weibel instability of plasma with thermal anisotropy L.V. Borodachev , D.O. Kolomiets"
Department of Mathematics, Faculty of Physics, M. V. Lomonosov Moscow State University, Moscow 119991, Russia.
E-mail: a [email protected], b [email protected].
The general picture of the thermal-anisotropy-driven Weibel instability obtained using a discrete numerical simulation of self-consistent radiationless plasma. The dependencies of the key instability parameters (the characteristic time and the maximum magnetic field energy) of the initial anisotropy of the medium are investigated.
Keywords: Weibel instability, particle-in-cell, PIC, Vlasov-Darwin model. PACS: 52.35.Qz. Received 12 October 2009.
English version: Moscow University Physics Bulletin 2(2010).
Сведения об авторах
1. Бородачёв Леонид Васильевич — канд. физ.-мат. наук, доцент; тел.: (495) 939-39-47, e-mail: [email protected].
2. Коломиец Дмитрий Олегович — аспирант; e-mail: [email protected].
