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

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

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

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

Добавлен: 11.04.2024

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

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

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

§ 3. Численные методы

117

щуюся в точке В. Точки между С А и В D становятся не­ нужными и их прибавляют к списку «свободных мест».

В случае фиг. 9, б точки двойной кривой левее А, С и правее D, В можно включить в список «свободных мест». Кроме этого, необходимо образовать замкнутую кривую, точки которой должны располагаться друг за другом по часовой стрелке. Для этого

Ф и г. 9. Три случая срастания контуров.

звено кривой между точками А жВ нужно перевернуть. В програм­ ме оператор CALL REVLINK (ІА, IB) изменяет направление звена между точками А и В соответственно с серийными номерами ІА и ІВ; после этого можно соединить точки А я С, а также D и В.

Остается еще одна трудность, так как значения координат меж­ ду А и С и между D и В могут отличаться на кратное длине пери­ ода L. В этом случае вызов подпрограммы MASCUR (ІС) исправит положение, и координаты точек на вновь образованной кривой ІС будут изменяться непрерывным образом.

На фиг. 9, в показан третий случай, когда двойная кривая разделяет две области 1 я 3 с разными значениями /. Вначале эти

118

Гл. 3. Модель «водяного мешка»

области не были смежными. В этом случае мы можем отбросить участок С D, соединить точки С и А, а также В и D, соответ­ ствующим образом изменить направление присоединенных участков и определить кривые нового типа (тип 2), которые соединяют точки А и В двух других контуров и для которых Д/ — / з — f 1. Эта возможность еще не использовалась.

5. Расчеты средней функции распределения

Большой интерес часто представляет средняя функция распре­ деления, которую можно рассчитать по формуле

І(ѵ) = Ь~1

J /(ж,

v)dx = L~1 j d x ^ A/j-Ѳ[к Vj(x)\.

(41)

 

 

i

 

Практически

более

эффективным оказывается вычисление

df(v)/dv и последующее интегрирование этой производной по ѵ. Из уравнения (41) мы имеем

■%- = Ь“1 2 M l j dx б [г; — Vj (ж)] =

3

рdx(vj)

= L_1 2 M i ) d v j d v j - б (V — V j ) =

3

= b -l 2 b f > ) ^ 7 d v l . (42)

3

Чтобы вычислить эту сумму, разделим фазовую плоскость на 2N горизонтальных эйлеровых полос шириной At;, как показано на фиг. 10. Для каждого контура Сj «средний наклон»

связан с п-й горизонтальной полосой, центр которой лежит при i/г = (п + Ѵ2) Ак, (—N ^ п ^ N — 1), следующим . соотно­

шением:

dx \ i

{ Ах \ .

2

-^ -[« (р )—*(«)].

(43)

по всем подсегментам

где Xj (ß) и Xj (а) — координаты точек контура Сj, в которых он пересекается с границами полосы. Поскольку контур просма­ тривается в направлении соединяющего отрезка, а встречается раньше ß.


§ 3. Численные методы

И

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

и/йѵ

Ф и г . 10. Аппроксимация формы контура для расчета F (ѵ).

-х', лежала в той же самой полосе, то разность координат (х" х ’) прибавлялась к текущему значению Ах этой полосы. Если эти точки лежали в соседних полосах, то интерполированием находи­ лась точка пересечения x'" и разности {х"' х') и (х" х'") прибавлялись к текущему значению Ах соответствующих полос. Кроме этого, программа могла анализировать более общие слу­ чаи, когда пару лагранжевских точек разделяло более одной горизонтальной границы или когда точки лежали за пределами сетки скоростей.

