Интегральные уравнения фотонно-кристаллических волноводов Текст научной статьи по специальности «Математика»

ФИЗИКА

УДК 537.8:621.371

ИНТЕГРАЛЬНЫЕ УРАВНЕНИЯ ФОТОННО-КРИСТАЛЛИЧЕСКИХ ВОЛНОВОДОВ

М.В. Давидович, Ю.В. Стефюк

Саратовский государственный университет E-mail: DavidovichMV@info.sgu.ru

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

Ключевые слова: фотонные кристаллы, диэлектрические волноводы, фотонно-кристаллические волноводы, интегральные уравнения, дисперсия.

INTEGRAL EQUATIONS FOR PHOTONIC CRYSTAL FIBERS

M.V. Davidovich, J.V. Stephuk

The integral and integrodifferential equations for photonic crystal waveguides (fibers) have been obtained both for finite and infinite dimensions of quasiperiodic dielectric coverings. Previously the equations for two-dimensional periodic photonic crystals with magnetodielectric and metallic periodic inclusions have been considered. The corresponding numerical results are presented.

Key words: photonic crystals, dielectric waveguides, photonic crystal fibers, integral equations, dispersion.

Введение

В последнее время широко используются и исследуются волноводы с оболочкой из метаматериалов в виде фотонных кристаллов (ФК) - фотонно-кристаллические волноводы (ФКВ) [1-3]. Имеется несколько типов таких волноводов, один из которых - волновод с полой сердцевиной, соответствующий диэлектрическим волноводам (ДВ) типа «полый канал в диэлектрике» (ПКД) [4-9]. Принцип работы таких ДВ и ФКВ основан на распространении в вакууме (сердцевине или канале) под малыми углами к оболочке при хорошем отражении от нее лучей, соответствующих модам. Указанные волноводы всегда являются многомодовыми с основной модой ИБц, вырожденной по поляризации, при этом с ростом частоты угол скольжения уменьшается, модуль коэффициента отражения от оболочки увеличивается, потери падают за счет преимущественного распространения в вакууме, а число распространяющихся мод растет. Все такие моды являются вытекающими, а потери основной и нескольких первых мод обусловлены в основном вытеканием и слабо зависят от потерь в толстой оболочке, в которую они почти не проникают. Для ФКВ предполагается наличие зоны bandgap (полного непропускания) соответствующего двумерно-периодического

МАУЧНЫЙ

ОТДЕЛ

© h.B. Давидович, Ю.В. Стефюк, 2009

ФК. Если бы полоса основной моды соответствовала этой зоне, а все высшие моды имели более высокочастотные отсечки, то это (с точностью до поляризационной ориентации) означало бы одномодовость ФКВ для моды ИЕП. Однако возмущение (полая сердцевина) приводит к нарушению зонной структуры бесконечно-периодического ФК и к разрушению зон bandgap (непропускания). Так, в такой зоне может появиться разрешенная полоса (дисперсионная ветвь). Наличие потерь в диэлектрике бесконечного ФК также приводит к разрушению зон bandgap [10], при этом дисперсионные ветви соответствующих прямых и обратных волн соседних зон Бриллюэна смыкаются, а в полосе bandgap возможно распространение волн с большими потерями (для двумерных ФК удобнее говорить о дисперсионных поверхностях). К таким же эффектам приводит и конечность ФК [10], т.е. наличие нескольких или нескольких десятков периодически повторяющихся слоев оболочки, что всегда имеет место в реальных ФКВ. Все указанные эффекты связаны с нарушением периодичности. Так, потери в материале нарушают периодичность из-за затухающего характера волн, а конечный ФК представляет собой, по сути, многополосовой фильтр, заграждение в полосах режекции которого всегда конечно. В конечном ФКВ также возможны и медленные поверхностные моды оболочки. Рассмотренный выше тип ФКВ, когда из сердцевины ФК удаляется материал, можно характеризовать как «акцепторный». Используются также и ФКВ «донорного» типа, когда в удаленную полость помещается другой диэлектрик (обычно оптически более плотный). Для анализа мод в ФКВ часто используют метод сверхрешетки (т.е. двумернопериодическое продолжение структуры реального ФКВ в поперечных направлениях) в сочетании с методом плоских волн [11]. Указанный подход является приближенным, поскольку не учитывает конечность оболочки. Он, в частности, не позволяет рассчитать потери ФКВ - как связанные с потерями в диэлектрике, так и с вытеканием.

Несколько десятилетий назад были предложены ДВ типа полый канал с многослойной квазипериодической оболочкой, об-

ладающие сходными с ФКВ свойствами [5]. Их анализ возможен, например, на основе метода матриц передачи [8, 9] или на основе методов типа прогонки [9]. Подобные волноводы с многослойной непериодической оболочкой могут демонстрировать улучшение селекции мод по потерям [8]. Их удобно использовать в миллиметровом, субмиллимет-ровом и квазиоптическом диапазонах. Для выбора направления линейной поляризации основной волны в таких ДВ используется поляризационная решетка из тонких параллельных металлических проволочек, ориентированных вдоль возбуждающего электрического поля (в КВЧ диапазоне это обычно прямоугольный волновод). Две повернутые на некоторый угол такие решетки представляют собой аттенюатор.

В данной работе предложены интегральные (ИУ) и интегродифференциальные уравнения (ИДУ) для ФКВ с неограниченной и ограниченной оболочками и приведен ряд численных результатов. Рассмотрение ведется на основе функций Грина (ФГ) двумернопериодически расположенных источников (2-D периодических ФГ), а также на основе квазипериодических ФГ [10]. Как частный случай, указанные уравнения применимы для многослойных ДВ класса ПКД.

1. Функция Грина 2-D ФК и интегральные уравнения

