Файл: Вычислительные методы в физике плазмы..pdf

ВУЗ: Не указан

Категория: Не указан

Дисциплина: Не указана

Добавлен: 11.04.2024

Просмотров: 343

Скачиваний: 6

ВНИМАНИЕ! Если данный файл нарушает Ваши авторские права, то обязательно сообщите нам.

§ 3. Модели плазмы, с малым ß

283

Тогда уравнение (78) принимает вид

 

d2i|)

1 dip

»m*ъ2 IP:

 

dr2

dr

r 2

 

 

 

 

r-\-a

re+ (Я, г) R dR

 

 

- 4леС[т J

- лг_ (r, o ] - (84)

 

[4Д2г2__(Д2_(_г2 _ а2)2]

 

 

 

i/2

Магнитное и электрическое поля того типа, который мы рассматриваем, приводят к появлению компонент дрейфовых ско­ ростей в направлении z, однако эти компоненты не включаются в уравнения. Из формулы (81) можно получить выражения для /•-компоненты дрейфовой скорости ионов и ее дивергенции:

== - -Т Ж Ь 5+ 4

Ѵ2Ф] *ітф.

(85)

divv+= —

(-§§•) +

Ѵ2ф]

(86)

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

сать в виде Ѵ+ = фо Ѵ+(г),

где

 

 

У+ (г)

= г (Qm +

Qe).

(87)

Предполагается, что функция распределения энергии ионов имеет вид 6-функции. Частота прецессии из-за градиента магнитного ноля равна

о _ сГ B z dB

,q„.

Ü M - ~ ^ W ~ d F ’

(88>

где Т — энергия ионов. Электрическая составляющая частоты пре­ цессии, вызываемой электрическим полем нулевого порядка, равна

Г

 

ec (1 — Г) j r\+(r')r' dr',

(89)

где ц+ (г) — невозмущенная плотность ионов.

В формуле (89)

использовано предположение, что ионная и электронная плотности имеют одинаковую пространственную зависимость, т. е. г]_ (г) =

= Гг]+ (г), где

Г — постоянная

величина. С

учетом

соотноше­

ний (85)—(89) уравнение (79) принимает вид

 

 

 

дп_і_

ІСП:

i f f

^

+

~ Ѵ3ФJ =

О,

(90)

 

dt

где

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

G(r) = ^ - V +(r)

■m

( & м Ң- £2_g),

 

 

 

 

 

B z

dN+

-N+

d ( BzY]

 

 

 

 

В 2

dr

 

dr \ B 2 } \ '

 

При решении

уравнения (90)

мы подставляем Д2ф из уравнения

 

Ѵ2ф = — 4лес I~п+п_ -f-

V2/i+J

 


284

Гл. 7. Модели плазмы без столкновений

Это уравнение получается подстановкой в уравнение Пуассона выражения п++ (а2/4 )Ѵ2п+ для плотности ионов, которое содер­ жит слагаемое второго порядка по ларморовскому радиусу ионов а Уравнение (90) можно теперь переписать в виде

+ iGn++ іН

леса2 ^ п+—

Ѵ2п+ ) J =

0.

(91)

Уравнение (91) имеет то преимущество,

что оно сходно

по

форме

с уравнением диффузии

относительно

n+(r, t) и для

его

реше­

ния можно использовать неявную разностную схему. Аналогичным образом нетрудно убедиться, что уравнение (80)

примет вид

 

 

 

 

 

о,

(92)

где

 

 

 

P(r) = ^ V - ( r ) = mQE, Q(r)

r

L B l dr

dr \ В 2 } j

 

Уравнения (84), (91) и (92) являются основными в нашей модели. Они образуют систему уравнений с частными производны­

ми для комплексных функций ф, щ

и п_.

Если положить

 

ф = V + iW, п+ = р +

іу,

= Ѳ + гб,

(93)

и приравнять действительные и мнимые части уравнений, то полу­ чится система действительных уравнений

d2V

1 ' dV

 

 

г-fa

р (R,

t) R dR

 

V = 4лес^~

j

 

 

 

 

— Ѳ(М)] ,

dr2

1 r dr

т а 2 _ Т /

 

 

9

'

 

