Асимптотические решения уравнения Гельмгольца для псевдопериодических структур Текст научной статьи по специальности «Математика»

АСИМПТОТИЧЕСКИЕ РЕШЕНИЯ УРАВНЕНИЯ ГЕЛЬМГОЛЬЦА ДЛЯ ПСЕВДОПЕРИОДИЧЕСКИХ СТРУКТУР

Досколович Л.Л., Харитонов С.И., Казанский Н.Л., Тулупова Е.А., Скуратов С.А. Институт систем обработки изображений РАН, Самарский государственный аэрокосмический университет

Аннотация

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

Введение

В общем случае решение задач дифракции представляет собой сложную математическую задачу и требует специальных численных методов и значительных вычислительных затрат [1]-[6]. Асимптотические методы относятся к категории аналитических методов и играют большую роль [1]-[5]. В ряде случаев они позволяют быстро и эффективно оценить результат, не прибегая к громоздким численным расчётам. Для развития асимптотических подходов для решения уравнений Максвелла целесообразно рассмотреть более простую задачу решения уравнения Гельмгольца. Асимптотические методы решения уравнения Гельмгольца представлены в работе [4]. Приближение геометрической оптики, в рамках которого решается уравнение Гельмгольца справедливо только с случае, когда диэлектрическая проницаемость меняется медленно. Во многих случаях, например, дифракции света на дифракционной решетке, период которой составляет несколько длин волн, это условие не выполняется. Для решения задач дифракции в случае, когда объект представляет периодическую структуру, обычно используется дифференциальный метод [8]. Однако этот метод нельзя применить для решения задачи дифракции на непериодических структурах, например, бинарной линзе Френеля. Асимптотические методы решения уравнения Гельмгольца в случае, когда объект представляет собой структуру со слабоизменяющимся пространственным периодом, были рассмотрены в работе [7]. Однако метод решения, изложенный в статье, не давал результатов по сравнению с использованием приближения геометрической оптики при расчёте поля в плоскости непосредственно прилегающей к оптическому элементу. В данной работе метод был модернизирован, и на его основе получен асимптотический метод решения трехмерной задачи расчёта поля на выходе дифракционного оптического элемента.

1. Метод локальной аппроксимации

В данной работе предлагается метод, основанный на представлении объекта в окрестности точки х0 псевдопериодической структуры дифракционной решеткой. Для расчета поля в окрестности произвольной точки х0 используется дифференциальный метод [8].

Рассмотрим дифракцию скалярного волнового поля на псевдопериодическом одномерном слое. В оптике этот случай описывает дифракцию волны, у которой вектор электрического поля лежит в плоскости

дифракционного оптического элемента (ДОЭ). ДОЭ предполагаем расположенным в области 0<х<а, где 0 и а минимальная и максимальная высота дифракционного микрорельефа. Рассмотрим вычисление волнового поля в окрестности произвольной точки х0, расположенной на ДОЭ.

Волновое поле описывается уравнением Гельм-гольца

д2Е д2Е Е 0

2- + —г + к2 е Е = 0,

дг2 дг2 '

, 2п

где к =--волновое число, Я - длина волны.

В окрестности х0 распределение диэлектриче ской проницаемости описывается выражением

2тп(х - х0)

• (х, г) = £ е п (х0)ехр

(2)

где й - размер области, в которой производится аппроксимация.

Решение уравнения Гельмгольца в области ДОЭ, в области перед ДОЭ и в области за ДОЭ представим в виде

Е(х, г) = ехр(1ках0) ^Е" (х0, г) ехр[кап (х - х0)] (3)

п

Е (х, г) =

= ехр (ках0) ^ I" (х0) ехр [¡кап (х - х0) + 1крпг] + (4)

п

+ ехр (ках0) ^ Яп (х0) ехр [кап (х - х0) - ¡квпг]

Е (х, г) = ехр (¡ках0) х

