Управление излучением и рассеянием света при помощи плазмонных наноструктур тема диссертации и автореферата по ВАК РФ 01.04.05, доктор наук Гинзбург Павел Борисович

• Гинзбург Павел Борисович
• доктор наук
• 2020

Оглавление диссертации доктор наук Гинзбург Павел Борисович

РЕФЕРАТ

SYNOPSIS

1. LOCALIZED PLASMON RESONANCES

1.1 Localized plasmon resonances as a topological problem

1.2 Nonlinear and nonlocal description of conduction electron plasma

1.3 Quantum effects, inspired by plasmonic particles

1.4 Other nonlinear effects

2. TAILORING LIGHT-MATTER INTERACTIONS WITH HYPERBOLIC METAMATERIALS

2.1 Optical properties of hyperbolic metamaterials

2.2 Light emission in the vicinity of metamaterials

2.3 Energy transfer assisted by metamaterials

2.3. Linear and nonlinear scattering, tailored by hyperbolic metamaterials

3. OPTO-MECHANICAL MANIPULATION WITH AUXILIARY STRUCTURED MEDIA

3.1 Opto-mechanical metamaterials

3.2 Opto-mechanical metasurfaces

3.3 Flat optical elements for opto-mechanical manipulation

3.4 Quantum opto-mechanical effects

4. EMULATION OF COMPLEX OPTICAL PHENOMENA WITH RF CIRCUITRY

4.1 Artificial magnon resonance

4.2 Hyperbolic metamaterials, emulated with spoof plasmons

4.3 Huygens elements and asymmetric electromagnetic responses

4.4. Scattering suppression in metamaterials

Summary and conclusions

Bibliography

APPENDIX A. MAIN JOURNAL PAPERS

РЕФЕРАТ

Введение диссертации (часть автореферата) на тему «Управление излучением и рассеянием света при помощи плазмонных наноструктур»

Общая характеристика работы

Взаимодействие электромагнитного поля и вещества является одним из самых фундаментальных явлений, которое позволяет описывать множество природных эффектов, окружающих нас ежедневно. Их классическое описание основано на макроскопических уравнениях Максвелла, дополненных материальными соотношениями, учитывающими влияние материалов. Эффекты хроматической и пространственной дисперсий могут быть введены с помощью классических феноменологических моделей. Более подробное описание материальных откликов требует введения нерелятивистских квантово-механических инструментов. Строгое квантово-механическое описание требует квантования вещества, полей и взаимодействия между ними. Достаточно интересным фактом является то, что большинство явлений, классических или квантовых, можно описать с помощью электромагнитных функций Грина. В то время как этот инструмент довольно широко используется в классическом электромагнетизме, он также может быть использован в макроскопической квантовой механике путем введения операторов квантового шума. Поэтому большая часть эффектов в этой работе может быть проанализирована с помощью функций Грина, которые позволяют ввести ряд хорошо известных и широко используемых параметров, таких как (1) плотность электромагнитных состояний, (п) эффекты Парселла и усиление излучения, (ш) поляризуемость частиц, (IV) оптомеханические силы, и многих других. Следовательно, способность модифицировать функции Грина с помощью структурированной среды открывает широкий спектр возможностей, как в фундаментальной, так и в прикладной науке. В частности, здесь будет осуществляться модификация функций Грина с помощью малых субволновых резонаторов (плазмонных антенн) и их массивов, образующих искусственно созданные среды (метаматериалы).

В данной работе основное внимание будет уделено исследованию оптических свойств малых частиц (наноантенн), массивов частиц (метаматериалов и метаповерхностей) и их влияния на ряд классических и квантовых взаимодействий, которые включают: (1) рассеяние волн, (п) спонтанное излучение, (ш) спонтанное двухфотонное излучение, (IV) нелинейную генерацию гармоник, (V) оптомеханические взаимодействия и некоторые другие. Анализ этих совершенно разных физических явлений имеет схожую теоретическую основу, а именно - анализ функций

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

Теоретическая и практическая значимость работы

С фундаментальной точки зрения, основными достижениями данной работы являются теоретические и экспериментальные исследования локализованных плазмонных резонансов, где показано, что локальная геометрия играет ключевую роль в спектральном сдвиге откликов в инфракрасный спектральный диапазон. В то время как предыдущие исследования в этой области были сосредоточены на изучении удлиненных частиц, автором было продемонстрировано, что использование вогнутых геометрий позволяет сохранять компактность плазмонных частиц и, в то же время, позволить им резонировать на относительно больших длинах волн. На практике этот результат очень важен для целого ряда биологических и биомедицинских исследований, где плазмонные частицы используются как пассивные метки и даже в качестве платформы для лечения рака. Особое значение здесь имеет инфракрасный спектральный диапазон, так как подавляющее большинство биологических тканей относительно прозрачны для этих частот (так называемое окно прозрачности ближнего инфракрасного диапазона 650-1300 нм). Тщательное проектирование плазмонных резонансов также позволяет повысить эффективность эмиссии излучателей, находящихся вблизи — это эффект Парселла, первоначально разработанный для анализа электромагнитных взаимодействий с веществом в закрытых резонаторах. Задача применения квантового формализма для открытия дисперсионных и поглощающих плазмонных полостей уходит корнями в основы квантовой механики, где эрмитовы операторы описывают измеримые величины (например, гамильтонианы систем). На основе последних статей и предыдущих изысканий данная диссертация развивает новый формализм, который способен описывать механизмы взаимодействия излучателей с дисперсными плазмонными частицами, обладающими потерями. Операторы шума Ланжевена введены для этой цели и позволяют построить квантовый формализм, основанный на классических электромагнитных функциях Грина, которые

используются в этой работе. Являясь универсальным инструментом, разработанный формализм может быть применен в случаях, в которых задействованы массивы плазмонных частиц. Резонансы частиц в этих массивах могут гибридизироваться, и формировать коллективные моды, напоминающие свойства непрерывных сред: это ключевая идея метаматериалов - искусственных структур с новыми и необычными (в данном случае электромагнитными) характеристиками. В некоторых структурах, таких как гиперболические метаматериалы, усиление Парселла, рассчитанное с использованием формализма функций Грина, может быть чрезвычайно высоким, что теоретически приводит к нефизическим результатам. Одна из фундаментальных целей этой работы состоит в том, чтобы выявить практические ограничения скорости излучения и показать, что она зависит от степени гранулярности (метод реализации конкретного метаматериала) метаматериала, а именно от размера его элементарной ячейки по отношению к рабочей длине волны. С прикладной точки зрения, повышение эффективности излучения с помощью наноструктурированных сред очень важно для ряда оптоэлектронных устройств, таких как светоизлучающие диоды, квантовые источники одиночных фотонов и запутанных фотонных пар, и ряда других приложений. В частности, в этой работе показаны несколько примеров, относящихся к вышеупомянутым областям. Первый пример — это усиление процесса спонтанного двухфотонного излучения, которое генерирует запутанные пары с переменными энергии-времени. Фундаментальный процесс спонтанного двухфотонного излучения из полупроводников также является одним из открытий, описанных в этой диссертации. Второй пример — это демонстрация рекордного значения коллективного усиления Парселла от молекул красителя, помещенных в гиперболический метаматериал. Этот подход даёт новые возможности в области разработки новых типов органических светодиодов и других светоизлучающих устройств. Применения в биологии также могут быть кратко упомянуты в контексте гиперболических метаматериалов. В частности, процесс резонансной передачи энергии между молекулами (также называемого биологической линейкой, поскольку он позволяет измерять наноразмерные расстояния между различными молекулярными комплексами) может эффективно контролироваться с помощью метаматериалов, которые позволяют значительно расширить соответствующие расстояния измерения, как показано в данной диссертации.

Другое значимое приложение, которое можно проанализировать с помощью разработанных инструментов (функций Грина), — это оптомеханическое манипулирование. Основываясь на ряде фундаментальных эффектов взаимодействия света с веществом, оптомеханический контроль стал мощным инструментом для ряда биологических исследований в области кинетики

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

Последняя глава диссертации посвящена разработке инструментов эмуляции, позволяющих изучать сложные явления взаимодействия света с веществом в структурированной среде. Весьма примечательно, что классические электромагнитные функции Грина и связанная с ними локальная плотность состояний могут быть непосредственно измерены в гигагерцовом диапазоне (и в диапазоне более низких частот) с помощью, так называемых S-параметров антенны. Здесь небольшие излучающие диполи (электрические или магнитные) эмулируют квантовые излучатели. Диссертация развивает несколько концепций локализованных резонансов в субволновой геометрии, повторяя феномен локализованных плазмонных резонансов. С практической точки зрения, эти небольшие резонирующие объекты могут использоваться в качестве компактных антенн и рассеивающих устройств, которые чрезвычайно востребованы в ряде современных технологий беспроводной связи. Другим важным примером радиочастотной эмуляции, продемонстрированным в диссертации, является реплика гиперболических метаматериалов в гигагерцовом диапазоне. В подавляющем большинстве экспериментально продемонстрированных гиперболических метаматериалов оптического диапазона используются плазмонные материалы с отрицательной диэлектрической проницаемостью. Эти металлические компоненты являются практически идеальными электрическими проводниками на более низких частотах, что не позволяет достичь желаемых гиперболических откликов. Несмотря на многочисленные попытки реализовать гиперболические метаматериалы в гигагерцовом диапазоне, эта проблема оставалась нерешенной до недавнего времени. Наша работа демонстрирует новый подход к решению этой проблемы с использованием, так называемых спуф-плазмонов — волн, которые распространяются на гофрированных поверхностях. С прикладной точки зрения, гиперболические метаматериалы в гигагерцовом диапазоне открывают новые возможности для реализации компактных антенн и соответствующих сверхширокополосных устройств для беспроводной связи.

Актуальность и степень разработанности темы

Контроль взаимодействия света и вещества с помощью наноструктурированной электромагнитной среды. Возможность контролировать световые взаимодействия необходима для получения различных функциональных возможностей в области оптоэлектронных технологий. Структуры на основе благородных металлов, рассматриваемые в этой работе, находят применение в качестве вспомогательных компонентов в широком спектре приложений, среди которых биомедицина [1], фотоэлектрические приборы [2], тонкопленочные оптические устройства [3], датчики [4] и многие другие. Полный потенциал этих элементов еще не до конца изучен, особенно в контексте квантовых приложений [5-8]. Одним из наиболее важных эффектов, которого требуется достичь для работы ряда нанооптических приборов, является контроль распространения электромагнитной энергии в наноразмерных объемах.

Одним из очень многообещающих и уже проверенных подходов к управлению электромагнитными полями на наноразмерном уровне является использование благородных (плазмонных) металлов, таких как серебро и золото. Ключевым оптическим свойством этих материалов является их отрицательная диэлектрическая проницаемость в видимой и инфракрасной областях спектра. Наноструктурированные материалы с отрицательной диэлектрической проницаемостью обладают уникальным свойством поддерживать локализованные плазмонные резонансы (LPR), приводящим к локализации энергии в наномасштабе [9]. Способность плазмонных структур концентрировать свет за гранью классического дифракционного предела [10-12] имеет ключевое преимущество для управления взаимодействием света и вещества, открывая путь к реализации новых явлений квантовой электродинамики в классических резонаторах (CQED) и плазмонных CQED (PCQED). Основное принципиальное отличие между традиционной CQED и PCQED заключается в том, как сила взаимодействия контролируется структурируемой электромагнитной средой. Наиболее ярким примером, подчеркивающим различия между этими подходами, является изменение скорости спонтанного излучения излучателя относительно свободного пространства (эффект Парселла [13]). Флуктуации электромагнитных мод вакуума в свободном пространстве зависят исключительно от фундаментальных физических констант. Несмотря на это, плотность состояний можно изменить, формируя окружающую материальную среду, сохраняя при этом общую интегральную плотность состояний, доступных для излучения [14]. Фактически, электромагнитное структурирование приводит к спектральному и пространственному перераспределению вакуумных флуктуаций. Соответствующей величиной, описывающей данный эффект, является локальная плотность

состояний (LDOS) [9] — константа пропорциональности, модифицирующая гамильтонианы взаимодействия и, следовательно, результирующие скорости излучения. В стандартной квантовой электродинамике коэффициент Парселла (с точностью до константы) пропорционален отношению добротности резонатора к его объему ^ / V) и может зависеть от обеих величин. Традиционный CQED влияет на динамику взаимодействия посредством модификации добротности, в то время как PCQED действует в основном через модальные объемы. Это принципиальное различие связано с оптическими свойствами компонентов и геометрий материалов. В то время как диэлектрические резонаторы могут иметь достаточно большие добротности [15,16], их модальные объемы ограничены снизу классическим дифракционным пределом (примерно равным кубу длины волны в материале). С другой стороны, плазмонные структуры могут концентрировать энергию в небольших областях, в то время как их добротности обычно не превышают 100 [17], с небольшими исключениями [18]. С теоретической точки зрения, разница между CQED и PCQED заключается в способе квантования электромагнитного поля. Кроме того, стандартные определения факторов Парселла в PCQED могут привести к неоднозначным результатам и могут быть использованы только для предварительных качественных оценок [19]. Основными причинами, приводящими к противоречию с классической теорией Парселла, являются наличие дополнительных каналов безызлучательного распада и открытая природа резонаторов, усложняющая строгое определение модальных объемов. Эти аспекты были рассмотрены в серии недавних работ (например, [20,21]), а пересмотренный формализм фактора Парселла был разработан в [22]. Некоторые из вышеупомянутых аспектов наряду с другими воздействиями поглощающей и дисперсионной природы плазмонных материалов на динамику взаимодействия света и материи будут подробно обсуждаться ниже.

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

1 Также была введена концепция квантового метаматериала с элементарными ячейками на основе полупроводников с размерным квантованием (например, [340-342]).

увеличение скорости спонтанного излучения [24], что получило несколько экспериментальных подтверждений (например, [25,26]).

