2007
НАУЧНЫЙ ВЕСТНИК МГТУ ГА серия Математика и физика
№ 114
УДК 533.9:51-7
ПРИМЕНЕНИЕ МЕТОДА КРУПНЫХ ЧАСТИЦ ДЛЯ АНАЛИЗА ПОВЕДЕНИЯ ДВУХКОМПОНЕНТНОЙ ПЛАЗМЫ С УЧЕТОМ СТОЛКНОВЕНИЙ МЕЖДУ ЗАРЯЖЕННЫМИ ЧАСТИЦАМИ
А.В. ПАНТЕЛЕЕВ, И.А. КУДРЯВЦЕВА
Рассмотрена зондовая задача применительно к двухкомпонентной плазме с учетом столкновений между заряженными частицами в случае плоского зонда. Сформирована математическая модель, представляющая собой самосогласованную систему уравнений, образованную уравнением Фоккера-Планка и уравнением Пуассона. Разработан алгоритм и программное обеспечение для решения поставленной задачи. Получены результаты численного моделирования поведения двухкомпонентной плазмы с учетом столкновений между заряженными частицами.
1. Постановка задачи
Рассматривается следующая физическая постановка зондовой задачи [1]. В невозмущенную бесконечно-протяженную плазму, состоящую из электронов и однозарядных ионов, внесена большая плоскость, заряженная до потенциала (рр. Плоскость расположена поперек потока
плазмы. Плоскость является идеально поглощающей для электронов, ионы при ударе о стенку нейтрализуются. Концентрации ионов п¥ и электронов пе¥, причем п¥ = пе¥ = п¥, а также температуры данных частиц Т¥ Те¥ в невозмущенной плазме, заданы. Каждая частица движется в
макроскопическом электрическом поле, магнитное поле отсутствует. Начальные распределения для всех типов частиц предполагаются максвелловскими.
Требуется найти с учетом столкновений между заряженными частицами самосогласованное электрическое поле Е (г, *) и функции распределения электронов /е (г, V, *) и однозарядных
ионов / (г, V, *), а также их моменты (плотности токов ионов и электронов £(г, * ) = q • | / (г, V, ёу , ]е (г, * ) = е • | /е (г, V, ёу , где q - заряд иона; е - заряд электрона; концентрации ионов и электронов п (г, * ) = | / (г, V, *) ёу , пе (г, * ) = | /е (г, V, *) ёу ). Поведение частиц во времени * характеризуется радиус-вектором г и вектором скорости V .
Математическая модель, соответствующая данной физической постановке, имеет вид:
д/а(r,V,і) + ^ д/а(г,V,і) + Ра(г,і) д/а(г,V,і) _fд/а(г,V,і)^
- + V---------------------+
а
ді дг т„ дv
е
ді
+ ^а(Т, V, і) ,
(1)
А ^(г, і) _----------(пі (г, і)-пе (г,і)), Е(г,і) _-У(р(г, і),
ео
где первое из уравнений системы (1) - уравнение Фоккера-Планка для частиц сорта а (а = г,е ); второе - уравнение Пуассона для самосогласованного электрического поля; /а (г, V, *)
э/а( ^ ^ *)'
- функция распределения частиц сорта а;
- оператор столкновения Фоккера-
э*
\ У С
Планка; функция £а(г,V,t) описывает источники или стоки частиц; ¥а(г,*) = qaЕ(г,*), где
с
г / V \ I е, а = е,
Е (г, *) - напряженность самосогласованного электрического поля, qa= < ; р( г, *) -
|^, а = ',
потенциал самосогласованного электрического поля; па(г, *), а = г, е - концентрация частиц сорта а; та - масса частицы сорта а; £0 - электрическая постоянная.
Оператор столкновения Фоккера-Планка для ионов имеет вид [2]:
1Г /
Г І ді
_ 2 (?.?.&( V V, ‘)): (V„Vу/; ) + VуC, (г, V, і )■?,/, + 4л Г Це +Л] /,,
где VVVVgi (г, V, і) - ковариантная тензорная производная второго ранга, символ двоеточия (:) обозначает операцию двойного суммирования
, Сі (г, V, і )_(1^-Г
гУ ’ 12 •> |( - V '|
g, (r, ^ і)_ | /г (r, V ', і ) |( - V 1 ^ ' + 12 | Ле (r, * ', і) |( - ( 1 & ', 1 _ ^ , Г =
1 е
Г кТ ^2
Є¥
2
ГЛПЄ¥Є J
й\' тг
т„
Оператор столкновения Фоккера-Планка для электронов имеет вид [2]:
^Г Л] _2(^. (V.V.‘)): ) + V,C, (г,і)-V,/, + 4л Г^ +
Г Г ді
Л Л,
3 кТ
(
2е
кТ
Г лпа е ,
\ е« у
е 2
те
У Л IV - V
^е (г, 1?,*) = | /е (г,V',*) [V -V'| ё.'+ 221/(г,V',*) [V -V'| ё.'.
К системе уравнений (1) необходимо добавить начальные и краевые условия:
* = 0: £а ( . . 0 ) = .С, а = ^ e,
г еА р : /а( г, V, * )| = 0, а = ', е,
р ¿а\ ’ ’ '\геПр ’ ’ ’
р( .' )| Геа р =рр ■
Г еЙ- : £а( г Г * )| ?еП = /ат ^ а = ^ ^
р(.,*)|. _ = 0
гт
где /а _ Па
Г т Л ,па /2 Г та 2 Л
ехр V - V«
Г2клТ у а¥ Г 2кТ а¥ у
(2)
а _ г, е ; ^ - множество радиус- векто-
ров частиц, концы которых принадлежат плоскости зонда и границе возмущенной зоны соответственно.
2. Метод решения уравнения Фоккера - Планка
с
На рис. 1 показано расположение заряженной плоскости относительно набегающего потока плазмы в прямоугольной декартовой системе координат.
Рис. 1. Расположение плоскости в плазме
Из рис. 1 вытекает, что положение частицы в пространстве будет определяться координатами х, у, 2, а скорость - координатами пх, уу, Уг. В силу того, что плоскость является бесконечно большой в сравнении с характерным размером задачи, функции распределения частиц в этом случае будут зависеть только от переменных у , ^, і. Тогда модель (1) - (2) будет выглядеть следующим образом:
дЛ
д/а , ГУ дїа _ 1 д2/а д2ga , д/а дС(
•Л (Л У >Л (Л >Л (Л *—■’ (Л V (Л с
дУ та^у 2 д [^ ]2 д [^ ]2 дvy д^
+На, а_ г,е,
где На
д ф е і \
^ = --( Пг - Пе ) ,
ду е0
і_0: Ла(у,vy,0)_/а
у _ 0: Л а (0, Vy, і)_ 0,
Ф( 0, і )_ ф, у _ у« : Ла (у«, V, і )_ Га ф( у«, і)_ 0,
Е _-дФ Е _ ду '
, а_ г, е, а _ г, е,
а_ г, е,
4л I
Л 12
4л
12 Л ]
+ £ £, а_е.
г У У
(3)
Для решения поставленной задачи необходимо обезразмерить величины, входящие в систему уравнений (1) - (2). Для этого необходимо воспользоваться следующим соотношением. Пусть Мх - масштаб размерной величины X, X - безразмерная величина X, тогда
X = Мх • X . Таким образом, используется следующая система масштабов:
Мь = ГП
(ад»)7
М
еп
М„ = п
МФ М*~ Ф
М
м} = емпм:, (4)
ЬТ
м9=-^
Ма = Мп
а = і, е.
е М (ма)
За характерные масштабы в данном случае были приняты: радиус Дебая, скорость теплового движения частиц, концентрация частиц в невозмущенной плазме, потенциал, возникающий при разделении зарядов в дебаевской сфере и производные от них величины.
В результате была получена следующая система безразмерных уравнений:
Э/а. А Э/а. В Э/
Э* “ Эу
Э 2ф
Э£„
7 К Э/ „
+ Кад-Т- + Н
Э [V, ]2
а = і, е,
Э,
-( пі - пе ) ,
Е =-Э^
= ЭР’
(5)
(6)
с начальными и краевыми условиями:
* = 0: /а (у у,,0 ) = }а, а = ие,
У = 0: /а(0, ^ *) = 0, а = ^ е,
Ф ( 0, *) = фр,
(7)
у=у¥: а у» , у,, * )=/а
а = і, е,
р (.У-, * ) = о,
где безразмерные коэффициенты Д^, 5а, 2а, Ка, На, а = I,е представлены в [3].
Для решения уравнения Фоккера-Планка (5) воспользуемся модификацией метода расщепления [4]. Согласно данной модификации оператор правой части уравнения необходимо разбить на сумму двух операторов, описывающих два физических процесса (перенос и столкновение частиц):
С 2 ~ ~ >
в э а + К„-^ +На
О/а =-
А„^/а + Ва-Э* ЭУ
Эу
а/а
у у
э [V, ]2
Эу„
Тогда уравнение (5) в операторной форме запишется в виде:
(8)
Далее зададим сетку у = I • Иу, у1 = у • Ъу
■п • т, 0 < I < Ь, 0 < у < J, 0 < п < #, где
Ь, 3, N - количество шагов по пространственной переменной, скорости и времени.
Решение задачи (8) - (7) находится в результате последовательного решения двух вспомо-
гательных задач.
Первая задача:
2
п
^а( Я Уу, *) „ і л л
-----Э*-= О^а^ ^ *) , а = l, е,
Э2Ф (* * ) * дф
зуг=-( «і-«.). "*, =-эу.
^а( У, Уу, *« )= }а( Я Уу, *« ) , « = ^.^ Я-1,
У = 0: ™а(0, Vу, *) = 0 а = ^ е (9)
ф ( 0,*) = ф, у=у»: ^а(у», *, *)=/а, а=l, е, ф* (.V», *) = °.
Данная задача представляет собой самосогласованную систему безразмерных уравнений Власова - Пуассона. Для ее решения применяется метод крупных частиц [5]. Сущность данного метода состоит в следующем. Решение рассматриваемой системы осуществляется путем расщепления на два этапа: на первом этапе в системе не учитываются конвективные члены и решение получается обычным интегрированием на неподвижной эйлеровой сетке, а на втором
этапе рассматривается система, которая описывает перенос частиц в лагранжевой системе ко-
ординат.
Вторая задача:
дча( Я ¿у, 1) ( _ п
-----^Яа(У,V,,, 1), а = г,е, , ч
у’ (10)
Яа( У *У, ) = ™а( у ¿У, ^ ) , П = 0,..., N - 1.
Для решения второй задачи необходимо воспользоваться ее конечно-разностной аппроксимацией и к полученной системе разностных уравнений применить один из классических методов решения систем линейных уравнений, например, метод Гаусса [4]. Поставленные задачи решаются последовательно. Решая первую задачу, можно найти функцию
м/а( У, ¿У, 1п),п = 0,..., N, которая определяет начальное условие для второй задачи. Решая вторую задачу, находим функцию да( у, V,, 1п ) = /а( у, V),, 1п), п = 1,..., N, а = г, е, которая определит решение уравнения Фоккера-Планка для рассматриваемых моментов времени.
Для нахождения значений потенциала электрического поля ф из уравнения (6) необходимо использовать метод сеток [3], а для нахождения напряженности Е - методы численного дифференцирования [4]. Моменты функций распределения /а,а = г, е находятся с помощью применения методов численного интегрирования, например, метода трапеций [4].
3. Алгоритм решения задачи
1. Задать Ь, 3, N - количество шагов по пространственной переменной у, скорости V, и времени 1, ширину возмущенной зоны утах , верхний предел скорости ¿тах, правую границу Сах промежутка времени, на котором ищется решение, величину скорости .
2. Найти значения шагов по пространственной переменной, скорости и времени:
К
Ут
К
2у„
т
Ь ' Уу 3 ' N
3. Задать разностную сетку и положить п = 0 :
Уг = I' Иу, Vу = ] • И , )п = п -т, 0 < I < Ь, 0 < ] < 3, 0 < п < N.
4. Вычислить сеточное представление начальных распределений :
Г = Р/2ехр
V Т» те у
5. Решить уравнение Пуассона
Э ф
эр2
:-( «і - ^*е ) ,
согласно алгоритму, предложенному в [3], вычислив Фр«, 0 < I < Ь. 6. Вычислить напряженность электрического поля:
-2-(-3ф*п + 4*п -рп), ё; =-
пп ^1+1 - ^1-1
2К
0 < I < Ь,
(11)
ё*;
= - 2Т (фЬ-2 - 4 ф*Ь-1 +3 ф*Ь).
7. Решить задачу (5) - (7) предложенным выше методом решения уравнения Фоккера-Планка:
Перейти от задачи (5) - (7) к задачам (9) -(10).
Для решения задачи (9) необходимо в первом уравнении системы (уравнении Власова) перейти к уравнениям характеристик:
ЖУ = А Ж* А
Ж*а( У, *, * )
ж Ва
Аа = лЙЛ, Ва=4ё~а^~Ёу, $,=“,
2еа
Ж
еа =
= 0,
Т
а»
Т
та
Ма=~
да.
а = і, е.
Применить к полученной системе обыкновенных дифференциальных уравнений один из методов из семейства методов Рунге-Кутты [4] на п -м временном отрезке. Получить значения функции м>а (у, ¿у,)п+1), п = 0,..., N -1, а = г, е, которые будут использоваться в качестве начальных условий для задачи (10).
Получить разностную аппроксимацию задачи (10).
п
0
Решить задачу (10) после получения разностной аппроксимации методом исключения Гаусса [4]. Получить да(у,уу, ¡п+ ) = /а(у,уу,гл+1), а = г,е.
8. Вычислить моменты функции распределения Та(у,7, Т ), а = г,е для ^ ^
С 7 п+1 _ ~ 1-1 + Тп+1. -т 1
«I, ¡-у ' У ^ а
п+1
пп+1 =
=I
1=0
"е, 1-1
л1,1
■к, .с=да-
?■
I
1=0
■ V
“1,у у
■к.
0 < I < Ь,
2
2
£ Т т
$а=е, е,= « , т« = т«,
т« Тг¥ т
а = г, е.
9. Если п = N, вычисления завершить. Иначе положить п = п +1 и перейти к п.5.
4. Результаты численного моделирования
На основании предложенного алгоритма создано программное обеспечение, с помощью которого были получены следующие результаты:
Картина поведения функций распределения ионов (ФРИ) и электронов (ФРЭ).
Эволюция токов ионов и электронов.
Изменение во времени концентраций ионов и электронов.
Эволюция самосогласованного электрического поля.
Некоторые из этих численных результатов представлены на рис. 1-3.
Рис. 1. Сечение ФРИ
-1 0 1 скорость
Рис. 2. Изолинии ФРИ
Рис. 3. Эволюция тока электронов Из анализа полученных результатов вытекают следующие выводы.
Сечения ФРИ и ФРЭ имеют вид, представленный на рис.1, из-за того, что плоскость, помещенная в плазму, является идеально поглощающей для электронов, а ионы при ударе об эту плоскость нейтрализуются. Вероятность наличия ионов и электронов со значениями скорости, достаточной, чтобы попасть на заряженную плоскость, близка к нулю. Это отчетливо видно на рис. 2 , на котором представлены изолинии ФРИ. Изолинии представляют собой "обрезанные" концентрические окружности в сравнении с изолиниями начальных распределений ионов и электронов.
Из анализа эволюции тока электронов, представленной на рис. 3, следует, что ток достигает стационарного значения, имея небольшой переходный режим. Данное поведение тока характерно для случая малости членов в уравнении Фоккера-Планка, отвечающих за столкновения. При малых значениях данных членов решение уравнения Фоккера-Планка стремится к решению уравнения Власова. Данное поведение тока полностью согласуется с поведением, полученным в [1].
Таким образом, в данной работе построена математическая модель, описывающая поведение двухкомпонентной плазмы с учетом столкновений между заряженными частицами. Разработан алгоритм решения зондовой задачи с применением метода крупных частиц. Получены и проанализированы результаты численного моделирования поведения двухкомпонентной плазмы.
ЛИТЕРАТУРА
1. Гурина Т.А., Демков В.П., Котельников В.А., Попов Г.А. Математическое моделирование электродинамики летательного аппарата в разреженной плазме. - М.: Изд. Нац. акад. прикл. наук, 1999.
2. Олдер Б. Вычислительные методы в физике плазмы. - М.: Мир, 1974.
3. Пантелеев А.В., Кудрявцева И.А. Формирование математической модели двухкомпонентной плазмы с учетом столкновений заряженных частиц в случае плоского зонда//Теоретические вопросы вычислительной техники и программного обеспечения: Межвузовский сборник научных трудов. - М.:МИРЭА, 2006.
4. Киреев В.И., Пантелеев А.В. Численные методы в примерах и задачах. - М.: Высшая школа, 2006.
5. Белоцерковский О.М., Давыдов Ю.М. Метод крупных частиц в газовой динамике. Вычислительный эксперимент. - М.: Наука, Физматгиз, 1982.
THE APPLICATION OF THE PARTICLE-IN-CELL METHOD TO BEHAVIORAL ANALYSIS OF THE TWO-COMPONENT PLASMA WITH CHARGED PARTICLES COLLISIONS
Panteleyev A.V., Kudryavtseva I.A.
The mathematical model describing the two-component plasma behavior is developed as well as accompanying algorithm and software. For two-component plasma with charged particles collisions, results of numerical modeling are obtained. These results agree with the previous research.
Сведения об авторах
Пантелеев Андрей Владимирович, 1955 г.р., окончил МВТУ им. Н.Э. Баумана (1978), доктор физико-математических наук, профессор, заведующий кафедрой прикладной математики и механики МАИ, автор более 90 научных работ, область научных интересов - синтез оптимальных систем управления, методы оптимизации.
Кудрявцева Ирина Анатольевна, окончила МАИ (2005), аспирант кафедры математической кибернетики МАИ, автор 4 научных работ, область научных интересов - физика плазмы, зондовые методы исследования.