Рассмотрим бесконечный однородный и изотропный магнитодиэлектрик с диэлектрической проницаемостью S и магнитной проницаемостью в который периодически в плоскости (х, у) внедрены одинаковые металлические и/или магнитодиэлектрические объекты, однородные по оси z. Зависимость волн в такой структуре будем искать в виде exp(jwt-jyz), где у =kz соответствует z-компоненте волнового вектора k . Пусть металлические объекты в ячейке периодичности с номером 0 характеризуются совокупной двумерной областью с площадью поперечного сечения £м G S0 (и совокупным ограничивающим контуром LM ), расположенной в элементарной ячейке S0. Область £м может

быть многосвязной. Принадлежность объекта к металлическому типу означает выполнение

на его поверхности импедансного граничного условия V х Е = Z(у х(т7 х Н)), где 2 -поверхностный импеданс, а V - внешняя нормаль в области Ем к контуру Ьм . Принадлежность объекта к магнитодиэлектрическому типу означает, что он электродинамически эквивалентен объемным электрическим и магнитным токам поляризации с плотностями соответственно

УР {г ) = &ео (5{г ) - 81 ^{г ),

УтР {)=7®/^о ) - $Н{г),

расположенными в области е £0. Здесь

8, /л - тензоры макроскопических диэлектрической и магнитной проницаемостей соответствующих периодических включений,

I - единичный тензор.

Периодическое внедрение означает трансляционную симметрию, т. е. наличие прими-

тивных трансляционных векторов аі (периодов решетки) и векторов периода р = п1а1 + п2а2 = Ап . Здесь п - вектор сдвига (нумерации) ячеек с целочисленными координатами пі, і = 1,2, а А - трансляционная матрица (составленная из векторов-столбцов аі). ФК разбивается на ячейки , нумеруемые вектором п или мультииндексом N = (п1, п2) с одинаковыми площадями

поверхностей, для которых имеет место

соотношение |£^| = |£01 = |а1 х а21 = £0.

Для решения подобных задач наиболее удобен метод функций Грина (ФГ) периодически расположенных источников [10,12-16], приводящий к ИУ в одной (обычно в нулевой) ячейке £0, а точнее в области и на контуре Ьм . Такая скалярная ФГ имеет вид [10]

О (к 0^8,и, к±,

г - Г 1 =

) = ) ехр (± і (к + хпХ7 - Г')-, ¡к: ( - г'))

