Nanoparticle above a Dielectric Interface - American Chemical Society

Nov 20, 2014 - Jean Lermé*. Institut Lumière Matière, UMR5306 Université Lyon 1-CNRS, Université de Lyon Bât. A. Kastler, 43 Bld du 11 Novembre 1918,...
2 downloads 0 Views 2MB Size
Article pubs.acs.org/JPCC

Nanoparticle above a Dielectric Interface: Plasmon Hybridization Model, Comparison with the Dimer System, and against Exact Electrodynamics Calculations Jean Lermé* Institut Lumière Matière, UMR5306 Université Lyon 1-CNRS, Université de Lyon Bât. A. Kastler, 43 Bld du 11 Novembre 1918, F-69622 Villeurbanne Cedex, France S Supporting Information *

ABSTRACT: In the presence of a nearby highly reflecting dielectric substrate the optical response of a metal nanoparticle exhibits a rich multipolar plasmonic structure. For interpreting and quantifying the spectral patterns, the quasistatic hydrodynamic plasmon hybridization (PH) method developed by the Nordlander group has been applied to the particle/interface system in the simple generic case of a jellium sphere (PI−PH model). In the quasistatic limit, the method of images allows the PI−PH model to be tightly related to the two-sphere dimer PH model (D−PH model). However, because the dynamics of the electron fluid in the mirror particle is enslaved to that of the particle above the interface, strong differences between both PH models result from the reduction by a factor two of the number of dynamical variables describing the electron fluids. In particular the PI−PH model involves only bonding bright modes, whereas bonding and antibonding, bright and dark modes are involved in the homodimer case. A very detailed description of both PH models is provided, with particular emphasis on the coupling with an external field and the transfer of the oscillator strength of the dipole excitations of the isolated particle(s) to high-order plasmon modes of the particle/interface or dimer systems. The predictions of both PH models are also compared, through absorption and extinction cross-section spectra, to exact results computed in the frame of the multipole expansion method, allowing the range of applicability and limits of both models to be determined. Some comments about heterodimers are also provided.

1. INTRODUCTION The optical properties of metal nanoparticles are of major importance in many technological, imaging, medical and sensing applications at the nanoscale.1−7 The unceasing growing interest for metallic species involving free conduction electrons is rooted in their plasmonic spectral response which can be tuned continuously over the entire near UV−visible-near IR spectral range through size, shape and local environment effects.3,4,7−11 From a fundamental point of view, these plasmonic properties, which manifest through narrow spectral ranges of magnified absorption/scattering peaks in optical spectra, as well as concomitant enhanced local fields, rest on the collective excitation modes (i.e., plasmons) of the confined conduction electrons in the metal. A major issue in this research field concerns the plasmonic properties of aggregates of particles, namely, groups of nanoparticles in very close mutual proximity.3,12−29 For such composite nanosystems, the cluster geometry, as well as the interparticle distances, are the main parameters ruling the overall optical spectral response of the entire system. In many cases the intrinsic properties of the individual subunits play a minor role and even can be absent from the optical spectra of the particle aggregate for specific incoming field polarization directions. For instance, in the case of dimers consisting of two gold spheres near the conducting contact limit, the longitudinal (electric field polarization parallel © 2014 American Chemical Society

to the dimer axis) excitation spectrum exhibits a very large dipolar band in the red or near IR spectral range, whereas the damped dipolar plasmon band of each individual sphere is located around 530 nm, below the interband transition threshold. Both experiment and simulation show that the spectral position of this band critically depends on the interparticle distance.21,30 Moreover, additional satellite bands in the visible spectral range emerge sequentially in the spectra when the interparticle distance decreases. All these features are straightforwardly rationalized in the frame of the plasmon hybridization (PH) formalism applied to metal sphere dimers (D−PH model).31 For these complex composite systems the optical spectra are, in general, very intricate, and the theoretical analysis requires disentangling the specific features arising from the optical properties of the individual subunits from those rooted in the aggregate geometry and particle−particle couplings. Indeed, numerous powerful analytical and numerical methods are available for computing the absorption/scattering properties of such systems.32−36 Most of them consist, basically, in directly solving Maxwell’s equations, discretized in a large computaReceived: September 12, 2014 Revised: October 31, 2014 Published: November 20, 2014 28118

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

ment, may induce strong changes of the plasmon frequencies and, consequently, of the optical response that is expected for the isolated nanostructure. A neighboring high-index dielectric or metallic substrate has thus imperatively to be taken into account in order to reliably model and analyze the measured absorption/scattering spectra of supported plasmonic nanosystems. Numerous works have been devoted to this issue this past decade, emphasizing that the substrate effects are of main importance in many cases. A (nonexhaustive) list of relevant works on this subject, as well as an introductive discussion of the substrate effects and of the simple models generally assumed, can be found in our previous work ref 39. In this respect, it would be very valuable to extend the PH picture to handle nanoparticle/surface systems. The case of metallic surfaces has been theoretically investigated in several papers, and addressed experimentally, showing that the hybridization between the localized plasmons of the nanostructure and the propagating plasmons of the extended metal surface can produce either a red-shift or a blue-shift of the nanostructure plasmon frequencies, depending of the relative energies of the nanostructure and surface plasmons.40,41 To our knowledge, the extension of the PH approach in the case of a nearby dielectric substrate (real refractive index), which does not sustain intrinsic plasmon modes, was never reported on in the literature, although PH-related concepts have been invoked by several authors to discuss and analyze qualitatively experimental and/or theoretical results.39,42,43 Actually, as stated in the section 3.5 of our previous work,39 the PH model relevant to single particle/interface systems is expected to differ in many respects from the D−PH model. This work aims first at developing the PH formalism suitable for this kind of hybrid system. For simplicity, the PH model will be detailed in the case of the simplest generic system, namely, a jellium-sphere located above a dielectric interface of refractive index Nsub (referred to as the PI−PH model). The refractive indexes of both the embedding medium surrounding the sphere and the spherical ionic background are thus assumed to be equal to 1. Extensions of the PI−PH model to deal with more complex cases (several neighboring particles above the interface, other particle shapes, dielectric effects from the embedding medium and the underlying ionic backgrounds) can be derived in following the same lines as those involved for handling this very simple generic system, but obviously would result in much more complex formalisms. The second purpose of this work is to compare the predictions of the quasistatic D−PH and PI− PH models with exact multipole expansion method results. For the two-sphere dimer case, only a brief outline of the D−PH model has been provided in the pioneering paper by Nordlander et al.31 Since the PI−PH model involves similar ingredients, a very detailed description of the D−PH formalism is provided in the Supporting Information (atomic units are used throughout this work).

tional domain, in taking into account appropriate boundary conditions. However, except for the simplest systems, for instance dimers or linear chains of spheres, all these methods result in heavy and very time-consuming computations. Some of these approaches, especially the (partially) analytical Tmatrix approach33 (often referred to as the “generalized Mie theory” or “multipole expansion method” in the literature), lead to results of very high accuracy that provide invaluable benchmarks for assessing the accuracy of the numerical approaches. Nevertheless, all these methods generate raw optical spectra, which in most cases give no straightforward information about the “physical meaning” of the spectral response of plasmonic systems. The appealing and versatile PH framework developed by the Nordlander group31,37,38 fills this shortcoming and is now of current use for analyzing spectral responses of complex single or multiparticle nanosystems.3,12−14,17,20,22 In contrast with the previous approaches, which consist in solving the relevant electrodynamics equations (that is, solving the Maxwell equations in the presence of the external exciting electromagnetic field, by using appropriate dielectric functions for describing the polarization/absorption properties of the involved materials), the PH approach consists in searching the (global) excitation modes ωk of the conduction electron fluid filling the metallic subunits of the entire composite system. These excitations are characterized by specific small-amplitude surface deformations of the conduction electron fluid. Because the electron charge density is neutralized in each individual subunit by a rigid positively charged ionic background, the surface deformations of the electron fluid give rise to a net superficial charge density on the metallic surfaces. The surface deformation eigen modes, that is, the collective plasmon resonances, are determined by solving the equations of motion derived from the Lagrangian L = T − V of the entire composite system expressed in terms of the dynamical variables associated with the net surface charge density (T and V are the kinetic energy and electrostatic energy of the system, respectively). Actually, for each system, the corresponding PH modeling can be viewed basically as a mere mechanical/electrostatic problem. The coupling to an external electromagnetic field is then introduced by including in the Lagrangian the corresponding nanosystem/electromagnetic field interaction energy term. Absorption can be achieved, for instance, in introducing dissipative frictional terms in the equations of motion, proportional to the velocities (let us point out that this seemingly ad-hoc procedure amounts to adding an imaginary component to the real frequency ω in the analytical expression of the imaginary part of the dynamical polarizability, as usually done). The PH formalism provides a physically intuitive description of coupled plasmons in complex single or multiparticle nanostructures, in describing the plasmon excitations of the entire system in terms of interactions between the plasmon modes of the individual subunits. Since its introduction 10 years ago in the case of layered concentric nanoshells,37,38 further developments have extended its applicability to a large variety of plasmonic nanosystems,3,12,31 in particular, aggregates of neighboring metal particles, explaining the widespread and almost systematic use of the PH picture for analyzing plasmonic optical spectra. In fact, most experiments on plasmonic nanostructures involve a nearby dielectric, or sometimes metallic, supporting surface. This unavoidable experimental condition, which breaks the symmetry (and homogeneity) of the nanosystem environ-

2. PLASMON HYBRIDIZATION (PH) MODEL FOR A JELLIUM SPHERE DIMER 2.1. Lagrangian, Equations of Motion and EigenModes. In the PH formalism the plasmon excitations are selfsustained small-amplitude deformations of the conduction electron fluid, assumed incompressible and irrotational. Since the electron fluid is incompressible, and the discrete ionic background modeled by a rigid homogeneous positively charged jellium, the deformation of the electron fluid is reflected through a departure from the charge neutrality on the 28119

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

≥ 0; S i ≥ max(1,m)] are given by eq 26 in the Supporting Information. At this stage it is valuable to introduce linear combinations of the dynamical variables SS(,im) and SS(,i−) m, allowing the singlesphere Lagrangian L(i) = T(i) − V(i) self to be expressed directly as the sum of independent oscillators, the frequencies of which coincide with those predicted by the classical electrodynamics theory (frequencies ωS ). Moreover, the new dipolar modes (S = 1) will be specifically excited by electric fields linearly polarized along the three axes of the Cartesian reference frame, facilitating therefore the intuition and the discussion when the interaction with a plane wave linearly polarized is addressed. We introduce below the convenient independent dynamical variables that are used throughout this work