Оптомеханические силы и вспомогательные структуры. Электромагнитное излучение передает диэлектрическим структурам механический импульс через механизмы рассеяния, поглощения и передачи углового момента. Самосогласованное электромагнитное поле взаимодействует с плотностями наведенного поляризационного заряда в объекте, в результате чего возникают макроскопические Ньютоновские силы, действующие на объект. Поскольку импульс фотона пропорционален частоте излучения, оптические частоты являются наиболее подходящим спектральным диапазоном для манипуляции. Оптические силы достаточно интересны с фундаментальной и с прикладной точек зрения, и уже доказали свою значимость в различных областях науки от астрономии до биологии. Эксперименты по оптическому захвату объектов фокусированными лазерными лучами (оптический пинцет [27], Нобелевская премия по физике 2018 года) положили начало эре оптических манипуляций [28]. В настоящее время оптический захват является хорошо отработанным методом и одним из наиболее часто используемых инструментов в биофизических исследованиях [29]. Это также связано с появлением систем независимых ловушек, достигаемых с помощью так называемого голографического пинцета [30]. Основным физическим механизмом оптомеханики является взаимодействие между силой рассеяния (в направлении распространения света) и градиентной силой (вдоль пространственного градиента света). Для стабильного захвата в трехмерном пространстве аксиальная компонента градиента силы, притягивающей частицу к фокусу луча, должна доминировать над компонентой рассеяния (светового давления) силы, выталкивающей ее из этой области. Возможны также другие типы движения, такие как передача момента вращения двулучепреломляющей частице [31], с помощью конфигурированных пучков, переносящих импульс [32], и другие [33]. С точки зрения «статических» характеристик, пространственная локализация захваченного объекта ограничена как стохастическим броуновским движением, так и фокусным объемом объектива. Одним из многообещающих подходов к преодолению дифракционного предела улавливающего луча является использование вспомогательных плазмонных структур, которые, как известно, позволяют в ряде случаев преодолеть дифракционный предел [34] (почти такая же идея будет описана в следующем разделе в контексте усиления взаимодействия света с веществом). Достигнутый прогресс в работе над плазмонным пинцетом [35] демонстрирует возможности достижения более высоких точности и стабильности захвата, а также снижения мощности лазера, используемого для захвата, и позволяет достичь многих других характеристик - все благодаря возможности

формировать геометрию ближних электромагнитных полей. Диэлектрические наноразмерные частицы труднее улавливать (по сравнению с плазмонными) из-за их малых поляризуемостей. Объекты с более высоким показателем преломления обычно легче улавливать. Особым преимуществом использования плазмонных структур [36] является их способность поддерживать резонансные явления, даже если они имеют субволновой размер. Что касается достижимых параметров, теория предсказывает, что сферы с радиусом 5 нм могут быть захвачены лазером с 6,5 мВт и с точностью 5 нм с помощью вспомогательных золотых иголок [37], либо отверстиями в непрозрачных пленках, которые обеспечивают аналогичные характеристики [38]. В настоящее время захват плазмонных частиц с радиусом 20 нм может быть осуществлен относительно просто. С точки зрения точности измерения, суб-пиконьютоновая сила может быть измерена с помощью правильно сконструированных оптических систем [39]. Стоит отметить дополнительные возможности различных типов оптических пинцетов, которые используют следующие вспомогательные эффекты: передачу силы от поверхностных плазмонных волн частицам [40], использование нелинейной генерации с оптически захваченными наностержнями [41], совмещение оптического улавливания с микрофотолюминесценцией [42], темнопольную спектроскопию захваченных обездвиженных частиц [43], плазмонные голограммы для пинцетов [44], построение трехмерных структур с двумя кольцевыми лазерными интерференционными лучами [45], силовую микроскопию с помощью кремниевых наностержней [46], поверхностное сканирование с помощью голографически захваченных зондов [46], исследование квазикристаллов с условно запрещенными кристаллографическими пространственными точечными группами [47], захватывание вблизи поверхностей [48] и др. В настоящее время возможны оптические манипуляции с частицами на наномасштабе, с нанометровой точностью и измерениями силы в области пиконьютонов [49].

Радиофизика и основы эмуляции для изучения сложных волновых явлений. Радиофизика— это устоявшаяся область с долгой историей, развивающейся с начала 20-го века. Большинство электромагнитных явлений может быть описано с точки зрения уравнений Максвелла. Возможность их решения либо аналитически (ограниченный набор задач имеет полное аналитическое решение), либо численно (множество методов было разработано рядом исследовательских групп и коммерческих компаний) является ключом к пониманию и разработке гигагерцовых электромагнитных приложений. Тщательно разработанные теоретические, численные и экспериментальные платформы привели к появлению множества коммерчески доступных устройств, таких как антенны, линии электропередачи и сложные электронные устройства [50]. Помимо этих разработок, строго ориентированных на практическое применение,

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

Метаматериалы - новый подход к контролю электромагнитного рассеяния. Управление рассеянием является важной задачей в области прикладного электромагнетизма. В антенной технике были разработаны и использованы различные подходы, в том числе для увеличения сечения рассеяния и направленности, являющихся наиболее востребованными характеристиками [50]. Хотя усиление сечения рассеяния имеет множество применений, обратная задача, а именно -его сокращение, также представляет потенциальный интерес. Концепция «маскирующего устройства» была введена около десяти лет назад [52,53], и уже привлекла к себе значительное внимание из-за потребности в достижении невидимости для радиолокационных волн [54] и видимого света [55-57]. Общая процедура создания «маскирующего устройства» основывается либо на концепциях трансформационной оптики [52], либо на конформных трансформациях [53]. Для обеспечения правильного действия маскирующих объектов необходимо использовать материалы с анизотропными, пространственно изменяющимися электрической и магнитной восприимчивостями. В то время как материалы, доступные в природе, не могут удовлетворить этим требованиям, искусственно созданные периодические композитные структуры с нетривиально структурированными элементарными субволновыми ячейками могут демонстрировать необычные свойства. Например, композитные структуры с отрицательным показателем преломления (в, ц <0) были созданы с использованием массивов кольцевых дипольных резонаторов [58-60]. Кроме того, плазмонные материалы (в <0) могут эмулироваться средой из проводов [61] с учетом сильных эффектов пространственной дисперсии [62]. Кроме того, плазмонный тип маскировки был продемонстрирован в радиоволновом диапазоне [63]. Среди многих других применений (в видимом диапазоне) могут быть кратко упомянуты визуализация с суперразрешением [64], резонаторы на основе метаматериалов [65], компактные резонаторы [66] и многие другие.

Свойства антенн и рассеивателей сильно подвержены влиянию электромагнитной среды, в которую они помещены. На языке радиотехники данные эффекты могут описаны с помощью функций Грина, в то время как в оптике более принятой терминологией является плотность состояний и связанный с этим эффект Парселла [13]. Воздействие окружающих структур довольно часто наблюдаются в радиоволновом диапазоне, где различные типы наведенных помех и другие паразитные эффекты должны быть приняты во внимание для разработки правильного дизайна устройства. Например, характеристики направленности антенн основаны на тщательном

проектировании взаимодействия между источником излучения и набором вспомогательных резонансных рассеивателей. Применение метаматериалов для контроля свойств излучения антенн также может быть очень полезным и показать немало преимуществ по сравнению с традиционными подходами к дизайну, например, с точки зрения направленности [67,68]. Значительные усилия были посвящены исследованию факторов Парселла и усилению исчезающих полей внутри структурированных сред из проводов [69,70].

Объекты исследования

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

Список литературы диссертационного исследования доктор наук Гинзбург Павел Борисович, 2020 год

o -10 -

-12 -

Nanospiral, a = 1.03-3^ Nanospiral, a = 1.03-3^, no QP Nanospiral, a = 5/2n Nanorod, d = 200 nm Nanorod, d = 12 nm

1 2 3 4 5 6

Photon energy (eV)

Figure 3 | Resonant nonlinear response. (a) Schematic of the nanospiral indicating the angle a, defining the positions of the nanospiral resonances. (b) Spectra of the extinction cross-sections of the nanospirals with a = 5/2p and a = 1.03 ■ 3p normalized to their geometrical area (s = 70 nm and w = 12 nm). (c) Nonlinear scattering spectra of the nanospirals with a = 1.03 ■ 3p (the resonance at the fundamental frequency) and a = 5/2p and the nanorods with d = 200 nm and d = 12nm.

b

a

c

0

nanostructures of 10 pm V " 1 (ref. 20) and 3.2 pmV "1 (ref. 23), taking into account the local field enhancements and the surface areas. The SHG enhancement is robust with respect to geometrical scaling under the resonant excitation of the fundamental nanospiral mode which size-dependent spectral position can be matched to the fundamental wavelength by varying a (Fig. 4a). The SHG intensity from the twice larger nanospiral increases by the factor of ~ 1.4 in accordance with a similar increase of the surface area. In the nonresonant case of the nanorods much smaller than the wavelength and without nonlocality, the SHG intensity scales with size as / d4 (ref. 24), which is in excellent agreement with 5 orders of magnitude difference in the SHG intensities from 12 and 200 nm nanorods (dashed lines in Fig. 2). Such behaviour is the result of much smaller dipolar (which is retardation-related) and quadrupolar moments excited at the SH frequency for a smaller nanorod.

To demonstrate robustness of the effect for a wide range of experimentally relevant scenarios, the three-dimensional spirals of different thicknesses under various illumination conditions were simulated (Fig. 4b-d). The high-Q factor nanospiral resonance, a key to the observed high nonlinear susceptibilities, was confirmed for the nanospiral thickness h ranging from infinity (two-dimensional case considered above) down to the deep-subwavelength thicknesses (Fig. 4c). With the decrease of the thickness of the spiral, the fundamental resonance experiences a red shift. However, the nature of the resonance, as can be seen

in the field maps in Fig. 4b, remains the same with the similar field distributions and field enhancement values. Thus, by adjusting the angle a, the resonance position can be kept at the same fundamental wavelength for spirals of different finite thicknesses, so that the efficiency of the nonlinear processes is similar. Furthermore, this behaviour remains the same when the light is obliquely incident at a nanospiral (Fig. 4b,d).

Discussion

The interplay between nonlocality and nonlinearity on metal nanostructures was investigated. It was shown that the quantum pressure, the manifestation of collective many body dynamics of electron plasma in metals, together with structure topology plays the key role in the process of nonlinear harmonic generation. The appearance of high harmonics (up to 6th) and broadband white light generation from spiral shaped nanostructures is the result of the interplay between local geometry and fractional harmonic generation by the nonlocal quantum pressure term. Fractional harmonics were shown to mix efficiently with natural integer ones via nonlinear interaction processes to give rise to strong supercontinuum generation. The interplay between integer and fractional harmonics is unique for plasmonic systems, favouring them over existent nonlinear materials particularly on the nanoscale.

Coupled nonlinearities, macroscopic and microscopic effects (nonlinear mesoscopic phenomena), being one of the hardest

I El /1 En|

20

10

1 2 3 4 5

Photon energy (eV)

2D, inf. h = 500 nm h = 50 nm h = 20 nm

2.5x10

0.8 1.0 Photon energy (eV)

0.6 0.8 1.0 Photon energy (eV)

Figure 4 | Dependence of the resonant properties on nanospiral geometric parameters. (a) Nonlinear scattering spectra of the nanospirals resonant at the fundamental frequency for the designs with s = 70 nm, w = 12 nm, and a = 1.03 • 3p (orange line) and s = 140 nm, w = 24 nm, and a = 2.69p (black line). (b) Electric field distributions |E|/|E01 at the fundamental nanospiral resonance for (top) different nanospiral thicknesses at normal incidence (the spectra are shown in c) and (bottom) different incident angles for nanospiral thickness h = 20 nm (the spectra are shown in d). The field maps were taken across the nanospiral midplane marked in the insets in c and d by a dashed red curve. (c) Dependence of the nanospiral extinction spectra on the thickness at normal incidence. (d) Dependence of the extinction spectra of a 20-nm-thick nanospiral on the angle of incidence. The other nanospiral parameters in b-d are s = 70 nm, w = 12 nm, and a = 1.03 • 3p.

0

c

tasks for analytical and numerical solutions, have been solved here for the first time. The implemented semi-phenomenological method enables addressing basic collective and, as a result, nonlocal effects of the electron plasma and is applicable in a range of validity of the hydrodynamic model. The latter neglects a number of quantum and classical phenomena, such as the electron spill out, tunnelling across small gaps, temperature effects, ultra-fast none-quilibrium dynamics and few others. Nevertheless, the majority of the effects, beyond the scope of the model, provide higher-order corrections. As a result, the hydrodynamic model is known to provide qualitative results in a good agreement with the majority of experimental observations. Its integrability with large scale, geometry-invariant structures makes the approach a universal tool of nonlinear analysis. From the application stand point, the developed approach provides a guideline for designing nanoscale nonlinear devices important in modern photonic technologies.

Methods

Modelling of the transient nonlinear optical response. The nonlinear optical response of the plasmonic nanostructures was studied using time-domain finite element method simulation of an electromagnetic problem defined by a set of Maxwell's equations coupled in a self-consistent way to additional partial differential equations implementing the hydrodynamic model in a framework of Comsol Multiphysics software. The Maxwell's equations were expressed in terms of a vector potential A

Vxm-1 (VxA) + moS@At + mo@t(eo@At -P) = 0; (5)

while the hydrodynamic description of nonlinear transient response of plasmonic nanostructures given by equations 2-4 was introduced through coefficient form (equations (2) and (3)) and general form (equation (4)) of partial differential equations. The material for nanostructures was chosen to be gold. The constants for gold permittivity were taken to be no = 5.98 x 1028m - 3, g = 1.075 x 1014s - 1, Op = 13.8 x 1015s - 1, relying on available and widely used tabulated data. The simulation domain size was set to be 6 x 6 mm2 to ensure that the outer domain boundaries had no effect on the simulation results, which was checked. The simulation time span T = 7t and offset fo = 3t were chosen so that both pumping and scattered light pulses (the latter containing higher harmonics) entirely propagated across the simulation domain, ensuring the complete and time-interval-independent modelling of nonlinear effects.

9. Bloembergen, N., Chang, R., Jha, S. & Lee, C. Optical second-harmonic generation in reflection from media with inversion symmetry. Phys. Rev. 174, 813-822 (1968).

10. Sipe, J., So, V., Fukui, M. & Stegeman, G. Analysis of second-harmonic generation at metal surfaces. Phys. Rev. B 21, 4389-4402 (1980).

11. Ginzburg, P., Hayat, A., Berkovitch, N. & Orenstein, M. Nonlocal ponderomotive nonlinearity in plasmonics. Opt. Lett. 35, 1551-1553 (2010).

12. Scalora, M. et al. Second- and third-harmonic generation in metal-based structures. Phys. Rev. A 82, 043828 (2010).

