Минимизация численных погрешностей в динамических моделях крупных частиц

МИНИМИЗАЦИЯ ЧИСЛЕННЫХ ПОГРЕШНОСТЕЙ В ДИНАМИЧЕСКИХ МОДЕЛЯХ КРУПНЫХ ЧАСТИЦ

MINIMIZATION OF NUMERICAL ERRORS IN DINAMIC MODELS OF PARTICLE-IN-CELL METHOD

И.И. Воротынцева, A.B. Аникин I.I. Vorotyntseva, A.V. Anikin

ГОУ ВПО МГСУ, ОАО «ВЫМПЕЛКОМ»

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

Criteria are developed and carried out an empirical assessment of the errors of modeling processes in the crossed fields. Was selected the optimal combination of numerical methods for solving the mathematical model.

В литературе по моделированию методом макрочастиц, наряду с разработкой адекватной математической модели и построением вычислительных алгоритмов, большое внимание уделяется исследованию погрешностей численных методов [1,2,4], поскольку ключевым моментом в численном эксперименте является доказательство физично-сти получаемых результатов и выделение из них счетных эффектов. Для многих численных схем решения уравнений, описывающих физику процессов, получены теоретические оценки погрешностей и величин, характеризующих такие понятия как согласованность, точность, устойчивость, эффективность, дисперсионные соотношения сеточных схем и т.п. [2].

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

Система уравнений тестовой модели состоит из уравнений движения частицы, уравнения Пуассона и способов вычисления сеточной плотности пространственного заряда и дифференцирования сеточного потенциала. Рассматривается цилиндрический диодный промежуток в магнитном поле. В этом промежутке задается равномерная плотность пространственного заряда с известным распределением потенциала. Это равномерное распределение плотности моделируется с помощью датчика случайных чисел некоторым заданным числом макрочастиц Ne. Затем с катода «выпускается» частица и рассчитывается ее траектория в течение нескольких циклотронных периодов

4./2011 ВЕСТНИК _4/20|Т_МГСУ

(среднее время пролета пространства взаимодействия в генерирующем приборе). При этом «фон» пространственного заряда остается неподвижным, а частица испытывает отклонения от регулярной траектории, обусловленные суммарной погрешностью.

Уравнение Пуассона для потенциала поля пространственного заряда решается сверхбыстрым методом Хокни [5]. По существу этот метод представляет собой сочетание БПФ с циклической редукцией - своеобразной разновидностью способа исключения, и относится к классу прямых не итерационных, а, следовательно, и наиболее быстрых с гарантированной устойчивостью и экономным использованием памяти.

Моделируемая система является консервативной, т.е. в любой момент времени должен выполняться закон сохранения энергии + — С, где ШК,ШП — кинетическая и потенциальная энергия частицы на каждом шаге по времени, С - константа, определяемая начальными условиями вылета с катода. Невыполнение закона сохранения свидетельствует о наличии погрешности. В относительных единицах погрешность

где вычисляется аналитически, а — численно.

В качестве погрешности можно принять также «разогрев» частицы при движении.

Рассмотрим некоторую окрестность вокруг частицы. Мгновенная скорость частицы V отличается от средней скорости частиц у0 этой окрестности. Независимо от причин отклонения удобно разность скоростей назвать тепловой скоростью ут — V — г?0, а кинетическую энергию теплового движения И^ = 0,5ту£ - тепловой энергией (т-масса модельной частицы). Суммарную погрешность численных методов можно рассматривать как «температуру» частицы Т, которая определяется из соотношений

Щ = Еп=1 Щп, ШТп = ± ту^ уТп =уп- у^

T = 1Xn=iTn, WTn = 1-kbTn,

где vn- составляющая мгновенной скорости частицы, v0n- составляющая скорости «точного» решения в тестовой программе или средней скорости в модели прибора, d-размерность модели. «Точное» решение находилось численно с аналитическим определением силы, действующей на частицу, методом Рунге-Кутта-Фельдберга четверто-го-пятого порядка [3] с автоматическим контролем погрешности и коррекцией шага интегрирования.

Среди испытывавшихся методов NGP, CIC и TSC с соответствующими схемами вычисления силы [1,2], самым эффективным (в смысле соотношения точность - затраты времени) оказался CIC метод в сочетании со сглаживающим дифференцированием сеточного потенциала. С этой комбинацией исследовались также методы интегрирования уравнений движения.

На рис. 1-6, представляющих некоторые результаты испытаний, приняты следующие обозначения: AW- невыполнение закона сохранения энергии в относительных единицах (1), Т - температура (2), tc - время в циклотронных периодах, Nc - число шагов в циклотронном периоде, Ме - коэффициент укрупнения заряда частицы. Цифрами обозначены кривые соответствующие модифицированному методу Эйлера (1), методу Рунге-Кутта четвертого порядка (2), методу, сочетающему экстраполяционные формулы Адамса для скоростей и степенные ряды для координат (ЭФАСР) (3), прогноза-коррекции второго порядка (4), методу циклоид (5).