surface of the particle.31,37,38 A very detailed description of the PH model appropriate to jellium-sphere dimers (same metal) is reported on in the Supporting Information (section A). This report provides, on the one hand, explicit formula allowing straightforward numerical implementation, and on the other hand, theoretical extensions regarding the coupling of the dimer system to an external electromagnetic field. In this section, only the ingredients that are necessary for facilitating the reading, and essential for elaborating the PH model appropriate to the sphere/dielectric interface system, are duplicated. Figure 1 displays the dimer geometry and the coordinate systems used.

ZS(,0i)(t ) = SS(,0i)(t )

(m = 0)

(5)

1 [SS(,im) (t ) + ( −1)m SS(,i−) m(t )] 2

XS(,im) (t ) = −

(m > 0) (6)

1 i[SS(,im) (t ) − ( −1)m SS(,i−) m(t )] 2

Y S(,im) (t ) = − Figure 1. Geometry of the jellium-sphere dimer, notations and coordinate systems used; d is the interparticle distance.

(7)

The selected normalization constants (−1/√2 and −i/√2) will allow the three coordinates x, y and z to be treated on an equal footing when the interaction with a linearly polarized electric field will be addressed. In terms of these variables, the surface charge density and Lagrangian of each sphere, and the electrostatic interaction Vint write as (eqs 20, 21, and 25 in the Supporting Information)

The surface charge density on each sphere 1 and 2 (indicated by the superscript (i) with i = 1, 2) is expanded upon the complete set of scalar spherical harmonics (standard definition of these complex angular functions are assumed)44 ⎛ S ⎞1/2 σ (θi , ϕi , t ) = −n0e ∑ ∑ ⎜ 3 ⎟ SS(,im) (t )Y Sm(θi , ϕi) R S = 1 m =−S ⎝ i ⎠ S



1/2 ∞ ⎛ S ⎞ σ (i)(θi , ϕi , t ) = −n0e ∑ ⎜ 3 ⎟ f S0 (θi)ZS(,0i)(t ) R S=1 ⎝ i ⎠

(i)

(1)



where −n0e (e > 0) is the homogeneous charge density of the incompressible electron fluid and Ri is the radius of the sphere i. The deformation amplitudes SS(,im) (unit: length5/2) are the dynamical variables in the PH model (in eq 1, specific normalization constants have been introduced for convenience, as in ref 31). The kinetic energy (T(i)(t)) and electrostatic energy (V(i) self(t)) of each sphere, and the electrostatic interaction between the surface charges of both spheres (Vint(t)) can be easily expressed in a compact form in terms of these dynamical variables (eqs 7, 10, and 14 in the Supporting Information) n0me 2

T (i)(t ) =

(i) V self (t ) =

n0me 2

(m > 0)

m=1

+

Y S(,im) (t )

L(i) =

⎛ 2S ⎞1/2 ⎟ f Sm (θi){XS(,im) (t ) cos(mϕi) 3 R ⎝ i ⎠ S=m ∞

+ n 0e ∑

∑⎜

sin(mϕi)}

(8)

∞ n0me ⎧ ⎨∑ [ŻS(,0i)2 − ωS2ZS(,0i)2] 2 ⎩ S=1 ⎪





+







∑ ∑ [(XṠ (,im)2 + YṠ (,mi)2) − ωS2(XS(,im)2 + Y S(,im)2)]⎬





m=1 S=m

S



∑∑

(i)

(9)

(i)

( −1)m SṠ , m(t )SṠ , −m(t )

S = 1 m =−S

(2)

Vint =

S



∑∑

( −1)m ωS2SS(,im) (t )SS(,i−) m(t )

S = 1 m =−S ∞

nm Vint(t ) = 0 e ωP2 ∑ 2 S1= 1





(3)

+

L



∑ ∑ m = 1 S1, S2 = m

∑ ∑

V (S1, m ,

X (2) S2)[X S(1) 1, m S2, m

+

⎫ ⎪

Y S(1) Y (2) ]⎬ 1, m S2, m ⎪ ⎭

(10)

S2 = 1 m =−L

× V (S1, |m| , S2)( −1)m S S(1) (t )S S(2) (t ) 1, m 2, −m

⎧ ∞ n0me 2⎪ Z (2) ωP ⎨ ∑ V (S1, 0, S2)Z S(1) 1,0 S2,0 ⎪ 2 ⎩ S1, S2 = 1

Note that different symmetries correspond to the Z-, X-, and Yterms in the expansion of the surface density (eq 8). The functions f Sm (θ) entering the surface charge density and coupling integral V(S 1,m,S 2) expressions are related to the spherical harmonics through the relation Y Sm (θ,ϕ) = f Sm (θ)eimϕ (m ≥ 0).

(4)