13. Sundararaman, R., Narang, P., Jermyn, A. S., Goddard, III W. A. & Atwater, H. A. Theoretical predictions for hot-carrier generation from surface plasmon decay. Nat. Commun. 5, 5788 (2014).

14. Ginzburg, P., Krasavin, A. V., Wurtz, G. A. & Zayats, A. V. Non-perturbative hydrodynamic model for multiple harmonics generation in metallic nanostructures. ACS Photon. 2, 8-13 (2014).

15. Higuchi, T., Stockman, M. I. & Hommelhoff, P. Strong-field perspective on high-harmonic radiation from bulk solids. Phys. Rev. Lett. 113, 213901 (2014).

16. Ginzburg, P. & Zayats, A. V. Localized surface plasmon resonances in spatially dispersive nano-objects: phenomenological treatise. ACS Nano 7, 4334-4342 (2013).

17. Mortensen, N. A., Raza, S., Wubs, M., Sondergaard, T. & Bozhevolnyi, S. I. A generalized non-local optical response theory for plasmonic nanostructures. Nat. Commun. 5, 3809 (2014).

18. Raza, S., Bozhevolnyi, S. I., Wubs, M. & Mortensen, N. A. Nonlocal optical response in metallic nanostructures. J. Phys. Condens. Matter 27, 183204 (2015).

19. Plech, A., Kotaidis, V., Lorenc, M. & Boneberg, J. Femtosecond laser near-field ablation from gold nanoparticles. Nat. Phys. 2, 44-47 (2005).

20. Davidson, II R. B. et al. Eflcient forward second-harmonic generation from planar archimedean nanospirals. Nanophotonics 4, 108-113 (2015).

21. Miyazaki, H. T. & Kurokawa, Y. Squeezing visible lightwaves into

a 3-nm-thick and 55-nm-long plasmon cavity. Phys. Rev. Lett. 96, 097401 (2006).

22. Nikogosyan, D. N. Nonlinear Optical Crystals: A Complete Survey (Springer Science + Business Media, Inc., 2005).

23. Zhang, Y., Grady, N. K., Ayala-Orozco, C. & Halas, N. J. Three-dimensional nanostructures as highly efficient generators of second harmonic light. Nano Lett. 11, 5519-5523 (2011).

24. Dadap, J. I. Optical second-harmonic scattering from cylindrical particles. Phys. Rev. B 78, 205322 (2008).

Modelling of resonant properties of plasmonic nanostructures. The spectral response of the plasmonic nanostructures was simulated using frequency-domain finite element method in scattered-field formulation. The size of the simulation domain and the linear optical parameters of gold were set to be consistent with the time-domain simulations. The nanostructures were illuminated with a plane wave with parametrically changed frequency, while perfectly matched layers were implemented around the simulation domain to ensure the absence of reflection of the scattered waves from the outer boundaries. The nanostructure's extinction spectrum was calculated through the sum of the scattering and absorption cross-sections. The latter two parameters were calculated via integrating the incoming total power flow (for the absorption) and outcoming scattered power flow (for the spattering) over a nanoscale cylindrical region around the nanostructure and normalizing the obtained value to the power flow incident on the nanostructure's geometrical cross-section.

Acknowledgements

This work was supported, in part, by EPSRC (UK) and US Army Research Office (W911NF-12-1-0533). A.V.Z. acknowledges support from the Royal Society and the Wolfson Foundation. G.A.W. acknowledges support from the EC FP7 project 304179 (Marie Curie Actions). P.G. acknowledges Newton International Fellowship follow up programme (Royal Society). All data supporting this research are provided in full in the results section.

Author contributions

A.V.K. performed numerical simulations and analysed the results, A.V.K. and P.G. developed the theoretical framework, numerical model and data analysis methods. G.A.W. and A.V.Z. guided the work and analysis. All authors discussed the results and wrote the manuscript.

References

1. Barnes, W. L., Dereux, A. & Ebbesen, T. W. Surface plasmon subwavelength optics. Nature 424, 824-830 (2003).

2. Boltasseva, A. & Atwater, H. A. Low-loss plasmonic metamaterials. Science 331, 290-291 (2011).

3. Gramotnev, D. K. & Bozhevolnyi, S. I. Plasmonics beyond the diffraction limit. Nat. Photon. 4, 83-91 (2010).

4. Cai, W., Vasudev, A. P. & Brongersma, M. L. Electrically controlled nonlinear generation of light with plasmonics. Science 333, 1720-1723 (2011).

5. Park, I.-Y. et al. Plasmonic generation of ultrashort extreme-ultraviolet light pulses. Nat. Photon. 5, 677-681 (2011).

6. Kauranen, M. & Zayats, A. V. Nonlinear plasmonics. Nat. Photon. 6, 737-748 (2012).

7. Berkovitch, N., Ginzburg, P. & Orenstein, M. Nano-plasmonic antennas in the near infrared regime. J. Phys. Condens. Matter 24, 073202 (2012).

8. Ginzburg, P. et al. Nonlinearly coupled localized plasmon resonances: Resonant second-harmonic generation. Phys. Rev. B 86, 085422 (2012).

Additional information

Competing financial interests: The authors declare no competing financial interests.

Reprints and permission information is available online at http://npg.nature.com/ reprintsandpermissions/

How to cite this article: Krasavin, A. V. et al. Nonlocality-driven supercontinuum white light generation in plasmonic nanostructures. Nat. Commun. 7:11497 doi: 10.1038/ncomms11497 (2016).

This work is licensed under a Creative Commons Attribution 4.0 International License. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in the credit line; if the material is not included under the Creative Commons license, users will need to obtain permission from the license holder to reproduce the material. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/

SCIENTIFIC REPORTS

OPEN Optical forces in nanorod metamaterial

Received: 24 July 2015 Accepted: 29 September 2015 Published: 30 October 2015

Andrey A. Bogdanov1'2'3, Alexander S. Shalin1 & Pavel Ginzburg4

Optomechanical manipulation of micro and nano-scale objects with laser beams finds use in a large span of multidisciplinary applications. Auxiliary nanostructuring could substantially improve performances of classical optical tweezers by means of spatial localization of objects and intensity required for trapping. Here we investigate a three-dimensional nanorod metamaterial platform, serving as an auxiliary tool for the optical manipulation, able to support and control near-field interactions and generate both steep and flat optical potential profiles. It was shown that the 'topological transition' from the elliptic to hyperbolic dispersion regime of the metamaterial, usually having a significant impact on various light-matter interaction processes, does not strongly affect the distribution of optical forces in the metamaterial. This effect is explained by the predominant near-fields contributions of the nanostructure to optomechanical interactions. Semi-analytical model, approximating the finite size nanoparticle by a point dipole and neglecting the mutual re-scattering between the particle and nanorod array, was found to be in a good agreement with full-wave numerical simulation. In-plane (perpendicular to the rods) trapping regime, saddle equilibrium points and optical puling forces (directed along the rods towards the light source), acting on a particle situated inside or at the nearby the metamaterial, were found.

The ability to control mechanical motion of micro- and nano-scale particles with focused laser beams is an essential tool, being a paramount for a wide range of applications, related to bio-physics, micro-fluidics, optomechanical devices and more1-4. Being first proposed and demonstrated by A. Ashkin5, the classical optical tweezers are nowadays a rapidly developing area of fundamental and applied research.

One of the promising and already conceptually proven approaches for improving performances of the optical manipulation schemes is to employ various auxiliary nanostructures, especially plasmonic ones1,6. The key idea of the plasmonic tweezers is to utilize strong light-matter interaction between nano-structured metals and focused laser beams. Noble metals, having a negative permittivity in the optical and infrared spectral ranges, support localized plasmon resonances enabling enhancement and control of near-fields at their vicinity7,8. In particular, the creation of strong intensity gradients is beneficial for obtaining substantial optical forces, which is important, for example, for achieving molecular manipula-tion9. Resonant amplification of optical forces, exceeding the amplification in plasmonic structures, can be achieved in all-dielectric planar metamaterials due to Fano resonances10.

Plasmonic nanostructures with subwavelength light concentration could be employed for obtaining new optomechanical effects, i.e. accelerating nanoparticles in an arbitrary direction (in relation to the light propagation direction)11, or for creating nano-modulators of plasmonic signals12.

Arrays of antennas and their integrations in photonic circuitry1,13, employed as auxiliary tools for optical trapping, were shown to outperform classical schemes (focused lasers in homogeneous media, e.g. liquid solutions) both in terms of spatial localization and optical power required per trapped particle. Antenna arrays were further extended for multifunctional platforms, enabling trapping, stacking, and sorting14. However, isolated plasmonic structures create limited number of hot spots (local enhancement of intensity) and are usually restricted to two-dimensional geometries. These constrains set significant

1ITMO University, St. Petersburg, 197101, Russia. 2Ioffe Institute, St. Petersburg, 194021, Russia. 3Peter the Great

St. Petersburg Polytechnic University, St. Petersburg, 195251, Russia. 4School of Electrical Engineering, Tel Aviv University, Tel Aviv, 69978, Israel. Correspondence and requests for materials should be addressed to A.A.B. (email: [e-mail])

Figure 1. (a) Schematics of the nanorod metamaterial with a spherical nanoparticle inside. Radius of the nanoparticle is R = 5 nm. Height and radius of the nanorods are H = 350 nm and r = 15 nm, respectively, the period is h = 60 nm. (b) Frequency dependence of the effective tensor components (real parts) of the nanorod metamaterial. Dashed line shows the transition between elliptic and hyperbolic dispersion regimes of the metamaterial - epsilon-near-zero (ENZ) point. Insets show characteristic shapes of iso-frequency surfaces, corresponding to the dispersion regimes.

limitations on the flexibility of optical manipulation by reducing potential degrees of freedom, available for optomechanical control. On the other hand, three-dimensional artificially created nanostructures or metamaterials could provide additional benefits and flexibility by configuring near-field interactions in large volumes15.

Hyperbolic metamaterials16 are one class of artificially created electromagnetic structures, capable to enhance efficiencies of various light-matter interaction processes owning to the unusual dispersion regimes of eigenmodes, supported by the structure - namely its hyperbolic dispersion. Among various possible designs of this type of metamaterials it is worth mentioning composites made of vertically aligned nanorods17,18, periodic metal dielectric layers19 and semiconductor quantum structures20,21. While the far-field interactions of waves with hyperbolic composites were proven to be well characterized in terms of the effective medium approximation, this description could be questioned if near-field mediated processes are involved22.

The general criterion for estimating impact of near-field contributions to an interaction is based on comparison of k-vector spectra with reciprocal vector of the metamaterial lattice. For example, scattering from objects within hyperbolic metamaterials involves consideration of the near-field effects23.

Analysis of near- and far-field contributions to optical force, acting on objects embedded in the nanorod metamaterial is the central topic of the manuscript. In particular, optical forces, acting on nano-sized spherical particle embedded inside the metamaterial assembly, are investigated both numerically and by using a semi-analytical approach, considering the finite size nanoparticle as a point dipole and neglecting re-scattering between the particle and nanorod array. The impact of the finite structure of the metamaterial unit cell and the relative arrangement of the particle in respect to it was analyzed as a function of the system's geometry and frequency of incident illumination. A semi-analytical approach based on dipole near-field interaction is developed and shown to be in a good agreement with results of the full-wave numerical analysis. The interplay between near- and far-field effects in the context of effective medium approximation is discussed.

Nanorod Metamaterial: far-field characteristics

The geometry under investigation is schematically represented in Fig. 1a - it shows an array of vertically aligned gold nanorods and a gold nanoparticle placed inside it. Material parameters of the constitutive elements were taken from widely used sources24. The parameters of the structure are indicated in the figure caption. While the nanorods in this model are situated in vacuum, substrate effects and host material filling the space between the rods could be taken into account straightforwardly. Similar structures have already found use in various multidisciplinary applications, among them bio-sensing25, enhancement of nonlinearities26, acoustic waves detection27, thin optical elements18 and others. The key properties of this auxiliary nanostructure leading to enhanced performance are large surface area and unusual collective optical response of the system, enabling control over both far- and near-field interactions. Hence, investigation of optical forces, mediated by nanorod metamaterials, has a profound potential interest.

Far-field interactions between electromagnetic waves and metamaterials, under certain circumstances, can be described within the effective medium approximation. The main idea of this homogenization procedure is to average the electromagnetic field over a unit cell of a structure. Therefore, the field inside the structure is assumed to be uniform. In the context of optical forces, as it will appear in the next section, the non-uniformities play a major role and, in fact, predefine the spatial structure of optical potentials.

Recently, a phenomenological approach taking into account the finite size of the metamaterial unit cell was proposed28. Inclusion of a depolarization volume around optically manipulated particles enabled investigations of far-field contributions to optical forces. However, near-field interactions, being strongly dependent on a specific metamaterial design, were not included explicitly.

One of the key properties of hyperbolic metamaterials, making them attractive for electromagnetic applications, is their unique ability to support an unusual regime of dispersion caused by having permittivity tensor components of opposite signs < 0). An immediate implication of this hyperbolic dispersion regime is the high density of photonic states, available for both emission and scattering29,30.

The effective permittivity tensor of nanorod metamaterial is given by:

s± 0 0

s = 0 s± 0

0 0 s

Here s± and £|| are effective permittivities perpendicular and along the wires, respectively.

Dispersion of the tensor components for the structure under consideration was calculated with the approach developed in31. The transition between elliptic and hyperbolic dispersion regimes occurs at the wavelength around 523 nm (see Fig. 1b). The transition point is called epsilon-near-zero (ENZ) regime, as the real part of the permittivity along the rods is vanishing, if the spatial dispersion effects are ignored. The high density of photonic states as well as the strong scattering emerges in the hyperbolic and ENZ regimes. The wavelength of the external illumination, exploited for optical manipulation in the subsequent investigations, is chosen around this ENZ point in order to distinguish between various dispersion regimes and their impact on optical forces.

Optical force distribution

Numerical model. The distribution of the optical forces, acting on the gold nanoparticle placed inside the wire medium is analyzed hereafter. The hyperbolic, ENZ, and elliptic dispersion regimes of the bulk metamaterial and their impacts on optical forces are compared and discussed.

In the first case, normal incidence scenario is considered - the illumination is chosen to be linearly polarized along the y-axis and it propagates along the z-axis (p = 0°, see Fig. 1a). Oblique incidence with p = 45° will be investigated hereafter. Full 3D numerical analysis, based on finite elements method32, is performed in order to calculate self-consistent electromagnetic fields in the system. Consequently, optical forces acting on the nanoparticle are calculated by integrating the Maxwell's stress tensor components over an imaginary spherical surface surrounding the nanoparticle.

