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

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

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

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

Добавлен: 11.04.2024

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

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

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

 

 

 

 

§ 2.

Прямые методы

 

 

 

 

 

147

 

 

 

1.

Описание программы

 

 

 

 

 

 

Вычисления проводятся на прямоугольной сетке х, у (NX Ф

Ф NY), причем шаги сетки могут

быть неравными (Н Х Ф HY).

В расчетах

может

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

N — 1,

N или N + 1

 

узлов

в зависимости

от рода граничных условий, но в любом случае

 

 

Заданный

граничный потенциал (IBCY=1)

 

 

 

 

NY=8 Г-

 

 

 

 

 

 

 

 

 

 

 

 

 

7 +

+

+

 

 

 

 

 

 

 

 

 

 

 

I

 

 

 

 

 

 

 

 

 

 

 

 

 

і

+

+

 

 

 

 

 

 

 

 

 

;

Ö

ß +

 

 

 

 

 

 

 

 

 

I

 

 

 

 

 

 

 

 

 

 

 

r

 

+

+

 

 

 

HX

 

 

 

 

!

I

5 +

 

 

 

 

 

 

 

 

 

 

+ r

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4 i

+

+

 

 

 

HY

 

 

 

 

 

1 ^

 

 

 

„ -L

+

+

+

+

 

 

 

 

 

 

+ ^

I

 

 

 

 

 

 

 

 

 

 

 

!

I

3 t

+ +

+

+

+

-*- +

+ + +

+

+

+

+

1

Qj

!

§

 

 

 

Заданный