In the above equations, me is the electron mass, ωS = ωP(S /(2S / + 1))1/2 is the single sphere plasmon frequency of multipolar order S (ωP is the bulk plasma frequency) and L = min(S 1,S 2). The dimensionless coupling integrals V(S 1,m,S 2) [m 28120

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

⎡ U (1) ⎤ S, m [Vk , m] = [M ]⎢ (2) ⎥ ⎢⎣ US ′ , m ⎥⎦

The Euler−Lagrange equations derived from the Lagrangian L = L(1) + L(2) − Vint lead to the following equations of motion for the Z, X, and Y dynamical variables of both spheres (1)

Z̈ S1,0 + ω S12Z S(1) + 1,0

1 2 ωP 2



with

∑ V (S1, 0, S2)Z S(2),0 = 0 2

S2 = 1

(m = 0; S1 ≥ 1)

⎡ U (1) ⎤ ⎢ S, m ⎥ = [U1,(1)m , U2,(1)m , ..., US(1), m , U1,(2)m , U2,(2)m , ..., US(2), m]T max max ⎢⎣ US(2) ⎥ ′,m ⎦

(11)

(16)



(2)

Z̈ S2,0 + ω S22 Z S(2) + 2,0

1 2 ωP ∑ V (S1, 0, S2)Z S(1) =0 1,0 2 S1= 1

(m = 0; S2 ≥ 1)

Ẍ S(1) + ω S12X S(1) + 1, m 1, m

The (uncoupled) equations of motion of the Vk,ms take the compact matrix form (12)

1 2 ωP 2

[Vk̈ , m] + [ωk2, m][Vk , m] = 0

where [ω2k,m] is a diagonal square matrix (of dimension 2S max) the diagonal coefficients of which (eigen-values of the matrix problem) are the squares of the eigen-mode frequencies of the dimer system. 2.2. Coupling to a Monochromatic Electromagnetic Field within the Quasistatic Approximation. In this specific case, the electric field is spatially homogeneous [E(t)] and the interaction can be derived from the scalar potential

∑ V (S1, m , S2)X S(2),m = 0 2

S2 = m

(13)

1 2 ωP 2

(17)



(m > 0; S1 ≥ m) Ẍ S(2) + ω S22 X S(2) + 2, m 2, m

(15)



∑ V (S1, m , S2)X S(1),m = 0 1

S1= m

(m > 0; S2 ≥ m)

(14)

V (r, t ) = −r·E(t ) = −xEx(t ) − yEy(t ) − zEz(t )

(i)

The equations of motion for the Y s (m > 0) are similar to eqs 13 and 14. As expected (axial symmetry of the dimer), these equations show that the hybridization mechanism, that is, the coupling of the plasmon modes of both spheres, is diagonal with respect to the m-value (m ≥ 0), and, for m > 0, that the X and Y variables are uncoupled and will correspond to identical plasmon modes (same eigen-values and eigen-vectors in the resulting matrix problem described in detail in the Supporting Information). The two-fold degeneracy of the m > 0 modes results from the equivalence of the x- and y-axes. As shown hereafter, in the presence of a homogeneous electromagnetic field, the m = 0 (variables Z) modes, the m = 1(variables X) modes, and the m = 1 (variables Y) modes are optically excited for, respectively, z-, x-, and y-polarizations. Each independent subset of the above infinite set of linear time-dependent differential equations can be handled through matrix techniques by prescribing a maximum S -value (S max), chosen large enough for ensuring convergence of the calculations and results in an eigenvalue problem. Solving the eigenvalue problems corresponding to each m-value allow the eigen-frequencies ωk,m and corresponding eigen-vectors Vk,m of the dimer to be determined (see the Supporting Information for details). For simplicity, only the values m = 0 and 1 are considered in the following. For m = 0 (variables Z) and m = 1 (variables X or Y), 2S max modes (i.e., independent onedimensional oscillators) are obtained. In particular, each eigen-vector set {Vk,m; k = 1, 2, ..., 2S max} can be chosen real and orthonormalized. The Vk,ms are thus related to the initial (2) dynamical variables {US(1) , m } and {U S , m } (U stands for Z if m = 0; ′ U stands for X or Y if m = 1, depending on the eigenvalue problem in hand) through a real orthogonal matrix M (square matrix of dimension 2S max, which depends on m) in which the k-th row consists of the expansion coefficients of Vk,m as a function of the 2S max initial dynamical variables, namely (hereafter [Vk,m], is a column vector containing the 2S max uncoupled dynamical variables)

(18)

In the presence of an external uniform electric field, the Lagrangian becomes (1) (2) L = L(1) + L(2) − Vint − V ext − V ext

(19)

with (omitting time dependence for clarity) ⎛ 4π ⎞1/2 (i) (i) (i) (i) = n0eR i3/2⎜ ⎟ [EzZ1,0 + ExX1,1 + EyY1,1 V ext ] ⎝ 3 ⎠ (i) (i) (i) = a(i)[EzZ1,0 + ExX1,1 + EyY1,1 ]

(20)

In the quasistatic approximation, only the dipolar modes of both spheres are directly coupled to the electromagnetic field, and therefore, only the equations of motion of the six (i) (i) dynamical variables X(i) 1,1, Y1,1, and Z1,0 (i = 1 and 2) are modified (eqs 58−61 in the Supporting Information). Solving the whole set of time-dependent differential equations is not straightforward because all other dynamical variables with m = 0 and 1 are coupled, through Vint, to one of the six dynamical variables directly coupled to the electromagnetic field. As shown below, this coupling will result in the transfer and sharing of the oscillator strengths of the S = 1 dipole excitations of the isolated spheres to all the m = 0 and m = 1 hybridized modes of the dimer system, most of them lying in remote spectral ranges. The smaller the intersphere distance d, the larger will be the oscillator strength transfer and the broadening of the excited plasmonic spectral range. The problem can be solved easily in expressing the equations of motion of the dynamical variables Vk, which remain uncoupled in the presence of the external field. Since the problem is diagonal with respect to the m-value and the symmetry involved (X- or Y-character for m > 0), we report on explicitly only the resulting equations corresponding to the m = 1 (variables X) case. All the other cases can be treated in an obvious similar way. One obtains (the label m is omitted in the following equations) 28121

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

⎡ a(1) ⎤ a(2) Vk̈ + ΓVk̇ + ωk2Vk = −⎢ Mk ,1 + Mk , Smax + 1⎥Ex(t ) n0me ⎣ n0me ⎦ = βk Ex(t )

only for small enough damping parameter Γ relative to the eigen-mode spacings. Strictly identical results are obtained for the m = 1 (variables Y) case, since the x- and y-axes are equivalent. The m = 0 (variables Z) case leads to the same formal equations, but the M-matrix coefficients that are involved are now different. 2.3. Comparison with “Exact” Calculations. In this section we report numerical tests demonstrating the relevance and the limits of the D−PH model for interpreting the optical absorption spectra of jellium sphere dimers, with respect to both the eigen-mode frequencies and the associated oscillator strengths. For this purpose, absorption cross sections computed in the frame of the D−PH model (eq 25) are compared against those computed using a home-developed computer code based on the exact multipole expansion method (MEM)33 that was previously applied to study the optical properties of gold dimers near the conductive contact limit (strong coupling regime).21 This method consists in expanding the involved electromagnetic fields (applied external field, scattered and internal fields of both spheres) in terms of vector spherical harmonics (VSH) and determining the unknown expansion coefficients in imposing the boundary conditions on the dielectric interfaces. The MEM is a very fast and accurate computational tool. The computational time depends on various parameters, in particular, the maximum S -order value selected for expanding the fields on the VSH. This time is typically on the order of a few seconds or a few tens of seconds for computing a spectrum involving 1000 wavelengths. A detailed account of this method is out the scope of this paper. It is sufficient for our purpose to recall that, on the condition that the z-axis coincides with the dimer axis (this axis will specify the azimuthal “quantum number” m of the VSH’s), the problem is diagonal with respect to the m-value. Consequently, the scattering, absorption and extinction cross sections are the sums of the various mcontributions for any incoming applied external electromagnetic field. Therefore, each m ≥ 0 mode frequency and associated oscillator strength sets in the D−PH model can be selectively investigated through the corresponding |m| contribution to the exact MEM absorption spectrum (sum of the m and −m contributions for m ≠ 0). For consistency, the MEM computations have been carried out in assuming a Drude dielectric function ε(ω) = 1 − ω2P/[ω(ω + iΓ)] for describing the dielectric properties of both metal spheres. Moreover, a small damping parameter ℏΓ (0.002 eV) has been selected in order to clearly separate the various hybridized plasmon excitation lines in the spectra. First, since retardation effects are disregarded in the electrostatic D−PH model, several numerical tests involving very tiny jellium sphere dimers have been carried out (see the note, ref 45). All these investigations have proved that the D− PH model is indeed a very accurate numerical tool in the quasistatic limit. Figure 2 illustrates this study in the case of a nanosized homodimer (radii R1 = R2 = 1 nm; interparticle distance d = 0.1 nm), for p-polarized light excitation at the oblique incidence θ = 45° with respect to the dimer axis [the electric field amplitude of the incoming wave Ei = E0 exp(ik·r − ωt) is (E0/√2) (ex + ez)]. The lower part of Figure 2 displays the m = 0 MEM absorption cross-section contribution. The D− PH model m = 0 absorption spectrum (upper part of Figure 2) is computed thanks to the analytical formula eq 25 divided by 2 for allowing a direct quantitative comparison with the MEM spectrum. Perfect agreement is obtained regarding the spectral location of the plasmon modes as well as the associated

(21)

where Mk,1 and Mk, Smax + 1 are the M-matrix coefficients of the Vk expansion on the two dynamical variables (m = 1; X-character) that are directly coupled to the external field, namely, X(1) 1,1 and X(2) 1,1 . In the above equation a “friction force” has been added, as in the Drude-Sommerfeld dielectric function model. It should be pointed out that the time-averaged light absorption is equal to zero if no oscillator damping is present. This ingredient, which can be suitably disregarded for determining the dimer eigen-modes, has necessarily to be introduced for describing the energy exchanges between the dimer and the applied electromagnetic field. For convenience, Γ is assumed modeindependent. From the steady-state solution of eq 21, the total induced dipole moment can be easily determined 2Smax

Dx (t ) =

∑ βk′ k=1

1 Ex(t ) = α(ω)Ex(t ) (ωk2 − ω 2) − iω Γ (22)

with βk′

⎫ R3 R 3/2 1 3 2⎧ 2 R1 ωP ⎨Mk ,1 + 23 Mk2, Smax + 1 + 2 23/2 Mk ,1Mk , Smax + 1⎬ 3 R1 R1 ⎭ ⎩ 1 3 2 = R1 ωP Ck′ (23) 3 =









The dimensionless coefficients C′k will be referred to as the “normalized oscillator strengths” in the manuscript. Since the matrix M is orthogonal, we have 2Smax

∑ βk′ = k=1

1 2 3 ωP [R1 + R 23] 3

(24)

The physical meaning of eqs 21−23 is obvious. When the two spheres are strongly coupled (small interparticle distances, more precisely small d/Ri values), the modes (S i = 1, 2, 3, ...; m = 1; X-character) of the individual spheres hybridize. In particular, each resulting mode k of the dimer system acquires a dipolar component and therefore will be coupled to the external monochromatic electromagnetic field when ω ≈ ωk. Actually, eq 24 is nothing else but the sum-rule describing the transfer (and sharing) of the oscillator strengths of the x-axis related dipole excitations of the two isolated spheres to the 2S max eigen-modes of the dimer system (modes with m = 1; Xcharacter). Since these eigen-modes are located in distinct spectral ranges, the absorption spectra will exhibit a broad plasmonic spectral pattern in the strong coupling regime. The absorption cross-section is given by σ(ω) = =

4πω Im[α(ω)] c 4πω c

2Smax

∑ βk′ k=1

ωΓ (ωk2

2

− ω 2)2 + ω Γ 2

(25)

In the case of a dimer, for small enough intersphere separations, a multipeak plasmonic structure will be observed in spectra recorded by varying the frequency of the applied external field. Obviously, a conspicuous multipeak pattern will be exhibited 28122

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

frequencies and eigen-vectors resulting from the free equations of motion (eqs 15 and 17). In addition, these computations show that involving the parameters (ωP, Γ) in the “mechanical/ electrostatic” D−PH model amounts to assuming a DrudeSommerfeld dielectric function with identical parameters for describing the optical properties of the conduction electrons in more conventional electrodynamics methods. For completeness, one has to mention that an analytical electrostatic solution of the Laplace equation can be obtained for two-sphere dimers by separation of variables in bispherical coordinates from which the dimer plasmon eigenmodes can be computed.46,47 Although the size of the dimer is very small, retardation effects are evidenced in the strongly magnified MEM spectrum (the high energy part of Figure 2 is reported in Figure 3). The

Figure 2. (a) m = 0 contribution to the absorption cross-section σabs of a jellium-sphere dimer (radii R1 = R2 = 1 nm; interparticle distance d = 0.1 nm) in vacuum, for p-polarized light excitation at the oblique incidence θ = 45° with respect to the dimer axis, computed within the D−PH model (eq 25). A small damping parameter ℏΓ (0.002 eV) has been selected in order to separate the various hybridized plasmon excitation lines. The electron density is characterized by rs ≈ 3.015 bohr (electron Wigner Seitz radius) and corresponds to a bulk plasma energy ℏωP equal to 9 eV (the plasma frequency corresponds to bulk silver and gold, but it should be emphasized that the spectra cannot be compared with experiment since the screening by the bound core electrons is disregarded throughout this work). In order to distinguish the high-order modes (m = 0, large k) that are weakly hybridized with the dipolar mode (m = 0, k = 1), the high energy part of the spectrum has been magnified (blue curve: σabs × 100 + 50; red curve: σabs × 4000 + 100; green curve: σabs × 16000 + 150). (b) m = 0 contribution to the absorption cross-section σabs computed within the exact multipole expansion method33 (identical parameters are involved). In this calculation, a Drude-Sommerfeld parametrization of the dielectric function of both jellium spheres has been taken (ℏωP = 9 eV and ℏΓ = 0.002 eV).

Figure 3. High energy part of Figure 2, displaying the emergence of the dark modes (indicated by arrows in the lower panel (b)) in the MEM spectrum (the retardation effects are taken into account). For homodimers the oscillator strengths of these modes are strictly vanishing in the D−PH model (upper panel (a)).

oscillator strengths. For the axial polarization we can see that, due to the sphere−sphere electrostatic coupling, a large part of the oscillator strength of the dipole excitations of both isolated spheres has been transferred to the bright quadrupole and octupolar modes (peaks labeled k = 2 and 3 in the spectra), as expected from the Figure 2(SI) in the Supporting Information. The high energy part of the spectra has been magnified in order to make visible to the eye high-order modes (m = 0, large k) that are weakly hybridized with the dipolar excitations (m = 0, k ≡ S = 1). The characteristics of these high-order dimer plasmon excitations are found also in perfect agreement, proving the accuracy of the D−PH model, despite both formalisms are of different nature, as emphasized in the Introduction. Actually all the features observed in the absorption spectra reflect intrinsic properties of the free jellium-sphere dimer (ie the dimer in the absence of any source of excitation) since the spectral locations and intensities of the excitation lines depend only on the eigen-

small bumps indicated by arrows in Figure 3 are the dark modes of the D−PH model, the oscillator strengths of which are strictly vanishing in this electrostatic approach. Indeed, the very high accuracy of the MEM computations allows revealing all the eigen-modes predicted by the D−PH model, permitting thus a complete checking of this appealing approach. In Figure 4 are plotted the D−PH and MEM spectra corresponding to the transverse ex polarization (m = 1 related eigen-modes are evidenced). As compared to the axial polarization ez (Figure 2), only the dipolar bright mode k = 1 of the dimer, the frequency of which is close to that of the isolated spheres (ωS =1 = ωP/√3) is noticeably excited, as expected from the Figure 2(SI) in the Supporting Information. As emphasized at the end of the Supporting Information, this is due to the weakness of the intersphere electrostatic coupling 28123

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

Figure 4. Same as Figure 2 for the |m| = 1 contribution (blue curves: σabs × 20 + 100; red curves: σabs × 80 + 350).

when the electron fluids oscillate along the transverse directions. In contrast to the m = 0 case, a lower magnification is sufficient for making visible to the eye several high-order modes. This is consistent with the fact that the m = 1 bright branches and, concomitantly, the curves describing the oscillator strength transfer, are flat, whatever the interparticle distance is (see Figures 1(SI) and 2(SI) in the Supporting Information). For a transverse polarization, the sharing of the dipole oscillator strengths of the isolated spheres among the various eigen-modes of the dimer is more homogeneous and does not depend very strongly on the distance d because of the weakness of the electrostatic coupling. Second, we have carried out calculations involving larger homodimers (R1 = R2 = R) in order to scrutinize the range of applicability of the electrostatic D−PH model. For addressing this objective it is sufficient to test the scaling properties of the D−PH model. In this approach the plasmon eigen-modes and eigen-vectors are unchanged if the geometrical parameters, namely, R1, R2, and d, are multiplied by a common scaling factor (see the Supporting Information) and, on the other hand, the absorption cross-section scales as the particle volume (see eqs 23 and 25). So, keeping constant the ratio d/R, we have computed a series of MEM extinction (sum of the energy absorption and scattering rates) cross-section spectra for increasing radius R. By way of illustration, we show in Figure 5 MEM spectra for homodimers with d/R = 0.2 and R ranging from 0.1 nm (lowermost spectrum) up to 15 nm (uppermost spectrum). Except for the length parameters, the same ingredients as those involved previously (Γ, ωP, irradiation conditions) have been assumed. In this figure the various plasmon excitation lines, as well as the shoulders in the broad bands observed for the largest radii, have been labeled by their associated |m| character.

Figure 5. Extinction cross sections σext computed using the multipole expansion method33 for jellium-sphere dimers, for p-polarized light excitation at the oblique incidence θ = 45° with respect to the dimer axis. The radius R of both spheres is ranging from 0.1 to 15 nm, in keeping constant the ratio d/R, where d is the sphere−sphere separation (d/R = 0.2). In these calculations, a Drude−Sommerfeld parametrization of the dielectric function of the jellium spheres has been taken (ℏωP = 9 eV and ℏΓ = 0.002 eV). The |m|-values corresponding to the observed plasmon excitation lines are given in the four uppermost spectra (R = 2.5 nm → 15 nm). Due to the axial symmetry of the physical system, the mode assignment can be determined straightforwardly in plotting the various |m| contributions to the overall extinction. The label “d” indicates that the line corresponds to a D−PH model “dark mode” in the limit of vanishing R values.

In Figure 5 the label “d” means that the peak corresponds to a dark mode in the electrostatic D−PH model. For the three smallest R values, the extinction is exhausted by the absorption and the spectra scale as R3, as predicted by the D−PH model. For R ≥ 5 nm several additional features can be observed. First, very small peaks, corresponding to dark modes in the D−PH model, begin to emerge. Second, a red-shift of the plasmon modes becomes noticeable. Basically, both features result from retardation effects and the concomitant size-dependence of the dynamic depolarization fields48,49 inside each sphere generated by the induced surface charge densities (the overall size of the system, i.e., 4R + d, is involved for quantifying the retardation effects in the strong coupling regime). Third, the various lines, more especially the dipolar ones, broaden due to radiation damping,48,49 and concomitantly, the scattering loss channel prevails over absorption. All these features, not accounted for by the electrostatic D−PH model, are magnified in the two uppermost spectra. Actually, for large sizes, the line widths are 28124

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

It should be emphasized that a very small damping parameter Γ has been taken for convenience in order to spectrally separate the various modes, explaining the strong prevalence of the scattering for the largest radii investigated. This small value is responsible for the very striking increase of the line widths with increasing radius, trend which is usually unexpected in this size range. Calculations involving more realistic Γ-values (roughly a few tens of meV) indicate that the absorption remains an efficient loss channel in the size range studied but confirm that the red-shifts of the modes do not depend on the selected Γvalue (see the Figure 4(SI) in the Supporting Information). This is consistent with the fact that the spectral shifts induced by the dynamic depolarization field depend only on the dimer size-to-wavelength ratio. Larger Γ-values result essentially in a smoothing off of the spectra (this affects essentially the very sharp structures related to D−PH model dark and high-order modes) and a change of the respective scattering and absorption contributions to the extinction cross-section (see the Figure 4(SI) in the Supporting Information). Moreover, in plotting separately the various |m| contributions to the spectra, I have checked that the sharing of the oscillator strength among the various plasmon modes do not depend on the Γ-value. All these calculations point out that the spectral mixing of modes of different S -values in narrow spectral ranges and their eventual reordering are general features in the strong coupling regime for large dimer sizes. One can also argue that, since the dielectric/screening effects stemming from both the embedding medium and ionic backgrounds have been disregarded, short wavelengths correspond to the plasmonic frequency range (near UV), resulting in a restricted range of validity of the quasistatic approximation. However, in including dielectric/ screening effects, the plasmon modes, as well as their energy spacings, are lowered, stressing that the mode mixing and reordering occur probably in this case, too, as soon as the dimer sizes exceed a few tens of nm. In conclusion, a thorough analysis is required to interpret and reliably label spectral structures in experimental spectra (extinction or scattering, depending on the experiment) beyond the quasistatic approximation. Throughout this work, all the illustrating examples (on both the D−PH and PI−PH models) involve dimers consisting of metal jellium spheres characterized by the plasma energy ℏωP = 9 eV, which corresponds to the conduction electron density of bulk silver and gold. Nevertheless, it should be emphasized that all the conclusions drawn previously for dimers, and hereafter in section 3.3, for the particle/interface system, concerning the range of applicability and limits of both models, apply also to metals characterized by, or a lower density, as sodium (rs ≈ 4a0), either a larger density, as aluminum (rs ≈ 2a0). However, as reminded previously, the dielectric/screening effects by the core electrons (related to the interband transitions) have been disregarded throughout this work, and so the displayed spectra cannot be compared with experimental spectra. Introduction of the core−electron contribution to the overall metal dielectric function is indeed necessary for predicting the correct plasmon energies and, in particular, explaining the differences observed in the optical properties of silver and gold. It should be stressed that the taking into account of the dielectric/screening effects from core-electrons in the D−PH and PI−PH models results in a much more complex formalism. These effects have been disregarded for convenience in this introductory work (a forthcoming paper, involving screening effects from the core electrons and the embedding matrix, is in preparation).

dominated by the radiation damping contribution [compare, for instance, the lowermost (line width = 0.002 eV) and uppermost spectra]. One can notice that the line widths of the dark and k ≡ S > 1 modes remain narrow as compared to those of the dipolar modes. This feature can be interpreted in terms of quasi-destructive interferences between several radiating dipoles (inducing thus weak losses), since, for these two cases, the overall electron fluid in the dimer consists of successive regions oscillating in opposition of phase (the sizes of each sphere and of the dimer are much lower than the wavelength). The radiation damping-induced broadening results in a strong slowing down of the oscillator strength R3-scaling law for increasing radius. In fact, when the radius is large, the spectra may become very intricate (see, for instance, the uppermost spectrum corresponding to R = 15 nm) and the labeling of the observed spectral structures requires to plot separately each |m| contribution to the MEM spectrum (see Figure 6, where the

Figure 6. m = 0 (black curve), |m| = 1 (blue curve) and |m| = 2 (red curve) contributions to the uppermost extinction cross-section MEM spectrum plotted in Figure 5 (R = 15 nm, d = 3 nm). For convenience, the |m| = 1 and |m| = 2 contributions have been vertically shifted (shift equal to 30000 nm2 and 45000 nm2, respectively). In the |m| = 0 and | m| = 1 spectra, the arrows indicate collective excitations that are related to D−PH model “dark modes” in the limit of vanishing R values [in this limit, these modes are located at 5.457 eV (m = 0), 4.76 eV (m = 1), and 5.357 eV (m = 1)].

arrows indicate dark modes of the D−PH model). This thorough analysis points out that several plasmon modes can contribute noticeably in a narrow spectral range (especially the range 4.9 eV ↔ 5.4 eV in the present example, in which a large part of the oscillator strength is concentrated). In all studied cases, only the lowest energy peak that corresponds to the bright m = 0 mode, does not suffer from this mode mixing on the frequency scale. In addition, because of the mode- and sizedependences of the red-shift induced by retardation effects, the order of the modes may change. For instance, the quadrupolar (k ≡ S = 2) m = 0 mode near 5.3 eV is red-shifted relative to the dipolar |m| = 1 mode for R ≤ 5 nm, becomes degenerate with the latter for R = 10 nm (see Figure 5) and the order of both modes is inverted for R = 15 nm (see Figure 6). 28125

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

Moreover, only homodimers have been considered because the main goal was to extend, in the framework of a very simple generic system, the PH approach to the particle−interface system (the radii of the jellium sphere and of its image are identical). For completeness, some relevant comments about heterodimers involving slightly different radii are reported on hereafter. 2.4. Additional Comments about Heterodimers. In the case of homodimers, the nanosystem is symmetric with respect to the transverse x−y plane cutting the z-axis at equal distance from the sphere centers (see Figure 1). In such symmetric systems, the decreasing/bonding and increasing/antibonding plasmon energy branches evolve steadily when the interparticle distance decreases and cross at specific separations (see Figure 1(SI) in the Supporting Information). The exact degeneracy of the plasmon modes at the crossings is actually a direct signature of the above-mentioned exact symmetry. For heterodimers (R1 ≠ R2), the symmetry is broken, resulting in the repelling of nearby plasmon frequency branches and avoided crossings are exhibited in the overall energy branch pattern.31 These features, which are directly linked to the symmetry breaking, are quite general phenomena and are actually met in many physical areas. By way of illustration, two figures are provided in the Supporting Information [Figures 5(SI) and 6(SI)], in the case of a heterodimer involving slightly different radii (R1 = 10 nm and R2 = 9 nm). The first Figure shows the distance dependence of the lowest plasmon energy branches for m = 0 and m = 1, which display conspicuous avoided crossings when two branches approach each other. The second one shows the striking evolution of the associated normalized oscillator strengths (defined in eq 23) at the avoided crossings. Indeed, when two branches k and k′ approach and repel from each other, their physical characteristics are exchanged. More precisely, the expansion coefficients Mk,S and Mk ′ ,S of their corresponding eigen-vectors Vk and Vk′ (see eq 15) are, to some extent, exchanged, in particular, the bright and dark characters (expansion coefficients with S = 1). Indeed, the crossing branches in the case of homodimers and the noncrossing branches in the case of heterodimers can be paralleled with the diabatic and adiabatic potential energy curves of a diatomic molecule plotted as a function of the internucleus distance, respectively. The “interaction” between the dipolar m = 1 antibonding branch and the quadrupolar m = 1 bonding branch is chosen for exemplifying the exchange of the bright and dark characters at the avoided crossings (located in this case near d = 2 nm; see Figure 5(SI) in the Supporting Information). In Figure 7 are plotted m = 1 absorption spectra (electric field polarized along the x-axis) computed using eq 25 for different interparticle separations. For large separations, these two modes are bright and dark, respectively. When the distance d is decreased, the bright and dark characters are steadily exchanged between the two modes, in keeping unchanged the total oscillator strength. Keeping in mind this feature is of particular importance since experiments involve, in general, heterodimers. This figure points out that, in addition to the particle radii, a perfect knowledge of the interparticle separation is required for analyzing quantitatively and labeling reliably the absorption spectrum. Moreover, for realistic damping rates, Γ, the spectral overlap of both modes involved at each avoided crossing may result in a broad and unexpected line width, as well as a pronounced asymmetry of the line shape.

Figure 7. Absorption spectra illustrating the exchange of the bright/ dark character at the avoided crossings between two plasmon frequency branches of heterodimers. The present example corresponds to an heterodimer in air (refractive index Nm = 1) involving slightly different sphere radii (R1 = 10 nm and R2 = 9 nm). The electron density is characterized by rs ≈ 3.015 bohr (electron Wigner Seitz radius) and corresponds to a bulk plasma energy ℏωP equal to 9 eV. This specific example correspond to the avoided crossing between the quadrupolar bonding transverse mode (m = 1) and the dipolar antibonding transverse mode (m = 1) for intersphere distances d close to 2 nm. These two modes correspond to the red and blue frequency branches plotted in the Figure 5(SI) [lower panel (b)] in the Supporting Information. When d is decreased the bright character of the antibonding dipolar mode is transferred to the bonding quadrupolar mode, which has a dark character for large distances (see Figure 6(SI) in the Supporting Information). The computations have been carried out within the D−PH model using eq 25. For clarity, the successive spectra have been shifted along the vertical coordinate.

3. PLASMON HYBRIDIZATION MODEL FOR A JELLIUM SPHERE ABOVE A DIELECTRIC INTERFACE 3.1. Lagrangian, Equations of Motion, and EigenModes. The extension of the D−PH model described in the previous section to the system consisting of a jellium sphere in vacuum above a dielectric interface of real index Nsub is straightforward. For simplicity, the index Nsub is assumed to be frequency-independent in the spectral range of interest. As compared to a plane metallic surface, a dielectric nonabsorbing medium contains only tightly bound core electrons and so does not involve self-sustained electron fluid excitations. Actually, in the absence of an external applied field, the dielectric medium acts on the jellium sphere only through screening effects, which are generated by the polarization charge distribution on the plane interface that is directly induced by the deformations of the electron fluid in the jellium sphere. This preliminary comment points out that no extra dynamical variables have to be introduced for characterizing the dielectric medium since the dynamics of the surface polarization charge distribution is entirely governed by that of the electron fluid in the jellium sphere. Indeed, only the equations of motion of the dynamical variables associated with the electron fluid deformations in the jellium sphere are required for elaborating this specific PH model. In particular, providing the explicit expression of the Lagrangian associated with the dielectric medium is unnecessary for determining the eigen-modes of the particle/interface system.50 Besides, in the presence of an incoming electromagnetic plane wave, the effects on the jellium sphere of the 28126

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

additional incoming field-induced polarization charge distribution on the plane interface can be taken into account merely in just adding the reflected electromagnetic field to the incoming one for defining the effective homogeneous electric field E(t) interacting with the jellium sphere dipole (linearity of the Maxwell equations). Figure 8 displays the jellium sphere/

Vint(t ) =

1 2

)σ ′(r′, t ) 2 2 d S d S′ ∫ σ(r, |tMM ′|

(28)

In fact, the total electrostatic energy, Vtot = Vself + Vint, could be deduced directly from the quite general equation of electrostatics (suitable on the condition that the underlying dielectric media are linear)51 1 2

Vtot =

∫ ρ(r)V (r)d3r

(29)

where ρ(r) is the free charge density and V(r) is the total electrostatic potential (generated by both the free and polarization charges). One has now to express the image density σ′(r′,t) as a function of the dynamical variables SS, m(t) describing the surface charge density σ(r,t) on the jellium sphere surface. One has σ ′(r′ ≡ (R , θ′, ϕ′), t ) ∞

= − n 0e ∑

Figure 8. Geometry of the jellium sphere/dielectric interface system, notations and coordinate systems used. The vertical red line (at equal distance between the two spheres) symbolizes the plane interface; 2d is the distance between the jellium sphere (left sphere) and its image (right sphere).

S



S = 1 m =−S

1/2

⎛ S ⎞ ⎜ ⎟ ⎝ R3 ⎠

SS′, m(t )Y Sm(θ′, ϕ′)

(30)

Since Y Sm (π − θ, ϕ) = ( − 1)S + m Y Sm(θ,ϕ), one obtains immediately from eq 27 SS′, m(t ) = −ξ( −1)S + m SS, m(t )

interface system geometry and the relevant coordinate systems used (note that d is defined as the distance between the jellium sphere and the interface). The surface charge density in the jellium sphere, that is σ(r), is expanded according to eq 1. In the absence of an applied external field, the Lagrangian of the jellium sphere/interface system writes as (retardation effects neglected, as assumed in the jellium sphere dimer case) L = T − Vself − Vint

After introducing the convenient dynamical variables eqs 5−7, the explicit expression of the Lagrangian is ∞ n0me ⎧ 2 ⎨ [ZṠ ,0 − ωS2ZS2,0] ∑ 2 ⎪ ⎩ S=1 ⎪

L=



+







∑ ∑ [(XṠ 2, m + YṠ 2, m) − ωS2(XS2, m + YS2, m)]⎬ − Vint ⎪



m=1 S=m

(26)

(32)

In terms of the dynamical variables SS, m(t), the kinetic energy and the electrostatic self-energy of the electron fluid in the jellium sphere, T and Vself, respectively, are given by eqs 2 and 3. Vint is the electrostatic interaction of the free charge density σ(r) with the polarization charge distribution that is induced by σ(r,t) on the plane interface. Actually, in the quasistatic limit, the sphere/dielectric medium interaction Vint can be computed, thanks to the method of images,51 as the interaction of σ(r) with a spherical charge density σ′(r′,t) located below the plane interface. The image charge density is related to σ(r,t) through the equation

with Vint =

1 n0me 2 ω P ( − ξ) 2 2

∞ ⎧ ⎫ ⎪ ⎪ × ⎨ ∑ V (S1, 0, S2)( −1)S1 Z S1,0Z S2,0 + M ⎬ ⎪ ⎪ ⎭ ⎩ S1, S2 = 1 ∞

M=

(33)



∑ ∑

V (S1, m , S2)( −1)S1+ m [X S1, mX S2, m + YS1, mYS2, m]

m = 1 S1, S2 = m

(34)

σ ′(r′ ≡ (R , θ′, ϕ′), t ) = −ξσ(r ≡ (R , θ = π − θ′, ϕ = ϕ′), t )

(31)

The dimensionless coupling integrals V(S 1,m,S 2) [m ≥ 0; S i ≥ max(1,m)] are the same as those involved in the D−PH model (with R1 = R2; eq 26 in the Supporting Information). Since the particle and its image have identical radii, the coupling integrals obey the following relationship

(27)

with ξ = (N2sub − 1)/(N2sub + 1). The origin O′ of the translated Cartesian reference frame in which is defined σ′(r′,t) is the mirror image of O (center of the jellium sphere) with respect to the plane interface (see Figure 8). As emphasized previously, the mirror density does not correspond to an additional, that is, independent, charge distribution (no additional dynamical variables are associated with this density because of the constraint eq 27). Moreover, since the building up of the mirror density adiabatically follows the building up of the real charge density from zero up to its actual value σ(r,t), the interaction energy expresses as (compare with eq 12 in the Supporting Information)

V (S2, m , S1) = ( −1)S1+ S2 V (S1, m , S2)

(35)

The R-independence of the plasmon modes ωS of the isolated sphere and the scale invariance of the coupling integrals V(S 1,m,S 2) point out that the plasmon frequencies of the particle/interface system will be unchanged if the geometrical parameters, that is R and d, are multiplied by a common scaling factor. The equations of motion deduced from the Lagrangian eq 32 express as 28127

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C ZS̈ ,0 + ωS2ZS,0 +

( − ξ) 2 ωP 2



=0

the squares of the eigen-mode frequencies of the particle/ dielectric interface system. As for the two-sphere dimer system, it should be pointed out that the electrostatic problem for the particle/interface system has been solved in a pioneering paper by Ruppin, by separation of variables in bispherical coordinates.52 During the course of this work we have checked that this alternative approach leads to the same plasmon frequencies as the present PH model. 3.2. Coupling to a Homogeneous Monochromatic Electromagnetic Field. Within the quasistatic approximation, the external electric field is assumed to be spatially homogeneous [E(t)] over the volume of the jellium sphere. The coupling of the jellium sphere to the external field writes as (the time dependence is omitted for lightening the equations; see eq 20)

∑ (−1)S′V (S′, 0, S)ZS′ ,0 = 0 S ′= 1

(m = 0; S ≥ 1) XS̈ , m + ωS2XS, m +

Article

(36)

( − ξ) 2 ωP 2





( −1)S ′+ m V (S′, m , S)XS ′ , m

S ′= m

(m ≥ 1; S ≥ m)

(37)

The equation of motion for the YS, ms (m ≥ 1) is similar to eq 37. The construction of the PI−PH model follows the same general lines as those involved in the dimer case. The above equations eqs 36 and 37 show that the hybridization mechanism, that is, the coupling of the plasmon modes, is diagonal with respect to the m-value (m ≥ 0), and, for m > 0, that the Xs and Ys variables are uncoupled and will correspond to identical plasmon modes (axial symmetry and equivalence of the x- and y-axes). As in the dimer case, the m = 0 (variables Z) modes, the m > 0 (variables X) modes and the m > 0 (variables Y) modes will be optically excited for z-, x-, and y-polarized homogeneous electromagnetic field, respectively. Each independent subset of the above infinite set of linear timedependent differential equations can be handled through matrix techniques in prescribing a maximum S -value (Smax ), chosen large enough for ensuring convergence of the calculations. Note that, for a given Smax value, the dimensions of the column vectors and of the matrix C (defined below) are two times smaller as compared to the dimer case (eq 34 in the Supporting Information). By way of illustration, we detail explicitly only the specific matrix problem m = 1 (x-axis; variables XS,1). The equations of motion are cast in the following matrix form [XS̈ ,1] + [C ][XS,1] = 0

Vext = a[EzZ1,0 + ExX1,1 + EyY1,1]

with

( − ξ) 2 ωP ( −1)S ′+ 1V (S′, 1, S) 2

(38)

(39)

The matrix C is real and symmetric (this property is ensured thanks to eq 35) and can thus be diagonalized thanks to standard numerical algorithms. In particular, the S max eigenvalues ωk2 (k = 1, 2, ..., S max) are real and the eigenvector-set {Vk;k = 1, 2, ..., S max} can be chosen real and orthonormalized. The Vks are thus related to the initial dynamical variables XS,1 through a real orthogonal matrix M (square matrix of dimension S max) in which the kth row consists of the expansion coefficients of Vk as a function of the S max initial dynamical variables, namely, ([Vk] is a column vector containing the S max uncoupled dynamical variables) [Vk] = [M ][XS,1]

⎡ a ⎤ Vk̈ + ΓVk̇ + ωk2Vk = −⎢ Mk ,1⎥Ex(t ) = βk Ex(t ) ⎣ n0me ⎦

(43)

The steady-state solution of this equation is (with E(t) = E0e−iωtex = Ex(t)ex) Vk(t ) = βk

(40)

(ωk2

1 Ex(t ) − ω 2) − iω Γ

(44)

The x-component of the induced dipole moment D is

The (uncoupled) equations of motion of the Vks take the compact matrix form [Vk̈ ] + [M ][C ][M ]−1 [Vk] = [Vk̈ ] + [ωk2][Vk] = 0

⎛ 4π ⎞1/2 a = n0eR3/2⎜ ⎟ ⎝ 3 ⎠

In the presence of the external field, the Lagrangian is L = T − Vself − Vint − Vext. Only the equations of motion of the three dynamical variables X1,1, Y1,1, and Z1,0 are modified since only the dipole modes (S = 1; m = 0, 1) are directly coupled to the external electric field (see eq 42). As in the dimer case, the dynamics of all other dynamical variables with m = 0 and m = 1 is, however, modified, because these variables are coupled, through Vint, to one of the three dipole variables directly coupled to the external electric field. The coupling will result in the transfer and sharing of the oscillator strengths of the S = 1 dipole excitations of the isolated sphere to higher-order (k ≡ S > 1) plasmon excitations of the particle/dielectric interface system lying in remote spectral ranges. The smaller the intersphere distance, 2d, the larger will be the oscillator strength transfer and the broadening of the excited plasmonic spectral range. As in the dimer case, the problem can be straightforwardly solved thanks to matrix techniques in expressing the equations of motion of the dynamical variables Vk, which remain uncoupled in the presence of the external field. Since the problem is diagonal with respect to the m-value and the symmetry involved (X or Y character for m > 0), we report on explicitly only the resulting equations of the m = 1 (variables X) problem (counterparts of eqs 21−25). The equation of motion of each mode Vk(t) expresses as

The column vector (Smax coefficients) consists of the set of the Smax coupled dynamical variables XS, m = 1, that is, [X1,1,X2,1, ..., X Smax,1]T. C is a square matrix of dimension S max, the elements of which are CS, S ′ = ωS2δS, S ′ +

(42)

Smax

Dx (t ) =

(41)

∑ βk′ k=1

1 Ex(t ) = α(ω)Ex(t ) (ωk2 − ω 2) − iω Γ (45)

[ω2k ]

where is a diagonal matrix (dimension S max), the diagonal coefficients of which (eigen-values of the matrix problem) are

with 28128

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C βk′ =

Article

1 3 2 2 1 R ωP Mk ,1 = R3ωP2Ck′ 3 3

(46)

Since the matrix M is orthogonal, we have Smax

∑ βk′ = k=1

1 3 2 R ωP 3

(47)

The absorption cross-section is given by σ(ω) = =

4πω Im[α(ω)] c 4πω c

Smax

∑ βk′ k=1

ωΓ (ωk2

2

− ω 2)2 + ω Γ 2

(48)

The physical meaning implied from eqs 46−48 is similar to that previously discussed in section 2.2 (eqs 23−25). Through the electrostatic coupling Vint each mode k of the jellium sphere/ dielectric interface system acquires a dipolar component and therefore will be coupled to the external field. In the present specific case, this results in the transfer/sharing of the oscillator strength of the x-axis related dipole mode of the isolated jellium sphere to the eigen-modes of the jellium sphere/dielectric interface system (modes with k ≡ S ≥ 1; m = 1; X-character). Since these modes are located in different spectral ranges, the absorption spectra will exhibit a broad plasmonic spectral pattern in the strong coupling regime (large index Nsub and small d/R ratio). Strictly identical formulas are obtained for the m = 1 (variables Y) case, since the x- and y-axes are equivalent. The m = 0 (variables Z) case leads to the same formal equations, but the coefficients of the M-matrix that are involved are different. Figures 9 and 10 illustrate the PI−PH model for two values of the dielectric index Nsub. For the sake of direct comparison with the dimer case, identical sphere radius R and plasma frequency ωP have been taken for computing the eigenmodes (Figure 9) and the normalized oscillator strengths C′k (defined in eq 46; Figure 10) as a function of the jellium sphere-image sphere distance 2d. In both figures the two uppermost panels correspond to the value Nsub = 2.4, the two lowermost ones to the value Nsub = 1.5. As emphasized in our previous work,39 the major difference with the D−PH model is rooted in the fact that the dynamics of the image charge density is enslaved to that of the jellium sphere (see eq 27), resulting in the concomitant reduction by a factor 2 of the number of dynamical variables. In consequence, for a given m-value and symmetry, each mode of the isolated jellium sphere gives rise to a single bonding (i.e., red-shifted) bright frequency branch. As expected, the larger the medium index Nsub and R/d ratio, the stronger will be the hybridization mechanism, regarding both the redshifts of the eigen-frequencies and the oscillator strength transfers. Indeed, when Nsub approaches the dielectric index of the embedding medium (ξ ≈ 0), the eigen-mode frequencies tend toward those of the isolated embedded sphere. Comparison between Figure 9 and the Figure 1(SI) in the Supporting Information indicates that, for high dielectric indexes Nsub, the redshifts of the bonding modes for both systems are approximately of the same order of magnitude. As expected, the distance-dependence of the m = 0 modes is found more pronounced than that of the m = 1 modes one, because of stronger electrostatic coupling (see the Figure 3(SI) in the Supporting Information for guiding the intuition). By contrast to the dimer case, the transfer of the oscillator strengths

Figure 9. Evolution of the plasmon energies of a jellium sphere of radius R = 10 nm located above a dielectric interface (refractive index Nsub) as a function of the sphere−interface distance d (2d is the separation between the jellium sphere and its image). The electron density is characterized by rs ≈ 3.015 bohr (electron Wigner Seitz radius) and corresponds to a bulk plasma energy ℏωP equal to 9 eV. For very large distances, the plasmon energies tend toward the multipolar excitation energies of the isolated sphere, namely, ℏωS = ℏωP(S /(2S + 1))1/2 with S = k. (a) Four lowest plasmon modes ℏωk,m for m = 0 and Nsub = 2.4. (b) Same as (a) for m = 1. (c and d) Same as (a) and (b), respectively, for Nsub = 1.5.

associated with the dipolar m = 1 and m = 0 modes to higher S -order modes are quite similar (compare Figure 10 with the Figure 2(SI) in the Supporting Information). This difference is directly correlated to the fact that, unlike the dimer case, both bright branches of the sphere/interface system are bonding branches exhibiting quite similar hybridization effects. 3.3. Comparison with “Exact” Calculations. As for the dimer system (section 2.3), numerical tests have been carried out, first in order to assess the suitability of the PI−PH model within the quasistatic approximation, and second to determine its range of applicability for larger sizes. To achieve these objectives, absorption cross sections computed in the frame of the PI−PH model (eq 48) have been compared against results computed using a home-developed code (see the Supporting 28129

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

modulus of the wave vector, λ the wavelength, and c the vacuum velocity of light]. In our numerical MEM code39 the absorption and scattering cross sections are defined conveniently as the time-averaged absorbed/scattered powers Pabs,sca divided by the incoming electromagnetic field intensity I0 = cE20/(8π) (atomic units), although the effective field irradiating the jellium sphere is the sum of the incoming and reflected fields (the total irradiating field, Ei + ER, is strongly inhomogeneous when the index Nsub is large, explaining why a conventional choice is required for defining cross sections in this context). As pointed out at the beginning of section 3.1 the effects on the jellium sphere of the additional polarization charge distribution, that is induced on the plane interface by the incoming plane wave, can be taken into account in the PI−PH model in a very simple way, namely in adding the reflected field to the incoming electric field. Consequently, for allowing a direct quantitative comparison, the m = 0 and |m|=1 contributions to the MEM spectra have to be compared with the corresponding PI−PH model formula eq 48 multiplied by, respectively, the square z- and x-components of the total irradiating electric field Ei + ER [assumed homogeneous (see eqs 42−45) since the electrostatic PI−PH model is grounded on the quasistatic approximation] divided by |E0|2 (counterpart of the factor 1/2 involved in the dimer case, in section 2.3, for comparing MEM and D−PH model spectra). The z- and xcomponents of Ei + ER at the center of the jellium sphere have been selected for defining the weighting multiplicative factors [the analytical formula allowing computing the axial and transverse components can be found in the Supporting Information of ref 39 (eqs 11 and 12)]. In the case of ppolarized light excitation, the weighting factors az and ax for the axial and transverse excitations write as az = [sin(θ )]2 {1 + r p2(θ ) + 2rp(θ )cos[2k(d + R )cos(θ)]} (49) 2

ax = [cos(θ )] {1 +

r p2(θ )

− 2rp(θ )cos[2k(d + R )cos(θ)]} (50)

Figure 10. Evolution of the normalized oscillator strengths Ck′ (eq 46) of the plasmon modes of the jellium sphere/interface system as a function of the jellium sphere-image sphere distance 2d (see the caption of Figure 9 for details). When the distance d is decreased the oscillator strength of the degenerate dipolar excitations (S = 1, m = 0, ±1) of the isolated jellium sphere is progressively transferred to highorder multipolar modes S = k > 1 (for each m value) of the sphere/ interface system, through the electrostatic coupling of the charge distributions located on the surfaces of the jellium sphere and of its image.

where rp(θ) is the Fresnel reflection coefficient appropriate for p-polarized light.51 Typical examples (with Nsub = 2.4) illustrating these investigations are provided in the Supporting Information (Figures 7(SI)−9(SI)). For allowing valuable comparison with the spectra illustrating the dimer system (Figures 2, 4, and 5), the same sphere radii and intersphere distance-to-radius ratios have been selected. Obviously one has to keep in mind that, in the present system, the hybridization mechanism (plasmon mode frequencies and oscillator strength transfers) depends on the value of the medium index Nsub, in addition to the geometrical parameters. First, the predictions of the PI−PH model described in this work are found to agree perfectly with the exact MEM results within the quasistatic approximation, with regard to both the plasmon frequencies and the peak intensities (see Figures 7(SI) and 8(SI) in the Supporting Information). As in the case of the dimer system the determination of the range of applicability of the PI−PH model has been addressed in testing its scaling properties, in enlarging progressively the sphere radius R and the sphere− interface distance d, in keeping constant the ratio 2d/R = 0.2 (see Figure 9(SI) in the Supporting Information). Roughly, the same conclusions as those stated in the case of the dimer system apply. For very small radii the absorption prevails and the extinction spectra scale as the jellium sphere volume. For

Information in ref 39),39 based on an exact multipole expansion method.53,54 As before, a Drude-Sommerfeld dielectric function has been taken in the MEM calculations for describing the dielectric properties of the jellium sphere, and a very small damping parameter Γ has been chosen in order to separate the plasmon excitation lines on the energy scale. A monochromatic (frequency ω/2π) p-polarized light excitation at the oblique incidence θ = 45° with respect to the z-axis (see Figure 8) has been assumed to ensure the simultaneous excitation of the m = 0 modes (longitudinal/axial modes; ez polarization) and m = 1 modes (transverse modes, ex polarization). Some technical remarks, of major importance in order to compare MEM and PH model spectra, deserve to be emphasized. The amplitude of the incoming electric field Ei = E0 exp(ik·r − ωt) is (E0/√2)(ex + ez) [|k| = ω/c = 2π/λ is the 28130

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

effects fully included). This comparative analysis demonstrates that the appealing electrostatic D−PH model, although of different nature compared to more conventional electrodynamic formalisms, is a perfect quantitative approach within the quasistatic limit with regard to both the plasmon mode energies and the sharing of the dipole oscillator strengths of the isolated spheres among the higher-order modes of the dimer. However, for large overall dimer sizes, roughly above 40−50 nm, quantitative and qualitative discrepancies are noticed, due to retardation effects. These effects, disregarded in the D−PH model, induce large modifications of the line-shapes (intensity and width) and red-shifts of the modes, with an eventual spectral mixing and reordering of the dimer mode spectrum on the frequency scale. In particular, the D−PH model dark modes can contribute noticeably to the spectra. This study has pointed out that, for large particles, caution has to be taken for analyzing and labeling spectral structures on the basis of the PH model predictions only. Third, the PH approach has been applied for dealing with the sphere/dielectric interface system. Within the quasistatic limit the method of images allows the jellium sphere-plane dielectric interface interaction to be expressed as the interaction of the jellium sphere with a fictitious mirror jellium sphere located below the interface, bringing into light the tight relationship between the dimer and particle/interface systems. However, because the dynamics of the mirror electron fluid is enslaved to that of the jellium sphere, two times less dynamical variables are involved, which results in large differences between both systems. In particular, each plasmon mode of the isolated jellium sphere gives birth to only a single bonding and bright - frequency branch. Finally, the predictions of the PI−PH model have been compared with those computed in the frame of the multipole expansion method, showing that the range of applicability and limits of the PH models appropriate to both systems are rather similar. It is clear that the approach described in the present work can be applied to more complex supported nanostructures involving several subunits and other shapes, as well as underlying dielectric backgrounds (embedding medium and polarizable jellium).

