УДК 537.876.4
ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ИЗГИБА ПЛАНАРНЫХ БРЭГГОВСКИХ ВОЛНОВОДОВ
Д. В. Прокопович1, А. В. Попов, А. В. Виноградов
Найдена аналитическая форма и осуществлена численная аппроксимация граничных условий прозрачности для параболического волнового уравнения в криволинейных координатах. Продемонстрировано соответствие решений, получаемых методом параболического уравнения, решениям спектральной задачи, определяющей моды брэг-говского волновода. Выполнено численное моделирование амплитуды поля и изгибных потерь в зависимости от радиуса кривизны и параметров брэгговского волновода.
1. Введение. В последнее время в оптике интенсивно развивается новое научное направление - фотонные кристаллы [1]. Фотонные кристаллы - это структуры с периодически изменяющейся диэлектрической проницаемостью. Периодичность приводит к дисперсии распространяющихся волн, в частности, к появлению фотонных запрещенных зон. Протяженные дефекты (отклонения от симметрии) в фотонных кристаллах приводят к образованию волноведущих структур - световодов на основе фотонных кристаллов. Механизм распространения излучения в таких световодах принципиально отличается от полного внутреннего отражения и связан с интерференцией, возникающей при отражении от периодической оболочки, окружающей сердцевину (дефект). Простейшими одномерными фотонными кристаллами являются цилиндрический брэг-говский световод и планарный брэгговский волновод в случае плоской симметрии.
Впервые идея использования резонансного отражения в волоконных световодах была предложена в 1978 г. [2]. В дальнейшем широкое распространение брэгговские световоды получили благодаря различным применениям [3, 4], которые связаны с их необычными свойствами по сравнению с традиционными диэлектрическими волокнами. Хотя в идеализированных моделях могут быть получены очень малые потери для основной моды брэгговского световода [5], при практическом изготовлении и использовании возникают проблемы устойчивости к малым изменениям: неточностям изготовления профиля, локальным неоднородностям и т.д. Наиболее часто встречающееся отклонение от идеальной формы - это изгиб световода, приводящий к искажению поля и дополнительным потерям [6]. Влияние изгиба световодов начали исследовать на заре волоконной оптики. Точные методы решения этой задачи довольно сложны. Тем не менее, для обычных оптических волокон получены надежные асимптотические оценки [6]. Изгиб брэгговских световодов изучен гораздо меньше, получены лишь качественные оценки дополнительных потерь на излучение [7, 8]. Большинство численных методов основано на использовании модели эффективного показателя преломления, в которой к исходному профилю добавляется линейная поправка, учитывающая кривизну [6]. Пока не изученными для брэгговских световодов с изгибом остаются интересные эффекты, такие как изменение структуры поля вдоль волокна, интерференция волн в периодической оболочке, межмодовое перетекание энергии при уменьшении радиуса кривизны.
Перспективы широкого применения брэгговских световодов и необходимость точно моделировать их оптические свойства при изгибе определили тему статьи. В данной работе развивается численный метод моделирования процессов установления поля и вытекания излучения в изогнутых брэгговских волноводах на основе параболического уравнения.
2. Метод параболического уравнения. Метод параболического уравнения (МПУ) является мощным вычислительным средством современной теории дифракции. Он был предложен Леонтовичем и Фоком для задач распространения радиоволн вдоль поверхности Земли [9]. Параболическое волновое уравнение аналогично нестационарному уравнению Шредингера в квантовой механике [10]. МПУ широко используется для вычисления волновых полей в сложных реальных сред [11, 12].
Поясним более подробно приближение МПУ. Пусть электромагнитная волна, распространяющаяся преимущественно вдоль оси г, имеет только одну компоненту электрического поля Еу = Е(х, г). Из уравнений Максвелла получаем волновое уравнение
для поля Е:
д2Е д2Е j2 „ п
Мы будем рассматривать среды с малым контрастом показателя преломления: п2(х) = ñ2[l + a(x,z)], где |а| <С 1 и а —► 0 при х —> оо. При а = 0 уравнение (1) имеет точное решение в виде плоской волны Е(х, z) = А ехр[г&п(х sin в + z cos 0)], где А - амплитуда, в - угол между волновым вектором к (к = 2тг/\) и осью z. Если волна распространяется под небольшим углом в <С 1 к оси z, то решение уравнения (1) выглядит как E(x,z) и Aexp(¿&ñz)exp[ífcñ(x0 — z02/2) = exp(ikñz)u(x, z). Легко видеть, что плавно меняющаяся волновая амплитуда u(x,z) удовлетворяет параболическому уравнению „., _ ди д2и
2гкп——|- —— = 0, справедливому как для отдельной плоской волны с малым углом oz ох2
скольжения 0, так и для произвольного пакета таких волн. В случае малого контраста методом двухмасштабных разложений [13, 14] получается более общее параболическое уравнение
ди д2и ,о9,ч /„ч
2гЬг— + — + к2п2а(х, г) = 0. (2)
Искомая амплитуда волны u(x,z) медленно меняется на расстояниях порядка А и связана с напряженностью электрического поля соотношением
Е = и(х, z) exp(ikz), к = kñ. (3)
Основными преимуществами приближения параболического уравнения являются:
1. Однонаправленность волнового процесса: вместо краевой задачи решается задача Коши с условием на начальном фронте и|г=0-
2. Рассмотрение медленно меняющейся комплексной волновой амплитуды u(x,z): подстановка (3) избавляет от основного осциллирующего множителя exp(ikz). Поэтому МПУ позволяет моделировать пространственное распределение поля в протяженных областях с небольшими вычислительными затратами.
3. Возможность аппроксимировать параболическое уравнение с помощью абсолютно устойчивой разностной схемы. Этими качествами определяется выбор МПУ для описания брэгговских волноводов. Поскольку основными модами брэгговского волновода являются вытекающие моды, выбранный вычислительный метод должен корректно описывать эффекты их излучения в свободное пространство. Важной особенностью МПУ является возможность постановки точных условий прозрачности на границе моделируемой области, где среда считается однородной. Такие условия получаются из сшивки
х
\ \
" \
-В
\ \ \ \ \ V \ \ \
р
Рис. 1. Схематическое изображение изогнутого планарного волновода.
решения параболического уравнения на границе вычислительной области с аналитическим решением в однородной среде [15]. В данной работе развивается этот метод в применении к задачам установления и вытекания мод в брэгговских волноводах с изгибом.
3. Аппроксимация граничных условий прозрачности. Согласно [16], параболическое уравнение в криволинейных координатах (х, г), где г - длина дуги кривой с кривизной X (радиус кривизны р = 1 /х), а ж - расстояние по нормали к этой кривой (рис. 1), записывается как
где к = кп. Брэгговский волновод, заданный показателем преломления п = п(х), ограничен в поперечном направлении |ж| < В. При |ж| > В будем считать среду однородной с п — п. Граничное условие прозрачности, описывающее свободное вытекание волн для дальней по отношению к центру кривизны границе х = В было найдено в [17]:
+ 2к -Аг1- l + 2Xx и = 0 nz
(4)
¿4/3,0 = д dv di
о
1
2тг7
—оо+гс
4 6 8 10
г, шш
Рис. 2. Амплитуда поля в брэгговском волноводе при А = 1 мкм. Рис. 3. Амплитуда поля в брэгговском волноводе при А = 0.5 мкм.
где введены безразмерные координаты £ = (кх2/2)1/3г, т/ = (2к2хУ^х и параметр /3 = (2к2Ху/3В. Граничное условие на вогнутой границе х = —В имеет вид
0 — оо+«с
здесь = = Л/5г[А1(<) + гВ1(<)] - функции Эйри-Фока [16, 18]. Точ-
ные аналитические формулы (5), (6) справедливы при любых радиусах кривизны, даже при предельно малых, когда нарушается сама применимость МПУ, однако они сложны для численной реализации. Получим более простые граничные условия прозрачности, справедливые в приближении малой кривизны \р\ В. Рассмотрим интегральное ядро входящее в граничное условие (5). Корни функции лежат на луче ехр(г7г/3)
[18], они же и точка £ = 0 являются полюсами подынтегрального выражения для регулярного в нижней полуплоскости. Сдвигая контур интегрирования на —гос и используя приближенные формулы шг(г) ~ ехр(2г3/2/3)/*1/4, ы'^) ~ *1/4 ехр(2г3/2/3) [18],
+оо-гс _^
получаем асимптотику для ядра А'+(£) «--/ ехр(г££)\Л — Р—■ Этот инте-
2ттг J * £
— 00+1С
грал разбивается на вычет в точке ! = 0 и интеграл по контуру, окружающему луч атg(t — /3) — 7г/3. Опуская выкладки, приводим конечный результат: +(£) =
ехр(—¿7Г/4)
Vñ
ext
. Таким образом, упрощенное граничное
I ехр(гй2)с?5 о
условие прозрачности (5) на выпуклой границе вычислительной области имеет вид
ди(В,г) [ ди(В, г — г')
dz
dz
-Q+(z')dz',
2 k
Q+(z') = — exp(—¿7r/4)y -
7TZ
_ Vx'kBz'
p(iXkBz') - 2i\JXkBz' j exp (it2)dt
(7)
Условие для вогнутой границы получается из выражений (7) заменой В —> —В. При ХВ —> 0 оба они сводятся к известному результату [15] для прямолинейной границы.
Разностная аппроксимация параболического уравнения (4) выполнялась с помощью абсолютно устойчивой неявной шеститочечной схемы Кранка-Николсона [19], затем система разностных уравнений решалась методом прогонки.
4■ Установление моды и потери в брэгговском волноводе. Рассмотрим прямолинейный плоский брэгговский волновод с тремя периодическими слоями; сердцевина волновода имеет ширину 2а = 20 мкм и показатель преломления п0 = ñ = 1.449. Периодические слои 1\ = 1.19 мкм, щ = ñ + 0.015 и /2 = 10.0 мкм, n2 = ñ оптимизированы под длину волны А = 1 мкм [20]. Волновод окружен однородной средой с показателем преломления п. Начальная амплитуда пучка выбрана в виде косинуса u|2=o = cos(7rx/2а) при |х| < а и u|z=o = 0 при |х| > а. На рис. 2, 3 показано рассчитанное численно распределение амплитуды поля \u(x,z)\ в продольном сечении волновода (рисунки сжаты в продольном направлении в 100 раз). Эти рисунки наглядно показывают роль интерференционных эффектов в формировании волноводной моды. Несмотря на малый оптический контраст, отражение скользящих волн от плотных слоев оболочки формирует в сердцевине волновода распространяющуюся моду. Минимум радиационных потерь достигается при А = 1 мкм (рис. 2), когда брэгговские слои имеют эффективную оптическую толщину А/4 (где А - "поперечная длина волны"). При удалении от оптимума потери растут, и при А = 0.5 мкм брэгговские слои уже имеют оптическую толщину, кратную А/2, и поле полностью высвечивается через периодическую оболочку (рис. 3).
2, ЦЛ1-10"
Рис. 4. Зависимость потерь в брэгговском волноводе от расстояния вдоль волновода г, штриховая линия соответствует потерям, полученным конечно-разностным методом, А = 1 мкм.
в
Несложно найти полный поток энергии через сечение z = const : P(z) = J \u(x, z)\2dx,
-в
R" 1 dP
и вычислить производную p =--стремящуюся в пределе к потерям уста-
2Р az
г ж л ГдБ1 20 •109
новившеися моды. I рафик потерь, точнее, А
км
log 10
/3"[мкм ] в зависимо-
сти от расстояния г представлен на рис. 4. Более детально изображены потери при 2 ~ 7,44,149 см. Для малых г ~ 0 — 10 см количество возбужденных в волноводе мод велико, а так как каждая имеет свою постоянную распространения /3 = /?' 4- г/3", то амплитуда результирующей волны \и(х,г)\ осциллирует и, как следствие, потери тоже. При г ~ 30 — 50 см интерферируют две слабозатухающие моды - колебания потерь имеют вид гармонической функции; далее при увеличении 2 амплитуда осцилляций потерь уменьшается - остается практически одна ТЕ\ мода с наименьшим затуханием. Потери, вычисленные решением точной спектральной задачи [20] для ТЕ\ моды на длине волны А = 1 мкм: А = 41.37 дБ/км, а с помощью МПУ: А = 41.5 ± 0.3 дБ/км. Это демонстрирует хорошую точность МПУ и возможность применения метода для оценки потерь в изогнутых брэгговских волноводах.
5. Моделирование брэгговского волновода с изгибом. Физический механизм увеличения потерь при изгибе ступенчатых диэлектрических волноводов известен давно [6]. В таком волокне захваченная мода аналогична квантово-механической частице в прямоугольной потенциальной яме, для нее выполняется условие п0болочки 5: ,/5'/к < тгСердцевины-Когда волновод изогнут, то к исходному ступенчатому показателю преломления добавляется линейная поправка, пропорциональная кривизне [6]. При увеличении расстояния от центра кривизны возникает ситуация, при котором пЭфф.об0лочки > Р'/к, таким образом мода начинает вытекать в оболочку волновода за счет туннельного эффекта. Для брэгговских волноводов также существует нижняя граница эффективного модового показателя преломления ¡3'/к, только она связана с шириной фотонной запрещенной зоны [21]. При изгибе брэгговского волновода зонная структура изменяется так, что мода может выйти из запрещенной зоны, пЭфф.заПр.зОНЫ > /З'/к, что приводит к резкому уменьшению коэффициента отражения брэгговской оболочки и увеличению излучательных потерь. Мы используем МПУ для количественного расчета этого эффекта.
50
г, шш х, цт
Рис. 5. Амплитуда поля в брэгговском волноводе с радиусом изгиба р = 8 см, А = 1 мкм.
Рис. 6. Поперечное распределение амплитуды и(х,г) в брэгговском волноводе с радиусом изгиба р = 8 см, А = 1 мкм.
На рис. 5 представлено вычисленное с помощью МПУ распределение амплитуды поля и(х, г) в брэгговском волноводе с радиусом изгиба /9 = 8 см. Интенсивное вытекание излучения происходит через верхнюю по отношению к центру кривизны оболочку. Это хорошо видно по асимметричному поперечному распределению поля (рис. 6).
При уменьшении радиуса изгиба волновода возрастают потери. Как показано в работе [21], критический радиус изгиба, т.е. радиус, при котором значительно увеличиваются излучательные потери, можно оценить по формуле рс = 4тгАп2/ |(/ЗД)2
~ пэфф.запр.зоны • Интересно сравнить эту приближенную формулу с нашими численными результатами. По формуле получаем, что для ТЕ\ моды (для волновода из раздела 4) на длине волны А = 1 мкм критический радиус рс ~ 50 см. Это качественно согласуется с рассчитанной зависимостью потерь от радиуса изгиба (рис. 7).
Рис. 7. Зависимости потерь в брэгговском волноводе от радиуса изгиба р, штриховая линия соответствует потерям, полученным конечно-разностным методом, А = 1 мкм.
Рис. 8. Зависимости потерь в брэгговском волноводе от расстояния вдоль волновода z при различных радиусах изгиба, штриховая линия соответствует потерям, полученным конечно-разностным методом, А = 1 мкм.
Действительно, при р рс потери почти не зависят от р: А{р — 500 см) = 41.8 дБ/км, А(р = 300 см) = 42 дБ/км, тогда как при р < рс потери резко увеличиваются по сравнению с прямым волноводом: А(р = 30 см) = 89.3 дБ/км, А(р = 20 см) = 179 дБ/км. Если еще увеличить изгиб (рис. 8), то при р ~ 7 см произойдет срыв основной ТЕ\ моды и в волноводе останутся (с малой амплитудой) только захваченные оболочкой моды, локализованные в слоях с повышенным показателем преломления (рис. 9, 10). Процесс такой трансформации мод при увеличении изгиба брэгговского
<я
50 0
-50 -100
N„ -150
vit а -200
-250
п(х)
( 1 1 1 I 1 [j ! !
V Л u(x,z) -----1-.-j-L
-50 -40 -30 -20 -10 0 10 20 30 40 50
х, цт
-50 -30 -10 0 10 20 30 40 50
х, цт
Рис. 9. Поперечное распределение амплитуды u(x,z) в брэгговском волноводе при z = 100 см, с радиусом изгиба р — 7 см, А = 1 мкм.
Рис. 10. Поперечное распределение амплитуды u(x,z) в брэгговском волноводе при z = 150 см, с радиусом изгиба р — 7 см, А = 1 мкм.
волновода имеет весьма сложный характер и вряд ли может быть количественно описан без детальных численных расчетов.
Итак, в работе предложен и численно реализован способ моделирования изгиба в слабоконтрастных брэгговских волноводах, основанный на использовании метода параболического уравнения с граничными условиями прозрачности в криволинейных координатах. Полученные асимптотические формулы позволяют существенно упростить известные аналитические выражения для граничного условия прозрачности, сохраняя достаточную точность для расчета радиационных потерь реальных волноводов.
Численно рассчитаны распределение амплитуды поля и потери на вытекание в прямом брэгговском волноводе с малым контрастом показателя преломления. Продемонстрировано соответствие решений, получаемых с помощью метода параболического уравнения, решениям спектральной модовой задачи. Показано, что метод параболического уравнения позволяет проследить весь процесс возбуждения, взаимодействия и вытекания мод.
Выполнено численное моделирование амплитуды поля и расчет изгибных потерь в зависимости от радиуса кривизны брэгговского волновода. Показано, что при малых радиусах кривизны происходит срыв основной ТЕ\ моды и перетекание оставшейся в волноводе энергии в захваченные оболочкой моды.
Работа выполнена при поддержке РФФИ, гранты NN 07-02-01244, 07-02-01177.
ЛИТЕРАТУРА
[1] J. Joannopouls, S. Johnson, J. Winn, et al., Photonic Crystals (Princeton University Press, New Jersey, 2008).
[2] P. Yeh, A. Yariv, and E. Marom, J. Opt. Soc. Am. 68, 1196 (1978).
[3] S. Johnson, M. Ibanescu, M. Skorobogatiy, et al., Opt. Express 9, 748 (2001).
[4] T. Engeness, M. Ibanescu, S. Johnson, et al., Opt. Express 11, 1175 (2003).
[5] B. Pal, S. Dasgupta, and M. Shenoy, Opt. Express 13, 621 (2005).
[6] А. Снайдер, Дж. Лав, Теория оптических волноводов (Радио и связь, Москва, 1987).
[7] Т. A. Birks, J. С. Knight, and P. St. J. Russell, Opt. Lett. 22, 961 (1997).
[8] J. C. Baggett, Т. M. Monro, K. Furusawa, and D. J. Richardson, Optics Comm. 227, 317 (2003).
[9] M. A. Leontovich and V. A. Fock, J. Phys. 10, 13 (1946).
[10] Г. Д. Малюжинец, Успехи физических наук 69, 320 (1959).
[11] Ф. Д. Тапперт, Метод параболического уравнения. В кн: Распространение волн и подводная акустика (Мир, Москва, 1980).
[12] М. Levy, Parabolic equation methods for electromagnetic wave propagation (IEE, London, 2000).
[13] Дж. Коул, Методы возмущений в прикладной математике (Мир, Москва, 1972).
[14] V. A. Baranov and А. V. Popov, Wave Motion 17, 337 (1993).
[15] V. A. Baskakov and A. V. Popov, Wave Motion 14, 123 (1991).
[16] В. А. Фок, Проблемы дифракции и распространения электромагнитных волн (Советское радио, Москва, 1970).
[17] А. V. Popov, Radio Science 31, 1781 (1996).
[18] В. М. Бабич, В. С. Булдырев, Асимптотические методы в задачах дифракции коротких волн (Наука, Москва, 1972).
[19] Г. И. Марчук, Методы вычислительной математики (Наука, Москва, 1977).
[20] Д. В. Прокопович, А. В. Попов, А. В. Виноградов, Квант, электроника 37, 873 (2007).
[21] Т. А. Birks, F. Luán, G. J. Pearce, et al., Opt. Express 14, 5688 (2006).
Поступила в редакцию 18 ноября 2008 г.