CC BY
39
18 января 2012 r.
I'll
ТЕХНОЛОГИИ ИНФОРМАЦИОННОГО ОБЩЕСТВА
Применение модифицированного метода нулевого поля к задаче дифракции на многорядной решетке
С помощью модифицированного метода нулевого поля решена двумерная задача рассеяния на многорядоой решетке, состоящей из имледансных или диэлектрических элементов. Выведена система интегральных уравнений и получены зависимости коэффициентов отражения и прохождения для различных геометрий элементов решетки. Приведены результаты расчетов для задачи дифракции на решетке из круговых цилиндров в низкочастотном приближении.
Кюркчан А.Г.,
Зав, кафедрой ТВ и ПАИ МТУСИ, профессор
Маненков С.А.,
Доцент кафедры маг. анализа
Смирнов В.И.,
Аспирант кафедрой ТВ и ПМ МТУСИ
Введем ие
В работе рассмотрена дифракция плоской волны на двумерной многорялной решетке, состоящей III импс-лансных или диэлектрических цилиндрических тел. Данная шлама представляет большой научный и практический интерес. В частности, с помощью многорялиых решеток моделируют так называемые фотонные кристаллы 11). Фотонные кристаллы обычно формируются из решеток. составленных либо из металлических, либо диэлектрических иилннлров. причем показатель преломления материала элементов решетки больше показателя преломления окружающей решетку среды.
Для решения задачи в работе использовался модифицированный метод нулевого поля (ММН11). предложенный и успешно апробированный в работах (2-5]. Метод нулевого ноля (МНИ), часто называемый в литературе также методом Т-матриц был впервые предложен Уотер-меном [б). В основе этого метода лежит некоторое соотношение (см. ниже), выполняющееся всюду внутри рассеивателя. Требование выполнения этого соотношения на некоторой замкнутой поверхности внутри рассеивателя позволяет свести краевую задачу дифракции к решению интегрального уравнения Фредгольма 1-го рода с гладким ядром. В работах (2-5) было показано, что интегральное уравнение МНИ имеет решение, соответствующее краевой задаче, лишь в том и только в том случае, если вспомогательная поверхность, на которой ставится условие нулевого поля, охватывает множество особенностей аналитического продолжения дифракционного поля внутрь рассеивателя. Кроме того, в этих работах покатано, что для получения наиболее быстродействующих и устойчивых алгоритмов вспомогательную поверхность целесообразно строить при помощи аналитической деформации границы рассеивателя |7|.
Отметим также, что в рассматриваемой ниже задаче рассеяния волн на периодической решетке используется периодическая функция Грина (ФГ). вычисление которой представляет определенные трудности. В работе эта функция рассчитывалась двумя способами 18-10]. В случае большого расстояния (вдоль координаты, перпендикулярной оси решетки) между точкой источника и точкой наблюдения ФГ рагзатазась в ряд. получаемый с помощью формулы Пуассона. В случае, когда указанное расстояние маю. возможно использование разложения ФГ в ряд по цилиндрическим гармоникам. Преимущество данного метода расчета ФГ решетки состоит в том. что он позволяет шбежать трудности, связанной с медленной
сходимостью ряда для ФГ, который получается е помощью формулы Пуассона.
I. Мос книжка задачи
Рассмотрим периодическую решетку, состоящую из нескольких рядов цилиндрических тел, образующие которых параллельны оси г (см. рис. I). Число рядов решетки обозначим через А/. Предполагаем, что расстояние между центрами тел у'-го ряда решетки равно (// (периоды решетки). где ] = |...I/ . При этом тела в разных рядах мо-
гут быть различны. Пусть Я0/ - контур центрального элемента у-го ряда решетки. Введем также локальные координаты (х,, ^у). связанные с центральным элементом у-го ряда. 11рн этом пусть х = х,, у = у{ ■ Обозначим через г0 = • радиус-вектор из начата координат в цен-
тра! ыюм элементе у - го ряда в начато координат цен-тратьного элемента /-го ряда (см. рис. I). Считаем, что структура облучается плоской волной:
и° =схр(-Мгсо5(<о-<^)). (1)
где (г,<р) - цилиндрические координаты, к - волновое число, - угол скольжения волны. Дифракционное иоле и'(х.у) вне решетки удовлетворяет однородному уравнению Гельмгольца:
Дм' + *:«'=0. (2)
На бесконечности вторичное поле удовлетворяет условию нхзучения:
(3)
ч\х,у) » £ ехр(-тл*,ть>»л).
ы »—«
t ... _ // ♦ ^ B U: _ и.: . причем знак кватратного
d, ' h h
корпя выоирастся из условия неположительности его мнимой части. = kd, cos<р0 ■ Знак плюс в формуле (3) относится к верхнему полупространству (нал решеткой). Рассмотрим два вида граничных условий на поверхности ненгразьных элементов каждого ряда решетки. В случае импслансного краевою условия имеем:
j= I..М■ И)
где с fen - производная по направлению внешней нормали к S0/. Zf ' импеданс. В случае диэлектрических элементов решетки будем иметь:
//) / ft, Е- пол яри »ация.
С и _ I <V *
СП Г). СП
е\ / гг. // - поляризация.
(5)
где и - поле внутри нентратыюго элемента у-го ряда
решетки, с, /л - соответственно диэлектрические
и магнитные пропииасмостн среды вис и внутри элементов решегки.
60
T-Comm, #11-2011
ТЕХНОЛОГИИ ИНФОРМАЦИОННОГО ОБЩЕСТВА
2. Вывод основных соотношений.
Будем решать поставленную задачу с помощью ММНИ. Представим вторичное поле в обдаст вне решетки в виде:
(6)
с
ж;
О)
и v а (J!l Выражение (7) представляет собой поле, рас' о»;
сеяное /—м рядом решетки, tan пса ни ос в системе координат. связанной с центральным элементом этого ряда. Функция - ФГ / -го ряда. Она имеет вид:
C/W.'D=j X п'>' * О/-у',у )«p(-/.v/,)
(8)
В случае диэлектрических элементов решетки поле в области внутри нейтрального элемента j -го ряда имеет вид:
V <9>
г сС'(г.г)
,<?,) = - J Tv------
с/('у.'/>■ ^ (*; 7(х/ - +^ ^• (|0)
причем у/ _ и - волновое число внутри элементов ' *>;
У -го ряда решетки.
Заладим далее конту ры 50 в локальной полярной системе координат:
х, = )<•<««>,. У, • Р, («>, )5>п«>,- (II)
Вспомогательные контуры (внутри и вне контура Я01). на которых ставится условие нулевого поля, будем строить, как уже отмечаюсь, при помощи анаштической
уравнения:
х'=р]у = Р] яла). (12)
где
а; =arg^/’(r). р) = |*/(г)|. (г) = р/г ±) схр( /г +й';).
гМ0.2я]. (|3)
Верхний так в формулах (13) соответствует контуру
£оу . Величины 3* и 8~ - положительные параметры,
отвечающие за степень деформаиии контура тела. Выбор этих параметров подробно описан в работах (2-5). (7). В соответствии с ММНП потребуем, чтобы выполнялись следующие равенства:
X |[^ '* -с/(?,-'п|>(• - -»•(*>)*
г,еЦ,. /• 1.2,...,.V. (14)
в случае импсдансных краевых условий на контурах элементов решетки, и
II
Л/Дг.г)
“/('П—т4------------v,(/=v;,(rif)
дп,
л; = -«“(г, >•
(15)
I
сЧ/'(г,г)
г, €ІІ,. у = 1.2..А/. (16)
в случае диэлектрических элементов решетки. Таким образом. залача сводится к системе интегральных уравнений первого рода относительно неизвестных токов на контурах центратьных элементов каждого ряда решетки. В формулах (14) и (15) функция С,(ггг{) выражается следующим образом:
см,-О=^ ^«іа(*7(*, -*іі +С>'/ і
(17)
где (р,гдІ.) - координаты «центра» у-го ряда в системе координат, связшшой с /-м рядом. Очевидно. что Рц = ~Р/г Ял в ~Яі, •
Система интегральных уравнений (14) или (15), (16) решается численно с помощью метода коллокации. описанного в работах (7), (10). Рассмотрим вопрос о вычислении периодической ФГ решетки. С использованием формулы Пуассона, получим
б»г>=—У—«р( * Рц»-'п, і У/ -у, *чц і)
(18)
/1а и пая формула пригодна дія вычислении ФГ при больших значениях величины к\у{-у’+ |. Исли ука-
занная величина мата то можно применить следующее преобразование ФГ. Имеем (8):
С,(г,.гГ) = ^| £к:'(ф < + /> >ї+9і,Г )«РН»й)+
1—1' ісчі,",ф 1 )•
(19)
(20) (21) (22)
Ч, = X H-:‘(khi' )«ЧЧ**/Г, -«/,)•
хгр.1
«»#, - */< - Рц- ч, si" А’, ■ <?/,•
Я„ С05ф() -X,-X,. R,, 5ІПФ,, = yj - у; ■
С помощью формулы
, ЖІ2+4
Н[:\кг)с\ЫЬнр) = -— f ехр(Tikrcos(у + <р)~ iny)dy'
n J
справедливой соответственно при условии cos^>0 (в
этом случае берется верхний знак в (23)) или cos^xO,
суммируем ряды (20):
ехр( -/(Q +1 )гг,* + / Д,. - iny)
--------<iy +
(24)
. f г «
'=7 J “
7І
I - exp(-/<r/) f" exp( -/(Q ♦ 1 )(7 -/А,, - m/)
-rm-gtl * ^
I -ехр(-/<т/ )
где Д;/ s kp,, cosY + kqh sin y. a,* = kd, cos y±Zr ^MTcrPa‘ лы в формуле (24) находились с помощью перехода к интегрированию по контуру, составляющему в окрестности точки / = 0 угол 45 с действительной осью у [II]. Таким образом. ФГ решетки записывается в виде (19), причем члены бесконечного ряда в этой формуле выражаются через интегралы Зоммсрфсльда . Заметим, что рал в формуле (19) сходится абсолютно при условии
T-Comm, #11-2011
61
ТЕХНОЛОГИИ ИНФОРМАЦИОННОГО ОБЩЕСТВА
-< I*
.. (25)
(СММ-1а,1
На практике номер {) выбирается так, чтобы отношение в левой части формулы (25) не превосходило 0,5 для тех значении аргументов ФГ. для которых используется представление (19).
После нахождения неизвестных токов можно определить коэффициенты отражения и прохождения плоской волны. В этом случае вторичное поле в верхнем и нижнем полупространстве выражается по формуле (3). в которой:
Л» - . *о/ (и/. • W) ’
(26)
(29)
£»/(»’/,.%) = goi(±aTCCOS(wh ,к))- <27)
а функция Rw(<p,,ip0) представляет собой диафамму рассеяния центрального элемента /-го ряда решетки и выражается по формуле:
^ J| И, )^Т - V, |jcxp(Лг/соЛа - *>,'))Л,’
/ = 1.2 .W ■ (28)
В формуле (28) явно указана зависимость диафаммы от угла падения волны. Коэффициенты отражения и прохождения плоской волны для многорядной решетки выражаются по формулам: и
*0 exp(/*r„ СОЯ(^>,, +Й,)).
«
н
7» = 1 + X 4*о *xp(«r„ cos(«n, - %)) где ru cos <ри = ри, ru sin (рц =
3. Численные результаты
Рассмофим результаты расчетов для различной геометрии задачи. Все дальнейшие результаты получены для случая Е-поляризации падающего поля. Ьыли рассчитаны зависимости коэффинисчгтов отражения и прохождения плоской волны от параметра d / А, (А - длина волны. d - период решетки) при падении ее на решетку, составленную из одинаковых круглых диэлектрических цилиндров. На рис.2 изображена зависимость коэффициента отражения по мощности двухрядной решетки для следующих параметров задачи: относительный радиус цилиндров ald = 0.3, расстояние между рядами hid - qK Id - 0.7. диэлектрическая проницаемость материала внуфи цилиндров е. =2. фп =90 (нормальное падение плоской волны). В силу того, что в рассматриваемом случае имеется одномодовый режим, коэффициент прохождения по мощности равен Tr - I - Rr и поэтому не
приводится. Пунктирной кривой на рисунке изображена зависимость коэффициента для решетки, составленной
из цилиндров с сечением в Blue суперэллипсов:
V* f ,.\г«
s I
(f) -(f) -
Параметры суперэллипсов имели значения aid = 0.3. b = а и q = 10. В силу того, что q достаточно велико, геометрия решетки из суперэллиптичсских цилиндров близка к гсомсфии решетки из цилиндров квадратного сечения. На рис. 3 показана зависимость коэффициента R^ от параметра d/л для шести рядной решетки
из одинаковых круговых цилиндров и цилиндров с супер*
эллиптическим сечением. Все параметры задачи имели те же значения, что и в предыдущем случае. Отметим появление нескольких значений волнового параметра dIА. соответствующих критическим частотам мод Флокс, рас-просфаияющнхся в рассмафивасмой структуре. Полученные зависимости для круговых элементов решетки полностью совпадают с аналогичными кривыми, приведенными в работе |12). Отмстим, что для приведенных параметров задачи была проведена проверка выполнения закона сохранения энергии. Расчеты показа!и. что закон сохранения выполняется с точностью 10 * для цилиндров кругового сечения и 10 ' для цилиндров с сечением в виде су нерэллипса.
На рис. 4 и 5 изображены зависимости коэффициента офаження и прохождения Тр плоской волны XIя гексагональной шсстирялной решетки от величины «зазора» kl, между элементами решетки. Геомсфия решетки такова, что горизонтальная координата нечетных слоев равна половине периода решетки. Элементы решетки идеально офажающис. Сплошные кривые на рисунке соответствуют зависимостям коэффициента офаження и прохождения для решетки, составленной из круговых цилиндров, а пунктирными кривыми показаны зависимости для решетки из цилиндров супсрэллиитичсского сечения, с помощью которых моделировались элементы квадратного сечения. Парамсфы геометрии задачи имели значения: 2а = 0.15л. Л«3>/ЗД/8, где 2а - диаметр круговых цилиндров. либо ширина цилиндров супсрэллиитичсского сечения. Рассмафнвалось нормальное падение плоской волны. На рис. 6 и 7 показаны аналогичные зависимости коэффициентов офаження и прохождения от параметра kL для гексагональной шсстирялной решетки из диэлск-фических цилиндров. Все размеры элементов и обозначения кривых те же, что и в предыдущем случае. Диэлск-фичсская проницаемость материала элементов решетки составляла g = 2.25 • Как видно из приведенных рисунков в случае идеальных элементов решетки при малых значениях параметра kL происходит практически полное отражение от решетки. С увеличением величины расстояния между элементами решетки коэффициент отражения по энерши резко уменьшается. В случае диэлектрических элементов решетки картина совершенно другая. Зависимости коэффициента отражения R имеют характерный
резкий максимум при некотором достаточно большом значении параметра kl.. В остальной части диапазона изменения кL коэффициент отражения не превосходит 10 процентов.
Представляет интерес сравнение строгого и асимптотическое решений задачи в предположении, что размеры элементов решетки много меньше длины волны. 11а рис. 8 показаны частотные зависимости коэффициентов отражения и прохождения для десятирядной решетки из одинаковых идеально офажаюших круговых цилиндров. Параметры гсомсфии имели следующие значения aid- 0.01» hid- 0.03. угол падения волны 0ц=45‘. Пунктирная кривая на рисунке иллюстрирует ас p5Q V от н чес кос решение задачи, получаемое с помощью метода диаграммных уравнений, который описан в работах (13, 14], а сплошная кривая соответствует строгому решению, полученному при помощи MMHI1. На рис. 9 изображены зависимости коэффициентов офаження и прохождения плоской волны от угла паления плоской волны для тех же параметров модели. Обозначения кривых те же. что и на предыдущем рисунке. Кривые построены для значения волнового па-
62
T-Comm, #11-2011