заряд (Неизвестный потенциалj

 

 

 

 

 

 

 

 

L .

 

\

 

 

 

 

Ф и г . 1. Примерный вид задачи, решаемой с помощью подпрограммы РОТ1.

Сетка периодически продолжена в х-направлении (І В С Х = 3), а на границе по коорди­ нате у потенциал принимает заданные значения (1BCY — 1). Число узлов по направле­

ниям ж и у не совпадает (IQX

= 4,

I Q Y = 3), так что

основная ячейка сетки имеет

форму прямоугольника ( НХ ф

HY).

Вначале значения

потенциала заданы в обозна­

ченных точками узлах сетки, а значения зарядов — в узлах, обозначенных крестиками.

Подпрограмма

заменяет

эти значения зарядов

на вычисленные

значения потенциала

 

 

в тех же узлах.

 

 

число N

должно

быть степенью

основания

2

(NX = 2ІС?Х,

N Y = 2z<2y). Это

ограничение упрощает программу

разложения

Фурье, в особенности при использовании метода быстрого преоб­ разования Фурье [7—9]. Оно упрощает также процесс рекурсив­ ной циклической редукции на шаге «в» в подпрограмме решения уравнения Пуассона. По координате х программа допускает три различных типа граничных условий, задаваемых параметром ІВСХ, и аналогично три различных типа граничных условий, определяемых параметром IBCY, допускается по координате у. Возможны все девять комбинаций этих граничных условий. При

10*


148

Гл. 4. Методы расчета потенциала

этом

допустимые граничные условия задаются следующим

образом.

IВ СХ = 1. Потенциал на границе принимает заданные зна­

чения, сведенные

в массивы ф0, t

и Фіѵх.ю перед входом в под­

программу РОТ1.

По зарядам qlt t ■■■

<> заданным в узлах

1 . . . (NX — 1), программа вычисляет потенциал в этих

 

точках.

ІВСХ = 2. Поле равно нулю

при s = 0 и s = NX.

 

Выпол­

нение этого

условия достигается

предположением

о симметрии

Фі, t = ф_і, і

и фіѵх+і, t = ^ nx-i, t-

Заданное сначала

в

узлах

О . . . N X

распределение зарядов

q0 . .

. qnx. t программа

преоб­

разует в значения потенциала.

 

 

 

 

 

 

ІВСХ =

3. Все функции периодичны по х, т. е. ф_8, t — q>Nx-s, t-

Заданные перед входом в программу

в узлах 0 . . . N X — 1

значения

зарядов

q0, t . . . qNx-\

t

преобразуются

в

значения

потенциала.

 

 

 

 

 

 

 

 

Аналогично определены возможности IBCY = 1, 2, 3.

Для обращения к этой программе необходимо задать вели­ чины IQX, IQY, HX, HY, ІВСХ, IBCY и распределение заря­ дов, а в случае ІВС = 1 — граничные значения потенциала. Для примера на фиг. 1 показана сетка задачи при IQX = 4,

IQY = 3, ІВСХ = 3, IBCY = 1.2

2.Программа разложения Фурье

Все, за единственным исключением, прямые методы вычисле­ ния потенциала основаны на преобразовании Фурье с конечным числом членов. Поэтому мы сначала опишем подпрограмму FOUR67, предназначенную для выполнения вещественного фурьепреобразования заданной в п узлах вещественной сеточной функции.

Идея метода восходит к Рунге [10] и заключается в том, что члены с общим множителем группируются прежде, чем выпол­ няется умножение. Процесс группировки содержит только опера­ ции сложения и вычитания, и его можно интерпретировать как перегибание пополам графика заданной сеточной функции, после чего совместившиеся значения функции складываются или вычи­ таются. Хокни [11] систематизировал этот процесс. В прило­ жении А дана программа FOUR67 на языке Фортран IV. Пере­ гибание сеточной функции осуществляется подпрограммами TFOLD и TFOLD1. Такое перегибание продолжается до тех пор, пока синусоидальные функции сохраняют свойство симметрии, после чего коэффициенты Фурье выражаются существенно уко­ роченными суммами. В этих укороченных суммах уже нельзя избежать перемножений, и они вычисляются подпрограммой KFOLD с помощью рекурсивной процедуры, аналогичной алго­ ритму быстрого преобразования Фурье. Подпрограмма FOUR67


§ 2. Прямые методы

149

имеет два входных параметра IQ и ІВС, которые задают коли­ чество узловых точек задачи и граничные условия аналогично тому, как это определялось параметрами подпрограммы РОТ1.

ІВС — 1. Прямое или обратное синус-преобразование опре­ деляется формулой

71=1

 

Y k = ^ Zss i n ( ~ ) , l < Ä < r a — 1.

(4)

S = 1

Переход к фурье-образу и обратный переход по этой формуле приводят к умножению исходной сеточной функции на га/2.

ІВС = 2. Прямое или обратное косинус-преобразование опре­ деляется формулой

 

П

 

 

 

П =

2 E(s,

ra)2Zsc o s ( - ^ ) , OCfcCra,

(5)

где

s—О

 

 

 

 

 

 

 

 

Г 1/

У 2 ,

если і = О,

 

Е (і,

і) = \( . 1 .

в остальных случаях.

 

Последовательный переход к фурье-образу и обратно по формуле

(5) приводит к умножению исходной функции на га/2.

ІВС = 3. Разложение в ряд Фурье определяется формулами

Y h = 2 2 Е (к, I ) Zsc o s ( - ^ ) ,

OCAC - J - ,

(6)

Y h = 2 2 Z,sin

+ l < f c < n - l .

 

s=0

 

 

ІВС = 4. Переход от коэффициентов Фурье к исходной функ­ ции определяется формулой

 

п/2

 

 

cos I12nsk j 1

 

У ,=

2 2 Е

(»■ - г ) 1

 

71-1

k=Q

 

 

 

 

 

 

Г2ns (к—■п/2) П

0</s<;ra— 1.

 

2

 

(7)

ZhSin L

п

 

ft=n/2+l

 

 

 

 

 

Последовательное применение операций разложения в ряд Фурье

(6) и суммирования ряда (7) приводит к умножению исходной функции на 2га.


150

Гл. 4. Методы расчета потенциала

Перегибание сеточной функции осуществляется подпрограм­ мой TFOLD в соответствии с формулами

Л /с

__ у(0

у(^)

y(t+1)_

у'(0

+ z^_

L i о

----- z > m —

Лтп—S

^

 

 

=0, 1,

m

 

(8)

 

 

1 Г

 

Выполнение перегибания сеточной функции, состоящей из чи­ стых гармоник (а в виде линейной комбинации последних пред­ ставима любая функция), показано на фиг. 2. Мы рассматриваем

Ф и г . 2. Разделение гармоник с помощью перегибания.

Показано изменение содержания массива Z при разложении в ряд по 16 точкам. Во избе­ жание путаницы на каждом уровне изображено только несколько первых гармоник.

L — размер складки.

периодический случай (ІВС = 3). Первое перегибание отделяет синусоидальные составляющие от косинусоидальных. Последую-

§ 2. Прямые методы

151

щие перегибания и тех и других компонент постепенно разделяют гармоники на группы. Первая группа состоит из гармоник с нечет­ ными номерами, вторая — из гармоник с номерами, образован­ ными умножением нечетного числа на два, затем на четыре, на восемь и т. д. Когда гармоника уже не допускает дальнейшего перегибания, она принимает одну из двух форм:

/і-і

( я ks \

 

/і-і

 

 

t

ks \

2

ИЛИ

■ чгі

/ J

а*81. П( т і Г )

 

 

s—0

 

 

Zi

 

COS[ - r — ) ’

 

 

s—О

 

 

 

 

 

О <Ск

(нечетное) <

2h.

 

(9)

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

Такие укороченные суммы вычисляются с помощью рекурсив­ ной процедуры, сходной по замыслу и скорости с алгоритмом быстрого преобразования Фурье. Однако эта процедура имеет некоторые особенности, так как она была разработана для непо­ средственного вычисления вещественных сумм (9) по интервалам, измеряемым полуцелыми значениями я, минуя вычисление ком­ плексных сумм по интервалам, кратным 2я, которые входят

валгоритм быстрого преобразования Фурье. Введем рекуррентные соотношения

S?+i) (к) = S ? (к) + {cos(Ѳ) S%hK2m) (к) + sin (Ѳ) C%h/(2m) (к)},

Sf+m2m) (4/n — к) = — Sf'* (к) + (аналогично),

 

C(rt+1) (к) =

(к) + (cos(0) Cf+hK2n) {к) -

(10)

 

 

— sin (Ѳ) <Sf+/i/(2m) (Щ,

 

Сг+л/(2т) №т к) =

(к) — (аналогично)

 

для t = 0, 1 . . . Т — 1, где

Т — log2 h.

t.

Параметры

m

и Ѳ зависят

от

уровня рекурсии

 

 

in — 2l,

Ѳ

як

 

 

 

Тѣ

 

 

 

 

 

 

 

соотношения

(10)

применимы

для

 

 

г = 0,

1,

1,

1<;/Ь (нечетное)< ;— 1.

Начальные условия задаются значениями (1) и С® (1). Если S™ = IS'“’ = 0, то результатом применения рекуррентных


152 Гл. 4. Методы расчета потенциала

соотношений

будет

 

 

 

 

 

Л-1

 

 

 

s p

№ =

2

. (1) + CT (1)) sin ( - f £ ) ,

 

 

 

