Диссертация (1149802), страница 7
Текст из файла (страница 7)
. . − α1 , . . . , cNp ).∂c k−1kcNpcNp(2.28)45Здесь α ∈ (0, 1). Для остановки итераций метода Ньютона можно использоватьследующий критерий: k I − I k−1 < εI ,Ik(2.29)где I k — ток эмиссии на k-й итерации; εI — заранее определенная точность.Выбор метода решения обусловливается тем, что для вычисления производнойфункции (2.27) требуется большой объем вычислений, при этом модифицированный метод Ньютона позволяет использовать только производную, рассчитанную на первом шаге.В соответствии с вышесказанным приведем оптимизационный алгоритм нахождения тока эмиссии:1.
Выбирается равномерное начальное распределение плотности тока эмиссии, итерационным методом решается задача (1.18)–(1.23).2. Вычисляются матрица Якоби функции (2.27) c использованием каких-либометодов численного дифференцирования.3. Определяется новое приближение тока эмиссии с помощью (2.28).4.
Решается задача (1.18)–(1.23).5. Пункты 3-4 повторяются до достижения выполнения условия (2.29).2.4.3 Модель эмиссии Гаусса в итерационном методеАлгоритмы для нахождения тока, ограниченного пространственным зарядом,основанные на применении закона Гаусса для ячеек сетки, находящихся вблизиэмиттера, используются во множестве программ,реализующих метод частиц вячейках [87]. Применение этого подхода в методе частиц в ячейках не составляет труда и дает достаточно точный результат. Далее в предлагается алгоритмприменения модели эмиссии на основе закона Гаусса к итерационному методу [58].Создадим дополнительную сетку, состоящую из прямоугольных гексаэдральных ячеек ("эмиттирующих ячеек"), таким образом, что одно ребро каждой из46ячеек лежит на эмиттирующей поверхности.
В дальнейшем рассмотрим алгоритм для двумерного случая (рис. 2.5). Для трехмерного случая алгоритм строится аналогично. Обозначим каждую эмиттирующую ячейку как Cj , j = 1 . . . NCи длины их ребер как Lem и Hem .Длякаждойэмиттирующейячееки с номером j будет выполняться закон Гаусса∫En dS = Qi /ε0 .(2.30)∂CiЗдесь ∂Ci — поверхность Ci ячейки, En — нормальная компонентаэлектрического поля к поверхноРиcунок 2.5: К модели эмиссии частиц.сти ячеки, Qi — пространственныйПунктирными линиями обозначенызаряд, находящийся в Ci ячейкеэмиссионные ячейки, сплошной жирнойсетки.
Для задач с током, ограни-линией обозначена поверхность эмиссии,ченным пространственным заря-сплошными линиями обозначены ячейкидом, поток вектора электрическо-расчетной сеткиго поля через поверхность эмиссии будет равен нулю. Таким об-разом, численно в двумерном случае закон Гаусса (2.30) для Ci ячейки можетбыть представлен как3∑fjs = Qj /ε0 .(2.31)s=1Здесь fij — поток вектора напряженности через расположенные в вакууме ребраj эмиттирующей ячейки. Эти значения известны из решения уравнения поля накаждой итерации метода. В трехмерном случае интегрирование закона Гауссабудет производиться аналогично: необходимо использовать потоки вектора напряженности через расположенные в вакууме грани эмиттирующих ячеек.
Вве-47дем новую индексацию для макрочастиц: сопоставим номер jp макрочастицам,эмиттированным внутри ячейки Cj p = 1, . . . , Njp . Обозначим суммарный токвсех частиц эмиттированных внутри одной ячейки Cj как ICj . Значение ICj можно вычислить следующим образом:ICj =Njp∑Ijp .k=1Общий заряд находящийся внутри некоторой эмиттирующей ячейки с номером i может быть представлен в виде линейной комбинации суммарных токовICjQi =NC∑aij ICj .(2.32)j=1Здесь aij — коэффициент, определяющий вклад тока ICj в суммарный заряд Ciячейки.Рассмотрим алгоритм вычисления коэффициентов aij . Сопоставим числоKjp = Ijp /ICj каждой макрочастице с номером jp . В процессе расчета необходимо вычислять пространственный заряд, вносимый траекториями всех макрочастиц в ячейки сетки следуя приведенному выше алгоритму.
Каждый k-й участокk−1/2[r jpk+1/2, r jp] траектории макрочастицы с номером jp вычисленный с помощьюсхемы с перешагиванием, или Бориса, полностью расположенный внутри некоk−1/2торой эмиттирующей ячейки Cp (r jpk+1/2∈ Cp и r jp∈ Cp ) вносит заряд ∆tjp Ijpв нее (∆tjp временной шаг интегрирования для макрочастицы с номером jp ).Введем функцию F (r1 , r 2 , C, jp )F (r1 , r2 , C, jp ) =∆tj /c,p0,еслиr1 ∈ C and r2 ∈ C,иначе.Таким образом, мы может вычислить суммарный заряд, вносимый траекториямивсех макрочастиц в эмиттирующую ячейку Ci :48Qi ==∑∑∑∑pjICjk+1/2, r jp, Ci , jp )Ijp =k∑∑pjk−1/2F (r jpk−1/2F (r jpk+1/2, r jp, Ci , jp )Kjp .kИспользуя это выражение, мы можем вычислить коэффиценты aij в уравнении(2.32) следующим образомaij =∑∑pk−1/2F (r jpk+1/2, r jp, Ci , jp )Kjp .(2.33)kПодставив (2.32) в (2.31) можно получить следующую систему линейных уравнений:AI = ε0 F.Здесь A — матрица коэффициентов aij , F ={j=4∑j=1f0j , .
. . ,j=4∑}fNC j(2.34)— вектор,j=1содержащий суммарные потоки вектора напряженности электрического поля через расположенные в вакууме грани эмиттирующих ячеек, I = {IC0 , . . . , INC 0 } —вектор токов, проходящих через эмиттирующие ячейки. Таким образом, мы можем представить алгоритм вычисления распределения плотности тока эмиссии,ограниченного пространственным зарядом:1. после задания начальных данных, необходимо поставить в соответствиекоэффициент Kjk каждой частице;2. в процессе расчета траекторий частиц, вычисляются элементы aij матрицыA используя уравнение (2.35);3.
после решения уравнения Пуассона, вычисляются элементы fNC j вектораF;4. на следующей итерации используется приближение к плотности тока эмиссии, полученное из решения системы линейных уравнений (2.36)I = ε0 A−1 F.Рассмотрим алгоритм вычисления коэффициентов aij . Сопоставим числоKjp = Ijp /ICj каждой макрочастице с номером jp . В процессе расчета мы долж-49ны вычислять пространственный заряд, вносимый траекториями всех макрочастиц в ячейки сетки следуя приведенному выше алгоритму. Каждый k-й участокk−1/2[r jpk+1/2, r jp] траектории макрочастицы с номером jp вычисленный с помощьюсхемы с перешагиванием, или Бориса, полностью расположенный внутри некоk−1/2торой эмиттирующей ячейки Cp (r jpk+1/2∈ Cp и r jp∈ Cp ) вносит заряд ∆tjp Ijpв нее (∆tjp временной шаг интегрирования для макрочастицы с номером jp ).Введем функцию F (r1 , r 2 , C, jp )∆tj ,pF (r1 , r2 , C, jp ) =0,r1 ∈ C and r2 ∈ C,еслииначе.Таким образом, мы может вычислить суммарный заряд, вносимый траекториямивсех макрочастиц в эмиттирующую ячейку Ci :Qi ==∑∑∑∑pjICjk+1/2, r jp, Ci , jp )Ijp =k∑∑pjk−1/2F (r jpk−1/2F (r jpk+1/2, r jp, Ci , jp )Kjp .kИспользуя это выражение, мы можем вычислить коэффиценты aij в уравнении(2.32) следующим образом∑∑k−1/2k+1/2aij =F (r jp , r jp , Ci , jp )Kjp .p(2.35)kПодставив (2.32) в (2.31) можно получить следующую систему линейных уравнений:AI = ε0 F.Здесь A — матрица коэффициентов aij , F ={j=4∑j=1f0j , .
. . ,j=4∑}fNC j(2.36)— вектор,j=1содержащий суммарные потоки вектора напряженности электрического поля через расположенные в вакууме грани эмиттирующих ячеек, I = {IC0 , . . . , INC 0 } —вектор токов, проходящих через эмиттирующие ячейки. Таким образом, мы можем представить алгоритм вычисления распределения плотности тока эмиссии,ограниченного пространственным зарядом:501. после задания начальных данных, необходимо поставить в соответствиекоэффициент Kjk каждой частице;2. в процессе расчета траекторий частиц, вычисляются элементы aij матрицыA используя уравнение (2.35);3. после решения уравнения Пуассона, вычисляются элементы fNC j вектораF;4. на следующей итерации используется приближение к плотности тока эмиссии, полученное из решения системы линейных уравнений (2.36)I = ε0 A−1 F.51Глава 3Комплекс программ для моделирования динамикипучков заряженных частиц в электростатическомприближении3.1 Основные характеристики комплекса программРазработанныйкомплекспрограммDAISI(DesignofAccelerators,optImizations and SImulations) полностью реализован с использованиемязыка программирования С++ [57].
При выборе дополнительных библиотеки программных решений, приоритет отдавался кроссплатформенному программному обеспечению (ПО) распространяемому бесплатно или по лицензиямGNU на свободное ПО. Для создания графического интерфейса пользователяиспользовался инструментарий разработки ПО Qt [88, 89]. Для визуализацииданных использовалась открытая библиотека Visualization Toolkit (VTK) [90],легко интегрируемая с виджетами Qt. Дополнительно, для решения отдельныхзадач использовались открытые библиотеки Boost [91] и Armadillo [92]. Дляраспараллеливания вычислений на персональном компьютере использовалсястандарт OpenMP [93, 94].Для компиляции использовался оптимизирующий компилятор Intel C++compiler (академическая лицензия СПбГУ).
При компиляции и проведении расчетов использовалась платформа конфигурации Win64, однако комплекс программ может быть скомпилирован для любой платформы в силу кроссплатформенности используемых инструментариев Qt, Boost, Armadillo, VTK, OpenMP.523.1.1 Формат входных данныхЕдинственными входными файлами, которые необходимы для проведениярасчетов, являются файлы, содержащие описание геометрии границ расчетнойобласти в специальном формате и файл конфигурации расчетной сетки. Например, для создания файлов геометрии можно использовать открытую интегрируемая платформу для численного моделирования Salome [95].
Подобный подходк конфигурации расчетов реализован также например в платформе численногомоделирования задач механики сплошных сред OpenFOAM [96]. Самостоятельное задание геометрии в рамках разработанного комплекса программ невозможно, однако доступна возможность визуализации заданной с помощью внешнихсредств геометрии.3.1.2 Варианты использованияКомплекс прикладных программ, предназначенный для моделирования динамики потоков заряженных частиц должен удовлетворять следующим требованиям:• Поддерживать возможность использования различных систем координатпри расчетах, таких как декартовые двумерная и трехмерная, цилиндрическая с учетом осевой симметрии и цилиндрическая трехмерная, полярная.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