[4 R г г і _ ( R

Г 2

 

 

2 _|_ Г 2 _ а 2) -.2f i 2

 

 

 

 

 

 

(94)

d2W , 1 dW dr2 1 r dr

m 2 -

W:

r 2

г-\-а

 

у (В, t) R dR

 

 

 

= ~ 4лес [т J [4Л2r2 —( № + r2 —a2)2]1/2

 

= Gy + H W - n e c a 2H [ y - ö - ' r

a2

j dfy2y_ .

1

dy

 

 

l4 dr2

r dr

dy_

 

2 / д2р

, 1

dp

dt

■= — Gp — HV + леса2Н £p — Ѳ+

 

^

r

dr

dddt = P8 + QW,

db = — PQ— QV.

dt

i(r, o ] ,

(95)

4t ) ] . (96>

^P ) ] . (97) (98)

(99)


§ 3. Модели плазмы с малым ß

285

Мы хотим найти решения уравнений (94)—(99), удовлетворяю­

щие следующим начальным условиям: при t

= 0 действительные

части возмущений плотности

электронов

и ионов имеют вид

 

 

 

р (г) =

Ѳ(г) =

го

 

для

г <

г0,

 

 

 

р(г) =

 

 

 

 

г > г 0,

 

 

 

 

 

Ѳ(г) = 0 для

 

 

где

гп = гмакс — а и а — ларморовский

радиус

ионов. Гранич­

ные

условия

при

г = 0:

 

dy _ dQ

 

 

 

 

 

dV _ dW

 

_ dp

q

 

 

 

 

dr

dr

 

dr

 

dr

 

dr

dr

 

 

Граничные условия при г =

гмакс;

 

 

 

 

 

 

 

 

 

F = H / = p = y = Ѳ = 8 = 0.

 

 

 

 

 

 

б.

Разностные методы

 

 

Итак, мы будем решать задачу с начальными условиями для

уравнений

(94)—(99)

в

области

 

 

 

 

 

 

 

 

 

 

0 < Г <

Гмакс,

t >

0.

 

 

Рассмотрим

в этой области конечно-разностную

сетку вида

 

 

 

Г) = /Аг, / = 0, 1, 2, 3,

...,

/ ,

 

тде

г,/ = гмакс =

/Д г и tn =

nlS.t,

п = 0,

1, 2,

3,

... . Мы исполь­

зуем

обычные обозначения, т. е.

 

 

 

 

 

 

 

 

 

 

Ѵ? = Ѵ (rj, tn) и т. д.,

 

 

и следующие

разностные аппроксимации:

 

 

 

 

/ ЗУ

у

V l

-

VU

(

dW

у

^ yn+ i - 2 V n + VU

 

\

dr

 

 

2Дг

 

’ \

дг2

/;'

 

(Дг)2

 

Для уравнений (94) и (95) можно записать разностную аппроксима­ цию вида

-

(ijVrj +

bjV y1-

CjVy_l =

dry \

 

 

(100)

-

cijWyi -'r bjW y 1 -

CjWyi = knj +1,

 

(101)

-где

1

 

 

 

 

 

 

 

1

h

2

m2

 

1

1

 

a j ~ (Ar)2

2rjAr

(Дг)2 ~ ~ r f '

C}~

(Ar)2 _

dj+1= к

r+a

 

p (fl) R dR________ y + i

 

ѲГ1]

,

 

 

 

 

[4fl2r2 _ (fl2 +

r2 _ a2)2]V2 ) j

 

r —a

jfcf1 = к

г 2 / " г 0 _______ у (fl)

L i r l J [4fl2r2 — (fl2

 

r—a

R dR________

\ n + 1 _ g n + i l

r2 — a2)2]1/2

/ j

J J ’


286 Гл. 7. Модели плазмы без столкновений

k = 4яес, а — ионный ларморовский радиус. Интегралы в форму­

лах находятся численно в точках с использованием значений

р"+ 1

ау ”+1, / = 0 , 1 , 2 Подынтегральные функции имеют

осо­

бенности в граничных точках, но интегрируемы. В окрестности

граничных точек

использовалась

теорема

о

главном

значении,

а в оставшейся

области, г — а +

Д ^ 7і