В результате такого просматривания определяется (Ах!АѵѴ

 

 

 

Ѵп+і/г

для контура любой

формы. Затем среднее распределение

/ (ѵп)

можно вычислить как сумму

 

 

2

а /, 5 ' ( ■ £ ) ; « „

- N < n < N .

(44)

 

 

Іp = - N

120

Гл. 3. Модель «водяного мешка»

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

§ 4. У ст о й ч и во ст ь м ет ода «с п ер еш а ги в а н и ем »

Лагранжевские точки на контурах удовлетворяют уравнениям движения частиц (2). На первый взгляд могло бы показаться, что решение уравнения Власова с помощью вычисления движения этих точек аналогично расчету движений индивидуальных частиц по программе, в которой вычисляются траектории большого числа частиц. В любом случае выгодно слегка поварьировать координа­ ты, скорости и ускорения, чтобы при интегрировании по времени добиться точности до величины второго порядка малости. Таким образом, в случае, когда силы не зависят от скорости, ускорение в момент t = геДг можно использовать для вычисления прира­ щения скорости от момента (п — V2) At к моменту (п + 1/2) Аt,

а эту скорость — для расчета приращения координат от момента

пМ к моменту

(n + 1) Аt,

т. е.

разностная

схема для частиц

имеет вид

 

[ ( и ~ т )

 

(

ДО»

«ДО ДО

(45)

Ѵі

[ ( ” + т )

=Ѵі

 

 

 

 

 

п

 

 

 

 

 

х і [(« + 1 ) ДО = % і [ п

Аі] +

Ѵ і

+

y

)

Д*] ДО

 

Эту схему можно использовать для вычислений, так как плотность заряда, которая нужна для расчета Е в момент t, зависит только от координат x t (t) и не зависит от скоростей, которые, следова­ тельно, не обязательно должны быть определены к этому моменту времени. Однако если используется модель «водяного мешка», то плотность заряда определяется формой контура, которая в каж­ дый момент зависит как от координат, так и от скоростей в тот же самый момент времени. Поэтому требуется более разработан­ ная разностная схема *).

г) Как недавно сообщил нам К. Саймон, интеграл j fdv, определяющий

плотность заряда, может быть вычислен, если использовать значения ско­ ростей в момент времени, смещенный на Аі/2, так как, согласно первому уравнению системы (45), все контуры при фиксированном х смещаются в этом случае как единое целое в вертикальном направлении. Из этого следует, что уравнения (45) можно использовать для расчетов с точностью до членов второго порядка по Д£ и в случае модели «водяного мешка». В период написа­ ния программы мы не обратили внимания на это обстоятельство и использо­ вали метод «с перешагиванием», который описывается в данном параграфе.


§ 4. Устойчивость метода «с перешагиваниемъ

121

1. Нечетные и четные фазовые пространства

Метод, позволяющий проводить вычисления с точностью дочленов второго порядка,— это схема «с перешагиванием». В этом методе используются два фазовых пространства: нечетное про­ странство 5 0, в котором координаты и скорость определены в «не­ четные» моменты времени (2п + 1)А4 и четное пространство S e, в котором эти величины определяются в «четные» моменты времения 2nAt (п — 0, 1, 2, ...). Четное пространство определяет движе­ ние контуров в нечетном пространстве, и наоборот. Разностныеуравнения имеют вид

х0 [(2п +

1)А4 =

х 0 [(2п — 1)А4 + 2ѵе (2nAt)At,

Vo \(2n +

1) A4 =

v0 l(2n — 1)

A4 -f 2Ee (2nAt)

At,

 

 

 

 

 

(46).

xe [(2n +

2) A4 =

xe (2nAt) +

2v0 [(2n +

1) A4

At,

ve t(2n +

2) A4 =

ve (2nAt) +

2E 0 [(2n +

1) A4 At,

где индексы о и е обозначают соответэтвенно нечетные и четные значения.

2. Неустойчивость расчета

Трудность, с которой приходится сталкиваться при использо­ вании различных схем «с перешагиванием», состоит в том, что могут возникать ложные движения точек. Эти движения возникают

Пространство S0

Пространство Se

Ф и г . 11. Сопряженные нечетное и четное пространства.

Форма контура в Se используется для расчета изменения формы в S 0 , и наоборот.

потому, что в этом случае при численном интегрировании степеней свободы в 2 раза больше, чем в реальной физической системе. В случае нашей модели каждый физический контур определяется двумя сопряженными кривыми С0 и Се, как показано на фиг. 11. Каждой точке (х0, ѵ0) соответствует сопряженная точка (хе, ѵе) на Се. Конфигурация контуров С0, е позволяет рассчитать поля Е 0і е


122

Гл. 3. Модель «водяного мешка»

в лагранжевых

точках х0і е, а затем, используя ѵ0, е и Е 0, е,—

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

В разностной схеме с конечным At нельзя ни определить конту­ ры С0 и Се точно в один и тот же момент времени, ни вычислить их истинные движения. Следовательно, нельзя требовать, чтобы они все время оставались согласованными. Асинхронные движения обязательно возникнут. Поэтому, чтобы предотвратить возник­ новение больших расхождений между двумя системами кривых, их необходимо периодически синхронизовать, усредняя сопря­ женные координаты.

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

Проанализируем теперь асинхронную неустойчивость с коли­ чественной стороны и опишем используемый в настоящее время

метод синхронизации.

Более

детальное рассмотрение

приведено

в приложении 2. Так как в

пределе Д£—ѵО неустойчивость не

исчезает,

то мы рассмотрим сначала этот более простой случай,

который

описывается

системой дифференциальных

уравнений.

3. Синхронная и антисинхронная моды

Рассмотрим равновесное состояние, в котором функция рас­ пределения представляется N горизонтальными контурами Cj, определяемыми уравнениями vj (х) = Vj = const. Скачок в ве­ личине/на каждом контуре равен Afj = fj fj-i- В случае равно­ весия нечетные и четные контуры идеально синхронизованы.

Линеаризованные уравнения движения, описывающие возму­ щение равновесного состояния, имеют вид

(47)

(48)

(49)