larger sizes a strong departure from this scaling law is observed, with regard to both the spectral location of the plasmon modes and the characteristics of the line-shapes (intensity and width). For this system, the spectra are however less intricate, as compared to the dimer case, because the number of modes is two times smaller, and the evolution of the spectra for increasing system size seems more regular. The Figures 7(SI)−9(SI) in the Supporting Information clearly show that the optical response of an isolated jellium sphere in air differs strongly from that of the same sphere located above a highly reflecting dielectric substrate since, in the first case a single peak at ωp/√3 is expected, whereas in the second case a rich multimodal spectral response, strongly dependent on the particle−interface distance, is observed. Numerous examples have been reported on in our previous work, ref 39, on silver and gold particles. For completeness absorption spectra of a jellium sphere of radius 5 nm, subject to a p-polarized light excitation at the oblique incidence 45°, computed within the multipole expansion method39 for various particle−interface distances, are displayed in Figure 10(SI) in the Supporting Information. In the displayed example, the dielectric index of the substrate is Nsub = 2.4. For large distances, the optical response consists of a single peak converging toward the response expected for the isolated jellium sphere.55 When the distance d is decreased, the peak red-shifts due to screening effects and concomitantly broadens owing to the substrate-induced removal of the transverse (m = 1; polarization along the x-axis) and longitudinal (m = 0; polarization along the z-axis) dipolar mode degeneracy. For a small enough distance (d = 2 nm) the splitting of the initial peak into two individual underlying excitations is evidenced. When the distance d is decreased further, a rich multimodal spectral pattern, extending over a very large energy range, develops. In this strong coupling regime the oscillator strengths of the transverse and longitudinal dipolar excitations of the isolated sphere are progressively and noticeably transferred to the higher S -order modes (belonging to the subsets m = 1 and m = 0) of the particle/interface system due to the hybridization mechanism. Taking into account the scaling properties of the quasistatic PI−PH model the spectral locations of each mode in Figure 10(SI) (MEM calculations) can be easily correlated to the plasmon energies plotted in Figure 9a,b, except for a very tiny spectral red-shift of the MEM peaks because of the retardation effects. Finally it should be emphasized that much lesser dramatic differences between spectra from an isolated sphere and the same sphere above a dielectric interface are expected in the case of substrate of low reflectance, as thin dielectric films.