х£ Тп (х0) ехр [¡кап (х - х0) +1кРпг\.

(5)

Коэффициенты 1п определяются падающим пучком. Коэффициенты Яп, Тп и функции Еп(х, г) подлежат определению: ат = а + 2л1т!Ы,

вт = V1 -ат2 . Подставляя выражение (3) в (1) получаем систему обыкновенных дифференциальных уравнений для определения функций Еп(х, г):

д2 Ет (г)

дг2

-к2 (т-,-а)8Г)) (г) = 0.

(6)

Решение системы дифференциальных уравнений и нахождение коэффициентов Кп, Тп рассмотрено в приложении 1. Поле в точке х0 непосредственно за оптическим элементом имеет вид:

п

Т1

Т1

Е(х, 2) = ехр(¡ках0) X Тп(х0) ехр[квпа]. (7)

п

Волновое поле на некотором расстоянии за оптическим элементом можно вычислить с помощью интеграла Кирхгофа ([1]-[2]). Предлагаемый подход позволяет вычислить поле от оптического элемента с большой апертурой.

2. Дифракция на оптическом элементе, обладающем зонной структурой

Рассмотрим дифракцию на ДОЭ, обладающем зонной структурой. В случае дифракционного оптического элемента распределение диэлектрической проницаемости описывается выражением

е(х, 2) = Xё"(2)exp[ikng(х)], (8)

п

где кё(х) - фазовая функция дифракционного оптического элемента. Границы зон определяются по формуле ё (хт) = ш!.

В этом случае

%

X1 I ёш (х0 +%)>

ш а -у

