Аналитическая техника при решении трехмерных задач электронной оптики Текст научной статьи по специальности «Физика»

ISSN 0868-5886

НАУЧНОЕ ПРИБОРОСТРОЕНИЕ, 2014, том 24, № 1, c. 96-103

РАБОТЫ ШКОЛЫ ПРОФ. Ю.К. ГОЛИКОВА: -

РАБОТЫ, ПОСВЯЩЕННЫЕ ПАМЯТИ Ю.К. ГОЛИКОВА

УДК 621.38; 517 © В. Я. Иванов

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

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

Кл. сл.: электронная оптика, расчеты полей, метод граничных элементов

ВВЕДЕНИЕ

Метод граничных элементов (МГЭ) является одним из наиболее мощных и гибких методов решения смешанных краевых задач математической физики для областей сложной формы [1-2]. В данной работе нас будут интересовать лишь задачи эллиптического типа (Лапласа и Пуассона), к которым сводятся расчеты электрических и магнитных полей, которые входят в правые части уравнений движения заряженных частиц. В таком случае интегральное представление потенциала ф электрического поля дается формулой

(р(г) = >1)^ +|у(г2) G(r, г2^2 +

S1 S2 СП

+ | G (г, гз)р(г^,

V

Г е S1,г2 е S2,г3 е V.

Здесь первый член описывает потенциал поверхностных зарядов о на поверхности Sl, второй — поверхностных диполей V на поверхности S2, а третий — поле заряда плотностью р пучка частиц в объеме V; G — функция Грина точечного источника поля, а п — вектор нормали к поверхности в точке источника. Структура интегрального представления для магнитного скалярного и векторного потенциалов магнитного поля аналогичны представленной выше, поэтому для иллюстрации развитой нами техники мы рассмотрим лишь электростатический случай и опустим дипольный член. Вид функции Грина для точечного заряда зависит от системы координат. В декартовых координатах в трехмерном случае эта функция имеет вид G(r, г') = !/(2яе0К), где R — расстояние между

точкой наблюдения г и точкой источника г', е0 — диэлектическая проницаемость вакуума.

Вычисление производных потенциала любого порядка осуществляется аналитическим дифференцированием функции Грина, что является несомненным достоинством метода интегральных представлений. Удовлетворяя заданным краевым условиям F(S) на поверхности границы S = S1 и S2, получим интегральные уравнения

Фредгольма СХ = F относительно плотности источников поля X = а. При численном решении требуется аккуратное выделение особенностей ядра С при г ^ г' и степенной особенности в окрестности нерегулярных точек границы. Техника выделения этих особенностей достаточно хорошо развита [3]. Один из возможных алгоритмов решения уравнения Фредгольма — метод интерполяций и коллокаций состоит в следующем. Представим границу области параметрическими уравнениями

N

S = и Sk, х = хк (т,Ц), У = Ук (т,ц),

к =1

г = ^ (т,ц), dS = Jk (т,ц)dтdц,

а для плотности источников предложим билинейную интерполяцию вида

а(т,ц) =

= {[а (т ~ Т-1 ) + а.'-1,1 (Т - т)] Ц- Ц-1 ) + + [а, 1-1(т - т,-1)+а -1,1 -1 (т - т)] (ц- -п)\1КК,,

где а = а(т,ц), те[т-1,т], це[ц_ьц], Нт = т- т-1, Нц = ц - ц--1. После дискретизации получим сис-

тему линеиных уравнении с матричными элементами

а =

^ £ Я

¥ш (т,П) г - г' (т,))

Здесь ут — моменты билинеИноИ формы в аппроксимации источника поля. Использование квадратурнои формулы

N N

И F = (т )

¿=1 ] =1

порядка N с весами Щ Щ- и узлами тг-, позволит

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

Сложность вычислительных процедур при решении задач электронной оптики заключается в том, что система линейных уравнений GX=F решается один раз, в то время как вычисление характеристик поля требуется в каждой точке интегрирования уравнений движения заряженных частиц. При этом объем вычислений значительно возрастает при выделении особенностей ядра а, а это имеет место не только в окрестности границы но и объеме V, занятом зарядами плотности р. Ситуация значительно упрощается, если использовать аналитику при интегрировании по поверхностным и объемным элементам

ВЫЧИСЛЕНИЕ ИНТЕГРАЛОВ ПО ПОВЕРХНОСТНЫМ ЭЛЕМЕНТАМ

В работах Е. Окона и Р. Харрингтона [4-6] развита аналитическая техника вычисления интегралов типа потенциала простого слоя для кусочно-постоянного, линейного и квадратичного распределений плотности источника на поверхности треугольника. В случае линейного распределения а = s1т + s2) формула для потенциала имеет вид

р( х, у, z) =

1

4яеп

<Я

^т + s2)) dт d)