The presence of a single nanoparticle breaks the inherent translation symmetry of the initial metamaterial geometry. In order to overcome the computation complexity of large systems modeling, Floquet periodical boundary conditions were imposed on finite size geometries. This type of model corresponds to a periodic system with variable unit cell, which consists of a square array of nanorods and the nan-oparticle. If the electromagnetic coupling between the particles in adjacent cells is minor, this type of analysis recovers the behavior of the infinite system with a single particle.

The numerical procedure is as follows: the number of rods in the unit cell is increased gradually and the convergence of a certain quantity (optical forces in our case) is checked. Recently, a similar approach was applied in studies of the Purcell effect in nanorod33 and wire34 metamaterials. Square unit cells containing 4, 9, and 16 nanorods were considered and the convergence of optical force values at different points of the metamaterial volume was checked. A unit cell of 4 nanorods (the smallest one) was shown to predict the behavior of an infinite array within the accuracy of several percent. All the subsequent results were obtained for this size of the unit cell. The direct consequence of this calculation is that (i) only nearest neighbor rods define the value of optical force and (ii) nanoparticles in different unit cells almostly do not interact with each other.

It should be noted, however, that the collective macroscopic behavior of the array is taken into account by imposing periodical Floquet boundary conditions.

Lateral force component. All the subsequent calculations were done for a particle of 10 nm in diameter. The optical force F, in the most general case, has three non-zero components (Fx, Fy, Fz). The lateral force F± = (Fx, Fy, 0) will be analyzed first. Values of optical forces are normalized to the intensity of the incident wave and volume of the particle in order to perform direct comparisons with other optical manipulation schemes.

The resulting normalized forces at the cut-plane z = 10 nm, calculated for wavelengths A = 450 and 600 nm, corresponding to the elliptic and hyperbolic dispersion regimes respectively, are shown in Fig. 2(a). Numerical simulation shows that spatial distribution of both electric field and optical force for different wavelengths of excitation (different dispersion regimes) are qualitatively similar and the quantitative difference is shown with different color bars on top panel in Fig. 2. The qualitative explanation of such a behavior is two fold: (i) small size of the sphere and nanorods (in comparison with wavelength) enables considering those structures as point dipoles; (ii) the same material of the sphere and nanorods providing similar frequency behaviour of their polarizability.

Figure 2. (a) Lateral optical force distribution at the cut-plane z = 10 nm [see Fig. 1(a)]. Magnitude of the

lateral force F = (Fx, Fy, 0) is shown with a color scale and the direction with arrows. (b) Distribution of

II v x y J

the electric field magnitude in the cut-plane z = 10 nm. Magnitude of the electric field is shown with the color scale. The arrow lines are electric field lines. Panel (c) shows the distribution of the force calculated with a semi-analytical dipolar approach [see Eq. (2)], i.e. where the perturbation of the field by the particle is neglected. Panel (d) shows the numerical simulation of the electric field magnitude distribution without the particle. Panel (a) should be compared with (c), while (b) with (d). Dark blue shells around nanorods on panels (a), (c) have the width equal to particle's radius. Optical forces are not calculated at those areas, as the nanoparticle's center cannot approach the nanorods that close. Upper and lower scales of the color bars correspond to the hyperbolic (A = 600 nm) and elliptic (A = 450 nm) dispersion regimes of the metamaterial, respectively. Electric field amplitude of the incident wave is 1 V/m. Electric field E of the incident wave is parallel to the y-axis. Wavevector of the incident wave k = (0, 0, k0).

Force maps at various cut-planes (with different z-coordinate) show qualitatively similar behavior too. The reasoning is as follows: at normal incidence, lateral optical force is mainly determined by first (gradient) term in Eq. (2). Within a good approximation |E|2 can be factorized as f(z)g(x, y), where f(z) is Fabry-Perot envelope function and g(x, y) is the in-plane field distribution. Therefore, the lateral force distribution is similar for different cut-planes while its absolute value is modulated by Fabry-Perot envelope function.

It should be noted, that the values and directions of forces are attributed to the geometrical center of the particle, hence certain regions (dark blue shells around nanowires with thickness equal to the radius of nanoparticle) on Fig. 2(a,c) are blank, as this center cannot approach the boundaries of the rods.

It can be seen from Fig. 2(a) that the force distribution has saddle points at the center of the unit cell and at its edges between the rods. These places correspond to the the saddle points of electromagnetic field magnitude distribution [Fig. 2(d)], and, consequently, to the unstable equilibrium positions of the nanoparticle. Some peculiarities in the optical force distribution appear on cut-planes near the edges of the nanorods (z = 350 nm and z = 0 nm), but they are attributed to the longitudinal (z-component) force component and will be discussed further.

The similarity of the spatial distribution of the forces at hyperbolic and elliptic dispersion regimes results from the dominating near-field coupling between the nanorods and the particle. Figure 2(b)

