CC BY
52
РАСЧЕТ ПРОСТРАНСТВЕННЫХ МОД ОПТИЧЕСКИХ ВОЛНОВОДОВ С НЕОДНОРОДНЫМ ПОПЕРЕЧНЫМ СЕЧЕНИЕМ МЕТОДОМ СОГЛАСОВАННЫХ СИНУСОИДАЛЬНЫХ МОД
В.В. Котляр, Я.О. Шуюпова* Институт систем обработки изображений РАН, *Самарский государственный аэрокосмический университет
Аннотация
В данной работе предложено некоторое усовершенствование метода согласованных синусоидальных мод, основанное на применении матричной формы записи характеристического уравнения. Метод реализован с помощью программного пакета МаНаЪ 6.0 и применен для анализа в рамках скалярной теории собственных мод волновода, выполненного в виде фотонного кристалла.
Введение
Для расчета векторных или скалярных пространственных мод оптических волноводов и волокон с неоднородным поперечным сечением применяются несколько методов: техника поперечного резонанса [1], метод согласованных синусоидальных мод [2, 3], метод линий [4], метод распространения собственных мод [5].
Также расчет электромагнитных полей в неоднородных по сечению волноводах можно проводить с помощью более универсальных расчетных методов: конечных элементов [6, 7], разностного решения уравнений Максвелла [8, 9] или распространения пучка [10, 11].
Поскольку результат расчета электромагнитных полей в волноводах зависит от возбуждающего поля, то для характеристики самого волновода требуется найти весь спектр возможных пространственных мод. Поэтому в данной работе для расчета мод неоднородных оптических волноводов выбран метод согласованных синусоидальных мод (ССМ) [2, 3], который наиболее органично подходит для анализа диэлектрических волноводов с кусочно-постоянным заполнением поперечного сечения.
ССМ-метод основан на представлении решения для пространственной моды в виде суперпозиции локальных синусоидальных мод, которые являются собственными модами однородных частей волновода с прямоугольным сечением.
1. Описание метода согласованных синусоидальных мод
Рассмотрим метод отыскания поля в волноводе неоднородного поперечного сечения. Для иллюстрации будем использовать волновод, изображенный на рис. 1. Предположим, что сечение волновода может быть разбито на N строк и М столбцов так, чтобы ни одна из прямоугольных ячеек этого разбиения не содержала неоднородностей (рис. 1).
Таким образом, каждой ячейке, находящейся на пресечении п - ой строки и т - го столбца, можно поставить в соответствие некоторое значение диэлектрической проницаемости г(т,п, постоянное для данной ячейки. Пусть оси координат располагаются, как показано на рис. 1. Тогда толщина п - ой строки:
4n) = y(n+1) - y(n), (1)
где y(n) - координата плоскости, разделяющей n - 1 и n строки. Аналогично, толщина m - го столбца
d(m) = х(m+l) - x(m), (2)
где х(m) - координата плоскости, разделяющей m -1 и m столбцы.
Í 9\—--------
7---------
6---------
5---------
4---------
3---------
2---------
1---------
0 123456789
X, мкм
Рис. 1. Схема поперечного сечения волновода, выполненного в виде фотонного кристалла; темным цветом показаны области со значением диэлектрической проницаемости е(т,п) = 3 ,
(т п) л
светлым - с ек ' ' = 1
Таким образом, диэлектрическая проницаемость в т - ом столбце есть функция
е(т) (у) = е(т,п) для у(п) < у < у(п+1). (3)
На краях сечения подразумевается наличие электрических или магнитных стенок, обеспечивающих равенство нулю функции поля или ее производной при X = 0, у = 0, X = X(М+1), у = у ^+1).
Рассмотрим уравнение Гельмгольца для монохроматического света:
V 2 Е( х, у, х) + к 02 е( х, у)Е( х, у, х) = 0, (4)
здесь Е(х, у, х) - проекция на продольную ось г вектора электрической напряженности электромаг-
2п
нитного поля; к 0 = —, где X0 - длина волны в ва-
X о
кууме; s(x, у) - диэлектрическая проницаемость среды, зависящая от пространственных координат, в общем случае есть зависимость и от координаты z , но мы считаем сечение волновода постоянным. Решение уравнения (4) будем искать в виде
E = и(х)ф(y)e. Подставим (5) в уравнение (4):
-ikzz \ . 7 2
(5)
V 2(и (х)ф( у )е) + к02е( х, у )и( х)ф( у)е" = 0, У2(м(х)ф(у))е-Кг - к]и{х)ф(у)е~'к-2 + +к02е( х, у)и( х)ф(у)е-к'2 = 0, V2(u( х)ф( у)) - к2и( х)ф( у) + +ко2е( х, у )и( х)ф( у) = 0, и (х)ф( у) + ф (у)и( х) + +(ко2е(х, у) - к2 )и(х)ф(у) = 0. Далее будем опускать аргументы у, и и ф , чтобы записи были более короткими. Итак, мы получили
иф + фи + (к0 е(x, у) - к,, )иф = 0 .
(6)
Проведем разделение переменных, поделим (6) на иф :
ф +1 + (к02е(x,у) - к,) = 0 ф и
или
ф + (к(2 s(x, у) - к22) = - -фи
Положим
и =-к2.
теперь (8) можно разделить на два уравнения
(ф+(к02 в(x, у) - к, - кх2)ф=0,
[и + к2 и = 0. Обозначим
к, + к, = к2, (11) тогда
(7)
(8)
(9) (10)
(ф + к2е( x, у) - к2)ф = 0, [и + (кк2 -к,)и = 0.
(12)
Первое уравнение в (10) и (12) будет далее решаться только при е(х, у) = е(т)(у) (17). Введем еще одно обозначение
к0 е( x, у) - кк = к2,
(13)
то есть теперь
к02е( x, у) = к2х + к2 + к2.
к0 е(X, у) = кх + ку + к2 . (14)
Значит кх, ку, к2 - проекции на соответствующие оси координат волнового вектора к , модуль которого равен | к |= ^к0е(х, у) .
Решение системы уравнений (12) зависит от функции е( х, у) и не может быть получено в общем
виде в алгебраической форме. Поэтому далее рассмотрим способ приближенного решения.
В скалярном случае поле в т - м столбце можно представить в виде
x, у) = £ икт)( x)ф кт)( у),
к=1
ikzZ
(15)
здесь опущен множитель е 2 . Каждая из множества функций и к(т) ( х) удовлетворяет второму, а каждая из фкт)( у) - первому уравнению системы (12) в столбце т, то есть для х(т) < х < х(т+1) и у(1) < у < у(м+1), причем в данной области все эти функции являются непрерывными вместе с первыми производными. Будем называть и(^п)(х) и фкт)(у) локальными модами, а индекс к определяет номер локальной моды и напрямую связан с введенной соотношением (11) величиной кк , их взаимосвязь будет показана далее.
Вычисление локальных синусоидальных мод Рассмотрим уравнение
ф + (к02е(x,у) - к2)ф= 0.
(16)
Так как е(х, у) = е(т)(у), то (16) можно решать в каждом столбце отдельно, независимо от других столбцов. Итак, в некотором столбце (индекс т опускаем) (16) примет вид:
ф + (к02е(у) - к2)ф= 0 .
(17)
здесь е(у) = 8(n) = const для у(n) < у < у(n+!). Решение (17) в n - ой строке между у(n) и у (n+1) будет иметь вид:
ф( у) =
ФSn,'- у(n))] + ^(-sm[к{"\у - у(n))], (18)
к (П) "у
Лу
где куп =^8(n)к2 - к2
(19)
фП) = ф(у(п) + 0) (20а)
- нижнее или левое значение функции ф(у) в строке п ;
ф^) =ф (у(п) + 0) (20б)
- нижнее или левое значение производной ф (у) в строке п .
Следует отметить, что поскольку куп) может
принимать только действительные или чисто мнимые значения, то функция ф(у), определяемая (18), будет вещественной для любых действительных ф <",!) и ф (ап,г). Аналогично (20а) и (20б) можно оп-
и
ределить правые значения ф(у) и ее производной ф(у) в строке п :
Ф [п,г) =ф(у(п+1) - 0), (20в)
Ф(а,Г) =Ф(у(п+1) -0). (20г)
Иллюстрация расположения ф [п,1), Ф[п'г), фаП'г) приводится на рисунке 2.
ф
(п,1)
У,
(N+1) У
ут 1 № Ф, < г Ж №Т>
II гф.
р) У Л М Ф* <
г -к <''т> ф,
т У
Рис. 2. Схема обозначения полей на границах разбиения сечения волновода Таким образом, чтобы решить (17) нужно определить значение кк и 2N - 2 - констант ф [п'г) и ф ап1), так чтобы функция ф( у) удовлетворяла (17) и
краевым условиям:
Гф (у(1)) = 0
|ф( у(1)) = 0
или
[ф (у(N+1)) = 0:
[ф( у(N+1)) = 0 соответственно для электрических
или магнитных стенок.
Требование непрерывности решения (17) и его производной приводит к следующим соотношениям:
ф[п,г) =ф [п+и), (21а)
ф!Г) =ф (Ги). (21 б)
Кроме того, для ф[п,'), ф(п1), ф[п'г), ф^) справедливы равенства:
ф[п'г) = ф(п')0С8[куп)й(п)] + 1^—8т[к<п)й(п)], (22а)
к(п)
«у
ф^ = ф^ соБкуп)й(п)]-ф^п япкуп)й(п)], (22б)
ф(п,Г )
ф(п-') = ф^) со8[куп)й(п)] -фапг81п[к^п) й(п)], (23а) ку
ф(п,1) = фС^ соък^й(п)] +ф[п,')к^ч) ыпк/й(п)]. (23 б)
Введем в рассмотрение матрицы, детерминант которых равен единице
р (;) =
е(° =
со^й(г)]
¡ип[ку°й(г)]/ ку
- зш[к уг)й (г)]к уг)
С08
[к«й (;)]
С08
[куг)й(г)] - ¡ап[ку0й(0]/ку'-[к^й(;)]
Яп[к«й (г)]к
С08
(24)
(25)
Очевидна справедливость равенств
ф[п,г) ф(ап,Г )
= р(п)
ф[пЛ ф(п,0
ф(п,') ф(;,')
=о(п)
ф[п,г) фап,г)
которые равносильны (22) и (23). Используя матрицы Р(;) и ((;), можно легко выражать значения
ф(;,г )' фУ )
через
фГ ) ф^)
рез
ф ™
ф аг)
ф [п,,г)' ф (п,г).
ф[п'+1,')" ф(.+1,')
; < у , и значения ; > у . Так, например,
ф Г'
ф у)
че-
= р(п,) р (п'-1)
••• Р1
ф [и) ф ?)
= ( (п'+1)( (п'+2) # (
( N )
ф^" )
фГ )
(26) (27)
Для вычисления кк, точнее квадрата кк , поскольку только квадрат фигурирует в приведенных выше формулах, воспользуемся следующим методом. Определим функцию
Д(п,)(кк2) = ф[п'г)фГМ) -ф(п,Г)фГ"), (28)
где п'е (1, N), значения ф
(п,г) ф(п,г)
фап,г) и ф
(п'+1,' )
ф
(п'+и )
вычисляются с использованием (26) и (27) и
потому зависят от к^ . Нули функции Д(п')(к1) и есть искомые значения. Поиск нулей (28) равносилен решению характеристического уравнения
ф Г) ф Ги)-ф (/,г) ф Ги) = 0,
(29)
причем п' может быть произвольным числом из интервала (1, N).
Можно убедиться, что зависимость решения
(29) от ненулевых значений
ф[и) ф(«и)
ф )
фГ)
от-
сутствует, поэтому, решая (29), в случае электриче-
ских
стенок
можно
считать
"ф^) ■ "0"
_фГ) _ 1
"ф^) ■ "1"
_фГ) _ 0
а в случае магнитных
"ф[и)" "0"
.ф(аи) _ 1
"1"
_ф(аи) _ 0
и
и
Различных значений кк , удовлетворяющих (29), бесконечно много, но все они находятся на интервале (-да, к^^ах], где их можно отсортировать по убыванию. Каждое из значений к2 определяет локальную моду ф кт)(у) в т - м столбце сечения
волновода. Точность построения поля по формуле (15) зависит от количества используемых локальных мод. Но как бы много мы их не взяли (как бы много корней уравнения (29) не нашли), это всегда будет некоторое конечное число К .
Когда все К значений кк найдены, нас интересуют константы ф^), фП), ф(",г), ф^), их можно найти из условия нормировки функции ф( у) на единицу: 1
у( N+1)
у(N+1)
Здесь
у(1) у(1)
{ф 2( у)йу = 1.
(30)
N У
I = { ф 2(у)ау = Х {ф 2Шу =
у(1)
1 , , (п I) I (и
2(к Г)2
с1(") ф'"-1)
) ф(".г)) + ¿((ф^ ))2 + (-фка___)2)].
2 ку
Когда кк2 известны, значение к(у") для каждой
строки может быть рассчитано по формуле (19). Принимая во внимание соотношения
ф ) ф (и,г)
= р("'')р("'-1) ... р1
ф 5и)
ф аи)
ф5",1)■ л ("- -1,г) "
ф (У) _ л ("фа -1,г )
" = 2, N ,
" = 1, N, (31)
(32)
значение интеграла I можно считать функцией от
фХ') или фа в зависимости от типа стенок. Тогда (30) решается стандартными методами и имеет ровно два действительных корня, отличающихся только знаком, любой из них может быть выбран в качестве значения ф(,1'1) (или фа'1"1), остальные константы рассчитываются с использованием (31), (32) и найденного значения.
Действуя по вышеописанному алгоритму, для каждого столбца сечения волновода можно найти К значений кк2 и для каждого из них сконструировать функцию у) , которая будет удовлетворять (17) при данном кк .
Сшивка локальных мод
Рассмотрим уравнение
икт)(X) + (ккт)2 -к2)«кт)(х) = 0
(33)
здесь мкт)(х) является к - ой модой в столбце т, соответствующей к((п) также как фкт)(у). Уравне-
ние (33) аналогично (17), соответственно его решение имеет вид
4т)( х) = ^[к^Ч X - X(т))] +
и (т,1)
-кктг51п[кХки)( X - х(т))];
кхк
где
кХт Чккт)2 - к2
и(т,1) = мкт)( х(т) + 0),
= г^кт)( х(т) + 0). Аналогичным образом вводятся
и
= и (т)( х(т+1)
= икт)(х(т+1) - 0),
и
(т,г) = и1<кт')(х(т+1) - 0).
(34)
(35)
(36)
(37)
(38)
(39)
Разница между уравнениями (33) и (17) состоит в том, что решением (33) является не одна непрерывная функция вида (34), а набор из М определенных и непрерывных на соответствующих интервалах функций. Формула (34) объединяет их в одну,
л. „.( т,1) т,г)
но разрывную функцию, поэтому для ик , ихк и
.(^) и(т,г) ак > ак
и)Г", и)Г' соотношения, аналогичные (21), не
выполняются, тем не менее, справедливы равенства
(т) и (т,1) к (т)и (т'г)
,, (т,1) 1ак
(т ,г)
кут) и ( кхк ик
кут )и ( кхк ик
1ап( ккт) d(хт)) ет* ккхт) ))
(40а)
к (т)и (т,1)
кхк ик
к (т)и (т,г)
кхк и*к (40б)
81И( ккхт) ахт)) 1аи( ккхт) ахт))
Поскольку (40) выполняются для всех к = 1, К удобно использовать матричную запись. Введем диагональные матрицы Т(т) и £(т) размерности К х К с диагональными элементами
Ткт = кхт)/1аи(ккхт) ахт)), (41а) ^ = к^/яп^"0ахт)). (41 б) Теперь (40)
(42а) (42б)
примут вид
тт(m,I) а = -Т (т)и (т,1) + ^ (т )тт(т,г)
тт(т,г) а = -Б (т)П (m,I) + т (т)-ц5т,г)
Здесь
ГТ (m,I) а =[иат,1) и (т,1) а2 , (т,1)]Т маК ]
тт(т,г) а =[иат,г) и (т,г ) а2 и(т,г)]т иаК ]
тт(m,I) = [иГ) и (т,1) и5 2 , 5т,|)]Т
тт(т,г ) =[и5т,г) и (m, г ) и$2 , (m, г)]Т
Рассмотрим теперь интеграл перекрытия
у (Ы+1)
</:!/ >=
1
y(N+1) - у(1) у(1)
{ /:(у) /2 (у)ау, (43)
и
"=1
составим матрицу из элементов вида (43), где в качестве / и f2 будут использованы ранее найденные
моды фкт)( у)
о^" =<фкт)| ф(рт,) >,
(44)
Ги™ = [0 0 - 0]Т , [и^) = [0 0 ... 0]Т.
Если стенки магнитные, то матрица Т(1) должна быть заменена на Т(1)* с диагональными элемен-
получится квадратная матрица 0(т,т) размера тами Т^* = -к^к ^п^Хк^Х1-1), а матрица Т(М) -К х К , обладающая очевидным свойством
(0(т,т1))Т = о(т',т) (45)
Кроме того, так как ортогональные функции
хк ) , а матрица т — на
матрицу Т(М)*, Т(М)* = -кХМ) ЪЩкХМ^Х?)).
Задачу, описываемую соотношениями (49), в общем виде можно представить как
фк )(у) нормированы, то матрица 0( , ) равна еди
т т')
ничной, и каждая из матриц О1 ' ' является унитарной.
0(
(т,т')о(т',т) =
1
(46)
Для конечных значений размерности К х К соотношение (46) выполняется довольно точно только для соседних столбцов, то есть если т' = т +1 или т' = т -1. Из ортонормированности функций фкт)( У) и требования непрерывности поля х, у) (вместе с производными первого порядка) следует выполнение условий
и
(т,к) _
А'к
= 10
(т,т+1) (т+1,1)
р=1
кр
и
(т,к) _
= 10
р=1
(т,т+1) (т+1,1)
кр
ар
Или в матричной записи
и (т, г ) = о(т,т+1)ц-(т+1,1) ц-(т, г ) = о(т,т+1)ц-(т+1,1)
(47а) (47б)
(48а) (48б)
Объединяя (42) и (48) и исключая и^1),
и(ат,г), ит), получим
(0(m, т-1)Т(т-1)0(т-1,т) + Т (т))ит) = О(т,т-1) 5(т-1)и (т-1,1) + 5 (т)о(т,т+1)и (т+1,1) (49а) для 2 < т < М ;
(0(2,1)Т (1)0(1Д) + Т (2))и^) = 5 (2)0(2,3)и^) (49б) для т = 2, когда М > 2 ;
(0(М, м-1)т (м-1)О(М-1,М ) + т (М ))и ) =
0(М^М-1)5(М-1)и(М-1,1) (49в)
для т = М .
В простейшем случае, когда М = 2 , получаем
(0(1,2)т(2) + Т(1)0(1,2))иР'1) = 0 . (49г)
Соотношения (49) выполняются в случае электрических стенок, когда краевые условия имеют вид
Л(кг )и = 0,
(50)
где матрица Л(к2) состоит из (М -1) х (М -1) блоков, размерности К х К каждый. Структура этой матрицы имеет вид:
А(2) с(2) О О О О О О
В(3) А(3) с(3) О О О О О
О В(4) А(4) с(4) О О О О
О О О О О в(М-1) А(М-1) С(м-1)
О О О О О О В(М) А(М)
Л(кг) =
Здесь
д(т) = о(т,т-1)т (т-1)о(т-1,т) + т (т)
В (т) = —о(т,т-1)5 (т-1) С (т) = —5 (т)о (т,т-1)
О - нулевая матрица размерности К х К . В (50) вектор
' и(2,') '
(51)
(52)
(53)
(54)
и =
и р
и.
(М ,1)
(55)
Мы получили задачу (50), которая имеет только тривиальное решение и = 0 , если det(Л(kz)) Ф 0. Так как матрица Л(к2) зависит от параметра, мы можем подобрать значение к2 таким образом, что хотя бы одно из собственных чисел Л(к2) будет нулевым, тогда, полагая и равным собственному вектору, соответствующему нулевому собственному числу, мы удовлетворим (50). Для поиска таких значений кг предлагается следующий метод.
Возьмем произвольный вектор V с ненулевыми компонентами, например, единичный. Будем решать уравнение
Л(к2 )и' = V (56)
относительно и'. Очевидно, для разных значений параметра кг будем получать разные решения уравнения (56). Определим функцию f(кг) = 1/и' , где и' есть р - ая компонента вектора и'. В окре-
1
1
0
0
1
стности искомых значений кг функция /(кг) является непрерывной функцией скалярного аргумента, и ее нули есть искомые к2.
Поиск нулей /(кг) осуществляется стандартными методами. Для каждого из найденных кг матрица Л(кг) становится числовой, определив собственный вектор и , соответствующий ее нулевому собственному числу, мы получаем некоторые из значений констант, необходимых для конструирования К функций вида (34). Остальные можно найти, воспользовавшись соотношениями (42) и (48), в следующем порядке.
Для т = М
и((М ,1) известно как часть и ;
1) Если стенки электрические, то
и(М,г) = [0 0 ••• 0]Т , тогда
и'М,1) =-Т(М)и(м,1) и и'М,г) (М)и(м,1);
2) Если стенки магнитные, то
и(М,г) = [0 0 ••• 0]Т
а L ,тогда
иМ ,1) = (т (М ))-1(-(т (М) )2 + (^ (М ))2)и(м ,1)
и(м,г) = (5(М))-1(иаМ,1) + Т(М)иМ)).
Для 1 < т < М
1) и|-т,г) = о(т,т+1)и (т+1,1);
2) и(т,г) = 0(т,т+1)и(т+1,1);
3) и(т,') = (Б (т))-1(-и(т,г) + Т (т)и,(т,г))
этот пункт выполняется только для т = 1 , так как для других значений т вектор U(sm,l) известен, как часть собственного вектора и ;
4) и(ат,') = -Т(т)и(т,1) + Б(т)и(т,г).
Таким образом, для каждого значения кг можно сконструировать К мод вида (34), объединив их с К модами вида (18) (они одни и те же для разных кг) по формуле (15), получаем поле в сечении волновода, распространяющееся в направлении оси г , с проекцией кг волнового вектора на эту ось.
2. Численные результаты
Ниже приводятся результаты расчетов, проведенных при помощи математического пакета Matlab 6.0, для волновода, сечение которого изображено на рис. 1.
Как следует из рис. 1, сечение волновода было разбито на М = 9 столбцов и N = 9 строк. В качестве граничного условия использовалось условие равенства нулю поля на краях сечения. Вычисления проводились для излучения с длиной волны в вакууме Х0 = 1,3 мкм. Для каждого столбца сечения было получено по К = 40 значений кк(т)2, удовлетворяющих характеристическому уравнению вида (29) (таблица 1). Для каж-
дого из этих значений была построена локальная мода Фкт)(У). Различных решений задачи (50) в данном случае оказалось более ста. Распределения интенсивности в сечении первых нескольких мод с наибольшими значениями параметра кг показаны на рис. 3.
Табл. 1. Корни характеристического уравнения (29) для трех различных типов столбцов сечения волновода
m k 1,3,7,9 2,4,6,8 5
1 64,2383 69,9582 69,1702
2 64,2343 69,5926 66,4472
3 64,2292 68,9834 64,2318
4 64,2251 68,1305 64,2317
5 47,6200 67,0339 61,9325
6 47,5711 65,6935 55,6687
7 47,5101 64,1095 47,7458
8 47,4602 62,2818 47,5406
9 25,9664 60,2104 47,5308
10 25,3814 57,8953 38,3052
11 24,5208 55,3365 27,9266
12 23,6360 52,5341 25,0060
13 15,6980 49,4879 24,8084
14 14,5517 46,1980 18,2359
15 12,5823 42,6645 15,0959
16 10,2083 38,8872 13,4231
17 8,0763 34,8663 10,9748
18 -0,6847 5,6230 1,3531
19 -3,7939 -0,1038 -3,9389
20 -7,7667 -6,0743 -7,1449
21 -12,2263 -12,2885 -12,4743
22 -16,8233 -18,7464 -17,0784
23 -22,1979 -25,4480 -22,5124
24 -27,5639 -32,3933 -28,6909
25 -33,0938 -39,5822 -33,4829
26 -38,4214 -47,0149 -40,7911
27 -47,9161 -54,6913 -49,0486
28 -53,3655 -62,6113 -55,0483
29 -59,8437 -70,7751 -61,6515
30 -66,6876 -79,1825 -68,6689
31 -73,1481 -87,8336 -75,7737
32 -81,8960 -96,7285 -84,6344
33 -89,5239 -105,8670 -92,5589
34 -97,4445 -115,2492 -100,1193
35 -105,4019 -124,8751 -109,8932
36 -115,8679 -134,7447 -119,1908
37 -123,7581 -144,8580 -127,5402
38 -132,6166 -155,2150 -136,7245
39 -141,7554 -165,8157 -145,7534
40 -150,5108 -176,6601 -155,6022
Из рис. 3 видно, что основная мода данного волновода, выполненного в виде фотонного кристалла, т. е. в виде периодически расположенных участков с высоким и низким показателями преломления, сосредоточена в осевом сердечнике размером 3 х 3 мкм.
123456789 X, мкм
б
123456789 X, мкм
в
Рис. 3. Распределение поперечной интенсивности
первых четырех мод неоднородного сечения волновода (рис. 1) с константой распространения: kI = 8,2786 мкм-1(а), kI = 8,1681 мкм-1(б), kz = 8,1513 мкм-1(в), ^ = 8,1373 мкм-1(г)
На рис. 4 показана в процентах доля энергии основной моды (рис. 3а), приходящаяся на ту или иную однородную область сечения, при этом в осевом сердечнике сосредоточено до 97%.
7 6 5 4 3 2 1
О
0 0 0 0 0 0 0 0 0
0 0 0 0,02 0,01 0,02 0 0 0
0 0 0,01 0,38 0,07 0,38 0,01 0 0
0 0,02 0,45 6,14 10,76 5,72 0,38 0,02 0
0 0,01 0,12 11,86 26,17 10,76 0,07 0,01 0
0 0,02 0,46 6,60 11,86 6,13 0,38 0,02 0
0 0 0,02 0,46 0,12 0,45 0,01 0 0
0 0 0 0,02 0,01 0,02 0 0 0
0 0 0 0 0 0 0 0 0
4
8 9 X, мкм
Рис. 4. Распределение энергии основной моды волновода (данные указаны в процентах) Очевидно, несимметричность распределения на рис. 4 связана с погрешностями в вычислении.
Поскольку данный метод является приближенным, возникает вопрос о его сходимости. В частности о том, какое количество локальных мод нужно использовать для построения поля по формуле (15), чтобы результат был более или менее точным. То есть насколько большим должно быть число K, чтобы результат и сложность его получения были приемлемы. На рис. 5 показана зависимость получаемого значения константы распространения kz
основной моды волновода (рис. 1) от количества локальных мод в столбцах K. Из рис. 5 видно, что при К>20 Kz рассчитывается с точностью до 4-го знака после после запятой. кг
8,2900 8,2880 8,2860 8,2840 8,2820 8,2800 8,2780
О 10 20 30 40 К
Рис. 5. Зависимость значения константы распространения kz (мкм.'1) основной моды от количества локальных мод в столбцах К
Заключение В работе получены следующие результаты. Предложена модернизация известного метода согласованных синусоидальных мод, применяемого для скалярного расчета модового состава неоднородных оптических волноводов, основанная на матричной форме записи характеристического уравнения.
Модернизированный метод согласованных синусоидальных мод реализован с помощью программного пакета МаНаЪ 6.0 и применен для расчета мод оптического волновода, выполненного в виде фотонного кристалла;
г
• показано, что почти вся энергия основной моды данного волновода сосредоточена в осевом сердечнике однородного сечения и не попадает в оболочку, которая представляет собой периодически расположенные однородные области с высоким и низким показателем преломления.
Благодарности
Работа поддержана российско-американской программой «Фундаментальные исследования и высшее образование» («BRHE»), а также президентским грантом РФ НШ-1007.2003.01.
Литература
1. Itoh T. Numerical technique for microwave and millimeter-wave passive structures // Wiley, New York, 1988.
2. Sudbo A.S. Film mode matching: a versatile method for mode film calculations in dielectric waveguides // Pure Appl. Opt., 1993. V. 2. P. 211-233.
3. Sudbo A.S. Improved formulation of the film mode matching method for mode film calculations in dielectric waveguides // Pure Appl. Opt., 1994. V. 3. P. 381-388.
4. Rogge U., Pregla R. Method of lines for the analysis of dielectric waveguides // J. Lightwave Techn., 1993. V. 11. P. 2015-2020.
5. Sztefka G., Nogling H.P. Bidirectional eigenmode propagation for large refractive index steps // IEEE Photonics Techn. Left., 1993. V. 5. P. 554-557.
6. Rahman B.M.A., Davies J.B. Finite-Elements Solution of integrated optical waveguides // J. Lightwave Techn., 1984. V. 2. P. 682-687.
7. Koshiba M., Maruyama S., Hirayama K. A vector finite element method with the high-order mixed-interpolation-type triangular elements for optical waveguiding problems // J. Lightwave Techn., 1994. V. 12. N. 3. P. 495-502.
8. Lusse P., Stuwe P., Schule J., Unger H. Analysis of vectorial mode fields in optical waveguides by a new finite difference method // J. Lighhtwave Techn., 1994. V. 12. N. 3. P. 487-493.
9. Hadley G.R., Smith R.E. Full-vector waveguide modeling using an iterative finite-difference method with transparent boundary conditions // J. Lighhtwave Techn., 1995. V. 13. N. 3. P. 465-469.
10. Lin P.-L., Li B.-J. Semivectorial Helmholtz beam propagations by Lanczos reduction // IEEE J. Quant. Electr., 1993. V. 29. N. 8. P. 2385-2389.
11. Lee P.-C., Voges E. Three dimensional semi-vectorial wide-angle beam propagation method // J. Lighhtwave Techn., 1994. V. 12. N. 2. P. 215-224.
