Диссертация (1097826), страница 11
Текст из файла (страница 11)
§6.3 данной главы). Вид системы алгебраических уравнений позволяетразвить рекуррентный способ построения матрицы рассеяния. Из рассчитаннойчисленно матрицы рассеяния сразу определяются компоненты отраженных ипрошедших волн.На основе методов RCWA и матрицы рассеяния может быть найдено распределение компонент электромагнитного поля внутри структуры. При этомдля нахождения поля в каждом слое, исходя из вида блоховских волн, найденных методом RCWA, вычисляются две матрицы рассеяния для двух подструк-63 тур, окружающих этот слой, и с их помощью определяются все компонентыэлектромагнитного поля.Метод расчета оптических свойств одномерно периодических многослойных структур, содержащих гиротропные среды, может быть расширен на случайструктур с двумерной периодичностью, что позволяет проводить расчеты электромагнитных свойств трехмерных структур.
Сложность обобщения методаRCWA на двумерно периодические структуры связана в первую очередь с тем,что в этом случае нельзя решать задачу рассеяния отдельно для двух (s- и p-)поляризаций. Это приводит к значительному усложнению итоговой системыматричных уравнений.Также из-за двумерной периодичности слоев необходимо производить неодномерное преобразование Фурье (как в случае одномерно периодическихструктур), а двумерное.
По этой причине значительно увеличиваются требуемые вычислительные мощности: вместо N (N – число гармоник в разложенииФурье) матричных элементов необходимо проводить операции с числом элементов порядка N2.Для получения точных значений коэффициентов отражения, пропускания,поглощения и других параметров исследуемой структуры при численном расчете необходимо использовать как можно большее, но конечное, значение N (точное решение вычисляется при бесконечном N, что невозможно в реальности).Так как при расчете используется лишь конечное число слагаемых Фурье-ряда,то результат моделирования отличается от истинного решения уравненийМаксвелла. Разность истинного и полученного численно решений зависит отколичества слагаемых ряда, использованного при расчете.Чтобы минимизировать эту разность при заданном числе N, применяютспециальные правила факторизации, улучшающие сходимость метода RCWA.Они заключаются в преобразовании материальных уравнений электромагнитного поля специальным образом для кусочно непрерывных функций и последую64 щим переводом их в пространство Фурье.
Это сокращает требуемое число Nпри заданной точности расчетов в 2-4 раза, в зависимости от параметров исследуемой структуры.6.3.Метод матрицы рассеянияВо избежание численного переполнения при решении задачи дифракцииметодом RCWA может быть использован подход с помощью матрицы рассеяния [155-159], позволяющий решить систему алгебраических уравнений длянахождения компонент электромагнитного поля отраженных и прошедших волнво всех дифракционных порядках гиротропных металло-диэлектрическихструктур.Матрица рассеяния (S-матрица) является адекватным физическим инструментом для описания электромагнитных свойств материалов, существенно неоднородных на субволновом масштабе, для которых стандартные оптическиеподходы, предполагающие однородность волнового фронта в латеральномнаправлении, не работают. К числу достоинств метода можно отнести возможность нахождения собственных колебаний структуры.
В данной работе методматрицы рассеяния применен для расчёта законов дисперсии плазмонных иволноводных мод различных плазмонных кристаллов.Свет при падении на периодическую структуру дифрагирует на ней, и возникающее при этом электромагнитное поле представимо в виде суперпозициипадающей, отраженных и прошедших плоских волн с продольными волновымичислами k x( m ) , отличающимися от волнового числа k x( i ) падающей волны на вектор обратной решетки mG.
В общем случае комплексные амплитуды отраженных и прошедших волн могут быть объединены в столбец A scat , а амплитуды65 падающих с обеих сторон на структуру волн во всех дифракционных порядках –в столбец A in .Комплексные амплитуды падающих и рассеянных волн связаны через матрицу рассеяния S:SAin A scat .(1.44)Размерность S-матрицы определяется размерностью столбцов амплитуд. Матрицу рассеяния рассчитывают численно.Собственные волны могут распространяться по структуре без внешнеговоздействия, поэтому им соответствует нетривиальное решение однородной задачи с нулевым столбцом амплитуд падающих волн Ain :S 1A scat 0 .(1.45)Для существования нетривиального решения необходимо, чтобыdet(S 1 ) 0(1.46)Из последней формулы определяются дисперсионные зависимости для собственных волн системы.Матрицу рассеяния для многослойных структур вычисляют рекуррентнымобразом.
В структуре выделяют N слоев, однородных вдоль оси z (перпендикулярной границам слоев), и N+1 соответствующих им интерфейсов. Комплексные амплитуды собственных волн на разных сторонах каждого слоя или интерфейса связаны через матрицу переноса Ti, а полная матрица переноса для всейструктуры имеет видT2 N 1Ti 1i(1.47)Вид матриц Ti определяется в результате нахождения собственных волндля каждого слоя и учета граничных условий на интерфейсе. Во внешней средевблизи поверхности крайнего слоя выделяют также виртуальный интерфейс,для которого матрица переноса T0 и матрица рассеяния S0 являются единичны66 ми. Обозначим Si матрицу рассеяния для той части структуры, которая заключена между виртуальным интерфейсом и i-м слоем или интерфейсом.
Для такихматриц имеют место рекуррентные соотношения, связывающие компонентыматриц Si с компонентами матриц Si-1 и Ti. Искомая матрица рассеяния для всейструктуры:S S 2 N 1 .(1.48)6.4. Метод конечных разностей во временной области (FDTD)В данной работе, в главе 7 для моделирования распространения плазмонных импульсов использован метод конечных разностей во временной области(FDTD - Finite-Difference Time-Domain). Этот метод в настоящее время являетсяодним из наиболее распространенных методов численного решения электродинамических задач. Метод FDTD предложен более 40 лет тому назад в работе[160], но получил большое распространение в последнее десятилетие, что связано со значительным прогрессом вычислительный технологий.
Метод FDTDуниверсален и может быть применен практически во всех задачах электродинамики, требующих численного решения. Это и задачи, связанные с расчетомближнего поля различных волноведущих и резонансных структур сложнойформы с неоднородностями, и задачи, требующие моделирования в дальнем поле излучающих структур, антенн, активных приборов СВЧ.Метод FDTD основан на дискретизации уравнений Максвелла, записанныхв дифференциальной форме. Сетки для электрического и магнитного полейсмещены по отношению друг к другу во времени и пространстве на половинушага дискретизации по каждой из переменных. Конечно-разностные уравненияпозволяют определить электрическое и магнитное поля в данный момент времени на основании известных значений полей в предыдущий момент времени.При этом необходимо задать начальные и граничные условия.67 В классическом виде FDTD метода используют эквидистантную ортогональную сетку – решетку Йее [161].
Однако для повышения устойчивости и эффективности метода используют также неэквидистантные и/или неортогональные сетки.Существенным ограничением применимости метода является ограничениешага по времени и требование больших вычислительных ресурсов.Рассмотрим явную разностную схему метода для случая распространенияТЕ-поляризованной волны, имеющей ненулевые компоненты поля H x , H z , E y .Будем считать, что задача двумерная, т.е.
структура однородна вдоль оси OY. Вслучае отсутствия токов и свободных зарядов из уравнений Максвелла следуетследующая система уравнений:H x E ytzE yH ztxE y 1 H x H z ztx (1.49)В области решения x 0, Lx , z 0, Lz необходимо ввести разностную сетку, содержащую N x N z узлов: r xi irx , z j jrz , где i 0,1,..., N x и j 0,1,..., N z .Временной промежуток, на котором решается задача, также развивается на Ntинтервалов, вводя временную сетку с узлами: tk k .
Значение поля в узле пространственно-временной сетки обозначим U ik, j U (irx , jrz , k ) , где U - компонентанапряженности электрического или магнитного поля. Компонента поля E y рассматривается в узлах сетки i, j , k , а компоненты магнитного поля H x , H z - вточках между узлами (i, j 12 , k 12 ) . Тогда пространственные и временные производные можно вычислить как68 U ik, jxU ik 1 , j U ik 1 , j2k1U ik, jt O rx2 2rxk1U i, j 2 U i, j 2rx.(1.50) O 2 Таким образом, уравнение (1.49) можно переписать в виде: H x i , j H x i , j k 12k 121212 H z i , j H z i , jk 12k 121212E Ey k 1y i, jkky i , j 1 Ey ki, jrzE Ey ky i 1, jki, jrx(1.51) H x i , j H x i , j k 12i, jE k 1212 rz H z i , j H z i , jk 1212k 121212 rxИз уравнения (1.51) можно получить выражения для компонент электромагнитного поля в следующий момент времени через значения в данный момент времени: H x i , j 2 1 ss E y i , j 1 E y i , j H z i , j H z i , j ss E y i 1, j E y i , j H x i , j s 12E s 1122s 12s 121212s 1 Ey y i, jsi, jrz rx (1.52) ssssH x i , j H x i , j H z i , j H z i , j . rx rz 1212121212121212Аналогичным образом, для случая ТМ поляризованной волны, имеющейкомпоненты поля Ex , Ez , H y , из уравнений Максвелла имеемEx1 H yt zEz 1 H yt xH y E E x ztx z.(1.53)Переходя к дискретным разностям, можно получить69 E x i , j E x i , j k 12k 121212H ky i , j 1Hy ki, j rz Ez i , j Ez i , j H y i 1, j H y i , jk 12kk 121212H k 1y i, jHy kk rx E x i , j E x i , j k 12i, j(1.54)k 1212 E z i , j E z i , jk 1212rzk 121212rx.Таким образом, разностная схема в этом случае имеет вид: Ex i , j 2 1 ss Hy H y ,1,j iji rz E z i , j E z i , j ss Hy H y i 1, ji, j rx E x i , j s 12s 1122s 12s 121212H s 1y i, j Hy si, j(1.55) ssssE x i , j E x i , j E z i , j E z i , j . rx rz 1212121212121212В методе FDTD важную роль играет правильный выбор граничных условий, позволяющий ограничить расчетную область и избежать отражения волнот ее границ.