Па рис.1 видно, что методы 1,2,4 неустойчивы и приводят к сильному росту погрешности со временем. Неустойчивость возникает уже при плотностях пространственного заряда выше, чем 0.1 от бриллюэновской, поэтому для сильноточных приборов не имеет смысла использовать обычные методы интегрирования дифференциальных уравнений (типа Рунге-Кутта или предиктор-корректор) и схемы высокого порядка точности (см. кривые 1,2,4).

Рис. 1 Накопление во времени относительной погрешности интегрирования уравнений движения

tc

Рис. 2 Зависимость температуры частицы на десятом циклотронном периоде от шага по времени

250 200 150

m

о юо н

50

о

4/2011

ВЕСТНИК

_МГСУ

Рис.3 Зависимость относительной погрешности на десятом циклотронном периоде от шага по времени

Метод 2 Метод 1

Метод 5

Метод 3

5

и

20 25 30 38 40 46 60 66 60

Ис

Рис.4 Зависимость температуры частицы от коэффициента укрупнения заряда

Метод 1 Метод 4

И

Метод 2

Метод 5

-- Метод 3

-*- ---♦----

Ш

о

0.4

0.8 1.2 1.6 Ме х 107 = 50, сетка 64x33)

ВЕСТНИК МГСУ

4/2011

Рис. 5 Зависимость относительной погрешности от величины шага сетки

1.0 0.9 0.8 0.7 0.6 0.6 0.4 0.3 0.2 0.1 О.О

\ = 50 N, = 20

N„ = 40

N„ = 30

6.0 4.5 4.0

3.5 3.0

2.6 2.0 1.5 1.0 0.5 0.0

11/

Рис. 6 Эффективность методов интегрирования уравнений движения

Метод 1

Метод 1

Мвто Д 5

Метод 3

30 35 40 45 50

Мс

На рис. 2 и 3 показано поведение погрешностей на десятом циклотронном периоде от шага по времени. Коэффициент укрупнения заряда Ме не может быть больше некоторой величины при заданной сетке и , как следует из рис.4. На рис.5 показывает удвоение числа шагов сетки по обеим координатам, т.е. точки 1,2,3,4,5 соответствуют сеткам 16x9, 32x17, 64x33, 128x65, 256x129. Из этого рисунка видно, что шаг сетки А. имеет оптимальное значение, которое зависит от временного шага. Кроме того, оно связано с числом модельных частиц (т.е. Ме). С увеличением к растет погрешность вычисления силы, действующей на частицу, а при уменьшении - снижается мо-

4/2011 ВЕСТНИК _4/2011_МГСУ

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

Эффективность ^^ выраженная как произведение нарушения сохранения энергии ДЖ на затрачиваемое процессорное время и нормированное на минимум кривой для метода ЭФАСР, показана на рис.6. Видно, что наиболее эффективным (т.е. быстрым и точным) является метод ЭФАСР, а шаг по времени имеет оптимальное значение в минимуме функции

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

Литература

1. Рошаль А.С. Моделирование заряженных пучков. М., Атомиздат, 1979.

2. Хокни Р., Иствуд Дж. Численное моделирование методом частиц. М., Мир, 1987.

3. Форсайд Дж., Малькольм М., Моулер К. Машинные методы математических вычислений. М., Мир, 1980.

4. Аникин А.В., Шадрин А.А. Численная динамическая модель макрочастиц генератора М-типа. Симферополь, СГУ, 1988.-33 с.

5. Hockney R.W. A Fast Direct Solution of Poisson's Equation using Fourier Analysis // J. Assoc. Comput. Math. 1965. Vol.12, pp.95-113.

Literature

1. Roshal A.S. Modelirovanie zaryazhennyh puchkov. M., Atomizdat, 1979.

2. Hokni R., Istvud Dzh. Chislennoe modelirovanie metodom chastits. M., Mir, 1987.

3. Forsa d Dzh., Malkolm M., Mouler K. Mashinnye metody matematicheskih vychislenii. M., Mir, 1980.

4. Forsa d Dzh., Malkolm M., Mouler K. Mashinnye metody matematicheskih vychislenii. M., Mir, 1980. 4. Anikin A.V., Shadrin A.A. Chislennaya dinamicheskaya model makrochastits generatora M-tipa. Simferopol, SGU, 1988.-33 c.

Ключевые слова: Численные методы, численные погрешности, метод крупных частиц, скрещенные поля

Key words: Numerical methods, numerical errors, particle-in-cell method, crossed field

129337, Москва, Ярославское шоссе, 4-4-364

(916) 650-02-56 ^п^ег@шЬох .ги

Рецензент: Попов Олег Николаевич, д.т.н., профессор Московского государственного института электроники и математики