УДК 621.372.2
ПОСТРОЕНИЕ КОНЕЧНО-ЭЛЕМЕНТНЫХ СХЕМ ВЫСОКОЙ ТОЧНОСТИ ДЛЯ ЦИЛИНДРИЧЕСКИХ волноводов
A. JI. Делицын, И. В. Степанов
(кафедра математики) E-mail: [email protected]; [email protected]
Разработан метод вычисления мод градиентных диэлектрических волноводов. Результаты проведенного тестирования для полого и двухслойного волноводов свидетельствуют о высокой точности метода. Предложенный метод позволяет решить проблему аппроксимации поля в нуле в цилиндрической системе координат.
Известно, что при применении метода конечных элементов к задаче расчета диэлектрического волновода в цилиндрической системе координат возникают две проблемы: аппроксимация поля в начале координат и аппроксимация нулевого собственного значения бесконечной кратности [1]. Одним из способов, позволяющих решить первую из них, может служить введение фиктивной дополнительной границы в окрестности начала цилиндрической системы координат, соответствующей идеально проводящему стержню малого радиуса. При стремлении этого радиуса к нулю собственные значения полученной задачи будут стремиться к собственным значениям исходной задачи. В случае задачи о собственных частотах мембраны это доказано в работе A.A. Самарского [2], что и положено в основу применяемого в работе метода.
Рассмотрим цилиндрический волновод с круговым поперечным сечением О единичного радиуса. В некоторую точку на его оси поместим цилиндрическую систему координат, ось Oz которой направим по оси цилиндра. Пусть волновод заполнен веществом с характеристикой:
e(r,cp,z) =e(r), ß(r,(p,z) = l,
е кусочно-непрерывна в О, стенки идеально проводящие.
Будем искать решения уравнений Максвелла в виде нормальных волн /.'. II ~ e-^t+%jz+%nnp _ дос. ле сокращения на экспоненциальный множитель приходим к задаче на отрезке [0,1], роль собственного значения играет параметр -у. Следуя [3], выберем из системы восьми уравнений Максвелла шесть основных для шести неизвестных следующим образом:
1 d Я
г dip
■ ijll r, + ikeEk = О,
i^Hj,
ld_ . г dr
d
-Я,
dr (re Er)
ike /'.'r, = 0,
(1)
ld_
r dip
(eEv) + iyeEz = 0,
1 d_
r dip
Ez — ijEv — //,://,. = 0,
d
ijEr —ikll,, = 0, dr
(2)
1 d . __ . 1 d „ -—(rHr) + -—Hv Kr dr r dip
ijHz= 0.
Введем обозначения X = (Нг, ЕХ)Т =
= (н±,Ег)т, у = (с/-;,..с/-;,.//,)г = (?/<; .//,)г.
Подставляя У из (2) в (1), приходим к уравнению на собственные значения:
/did г—---—г + к
dr г dr
. 1 d im——г г dr
kme
7 .dl
"er imr—--
dr r
-.2
u2
m
r
kern d
\
er
-iker
dr
d d d e о
гк—re —те—---m
dr dr dr r
X =
= 7 2r
(3)
где X принадлежит множеству векторов из С°°[0,1]. Уравнение (3) дополним граничными условиями
Яг(1)=0 и £г(1) = 0, |Х(0)|<оо (4)
и неиспользованным выше уравнением Максвелла
1 d 1
-//.:? /'.', =--—(rHv)--inHr
г dr г
(5)
Введем операторы
gradip(r) = ( ^-ip, %-^-ip
t \ i im rot ip(r) = ( —<p,
d_
dr
4>
T
T
1 (1 1171
(Ну II = ——/•//, Н--//.,.
г г
Распределение собственных значений в комплексной плоскости и их асимптотики получены в [3].
Выпишем вариационный функционал для нашей задачи. Для этого умножим (3) слева на произволь-
—Т
ный вектор Y и проинтегрируем по г в пределах от 0 до 1. Тогда получим 1 1
г dr div II div II + er dr (grad Ez, grad Ez)
dr^k2er (HrHr + HVHV) - kemHrEz
knie /'.'-. //,. — ikerHv^-Ez — ikerHv^-Ez\ = dr dr >
d_
' dr
= -72 r dr (HrHr + //r,//r, + eEzEz) . (6)
о
Поскольку задача (3)-(4) имеет особенность в нуле, ее численная реализация методом конечных элементов представляет определенные трудности. Заменим ее на задачу, в которой в окрестности нуля расположим проводящий цилиндр малого радиуса а, что соответствует с математической точки зрения рассмотрению задачи 1 1
г dr div // div // + er dr (grad Ez, grad Ez)
dr$yk2er (ffrffr + HVHV) - kemHrEz
knie /'.'-. //,. — ikerHv^-Ez — ikerHv^-Ez \ = dr dr >
d_
' dr
= —г2 j г dr (НГНГ + //г,//г, + ?/•;./%) (7)
а
на отрезке [а, 1] и постановке граничных условий вида
Яг(1)=0 и £г(1) = 0, |А"(о)|<оо, (8)
X принадлежит множеству векторов из С°°[а, 1], удовлетворяющих граничным условиям, и неиспользованному выше уравнению Максвелла
\ с1 1
—/А:?/','- = - — (гЯ„)--гтНг. (9)
г ¿г * г
Для упрощения вычислений производятся замены гНТ = II,.. ¡11 г, = Нк.
Решение будем искать методом смешанных конечных элементов [4]. Введем на отрезке г е [а, 1] сетку {ж.?}: 3 = 0, • • •, IV, хо = а, хм = 1, N — количество конечных элементов.
Разложим поля по базисным функциям с конечными носителями:
Щ
/ = ^ ^ стфт1 ' 1 ?
п=0
Щ — количество базисных функций, зависящее от количества конечных элементов и их порядка,
Cin — неизвестные коэффициенты. При этом для IIг, берутся разрывные элементы первого порядка
ф2п = <
О,
,0,
= <
о,
X ^ X j >
X j ^^ X ^^ X j > X ^ Xj+1-X ^ X j >
Ж ^ Xj+1-
Для и Ех берутся непрерывные элементы второго порядка
0,
<РЗ,4п = <
4(ж — xj)(x — Xj+1)
(Xj+1-Xj)2
0,
УЗ,4(п+1)
0,
2(х — xj)(x — 0.5 (xj
X ^ X j >
Ж ^ ^^ Ж ^^ Ж ^ > Ж ^
Ж ^ Ж j > •gj+l))
= <
(Sj + l-Sj)2
X j ^^ Ж ^^ X j >
2(ж — ®j+i)(a; — 0.5(®j + Xj+1))
(®i+i - xj)2
Xj+1 Ж
0,
При этом в точках xq я хм базисные функции должны удовлетворять граничным условиям исходной задачи.
Подставляя полученные разложения в вариационный функционал и производя интегрирование, получим задачу на собственные значения: А(к,т) = 2В, где А(к,т) = а(к,т)\, В = b\, I 0..... Л;| + Л;2 + /V;s- — квадратные разряженные матрицы.
Программа составлена для расчета цилиндрических волноводов с произвольным куеочно-непре-рывным заполнением диэлектриком. Тестирование проводилось для случаев полого волновода и двухслойного волновода, допускающих точное решение.
Таблица 1
Номер моды Номер 7 7
ТМ-волна ТЕ-волна
0 1 2 3 2.4048256 (2.4048261) 5.5200791 (5.5200780) 8.6537380 (8.6537277) 3.8317064 (3.8317061) 7.0155942 (7.0155858) 10.173511 (10.1734682)
1 1 2 3 3.8317062 (3.8317061) 7.0155895 (7.0155858) 10.173489 (10.1734682) 1.8409158 (1.8411838) 5.3298359 (5.3314429) 8.5322647 (8.5363167)
2 1 2 3 5.1356226 (5.13562231) 8.4172498 (8.41724411) 11.619875 (11.6198408) 3.0542370 (3.05423697) 6.7061354 (6.70613322) 9.9694866 (9.96946724)
6 ВМУ, физика, астрономия, №5
Таблица 2
Точные собственные значения Радиус внутреннего проводящего стержня а
а = 1(Г2 а = 1(Г3 а = 1(Г4 а = 10-5 а = 10"® а = 1(Г7
2.4048261 2.4052714 2.4048302 2.4048258 2.4048257 2.4048257 2.4048257
5.5200780 5.5224651 5.5201168 5.5200932 5.5200930 5.5200930 5.5200930
3.8317061 3.8328913 3.8317252 3.8317135 3.8317134 3.8317134 3.8317134
7.0155858 7.0195663 7.0157455 7.0157080 7.0157080 7.0157079 7.0157081
В случае полого волновода система волн представляет собой суперпозицию полей магнитного и электрического типа. Решение можно выписать через вектор Герца П = CJm(\r) cosrmp или П = С Jm (Ar) sin rmp.
Тогда уравнение на собственные значения для электрических волн имеет вид Пе = 0 при г = 1, т.е. Jm(А) =0, Хтп — п-й положительный корень. Для магнитных: = 0 при г = 1, т. е. J'm(А) = 0, Атп — п-я положительный корень.
В табл. 1 приведены вычисленные собственные значения и точные собственные значения (указаны в скобках) при 40 конечных элементах и радиусе внутреннего стержня а = 1СГ6.
В табл. 2 приведена зависимость собственных значений нулевой моды от внутреннего радиуса а при 20 конечных элементах.
Полученные результаты позволяют утверждать,
Ref. 121086420-2-4-Imy-
0 1 2,3 4 5
к
Рис. 1
Rey. 543210-1-2-3-Im у-
0 0.5 1 1.5 2
к
Рис. 2
что применение конечных элементов высокого порядка позволяет аппроксимировать поле в нуле и приводит к высокой точности на достаточно грубой сетке.
Рассмотрим теперь случай двухслойного волновода с оеееимметричным заполнением радиуса Ъ и диэлектрической проницаемостью е. Моды круглых слоистых волноводов обладают рядом особенностей, отсутствующих у волн однородно заполненных волноводов и волн с осевой симметрией. В частности, им присуще, например, явление аномальной дисперсии (d-y/dk) < 0. При определенных параметрах волновода эти волны даже при отсутствии диссипации энергии преобразуются в волны с комплексными постоянными распространения — комплексные волны.
Решение будем искать методом смешанных конечных элементов описанным выше образом. Для расчета берется 20 конечных элементов и радиус внутреннего стержня а = Ю-6.
Ниже приведены дисперсионные характеристики низших типов гибридных волн HЕп и /-.'//| | волновода при е = 10, 6 = 0.2 (рис. 1), 6 = 0.6 (рис. 2). Непрерывными линиями изображены зависимости Re7 от к, точками — зависимости I1117 от к. Результаты расчетов с графической точностью совпадают с точными значениями, приведенными в [5].
Как видно из рис. 1, в данном случае постоянные распространения волн низших типов могут быть либо чисто мнимыми, либо чисто действительными величинами, причем обе волны являются нормальными во всем частотном диапазоне.
Как видно из рис. 2, в данном случае присутствуют как нормальные, так и аномальные волны. При этом возникают комплексные значения 7, имеющие как вещественную, так и комплексную части.
Литература
1. Боголюбов А.Н., Делицын A.J1. // Вестн. Моск. ун-та. Физ. Астрон. 1996. № 1. С. 9 (Moscow University Phys. Bull. 1996. N 1. P. 7).
2. Самарский АЛ. 11 Докл. АН СССР. 1948. № 1. С. 631.
3. Боголюбов А.Н., Делицын A.J1., Малых М.Д. // ЖВМ и МФ. 2001. 41, №1. С. 126.
4. Siarlit P.J., Lions J.L. Handbook on numerical analysis. V. 2. Finite element method. North Holland, Amsterdam, 1991.
5. Веселовский Г.И., Раевский С.Б. Слоистые металло-диэлек-трические волноводы. М., 1988.
Поступила в редакцию 11.10.04