1

2

4же

(V + V2)

/ (Л) = (с + СЛ)(2аС - Ьс + ЪСЛ) + Ье2 ¿Л) -

с + СЛ„ к2 - е2(Ь2 + С 2)и. ... X (а, Ь; Л) +-, Р(а, Ь; Л) +

+

2е 2С V Ь2 + С2

(Ь' 2 + С' 2)СЛ + С (а Ъ ' + сС) - Ьк'

2С (Ъ '2 + С2)

X (а ', Ъ' ; Л) +

ек,

+ <(а, Ъ;Л) - 0(а', Ъ ;Л)] +

+

-[(ЪЪ' + С 2)к' 2 +

4С 2(Ъ '2 + С 2)3/2 + е2(Ъ '2 + С 2)(ЪЪ ' + С2) -2кк' (Ъ' 2 + С 2)]х х Р(а ', Ъ '; Л),

Р(а,Ъ; Л) = 21п [(Ъ2 + С 2)Л + аЪ + сС +

+^[(Ъ2 + С 2)Л+аЪ + сС]2 + (аС - Ъс)2 + е2(Ъ2 + С2)

2(а,Ъ;Л) = Parccos1еС| х

1

(а + ЪЛ)2 + (с + С Л)2

[(с + СЛ)2 + е2][(аС - Ъс)2 + е2(Ъ2 + С2)] [

а + ЪЛ + X (а, Ъ;Л)

¿(Л) = 1п — ,

а' + Ъ 'Л + X (а ', Ъ' ;Л)

X(а,Ъ;Л) = ^(а + ЪЛ)2 + (с + СЛ)2 + е2,

1 - г3|, р = -sign[(с + СЛ)к - Ъе2],

а = г - г

а' = (г1 -г3)(г3 -г)/а2, а = 1 + а',

Ъ' = (г1 -г3)(г2 - г)/а2, Ъ = Ъ' -1, с = п[(г - г3) х (г3 - г)]/а2, С = |(г1 -г3) х (г2 -г3)|/а2, е = |п(г3 - г)/ С, к = аС - Ъс, к' = а 'С - Ъ 'с,

(г - г3) х (г -1(г1 - г3) х (г2 - г3>|

, А =1 а2С. 2

где т, г) — локальные координаты треугольника с вершинами г1, г2, г3. Интеграл VI вычисляется по формуле V = 2А[ / (1) - / (0)]/а. Здесь

Аналогичные выражения можно написать для интеграла Выведенные авторами формулы вполне регулярны. К сожалению, в этих работах отсутствуют выражения для градиентов потенциала. Напряженность электрического поля Е = -Ур будем вычислять прямым дифференцированием формулы для потенциала. Эта формула содержит особенности в окрестностях вершин и ребер тре-

1

х

п

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

Уб = -

Р(А - В)

где

S = (а + ЬХ)2 + (с + dХ)2 + е2, Z = [(с + dХ)2 + е2]Г, Y = [(а + ЬХ)2 + е2(Ь2 + d2)],

А = d Уел/З + 4е [(а + ЬХ)Уа + (с + dХ)Vc + еУе],

В =

de^JS

Z

Л

{[(с + dХ)Vc + еУе^ + [(с + dХ)2 + е2 ] х

Уб = -

рУа

с + dХ

2 + е2^

2 + е2

Уб = -

рУс

(

4а2 + е2 ^

4 а2 + е2

х [(аd - bc)(d Уа - ЬУс) + еУе + (Ь2 + d2)]}.

Существуют два случая, когда выражения для У б приводят к неопределенности вида 0/0. Первый случай соответствует а ^ 0, т. е. г ^ г3 . Разлагая по малому параметру, имеем

Во втором случае с + dХ^ 0 отвечает ситуации, когда вектор г3 - г коллинеарен вектору нормали в углу г3 . В этом случае получим формулу

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

ВЫЧИСЛЕНИЕ ИНТЕГРАЛОВ ПО ОБЪЕМНЫМ ЭЛЕМЕНТАМ

Рассмотрим параллелепипедальную неравномерную сетку в декартовой системе координат

{хг} х {у- } х {гк } и предложим трилинейную аппроксимацию для плотности объемного заряда

Р(Х У, г) = {[( Р+и,к (Х - X ) + Рг,-,к (Х+1 - Х) ) Х

х(У-+1 - У) + (Рг+1,1,к (Х - хг) + Рг, 1,к (хг+1 - х)) Х

х(У - У1)] (гк+1 - г) + [(Рг+1,1 ,к+1(Х - хг) + +Рг, 1 ,к+1(Хг+1 - Х)) (У]+1 - У) + (Рц+1(Х - Хг) + +Рг, 1,к+1(Х+1 - Х)) (У - У1)] (г - гк )} х

х{(Хг+1 - Хг)(У1+1 - У1)(гк+1 - гк )} 1.

Тогда при вычислении потенциала объемного заряда

р( Х, У, г) =

Х1 +1 У]+ гк+1

= Х | | | Р(х', у', г')С(Х, У, г; х', у', z')dxdydz

г,1 ,к х, У, г,

нам понадобится вычислить 4 интеграла от моментов плотности заряда:

•/1=Ы К •/2=Ы К ^у^,

Jз = Ш ^ /1 = Ш

где К2 = (х - х')2 + (у - у')2 + (г - г')2. Остальные интегралы, необходимые для вычисления потенциала, получаются циклической заменой координат Х, У, г.

После замены переменных Х = х' - х, у = у' - у, г = г' - г и введения обозначения г2 = х2 + у2 + г2 необходимые нам интегралы вычисляются аналитически

3

/1 = У [гг + (х2 + У2) 1п| г + г| ] + 1п| у + г| +

+

2

г 1п| у + > + У 1п| г + г| - г + х tg 11 — I -

- Х «1

V хг

+

+

]_

36

6х2г - 2г3 + 6х3

«Ч ^ ]-«Ч г"

V хг

3 У ( У2 + 3х21п| г + г|)

е

е

е

а

е

2

х

х

„3

72= 12+24 (х2 + у2)[ .г+(х 2+у 2) 1п12 + г1 ]'

3

л=

3 15

= ху 1п|2 + г\ + у. 1п|X + г \ + 2Х 1п|у + г\ -

х^-1 ^ + у ^ -1 ^ + г ^ -1 хУ хг уг гг

дх dy 1 ху | | | |

-= -2 tg--+ у 1п х + г + х1п у + г ,

г .г

= Я

72. =Ц «Ыу = 1 [ уг + (х2 + . 2)1п| у + г| ],

7з. = Я ^ = ^ • ^4. =1/

ху d х dy гъ

г

3

75.=я ^=2

^..=Л

хугдх dy . 3/2

г ~ 3г .

Подобная техника для вычисления интегралов по ячейкам двумерной прямоугольной сетки в декартовой и цилиндрической системах координат при билинейной аппроксимации объемного заряда была разработана В. Ивановым [8] в 1989 г. Трехмерный случай для кусочно-постоянной аппроксимации рассмотрен в работе Дж. Куанга и др. [9]. Там используется только интеграл который выписан с ошибкой, к тому же в комплексной форме, хотя статический потенциал является вещественной функцией.

Для получения компонент электрического поля следует дифференцировать вышеприведенное представление потенциала. При этом понадобится вычисление новых интегралов. Например, для компоненты Е. нам понадобятся следующие интегралы:

Остальные компоненты напряженности поля вычисляются циклической заменой переменных х, у, Обратим особое внимание на проблему сингу-лярностей. Заряженный кубик имеет особенности в окрестностях вершин и ребер. На рис. 1 показано, что при кусочно-линейной аппроксимации объемного заряда в точке "Ь", лежащей на общих гранях соседних двух кубиков, эти особенности имеют искусственный характер, поскольку они взаимно компенсируются вкладами соседних ячеек сетки, в то время как в окрестности точки "а" — это реальные особенности. На самом деле плотность заряда на границе пучка плавно падает до нуля, поэтому при линейной аппроксимации особенности градиента потенциала обнуляются умножением на нулевое значение плотности заряда, соответствующего граничным узлам сетки. Зато при кусочно-постоянной аппроксимации все эти особенности проявляются не только в окрестности граничных ребер и вершин сетки, но и на внутренних границах ячеек сетки, поскольку плотность заряда на этих границах меняется скачком и компенсации особенностей не происходит. По этой причине кусочно-постоянная аппроксимация заряда в принципе неприемлема для задач электронной оптики, поскольку она дает высокую погрешность расчета компонент поля. Детальный анализ погрешностей для кусочно-постоянной и линейной аппроксимаций дан в работе [10].

Я И

ь

У

0 Ь a

Рис. 1. Ячейки сетки (слева) и алгоритмы аппроксимации объемного заряда (справа): кусочно-постоянный, кусочно-линейный и квадратичный

В общем случае для пучков высокой компрессии аппроксимация объемного заряда в ячейках прямоугольной формы может давать значительную погрешность. Развитую нами методику легко обобщить на случай криволинейной сетки, согласованной с границей пучка, что обеспечит гораздо более высокую точность расчета, если нам удастся построить непрерывное взаимно-однозначное отображение Т криволинейных координат (и, V, м) в декартовы (х, у, г). В этом случае интегрирование по криволинейным ячейкам сведется к рассмотренному выше случаю параллелепидальных ячеек с заменой Р(и, V, м^и dv dw на р(х, у, г)/dx dy дг, где / — якобиан преобразования координат Т. В большинстве практически важных случаев построение такого преобразования не составляет особой сложности. Для этого с поверхности эмиттера выпускается ряд траекторий, образующих контур пучка. Если интегрировать их с постоянным шагом по времени, то сами траектории образуют множество линий сетки в продольном направлении, а совокупность их координат для фиксированного момента времени дают возможность построения двух множеств линий сетки в поперечных направлениях. Таким образом можно построить кусочно-линейное представление криволинейной сетки или более гладкие представления, если использовать сплайны нужного порядка.

ПРИМЕР РЕШЕНИЯ ПРАКТИЧЕСКОЙ ЗАДАЧИ

В качестве практической задачи рассмотрим расчет пушки перспективного 80-мегаваттного клистрона с плоским пучком [11]. Параметры пушки: ток 250 А, ускоряющее напряжение 415 кВ, рабочая частота клистрона 91 ГГц. В традиционных клистронах с цилиндрическим пучком влияние объемного заряда создает серьезные препятствия для построения приборов с мощностью свыше 75 мегаватт, поэтому более мощные клистроны могут быть построены на новых принципах: многоструйные или с плоским пучком. В многоструйных клистронах при сведении отдельных пучков в единый жгут трудно получить пучок с малым поперечным фазовым объемом, поэтому получение ламинарных пучков имеет перспективы в схеме с плоским пучком. Геометрия пушки такого клистрона показана на рис. 2. Значительный ток пучка и высокая степень релятивизма приводят к задаче самосогласованного поля, в которой требуется учитывать кулоновское поле объемного заряда и собственное магнитное поле пучка. Полная система уравнений в этом случае включает уравнения поля (уравнение Пуассона), релятивистское уравнение движения частиц и уравнение непрерывности токов и зарядов:

Рис. 2. Геометрия пушки клистрона.

1 — катод; 2 — фокусирующий электрод; 3 — анод

Cathode

ssseff r • «fiS^S! ÎSgl*-k^MtL :

2

Рис. 3. Результаты траекторного анализа пушки клистрона с плоским пучком

s

p = e ^ E + -[v x B]j, E = -Vp, div pv = 0, p = my v, y = -

л/1 - (V / с)2 '

где е — заряд частицы, т — ее масса покоя, V — скорость, р — импульс, у — гамма-фактор, Е — напряженность электрического поля, В — индукция магнитного поля, с — скорость света в вакууме.

Эти уравнения должны быть дополнены краевыми условиями для поля, начальными условиями для координат, скоростей частиц и плотности тока эмиссии. Собственное магнитное поле пучка вычисляется интегрированием плотности тока почка ] = pv по занимаемому им объему V

B( ) 4 r[j(r') x r']

B(r ) = 4пИ-о p-R-1

r )xi ,

3 dV, r ' eV. R3

Здесь ^о — магнитная проницаемость вакуума. Поскольку полученная система уравнений нелинейна, она решается итерационным методом с релаксацией зарядов рп = юпрп-1 + (1 - оп)рп и токов эмиссии ]п = (оп]п- + (1 - соп)]п. Здесь верхний индекс соответствует номеру итерации, а величины со значком ... отвечают зарядам и токам, вычисляемым на текущей итерации п с помощью специальных процедур метода "трубок тока" для плотности заряда в ячейках сетки, закону Ленгмюра для термоэмиссии или Фаулера—Нордгейма для автоэмиссии. Правильный подбор коэффициентов

Рис. 4. Конструкция пушки, изготовленной компанией Calabasas Creek Inc.

Рис. 5. Прототип клистрона, разработанного отделом клистронов Стенфордского центра линейных ускорителей (SLAC)

релаксации 0 < rnn <1 позволяет значительно сократить необходимое число итераций до сходимости процесса [11].

Результаты расчета представлены на рис. 3. Здесь поверхности электродов аппроксимированы набором треугольников, измельчающихся в областях с малыми радиусами кривизны. Для пучка с током 250 А потребовалось 10 итераций, чтобы получить сходимость по току с точностью 0.5 %. На рис. 4 представлена пушка клистрона, а на рис. 5 дан общий вид клистрона с плоским пучком, впервые изготовленного и протестированного в SLAC [12].

СПИСОК ЛИТЕРАТУРЫ

1. Бреббия К., Теллес Ж., Вроубел Л. Методы граничных элементов. Пер. с англ. М.: Мир, 1987. 524 с.

2. Бенерджи П., Баттерфилд Р. Метод граничных элементов в прикладных науках. Пер. с англ. М.: Мир, 1984. 494 с.

3. Иванов В.Я. Методы автоматизированного проектирования приборов электроники. Ч. 1: Методы расчета физических полей. Новосибирск: Изд-во ин-та математики, 1986. 194 с.

4. Okon E.E., Harrington R.F. The potential due to a uniform source distribution over a triangular domain // Int. J. Numer. Methods Eng. 1982. V. 18, N 9. P. 14011411.

5. Okon E.E., Harrington R.F. The potential integral for a linear distribution over a triangular domain // Int. J. Numer. Methods Eng. 1982. V. 18, N 12. P. 18211828.

6. Okon E.E. Potential integrals associated with quadratic distribution in a triangular domain // Int. J. Numer. Methods Eng. 1982. V. 21, N 2. 197-209.

7. Ivanov V., Krasnykh A. Analytical computations in solution of 3D problems of physical electronics by BEM // ACES-2002, March 18-22, 2002, Monterey, CA.

8. Иванов В.Я. Метод анализа трехмерных нестационарных потоков заряженных частиц // "Численный анализ". Тр. Ин-та математики СОАН СССР, Т. 15. Новосибирск: Наука, 1989. C. 172-187.

9. Qiang J., Lidia S., Ryne R.D., Limborg-Deprey C. Three-dimensional quasistatic model for high brightness beam dynamic simulation // Phys. Rev. Special Topics - Accel. & Beams. 2007. V. 10. 129901(E) (2 pages).

10. Ivanov V. Green's function technique in forming of intensive beams // Int. J. of Modern Physics A. 2009. V. 24, N 5. P. 869-878.

11. Иванов В.Я. Методы автоматизированного проектирования приборов электроники. Ч. 2: Методы решения задач электронной оптики. Новосибирск: Изд-во ин-та математики, 1986. 197 с.

12. Read M.E., Miram G., Ives R.L., Krasnykh A., Ivanov V. An electron gun for a sheet beam klystron // Proc. of IEEE 29th International Conference on Plasma Physics, Piscataway, New Jersey. IEEE, 2002. P. 190.

Институт вычислительных технологий СО РАН, г. Новосибирск

Контакты: Иванов Валентин Яковлевич, vivanov.48@mai1.ru

Материал поступил в редакцию 7.01.2014

ANALYTICAL TECHNIQUE FOR THREE-DIMENSIONAL PROBLEMS OF ELECTRON OPTICS

V. Ya. Ivanov

Institute of Computational Technologies, Siberian Branch of RAS, Novosibirsk, RF

The analytical technique for three dimensional electrostatic fields of complex-shape electrode system and space charge of electron beams is presented in the paper. The singularity analysis of integral representations for the potential and field gradients is provides. The practical problem solution for an electron gun of power sheet-beam klystron is discussed.

Keywords: electron optics, field computations, boundary element method

REFERENСES

1. Brebbiya K., Telles Zh., Vroubel L. Metody granich-nykh elementov. Per. s angl. M.: Mir, 1987. 524 s. (in Russian).

2. Benerdzhi P., BatterfildR. Metod granichnykh elementov v prikladnykh naukakh. Per. s angl. M.: Mir, 1984. 494 s. (in Russian).

3. Ivanov V.Ya. Metody avtomatizirovannogo proyek-tirovaniya priborov elektroniki. Ch. 1: Metody rascheta fizicheskikh poley. Novosibirsk: Izd-vo in-ta matema-tiki, 1986. 194 s. (in Russian).

4. Okon E.E., Harrington R.F. The potential due to a uniform source distribution over a triangular domain // Int. J. Numer. Methods Eng. 1982. V. 18, N 9. P. 1401-1411.

5. Okon E.E., Harrington R.F. The potential integral for a linear distribution over a triangular domain // Int. J. Nu-mer. Methods Eng. 1982. V. 18, N 12. P. 1821-1828.

6. Okon E.E. Potential integrals associated with quadratic distribution in a triangular domain // Int. J. Numer. Methods Eng. 1982. V. 21, N 2. 197-209.

7. Ivanov V., Krasnykh A. Analytical computations in solu-

tion of 3D problems of physical electronics by BEM // ACES-2002, March 18-22, 2002, Monterey, CA.

8. Ivanov V.Ya. Metod analiza trekhmernykh nestatsio-narnykh potokov zaryazhennykh chastits // "Chislen-nyy analiz". Tr. In-ta matematiki SOAN SSSR, T. 15. Novosibirsk: Nauka, 1989. C. 172-187. (in Russian).

9. Qiang J., Lidia S., Ryne R.D., Limborg-Deprey C. Three-dimensional quasistatic model for high brightness beam dynamic simulation // Phys. Rev. Special Topics - Accel. & Beams. 2007. V. 10. 129901(E) (2 pages).

10. Ivanov V. Green's function technique in forming of intensive beams // Int. J. of Modern Physics A. 2009. V. 24, N 5. P. 869-878.

11. Ivanov V.Ya. Metody avtomatizirovannogo proyektiro-vaniya priborov elektroniki. Ch. 2: Metody resheniya zadach elektronnoy optiki. Novosibirsk: Izd-vo in-ta matematiki, 1986. 197 s. (in Russian).

12. Read M.E., Miram G., Ives R.L., Krasnykh A., Ivanov V. An electron gun for a sheet beam klystron // Proc. of IEEE 29th International Conference on Plasma Physics, Piscataway, New Jersey. IEEE, 2002. P. 190.