shows the magnitude distribution of total electric field |E| = (( + Ey + EZ) ^ ,while the arrows show its direction. One can see that the field map is formed by electrical dipoles induced on the rods and the particle by the incident wave. Orientations of the dipoles nearly coincide with the polarization of the incident wave. Minor deviations from the above description are related to the higher multipole contribution and the interaction between the particle and the nanorods.

Lateral force distribution analysis can be provided with the following semi-analytical approach. First, the total electric field distribution in the nanorod array under external incident wave without the particle is calculated numerically with the periodic boundary conditions applied. Results of the simulation are shown in Fig. 2(d). The knowledge of the spatial distribution of the electric field magnitude enables calculation of optical forces with two assumptions: (i) the nanoparticle is represented by a structureless point electric dipole with a moment / (ii) the dipole is assumed to act as a small perturbation to the fields of the standalone metamaterial. This means, that only collective scattering properties of the nano-rod array were taken into account, while the mutual re-scattering of the field between the particle and nanorods was neglected. Comparison between Fig. 2(b,d) verifies this approximation - both the structure and values of the field magnitude are similar.

The time averaged optical force acting on the point dipole is given by7:

<F> = E |2 + if Re[E x HI,

42 c V £0 (2)

where a = a' + ia'' is the complex particle's polarizability. The polarizability of the spherical particle is given by7:

a = 4nsnR %u - 1

+ 2 (3)

The resulting optical force map, calculated using the dipolar approximation [Eq. (2)] appears on Fig. 2(c). The arrows show the direction of the force at corresponding points. The color pattern corresponds to the absolute value of the force. The remarkable similarities between Maxwell's stress tensor calculations [Fig. 2(a)] and the approximate analytical model [Fig. 2(c)] suggest the validity of the dipolar model and highlights the impact of near-fields on the optical force. It should be noted, however, that overall values of optical forces, calculated within those approaches, have about 20% difference, which is related to the finite size of the particle and the mutual field re-scattering between the particle and nanorods.

Distribution of optical force in the case of oblique incidence (p = 45°, see Fig. 1a) obtained with full numerical simulations and semi-analytical approach is shown in Fig. 3. One can see that in contrast to the case of normal incidence electric field distribution around the nanorods is asymmetric [Fig. 3(b,d)]. It results in asymmetry of optical force distribution. The asymmetry can be explained by excitation of non-dipole modes in the nanorods at oblique illumination35.

Vertical force component. The distribution of the lateral optical force component (perpendicular to the nanorods) was analyzed in the previous section. Longitudinal force component Fz (parallel to the nanorods) is analyzed here.

As it was already mentioned, the homogenization procedure averages the near-fields over the unit cell, hence, it is inapplicable for estimation of gradient optical force in the lateral plane. Nevertheless, the field distribution along the z-axis can be roughly estimated considering the slab of the nanorod metamaterial as a Fabry-Perot resonator in z-direction33. Therefore, it is reasonable to expect a standing wave in the slab (along the z-axis) and maxima of the electric field resulting in in-plain trapping of the particle. The results of the numerical calculation suggest the validity of this hypothesis. Estimation of the field maxima position deeply inside the slab can be provided by the effective medium approximation36,37 but a more detailed analysis, that takes into account boundary effects, demands numerical simulation.

The profiles of the total electric field magnitude along the line parallel to the rods and passing through the point x = 30 nm and y = 5 nm [see inset in Fig. 4(a)] calculated without nanoparticle for the elliptic (A = 450 nm), ENZ (A = 523 nm), and hyperbolic (A = 600 nm) regimes are shown in Fig. 4(a). The insets show distribution of the total electric field magnitude in xz-plane passing through y = 5 nm. One can see that electric field distribution along the z-axis strongly depends on the wavelength of the incident wave. In the hyperbolic regime (A = 600 nm), three distinct field maxima are observed - one inside the slab and two in the vicinity of its boundaries. In the elliptic (A = 450 nm) and ENZ (A = 523 nm) regimes, field decays inside the metamaterial and weak oscillations do not possess sharp field maxima. Additional contribution to those differences (apart from the interplay of dispersion regime and geometry, namely Fabry-Perot conditions) comes from a strong wavelength dependence of losses in gold24:

X = 600nm 0 1.5 3.0 4.5 6.0

X = 600 nm 0.0 0.6 1.2 1.8 2.4 3.0

Figure 3. (a) Lateral optical force distribution at the cut-plane z = 10 nm see Fig. 1(a)]. Magnitude of the lateral force k = (0, k0/42, k0/42) is shown with a color scale and the direction with arrows. (b) Distribution of the electric field magnitude in the cut-plane z = 10 nm. Magnitude of the electric field is shown with the color scale. The arrow lines are electric field lines. Panel (c) shows the distribution of the force calculated with a semi-analytical dipolar approach [see Eq. (2)], i.e. where the perturbation of the field by the particle is neglected. Panel (d) shows the numerical simulation of the electric field magnitude distribution without the particle. Panel (a) should be compared with (c), while (b) with (d). Dark blue shells around nanorods on panels (a), (c) have the width equal to particle's radius. Optical forces are not calculated at those areas, as the nanoparticle's center cannot approach the nanorods that close. Upper and lower scales of the color bars correspond to the hyperbolic (A = 600 nm) and elliptic (A = 450 nm) dispersion regimes of the metamaterial, respectively. Electric field amplitude of the incident wave is 1 V/m. Electric field E of the incident wave lies in the zy-plane. Wavevector of the incident wave k _ (0,k0/^> k0/^).

Re(£au) ' 0.3 for A = 450nm;

Im (^Au) v 6.2 for A = 600nm. (4)

Optical losses cause the reduction in quality factors of the modes, smearing out the sharp peaks, as could be seen in the case of elliptic dispersion.

Distributions of the z-component of the optical force along the nanorod for elliptic (A = 450 nm), ENZ (A = 523 nm) and hyperbolic (A = 600 nm) regimes are shown in Fig. 4(b). Positions of the stable trapping in transverse planes are marked with arrows (note, that the force derivative should be negative in order to obtain a stable equilibrium). Shaded areas on the figure show the regions within the metamaterial, where the optical force component Fz is directed towards the light source (for A = 600 nm). Optical pulling forces or optical attraction gained considerable attention over the last decade, as it provides additional flexible degree of freedom in optical manipulation38.

Figure 4. Distribution of: (a) Total electric field magnitude |E| and (b) longitudinal optical force component Fz along the line passing through the point with coordinates x = 30 nm and y = 5 nm (see the inset with geometrical arrangement) and parallel to the nanorods for the elliptic (A = 450 nm), ENZ (X = 523 nm), and hyperbolic (A = 600 nm) regimes of the metamaterial. The insets in panel (a) show the distribution of the electric field magnitude in the zx-cut-plane passing through y = 5 nm for the elliptic (A = 450 nm), ENZ (A = 523 nm), and hyperbolic (A = 600 nm) regimes of the metamaterial. Shaded areas in panel (b) show the region where a pulling force emerges. Electric field amplitude of the incident wave is 1 V/m.

-6.4 -5.6 -4.8 -4.0 -3.2 -2.4 -1.6 -0.8 0 0.8

Figure 5. Distribution of the normalized optical force acting on the nanoparticle, situated at the lateral plane above the metamaterial. Distance between the center of the nanoparticle and top faces of the nanorods is 10 nm. Arrows show the direction of the optical force. Color map shows the distribution of the optical force component parallel to the nanorods (Fz). Black solid lines show the the geometrical edges of the nanorods.

The hyperbolic regime supports three regions of optical attraction, while the elliptic and ENZ have only one, as could be seen in Fig. 4(b). This occurrence could be understood as follows: in the elliptic regime both weak gradient of the electric field magnitude [see Fig. 4(a)] and high material losses of the particle result in the domination of radiation pressure [second term in Eq. (2)] over the gradient force. The radiation pressure is co-directional with the Poynting vector of the incident radiation, so the optical attraction cannot be obtained in this case. Nevertheless, the first term of Eq. (2) overcomes the second one in the vicinity of the nanorod's edge where strong gradient of the electric field intensity is observed [see Fig. 4(a)]. In the hyperbolic regime, on the other hand, there are several regions where the optical force component is directed to the light source - that's the result of high quality factor Fabry-Perot modes and dominating real part of the particle's polarizability.

For the hyperbolic regime (A = 600 nm), optical potential in the vicinity of nanorod's edge is stronger than in elliptical and ENZ regimes. For example, optical traping potential of 26 meV for 5nm radius particle can be achieved with electric field intensity ~107 V/m that corresponds to focusing of 1 W beam to 3 /im spot in diameter.

As a separate case, the particle situated over the metamaterial slab will be considered next. This scenario describes the case where the metamaterial is used as a substrate for advanced optical manipulation. Results of numerical studies appear in Fig. 5, showing the distribution of the vertical optical force acting on the nanoparticle in the lateral plane of z = 10 nm above the nanorods. It could be seen, that

the maximal attraction force on the particle emerges in the vicinity of nanorods edges (the sample is illuminated from below - see Fig. 5). Arrows indicate the direction of the optical force. Color pattern shows the distribution of the optical force component parallel to the nanorods (Fz). Solid white lines correspond to Fz = 0. Remarkable behaviour of forces above the metamaterial substrate could suggest the later as an auxiliary nanostructure or metasurface, providing additional flexibility in optical manipulation.

Conclusion

In this work, comprehensive analysis of the optical forces acting on a metal nanoparticle placed inside or in the vicinity of three-dimensional nanorod metamaterial slab was performed. Numerical simulations of finite size square unit cells with periodical Floquet boundary conditions enable to take into account all collective effects in the metamaterial and estimate optical forces on small particles. Unit cells containing 4, 9, and 16 nanorods were analyzed and the convergence of the optical forces for different positions of the particle was checked. It was shown that the smallest unit cell already reproduces the effect of optical forces on a particle, situated within the infinite metamaterial. Therefore, only four neighboring nanorods nearest to the particle make the dominant contribution to the optical forces. This statement has been confirmed with the developed semi-analytical model which neglects the particle's interior and the re-scattering effects between the particle and nanorods. Furthermore, it was shown that the 'topological transition' from the elliptic to hyperbolic dispersion regime of the metamaterial, usually having an impact on various light-matter interaction processes, is less important for optical forces.

In-plane optical trapping and optical pulling forces were observed. The comprehensive numerical modeling enables estimation of optical forces values, normalized to incident power and particle's volume. Values as high as 2.3 x 103 pN/W/nm for both lateral and optical pulling forces were predicted. Those results overcome other reported values39,40.

The remarkable structure of predicted optomechanical interactions (in particular pulling forces), mediated by the metamaterial, makes the later to be a promising platform for large span of multidiscipli-nary applications, involving demands for precise nanoscale mechanical manipulation, including trapping sorting, mixing and more.

References

1. Juan, M. L., Righini, M. & Quidant, R. Plasmon nano-optical tweezers. Nature Photonics 5, 349-356 (2011).

2. Dienerowitz, M., Mazilu, M. & Dholakia, K. Optical manipulation of nanoparticles: a review. Journal of Nanophotonics 2, 021875 (2008). URL http://link.aip.org/link/JNOACQ/v2/i1/p021875/s1&Agg=doi.

3. Marago, O. M., Jones, P. H., Gucciardi, P. G., Volpe, G. & Ferrari, A. C. Optical trapping and manipulation of nanostructures. Nature Nanotechnology 8, 807-819 (2013).

4. Grier, D. G. A revolution in optical manipulation. Nature 424, 810-816 (2003).

5. Ashkin, A. Acceleration and trapping of particles by radiation pressure. Physical Review Letters 24, 156 (1970).

6. Zhang, J., MacDonald, K. F. & Zheludev, N. I. Optical gecko toe: Optically controlled attractive near-field forces between plasmonic metamaterials and dielectric or metal surfaces. Physical Review B 85, 205123 (2012). URL http://link.aps.org/ doi/10.1103/PhysRevB.85.205123.

7. Novotny, L. & Hecht, B. Principles of nano-optics (Cambridge university press, 2012).

8. Berkovitch, N., Ginzburg, P. & Orenstein, M. Nano-plasmonic antennas in the near infrared regime. Journal of Physics: Condensed Matter 24, 073202 (2012).

9. Shoji, T. & Tsuboi, Y. Plasmonic optical tweezers toward molecular manipulation: Tailoring plasmonic nanostructure, light source, and resonant trapping. The Journal of Physical Chemistry Letters 5, 2957-2967 (2014).

10. Zhang, J., MacDonald, K. F. & Zheludev, N. I. Giant optical forces in planar dielectric photonic metamaterials. Optics Letters 39, 4883 (2014). URL http://www.ncbi.nlm.nih.gov/pubmed/25121899 https://www.osapublishing.org/ol/abstract.cfm?uri=ol-39-16-4883.

11. Shalin, A. S. & Sukhov, S. V. Plasmonic nanostructures as accelerators for nanoparticles: optical nanocannon. Plasmonics 8, 625-629 (2013).

12. Shalin, A. S., Ginzburg, P., Belov, P. A., Kivshar, Y. S. & Zayats, A. V. Nano-opto-mechanical effects in plasmonic waveguides. Laser & Photonics Reviews 8, 131-136 (2014).

13. Lin, P.-T., Chu, H.-Y., Lu, T.-W. & Lee, P.-T. Trapping particles using waveguide-coupled gold bowtie plasmonic tweezers. Lab on a Chip 14, 4647-4652 (2014).

14. Roxworthy, B. J. et al. Application of plasmonic bowtie nanoantenna arrays for optical trapping, stacking, and sorting, Nano Letters 12, 796-801 (2012).

15. Cai, W. & Shalaev, V. Optical metamaterials: fundamentals and applications (Springer, 2010).

16. Shekhar, P., Atkinson, J. & Jacob, Z. Hyperbolic metamaterials: fundamentals and applications. Nano Convergence 1, 1-17 (2014).

17. Atkinson, R. et al. Anisotropic optical properties of arrays of gold nanorods embedded in alumina. Physical Review B 73, 235402 (2006).

18. Ginzburg, P. et al. Manipulating polarization of light with ultrathin epsilon-near-zero metamaterials. Optics Express 21, 14907-14917 (2013).

19. Krishnamoorthy, H. N., Jacob, Z., Narimanov, E., Kretzschmar, I. & Menon, V. M. Topological transitions in metamaterials. Science 336, 205-209 (2012).

20. Hoffman, A. J. et al. Negative refraction in semiconductor metamaterials. Nature Materials 6, 946-950 (2007).

21. Ginzburg, P. & Orenstein, M. Nonmetallic left-handed material based on negative-positive anisotropy in low-dimensional quantum structures. Journal of Applied Physics 103, 083105 (2008).

22. Ginzburg, P. et al. Self-induced torque in hyperbolic metamaterials. Physical Review Letters 111, 036804 (2013).

23. Shalin, A. S. et al. Scattering suppression from arbitrary objects in spatially dispersive layered metamaterials. Phys. Rev. B 91, 125426 (2015). URL http://link.aps.org/doi/10.1103/PhysRevB.91.125426.

24. Johnson, P. B. & Christy, R. W. Optical Constants of the Noble Metals. Physical Review B 6, 4370 (1972). URL http://dx.doi. org/10.1103/PhysRevB.6.4370.

25. Kabashin, A. et al. Plasmonic nanorod metamaterials for biosensing. Nature Materials 8, 867-871 (2009).

26. Wurtz, G. A. et al. Designed ultrafast optical nonlinearity in a plasmonic nanorod metamaterial enhanced by nonlocality. Nature Nanotechnology 6, 107-111 (2011).

27. Yakovlev, V. V. et al. Ultrasensitive non-resonant detection of ultrasound with plasmonic metamaterials. Advanced Materials 25, 2351-2356 (2013).

28. Shalin, A., Sukhov, S., Bogdanov, A., Belov, P. & Ginzburg, P. Optical pulling forces in hyperbolic metamaterials. Phys. Rev. A (2015).

29. Jacob, Z. et al. Engineering photonic density of states using metamaterials. Applied Physics B 100, 215-218 (2010).

30. Iorsh, I. V., Poddubny, A. N., Ginzburg, P., Belov, P. A. & Kivshar, Y. S. Compton-like polariton scattering in hyperbolic metamaterials. Physical Review Letters 114, 185501 (2015).

31. Elser, J., Wangberg, R., Podolskiy, V. A. & Narimanov, E. E. Nanowire metamaterials with extreme optical anisotropy, Applied Physics Letters 89, 261102 (2006). URL http://scitation.aip.org/content/aip/journal/apl/89/26/10.1063/L2422893.

32. Pryor, R. W. Multiphysics Modeling Using COMSOL (Jones & Bartlett Learning, 2011).

33. Slobozhanyuk, A. P. et al. Purcell effect in hyperbolic metamaterial resonators, arXiv preprint arXiv:1504.06950 (2015).

34. Mirmoosa, M., Kosulnikov, S. Y. & Simovski, C. Topological phase transition in wire medium enables high purcell factor at infrared frequencies, arXiv preprint arXiv:1504.05936 (2015).

35. Velichko, E. A. & Nosich, A. I. Refractive-index sensitivities of hybrid surface-plasmon resonances for a core-shell circular silver nanotube sensor. Optics Letters 38, 4978 (2013).

36. Silveirinha, M. G. Nonlocal homogenization model for a periodic array of e-negative rods. Physical Review E 73, 046612 (2006). URL http://link.aps.org/doi/10.1103/PhysRevE.73.046612.

37. Tyshetskiy, Y. et al. Guided modes in a spatially dispersive wire medium slab. JOSA B 31, 1753-1760 (2014).

38. Chen, J., Ng, J., Lin, Z. & Chan, C. Optical pulling force. Nature Photonics 5, 531-534 (2011).

39. Yang, X., Liu, Y., Oulton, R. F., Yin, X. & Zhang, X. Optical forces in hybrid plasmonic waveguides. Nano Letters 11, 321-328 (2011).

40. Wang, K., Schonbrun, E. & Crozier, K. B. Propulsion of gold nanoparticles with surface plasmon polaritons: Evidence of enhanced optical force from near-field coupling between gold particle and gold film. Nano Letters 9, 2623-2629 (2009).

Acknowledgements

This work was partially supported by the Government of the Russian Federation (Grant 074-U01), by the Russian Foundation for Basic Research (No. 15-02-01344), by the Program on Fundamental Research in Nanotechnology and Nanomaterials of the Presidium of the Russian Academy of Sciences. The investigation of optical forces distributions has been supported by the Russian Science Foundation Grant No. 14-12-01227. A.B. thanks Russian Federation President support program of leading scientific schools (NSh-5062.2014.2) and RFBR (No. 14-02-01223).

Author Contributions

A.B. performed the numerical simulations. A.S. supervised the project. The paper was written by A.B. and P.G. All authors discussed the results and contributed to the manuscript.

Additional Information

Competing financial interests: The authors declare no competing financial interests.

How to cite this article: Bogdanov, A. A. et al. Optical forces in nanorod metamaterial. Sci. Rep. 5, 15846; doi: 10.1038/srep15846 (2015).

l/jjv 0 I This work is licensed under a Creative Commons Attribution 4.0 International License. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in the credit line; if the material is not included under the Creative Commons license, users will need to obtain permission from the license holder to reproduce the material. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/

Pfiotonics

pubs.acs.org/journal/apchd5

Nonperturbative Hydrodynamic Model for Multiple Harmonics Generation in Metallic Nanostructures

Pavel Ginzburg,*'^'§ Alexey V. Krasavm,^'§ Gregory A. Wurtz,^ and Anatoly V. Zayats^

^Department of Physics, King's College London, Strand, London WC2R 2LS, United Kingdom *ITMO University, St. Petersburg 197101, Russia

Letter

ft.

SHG fieldmap r = 100 nm Au rod

Hydrodynamic model

ABSTRACT: The electromagnetic response of a free-electron gas leads to the inherent nonlinear optical behavior of nanostructured plasmonic materials enabled through both strong local field enhancements and complex collective electron dynamics. Here, a time-domain implementation of the hydrodynamic model for conduction electrons in metals has been developed to enable nonperturbative studies of nonlinear coherent interactions between light and plasmonic nanostructures. The effects originating from the convective acceleration, the magnetic contribution of the Lorenz force, the quantum electron pressure, and the presence of the nanostructure's boundaries have been taken into account, leading to the appearance of both second and third harmonics. The proposed time-domain method enables obtaining a universal, self-consistent numerical solution, free from any approximations, allowing investigations of nonlinear optical interactions with arbitrary spatially and temporally shaped optical pulses, opening unique opportunities to approach description of realistic experimental scenarios.

KEYWORDS: plasmonics, nonlinear optics

Nonlinear optical interactions give rise to a variety of phenomena extensively used in numerous applications in lasers, classical and quantum optical information processing, bioimaging, and sensing. Consequently, tailoring, enhancing, and controlling nonlinear processes, which are inherently weak and require high light intensities, is a task of prime importance. One approach offering promising developments in this area employs nanostructured materials to tailor enhanced effective nonlinear susceptibilities through geometrical means by shaping local electromagnetic fields.

Plasmonic nanostructures deliver both true nanoscale light confinement and strong nonlinear response, as both rely on peculiarities of the response of conduction electrons to the electromagnetic field of incident light.1 Both second and third harmonics generation (SHG and THG), optical Kerr non-linearities, and the possibility to achieve solitonic waves have been demonstrated via the coupling to plasmonic resonances in metal films and metallic nanostructures.2-6 Typical theoretical approaches used in the description of the nonlinear response of plasmonic materials rely on a perturbative treatment of the material's polarization, either through a hydrodynamic model or phenomenologically, enabling the filtering of the nonlinear harmonic fields generated via a quasi-Fourier transformation.7

In this Letter, we develop a comprehensive and non-perturbative numerical model for the investigation of nonlinear interactions of light with plasmonic nanostructures. The complete time-domain analysis fully addresses the nonlinear dynamics of free electrons without any additional assumptions on the nature of the interaction, providing the opportunity to

explore the interplay between various nonlinear optical processes. The optical properties of metal particles are described with the help of a full hydrodynamic model,8 which enables accounting for both linear and nonlinear dynamics of the conduction electrons under visible and infrared light illumination, away from the spectral range of interband transitions. While the hydrodynamic equations straightforwardly reproduce the Drude model in the linear regime of interaction with electromagnetic fields, at higher intensities, convective acceleration, magnetic contributions from the Lorenz force, and quantum electron pressure terms, present in the hydrodynamic model, introduce strong nonlinear contributions to the optical response of the system. Timedomain studies give an opportunity to explore simultaneously both bulk and surface contributions to nonlinear generation processes, as well as the efficiency of sideband generation. All these effects have simultaneously been taken into account by coupling nonlinear hydrodynamic equations, describing the behavior of the electron plasma, with Maxwell's equations to model the response to electromagnetic fields. As a model example, we analyze the nonlinear scattering from a gold nanocylinder. Our theoretical approach allows the observation of third-harmonic radiation as well as boundary-enabled second-harmonic responses.9 As will be shown here, clear differences emerge when comparing the results obtained from the hydrodynamic formulation with the phenomenological

Received: October 2, 2014 Published: December 12, 2014

-^p-ACS Publications

© 2014 American Chemical Society

8

dx.doi.org/10.1021/ph500362y I ACS Photonics 2015, 2, 8-13

model, in both the properties of nonlinear harmonic generation and the associated radiation pattern. These differences are discussed and reinforce the importance of using the presented approach in the nonlinear optical study of nanoscale objects. Moreover, it will be shown that frequently used approximated models should be reconsidered. No approximations on the pulse shape and the time dependence of electromagnetic fields were made in the model, applicable unless the hydrodynamic description of the metal breaks down, e.g., for very short time scales not addressed here. The presented approach also enables studying carrier envelop phase effects and other ultrafast interaction dynamics.10

Self-Consistent Formulation of the Electromagnetic-Hydrodynamic Problem. The interaction of electromagnetic fields with objects made from arbitrary (nonmagnetic) materials is described in terms of the induced polarization P via the wave equation

V X V X E(r, t) + 44E(r, t) + r0dp(r, t) = 0 c

(1)

where E(r,t) is the electric field, c is the speed of electromagnetic waves in a vacuum, and ju0 is the vacuum permeability. In general, the coordinate-dependent polarization term contains all the information on both linear and nonlinear contributions, also including chromatic dispersion. In the framework of the hydrodynamic model, P is introduced via polarization currents, which are defined with the help of natural hydrodynamic variables: the macroscopic position-dependent electron density n(r,t) and velocity v(r,t). The basic set of hydrodynamic equations is then given by8

men(dtv + v-Vv) + ymenv = -en(E + v X H) — Vp dtn + V-(nv) = 0

(2)

where me and e are the electron mass and charge, respectively, y is the effective scattering rate, representing optical losses in a phenomenological way, and p = (3n2)2/i(h2/5me)ns/i is the quantum pressure evaluated within the Thomas-Fermi theory of an ideal Fermionic gas. The v-Vv term is the convective acceleration (in analogy to fluid dynamics) and is one of the key contributors to the nonlinear harmonic generation process. The electromagnetic (eq 1) and hydrodynamic (eq 2) sets of equations are coupled via the microscopic polarization term

TM plane waves by an infinitely long cylinder of nanoscale diameter is considered (Figure 1). The simulation domain

dp = -env

(3)

Equations 1-3 provide a self-consistent formulation of nonlinear optical processes originating from free conduction electrons in plasmonic systems. The effects of both surface nonlinearities and nonlocality are taken into account via the boundary conditions imposed by Maxwell's equations and vanishing current perpendicular to the boundaries. The proposed approach is a new nonperturbative description of free-electron nonlinearities allowing accounting for all hydro-dynamic processes of the electron plasma. It is worth noting that, to the best of our knowledge, all previous studies addressed only linear, single-frequency regimes, while nonlinear contributions were taken into account perturbatively via quasiFourier transform techniques.

Details of the Numerical Model. While the proposed method is universal and enables addressing any geometry, here, as a particular example, the problem of nonlinear scattering of

Figure 1. Simulation layout of nonlinear optical generation from a gold cylinder illuminated by a short pulse.

includes both near- and far-field regions, thus taking into account possible field enhancements in the near-field of the nanoparticle. This will allow for a comparison of our results to other known approaches developed for metal nanoparticles. Tabulated physical constants were used for the implementations of eqs 1-3, and the constants for gold were taken to be n0 = 5.98 X 1028 m-3, y = 1.075 X 1014 s-1, and mp = 9.0834 eV. In the linear regime, eqs 2 and 3 provide the general Drude-like response of the electron gas. This was numerically tested by comparing the low-intensity linear scattering in the full hydrodynamic model with the linear scattering simulated as when the Drude response is conventionally introduced via the dielectric permittivity. The nonlinear terms in the full hydrodynamic model become significant only under high-intensity excitation.

The set of eqs 1-3 was numerically solved by employing the finite element method. A Gaussian pulse of the form Em(y,t) = (0,E®) exp[- y2/(2w2)] exp[-(t - t0)2/(2T2)] cos[mt] was considered at the fundamental driving frequency of m = 1.257 X 1015 rad/s, corresponding to a free-space wavelength of X = 1500 nm («0.83 eV), with a temporal width t and a spatial width w = X/2. The pulse has a linear polarization in the y-direction and is incident on the metal cylinder of radius r = 100 nm along the x-direction (Figure 1). The time offset t0 = -3t and simulation time span T = 7t were chosen so that the scattered light pulse containing higher harmonics is able to entirely propagate across the 6 X 6 ^m2 simulation domain. Additionally, it was checked that any further increase of the simulation time span does not affect the results. A maximum peak intensity of I¡¡ = 9 X 1018 W/m2, corresponding to a typical intensity the particle can still withstand due to surface-induced ablation effects11 (with the correction for the incident light wavelength), was used to ensure convergence of the model, for any experimentally relevant intensity.

Nonlinear Spectrum. The time dependence of the nonlinear scattered signal obtained after subtracting the linear scattered excitation field was probed in the near and far field regions. Subsequently, it was Fourier transformed in the frequency domain exhibiting typical multiple harmonics spectra (Figure 2a). For a very short pulse (t =5 fs), the separation between harmonics is comparable to their spectral width and shows the overlap between second and third harmonics (SH

Figure 2. (a) Nonlinear scattering spectra from the infinitely long Au cylinder of 100 nm radius simulated for the excitation pulse of T = 5 fs (top) and T = 20 fs (bottom). (b, c) Integrated flux F2q for (b) SHG and (c) THG normalized to the maximum incident fundamental flux Fq (corresponding to the peak intensity IQ) as a function of the fundamental light intensity.

and TH) spectra (Figure 2a, top panel). For longer excitation pulses (t = 20 fs), nonlinear harmonics are well separated (Figure 2a, bottom panel). The spectral widths of the SH and TH signals are larger than those of the fundamental pulse, as the temporal nonlinear pulses are shorter due to the nonlinear dependence on the excitation field. For further comparison, the pump pulse energy at both pulse widths was normalized. As a result, the peak intensity for the t = 5 fs pulse is larger than that observed for t = 20 fs pulse, leading to an increased ratio between TH and SH intensities for the shorter pulse. In the short-pulse excitation regime, a small (several tens of nm) blue shift of both SH and TH peaks is observed with respect to the respective harmonics of the central fundamental frequency. This is due to the dispersion of the imaginary part of the frequency-dependent permittivity of Au.

The power dependence of the nonlinear harmonic generation was determined via the spectral dependence of the far-field intensity LE(®)I2 integrated over the harmonic spectra to extract the relative intensities of the nonlinear signals. The latter were then integrated over a cylindrical surface of radius 2A centered on the rod to obtain the total fluxes of the generated harmonics. As a crucial proof of the validity of the model, it was found that those dependencies follow the well-defined powers of the pump, quadratic and cubic for SH and TH, respectively, as expected for nonlinear optical processes (Figure 2b,c). Since the THG intensity grows with the third power of the fundamental intensity and SHG with the second, the former shows comparably faster growth with the excitation intensity (Figure 2a). The same results are obtained for different integration cylinder radii of A and 3/2A, confirming that the integration was performed in the radiation zone. From the slopes of the dependences one can evaluate the effective nonlinear susceptibilities of the Au nanorod to be x^ « 10-13

m/V and « 10 26 m2/V2. The value is consistent with the nonlinear susceptibility of Au examined in ref 12 taking into account fUndamental-pulse duration and wavelength, while the second-order susceptibility, although strongly geometry-dependent, is of a value that is also typical for nanostructured objects.

Second Harmonic Generation. A frequency domain analysis of the hydrodynamic equations allows deriving surface polarizabilities9 for a particle under the undepleted pump assumption. The effective nonlinear surface polarizability can also be introduced phenomenologically and related to the experimental data. The SHG from nanoparticles was intensively studied in the quasistatic limit13'14 as well as relying on the extended Mie theory15'16 with one of the central points being to account for nonlocal and retardation effects responsible for the radiation pattern formation. Advanced numerical modeling carried out in the frequency domain enables addressing arbitrary particle geometries,17 but restricted to an a priori chosen model for the nonlinear response.

In order to compare our time-domain approach to existing models, we have studied the second-harmonic radiation pattern calculated with the help of our full hydrodynamic description (Figure 3a). The SH intensity distribution was calculated by spectral frequency filtering and represents the intensity LE(®)I2 of the SH signal at each simulation domain point integrated over the spectral spread of the SH resonance. The emission diagram has two radiation lobes pointing predominantly in the vertical (y-) direction and showing the appearance of a modified dipole radiation. It is worth noting that in the chosen geometry, we can identify the main contribution to the SHG from the convective acceleration and the Lorenz force with comparable contributions, while the quantum pressure effects, expected to result in 5/3 harmonic generation, were not

Hydrodynamic Effective

(a) (b)

3 <x

1 |jm

(c) (d)

M

0.0 0.2 0.4 0.6 0.8 1.0 Intensity (a.u.)

Figure 3. Intensity distributions of the nonlinear scattering around the infinitely long Au cylinder of 100 nm radius simulated for the excitation pulse of T = 20 fs for (a, b) SHG and (c, d) THG in (a, c) microscopic hydrodynamic and (b, d) phenomenological effective models. The color scale is internal for each plot.

observed (the quantum pressure term enters eq 2 and is proportional to the carrier concentration in the power of 5/3; hence the appearance of a fractional harmonic might be expected). Note that the relative contribution of the different terms may differ depending on a particular geometry.

To compare this radiation pattern to that obtained from a phenomenological model, Figure 3b shows the SH radiation pattern evaluated using a two-step model, in which the fundamental field inside the particle is calculated in a first step and the nonlinear field distribution is subsequently derived in a second step using a source of the surface nonlinear polarization Pf = yXsurf,in(Epump,i)2, where 1 stands for the local normal to the surface. This assumes the nondepleted pump regime, relies on a quasi-Fourier separation of the harmonics, and makes the explicit restricting assumption on the interaction nature to originate from local boundary terms only. The latter follows from the centrosymmetric nature of considered particle's material, so that the dipolar SH radiation from the bulk is forbidden by the selection rules. In this case, one can clearly observe four lobes determined by a quadrupole-like emission added to a dipolar (two lobes) contribution. The phase relations between the two contributions determine the directivity of total SH emission. One can see the similarities of the nonperturbative and phenomenological models (cf. Figure 3a and b); however, the ratio between the dipolar and quadrupolar contribution is heavily distorted in favor of the latter in the phenomenological model. The phenomenologicaly defined surface polarization term also results in an unphysical singularity at the boundary, failing to describe the near-field distribution of the nonlinear source; conversely, the full timedomain method reproduces the physics of surface interactions with much better accuracy in both near- and far-field regions.

Third Harmonic Generation. Bulk third-harmonic generation in metals has been considered by the introduction of a nonlocal ponderomotive force, acting on electrons subjected to an electromagnetic field gradient. Higher-harmonic generation at boundaries (e.g., flat metal surfaces), based on the

Sommerfeld free-electron models with subsequent solution of the Schrodinger equation in the Kramers—Henneberger accelerating frame, were also developed19 and are in good agreement with experimental data.20 However, the impact of complex geometries on higher-harmonic generation is virtually not studied to date, primarily due to both experimental difficulties and the complex theoretical treatment required. The method developed here enables efficiently addressing this theoretical challenge.

The TH radiation pattern generated by the nanorod (Figure 3c) possesses strong beaming characteristics in the forward scattering direction as expected for large particles, for which retardation effects lead to constructive interferences in the forward direction only, lowering the backward scattered intensity.21 In order to perform a comparison between the time-domain approach and phenomenological model, the TH intensity distribution was simulated with the above-described two-step model using a nondispersive, bulk third-order susceptibility. The radiation pattern derived phenomenologi-cally has similar features to the microscopic model with a stronger backward THG (Figure 3d). The same as for the SHG pattern, strong differences are observed in the near field of the nanocylinder; while the phenomenological model assumes a third-order susceptibility to be homogeneous across the nanorod, the differences between models may be indicative of a position-sensitive effective third-order susceptibility arising from the hydrodynamic description.

For the geometry and the pump powers considered here, the signature of the fourth harmonic has been reliably observed, but higher harmonics are beyond the numerical accuracy of the simulations.

Resonance Effects. By tuning the nanoparticle geometry, the spectrum of its resonances can be engineered on demand.22,23 This leads to enhanced electromagnetic fields and scattering, and the combination of various resonant conditions at both fundamental and harmonic frequencies can be used to significantly modify the nonlinear scattering.1 Even at moderate illumination intensities, resonant effects start to be especially important in small metal nanostructures, as they enhanced the effective nonlinear surface response.24 In order to study the impact of localized surface plasmon (LSP) resonances at either fundamental or nonlinear harmonic frequencies, either the particle shape or material can be modified.25 We have chosen to change the electron concentration of the metal while keeping the particle geometry and dimensions the same, in order to obtain a plasmonic dipolar resonance in free space (£metai = —1) at either the fundamental, second, or third harmonic frequencies. This allows us to compare the effects of the resonant field enhancement while not significantly influencing the carrier distribution n(r) within the particle, which would be the case for varying geometry. It should be noted, however, that changing the electron density influences not only the plasma frequency but also the nonlinear susceptibilities (eq 1).17 The nonlinear generation efficiency—the parameter that depends on the _x(2) and _x(3) of the medium as well as the excitation conditions—was evaluated in all the studied cases via a numerical estimation of the slopes on power-dependent graphs similar to those in Figure 2b,c. The enhancements were calculated relative to the geometry considering the natural electron concentration described above.

Both second and third harmonics show the highest enhancements, 5 and 10 orders of magnitude, respectively, for the 20 fs excitation pulse, when the nanocylinder has a LSP

resonance at the fundamental frequency. For the LSPs at the SH frequency, the observed enhancement is 4 and 6 orders of magnitude, and at the TH frequency, it is 2 and 5 orders of magnitude for SHG and THG efficiencies, respectively. For all the LSP resonances, the ratio between TH and SH intensities increases, so that at the fundamental LSP frequency TH intensity is 1 order of magnitude higher than the SH intensity at the highest fundamental intensity studied, with comparable TH and SH intensities for other resonant cases, in contrast to the nonresonant situation where SHG dominates. This indicates the stronger importance of the excitation field enhancement over the increased scattering cross-section of the nonlinear harmonics. At the fundamental frequency resonance, the effective "enhanced" susceptibility corresponds to ^(2) « 10-10 m/V and ^(3) « 10-21 m2/V2. Taking into account that these were obtained with the reduced electron concentration (n/n0 = 60), the effective nonlinearity of Au nanoparticles in these resonant conditions may be comparable to the best nonlinear dielectrics, even under femtosecond excitation.

In conclusion, we presented a comprehensive time-domain treatment of the microscopic polarization of conduction electrons using a full hydrodynamic description allowing for a self-consistent modeling of both the linear and nonlinear response of plasmonic nanostructures including the generation of multiple harmonics. The effects originating from the convective acceleration, the magnetic contribution of the Lorenz force, the quantum electron pressure, and the presence of the nanostructure's boundaries have been taken into account, leading to the appearance of both second and third harmonics. The developed method provides an ultimate approach to investigate the nonlinear response of arbitrarily shaped complex nanoscale plasmonic structures and enables addressing their self-consistent nonlinear dynamics. While various approaches typically make restricting assumptions on the nonlinear dynamics, such as undepleted pump or uncoupled frequencies approximations, in order to simplify the solutions of the coupled nonlinear equations, the proposed method enables one to obtain a universal, self-consistent numerical solution, free from any approximations. Moreover, it should be emphasized that the straightforward perturbative hydrodynamic description (e.g., relying on phenomenological second-order polarizability) is inconsistent with the full model reported here, underlining the fact that the majority of previously employed approaches should be reconsidered. Furthermore, nanostructures with arbitrary geometry and resonant frequencies for excitation, nonlinear scattering, or both can be comprehensively studied using the proposed model. The developed nonperturbative model enables investigating a vast number of multidisciplinary problems, involving metal composites interacting with weak, moderate, or intense optical pulses. The role of nonlocal electromagnetic response,26 the effects accounting for the response of core electrons,27 and employing quantum approaches for electron exchange correlations28 can also be straightforwardly included.29 Furthermore, the developed formalism paves the way for investigating ultrafast dynamics in mesoscopic and nanoscopic systems with properties defined via microscopic degrees of freedom, which can be introduced in the permittivity model.

■ AUTHOR INFORMATION Corresponding Author

*E-mail: [e-mail].

Author Contributions

§P. Ginzburg and A. V. Krasavin contributed equally. Notes

The authors declare no competing financial interest.

■ ACKNOWLEDGMENTS

This work has been supported, in part, by EPSRC (UK). A.Z. acknowledges support from the Royal Society and the Wolfson Foundation. P.G. and A.Z.'s work was supported by the U.S. Army Research Office. G.W. is grateful for support from the People Programme (Marie Curie Actions) of the EC FP7 project 304179. The authors are thankful to Michael Scalora for the discussions.

■ REFERENCES

(1) Kauranen, M.; Zayats, A. V. Nonlinear plasmonics. Nat. Photonics 2012, 6, 737-748.

(2) Wokaun, A.; Bergman, J. G.; Heritage, J. P.; Glass, A. M.; Liao, P. F.; Olson, D. H. Surface second-harmonic generation from metal island films and microlithographic structures. Phys. Rev. B 1981, 24, 849-856.

(3) Kim, E. M.; Elovikov, S. S.; Murzina, T. V.; Nikulin, A. A.; Aktsipetrov, O. A.; Bader, M. A.; Marowsky, G. Surface-enhanced optical third-harmonic generation in Ag island films. Phys. Rev. Lett. 2005, 95, 227402.

(4) Zhang, Y.; Grady, N. K.; Ayala-Orozco, C.; Halas, N. J. Three-dimensional nanostructures as highly efficient generators of second harmonic light. Nano Lett. 2011, 11, 5519-5523.

(5) Fan, W.; Zhang, S.; Panoiu, N.-C.; Abdenour, A.; Krishna, S.; Osgood, R. M.; Malloy, K. J.; Brueck, S. R. J. Second harmonic generation from a nanopatterned isotropic nonlinear material. Nano Lett. 2006, 6, 1027.

(6) Ginzburg, P.; Krasavin, A.; Zayats, A. V. Cascaded second-order surface plasmon solitons due to intrinsic metal nonlinearity. New J.

Phys. 2013, 15, 013031.

(7) Boyd, R. W. Nonlinear Optics, 3rd ed.; Academic Press, 2008.

(8) Sipe, J. E.; Y.So, V. C.; Fukui, M.; Stegeman, G. I. Analysis of second-harmonic generation at metal surfaces. Phys. Rev. B 1980, 21 , 4389.

(9) Dadap, J. I.; Shan, J.; Eisenthal, K. B.; Heinz, T. F. Second-harmonic Rayleigh scattering from a sphere of centrosymmetric material. Phys. Rev. Lett. 1999, 83, 4045-4048.

(10) Krausz, F.; Ivanov, M. Attosecond physics. Rev. Mod. Phys. 2009, 81, 163.

(11) Pletch, A.; Kotaidis, V.; Lotenc, M.; Boneberg, J. Femtosecond laser near-field ablation fromgold nanoparticles. Nat. Phys. 2006, 2, 44.

(12) Boyd, R. W.; Shi, Z.; De Leon, I. The third-order nonlinear optical susceptibility of gold. Opt. Commun. 2014, 326, 74-79.

(13) Dadap, J. I.; Shan, J.; Heinz, T. F. Theory of optical second-harmonic generation from a sphere of centrosymmetric material: small-particle limit. J. Opt. Soc. Am. B 2004, 21, 1328-1347.

(14) Makeev, E. V.; Skipetrov, S. E. Second harmonic generation in suspensions of spherical particles. Opt. Commun. 2003, 224, 139-147.

(15) Pavlyukh, Y.; Hubner, W. Nonlinear Mie scattering from spherical particles. Phys. Rev. B 2004, 70, 245434.

(16) Capretti, A.; Forestiere, C.; Dal Negro, L.; Miano, G. Full-wave analytical solution of second-harmonic generation in metal nano-spheres. Plasmonics 2014, 9, 151-166.

(17) Bachelier, G.; Russier-Antoine, I.; Benichou, E.; Jonin, C.; Brevet, P.-F. Multipolar second-harmonic generation in noble metal nanoparticles. J. Opt. Soc. Am. B 2008, 25, 955-960.

(18) Ginzburg, P.; Hayat, A.; Berkovitch, N.; Orenstein, M. Nonlocal ponderomotive nonlinearity in plasmonics. Opt. Lett. 2010, 35, 15511553.

(19) Varro, S.; Ehlotzky, F. Higher-harmonic generation from a metal surface in a powerful laser field. Phys. Rev. A 1994, 49, 3106-3109.

(20) Farkas, Gy.; Toth, Cs.; Moustaizis, S. D.; Papadogiannis, N. A.; Fotakis, C. Observation of multiple-harmonic radiation induced from a gold surface by picosecond neodymium-doped yttrium aluminum garnet laser pulses. Phys. Rev. A 1992, 46, R3605.

(21) Bohren, C. F.; Huffman, D. R. Absorption and Scattering of Light by Small Particles; Wiley VCH, 1998.

(22) Ginzburg, P.; Berkovitch, N.; Nevet, A.; Shor, I.; Orenstein, M. Resonances on-demand for plasmonic nano-particles. Nano Lett. 2011, 11, 2329—2333.

(23) Prodan, E.; Radloff, C.; Halas, N. J.; Nordlander, P. A. Hybridization model for the plasmon response of complex nanostructures. Science 2003, 302, 419—422.

(24) Ginzburg, P.; Krasavin, A.; Sonnefraud, S.; Murphy, A.; Pollard, R. J.; Maier, S. A.; Zayats, A. V. Nonlinearly coupled localized plasmon resonances: resonant second-harmonic generation. Phys. Rev. B 2012, 86, 085422.

(25) Berkovitch, N.; Ginzburg, P.; Orenstein, M. Nano-plasmonic antennas in the near infrared regime. J. Phys.: Condens. Matter 2012, 24, 073202.

(26) Ginzburg, P.; Zayats, A. V. Localized surface plasmon resonances in spatially dispersive nano-objects: phenomenological treatise. ACS Nano 2013, 7, 4334—4342.

(27) Scalora, M.; Vincenti, M. A.; de Ceglia, D.; Roppo, V.; Centini, M.; Akozbek, N.; Bloemer, M. J. Second- and third-harmonic generation in metal-based structures. Phys. Rev. A 2010, 82, 043828.

(28) Crouseilles, N.; Hervieux, P.-A.; Manfredi, G. Quantum hydrodynamic model for the nonlinear electron dynamics in thin metal films. Phys. Rev. B 2008, 78, 155412.

(29) Hurst, J.; Haas, F.; Manfredi, G.; Hervieux, P.-A. High-harmonic generation by nonlinear resonant excitation of surface plasmon modes in metallic nanoparticles. Phys. Rev. B 2014, 89, 161111(R).

Scattering suppression from arbitrary objects in spatially dispersive layered metamaterials

Alexander S. Shalta,1,2,3^ Pavel Ginzburg,4,t Alexey A. Orlov,1 Ivan Iorsh,1 Pavel A. Belov,1 Yuri S. Kivshar,1,5 and Anatoly V. Zayats4 1ITMO University, St. Petersburg 197101, Russia 2 Kotel'nikov Institute of Radio Engineering and Electronics ofRAS (Ulyanovsk branch), Ulyanovsk 432011, Russia 3 Ulyanovsk State University, Ulyanovsk 432017, Russia 4 Department of Physics, King's College London, Strand, London WC2R 2LS, United Kingdom 5Nonlinear Physics Center and Center for Ultra-high Bandwidth Devices for Optics Systems (CUDOS), Australian National University,

Canberra ACT 0200, Australia (Received 6 June 2014; revised manuscript received 30 December 2014; published 19 March 2015)

Concealing objects by making them invisible to an external electromagnetic probe is coined by the term "cloaking." Cloaking devices, having numerous potential applications, are still facing challenges in realization, especially in the visible spectral range. In particular, inherent losses and extreme parameters of metamaterials required for the cloak implementation are the limiting factors. Here, we numerically demonstrate nearly perfect suppression of scattering from arbitrary-shaped objects in spatially dispersive metamaterial acting as an alignment-free concealing cover. We consider a realization of a metamaterial as a metal-dielectric multilayer and demonstrate suppression of scattering from an arbitrary object in forward and backward directions with perfectly preserved wave fronts and less than 10% absolute intensity change, despite spatial dispersion effects present in the composite metamaterial. Beyond the usual scattering suppression applications, the proposed configuration may be used for a simple realization of scattering-free detectors and sensors.

DOI: 10.1103/PhysRevB.91.125426

I. INTRODUCTION

Controlling scattering from an object illuminated by an external wave is one of the prominent important objectives of applied electromagnetism. Various approaches have been developed and employed in antenna engineering for the enhancement of scattering cross sections and directivities which are the most demanded characteristics [1]. While the main technological effort has been focused at the millimeter-and radio-frequencies ranges where the majority of applications are, considerable interest in the optical frequencies has recently emerged. Optical antennas have been already shown to manipulate radiation properties of single emitters as well as to enhance absorption cross sections [2].

The reduction of scattering cross sections in a specifically designed material environment leads to reduced detectability of objects, and may, ultimately, result in "invisible" objects if the scattering is absent. The concept of "cloaking" was introduced in [3,4] and gained considerable attention due to a continuous demand to achieve invisibility for radar waves [5] and visible light [6-8]. The general approaches for cloaking rely either on transformation optics concepts [3] or conformal optical mapping of complex electromagnetic potentials [4]. While the former generally results in the requirements of highly anisotropic (and sometimes singular) electric and magnetic susceptibilities of a medium of a cloak, the latter approach requires a position-dependent refractive index variation.

One of the main challenges in the development of practical cloaking devices is to minimize demands on permeabilities and anisotropy and reduce inherent material losses of required

'Corresponding author: [e-mail] t [e-mail]

PACS number(s): 78.67.Pt, 42.25.Fx, 42.79.Wc

materials. The so-called "carpet cloak" has been proposed and implemented to mitigate these factors by imposing certain geometrical restrictions and utilizing quasiconformal mapping [9,10]. A qualitatively different approach to cloaking utilizes e-near-zero (ENZ) metamaterials to suppress a dipolar scattering of a concealed object [11]. ENZ regime, where a real part of the permittivity is close to zero, can be achieved in anisotropic configurations where a wave with certain polarization does not have phase advancement. While a number of the ENZ metamaterial realizations exist in the radio-frequency range [12,13], in optics such an anisotropic response may be achieved through metal-dielectric layered structures [14], semiconductor heterostructures [15], or vertically aligned arrays of nanorods [16]. Recently, an idealized homogeneous uniaxial ENZ material was proposed for partial cloaking [17], but its realization was not specified. In particular, the presence of inherent strong spatial dispersion, associated with realistic plasmonic metamaterial geometry, was not addressed.

In this paper, we demonstrate suppression of scattering from arbitrary-shaped and large (not necessarily subwavelength) objects placed inside a layered metal-dielectric metamaterial. We show that a metamaterial realization has a major influence on the scattering phenomenon since the electromagnetic response of the plasmonic multilayers is substantially affected by spatial dispersion effects [18] but does not impede scattering suppression. Investigation of the exact numerical model, taking into account material losses and finite dimensions of the metamaterial realization, shows the possibility of almost perfect scattering suppression from arbitrary-shaped objects with minimal variations of the phase front of the transmitted/reflected optical wave and small amplitude modifications. The proposed scheme does not require extreme electric and magnetic susceptibilities and can operate in an alignment-free manner. This approach is qualitatively different from other proposals when the incident light does not interact with the

1098-0121/2015/91(12)/125426(7)

125426-1

©2015 American Physical Society

cloaked object but is bent around it by a material layer [3,5] or when the scattered field is suppressed with another antiparallel dipole [11]. In our approach, the light scattering by an object placed inside a metamaterial is strongly anisotropic and suppressed in the forward (transmission) and backward (reflection) directions.

The paper is organized as follows. First, various dispersion regimes in layered metamaterial composites are studied. This is followed by investigations of a dipole radiation pattern inside the metamaterial. Subsequently, concealing of objects of various shapes and dimensions is numerically analyzed in two-dimensional geometries. The scattering properties in a three-dimensional case are studied in the final section where both numerical modeling and fully analytical formulation of the scattering problem in the dipolar approximation underline advantages and limitations of the proposed approach.

II. RESULTS AND DISCUSSION

A. Dispersion of electromagnetic waves in plasmonic multilayer metamaterials

The effective permittivity tensor of metamaterials based on multilayers can be directly evaluated from their thickness, periodicity, and optical parameters of the constitutive materials [18]. While this approach can predict optical properties of the composites made of low-contrast, positive permittivity layers, it faces severe challenges once negative permittivity materials (plasmonic metals) are involved. Optical properties of such multilayered composites are influenced by guided surface electromagnetic modes on metal-dielectric interfaces, surface plasmon polaritons (SPPs), and spatial dispersion effects become especially significant in the parameter range where the effective medium theory (EMT) predicts vanishing values of permittivity [19].

SPP modes in planar layered composites, as well as their dispersion, can be found via semianalytical formulation using the transfer-matrix method, and the dispersion equation of eigenmodes in arbitrary periodic multilayered structures is given by [18]

2cos[ß (dd + dm)] =

cos (ßddd)cos(ßmdm)

- ßd sin (ßddd) sin (ßmdm ) Pm

(1)

where j is the Bloch wave number, dd and dm are the thicknesses of dielectric and metal layers, respectively, and jm and jd are the transverse parts of the wave vectors in metal and

dielectric materials, respectively, pmj = ^smj(m/c)2 - kx, where kx is the transverse wave vector in the direction parallel to the layers [18]). Solution of Eq. (1) enables one to obtain isofrequency contours of the modes, allowed to propagate in the layered medium for the frequency fixed at arbitrary value. Hereafter, generic metamaterial parameters were chosen to be dm = 20 nm for the Au [20] layer, and dd = 100 nm for air (note that any dielectric material could be used with proper recalculation of the cloak parameters).