^

г Jr а

А, приме­

нялась формула трапеций. Для обеспечения нужной точности интегрирования использовалось не менее десяти точек, так что

на отрезке а умещалось не менее пяти точек г;-. Для малых значений г, когда область интегрирования содержит меньше десяти точек гj, требуется интерполяция.

Для разностной аппроксимации уравнения (96) используется следующее неявное разностное уравнение:

р П + і ----р п

.

- J ж = Т ъ (ѵГ1+ чЪ + HjW] +

щн, [уГ-бГ1+-т(-алм+ълГ1-с№!)] +

+ka)Hj [у? — 6? + -^- ( — ajy% 4+ bjyf - crf$-i)] -

Аналогично аппроксимируется и уравнение (97):

7-1

— jG j (рГ 1 + Pi) -

HjV] -

4

kaW i [ p r 1 - ѳ г 1 +

At

 

+ 4

( - ^ р м + ^ рГ

- ^ р"-!1)] -

 

-

ka)Hj- i [ p

?

Ц- -

( aj- pѲ ? ? + bjp71 + ]+ - cjp^) ]

Написанные выше разностные уравнения можно представить в виде системы линейных алгебраических уравнений

 

Ajyr;: i + р

г - B j y ^ +

Cjy]_+i =

d -,

(102)

-

A j Р - 1 +

B jpJ+1+

тГ 1 -

Cjp™ = К Ь

(1оз>

где

 

 

 

 

 

 

 

 

Aj = 22 ka*HjdjAt,

 

 

Bj ==y GjAt -\-^-ka2jHjAt -f-^ kajHjbjAt,

 

 

Cj = — kajHjCjAt,

 

 

Dj = И;-у™+1 -f- p" -j- Bjy1}

+

HjW’jAt

g kaj (öf -(- 8j+1) At,.

K nj = Ajpl, -

Bjp] -f y? + Cjp?„, -

HjVjAt + 4

kä) (&) +

Ѳ"+1) At.


§ 3. Модели плазмы с малым ß

287

Уравнения (98) и (99) содержат только производные по времени, и мы аппроксимируем их разностными уравнениями

Ѳ Г = ѲГ1 + 2Аt (Pj6] + QjW?),

(104)

б г = 6Г 1 ~ 2Al (P ßf + QjV]).

(105)

Рассмотрим вопрос о решении системы уравнений (100)—(105). При фиксированном значении времени величины с верхним индек-. сом п и п 1 известны, а величины с верхним индексом п + 1 неизвестны. Первый шаг состоит в вычислении Ѳ”+1 и öf+1 для всех значений / по формулам (104) и (105). Следующим шагом служит вычисление р" +1 и y f+1 для всех значений /. Чтобы решить системы уравнений (102) и (103), запишем их в виде одной системы

~ ÄjV]:} + B j f r 1 - C j V ^ = ф?,

(106)

где

1

0

1

_

 

,

l — Aj

0 J 1

'р Г '

II

'D T

 

 

 

>

Фі =

 

 

 

B i = '' 1

- Я Д ,

— C; = ' 0

C f

J

0

?

J

_ — Cj

0_

 

 

 

Для решения системы (106) мы воспользуемся следующим алго­ ритмом:

ѴГ1= ËjV^l + J r \ / = = 0 , 1 , 2 , . . . , / - 1 ,

где матрицы Е }- и векторы / ? +1 определяются рекуррентными соот­ ношениями

Ei — ( E j С j E j - 1) 1 Aj,

/ Г = ( B i - CjËj-іУ1(ф? + C f f l ).

Из граничных условий при г = 0 следует

' 1 0

0

"

£o = 0 1 _ и / Г =

0 _

Вычисления выполняются в два этапа. На первом этапе вычисля­ ются все коэффициенты Ej и /™+1, на втором этапе вычисляются

значения F"+1, начиная с тех F"+1, которые задаются с помощью величин р (гмакс, 1) и у (гмакс, *). Вычисление величин Vf+1 и lFf+1 для всех значений / составляет третий шаг вычислительного цикла для момента времени t. Метод решения уравнений (100) и (101) такой же, как только что описанный метод решения уравне­ ний (102) и (103).