ASSOCIATED CONTENT

S Supporting Information *

Section A, Plasmon hybridization (PH) model for a jellium sphere dimer; Section B, Illustrative examples and qualitative analysis; Section C, Plasmon hybridization (PH) model for the sphere/interface system: tests of accuracy and range of applicability. This material is available free of charge via the Internet at http://pubs.acs.org.



4. CONCLUSION In this work the PH model developed by the Nordlander group for dealing with nanoparticle dimers has been generalized to the case of a particle above a plane dielectric interface, and both models have been exemplified in the simple case of jellium metal spheres. First, a detailed description of the PH model for jellium sphere dimers has been reported, as well as additional extensions regarding the coupling with an external electric field and the transfer of the oscillator strength of the dipole excitations of both isolated spheres to higher-order plasmon modes of the dimer. Second the predictions of the electrostatic D−PH model have been compared, through absorption and extinction cross-section spectra, to exact results computed in the frame of the multipole expansion method (retardation

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



REFERENCES

(1) Stewart, M. E.; Anderton, C. R.; Thompson, L. B.; Maria, J.; Gray, S. K.; Rogers, J. A.; Nuzzo, R. G. Nanostructured Plasmonic Sensors. Chem. Rev. 2008, 108, 494−521. (2) Giannini, V.; Fernandez-Dominguez, A. I.; Heck, S. C.; Maier, S. A. Plasmonic Nanoantennas: Fundamentals and Their Use in Controlling the Radiative Properties of Nanoemitters. Chem. Rev. 2011, 111, 3888−3912.

28131

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

(3) Halas, N. J.; Lal, S.; Chang, W. S.; Link, S.; Nordlander, P. Plasmons in Strongly Coupled Metallic Nanostructures. Chem. Rev. 2011, 111, 3913−3961. (4) Mayer, K. M.; Hafner, J. H. Localized Surface Plasmon Resonance Sensors. Chem. Rev. 2011, 111, 3828−3857. (5) Boisseau, P.; Houdy, Ph.; Lahmani, M. Nanoscience: Nanobiotechnology and Nanobiology; Springer: New York, 2010. (6) Jain, P. K.; El-Sayed, I. H.; El-Sayed, M. A. Au Nanoparticles Target Cancer. Nanotoday 2007, 2, 18−29. (7) Kelly, K. L.; Coronado, E.; Zhao, L. L.; Schatz, G. C. The Optical Properties of Metal Nanoparticles: The influence of Size, Shape, and Dielectric Environment. J. Phys. Chem. B 2003, 107, 668−677. (8) Kreibig, U.; Vollmer, M. Optical Properties of Metal Clusters; Springer: Berlin, 1995. (9) Bohren, C. F.; Huffman, D. P. Absorption and Scattering of Light by Small Particles; Wiley: New York, 1983. (10) Wiley, B. J.; Im, S. H.; Li, Z. Y.; McLellan, J.; Siekkinen, A.; Xia, Y. Maneuvering the Surface Plasmon Resonance of Silver Nanostructures through Shape-Controlled Synthesis. J. Phys. Chem. B 2006, 110, 15666−15675. (11) Liz-Marzan, L. M. Nanometals: Formation and Color. Mater. Today 2004, 7, 26−31. (12) Halas, N. J.; Lal, S.; Link, S.; Chang, W. S.; Natelson, D.; Hafner, J. H.; Nordlander, P. A Plethora of Plasmonics from the Laboratory for Nanophotonics at Rice University. Adv. Mater. 2012, 24, 4842−4877. (13) Brandl, D. W.; Mirin, N. A.; Nordlander, P. Plasmon Modes of Nanosphere Trimers and Quadrumers. J. Phys. Chem. B 2006, 110, 12302−12310. (14) Mirin, N. A.; Bao, K.; Nordlander, P. Fano Resonances in Plasmonic Nanoparticle Aggregates. J. Phys. Chem. A 2009, 113, 4028−4034. (15) Jain, P. K.; Eustis, S.; El-Sayed, M. A. Plasmon Coupling in Nanorod Assemblies: Optical Absorption, Discrete Dipole Approximation Simulation, and Exciton-Coupling Model. J. Phys. Chem. B 2006, 110, 18243−18253. (16) Grillet, N.; Manchon, D.; Bertorelle, F.; Bonnet, C.; Broyer, M.; Cottancin, E.; Lermé, J.; Hillenkamp, M.; Pellarin, M. Plasmon Coupling in Silver Nanocube Dimers: Resonance Splitting Induced by Edge Rounding. ACS Nano 2011, 5, 9450−9462. (17) Wang, H.; Brandl, D. W.; Nordlander, P.; Halas, N. J. Plasmonic Nanostructures: Artificial Molecules. Acc. Chem. Res. 2007, 40, 53−62. (18) Lombardi, A.; Grzelczak, M. P.; Crut, A.; Maioli, P.; PastorizaSantos, I.; Liz-Marzan, L. M.; Del Fatti, N.; Vallée, F. Optical Response of Individual Au-Ag@SiO2 Heterodimers. ACS Nano 2013, 7, 2522− 2531. (19) Lassiter, J. B.; Sobhani, H.; Fan, J. A.; Kundu, J.; Capasso, F.; Nordlander, P.; Halas, N. J. Fano Resonances in Plasmonic Nanoclusters: Geometrical and Chemical Tunability. Nano Lett. 2010, 10, 3184−3189. (20) Funston, A. M.; Novo, C.; Davis, T. J.; Mulvaney, P. Plasmon Coupling of Gold Nanorods at Short Distances and in Different Geometries. Nano Lett. 2009, 9, 1651−1658. (21) Marhaba, S.; Bachelier, G.; Bonnet, C.; Broyer, M.; Cottancin, E.; Grillet, N.; Lermé, J.; Vialle, J. L.; Pellarin, M. Surface Plasmon Resonance of Single Gold Nanodimers near the Conductive Contact Limit. J. Phys. Chem. C 2009, 113, 4349−4356. (22) Slaughter, L. S.; Wu, Y.; Willingham, B. A.; Nordlander, P.; Link, S. Effects of Symmetry Breaking and Conductive Contact on the Plasmon Coupling in Gold Nanorod Dimers. ACS Nano 2010, 4, 4657−4666. (23) Rechberger, W.; Hohenau, A.; Leitner, A.; Krenn, J. R.; Lamprecht, B.; Aussenegg, F. R. Optical Properties of Two Interacting Gold Nanoparticles. Opt. Commun. 2003, 220, 137−141. (24) Acimovic, S. S.; Kreuzer, M. P.; Gonzales, M. U.; Quidant, R. Plasmon Near-Field Coupling in Metal Dimers as a Step towards Single-Molecule Sensing. ACS Nano 2009, 3, 1231−1237. (25) Atay, T.; Song, J. H.; Nurmikko, A. V. Strongly Interacting Plasmon Nanoparticle Pairs: From Dipole-Dipole Interaction to Conductively Coupled Regime. Nano Lett. 2004, 4, 1627−1631.

(26) Jain, P. K.; El-Sayed, M. A. Universal Scaling of Plasmon Coupling in Metal Nanostructures: Extension from Particle Pairs to Nanoshells. Nano Lett. 2007, 7, 2854−2858. (27) Barrow, S. J.; Funston, A. M.; Gomez, D. E.; Davis, T. J.; Mulvaney, P. Surface Plasmon Resonances in Strongly Coupled Gold Nanosphere Chains from Monomer to Hexamer. Nano Lett. 2011, 11, 4180−4187. (28) Jain, P. K.; Huang, W.; El-Sayed, M. A. On the Universal Scaling Behavior of the Distance Decay of Plasmon Coupling in Metal Nanoparticle Pairs: A Plasmon Ruler Equation. Nano Lett. 2007, 7, 2080−2088. (29) Gunnarsson, L.; Rindzevicius, T.; Prikulis, J.; Kasemmo, B.; Käll, M.; Zou, S.; Schatz, G. C. Confined Plasmons in Nanofabricated Single Silver Particle Pairs: Experimental Observations of Strong Interparticle Interactions. J. Phys. Chem. B 2005, 109, 1079−1087. (30) Romero, I.; Aizpurua, J.; Bryant, G. W.; Garcia de Abajo, F. J. Plasmons in Nearly Touching Metallic Nanoparticles: Singular Response in the Limit of Touching Dimers. Opt. Express 2006, 14, 9988−9999. (31) Nordlander, P.; Oubre, C.; Prodan, E.; Li, K.; Stockman, M. I. Plasmon Hybridization in Nanoparticle Dimers. Nano Lett. 2004, 4, 899−903. (32) Mishchenko, M. I.; Hovenier, J. W.; Travis, L. D. Light Scattering by Nonspherical Particles: Theory, Measurements, and Applications; Academic Press: San Diego, 2000. (33) Mishchenko, M. I.; Travis, L. D.; Lacis, A. A. Scattering, Absorption, and Emission of Light by Small Particles; Cambridge University Press: Cambridge, 2002. (34) Kahnert, F. M. Numerical Methods in Electromagnetic Scattering Theory. J. Quant. Spectrosc. Radiat. Transfer 2003, 79−80, 775−824. (35) Doicu, A.; Wriedt, T.; Eremin, Y. A. Light Scattering by Systems of Particles: Null-field Method with Discrete Sources; Springer: Berlin, 2006. (36) Yurkin, M. A.; Hoekstra, A. G. The discrete dipole approximation: an overview and recent developments. J. Quant. Spectrosc. Radiat. Transfer 2007, 106, 558−589. (37) Prodan, E.; Radloff, C.; Halas, N. J.; Nordlander, P. A Hybridization Model for the Plasmon Response of Complex Nanostructures. Science 2003, 302, 419−422. (38) Prodan, E.; Nordlander, P. Plasmon Hybridization in Spherical Nanoparticles. J. Chem. Phys. 2004, 120, 5444−5454. (39) Lermé, J.; Bonnet, C.; Broyer, M.; Cottancin, E.; Manchon, D.; Pellarin, M. Optical Properties of a Particle above a Dielectric Interface: Cross Sections, Benchmark Calculations, and Analysis of the Intrinsic Substrate Effects. J. Phys. Chem. C 2013, 117, 6383−6398. (40) Nordlander, P.; Prodan, E. Plasmon Hybridization in Nanoparticles near Metallic Surfaces. Nano Lett. 2004, 4, 2209−2213. (41) Le, F.; Lwin, N. Z.; Steele, J. M.; Käll, M.; Halas, N. J.; Nordlander, P. Plasmons in the Metallic Nanoparticle-Film System as a Tunable Impurity Problem. Nano Lett. 2005, 5, 2009−2013. (42) Knight, M. W.; Wu, Y.; Lassiter, J. B.; Nordlander, P.; Halas, N. J. Substrates Matter: Influence of an Adjacent Dielectric on an Individual Plasmonic Nanoparticle. Nano Lett. 2009, 9, 2188−2192. (43) Wu, Y.; Nordlander, P. Finite-Difference Time-Domain Modeling of the Optical Properties of Nanoparticles near Dielectric Substrates. J. Phys. Chem. C 2010, 114, 7302−7307. (44) Arfken, G. B.; Weber, H. J. Mathematical Methods for Physicists, 4th ed.; Academic Press: San Diego, 1995. (45) The present numerical investigations aim only at testing the accuracy of the D-PH model in the frame of classical electrodynamics. So the D-PH and MEM approaches are applied, whatever the jellium radii are, irrespective of the fact that, for very small sizes, a selfconsistent quantum description including nonlocal screening effects is required for ensuring correct predictions. (46) Ruppin, R. Surface Modes of Two Spheres. Phys. Rev. B 1982, 26, 3440−3444. (47) Chu, P.; Mills, D. L. Electromagnetic Response of Nanosphere Pairs: Collective Plasmon Resonances, Enhanced Fields, And LaserInduced Forces. Phys. Rev. B 2008, 77, 045416 (1−10). 28132

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133

The Journal of Physical Chemistry C

Article

(48) Wokaun, A.; Gordan, J. P.; Liao, P. F. Radiation Damping in Surface-Enhanced Raman Scattering. Phys. Rev. Lett. 1982, 48, 957− 960. (49) Meier, M.; Wokaun, A. Enhanced Fields on Large Metal Particles: Dynamic Depolarization. Opt. Lett. 1983, 8, 581−583. (50) These preliminary comments apply likewise in the case of concentric metal nanoshells including underlying dielectric backgrounds (see, for instance, ref 38). The dielectric media, more precisely the induced surface polarization charge densities, are involved in the PH modeling only through their contributions to the total electrostatic potential that is experienced by the free charge densities on the metallic surfaces. (51) Jackson, J. D. Classical Electrodynamics, 3rd ed.; John Wiley & Sons, Inc.: New York, 1998. (52) Ruppin, R. Surface Modes and Optical Absorption of a Small Sphere above a Substrate. Surf. Sci. 1983, 127, 108−118. (53) Fucile, E.; Denti, P.; Borghese, F.; Saija, R.; Sindoni, O. I. Optical Properties of a Sphere in the Vicinity of a Plane Interface. J. Opt. Soc. Am. A 1997, 14, 1505−1514. (54) Wriedt, T.; Doicu, A. Light Scattering from a Particle on or near a Surface. Opt. Commun. 1998, 152, 376−384. (55) The uppermost spectrum (d = 25 nm) has to be divided by az + ax (eqs 49 and 50) to obtain the absolute cross section computed for an isolated sphere subject to a homogeneous plane wave irradiation.

28133

dx.doi.org/10.1021/jp509238j | J. Phys. Chem. C 2014, 118, 28118−28133