The isofrequency curves of an infinite metal-dielectric multilayer structure for the extraordinary, TM polarized,

(a)

0

(b)

(c) 1

0

(d)

(e)

■ ^ 0

-1 0 1

k

X

air

U k effective

medium

-2-1 0 1 2

FIG. 1. (Color online) (a,c,e) Isofrequency contours for the infinite, multilayered metamaterial simulated for (a) the wavelength of 600 nm, (c) the operational wavelength of 540 nm, and (e) 500 nm wavelength. The metamaterial consists of the air and Au layers of da = 100 nm and dg = 20 nm thickness, respectively. Red curves show propagating modes in the layered metamaterial: (bright) main modes; (pale) additional modes. Green circle represents light cone in air. Yellow curves are the isofrequency contours of the same metamaterial simulated in the EMT. (b),(d),(f) The electric field (y-component) distributions of the radiating y -polarized dipole situated in the center of the metamaterial block of 4.8 x 4 /m size at the wavelength corresponding to the isofrequency contours in (a), (c), and (e). Wave numbers kx and ky are normalized to the Brillouin zone boundaries n/(da + dg). White lines represent the power flow.

modes have strong dependence on the wavelength, and different characteristic dispersions can be observed [Figs. 1(a), 1(c), and 1(e)]. The distinctive feature of these plots is the simultaneous emergence of multiple bands at the same frequencies. In other words, for certain values of ky, there exist two nontrivial solutions for kx. This effect is the manifestation of strong spatial dispersion, inherently attributed to plasmonic excitations. The transition between various dispersion regimes occurs in the vicinity of the ENZ point, predicted by the EMT calculations (around 540 nm for considered multilayer parameters). According to the EMT, isofrequency surfaces at the wavelengths longer than the ENZ wavelength have pure hyperbolic shape [Fig. 1(a)], while at slightly shorter wavelengths the dispersion is purely elliptic [Fig. 1(e)].

