CC BY
95
ПРИМЕНЕНИЕ МЕТОДА РАССЕЯННОГО ПОЛЯ ДЛЯ FDTD-МОДЕЛИРОВАНИЯ ЭЛЕКТРОМАГНИТНЫХ ПОЛЕЙ ВБЛИЗИ ДИЭЛЕКТРИЧЕСКИХ И МЕТАЛЛИЧЕСКИХ НАНОЧАСТИЦ
Асенчик О.Д., Стародубцев Е.Г.
Гомельский государственный технический университет имени П. О. Сухого
Аннотация
На основе метода рассеянного поля для длины волны возбуждающего электромагнитного излучения 405 нм выполнено по алгоритму РБТО моделирование интенсивности локального поля вблизи поверхности сферических диэлектрических и металлических (серебряных) наночастиц радиуса 40 нм, размещённых в воздухе и вблизи плоской границы раздела двух диэлектриков. Получены оценки размеров области и величины максимального усиления ближнего поля, проведён анализ влияния границы раздела на распределение поля в моделируемой системе. Как отдельный случай исследованы нанокомпозиты на основе серебряных наночастиц и высококремнезёмных стёкол, полученные при использовании золь-гель технологий.
Ключевые слова: РБТБ-метод, метод рассеянного поля, диэлектрические и металлические наночастицы, ближнее электромагнитное поле.
Введение
Для моделирования проблем рассеяния электромагнитных плоских волн пространственно ограниченными объектами различных размеров широко используется метод FDTD (Finite-Difference TimeDomain [1, 2]). В рамках этого метода возбуждение плоских волн обычно реализуется по методам TF-SF (Total Field - Scattered Field) или SF (Scattered Field), см., например, [2 - 9] и ссылки в этих работах. При этом моделируются несколько типов электромагнитных полей: падающее - поле в отсутствие объектов, рассеивающих излучение (далее такие объекты для краткости будем называть «рассеивателями»); рассеянное - поле, возникающее в результате взаимодействия рассеивателей с падающим полем; полное поле
- векторная сумма падающего и рассеянного поля. В случае метода TF-SF расчётная область разбивается на подобласти полного и рассеянного поля, разделённые границей, генерирующей падающее поле и «проводящей» его через область полного поля. При использовании метода SF уравнения Максвелла преобразуются таким образом, что требуется расчёт FDTD-методом только рассеянных полей (в предположении, что падающее поле изначально известно).
Метод ТF-SF в классической формулировке имеет ряд ограничений в случае, если геометрия задачи предполагает наличие бесконечных в одном из направлений тел (например - плоская граница раздела двух сред) или пространственная форма возбуждения отличается от волны с плоским фронтом. Следует также отметить, что для полноценной реализации TFSF-метода или его модификации необходимо внести изменения непосредственно в исполняемый код программы, реализующий алгоритм FDTD, с последующей её отладкой и компиляцией.
Реализацию метода SF можно осуществить на уровне языка задания параметров, имеющегося в настоящее время у большинства из распространённых программных пакетов, использующих алгоритм FDTD, то есть усилиями непосредственно пользователя программного комплекса. Такой подход реали-
зован в настоящей работе при использовании распространённого пакета Meep [10], где в качестве входного языка подготовки заданий для выполнения FDTD-ядром использован язык функционального программирования Scheme.
Одним из ключевых аспектов применения любого численного метода является верификация его результатов и определение условий использования метода для решения конкретных прикладных задач. В настоящей работе FDTD SF-моделирование использовано в качестве основного метода исследования интенсивностей электромагнитного поля вблизи сферических диэлектрических и металлических наночастиц (НЧ) при условиях, близких к условиям реализации локализованных поверхностных плазмонных резонансов [11, 12]. Для простейших модельных систем (НЧ в воздухе, НЧ вблизи границы раздела двух диэлектриков) выполнено сравнение результатов, полученных по SF- и TF-SF-методам, с точным аналитическим решением согласно теории Ми [11, 12].
Как отдельный случай исследованы нанокомпозиты на основе металлических НЧ и высококремнезёмных стёкол, полученные при использовании золь-гель технологий [13-15]. Наноструктурированные стёкла, допированные НЧ восстановленных металлов, используются в качестве материалов, применяемых для управления параметрами лазерного излучения, оптической записи информации, создания сенсорных элементов [13]. Свойства этих материалов сильно зависят от формы, концентрации, взаимного расположения НЧ, свойств диэлектрической матрицы. В частности, установлено, что для золь-гель стёкол, содержащих НЧ восстановленных металлов, возможно создание условий, при которых наведённое локальное изменение показателя преломления остаётся стабильным в течение длительного времени [13-15], что представляет интерес для приложений оптики и на-ноплазмоники. Средний размер металлических НЧ (Ag, Cu, Ni, Co, Mn) в золь-гель стекле лежит в диапазоне от единиц до сотен нанометров. Экспериментальная реализация эмпирически подобранной струк-
туры для получения материала с заданными свойствами на основе таких нанокомпозитов является дорогостоящей и длительной. Компьютерное моделирование распределения электромагнитных полей (особенно в области ближнего поля) существенно сокращает время перебора оптимальных вариантов нано- и микроструктур и позволяет предложить новые материалы, формирующие заданный оптический отклик [16].
Таким образом, основными целями работы являются: во-первых, верификация и сравнение результатов РБТО-(по методикам 8Р, ТР-8Р) моделирования диэлектрических и металлических НЧ; во-вторых, исследование с помощью разработанных моделей характеристик локального оптического отклика (ближнего поля) базовых наноструктур для указанного класса композитных материалов, важных для оптических приложений.
Остальная часть статьи организована следующим образом. Основы метода 8Р приведены в разделе 1. Результаты численных экспериментов, графического анализа, а также сравнение данных (пространственных распределений интенсивности и характеристик усиления ближнего поля), полученных при РБТБ-(8Р, ТР-8Р) моделировании и использовании точного решения Ми для отдельной сферической диэлектрической или металлической (серебряной) НЧ, рассмотрены в разделе 2. Аналогичные данные, но при наличии плоской границы раздела двух диэлектриков (воздуха и золь-гель стекла) вблизи НЧ рассмотрены в разделе 3. В разделе 4 предложена простая аналитическая оценка влияния плоской границы раздела на ближнее поле НЧ (на основе результатов РБТБ-моделирования), которая согласуется с результатами раздела 3. В заключении сделаны выводы по полученным результатам с учётом их возможного использования при создании новых оптических материалов на основе золь-гель стёкол и НЧ.
1. Постановка задачи и основы реализации БЖ-метода
Рассмотрим постановку задачи моделирования при использовании 8Р-метода [5, 7 - 9]. Уравнения Максвелла в частотной области для линейной, изотропной и немагнитной среды можно записать в виде:
Ух Е(г, ю) = — Н(г, ю), (1)
с
УхН(г,ю) = -— е(г,ю)Е(г,ю), (2)
с
где /2 = -1, ю - частота, с - скорость света в вакууме, е (г, ю) - комплексная диэлектрическая проницаемость на частоте ю (мнимая часть е положительна для поглощающих сред), г - радиус-вектор точки наблюдения электромагнитного поля.
Вследствие линейности уравнений Максвелла электрическое и магнитное поле может быть представлено суммой падающего и рассеянного полей:
Е(г, ю) = Е1ПС (г, ю) + Е ^ (г, ю), (3)
Н(г, ю) = Ищс (г, ю) + Н ^ (г, ю), (4)
где Е1пс (г, ю) и Н1ИС (г, ю) - поля при отсутствии рассеивателей в моделируемой задаче (эти поля полагаются известными), Ес (г, ю) и Нха (г, ю) - поля, возникающие в результате взаимодействия падающих полей с рассеивателями (эти поля должны быть определены в результате моделирования).
Из уравнений (1 - 4) можно получить следующие
соотношения:
Ух ЕШс (г, ю) = — Нтс (г, ю), (5)
С
У х Нщс (г, ю) = -— &1пс (г, ю)ЕШс (г, ю), (6)
с
Ух Е _ (г, ю) = - Н _ (г, ю), (7)
с
У х Н^а (г, ю) = -— е(г, ю)Е5са (г, ю) + Б(г, ю), (8)
с
8(г, ю) = - — (е(г, ю) - Є і пс (г, ю))Е„с (г, ю). (9)
с
Мы предполагаем, что решения уравнений (5, 6) и вид зависимости е пс (г, ю) от частоты и координат известны. Выполняя операцию обратного преобразования Фурье Р -1 над каждым членом уравнений (7 - 9), получим набор уравнений Максвелла во временной области, которые требуется решить:
Ух Е _ (г, г) =1 Э, Н _ (г, г), (10)
с
Ух Н _ (г, г) =
= Р-11- ^е(г, ю^ (г, ю) | + 8(г, г)/ (11)
8(г, г) = Р-1 |- 1С (е(г, ю)-еШс (г, ю))Е пс (г, ю)|. (12)
В случае отсутствия у каждой из частей среды свойства частотной дисперсии первый член из правой части уравнения (11) и выражение (12), содержащие оператор Р-1, примут вид:
Р-1 <[- с e(Г, ю)Е"с° (Г, ^ = -1 £(Г)д(Е^са (г, г) ,
S(г, г) =- 1(е(г) -еШс (г))Э,Епс (г,1) .
с
В общем случае использование формул (10), (11) при реализации метода РБТБ в средах с частотной дисперсией требует задания аналитического выражения для зависимости диэлектрической проницаемости от частоты [2, 5]. В результате возможно использование метода 8Р и при учёте частотной дисперсии среды. Для расчёта вектора S (г, г) также необходимо
знать явный вид Einc (r, w). Возможна реализация метода SF при возбуждении: монохроматической волной; импульсом конечной длительности, определяемым зависимостью Einc (r, w). Для получения стационарного распределения поля на заданной частоте в первом случае достаточно установить время моделирования, существенно превышающее период возбуждающей волны, а во втором - сохранять результаты расчётов пространственного распределения полей в промежуточных временных точках и выполнить впоследствии соответствующее преобразование Фурье [2, 5]. Далее будем рассматривать случай стационарного монохроматического возбуждающего поля с
частотой w0 : Einc (r, t) = <c (r)eXp(-/'woi). В этом случае выражение (12) принимает явный вид:
S(г,t) = -^(e(r, Wo) - Einc(r))E0nc(r)exp(-iWot). (13)
c
Для расчёта распределения поля вблизи металлической НЧ выражение для частотной зависимости e (r, Wo) выбиралось в характерной для металлов форме Друде, параметры этой зависимости выражались через действительную и мнимую части диэлектрической проницаемости на частоте возбуждения в виде, предложенном в [17]. Путём непосредственного численного эксперимента нами установлено, что при таком выборе дисперсионной зависимости диэлектрической проницаемости итерационное решение уравнений (10) - (12), выполняемое методом FDTD, сходилось к некоторому стационарному решению.
Согласно формулам (7) - (9), уравнения Максвелла преобразованы таким образом, что FDTD-методом требуется определить только рассеянные поля (с учётом того, что изначально известны падающие поля). Диэлектрические проницаемости einc и e соответствуют геометриям моделируемой системы при отсутствии и наличии рассеивателей. Вектор-«источник» S (единственная величина в уравнениях (7), (8), зависящая от падающего поля) ненулевой только в областях модельного пространства, где значения einc и e отличаются. В частности, если рассеиватель размещён в областях с различными значениями einc (например, включает границу раздела двух сред для начальной постановки задачи), то такому рассеивателю соответствуют несколько источников (векторов S) в (8). Много источников также появляется при наличии многих дискретных рассеивателей в однородных областях модельного пространства.
Если поля Einc, Hinc, соответствующие начальной постановке задачи, характеризуются фазовым множителем плоской монохроматической волны exp [i (kr) - iwt], где k - волновой вектор (комплексный в поглощающих средах), то определение этих полей существенно упрощается для многих важных для приложений задач. В частности, для задач о рассеивании излучения объектами (например, НЧ), размещёнными вблизи плоской границы раздела двух однородных сред с локальными материальными параметрами, падающее поле имеет вид:
Епс (ю.г) = Кс <ютЕ> + г), (14)
[ ТЕпс(ю г),
где верхнее выражение соответствует области падающей на границу волны (первая среда), Е0пс (ю,г) , и волны, отражённой от границы, КЕ°пс (ю, г), а нижнее выражение - области прошедшей через границу волны (вторая среда); Я и Т - известные матрицы, характеризующие коэффициенты отражения и пропускания. Метод 8Р также достаточно просто обобщается на случай многослойных сред, включающих различные виды дискретных рассеивателей [7, 8].
Таким образом, при применении метода 8Р моделируемая задача может упрощаться за счёт использования известных точных аналитических решений граничных задач для начальной (без рассеивателей) постановки задачи. Согласно выражению (9), данные решения, а также геометрия задачи определяют один или несколько векторов-источников S. После задания величин S для моделирования рассеянного поля возможно использование стандартных пакетов прикладных программ, включая пакеты на основе РБТБ-алго-ритмов для численного решения уравнений (7), (8) конечно-разностными методами. Преимуществами 8Р подхода при РБТБ-моделировании также являются отсутствие ошибок дискретизации для падающего поля, возможности задания поглощающих граничных условий только для рассеянного поля и уменьшения размеров расчётной области при относительно малых рассеивающих объектах [5, 7 - 9].
2. Моделирование полей вблизи диэлектрической и металлической НЧ
Для оценки и сравнения точности РБТБ-моделирования методами 8Р (основной метод, используемый в данной работе) и ТР-8Р были проведены вычислительные эксперименты для исследования ближнего поля одной НЧ. На рис. 1 и 2 приведены распределения квадрата модуля амплитуды электрического поля (здесь и ниже будем использовать относительные единицы, | Е |2/| Е0 |2, где Е0 -векторная амплитуда падающей на НЧ волны) вблизи трёхмерных сферических диэлектрической (рис. 1) и металлической (рис. 2) НЧ радиусом Я = 40 нм. НЧ находятся в воздухе (е = 1), их центры совпадают с центром декартовой системы координат. На рис. 1, 2 также приведены для сравнения соответствующие кривые, полученные аналитическим методом (точное решение) на основе теории Ми [11, 12] (в областях вне НЧ). Для РБТБ-моделирования (8Р, ТР-8Р) использовался пакет Меер [10] и средства платформы удалённых параллельных вычислений папоНиВ [18, 19]. Различные кривые на рис. 1, 2 соответствуют изменению одной из координат при нулевых значениях двух других координат, т.е. рассматривается распределение поля в центральных сечениях сферической НЧ, параллельных координатным осям.
а) X, У,Ъ, нм
Рис. 1. Распределение нормированного квадрата модуля амплитуды электрического поля вблизи сферической диэлектрической НЧв воздухе для направлений вдоль одной из координатных осей и центрального сечения НЧ. Расчёты выполнены по методам SF (а), TF-SF (б), согласно теории Ми (для областей вне НЧ; на (а), (б) приведены для сравнения одинаковые кривые). На вставке схематически показаны направление волнового вектора (параллельно оси Z) и TE-поляризация падающей волны для данных на (а), (б). Здесь и на рисунках ниже на оси абсцисс выделена область, соответствующая расположению центрального сечения НЧ
а)
-100
-50
О
50
100 -100 -50 0 50 100
Х,Х2,нм б) Х,Х2,нм
Рис. 2. Распределение нормированного квадрата модуля амплитуды электрического поля вблизи сферической металлической (Ag) НЧ в воздухе для направлений вдоль одной из координатных осей и центрального сечения НЧ.
Остальные данные аналогичны приведённым для рис. 1
Возбуждение осуществлялось плоской, линейно поляризованной вдоль оси X волной с длиной 405 нм. Моделирование проводилось до момента времени установления стационарного распределения поля в расчётной области. Выбирались значения диэлектрической проницаемости материала НЧ eNP = 3,5191 + 0,5798i (для данных на рис. 1, 3) и eNP=-3,5191 + 0,5798i (для данных на рис. 2, 4, соответствует серебру [20], причём при указанных размерах НЧ и условиях возбуждения могут иметь место локализованные плазмонные резонансы НЧ [11, 12, 15]). Характеристики серебряных НЧ, материалов вне НЧ и условий возбуждения для данных рис. 2, 4 выбирались в областях значений параметров, реализуемых для нанокомпозитов на основе высококремнезёмных золь-гель стёкол [13 - 15].
Для FDTD-моделирования использовалась пространственная расчётная область в форме куба с ребром 240 нм. После установления стационарного распределения поля результаты моделирования практически не зависели от толщины согласующего слоя (Perfectly Matched Layer, PML [2, 10]) вокруг расчётной области. Поэтому для ускорения расчётов использовался тонкий PML с толщиной 10 нм. Шаг пространственной сетки соответствовал разрешению 600 точек на 1 мкм (около 1,67 нм) в направлениях всех координатных осей.
Данные РБТБ-моделирования на рис. 1б и 2б получены при реализации возбуждения плоской волной, заданием требуемого согласно методу ТР-8Р [2 - 4, 10, 21] распределения поля на границах расчётной области - на плоскостях х = -120 нм, у = -120 нм, г = -120 нм.
Согласно данным рис. 1, 2, при удалении от поверхности НЧ на расстояние порядка Я в направлении вдоль одной из координатных осей относительная амплитуда электрического поля быстро уменьшается (увеличивается) до 1 для направлений, параллельных (перпендикулярных) полю Е0. Усиление ближнего поля в окружающей НЧ среде имеет место как для металлической (рис. 2), так и для диэлектрической (рис. 1) НЧ (в основном, для направлений, параллельных плоскости поляризации возбуждающего излучения). При этом максимальные значения | Е |2/ | Е0 |2для серебряной НЧ на порядок и более превосходят аналогичные значения для диэлектрической НЧ.
Сравнение результатов моделирования с точным решением Ми (в области вне НЧ) показывает, что для диэлектрической НЧ (рис. 1) точность метода 8Р несколько выше, чем ТР-8Р. Для металлической НЧ (рис. 2) оба метода близки по точности при расчёте поля в направлении поляризации возбуждающей вол-
ны (вдоль оси X). Для направлений, перпендикулярных поляризации возбуждающей волны, результаты метода ТР-8Р менее точны (рис. 2б). Такая погрешность ТР-8Р моделирования обусловлена (в числе прочего) особенностями реализации ТР-8Р в пакете Меер (введением источников в виде плоскостей, а не в численной формулировке [2]).
\Е\2,отн. ед.
без границы, вдоль X с границей, вдаль X
\Е\,отн.ед.
а)
-100
-50
50
3. Моделирование полей вблизи НЧ: влияние плоской границы раздела Для исследования влияния границы раздела двух диэлектрических сред на ближнее поле НЧ было выполнено моделирование по методу 8Р для расчётных областей и условий возбуждения, схематически показанных на вставках к рис. 3, 4.
\Е\2,отн.ед.
с границей, вдоль У без границы, вдоль ¥ _
1оо
Х.У.нм
б)
Х,Хнм
Рис. 3. Распределение нормированного квадрата модуля амплитуды электрического поля вблизи сферической диэлектрической НЧ на границе раздела г = Я воздуха (область г < Я) и диэлектрика (г > Я, е2 = 2,4) для направлений вдоль одной из осей X, У и центрального сечения НЧ. Моделирование выполнено по методу БЕ для нормального (а) и наклонного падения (угол падения 45°) возбуждающей волны с ТН- (б) и ТЕ- (в) поляризацией. На вставках схематически показаны направления волнового вектора и поляризации падающей волны. Для сравнения также приведены соответствующие кривые, полученные при отсутствии границы раздела, в том числе по теории Ми для области вне НЧ (а)
\Е\2,отн.ед.
100
XXнм
б)
ХХнм
в)
-100
-50
О
50
100
Х,Хнм
Из первой среды (воздуха, е! = 1) на НЧ падает линейно поляризованная (ТН- или ТЕ-) волна с длиной 405 нм, которая отражается и преломляется на границе раздела. Для второй среды (в области г > Я) задана
Рис. 4. Распределение нормированного квадрата модуля амплитуды электрического поля вблизи сферической металлической (Ag) НЧ на границе раздела г = Я воздуха (г < Я, е1 = 1) и диэлектрика (г > Я, е2 = 2,4) для направлений вдоль осей X, У и центрального сечения НЧ. Остальные данные аналогичны приведённым для рис. 3
Предполагалось, что на границе раздела (плоскости г=Я) размещена диэлектрическая (рис. 3) или серебряная (рис. 4) НЧ с центром в начале координат и радиусом Я = 40 нм.
диэлектрическая проницаемость е2 = 2,4, соответствующая диапазону характеристик высококремнезёмных золь-гель стёкол, применяемых для создания на-нокомпозитных материалов [13 - 15]. Отметим, что
модельная структура «НЧ - плоская граница раздела диэлектриков» также может быть использована для исследования более сложных систем «НЧ - сферические микрочастицы кремнезёма (в составе аэросилов, силохромов)» [22 - 24]. Размеры таких микрочастиц могут существенно превышать 1 мкм, поэтому в области ближнего поля НЧ кривизна поверхности микрочастицы несущественна.
Исследовалось распределение ближнего поля на границе раздела вблизи НЧ для направлений вдоль одной из осей X, У. Таким образом, данные на рис. 3, 4 соответствуют распределению величины | Е |2/ | Е0 |2 при изменении координаты X (при у = 0, г=Я) и при изменении координаты У (при х = 0, г = Я) в присутствии или отсутствие границы раздела.
Рассмотрим характерные особенности изменения локального поля вблизи НЧ, обусловленные влиянием границы раздела, на основе сравнения данных рис. 1 и 3 (для диэлектрической НЧ), рис. 2 и 4 (для металлической НЧ). Видно, что (как при наличии, так и при отсутствии границы раздела) возмущения ближнего поля, существенно отличающиеся по величине от поля падающей волны, имеют место, в основном, на расстояниях порядка 3Я от центра НЧ. При этом направление наибольшего изменения локального поля определяется направлением поляризации возбуждающей волны. Также можно отметить (для обоих типов НЧ) уменьшение в 2 - 4 раза интенсивности локального поля вне НЧ при переходе от плоскости центрального сечения НЧ (данные рис. 1, 2) к плоскости, касающейся НЧ (рис. 3, 4).
Диэлектрическая НЧ. В этом случае (рис. 3) появление границы раздела ведёт к монотонному (в зависимости от расстояния) уменьшению ближнего поля в направлении поляризации падающей волны вблизи поверхности НЧ (на расстояниях порядка Я от поверхности НЧ). При этом значения величины | Е |2/ | Е0 |2 уменьшаются на 5-10% при нормальном падении (рис. 3а) и на 30-50% при наклонном падении (рис. 3б, 3в). Для наклонного падения и направлений, перпендикулярных поляризации падающей волны (X, У для данных рис. 3б и 3в соответственно), величина локального поля также уменьшается с появлением границы раздела, но в меньших пределах (порядка 10%).
Серебряная НЧ. В этом случае (рис. 4) появление границы раздела ведёт к изменениям ближнего поля НЧ, качественно похожим на случай диэлектрической НЧ. Однако области максимального возмущения ближнего поля НЧ становятся гораздо больше - до расстояний порядка 2Я от поверхности НЧ в направлениях поляризации возбуждающего излучения. Кроме того, в случае серебряной НЧ существенно возрастает отличие интенсивности локального поля для направлений вдоль и перпендикулярно поляризации падающей волны. Для направлений, перпендикулярных поляризации падающей волны, усиление локального поля практически отсутствует в области вне НЧ (рис. 4б) либо имеет место не усиление, а ослабление ближнего
поля вблизи НЧ по сравнению с падающей волной (рис. 4а, 4в).
При нормальном падении появление границы раздела может вести как к усилению локального поля вблизи поверхности НЧ, так и к его ослаблению с ростом расстояния от НЧ (рис. 3а, направление вдоль оси У, рис. 4а, направление вдоль оси X).
На рис. 5 приведена карта распределения величины |Е|2/|Е0|2 в плоскости поляризации возбуждающей волны, включающей центральное сечение серебряной НЧ и границу раздела «воздух - золь-гель стекло». Рассмотрены случай касания НЧ границы раздела (рис. 5а, эти данные дополняют результаты, приведённые на рис. 4а), а также случай наличия малого зазора Я/10 = 4 нм между поверхностью НЧ и границей раздела (рис. 5б). Геометрию системы на рис. 5б предлагается использовать для создания нанокомпо-зитных покрытий, являющихся эффективными плаз-монными поглотителями света [25] (при выполнении ряда дополнительных требований к параметрам системы, ведущих к усилению магнитных мод НЧ). Сравнение данных рис. 2 и рис. 4а, 5 показывает, что появление границы раздела вблизи поверхности серебряной НЧ (на расстояниях порядка нескольких нм) качественно модифицирует пространственное распределение ближнего поля НЧ.
Е
и
а) ^ ^ Рис. 5. Карта распределения значений величины \Е\2/\Е0\2 вблизи сферической серебряной НЧ на границе раздела воздуха (£]=1) и диэлектрика (е2=2,4) в плоскости поляризации падающей волны для центрального сечения
НЧ: НЧ касается границы раздела (а); зазор между поверхностью НЧ и границей раздела равен 4 нм (б). На вставках схематически показаны направление волнового вектора и поляризация падающей волны (соответственно перпендикулярно и параллельно границе раздела)
Максимальное усиление локального поля имеет место не вблизи концов отрезка - диаметра НЧ, параллельного плоскости поляризации возбуждающей волны (что характерно для отдельной НЧ, рис. 1, 2), а вблизи точки касания НЧ и границы раздела (рис. 4а, 5а) или зазора между НЧ и границей раздела (рис. 5б). При этом геометрия системы без зазора (рис. 5а) предпочтительнее для получения больших областей эффективного усиления ближнего поля, чем геометрия с зазором (рис. 5б).
Согласно данным рис. 5, наибольший линейный размер области (с центром в точке касания НЧ и границы раздела), соответствующей максимальному усилению локального поля (до значений |Е|2/|Е0 |2 @ 60 и более), можно оценить как (1^2) Я в
направлении поляризации возбуждающего излучения. В направлении, перпендикулярном плоскости поляризации падающей волны, область максимального усиления ближнего поля имеет меньшие размеры: порядка Я и менее (рис. 4, это также подтверждается данными, аналогичными рис. 5 для этого случая, которые здесь не приводятся). Можно предположить, что именно в этой области (причём как на границе раздела, так и вне её) целесообразно размещение различных фотоактивных центров (молекул, квантовых точек) для приложений нано-плазмоники, использующих металлические НЧ в качестве «генераторов» больших локальных полей вблизи границ раздела диэлектриков.
4. Оценка влияния наличия плоской границы раздела на поле вблизи НЧ Предположим, что сферическая НЧ радиуса Я с центром в начале координат и некоторые поглощающие (пробные) центры (например, фотоактивные центры, молекулы, квантовые точки) расположены на участке плоской поверхности Ж (границы раздела двух сред), характеризуемой уравнением 7 = Я в декартовой системе координат. При моделировании электромагнитных полей методом конечных разностей используется разбиение пространственной счётной области с помощью прямоугольной трёхмерной сетки. С учётом этого определим коэффициент К, позволяющий оценить «в среднем» влияние заданного участка поверхности Ж на интенсивность поля вблизи НЧ, следующим образом:
1 _ _ Е( х, у
К(2 = Я) = - X ------------------------
^ Iх,Му,\<а Еге^(х,
(15)
У!, Я)
Таблица. Значения коэффициента К в зависимости от параметров моделируемой системы
где х,, у,, 2 - координаты узла сетки, лежащего на поверхности Ж, Е (х,, у¡, 2) и Егеу(х,, у,, 2) - полученные в результате моделирования напряжённости электрического поля вблизи НЧ в присутствие и отсутствие границы раздела Ж соответственно. В выражении (15) в качестве участка поверхности Ж выбран квадрат со стороной, равной 2а, симметрично расположенный относительно НЧ, N - количество расчётных точек на этом участке поверхности. Согласно данным графического анализа, приведённым выше (рис. 3 - 5), размеры области, где рассеянное НЧ поле вносит основной вклад в полное поле вблизи поверхности Ж, можно оценить условием 2а < 6Я .
Таким образом, значение коэффициента К показывает, насколько интенсивность поля вблизи НЧ и границы раздела, определяющая вероятность поглощения излучения фотоактивными центрами, размещёнными вне НЧ на поверхности Ж, изменяется в среднем на заданной квадратной области при появлении плоской границы раздела двух однородных сред в области ближнего поля НЧ.
С использованием метода оценено влияние границы раздела двух диэлектрических сред на распределение поля вблизи диэлектрической и металлической НЧ. Результаты расчётов при выборе значений 2а = 6Я , 4Я , 3Я приведены в таблице ниже. При проведении расчётов использовались следующие значения параметров, не указанных в таблице: радиус НЧ Я = 40 нм, первая среда - воздух (е1 = 1), вторая среда -высококремнезёмное золь-гель стекло (е2=2,4), для всех случаев бралось Іт (є№) = 0,5798 (изменялась вещественная часть диэлектрической проницаемости материала НЧ, Яе (є№)), возбуждение осуществлялось плоской, линейно поляризованной волной с длиной 405 нм.
2
Тип НЧ Яе(е ^) Поляризация Угол падения Коэффициент влияния поверхности К
2а = 6Я 2а = 4 Я 2а = 3Я
Метал- личе- ская -2,5 - 0 1,1578 1,2450 1,4657
-3 - 0 1,1681 1,3207 1,5989
-3,51 - 0 1,1707 1,4777 1,9106
-5,5 - 0 1,1080 1,4731 1,9122
-3,51 ТЕ 45° 0,8548 1,0779 1,3870
-3,51 ТМ 45° 0,7862 0,9120 1,1063
Ди- элек- тричес- кая 2,5 - 0 0,8983 0,9206 0,9506
3 - 0 0,8992 0,9270 0,9642
3,51 - 0 0,9006 0,9338 0,9779
5,5 - 0 0,9069 0,9566 1,0211
3,51 ТЕ 45° 0,6962 0,7244 0,7617
3,51 ТМ 45° 0,8407 0,8449 0,8542
Как видно из таблицы, появление границы раздела может как усиливать (К > 1), так и ослаблять (К < 1) поле вблизи НЧ. Для диэлектрической НЧ изменение поля небольшое: |К-1| < 0,1 при нормальном падении, |К-1| < 0,3 при наклонном падении (угол падения 45°). Для такой НЧ при наклонном падении наличие границы раздела оказывает более существенное влияние и, по-видимому, в большинстве случаев ослабляет локальное поле на границе раздела вблизи НЧ. Также видно, что
уменьшение площади участка поверхности, влияние которого исследуется, от 36Я2 до 9Я2 ведёт к небольшому (порядка 10%) увеличению значений К.
Для металлической НЧ характерно существенное усиление влияния границы раздела с уменьшением площади участка поверхности (до значений |К-1| < 0,9 для данных таблицы). Видно, что и в этом случае усиление локального поля максимально при нормальном падении возбуждающего излучения.
Отметим, что данные таблицы учитывают вклад в усиление локального поля как участков поверхности, близких к центральному сечению НЧ (параллельному плоскости поляризации падающей волны), где значения К гораздо больше приведённых в таблице (см. данные рис. 4, 5), так и вклад других участков области поверхности с меньшими (близкими к 1) значениями К. Приведённые данные дополняют результаты графического анализа (раздел 4) и подтверждают сделанные оценки размеров области на границе раздела вблизи НЧ, соответствующей наиболее эффективному усилению поля: 2а < 3Я и 2а < Я для направлений соответственно вдоль и перпендикулярно поляризации возбуждающего излучения.
Заключение
Результаты выполненного моделирования и графического анализа позволяют уточнить ряд аспектов, связанных как с методикой РБТО-моделирования, так и с исследованием характеристик локального оптического отклика рассматриваемых наноструктур.
Сравнение данных моделирования, полученных по методикам ТБ-8Р и 8Б, с точным решением Ми (раздел 2) показывает, что точность методики 8Б для исследованных геометрий задач и значений параметров оказывается выше, чем ТБ-8Р (в частности, существенно выше для диэлектрических НЧ).
Исследовано изменение характеристик усиления ближнего поля диэлектрической и металлической НЧ в воздухе, обусловленное наличием границы раздела двух диэлектриков (воздух - высококремнезёмное золь-гель стекло) в области ближнего поля НЧ (разделы 2, 3). Показано, что (как при наличии, так и при отсутствии границы раздела) возмущения ближнего поля НЧ по сравнению с полем падающей волны имеют место, в основном, на расстояниях порядка 3Я от центра НЧ. При наличии границы раздела область максимального усиления локального поля находится вблизи точки касания НЧ и границы раздела или зазора между НЧ и границей раздела. Для серебряной НЧ размеры данной области на границе раздела можно оценить как (1 + 2)Я и Я для направлений параллельно и перпендикулярно плоскости поляризации возбуждающей волны соответственно.
В разделе 4 предложена простая методика оценки влияния границы раздела на характеристики ближнего поля НЧ и показано, что появление такой границы может, в среднем, как усиливать, так и ослаблять локальное поле НЧ. Данные изменения поля невелики для диэлектрических НЧ и могут быть существенны для металлических НЧ, особенно в случае нормального падения возбуждающего излучения на систему «НЧ - граница раздела двух диэлектриков».
Полученные результаты также могут быть использованы для оптимизации характеристик нанокомпозитов «высококремнезёмные золь-гель стекла - НЧ благородных металлов» для приложений оптики и наноплазмоники [13 - 15].
Благодарности Работа выполнена при поддержке программы
ГКПНИ «Функциональные и машиностроительные
материалы, наноматериалы» (задание 2.4.02.01) Республики Беларусь.
Литература
1. Yee, K.S. Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media / K.S. Yee // IEEE Trans. Antenn. Propag. - 1966. -Vol. 14, N 3. - P. 302-307.
2. Taflove, A. Computational Electrodynamics: the Finite-Difference Time-Domain method / A. Taflove, S.C. Hag-ness. - 2nd ed. - Artech House, 2000. - 866 p.
3. Umashankar, K.R. A novel method to analyse electromagnetic scattering of complex objects / K.R. Umashankar, A. Taflove // IEEE Trans. Electromagn. Compat. - 1982. -Vol. 24, N 4. - P. 397-405.
4. Schneider, J.B. Plane waves in FDTD simulations and a nearly perfect total-field/scattered-field boundary / J.B. Schneider // IEEE Trans. Antenn. Propag. - 2004. - Vol. 52, N 12. -P. 3280-3287.
5. Kunz, K.S. The Finite Difference Time Domain Method for Electromagnetics / K.S. Kunz, R.J. Luebbers. - CRC Press, 1993. - 461 p.
6. Котляр, В. В. Численное решение уравнений Максвелла в задачах дифракционной оптики / В.В. Котляр // Компьютерная оптика. - 2006. - Вып. 29. - С. 24-40.
7. Kong, S.-C. ADE-FDTD scattered-field formulation for dispersive materials / S.-C. Kong, J.J. Simpson, V. Backman // IEEE Microwave and Wireless Comp. Lett. - 2008. -Vol. 18, N 1. - P. 4-6.
8. Olkkonen, J. FDTD scattered field formulation for scatter-ers in stratified dispersive media / J. Olkkonen // Optics Express. - 2010. - Vol. 18, N 5. - P. 4380-4389.
9. Challener, W.A. Scattered field formulation of finite difference time domain for a focused light beam in dense media with lossy materials / W.A. Challener, I.K. Sendur // Optics Express. - 2003. - Vol. 11, N 23. - P. 3160-3170.
10. Oskooi, A.F. MEEP: A flexible free-software package for electromagnetic simulations by the FDTD method / A.F. Os-kooi, D. Roundy, M. Ibanescu, P. Bermel, J.D. Joannopou-los, S.G. Johnson // Computer Physics Communications. -
2010. - Vol. 181, Issue 3. - P. 687-702.
11. Bohren, C.F. Absorption and Scattering of Light by Small Particles / C.F. Bohren, D.R. Human. - NY.: Wiley, 1998. - 544 p.
12. Kreibig, U. Optical Properties of Metal Clusters / U. Krei-big, M. Vollmer. - Wiley, Berlin, Heidelberg: Springer, 1995. - 553 p.
13. Алексеенко, А. А. Функциональные материалы на основе диоксида кремния, получаемые золь-гель методом / А. А. Алексеенко, А. А. Бойко, Е.Н. Подденежный. - Гомель: ГГТУ им. П.О. Сухого, 2008. - 183 с.
14. Yeshchenko, O. Optical properties of sol-gel fabricated Co/SiO2 nanocomposites / O. Yeshchenko, I. Dmitruk,
A. Alexeenko, A. Dmitruk, V. Tinkov // Physica E. - 2008.
- Vol. 41. - P. 60-65.
15. Yeshchenko, O.A. Size-dependent surface-plasmon-enhan-ced photoluminescence from silver nanoparticles embedded in silica / O.A. Yeshchenko, I.M. Dmitruk, A.A. Alexeenko, M.Yu. Losytskyy, A.V. Kotko, A.O. Pinchuk // Phys. Rev.
B. - 2009. - Vol. 79. - P. 235438-1-8.
16. Курочка, К.С. Организация распределённых вычислений при моделировании электромагнитных полей в композитных материалах, содержащих наночастицы металлов /
К.С. Курочка, О.Д. Асенчик, Е.Г. Стародубцев // Между-нар. конгресс по информатике: информационные системы и технологии: матер. междунар. науч. конгр. Минск, 31 окт. - 3 нояб. 2011 г.: в 2 ч. Ч. 2 / редкол.: С.В. Абламейко (отв. ред.) [и др.]. - Минск: БГУ, 2011. - С. 112-116.
17. Johnson, P.B. Optical constants of the noble metals / P.B. Johnson, R.W. Christy // Phys. Rev. B. - 1972. - Vol. 6.
- P. 4370-4379.
18. http://nanohub.org/resources/Meep [Electronic resource]. -Date of access: 06.03.2013.
19. Кильдишев, А.В. Трансформационная оптика и метаматериалы / А.В. Кильдишев, В.М. Шалаев // Успехи физических наук. - 2011. - Т. 181, № 1. - C. 59-70.
20. Rakic, A.D. Optical properties of metallic films for vertical-cavity optoelectronic devices / A.D. Rakic, A.B. Djurisic, J.M. Elazar, M.L. Majewski // Applied Optics. - 1998. -Vol. 37, N 22.- P. 5271-5283.
21. http://www.mail-archive.com/[email protected]/info.html [Electronic resource]. - Date of access: 06.03.2013.
22. Муха, Ю.П. Усиление поглощения и флуоресценции молекул родамина 6Ж вблизи наночастиц золота в матрице SiO2 / Ю.П. Муха, А.М. Ерёменко, Н.П. Смирнова, М.Я. Валах, В.И. Джаган // Химия, физика и технология поверхности. - 2011. - Т. 2, № 3. - С. 284-288.
23. Брюханов, В.В. Взаимодействие поверхностных плазмонов наночастиц серебра на силохроме и шероховатых пленках серебра с электронно-возбужденными адсорба-тами молекул родамина 6Ж / В.В. Брюханов, Н.С. Тихомирова, Р.В. Горлов, В. А. Слежкин // Известия КГТУ. -
2011. - № 23. - С. 11-17.
24. Брюханов, В.В. Плазмонное усиление флуоресценции органолюминофоров в полимере и на поверхности кремнезема / В.В. Брюханов, В.А. Слежкин, Н.С. Тихомирова, Р.В. Горлов // Вестник Балтийского федерального университета им. И. Канта. - 2012. - Вып. 4. - С. 52-59.
25. Albooyeh, M. Huge local field enhancement in perfect
plasmonic absorbers / M. Albooyeh, C.R. Simovski // arXiv:1203.2100v3 [cond-mat.mtrl-sci] (2012) [Electronic resource]. - Mode of access: http://arxiv.org/pdf/
1203.2100v3.pdf - Date of access: 06.03.2013.
References
1. Yee, K.S. Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media / K.S. Yee // IEEE Trans. Antenn. Propag. - 1966. -Vol. 14, N 3. - P. 302-307.
2. Taflove, A. Computational Electrodynamics: the Finite-Difference Time-Domain method / A. Taflove, S.C. Hag-ness. - 2nd ed. - Artech House, 2000. - 866 p.
3. Umashankar, K.R. A novel method to analyse electromagnetic scattering of complex objects / K.R. Umashankar,
A. Taflove // IEEE Trans. Electromagn. Compat. - 1982. -Vol. 24, N 4. - P. 397-405.
4. Schneider, J.B. Plane waves in FDTD simulations and a nearly perfect total-field/scattered-field boundary / J.B. Schneider // IEEE Trans. Antenn. Propag. - 2004. - Vol. 52, N 12. -P. 3280-3287.
5. Kunz, K.S. The Finite Difference Time Domain Method for Electromagnetics / K.S. Kunz, R.J. Luebbers. - CRC Press, 1993. - 461 p.
6. Kotlyar, V.V. Numerical solution of Maxwell’s equations in diffraction optics problems / V.V. Kotlyar // Computer Optics. - 2006. - V. 29. - P. 24-40. - (In Russian).
7. Kong, S.-C. ADE-FDTD scattered-field formulation for dispersive materials / S.-C. Kong, J.J. Simpson, V. Backman //
IEEE Microwave and Wireless Comp. Lett. - 2008. -Vol. 18, N 1. - P. 4-6.
8. Olkkonen, J. FDTD scattered field formulation for scatter-ers in stratified dispersive media / J. Olkkonen // Optics Express. - 2010. - Vol. 18, N 5. - P. 4380-4389.
9. Challener, W.A. Scattered field formulation of finite difference time domain for a focused light beam in dense media with lossy materials / W.A. Challener, I.K. Sendur // Optics Express. - 2003. - Vol. 11, N 23. - P. 3160-3170.
10. Oskooi, A.F. MEEP: A flexible free-software package for electromagnetic simulations by the FDTD method / A.F. Oskooi, D. Roundy, M. Ibanescu, P. Bermel, J.D. Joannopou-los, S.G. Johnson // Computer Physics Communications. -2010. - Vol. 181, Issue 3. - P. 687-702.
11. Bohren, C.F. Absorption and Scattering of Light by Small Particles / C.F. Bohren, D.R. Human. - NY.: Wiley, 1998. -544 p.
12. Kreibig, U. Optical Properties of Metal Clusters / U. Krei-big, M. Vollmer. - Wiley, Berlin, Heidelberg: Springer, 1995. - 553 p.
13. Alexeenko, A.A. Functional Materials on Basis of Silicon Dioxide Obtained by Sol-gel Method / A.A. Alexeenko, A.A. Bojko, E.N. Poddenezhny. - Gomel: P.O. Sukhoi GSTU, 2008. - 183 p. - (In Russian).
14. Yeshchenko, O. Optical properties of sol-gel fabricated Co/SiO2 nanocomposites / O. Yeshchenko, I. Dmitruk,
A. Alexeenko, A. Dmitruk, V. Tinkov // Physica E. - 2008.
- Vol. 41. - P. 60-65.
15. Yeshchenko, O.A. Size-dependent surface-plasmon-enhan-ced photoluminescence from silver nanoparticles embedded in silica / O.A. Yeshchenko, I.M. Dmitruk, A.A. Alexeenko, M.Yu. Losytskyy, A.V. Kotko, A.O. Pinchuk // Phys. Rev.
B. - 2009. - Vol. 79. - P. 235438-1-8.
16. Kurochka, K.S. Organization of distributed calculations at electromagnetic field modeling in composite materials including metal nanoparticles / K.S. Kurochka, O.D. Asen-chik, E.G. Starodubtsev // Proc. of Intern. Congress on Informatics: Information Systems and Technologies. - Minsk, 31 Oct. - 3 Nov. 2011: Part 2 / Ed.: S.V. Ablamejko [et al.].
- Minsk: BSU, 2011. - P. 112-116. - (In Russian).
17. Johnson, P.B. Optical constants of the noble metals / P.B. Johnson, R.W. Christy // Phys. Rev. B. - 1972. - Vol. 6.
- P. 4370-4379.
18. http://nanohub.org/resources/Meep [Electronic resource]. -Date of access: 06.03.2013.
19. Kildishev, A.V. Transformation optics and metamaterials / A.V. Kildishev, V.M. Shalaev // Physics-Uspekhi. - 2011. -V. 181, N 1. - P. 59-70. - (In Russian).
20. Rakic, A.D. Optical properties of metallic films for vertical-cavity optoelectronic devices / A.D. Rakic, A.B. Djurisic, J.M. Elazar, M.L. Majewski // Applied Optics. - 1998. -Vol. 37, N 22. - P. 5271-5283.
21. http://www.mail-archive.com/[email protected]/info.html [Electronic resource]. - Date of access: 06.03.2013.
22. Mukha, I. The absorption and fluorescence enhancement of Rhodamine 6G molecules near gold nanoparticles in SiO2 matrix / I. Mukha, A. Eremenko, N. Smirnova, M. Valakh, V. Dzhagan // Chemistry, Physics and Surface Technology.
- 2011. - V. 2, N 3. - P. 284-288. - (In Russian).
23. Bryukhanov, V.V. Interaction of surface plasmons of silver nanoparticles on silochrom and rough silver films with electron-induced adsorbates of Rhodamine 6G molecules / V.V. Bryukhanov, N.S. Tikhomirova, R.V. Gorlov, V.A. Slezhkin // Izvestiya KSTU. - 2011. - N 23. - P. 11-17.
- (In Russian).
24. Bryukhanov, V.V. Plasmon enhancement of organic luminophors fluorescence in polymer and on silica surface / V.V. Bryukhanov, V.A. Slezhkin, N.S. Tikhomirova, R.V. Gorlov // Vestnik Baltij skogo Federal’nogo Universiteta imeni I. Kanta. - 2012. - N. 4. - P. 52-59. - (In Russian).
25. Albooyeh, M. Huge local field enhancement in perfect plasmonic absorbers / M. Albooyeh, C.R. Simovski // arXiv:1203.2100v3 [cond-mat.mtrl-sci] (2012) [Electronic resource]. - Mode of access: http://arxiv.org/pdf/
1203.2100v3.pdf - Date of access: 06.03.2013.
APPLICATION OF THE SCATTERED FIELD TECHNIQUE FOR FDTD MODELING OF ELECTROMAGNETIC FIELDS NEAR DIELECTRIC AND METALLIC NANOPARTICLES
O.D. Asenchik, E. G. Starodubtsev P. O. Sukhoi Gomel State Technical University
Abstract
Modeling of local electromagnetic field intensity near surfaces of spherical dielectric and metal (silver) nanoparticles of 40 nm radius (which are placed in air and near a two dielectrics plane boundary) is carried out using the scattered field technique and FDTD algorithm for exciting radiation wavelength 405 nm. Estimations of the region and the quantity of near field maximal enhancement are obtained, and effects of the boundary on the field distribution in the system are analyzed. As a special case, we investigate nanocomposites on the basis of silver nanoparticles and high-silicious glasses obtained by sol-gel technologies.
Key words: FDTD-method, scattered field technique, dielectric and metal nanoparticles, near electromagnetic field.
Сведения об авторах Асенчик Олег Даниилович, 1967 года рождения. В 1989 году окончил Белорусский ордена Трудового Красного Знамени государственный университет имени В.И. Ленина (ныне - Белорусский государственный университет) по специальности «Физика». Кандидат физико-математических наук (1998 год) по специальности «Оптика», доцент, работает первым проректором учреждения образования «Г омельский государственный технический университет имени П.О. Сухого» (ГГТУ им П.О. Сухого), доцентом кафедры информационных технологий ГГТУ им. П.О. Сухого. Область научных интересов: фотопроцессы в конденсированных средах, математическое моделирование, нанооптика и на-нофотоника.
E-mail: [email protected]. by .
Oleg Daniilovich Asenchik (b. 1967) graduated from Belarusian State University in 1989 with degree of "Physics". He received his Candidate in Physics & Maths (1998) from Belarusian State University, associate professor, works as the first vice-rector of P. Sukhoi Gomel State Technical University (GSTU), assistant professor of GSTU Information Technology Chair. Research interests are photoprocesses in condensed matter, mathematical modeling, nanooptics and nanophotonics.
Стародубцев Евгений Генрихович, 1967 года рождения. В 1991 году окончил Белорусский государственный университет по специальности «Физика». Кандидат физикоматематических наук (1996 год) по специальности «Оптика», доцент, работает доцентом кафедры информационных технологий ГГТУ им. П.О. Сухого. Научные интересы: электродинамика, фотоника, математическое моделирование.
E-mail: [email protected] .
Evgenii Genrichovich Starodubtsev (b. 1967) graduated from Belarusian State University in 1991 with degree of "Physics". He is Candidate in Physics & Maths (1996), associate professor, works as assistant professor of GSTU Information Technology Chair. Research interests are electrodynamics, photonics, mathematical modeling.
Поступила в редакцию 24 июня 2013 г.
