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

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

чаях, когда возмущения малы; сравнением с точными аналитическими решениями для задачи

свободного течения (Е = 0); это обсуждалось в § 1, п. 2; интегрированием по времени с последующим его обращением

(уравнения (1) и (2) точно обратимы во времени]; вычислением полной энергии, которая сохраняется согласно

уравнениям (1), (2) (полное число частиц точно сохраняется из-за самой формы преобразованных уравнений);

сравнением результатов численных расчетов с различными М, N и М.

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

Из линейной аналитической теории колебаний устойчивой плазмы следует, что Еп (t) ~ ехр [гсо (nk0)t], где время t доста­ точно большое, чтобы исчезли эффекты от начальных условий.

результаты известны из приближенных решений дисперсионного уравнения. Строго экспоненциальное поведение Е і (t) продемон­ стрировано на фиг. 10.

Уравнение свободного течения нейтрального газа [оно полу­ чается, если положить в (39) поле Е равным 0] описывает развитие больших градиентов скорости, которые будут обсуждаться в п. 2. Хотя уравнение свободного течения обладает свойствами, которые затрудняют численное интегрирование, его точное аналитическое решение получается легко для всех времен. Если использовать начальные условия типа (41), то точное решение для матрицы Zmn записывается в виде

Это решение сравнивается в табл. 1 и на фиг. 11с решением, кото­ рое получается в результате численного интегрирования уравне-

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

В недавнем расчете при неустойчивом начальном условии вида «горб на хвосте» [см. (43)] полная энергия сохранялась с точно­ стью до ІО“8 в течение времени 65 Юр1. Для сравнения укажем, что-

свободного течения для Z 0n (t), которое показывает величину ошибок округле­

ния и обрывания.

примерно 3 % от начальной кинетической энергии частиц перешло в энергию волн. Однако полная энергия относительно нечувстви­ тельна к существенным ошибкам в f (х, ѵ, t), которые не распро­ страняются на большие области фазового пространства. В упомя­ нутом расчете полная энергия сохранялась очень хорошо, несмо­ тря на то обстоятельство, что функция / (х, v, t) в малой области фазового пространства достигала небольших отрицательных зна­ чений. Тот факт, что функция / /х, v, t) в нескольких местах при­

нимала отрицательные значения, не был неожиданным в этом расчете, поскольку функция / (х , ѵ, 0) в начальный момент равня­ лась нулю на линии ѵ2 = 3; так что если ошибки при вычислении / случайны, то половина значений функции / (х , ѵ, 0), равная вначале нулю, смещается немного в область отрицательных значе­ ний. Участки отрицательных значений / едва ли влияют на точ­ ность необходимой информации об Еп (t).

Универсальной проверкой точности численного решения явля­ ется исследование его чувствительности к изменению размера используемых конечных разностей, в нашем случае At. В табл. 2 представлены итоги сравнения результатов для одного частного

расчета. Хотя оценки для ошибки обрывания в методах Рунге — Кутта неизбежно грубы, кажется, что отклонения в табл. 2 больше, чем можно было бы ожидать от случайно накопившихся ошибок обрывания; их величины соответствуют скорее прямому сложе­ нию ошибок. Даже в этом случае неопределенность в решениях, продемонстрированная в табл. 2, является вполне допустимой для любых известных приложений.

Дальнейшее обсуждение ошибок, особенно тех, которые воз­ никают из-за обрывания рядов Фурье — Эрмита, будет проведено в п. 2. В качестве последнего замечания к описанию машинных кодов нужно отметить, что наше внимание было сконцентрировано на информации, содержащейся в коэффициентах Z0,n (дает элек­ трическое поле) и Z2,о (дает кинетическую энергию частиц), по­ скольку эти величины представляют наибольший физический интерес. Можно было бы изучить другие матричные элементы и получить в общем подобные результаты, за исключением того, что для элементов вблизи границ матрицы ошибки будут значи­ тельно больше, чем указанные выше. К счастью, когда матрица достаточно большая, влияние граничных элементов на физически интересную информацию, содеряшщуюся в матрице, по-видимому, мало.

Легко видеть, что уравнения (39) и (40) не образуют замкнутой системы ни по индексу Фурье п, ни по индексу Эрмита т. Нелиней­ ный член, возникающий из выражения Е (df/dv), связывает каж­ дую моду Фурье с бесконечным набором других мод, а конвектив­ ный член, возникающий из выражения v (dfldx), связывает каж­ дую моду Эрмита с соседними. Рассмотрим вначале эффект обры­ вания рядов Фурье.

Требование разложимости функции / (х , v, t) и Е (х , t) в быстро сходящиеся ряды Фурье ограничивает класс начальных условий, для которых применим этот метод. В результате пригодными оказы­ ваются только те условия, для которых по крайней мере начальное состояние почти однородно, т. е. отношение

Zq,i (0) достаточно мало.

Zo, о(0)

Проблемы ударных волн или оболочек могли бы потребовать другого набора базисных функций для пространственного разло­ жения [38]. Обычным путем линеаризованное аналитическое реше­ ние уравнений (39) и (40) описывает эволюцию малых возмущений однородного начального состояния. Этот метод создает основу для изучения эффекта обрыва рядов Фурье в нелинейной задаче. Численно было обнаружено, что если в начальный момент в устой­ чивом состоянии (имеется затухание Ландау) возбуждена волна с п = 1 и амплитудой е, то потом нелинейно развивается п-я гар­ моника до амплитуды О (еп). Требуемая точность численных решений получается в результате выбора начального набора е и к0 и просчета задачи для нескольких возрастающих значений N, пока изменения (при фиксированном времени) величин, которые нас интересуют (обычно это Еп), станут меньше чем допустимые неточности. В табл. 3 сравниваются значения Z01 (t) и Z02 (t), полученные при одних и тех же устойчивых начальных условиях:

изменения Z0i (t) оказались меньше чем О (е3), а изменения Z02 (t) гораздо больше; сравнение затруднено, потому что введение третьей гармоники возмущает фазу второй. Если был а бы необ­ ходима более точная информация о второй гармонике, то требо­ вался бы счет с N = 4 для установления сходимости. Главным препятствием для расчетов с большим числом гармоник Фурье является время вычислений, которое возрастает как N 2.

При неустойчивых начальных условиях число N выбирается таким, чтобы необходимая сходимость получалась при ограни­ ченной амплитуде. Обычно, если N содержит две или более волн,