Математическое моделирование киральных волноведущих систем (1103709), страница 2
Текст из файла (страница 2)
Поле в волноводе определяется системойуравнений Максвелла, которая с учётом зависимости полей от времени в виде e − iω t ,записывается следующим образом:riω r 4π rrotH=−D+j,ccrr rot E = iω B,crdiv B = 0,rdiv D= 4πρ ,(1)r rгде E , B – это напряжённость электрического поля и индукция магнитного поляrr rrсоответственно, D = E + 4π P – вектор электрической индукции ( P – макроскопическийrr rrвектор поляризации среды), H = B − 4π M – вектор напряжённости магнитного поля ( M –rмакроскопический вектор намагниченности среды), а ρ , j – усреднённые плотности зарядаи тока свободных носителей зарядов.6Для замыкания системы уравнений (1) необходимо использовать дополнительныесоотношения, называемые материальными уравнениями среды.
Материальные уравнениякиральной среды имеют следующий вид:rrrD = ε E − iχ H ,rrrB = µ H + iχ E,(2)где ε, µ, χ – материальные константы. Параметр связи электрического и магнитного полейχ , называемый параметром киральности, пропорционален отношениюразмер частицы-элемента среды, а λ – длина волны. Приa, где a – линейныйλa→ 0 киральные свойства средыλисчезают. В искусственных киральных средах величину χ = Caудалось значительноλувеличить. В этом случае киральность уже не является малой поправкой и даже для волнСВЧ диапазона характеристики распространения электромагнитных волн в киральных инекиральных средах могут значительно отличаться.rrНаличие в D слагаемого, пропорционального H , означает, что ток, индуцируемыйпеременным магнитным полем в киральных элементах, образующих киральную среду,создаёт не только магнитный дипольный момент, но и электрический дипольный момент.Требование взаимности приводит к тому, что переменное электрическое поле индуцирует втаких элементах ток, который создает как электрический, так и магнитный дипольныеrrмоменты.
Это налагает ограничения на коэффициенты при H в первом уравнении и при Eво втором уравнении системы (2): они должны быть комплексно сопряжены. В дальнейшембудем рассматривать изотропные киральные среды, для которых ε , µ , χявляютсяскалярными величинами.Существует другая распространенная форма записи материальных уравнений:rrrD = ε c E + iξ c B,rr rH = iξ c E + B µ ,(3)cгде ε с = ε − χ2µ, µс = µ , ξ с = − χµ( ξ с называют киральным адмитансом). В такой формезаписи материальных уравнений ε c зависит от параметра киральности и в случае обращенияпараметра киральности в ноль переходит в диэлектрическую проницаемость для обычнойдиэлектрической среды.7В третьем разделе рассматриваются основные методы расчёта киральных волноведущихсистем.
Для систем со сложным анизотропным заполнением, в том числе и киральным,точные аналитические решения можно получить только для сильно ограниченного классаобъектов.Вбольшинствеслучаевнеобходимоприменятьчисленныеметодыикомпьютерное моделирование. В диссертации рассмотрено несколько методов, достаточнохорошо зарекомендовавших себя при моделировании волноведущих систем со сложныманизотропным заполнением. Для каждого метода дана краткая характеристика, а такжепример применения данного метода для расчёта конкретной задачи.1.
Метод конечных разностей. Применение метода конечных разностей к расчётуволноведущих систем с киральным заполнением проиллюстрировано на примерепрямоугольноговолновода,длякоторогопроанализированысвойствараспространяющихся мод и построены дисперсионные кривые для первых двухмод.2. Метод Бубнова-Галёркина. Его применение рассмотрено на примере расчётапостоянных распространения кирального прямоугольного волновода.3. Метод конечных элементов. Рассмотрено его применение для расчёта постоянныхраспространения электромагнитных волн в киральном волноводе с круглойгеометрией поперечного сечения, причём благодаря введению дополнительныхграничных условий на нормальные составляющие поля, удаётся избавиться отпоявления фиктивных решений, возникающих при применении метода конечныхэлементов к расчёту волноведущих систем, и создать весьма экономичныйалгоритм, гарантирующий хорошую точность вычислений.В четвёртом разделе рассмотрены аналитические методы, использующиеся для расчётакиральных волноводов:1.
Метод векторных цепей.2. Метод диадных функций Грина.Дано представление об основных положениях этих методов и их применении для решенияконкретных задач.Во второй главе строится эффективный алгоритм расчёта киральных волноведущихсистем. В первом разделе решается задача расчёта цилиндрического кирального волноводаметодом конечных элементов.
Поле в таком волноводе в отсутствие токов и зарядовописывается системой уравнений Максвелла8rr rotH = −ikD ,r r rotE = ikB,гдеr rH, E– вектора напряжённости магнитного и электрического поля,электрическая и магнитная индукции, k =(4)r rD, B–ω.cНа стенке волновода зададим однородное граничное условиеr rr E × n = 0 (n - нормаль к ∑ ) .∑(5)Материальные уравнения киральной среды запишем в наиболее общем виде:rrrD = a11 (ε , µ , χ ) E + a12 (ε , µ , χ ) H ,rrrB = a21 (ε , µ , χ ) E + a22 (ε , µ , χ ) H ,(6)где ε и µ – диэлектрическая и магнитная проницаемость соответственно, а χ – параметркиральности.
Матрица коэффициентов в (6) должна быть невырожденной, а самикоэффициенты должны удовлетворять специальным предельным соотношениям, чтобы приустремлении χ к нулю киральная среда вырождалась в диэлектрик или магнетик:lim a11 (ε , µ , χ ) = ε ,χ →0lim a12 (ε , µ , χ ) = 0,lim a21 (ε , µ , χ ) = 0,lim a22 (ε , µ , χ ) = µ .χ →0χ →0(7)χ →0Отметим, что материальные уравнения (6), на самом деле, описывают более общий класссред, называемых биизотропными. Однако, в случае когда коэффициенты a12 и a21комплексно сопряжены, а a11 и a22 некоторые скалярные величины, мы переходим кпредельному случаю киральной среды.В силу регулярности рассматриваемой волноводной системы вдоль оси z можно записать∂= i β , где β – постоянная распространения. Введём обозначение λ = -iβ, и тогда система∂zуравнений Максвелла в покомпонентной записи в цилиндрической системе координат будетиметь следующий вид:22ra 1 ∂ ( rEϕ ) ∂ 2 Er a21 1 ∂ ( rHϕ ) ∂ 2 H r −−−ir − i 11 2 i +2 2 k ∆ r 2 ∂ϕ∂rk ∆ r ∂ϕ∂r∂ϕ∂ϕ + λ Hϕ + ika11E r + ika12 H r = 0,ra ∂ 1 ∂( rEϕ ) ∂Er a21 ∂ 1 ∂( rH ϕ ) ∂H r iϕ − i 11 −++ − − i −k ∆ ∂r r ∂r∂ϕ k ∆ ∂r r ∂r∂ϕ − λ H r + ika11 Eϕ + ika12 H ϕ = 0,9(8)rir2a12 1 ∂ ( rEϕ ) ∂ 2 Er a22i− +i2 k ∆k ∆ r 2 ∂ϕ∂rϕ∂+ λ Eϕ − ika21 Er − ika22 H r = 0,riϕi 1 ∂ 2 ( rH ϕ ) ∂ 2 Hr− 22rϕ∂∂∂ϕ r a12 ∂ 1 ∂ (rEϕ ) ∂Er a22+ − + ik ∆ ∂r r ∂r∂ϕ k ∆ + ∂ 1 ∂( rHϕ ) ∂H r + − −∂r∂ϕ ∂r r (9)− λ Er − ika21 Eϕ − ika22 H ϕ = 0,где ∆ = a11a22 − a12a21 .
Система уравнений (8)-(9) дополняется граничным условием (5).Для решения системы уравнений (8)-(9) используется метод конечных элементов. Этотметод имеет такое неоспоримое преимущество, как широкая универсальность, что позволяетрассчитывать волноводы со сложной геометрией поперечного сечения. Кроме того, выборэтого метода объясняется простотой его физической интерпретации, а также ясностью ичёткостью численного алгоритма, что существенно облегчает программирование сложныхзадач математической физики. Введём в области Ω сетку ω = ω r × ωϕ :{ω r = ri = ihr ,i = 1, 2,...K ;ωϕ = ϕ j = jhϕ ;}hr > 0 ,j = 1, 2,..., M ; hϕ =2π ,Mразбив тем самым область Ω на подобласти Ωij. Осуществим триангуляцию области Ω,разделив каждый прямоугольник Ωij диагональю.
Основная идея метода конечных элементовзаключается в разложении искомой функции по системе выбранных базисных функций{N }. Базисную функцию Nij для конечного элемента, соответствующего узлу (i,j) введённойijнами сетки, выберем так, чтобы она была равна единице в узле (i,j) и нулю во всех остальныхузлах.
Выберем носитель базисной функции, как показано на рис. 1.10DijsuppN ijjiРис. 1. Носитель базисной функции в прямоугольных координатах (r, φ)Каждую из функций Nij можно выразить через «стандартную» функцию φ(s,t) вида1 − s,0 ≤ s ≤ 1,0 ≤ t ≤ s,0 ≤ s ≤ 1,s ≤ t ≤ 1,1 − t ,− 1 ≤ s ≤ 0,0 ≤ t ≤ s + 1,1 + s − t ,ϕ ( s, t ) = 1 + s,s ≤ t ≤ 0,− 1 ≤ s ≤ 0,1 + t ,− 1 ≤ s ≤ 0,− 1 ≤ t ≤ s,s − 1 ≤ t ≤ 0.0 ≤ s ≤ 1,1 − s + t ,Общий вид этой функции показан на рис.
2.1stРис. 2. Общий вид функции φ(s,t)Тогда базисные функции можно представить следующим образом: rϕN ij (r , ϕ ) = ϕ − i,− j.hhϕ rФункции такого вида часто называют функциями Куранта.11Согласно методу конечных элементов, решение системы уравнений (8)-(9) будем искать ввидеK MK MHϕ = ∑ ∑ H ijϕ N ij ( r, ϕ ),H r = ∑ ∑ H ijr N ij ( r , ϕ ),i =1 j =1K MEr = ∑ ∑i =1 j =1i =1 j =1K MEϕ = ∑ ∑Eijr N ij ( r ,ϕ ),i =1 j =1(10)Eijϕ N ij ( r, ϕ ),где Hijr,φ, Eijr,φ – коэффициенты, являющиеся ввиду выбора базисных функций значениямиискомых функций в узлах сетки.В результате подстановки разложения в исходную систему уравнений, приходим кобобщённой задаче на собственные значения для определения постоянной распространения.Алгоритм решения задачи на собственные значения (8)-(9) реализован на языке FORTRAN.Рис.
3. Результаты расчёта постоянной распространения цилиндрического киральноговолновода радиуса RНа рис. 3 приведено сравнение результатов расчёта постоянной распространенияцилиндрического волновода, полученных в диссертации, а также на основе другого подходав работе Сведина1. Обозначение HEmn использовано для моды, переходящей в TEmn волну,где n означает зависимость от ϕ в виде einϕ, а m обозначает, что мода имеет номер m, приупорядочивании по возрастанию частоты отсечки для данного n.
Видно, что результатырасчётов согласуются друг с другом.1J.A.M. Svedin. Propagation analysis of chirowaveguides // IEEE Trans. Microwave Theory Tech. – 1990. – Vol. 38. –№10. – P. 1488-1495.12Отметимнекоторыеособенностикиральноговолновода,выявленныеприегомоделировании:§Все моды являются гибридными или смешанными.§Каждая киральная мода разделяется на две ветви с различными постояннымираспространения, которые зависят от знака n, однако частоты отсечки для двухтаких ветвей одинаковые (явление бифуркации мод).§Существует «обратная область», где фазовые и групповые скорости имеютпротивоположный знак.Во втором разделе даны основные положения метода смешанных конечных элементов,описан разработанный алгоритм и представлены результаты решения задачи нахожденияпостоянных распространения прямоугольного кирального волновода методом смешанныхконечных элементов.Поле в волноводе в отсутствие токов и зарядов описывается системой уравненийМаксвелла (4) с граничным условием (5) и материальными уравнениями в виде (6).rrВ силу регулярности системы вдоль оси z, E ( x, y, z ) = E ( x, y )e γz , и систему уравнений(8)-(9) можно переписать следующим образомrrrrotγ H = −ika11 E − ika12 H ,(11)rrrrotγ E = ika21E + ika22 H ,(12)r ∂Xr ∂Xгде rotγ X = z − γ X y ex + γ X x − z∂x ∂y r ∂X y ∂X x e y + ∂x − ∂yr ez .rrИз уравнения (12) выразим H , подставим выражение для H в уравнение (11) и запишемполученное равенство в видеrotγrrr1a raa a rrotγ E − ikrotγ 21 E + ik 12 rotγ E − k 2a11E + k 2 12 21 E = 0.a22a22a22a22(13)r~rУмножим (13) скалярно справа на E * из того же пространства, что и E , проинтегрируем повсей области и учтём свойства скалярного произведения, дивергенции, ротора и граничныеусловия.