г (к1+Ф У+кг- кіи

-I-

ехРІ ± і

(1)

Здесь к0 = со / с - волновое число, трехмерное суммирование осуществляется по вектору п (мультииндексу N) в бесконечных пределах - да < пі < да , g = 2пА- -

тензор обратной решетки, к1 - приведенный волновой вектор, связанный с вектором фазовых сдвигов на ячейку: р = Ак1. Соответ-

ственно кL = А 1ф . В формуле (1) можно брать любой знак, поскольку в бесконечной периодической среде прямые и обратные волны неразличимы. Будем брать знак «минус». Если векторы аг ортогональны и направлены соответственно по осям X, у , то матрица А становится диагональной, а вид ФГ (1) упрощается:

1

да да

ех

2п1а2 _дап =-да п2=-да

Р(- і(к1 + ^п1 V - і(к2 + g2n2 )У -^Ь1^к =

2, 3

дада

2а а

пх--да п=-да

ехР

(1 + glnl )2 + (2 + g2n2 )2 + к32 - к02Є

і(1 + ^п1 )х - і(2 + g2n2 )у - ^ (2)

да

Здесь g = (п / а1) + (п / а2 )у0 - вектор обратной решетки. Индексы 1 и 2 будем сопоставлять х и у, если векторы а1, а2 колли-неарны орт-векторам Х0, у0. ИУ и ИДУ для

различных включений с ФГ трехмерно-периодической (3-0), двумерно-периодической (2-Б) и одномерно-периодической (1-0) ФГ сформулированы в работе [10]. Электромагнитное поле создается источниками в виде плотности поверхностного электрического

тока Ут на металлических поверхностях и

токами поляризации во включениях. Приведенные ФГ (1), (2) пригодны для любых зависимостей полей от г. Если же использовать экспоненциальную зависимость полей и тока с волновым числом кг = у, то интегрирование по г выполняется элементарно, при этом в ИУ и ИДУ появляется множитель exp(—}уг). Он, в частности, появляется в зависимости электрических и магнитных вектор-потенциалов:

(г) = |о(к0^в/й,к,г -г'рг(г')//'ёг' + іоє0|О(к0Л[е/,к,г -г')?^')-єі]Е(г')ъг', (3)

Ат (г) = і о,, | о(к^л[ей, к, г - г1')[^(гr,) - / ]н (г')ё Зг'

(4)

через которые стандартно выражаются поля. В зависимости (3) и (4) входят полные поля

волновода, при этом плотность тока УТ определена на поверхности цилиндра с конту-

ром Ьм . Сокращая на указанный экспоненциальный множитель, получим парные ги-персингулярные ИУ для электрического и магнитного поля, а также и тока:

£() = [(1 -і[)®(і -№0) + кІє/ |0/х (і -г1\К(1 ^'^у'

ІЛ

+

(5)

+ (і°є0єУ І0 гі -УН-іо/(і -і^0)х І0 г1 -^іМ'їОЖ'І ^Уу',

Цм J

н (гГ1)= іоє0 (у1 - і?г0 )х 10 є,(Хє, |Гі - ^іК^і'^ (г[)ёх'ёу' +

+

(1 - і[ )®(Уі - № )+ к 02єи]|0 |] - г[§Км{г[)Н ЙУ^у

(6)

В приведенных уравнениях величина Х8Л = д/к028// - у2 представляет собой поперечное волновое число в основе, тензоры К8{г]) = е{г])18-1 и км {г1) = $(гГ1)/¡Л- 1

соответствуют восприимчивостям включений соответственно для магнитодиэлектрической основы, ® - символ тензорного (ди-адного) произведения векторов, 2-0 ФГ G

є,

имеет вид

0 є,(0 1^1 - г1|)= ^ I Є( (+ і (кі)2+ '^п)(г 2 гі)) = I °І )(х0 (1 - К - Р|)еХР (- іп • р).(7)

£0 Г к1 + gnj + ” - к є// п

22 у - к0є/

Здесь введена двумерная скалярная ФГ [17], выражающаяся через функцию Ганкеля второго рода индекса 0:

о/'(хє, г |) = -(і / 4)н 02 Х Г1). (8)

Решать ИУ (5) и (6) следует одновременно. При этом, задавая к1 и к2, из решения можно определить дисперсию у = Г(ко ).

Объемно-поверхностные уравнения (5), (6) являются наиболее общими для продольных

Є

Б

волн в 2-0 ФК, при этом они решаются в «двумерном объеме» - на поверхности 50 и

на «одномерной поверхности» - контуре Ьм . Обычно желательно сформулировать уравнения только относительно электрического или магнитного поля. Для приведенных ИУ

это можно сделать, выразив Е и Н из уравнений Максвелла и подставив их под интегралы соответственно в уравнения (6) и (5). При этом получаются ИДУ. Для идеально проводящих включений с замкнутым контуром контурный интеграл можно преобразовать, выделив внеинтегральный член. Однако в высокочастотных диапазонах металлические включения обладают заметным импедансом. Вводя соотношение Ут = сЕт, мы не рассматриваем поля внутри контура Ьм , при этом учитывается как касательная к контуру в его плоскости компонента поля, так и Ег. В квазиоптическом и оптическом диапазоне следует учитывать электрофизические свойства металлов [18], при этом металл становится подобен плазме с соответствующими диэлектрическими свойствами, а контурный интеграл превращается в соответствующий

интеграл по поперечному сечению £ . Металлические включения в ФК оптического диапазона практически не используются, как не используется и магнитная основа, т. е. обычно л = 1. Магнитные включения, например, ферритовые, используются, но реже и в основном для управления магнитным полем дисперсионными свойствами. Таким образом, ФКВ обычно представляет собой диэлектрическую основу, в которую периодически в поперечных направлениях включены диэлектрические цилиндрические тела или полости. Заметим, что полость представляет собой дефицит тока поляризации основы, т.е. ИУ (5), (6) в этом случае описывают такие ФК, для которых решение ищется только в

полости, при этом к8 (т ) = -1)1,

кЛ (т ) = (л- -1)1. Часто полости заполняются другим диэлектриком, например, одноосным нематическим жидким кристаллом, что позволяет электрически управлять дисперсией [19]. В случае конечных ФК следует использовать ФГ вакуума (8 = л = 1). Далее будем рассматривать наиболее характерный для ФКВ случай л = 1, Л = I. В результате выражения (5) и (6) можно записать в виде

Е(і) = [і - і[)'®(і - і>^0) + к02/{] (] |] - (г1 )Е(г1]dХ 'ф', (9)

н (1)=

(Vl- і^?0)х Е(гД

°,0,

іоє0(і - і>?0)х 10є,{хє,\гі - г[\)(Ї)Е^'^ ’ёу’. (10)

Общий интеграл соотношения (10) показывает, что электрическое поле равно интегралу в правой части (10), умноженному на к0 плюс градиент некой произвольной функции: (і- і'?5 0) (). Если взять эту функцию в виде дивергенции от указанного интеграла

Р(Гі) = (у 1 - № )х

х I 0 , ( 1Г - ^іК (1' )Е (гІ' ^ ' ёУ' ,

то уравнения (9) и (10) совпадут. Используя известное векторное тождество и уравнение для ФГ, уравнение (9) можно переписать так:

є(Гі)Е(Гі)=(Уі - і^?0)х (Vl- утёс)х

х 10є,(Хє,\гі- ^іКОї^^'^ ' ёу'. (11)

^Б

Удобно использовать также ИДУ для магнитного поля:

н^іЬ ^1 - №0)х |0є,(Хє,|гі - г1$еє(гі)є-1 (г1')^1 - №0)х н2г', (12)

Б

Б

Б

Б

получающееся из уравнения (10) путем использования уравнения Максвелла в области

80 в виде ]&Е(г1 )=8- (г1 )(У1 - да )х

х Н (г1). Подстановкой этого соотношения в уравнения (9) или (11) получим ИДУ для магнитного поля. Такие уравнения использовать неудобно. Следует заметить, что ИУ (9) приводит к большому числу эквивалентных ИДУ для электрического поля, нагруженных контурными интегралами по контуру Ьв области и имеющих пониженную сингулярность ядер. Такие преобразования можно выполнить несколькими способами. Один из них - перенос однократный или двукратный

действия указанных операторов с ФГ на подынтегральную функцию. Именно функция Оє, логарифмически сингулярная. Действуя

на нее оператором V1 ® V1 = ((гаё)(йу) = = V1V1 •, получаем гиперсингулярное ядро с неинтегрируемой особенностью типа

г - г1\2. Поскольку V10' Лхє/Ті- г1|)=

г1- ги

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

і).

где штрих у опера-

VIV

-V|

• 10-Гі^є^'іЩ'іУМііу' = V110єм{хєм\гі -Гі'|)У 1 • [кє^Щг'^'ф'-

£б £б

10є,(хє,\гі - 10є,(%є,\Гі - •{кє^іЩг^'Ф' - (13)

ЦБ £Б

10є,(хє,\гі -- |і^гір/х/Ті -гЦ)уі •{кє^іЩг^Я'.

Здесь использована теорема о дивергенции (Гаусса) и теорема о градиенте. В уравнениях типа (12) следует использовать теорему о роторе. При двойном переносе ядро имеет слабую логарифмическую сингулярность, однако использовать подобные ИДУ проблематично из-за алгоритмических сложностей при вычислении контурных интегралов, особенно если область имеет сложную конфигурацию. Альтернативные методы получения ИДУ основываются на лемме Лоренца, формулах Стрэттона-Чу и формулах Грина. Понижение сингулярности возможно и путем непосредственного интегрирования по координатам точки наблюдения [20].

Другой подход основан на получении сингулярных ИУ. Пользуясь значениями функции Ганкеля при малых аргументах, особенность ядра в (13) представим в виде

0 8Л(%8ЛК)~ (2п) 1п(2/(Г0к0л18Лк)), где /0 =

= 1.7810724, Я = |К - К'|. Поэтому можно записать О є,(хє,

я)=ф(я НОМ) где

є,

Ф(я) = -(2п)-1 1п(к) - логарифмический потенциал, а 80 8Л - функция регулярная. Для

логарифмического потенциала обычным образом с использованием формулы Остроградского доказывается свойство дифференцирования его распределения с плотностью р [21]. Рассмотрим объемный двумерный потенциал

V(Ю = \р(г1)08Л(Х8л|Г1 - Г1|К'Ф'. (14)

Очевидно, у него существуют первые производные. Окружим точку г1 кругом 58

радиуса 8. Очевидно, в области 5» - 58 существуют все вторые производные потенциала (4). В области 58 представим функцию О8Л указанным выше образом. Стягивая 58 в точку и пользуясь теоремой о среднем [21], получим:

Б

Б

ііш

д2

V(г1)= 1іШ 5^0 V (г1)= 1іШ 5^0 | Р(і У°єи (-

(|К -П)у' = -£&)

2

+

+

8 —— 0 2 \ 1 / 8 —— 0

дх ду

' 5 8

Запишем уравнение (9) относительно поперечной и продольной частей поля:

Е1(г1 Ь - 1(к8(( )Е(г1 ))1 + )028л\ 0 8Л (I Т1- А1 |)(к8 (г1 )Е(г1 ))1¿х''¿у '

2 5»

|((1)Ей)-^к)((|Т1 -г1|)у ' + р.У. |((1)ЕЙ!))• ^Кф( -г^'ёу'

5В 5В

- ^¡У108Л(8М\Т1 - Г^К8((Е(г[))У

5В

Ег (г1) = -Уг1 Vl0 8Л(| Т1- Ах!)- ((г1)Е(г10)у' +

(8Л-У2 )) 08Л(Т1 - Г1|)(*8 (( )Е(Г1 ))г ¿х)'

+

(15)

(16)

Здесь р.у. означает главное значение интеграла, получаемое при выкалывании 8 -окрестности точки наблюдения. Все рассмотренные уравнения упрощаются для изотропных включений, для которых к8 (г1) скаляры, а также для однородных сред, для которых к8(П) - константы. При совпадении обоих случаев величина к8 может быть вынесена за знак интеграла.

2. Дисперсионные уравнения бесконечного 2^ ФК

Далее для простоты будем рассматривать ИУ и ИДУ только для электрического поля в ФК без магнитных и металлических включений. Дисперсионное уравнение (ДУ) для собственных волн определяет зависимость Р (к0, к1 ,у) = 0 или к0 = / (к1 ,у).

Указанные функции зависят также от параметров, определяющих конфигурацию ФК. Для получения ДУ следует решать приведенные уравнения в области 5» . Введем обычным образом гильбертово пространство вектор-функций в области 5 » со скалярным произведением

(Е1,р)5 =|ЕГ(г1)-Р2(г1 ^хФ ,

Пусть Е(П) - некоторая функция из указанного пространства. Пусть Е£б (г^) -

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

л(р йЛ Е ЙЛ кс кі ,г)=

= (Р ЕЮ- Е£б

(18)

зависящий от параметров к0, к1, у . Его стационарное значение 0 достигается на точном решении задачи Е0 (г1). При подстановке в

уравнение (12) точного решения ДУ принимает вид

л(Е Е0 (г1), к0, к1 ,г) = 0. (19)

Для приближенного получения указанного решения пробную и весовую функции необходимо представить в виде разложения по полной системе функций ит (г1) рассмотренного пространства соответственно с коэффициентами ат и вт, а затем наложить

условия экстремума (д / двт )л = 0 , т = 1,

2, ..., М (здесь Л - квадратичная форма, со-

Б

£

Б

ответствующая уравнению (19) при учете М функций). Тогда ДУ (19) приближенно имеет вид равенства нулю соответствующего определителя. При вычислении матричных элементов следует учитывать сингулярности ядер, которые во всех приведенных соотношениях являются интегрируемыми. В частности, при использовании ИУ (9) действие одного из операторов можно перенести на

функцию Р . Если потребовать, чтобы нормальная компонента этой функции обращалась в нуль на контуре Ь0, контурные интегралы не возникают. В этом смысле билинейный функционал удобнее, чем квадратичный, для которого Р = Е, и контурные интегралы при использовании ИУ (9) обычно появляются. В качестве базисов удобно использовать непрерывно-дифференцируемые двумерные векторные конечные элементы (КЭ), при этом достаточно легко аналитически определяются диагональные матричные элементы. Итак, считаем, что точное решение Е0 и точное значение к0 при заданных

значениях к1 и у определены. Если бесконечный 2-0 ФК имеет диэлектрические потери, величины к1 и у становятся комплексными. В этом случае положим к1 = 0, т.е. будем искать комплексные значения у = у’ - ]у” при у' , у” > 0 и периодических полях в поперечных направлениях. В случае комплексных у удобно задавать у” и у" и

для них искать к0 . Если такое значение не

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

3. Интегральные уравнения неограниченного квазипериодического ФКВ

Пусть бесконечный 2-0 ФКВ представляет собой дефект ФК в виде акцепторных или донорных включений в одной или нескольких соседних ячейках, центр которых мы расположим в начале координат. В качестве акцепторного включения может, например, выступать полая сердцевина, а в качестве донорного - включение с проницаемостью 8(г1), обычно (но не обязательно) оптически

более плотного, чем 8 . В частности, ~(г_)=I соответствует акцепторному включению. Пусть указанная область занимает двумерный объем , который состоит из минимального целого количества ячеек ФК, покрывающих неоднородность. Естественно, в ФКВ следует рассматривать лишь волны

вдоль оси г, т.е. полагать к1 = 0. При этом на больших расстояниях от неоднородности в поперечном направлении волна либо вытекает, либо носит поверхностный характер и экспоненциально затухает. В первом случае у всегда комплексное, даже при отсутствии потерь в диэлектрике. Соответственно неоднородность приводит к тому, что полное поле Е (А1) отличается от решения Е0 (г1) на

величину АЕ(Г1), т.е. Е(г1) = Е0(^ + ,

так же, как и изменяется постоянная распространения у = у + Ау при данном волновом числе к0. Соответственно для поправки аЕ (г1) следует получить уравнение. Если волна поверхностная, то поле аЕ (1) затухает в поперечном направлении. Пусть 50 -

область, вне которой полем АЕ(1) можно пренебречь. Она включает в себя область неоднородности, а вне ее поле можно считать периодическим. В этом случае сформулировать уравнения наиболее просто. Неоднородность приводит к дополнительному

току поляризации АУР (?1) = ]о808еЕ(), 88 = 8(^1)-8^1))/ 8 . Этот ток отличен от нуля только в области и только там в

ней, где 88 0 . Все приведенные выше ИУ

для собственных продольных волн ФК можно записать в виде

Е0 (г1)= (20)

(, г1,к 0, у)к8 ( )Е0 (()х '¿у ” ,

где К((,г",к0,у) - некое ядро, зависящее

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

торы дифференцирования (для ИДУ). Указанное ядро разложимо в двумерный ряд

ад ад

К (Д , Г1, к0,г)= X X Ктп (Г1. Г1. к0,Г)

т=-ад п=-ад

ад ад

(21)

X X к( ,г1 -тах -па2,к0,х).

т=-ад п=-ад

Этот ряд означает суммирование в уравнении (20) по всем источникам тока поляризации во всех ячейках. Учитывая дополнительный ток поляризации, запишем ИУ ФКВ:

АЕ (/г1 )= | К0 (/г1, ^, к0, г)ке (г1 )аЕ 0± Ух ¿у ” + {8К (Т1, , к0,у,у К (г1 )Е0 (/1)х '¿у ” +

+

18К~0 ((, Т1, к0, г, г)) 0± Т (/1 )х'Ф ”+18К ^, Т1, к0, г, г )е0 (1)+АЕ (1 )]^х '¿у'.

(22)

Здесь 8К~ (, г1, к0, / ) =К (/1, г1, к0,у ) -

- К (, Г1, индексы 0 и Ж у ядер в

уравнении (22) означают, что в соответствующих суммах типа (21) учтены только источники, попадающие в области 50 и 5Ж .

Итак, при условии, что к0, у и Е0 известны,

ИУ (22) позволяет определить у и АЕ, например, путем поиска экстремума функционала типа (18), а также используя итерационный процесс для ИУ (22) и соответствующего ему функционала. Заметим, что из функционалов можно непосредственно выразить постоянную распространения в виде

у = х¥(к00,у), где в правую часть у входит нелинейно. Если она вытекает, то у”> 0, и поле в радиальном направлении далеко от неоднородности возрастает по закону р~1/2exp(- у'кр), к" > 0, к" > 0, если потерь в диэлектрике нет. Указанное свойство обусловлено свойством функций Ганкеля. Если есть потери, например, в основе, то поглощение может превысить возрастание от вытекания, и тогда вышеприведенные уравнения применимы.

Рассмотрим ФКВ типа полый канал в

2-0 ФК без потерь. Очевидно, ИУ для АЕ следует из соотношения

Е (г1 )= Е0 (г1 )+ АЕ (г1 ) = |К (г1. .к0. Г ')к8 (/1)Е ^ Vх ' ¿у' - ( - - 1)) К „ (Т1 > > к 0, Г )[Е 0 (/1 )+ А Е (Г1 )]у'

и имеет вид

АЕ (/Г1 )= 188~ (г1 > Г1 >к0, Г, Г К (Г1 )[е0 (Г1) + АЕ (г1 )]^х '¿у' -

(23)

( -1 -1)| кж ( ’ Г1 ’ к0’ ^ )Е0 Й)+ аеТ (Г1 )]^х '¿у' +|К (Г1, Г1, к0, г)х-8 Й )аЕ 01 Vх '¿у' •

О

О

Ж

Ж

О

Ж

О

4. Интегральное уравнение ограниченного ФКВ

Ограниченность всегда означает нахождение структуры в вакууме, т.е. є=л = 1. Ограниченность также означает конечное число ячеек у ФК, что соответствует ФГ свободного пространства с конечным числом членов в суммах [10]. Обычно ФКВ имеют оболочки, близкие к круглым. В качестве примера рассмотрим наиболее простой для анализа квадратный ФКВ типа полый квадратный канал в квадратной оболочке (рис. 1).

Рис. 1. Прямоугольный ФКВ типа «полый канал в оболочке»

Пусть ФК состоит из 81 ячейки размера а1 = а2 = а с проницаемостью є, а каждая ячейка имеет в центре полость размера Ь х Ь, Ь < а. ФКВ получается из ограниченного ФК путем удаления девяти центральных ячеек. Таким образом, размер полого канала 3а х 3а , а оболочка имеет по три слоя в каж-

дом из поперечных направлений. В данном случае квазипериодическая ФГ представляет собой конечную сумму типа (7) по 72 ячейкам при k± = 0. Рассмотрим сначала бесконечный ФК и найдем его дисперсию y = y(k0) при kx = ky = 0. Выберем начало

системы координат в центре полого прямоугольного канала. Удобно взять ФГ для диэлектрической среды. Тогда область поиска решения есть полый канал — b/ 2 < x < b/ 2, - b / 2 < у < b / 2. В этой области электрическое поле соленоидально, удовлетворяет волновому уравнению и имеет все три компоненты. С учетом этого и симметрии решение можно искать в виде

Ex = jya sin (ax )sin (py )exp(- jyz ),

Ey = - jrP cos(ax )cos(ey )exp(- jyz ), (24)

Ez = (k02 - Y2 )cos(ax)sin(p)exp(- jy),

где a2 + p2 + y2 = k02. Можно рассмотреть и другие типы симметрии (поляризации). Весовые функции F удобно выбирать так, чтобы они обладали той же симметрией, что и искомое решение. Возьмем соленоидаль-

ную вектор-функцию F такую, чтобы ее нормальная компонента равнялась нулю на контуре. Очевидно

Fx = уу(2п/ b)sin (2пх / b)sin (пу / b)exp (- jyz), Fy = - jy(n / b )cos (2nx / b )cos (пу / b )exp(- jyz ), Fz = [(2п / b )2 + (п / b )2 ] cos (2nx / b )sin (пу / b )exp (- jyz ).

(25)

Для высших мод можно рассматривать соот- налы на основе ИУ (9), умножая его скаляр-

ношения типа (25) с соответствующим чис- но на Р1(х, у) и на ?0. Экспоненциальную

лом вариаций по осям. Построим функцио- зависимость от г будем опускать. Тогда

1i =

+

jFi'(fi)Ei(fi)dxdy -k0(1 -s)J JgЄ(ХЄ\А -rÍl)Fi*(rÍ)-E!(ri,)dx 'dy

SD SDSD

( -i)J/Vi • Fi'(ri)?iO (Iri -r$-¿i(r1)dx'dy' + (26)

SDSD

+ j7(l - 1)J JFi*(ri^ViGe(cs|ri - r$Ez(r1,)dx 'dy

+

snsn

I =

IРг*(г±К(г^хЛу-\к18-У2] -!)|Е{хе\г -Г±|И-

^л^л

-7У(1 -1)) \К{г±)У±& Е(е\г± -㱧■Е±{r±)dx'dy'.

(27)

Интеграл (26) можно разбить на сумму 1±= 4 + 1у . Заметим, что для действительных у < к0 (т.е. для быстрых волн в отсутствие потерь) а ив должны быть действительными. Однако это возможно только в пределе к0 ^ да . В этом случае отражение от стенок полное, а а ив соответствуют волноводу с идеальными стенками. Поиск решения будем производить методом итераций. Пусть при большом к0 заданы начальные значения а , в. Тем самым определено и начальное значение у. Соответственно можно вычислить интегралы (26), (27). Из последнего определим скорректированное значение у . С его использованием определяем из 1Х новое значение к0 и соответственно

значение а = -\/к0 - в - у1 . Затем тем же

путем из I определяем к0 и в=\к2 -О2 -^ . Наконец, из функционала I = 1± + 1г находим новое значение к0 . Процесс повторяем

до сходимости результатов. Заметим, что в силу определения в (26) можно заменить

У± • Р± (г±) на - ]у*Ег (г±). При итерациях мы последовательно изменяем как собственные значения (волновое число к0), так и собственную функцию (величины а,в,у). Коль скоро собственная функция найдена внутри полости, она с помощью интегральных соотношений определена везде. Когда сходимость для данного к0 достигнута, величина

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

Для ФКВ с бесконечной оболочкой следует использовать ИУ (23), решение которого следует искать в области полого канала (поскольку е ). Для этого также

можно использовать представления решений (34) и пробные функции типа (25) с заменой у на у и Ь на 3а. В качестве исходного интегрального оператора в (23) можно взять ИУ (9) или другое, например, (16), а также любое ИДУ, нагруженное контурными интегралами. Использование представлений полей типа (24) приводит к экономным алгоритмам, поскольку неизвестными являются только три в общем случае комплексных параметра. Применение разложения полей по каким-либо базисам (например, по двумерным конечным элементам) приводит к алгоритмам достаточно большой размерности.

Для ФКВ (см. рис. 1) примем соотношение Ь = а /3 и будем использовать кусочнопостоянные базисные КЭ. Тогда каждая из 72 ячеек разбивается на 8 конечных элементов с тремя значениями компонент поля, а размерность задачи становится равной 1728. Однако в силу симметрии достаточно учесть лишь четверть области ФКВ, и реальная размерность равна 432. В случае использования ИУ (16) в качестве весовых функций можно также взять кусочно-постоянные КЭ. Тогда на основе метода Галеркина задача сводится к векторному уравнению Е = К{к0,у)Е в 432мерном пространстве. Здесь Е - вектор из 432 коэффициентов разложения, а К - соответствующая квадратная матрица, зависящая от к0 и у. Интегралы в ее диагональных

элементах сходятся. Величина к0 стандартно

выражается через отношение двух квадратичных форм:

2 = (Е, Р(к„,г)Е) 0 (Е, 0[к„,Г)Е)

(28)

л

В них Р(к0,г) и б((0,Г) - определенные матрицы, в которые параметры к0 и у входят нелинейно. Находить комплексные корни определителя порядка 432 неэффективно. Более эффективно методом совместной итерации искать решения уравнения (28) и

уравнения Е = К(к0, у)Е . Использовать кусочно-постоянные весовые КЭ для ИУ (9) нельзя, поскольку они не дифференцируемые, а действие оператора У± ® У± на логарифмически сингулярное ядро приводит к неинтегрируемой особенности. Использование ИДУ требует применения дифференцируемых КЭ и приводит к громоздким контурным интегралам и алгоритмам. Для получения численных результатов воспользуемся дифференцируемыми весовыми КЭ типа (25), например, взяв Рутп(х, у) = -]у(2п - 1)п/с) х

х сов(2(т - 1)х / с)сов((2п -1) / с), с = 9а.

Указанные функции заданы на всей области ФКВ. Использовать их в качестве пробных КЭ нельзя, так как их нормальные компоненты непрерывны на границе раздела диэлектрик - вакуум. Взяв 1 < т, п < 12, получим 432 весовые функции. При интегрировании

по координатам точки наблюдения можно сначала вычислить интеграл по всему поперечному сечению ФКВ, а затем вычесть из него интегралы по светлым участкам (см. рис. 1). Итерационный подход с использованием квадратичных и билинейных форм типа (28) позволяет определить все волны ФКВ, включая и поверхностные. Для того чтобы метод сходился к требуемой моде, надо налагать соответствующие условия симметрии, нормировать решение на каждом шаге итерации и требовать ортогональности искомой моды к ранее найденным. В случае поиска вытекающих мод параметр у комплексный, поэтому сходимости к собственным поверхностным модам быть не может.

Рассмотрим, наконец, интегрирование для логарифмической особенности с кусочно-постоянными КЭ для диагональных матричных элементов. Соответствующие подынтегральные выражения можно записать так: [(д / дх) 1п(Я)Ех (х;, у;) + (д / ду)1п(Я) х

х Еу (х;, у;)] ¿х'¿у'. Рассмотрим интеграл для

первого члена, обозначив значение х-компо-ненты Ех (х; , у;) в центре КЭ с номером г

как Е.:

Ехг | (х х Мх х ) Лх'¿у' = Ех. | (х^2Хг ) Лх'¿у' +Ехг | -С°Р) Лх'¿у' . я Я я Я я Р

(29)

При замене прямоугольного КЭ Яг кругом бесконечно дифференцируемой функцией

эквивалентной площади последний интеграл в уравнении (29) исчезает. Интегрирование первого интеграла по точке наблюдения с

/ (х, у )= У± • Рг ( х, у) по элементу Яг путем перехода к полярным координатам дает

Е,

/г (х, у )х - х; ) ¿х ; ;

¿х' ¿у' ¿хЛу =Ехг + р СОБ (р) у; + р б1п (р))сО8 )(рУ1рЛф(1хх¿у' =

Ех № (-у ) + Их (х;, у; )рсоб (р)+ / (х;, у; (р)+ ...] СОБ (р)dpd рх1х '¿у'.

Очевидно, если область Я. является

кругом, интеграл от первых четырех членов в квадратных скобках равен нулю.

Для малой области рассмотренный интеграл мал, и им можно пренебречь.

5. Численные результаты

На рис. 1-6 приведены результаты ряда численных экспериментов, иллюстрирующие свойства ФКВ и волноводов типа ПКД.

На рис. 2 демонстрируется дисперсия бесконечного двумерно-периодического ФК при поперечном распространении волн вдоль одной из координат.

ка

у а

Рис. 2. Дисперсия кх (к0) в 2-Э ФК при ку = у = 0 : е= 12.0 (кривые 1) и е= 3.0 (кривые 2)

у а

ко а

к о а

б

Рис. 3. Дисперсия для 1М- моды (а) и замедление (б) в плоскопараллельном полом канале при следующих параметрах бесконечной оболочки: е= 3.0 (кривая 1); е= 12.0 (кривая 2); е= 3.0 -/12.0 (кривая 3); е= 1.0 - /104 (кривая 4)

к о а

Рис. 4. Нормированные потери Яе(у)а ¿М-моды. Обозначения кривых см. на рис. 3

у' а, п

к о а

Рис. 5. Замедление п = Яе(у) / к0 (штриховые кривые) и нормированные потери -1ш(у)а (сплошные кривые) волны НЕ11 для волновода класса ПКД (см. рис. 1) без использования фотонно-кристаллической структуры (без полостей): е= 3.0 (кривые 1); е= 3.0 -- /3.0 (кривые 2); е= 12.0 -/12.0 (кривые 3)

Рис. 3 и 4 соответствуют результатам для волновода типа полый плоскопараллельный канал ширины а в магнитодиэлектрике, ограниченный полубесконечными импеданс-ными стенками. Стенки предполагаются выполненными из полубесконечного магнито-диэлектрика, в результате чего на его границе раздела можно задать импеданс, а дисперсионные уравнения записываются явно ана-

а

п

У а, !п(л)

ко а

Рис. 6. Логарифм замедления (сплошные кривые) и нормированные потери (штриховые кривые) первой НЕ-моды полого канала размера а*а со сплошной бесконечной плоскопараллельной оболочкой (1), полого канала с плоскопараллельной оболочкой (2) и фотонно-кристаллического волновода (3). Во всех результатах е= 3.0 -- 0.05/, а характерные размеры каналов совпадают

литически [22], что важно для качественной оценки других подобных результатов. С практической точки зрения магнитодиэлектрическую пластину можно считать полубес-конечной, если все излучение из канала в оболочку затухает. Рассмотренный случай моделирует ряд задач. В частности, задавая плазменную диэлектрическую проницаемость, на низких частотах (существенно ниже плазменной частоты и частоты столкновений) получаем импеданс металла.

На рис. 5 приведены результаты для прямоугольного волновода типа ПКД, соответствующего ФКВ без полостей (см. рис. 1). В этом случае вместо двумерной суммы в ФГ взят один член с функцией Ганкеля. Заметим, что как для плоскопараллельного, так и для прямоугольного волновода типа ПКД имеет место переход быстрой вытекающей волны на высоких частотах в медленную вытекающую волну на очень низких частотах. Это связано с тем, что соответствующие строгие дисперсионные уравнения не содержат ветвей, исходящих из точки к0 = 0 , у = 0 . Для ФКВ (см. рис. 1) и соответствующего волно-

вода типа ПКД выше частоты отсечки, определяемой из условия п = 1, возможны и решения в виде медленных поверхностных мод. Однако на частотах ниже отсечки вытекающая волна сначала становится быстрой, а замедление п достигает минимума. Затем с уменьшением частоты замедление растет, становится больше 1 и продолжает рост (медленная вытекающая волна). Указанные моды следует искать в комплексной области, и в данной работе они не рассматривались.

На рис. 6 приведены результаты влияния квазипериодичности оболочки на замедление и потери для плоскопараллельного ФКВ с двумя одномерными квазипериодическими оболочками и ФКВ (см. рис. 1). Первый случай соответствует рис. 1, когда один из размеров полого канала бесконечный. Кроме этого жирными линиями (кривые 1) нанесены результаты для ДВ типа ПКД без полостей в оболочке. Видно, что имеется набор полос, в которых потери малы, а фазовая скорость близка к скорости света. Они разделены участками с более высокими потерями. Соответствующие полосы для ближайших высших мод обычно смещены, а затухание их выше, что при большой длине позволяет осуществить квазиодномодовый режим. Рассмотренный ФКВ имеет двукратное вырождение основной моды по поляризации.

Заключение

На основе метода ФГ периодических и квазипериодических источников сформулированы интегральные и интегродифференци-альные уравнения для фотонно-кристаллических волноводов. Предложенные уравнения использованы для построения дисперсионных уравнений в виде билинейных и квадратичных функционалов. Также предложены итерационные процедуры решения полученных уравнений, позволяющие проводить решение в комплексной области и выделять требуемую моду из всего спектра мод. Приведены результаты численных исследований дисперсионных уравнений, включая случаи волноводов типа ПКД, для которых имеются строгие аналитические соотношения. Указанные решения качественно согласуются с

приведенными решениями для двумерных структур. Результаты показывают, что рассмотренная толщина оболочки ФКВ недостаточна для хорошей селекции мод. Уравнения для полубесконечной периодической оболочкой приведены в [23]. Для плоскослоистого волновода типа ПКД их можно записать в явном аналитическом виде.

Список литературы

1. Briks T.A., Roberts P.J., Russell P.St.J., Atkin D.M., Shepherd T.J. Full 2-D photonic bandgap in silica/air structures // Electron. Lett. 1995. Vol.31, №22. P.1941-1943.

2. Cregan R.F., Mangan B.J., Knight J.C., Briks T.A., Russell P.St.J., Roberts P.J., Allan D.C. Single mode photonic bandgap guidanceof light ion air // Science. 1999. Vol.285, №5433.P.1537-1539.

3. Желтиков А.М. Дырчатые волноводы // УФН. 2000. Т.170, №11. С.1203.

4. Мелехин В.Н., Маненков А.Б. Диэлектрические трубы как волноводы с малыми потерями // ЖТФ. 1968. Т.38, №12. С.2113—2115.

5. Мелехин В.Н., Маненков А.Б. Диэлектрические трубы -открытые волноводы с малыми потерями и редким спектром // Электроника больших мощностей. М.: Наука, 1969. Вып.6. С.161-178.

6. Казанцев Ю.Н., Харлашкин О.А. Прямоугольные волноводы класса «полый канал в диэлектрике» // РЭ. 1978. Т.23, №10. С.2060-2068.

7. Chou R.C., Lee S.W. Modal attenuation in mu1tilayered coated waveguides // IEEE Trans. 1988. Vol.MTT-36, №7. P.1167-1176.

8. Давидович М.В. Улучшение самофильтрации мод в полых волноводах с многослойной магнитодиэлектрической оболочкой // РЭ. 1994. Т.39, №1. С.53-61.

9. Голант Е.И., Голант К.М. Новый метод расчета спектра и радиационных потерь вытекающих мод многослойных оптических волноводов // ЖТФ. 2006. Т.76, вып.8. С.99-106.

10. Давидович М.В. Фотонные кристаллы: функции Грина, интегродифференциальные уравнения, результаты моделирования // Изв. вузов. Радиофизика. 2006. Т.49, №2. С.150-163.

11. Benisty H. Modal analysis of optical guides with twodimensional photonic band-gap boundaries // J. Appl. Phys. 19965. Vol.79. P.7483-7492.

12. Kohn W, RostokerN. // Phys. Rev. 1954. Vol.94. P.1111.

13. Цидильковский И.М. Электроны и дырки в полупроводниках. М.: Наука, 1972. 640 с.

14. Иваненко Д.Л., Соколов А.А. Классическая теория поля. М.: Гостехиздат, 1951. 480 c.

15. Singh S., Richards W.F., Zinecker J.R., Wilton D.R. Accelerating the convergence of series representing the free space periodic Green's function // IEEE Trans. 1990. Vol.AP-38, №12. P.1958-1962.

16. Yasumoto K., Yoshitomi K. Efficient calculation of lattice sums for free-space periodic Green's function // IEEE Trans Antennas and Propagat. 1999. Vol.AP-47, №6. P.1050-1055.

17. Марков Г.Т., Чаплин А.Ф. Возбуждение электромагнитных волн. М.: Радио и связь, 1983. 296 с.

18. OrdalM.A., Long L.L., BellR.J., BellS.E., BellR.R., Alexander P.W., WordC.A. Optical properties of the metals Al, Co, Cu, Au, Fe, Pb, Ni, Pd, Pt, Ag, Ti, and W in the infrared and far infrared // Appl. Optics. 1983. Vol.22, №7. P.1099-1120.

19. Du F., Lu Y., Wu S. Electrically tunable liquid-crystal photonic crystal fiber // Appl. Phys. Lett. 2004. Vol.85, №12. P.2181-2183.

20. Давидович М.В. Численный метод для задач дифракции на полубесконечных структурах // Изв. вузов. Радиофизика. 2006. Т.49, №3. С.235-245.

21. Тихонов А.Н., Самарский А.А. Уравнения математической физики. М.: Наука, 1977. 736 с.

22. Вайнштейн Л.А. Открытые резонаторы и открытые волноводы. М.: Сов. радио, 1966. 476 с.

23. Davidovich M. V. Integral and integrodifferential equations for unbounded pseudoperiodic structures // Mathematical Methods in Electromagnetic Theory: Proc. of 12th Intern. Conf. (MMET’2008). Odesa, Ukraine. 2008. P.346-348.