> exP (((х0 + «.)) (. П)

Падающее поле имеет вид: Е( х) = К (х) exp (кр( х)) а/

1

1п = — | К(х0 + (к(р(х0 +£) -8£))х

а -а/2

х exp (-|

(9)

(10)

(11)

а )

др

где 8 =-, К(х) - амплитуда падающей волны,

дх0

<р(х) - эйконал падающей волны.

В случае, когда функция ё(х) плавно изменяется, выражение для коэффициентов Фурье имеет вид:

02

1

еп = Х 0 1 ёш (^^ ( (х0)ш )х

ш а -<У2

х exp ( Iк ¡а£ + -в- ) ш ) exp (-юп£)а^,

(12)

где юп = ■

2пп

ка

1п =

К (х0)

02

| exp (к(р(х0) )>

(13)

с exp ^ 1к У | exp (-гап£) а%.

После элементарных преобразований выражение для коэффициентов Фурье имеет вид:

еп = 0 Хёш (^^ ( (х0))

<exp

ш

( ( -Iк

V V

(аш -юп ) 2вш

2 Л |

//

( / \2\ 1крш юп-аш |

> | exp

-02

вш )

К (Xo)exp (¡кр( х0))

1п =-^-exp

а

где а =

дх0

С

-к

V 2У

2|

у '

д ё д !Р Г^ "

( /

1к У

V

д2 ё др

^ у = дх02 дх0 '

(15)

При а = м -2- и в = 0 приближенно можно по-

дх0

ложить

еш = ёш(Xo)exp(%(х0)ш), (16)

1п = exp(-1ках0)К(x0)exp(кр(х0))8%, (17)

где 88 - символ Кронекера.

В этом случае поле на выходе дифракционного оптического элемента имеет вид:

Е(х0) = ХТп\1кпё(х0)] . (18)

п

Это выражение совпадает с выражением, полученным в работе [7]. Оно приближенно описывает поле в случае, когда размер зоны на ДОЭ меняется слабо. В общем случае метод давал незначительные поправки по сравнению с использованием приближения геометрической оптики [4]. Тогда при расчете надо использовать для коэффициентов Фурье выражения (14), (15)

На рис. 1 - 3 представлены графики интенсивности поля в фокальной плоскости бинарной линзы Френеля, имеющей следующие параметры: апертура - 130 мкм, фокусное расстояние - 134 мкм, высота штрихов - 1 мкм, диэлектрическая проницаемость штрихов - 2.25, длина волны - 1мкм, угол падения -30°, диэлектрическая проницаемость перед линзой и за линзой - 1. Для расчета поля на рис. 1 был использован дифференциальный метод [8]. Для расчета поля на рис. 2 и 3 был использован метод локальной аппроксимации. Поле за ДОЭ рассчитывалось с использованием интеграла Кирхгофа. Отличие интенсивности в центре во всех рассмотренных случаях составляет несколько процентов, но при использовании метода локальной аппроксимации потребовалось значительно меньше ресурсов оперативной памяти и процессорного времени. Отсутствие высокочастотных биений на последних двух графиках объясняется фильтрующим свойством интеграла Кирхгофа.

0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0

1

■МЛ-

0 16.25 52.5 48.75 65 81.25 97.5 113.75 130 Рис. 1. Интенсивность в фокальной плоскости бинарной линзы, рассчитанная с использованием дифференциального метода, при числе учитываемых коэффициентов Фурье -160

1

0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0

—

О 16.25 32.5 48.75 65 81.25 97.5 113.75 130 Рис. 2. Интенсивность в фокальной плоскости бинарной линзы, рассчитанная с использованием метода локальной аппроксимации при числе локальных периодов - 10 и коэффициентов Фурье 60 в каждом периоде

1

0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 О

— —

О 16.25 32.5 48.75 65 81.25 97.5 113.75 130

Рис .3. Интенсивность в фокальной плоскости бинарной линзы, рассчитанная с использованием метода локальной аппроксимации при числе локальных периодов - 10 и коэффициентов Фурье 60 в каждом периоде.

3. Решение уравнения Гельмгольца в трёхмерном случае Рассмотрим дифракцию скалярного на объекте, представляющем собой псевдопериодический объект. Предполагаем расположенным в области 0 < х < а , где 0 и а - минимальная и максимальная высота слоя. Волновое поле описывается уравнением Гельмгольца:

^ + + + к2 еЕ = 0. (19)

дх дУ дг2

В окрестности (х0, у0), распределение проницаемости описывается выражением

2п/п1 (х - х0)

е( х, у, г) = Х еп (х0)ехр

п

Т.жт2( у - У0)

йх

х ехр

(20)

(22)

где йх й2 - периоды по различным осям.

Решение уравнения Гельмгольца в области ДОЭ, в области перед ДОЭ и в области за ДОЭ представим в виде: Е (х, у, г) =

= ехр (¡к (ах0 + вУ0)) X Е-2 ( у,, г) х (21)

х ехр[¡к (а (х - х0) + в (у - У0))] Е(х, у, г) = ехр (¡к (ах0 + Руа)) х

хХ 1пщ (ь, у0 )х

х ехр[¡к ((х - х0)+в (у- у0))+ис/п^г]+ + ехр (¡к (ах0 + Руа ))х хХ R"1"2 (х0,у0)х

пищ

х ехр [¡к ( (х - х0 ) ■+ РП1 (у - у0 )) - /kY"l"2 г]

Е(х, у, г) = ехр (¡к (ах0 + ву0)) X Т^2 (0, у0) х

"l,"2

х ехр [¡к (ащ (х - х0 ) + РПг (у - у0 )) + 1кГ„Л г].

Коэффициенты 1пПг определяются падающим пучком. Коэффициенты R"1"2, Т^2 и функции Е (х, у, г) подлежат определению.

ат = а + 2я7'т/Ы1 , вп = в + 2я7п/Ы2 ,

Гтп =^11 -ат2 -в„2 .

Подставляя выражение (21) в (19), получаем систему обыкновенных дифференциальных уравнений для определения функций Етт2 (х, у, г)

(23)

д2 Ет1т2 (г)

дг2

"к2 (ет1 -Ая

-(а2 +а2 )Е'1'2 (г) = 0.

(24)

Решение системы уравнений аналогично решению системы в пункте 1.

Поле непосредственно за оптическим элементом

Е( х, у, г) = ехр (¡к (ах0 + ву0 ))х хХ Тп1п (0,у0)ехра).

(25)

Существует случай, когда дифракция на ДОЭ сводится к дифракции на одномерной решетке:

e( x, y, z) = Xe" (z)exp

2nin

' x - x0 + У - y0 А

d1

d,

2

С помощью преобразования координат:

í x = x0 + u cos t - v sin t,

.(26)

(27)

cos t sin t Y x - x

- sin t cos t Д y - y0

(36)

Для дальнейших рассуждений разложим функцию ё(х, у) в ряд Тейлора.

[ y = y0 + u sin t + v cos t. Если угол равен t = -arctg (dl¡ d2)

g(x y) = g(x0. y0) +^—(x - x0 ) + т—(y - y0) +

dxx0 dy0

1 1

1

dj d ^ d

Выражение для диэлектрической проницаемости представляется в виде:

2ninu

1

+

2

'Й (x - x- )2 +0 (y - >0)

K9x0

+ (37)

д2 g dx0 дУч

s( x, z) = Xs" (z)exp

d

(28)

Выражение для поля имеет вид:

E(и, v, z) = XE"' (z)exp[ik(anu + ptv)], (29)

E(u, v, z) = XIn' (z)exp[ik(anu + p¡v) + ikynIz] +

(30)

+X R"' (z) exp [ik (anu + fi,v)- ikynlz]

n,l

E(u, v, z) = X Tnl (z)exp [ik (anu + p,v) + ikynlz], (31)

d2 Em1( z) dz2

(x- x0 )(y -y0 ).

Пусть t = аг^ё I -2- -2- |. В этом случае функ-

VдУo| дх0 )

ция ё(х,у) в окрестности точки (х0,у0) имеет вид:

ё(x, у) = ё(xo, У0)+\}((я(^ У0) )2 и + а(х У )и2 1 (38)

+ а( х02у0)и +1 (, у0)у2 + 2С( х0, у0)ш ).

Выражения для коэффициентов приведены в Приложении 3.

Рассмотрим случай когда с=0, что соответствует переходу к случаю радиальной симметрии. В этом случае функция ё(х, у) в окрестности точки (х0, у0) имеет вид:

- / / п п\ \ ,, ИМСС1 ьид.

+ k2 (-1 -( +af ) )E1(z) = 0.(32) f-^

g(xy) = g(x0>y0) +j(g(x0'y0)) u

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

4. Дифракция на двумерных зонных структурах

Пусть волна с комплексной амплитудой Е(х) = К (х,у)exp(кр(х, у)) падает на область с проницаемостью в виде:

е(х, у, 2) = X ёп (z)exp[ik"g (х, у)]. (33)

п

Выберем систему координат таким образом, чтобы ось совпадала по направлению с градиентом функции ё (х, у)

Уё (х, у)

-2 (a( x0' >0)u 2 + b( x0' >0)v 2 ).

(39)

Выражение для Фурье коэффициентов проницаемости принимают вид:

еп = X ёш (х0 , у> ) eXP (кё (х0 , у0)) х

с exp

m

' ' -ik

(hm -юп )

Л А

V V

д/

sn = D J exp

_ n /

2am

'

exp

( A -ik ^ 2bm

SnFl.

i - Л/

ikam ' „ (on - hm £--n-

2A

F =

1 r

D -Lexp

' ikbm (on

2 V bm

Л

(40)

(41)

(42)

V(g (x, у) )2

(34)

где еи, еу - базисные вектора новой системы координат. Преобразование координат имеет вид:

x ) Í x0 А ' cos t - sin t Y u

VУ) Vy0 ) V sin t cos t JVvy Обратное преобразование имеет вид:

(35)

Рассмотрим случай, когда волновое число стремится к бесконечности (k ), но kDj ^ const, kD2 ^ const. В этом случае Fl = 0 для ' ф 0 и отличны от нуля только S"0. В этом случае задача сводится к задаче конической дифракции, рассмотренной в пункте 3.

Рассмотрим теперь падающее поле. Падающее поле в окрестности точки (х0, _у0) имеет вид:

E(u, v) = K (х, y) exp (ikp(x, y)) = K x

<exi

:p (/((x0 + u cos t - v sin t, y0 + u sin t + v cos t)).

(43)

n

n

Разложим функцию <р(х, у) в окрестности точки (хс, ус) в ряд Тейлора:

<( х, у) = <( х^ у0) +

1 / 2 2 V (44)

+2 (((х0, у>)У + аД х0, у0)и + 2С1( х0, у0)му),

где <(х, у) = <( х^ у^ + ем + /у.

Выражения для коэффициентов приведены в Приложении 3.

Разложим функцию К ехр (¡к<(х, у)) в квазипериодический ряд Фурье: К ехр (¡к<(х, у)) =

= X ^ ехр [¡к ( м + вПг V) + ¡кущщ г] , (45)

где

/И1И2 =

х exp х exp

JJk (x» y)х

(i(x0, y0)v 2 + «1 (x0, y0)u2 + 2с1 (x0, У0 )uv)

ik

-ik

2nn,u 2nnV

v kD,

kD2

dudv,

2пп1 2пп2

а1=е+-щ> в = /+-щ-

В качестве примера рассмотрим дифракцию ра-диально-симметричной скалярной волны на ради-ально-симметричном объекте. В этом случае разложения имеют вид:

( ) () +М») + 1 д 2 g (Г0 ) 2 + 1 д 2 g (Г0 ) 2 (46)

g( х у) = gOo) + , м ,2 м --у (46)

дг0 2 дг0 2г0 дг0

ду (r0)

У x, у) = У (г0) + —u +

дг0

+1 дУ0 (Г0 ) u 2 + ±_ дУр (Г0 ) v2

(47)

2 дг:

2r0 дг0

где r0 =yl Х0 + У02 .

Рассмотрим теперь асимптотический случай k ^го , kD, ^ const, kD2 ^ const. В этом случае задача сводится к дифракции нескольких плоских волн на одномерной дифракционной решетке. Выражения для Фурье коэффициентов проницаемости и падающего поля имеют вид:

Sn = 77 X gm (r0 ) eXP (ikg (r0 ) ) х

D1 m

<exp

( (/ \2 \\ (am -wn )

-ik

v v % (

i eXP

2ftm ik в m

(48)

//

4-

Pm

2Л

d4,

In =

K (r>)exp (ikp(r0))

где a =

3g_ дг0

D1 % (

х J exp

P=4 дГ

0

exp

2A

-ik-

* f (^

_ ду _

дг0 '

2f

2 Л

d4,

д 2p дг02

'0 ^'0

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

Заключение

В работе предложен метод решения задачи дифракции в случае, когда объект представляет собой дифракционный оптический элемент, обладающий зонной структурой. Получено, что в случае, когда к ^го задача дифракции произвольной волны сводится к решению задачи конической дифракции на решетке. Это позволяет значительно уменьшить вычислительную сложность при численной реализации методов решения задач дифракции. Предложенный метод позволяет оценить параметры волнового поля без использования суперкомпьютеров.

В дальнейшем предполагается использовать аналогичный подход для решения задач в рамках векторной электромагнитной теории.

Приложение 1

Решение уравнения имеет вид:

Ет = Х(акехр(1кХкг) + Ькехр(-1кХк (г - а)))ет ,

к

где ек и Хк - собственные вектора и собственные

числа матрицы (ет-' - а2 8? ) ек = ек .

Поле внутри дифракционной решетки имеет вид:

Е ( х, г ) = ехр (¡ках0) X (акехр (кХкг ) +

п

Ькехр (-¡кХк (г - а))) е1ехр (¡кап (х - х0)).

Коэффициенты отражения и пропускания находятся из системы линейных уравнений

X( s;Rk (Х0)-e;ak -k

-enkexp(ikXkd)bk + 8"Jk (x0 )) = -In (x0 ),

X(-PnSkRk(x0)-Xke;ak -k

-Xkel exp(ikXkd)bk +b"kTk(x0)) = (

X(S"nRk (x0 ) + elexp(ikXkd)ak + enkbk -k

-S^exp(ikfikd)Tk (x0 )) = 0,

X (SlRk (x0 ) + Xkelexp (ikXkd) ak -

k

-Xke"kbk - fikexp(ik£kd)8;Tk ( )) = 0. Приложение 2

V ' / ? n j-y kl nl k

X (8kR - eka -

k

- e"k,exp (ikXkla )bkl + 8 nT kl ) = - Inl,

X(-Pi8lRl -h,e>kl -k

-X^exp (ikXka )bkl + 8"Jkl) = - PnlIkl,

X(8k"Rkl +e'klexp (ikXkla )akl +

k

+e"kbkl - 8nkexp (ikpua) Tkl) = 0,

X(0• 8;Rkl+ Xuenuexp(ikXkla)ak -k

-Xkle;,bkl - Pep (ikPlfl)8;Tkl) = 0.

Здесь суммирование по индексу l не производится, и система уравнений распадается на несколько систем уравнений меньшей размерности.

Приложение 3

Кх0 , Уо ) =

д2 g , д2 g . 2 д2 g —2-cos21 + —f sin21 +-— -

дхо дУо дхо дУо

sin2t,

it \ д2g . 2 д2 g 2 д2 g . „ b (x Уо ) = —fsin21 + —fcos2 t--— sin2t,

дУо2

дходуо

-(хо , Уо)

д г д г д g

-sin2t+—#- sin 2t + 2-2—cos2t,

дх„2

Чдхо

дg . -cos t+^- sin t

дУо

дУо

Л

дхо дУо

г

/=

дg . —cos t--sin t

ЧдУо дхо

Л

31 (хо , Уо )

д2 р дх1

cos2 2t

д2 р .

sin2 2t -

д2 р

, . д2 р . 2 д2 р 2 b1 (хо,Уо) = —^sin2 t^—-^cos2 t-

дходУо д2 р

sin2t,

дхо

дУо2

ч 52 р д2р . „ „ с1 (хо,Уо) = _ 2 sin2t+ _ 2 sin2t-2

dx¡ дУ2

дГодУо д 2р дхо дУо

sin2t,

cos 2t.

Благодарности

Работа выполнена при поддержке российско-американской программы «Фундаментальные исследования и высшее образование» («BRHE»), а так же грантов РФФИ №№ о3-о1-ооЮ9, о4-о7-9о149 и о4-о1-96517.

Литература

1. Борн М., Вольф Э. Основы оптики // Pergamon Press, 1986.

2. Виноградова М.Б., Руденко О.В., Сухоруков А.П. // Теория волн. М.: Наука, 1976.

3. Ильинский А.С., Кравцов В.В., Свешников А.Г. // Математические модели электродинамики. М.: Высшая школа, 1991.

4. Кравцов В.В., Орлов А.А. // Геометрическая оптика неоднородных сред. М.: Наука, 1979.

5. Федорюк М.В. Асимптотики, интегралы и ряды. // М.: Наука, 1987.

6. Golub M.A. , Doskolovich L.L. ,. Kazansky N.L, Khari-tonov S.I. , Soifer V.A. // Journal of Modern Optics. 1992. 39(6) 1245.

7. Досколович Л.Л., Харитонов С.И., Казанский Н.Л. // Асимптотические решения скалярного волнового уравнения.

8. Moharm M.G. Gaylord T.K. J.Opt. Sos. Am. A 3(11) 1986. 178о.