At the ENZ wavelength, the EMT-obtained isofrequency contours are strongly prolate spheroids [Fig. 1(c)] due to the finite losses in the metal layers (in the lossless case,

y

x

x

it degenerates to a single line). The exact solution of the dispersion relation [Eq. (1)] shows that the EMT description fails near the ENZ frequency. First, the finite period of the structure imposes limitations on the range of k vectors, in the direction normal to the layers. Thus, the hyperbolic shape of the dispersion will deviate close to the edges of the Brillouin zone, ky = n/(dm + dd). Furthermore, due to the presence of SPP modes on the multilayer interfaces, strong inherent spatial dispersion of the metamaterial will result in the coexistence of two TM polarized waves. Hence, there are two dispersion branches for the same polarization, contrary to conventional media [19,21]. The mode that follows predictions of the EMT is the main mode, while another is called "additional." The behavior of structural spatial dispersion in layered composites is similar to nanorod metamaterials [21,22]. It should be noted, however, that this type of spatial dispersion has a pure "structural" nature and is not related to nonlocalities in material components [23,24].

B. Dipole radiation in metal-dielectric multilayered metamaterial

Scattering of electromagnetic waves by an arbitrary object can be evaluated using multipolar decomposition of the field taking into account the geometry of the problem [25]. The most significant contribution to the far-field scattering usually emerges from the dipolar term. Following this rule of thumb, we first study the properties of the scattering by a point dipole inside the multilayer metamaterial, before considering a general case of concealment of an arbitrary object.