s — 0

 

 

 

 

 

Л -1

 

 

 

C<T>(*) =

2

( -

^ s(l) + C r ( l ) ) c o s ( - ^ ^ ) ,

(H)

 

 

s = 0

 

 

 

 

 

l<;/c

(нечетное) <12h—1.

 

Подходящим выбором начальных условий можно получить любые требуемые суммы. В случае периодических граничных условий

1

У ровень t

 

 

 

 

Ф ( * )

 

 

 

 

 

 

Ф ( к )

 

 

 

 

 

г

0

1

2

3

4

5

6

7

0

1

2

3

4

 

5

6

7

 

0

1

1

1

1

1

1

1

1

1

1

1

1

1

1

1

1

 

m

 

•— г —

ff

1tNi

 

•—г — ».------ h / 2 т------

 

1

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4 т - к

 

 

 

 

 

 

4 т ~ к

 

 

1

1

1

1

1

3

3

3

3

1

1

1

1

3

 

3

3

3

 

г

1

I

7

1

3

1

 

1

/

I

 

1

3

 

1

5

1

 

1

7I

3I

5 5

1/і

7 7;

 

31

5;

г

 

 

1

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4

3

/ 1 / 5 19

3 1/5 5 111 / 1 / 5 7 1 9 3 113 5 111

 

НС ум7нм ы

п о Нс и н у с а мн нС ум мны по к оНс и н уса мн

Ф и г . 3. Формирование и хранение сумм для различных значений к при рекурсивном обращении к соотношениям (1 0 ).

Целые числа в центральной части означают величину h в суммах на различных уровнях рекурсии, а стрелки указывают позиции, объединяемые в частичные суммы на каждом уровне.

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

Эти рекуррентные соотношения использованы в подпрограмме KFOLD, и на фиг. 3 показано, каким образом сочетается хра­ нение переменных в машинной памяти на разных уровнях рекур­

сии. На всех уровнях новые значения записываются на место старых.

Количество вещественных арифметических операций (т. е. операций сложения, вычитания, умножения и деления), необ­