It is well known, that the dipole radiation patterns are strongly affected by electromagnetic environment. For example, periodically structured media such as photonic crystals with engineered dispersion may enforce flat fronts of dipole radiation [26]. Furthermore, radiation patterns of an emitter situated inside a homogeneous, but anisotropic media may have quite remarkable shapes. The characteristic radiation pattern of a dipole in free space (a spherical wave front modulated by the cosine of the angle of the dipole axis orientation) becomes cross shaped if placed in a hyperbolic medium, with the opening angle being solely defined by the ratio of the permittivity tensor components (e.g., Ref. [27]). The details of this cross-shaped radiation pattern significantly depend on specific metamaterial realization.

We used the scattered-field formulation in the finite element modeling software [28] to simulate the field distributions inside the simulation domain which was surrounded with perfectly matched layers (the periodic boundary conditions were not used in the simulations). The radiation patterns of a vertically polarized (y-oriented) dipole, placed inside a multilayered metamaterial block have qualitatively different behaviors in the hyperbolic, elliptic, and ENZ regimes (Fig. 1). While all the patterns have a cross-shaped envelope, the wave fronts in the hyperbolic and elliptic regimes are significantly curved, whereas the dispersion near the ENZ frequency maintains the wave front flat. Moreover, the metamaterial was designed to be impedance matched with a free space at normal incidence to minimize reflections on the boundaries (the circle representing the dispersion of light in free space passes the vicinity of the degeneracy point of the metamaterial dispersion

for this choice of parameters [Fig. 1(c)]. As a result, multiple reflections inside the layered composite are mostly suppressed. The impedance matching design and the flat dipole radiation front inside the composite provide the necessary ingredients for maintaining a flat front of the incident wave after interaction with an object and low reflection from the metamaterial itself. The presence of the additional waves in the nonlocal regime inherent to any realistic realization of the metamaterial does not impede flatness of the forward and backward radiated wave fronts.

C. Scattering suppression from two-dimensional objects

A plane-wave-like radiation pattern of a point dipole directly implies that the major dipolar scattering of an object will not distort the front of an incident wave. In order to verify this, we will first consider scattering from perfect electric conductor (PEC) cylinders of different sizes. The cylinder is placed in the center of a layered metamaterial block of a finite volume (to simulate a realistic scenario) which is illuminated by a plane wave polarized perpendicular to the layers and with the wave vector along them. It is worth noting, that higherorder multipoles also contribute to the scattering, especially in the case of larger objects. Nevertheless, the scattering suppression still can be observed, as will be confirmed in numerical simulations.

The electric field distributions in and outside the metamaterial when the cylinder is placed inside and for the empty metamaterial cover show that the flat wave front of the incident wave is maintained in both cases, while the amplitude slightly decreases due to inherent material losses of the metal layers [Figs. 2(a) and 2(b)]. Since the metamaterial is impedance matched to free space, as was shown above, reflections from both boundaries are almost completely suppressed. Figure 2(c) shows the difference of the fields between the case of empty cloak and while cloaking the PEC cylinders of 200 nm diameter, demonstrating nearly perfect cancellation of the scattering from the object in the far field: The amplitude difference does not exceed 3%. It should be emphasized that the general detectable signature of an object is the wave-front distortion of an incident illuminating wave and/or intensity variations across the wave front (shadow areas). At the same time, absolute transmitted or reflected intensity is hardly quantifiable without reference measurements. Thus, preserving a phase front and an amplitude front undistorted, the proposed approach is perfect for simple and practical implementation of scattering suppression. The same considerations are applicable for observation in reflection with the impedance-matched metamaterial cover (Fig. 2).

Large PEC cylinders and arbitrary-shaped objects have also been tested [Figs. 2(d)-2(g)] as well as strongly absorbing and high-index transparent objects of different shapes (not presented, but showing very similar far-field distributions). Nearly perfect scattering suppression was observed in all scenarios with nearly the same performance [cf. Figs. 2(e) and 2(g)]. The amplitude difference observed for a large arbitrary-shaped object is about 7% [Fig. 2(g)], only slightly larger than for the big cylinder, showing about 5% variations [Fig. 2(e)]. Furthermore, increasing object sizes leads to the increased amplitude difference due to the stronger contribution from the

-4 -3 -2 -1 0 1 2 3 4 -4 -3 -2 -1 0 1 2 3 4

x (|jm) x (um)

x (|jm) x (|jm)

FIG. 2. (Color online) Electric field (y-component) distribution for the plane wave illumination of (a) the empty cloak, (b,d,f) PEC object inside the cloak: cylinder of (b) 200 nm and (d) 500 nm diameters and (f) arbitrary-shaped object. The field distribution differences between (a) and (b,d,f) are shown in (c,e,g), respectively. The metamaterial cloak parameters and the wavelength correspond to the dispersion in Fig. 1(b). White lines represent the power flow. The simulations are performed in 2D geometry.

multipole scattering terms. While the above considerations were presented for normally incident collimated beams, we also checked the angular dependence of the scattering. While scattering from the object is suppressed, the same as for normal incidence, the angular dependent impedance matching leads to a strong shadow at larger angles of incidence produced by the metamaterial cover itself.

As an intuitive explanation of the described effect, one can consider that the light scattering by an object inside the metamaterial is strongly anisotropic. In the elliptic and hyperbolic regimes, this results in convergent and divergent wave fronts, respectively, thus, leading to "shadows" from the objects located inside the metamaterial. In the ENZ regime, the dipole radiation and, thus, scattering from the object result in the plane wave fronts due to peculiarities of the dispersion. The presence of an additional wave in the composite metamaterial, absent in the homogeneous metamaterial, leads to smaller transmission but does not destroy plane wave fronts and cloaking. As a result, the presence of the object under such

specific conditions causes only the changes of the transmitted plane wave intensity with the rest of the energy being dissipated due to material loss.

For a better understanding of the key benefits of the scattering reduction inside the metamaterial, it is useful to compare a half-transparent glass block and the one made of the layered metamaterial (our design), assuming that both structures have the same transmission coefficient (less than 100%). Now, a scattering object is placed inside each one of those blocks. The glass cover will not be transparent at the objects' location, whereas the metamaterial composite will maintain the same optical properties, as if no scatterer is present inside. It should also be noted that the importance of the own shadow of the metamaterial cloak on its performance is not unambiguous and strongly depends on the final goal of the application. For perfect invisibility, the shadow should be suppressed; however, even without this, the metamaterial may have possible applications in, e.g., static camouflage cloaking.

(a)

1 0 -1

(c) "E

0

x (Mm)

0

x (Mm)

1.5 1

0.5 0

-0.5 -1

-1.5

Uj

FIG. 3. (Color online) Electric field (y-component) distributions for the plane wave illumination of the multilayered cloak in 3D geometry with (a) the infinite and (b),(c),(d) finite (200 nm) length cylinder of 500 nm diameter inside: (c) the electric field distribution in the x-y plane and (d) the electric field distribution in the x-z plane. Inset shows scattering on the cylinder without the cloak. The metamaterial parameters and the wavelength correspond to the dispersions in Fig. 1(b).

D. Scattering suppression from three-dimensional objects

For three-dimensional objects, the influence of the third dimension introduces additional scattering channels emerging from edge effects of finite length objects placed inside anisotropic metamaterial. In order to estimate the influence of these edge effects, we compared scattering inside the metamaterial from infinite and finite (200 nm length) cylinders. In these simulations, the metamaterial cover was considered with periodic boundary conditions imposed in the direction normal to the layers (y direction) in order to mitigate the requirements on computational resources.

The scattering amplitude distribution of an infinitely long cylinder [Fig. 3(a)] possesses identical behavior to the pure two-dimensional (2D) case [Fig. 2(d)] confirming that the cloaking effect is not attributed to properties of Maxwell's equations in low-dimensional space. In the three-dimensional (3D) case, the scattering suppression takes place in both x-z and x-y planes (the x-y plane is not shown but is identical to the 2D case) as illustrated by the electric field distributions. Consequently, an observer situated outside the metamaterial slab will not detect any wave-front distortions that take place without the cover (amplitude variations are about 5%). In the case of a finite-size cylinder [Fig. 3(b)], while the cloaking in the x-y plane is still preserved [Fig. 3(c)] with the field distribution and the energy flow lines are almost identical to the 2D case, the wave front in the x-z plane is slightly distorted [Fig. 3(d)]. This is due to the effect of the cylinder edges. The diffraction on the edges gives rise to the appearance of x and z components of the scattered field and to the related energy outflow in the direction normal to the propagation direction of the incident wave. It is worth noting that this effect is diminished when long enough (larger than wavelength) particles are considered.

E. Theory of a three-dimensional homogenized metamaterial cover

In order to understand concealing performance of the metamaterial cover, we have developed a theoretical model treating the scattering process, in the dipolar approximation, in a homogenized anisotropic medium with the permittivity tensor components corresponding to the effective medium permittivity of the metamaterial realization, described above.

Electromagnetic scattering can be considered in the framework of the Green's function formulation. We first describe a two-dimensional geometry with infinitely-long cylindrical scatterers situated in a uniaxial anisotropic medium, illuminated perpendicular to its axis with the plane wave having the electric field in the y direction, also perpendicular to the cylinder axis. Using the coordinate transformations, the y component of the electric field scattered by a 2D electric dipole oriented in the y direction is given by [29]

#0(1We)y2+ X2Syy[koreHo(1)(£ore) - Hi(1We)]

Ey = p

(1)/ 0

rd)/

(2)

where re = ^x2eyy + y2exx, p is the polarizability of the dipole, k0 = to/c, Hq \ are the Hankel functions of the first kind, and en are the diagonal components of the permittivity tensor of the medium. In the limit of exx ~ 0, Eq. (2) simplifies to

Ey

H(1)(ko lx |

H(1)(ko lx |

k(( lx l £yy

(3)

Equation (3) shows that the scattered field has a plane wave front, the same since the incident field. Therefore, such a

1

configuration provides the low scattering cross section and the invisibility effect in the far-field observations since the asymptotic Hankel functions for large \x | further simplify Eq. (3) to

Ev

i - 1

eikox^syy

.Jsyyw (kox ^ëyy)1/2 '

(4)

corresponding to an unperturbed plane wave.