1. Introduction

Periodic arrays of plasmonic nanospheres have been studied because of their interesting performance in several practical applications. For example, to achieve near-field enhancement and super-resolution at optical wavelengths [1–5]; to achieve directive radiation [6,7]; to perform microscopy and spectroscopy using SERS [8,9]; to design biosensors [10,11], and DNA detectors [12]. However, despite these innovative applications, it is still not clear how to associate determinate characteristics with a particular set of excitation wavelengths. Several difficulties are still encountered in the analysis and interpretation of mode propagation in one dimensional (1D-periodic, linear chains) and two dimensional (2D-periodic) arrays of metallic nanospheres when varying frequency, for example. Based on what it has been published so far it is somehow hard to discern if a mode is physical (excitable) or not. The scope of this paper is the characterization of the bound (non radiating) and leaky (radiating) modes with complex wavenumber in linear chains of plasmonic nanospheres (Fig. 1 ), accounting for metal losses, as we recently did in [13] for 2D-periodic arrays.

 

Fig. 1 Linear chain of plasmonic nanospheres embedded in a homogeneous medium with permittivity ${\epsilon }_h$. The radius of each nanosphere is r, and the array spatial period is d.

Download Full Size | PDF

Modes of linear chains of metallic nanospheres, object of this paper, have been studied in [7,14–26], assuming either lossless or lossy metal and using the single dipole approximation (SDA), described in [27,28]. In particular, in [14] it has been shown that the dispersion relation k-β (where k is the host medium wavenumber, and β is the wavenumber of the guided wave) for a linear chain of gold nanospheres, obtained from finite-difference time-domain (FDTD) calculations, closely agreed with that predicted by the SDA in case of real wavenumbers. In [15], the authors have calculated the dispersion relation for a linear chain of gold spherical particles in a uniform host in the case of particle dimensions much smaller than the wavelength, including all the multipolar terms within the quasistatic approximation, finding that the dipole approximation was reasonably accurate for a periodicity larger than three times the nanosphere radius, becoming increasingly inaccurate for smaller periodicity. In [16], the SDA with the expression for the nanosphere polarizability taken from Mie theory has been applied to compute the dispersion diagram (the real part of the wavenumber) of lossless chains of gold nanospheres for longitudinal and transverse polarizations, comparing the results with the ones reported in [14,15]. In [7,18–21,26], the authors used polylogarithmic functions to analytically continue the validity of the dispersion relation into the complex wavenumber plane, and they also accounted for metal losses. However, conditions to classify proper and improper modes (proper/improper modes are defined as those that vanish/grow getting away from the 1D array axis, and this concept will be better clarified in Sec. 2.2), as well as to distinguish the subset of physical waves allowed in the linear chain from the complete mathematical set of modes, have not been developed yet when using polylogarithmic functions. In particular in [26], the authors derived the 1D dispersion relation by considering all multipole field contributions, showing that in the case of a linear chain of plasmonic nanospheres the dispersion relation is well approximated by the one obtained through the SDA with the expression for the nanosphere polarizability taken from Mie theory, as we do in this paper, assuming small nanospheres with respect to the wavelength. Experimentally, modes have been analyzed in linear chains of gold nanoparticles comparing the dispersion diagram results from measurements (real part of the wavenumber) with corresponding ones obtained from theoretical fully retarded solutions and numerical FDTD simulations [22,23]. In [24], the effect of small random disorder (e.g., due to fabrication limitations) in linear chains of metallic nanoparticles has been analyzed theoretically. Specifically, guidance properties of these aperiodic arrays (in fact, the periodicity has been broken by a small disorder) have been analyzed, resulting into additional radiation losses for the guided mode. A review of the state of the art of studies performed for both 1D- and 2D-periodic arrays of plasmonic nanospheres, including mode analysis, fabrication possibilities, and possible applications, has been outlined in [25]. To conclude, we want to emphasize that all the aforementioned excellent referenced work did not provide yet a complete characterization and classification of the modes in a linear chain of plasmonic nanospheres, in terms of leaky/bound, proper/improper, forward/backward, including the conditions for their physical existence, which is the goal of the present paper. We want to stress that for a physical mode we mean a mode that can be excited (launched) into the array by a localized source in proximity of the array, or by a defect or by an array truncation. Other modes that are not classified as physical here may, however, be excited by much more complicated source distributions, though this study is not within the scope of this paper. The analysis in this paper is a necessary step towards a thorough understanding of modal description in periodic arrays. Several results and observations provided in this paper have not been provided elsewhere.

The novelty of the approach presented here is that it allows for the tracking and especially for the characterization of the evolution of all the modes varying frequency. Each nanosphere is modeled as a single dipole, by using the SDA [25,27,28] and the metal permittivity is described by the Drude model. The computation of the complex zeroes of the dispersion relation of the 1D-periodic array is performed by using a periodic Green’s function (GF), analytically continued into the complex wavenumber space by using the Ewald representation. The computation of the complex roots of complex valued functions is a difficult problem in its own rights. We have (i) implemented the Muller’s method [29,30] (which can be used to find complex zeroes of analytic functions) and (ii) used the built-in IMSL function ZANLY for Fortran [31] (which uses Muller’s method) for verification purposes. The two obtained solutions are in excellent agreement with each other. The use of the periodic GF in the computation of the complex zeroes will be explained in Sec. 2.1. As a consequence, full characterization of the modes with complex wavenumber in the linear chain is provided, in terms of guidance/radiance, proper/improper, physical/nonphysical, and propagation direction conditions. In this paper, a new numerical procedure that uses the Ewald representation for the periodic dyadic GF for linear arrays has been introduced, based on previous scalar developments [32–34]. Our method is alternative to the one used in [18,20,21,26] based on polylogarithms, and future studies will be devoted to relate these two techniques.

The structure of the paper is as follows. The mathematical description for the simulation model adopted to perform the analyses described in this paper is briefly summarized in Sec. 2, including a discussion regarding conditions of guidance, radiance and physical existence of modes and their constituent Floquet waves. Then, the characterization of the complex bound and leaky modes in linear chains of plasmonic nanospheres is reported in Sec. 3, considering nanospheres made of silver, and accounting for metal losses. Though our examples relate to an array with a precise set of representative parameters, the analysis, wave classification and results are general. The analysis of guidance and radiance is reported in Sec. 4. Parameterization of modes with respect to the geometrical parameters of the array is shown in Sec. 5. Criteria for physical existence (excitability) of modes are presented in Appendix A. The dyadic form of the Ewald representation of the periodic GF evaluation used in this paper is provided in Appendix B.

2. Theoretical background

The monochromatic time harmonic convention $\exp \left(-i\omega t\right)$ is assumed. Moreover, in the following equations, bold letters refer to vector quantities, a caret on top of a bold letter refers to unit vector quantities, and a bar under a bold letter refers to dyadic quantities.

2.1 Simulation model

As stated in Sec. 1, according to SDA, each plasmonic nanosphere of the linear chain, at optical frequencies, acts as a single electric dipole and can be described by the induced dipole moment

Consider now a linear chain of nanospheres as in Fig. 1, for which each nanosphere is placed at positions $z_n=nd$, with $n=0,\pm 1,\pm 2,\dots$, and d is the spatial period of the chain. Suppose that the array is immersed in a homogeneous background, with relative permittivity ${\epsilon }_h$, and that the array is then excited by a plane wave or by a quasi-periodic excitation with wavevector $k_{\text{B}}=k_z\widehat{z}$ [28]. Consequently, each spherical nanoparticle will have a dipole moment equal to $p_n=p_0\text{exp}\left(in{k}_{z}d\right)$, with $p_0$ defined as in Eq. (1) by using n = 0. Then, the local field acting on a nanosphere at position $z_0=0$ is given by

Mode analysis in the linear chain is then performed by computing the complex zeroes of the homogeneous linear system in Eq. (6), after having imposed no excitation source (i.e., $E^{\text{inc}}\left(k_z\right)=0$). In other words, the complex $k_z$-zeroes of the determinant of $\underset{\_}{A}\left(k_z\right)$ have to be determined (i.e., $k_z$ such that $\det \underset{\_}{A}\left(k_z\right)=0$) using the methods described in Sec. 1. The series inside the GF definition in Eq. (4) is slowly convergent, and it is even divergent for complex values of the wavenumber $k_z$, and some technique must be adopted to allow analytic continuation of the GF into the complex $k_z$ plane. The GF representation adopted in this work makes use of Ewald’s method for linear arrays, which, besides providing analytic continuation, it also provides series with Gaussian convergence and only a handful of terms are needed [32–34]. The dyadic form of the Ewald representation of the periodic GF evaluation used in this paper is provided in Appendix B.

2.2 Floquet waves representation

The electric field in the linear chain can be described as $E\left(r+nd\widehat{z},k_z\right)=E\left(r,k_z\right)\text{exp}\left(in{k}_{z}d\right)$, where d is the spatial period of the periodic structure, $n=0,\pm 1,\pm 2,\dots$, $k_z={\beta }_z+i{\alpha }_z$ is the wavenumber along the array axis, and $r=x\widehat{x}+y\widehat{y}+z\widehat{z}$. The field relative to a mode in the linear chain is expressed in term of a Fourier series expansion, and thus represented as the superposition of Floquet wave (FW) spatial harmonics as

2.3 Physical excitation conditions by a source, a defect or a truncation

Among all the mathematical wavenumber solutions of $\det \underset{\_}{A}\left(k_z\right)=0$ only a subset represents physical waves, i.e., waves that can be excited by a localized source, or defect, or array truncation, as described in Appendix A. The physical modes that can be excited in linear chains of plasmonic nanospheres, according to Appendix A and what follows, are summarized in Table 1 , similar to what was previously reported for the 2D-periodic case in [41,43], based on the complex wavenumber of the fundamental harmonic (with p = 0). Table 1 applies to the case of chain period shorter than half wavelength, since it implies that only wavenumbers in the fundamental BZ may represent fast waves.

Table 1. Classification of Physical Modes with Complex Wavenumbers, Based on the Wavenumber of the Mode-Harmonic in the Fundamental BZ

View Table

Assume that a point source is near the origin (Fig. 2 ) and an observer is located along the positive z axis of the array. The field at the observer is represented by a spectral path integral on the real axis as shown in the figure in Appendix A. Deformation onto the steepest descent paths encounters one or more periodic sets of poles that represent the wavenumbers of the physical FWs [41,44,45]. A mode is represented by an infinite set of wavenumbers, according to Eq. (7). Periodic sets of poles that are not captured in the path deformation do not contribute to the final field value and therefore represent modes that are not excited. In order to be captured in the path deformation in Appendix A, a pole representing a wave with complex wavenumber $k_z={\beta }_z+i{\alpha }_z$ should have ${\alpha }_z>0$. Therefore, considering the modal harmonic in the fundamental BZ (often the one with p = 0), a physical forward wave belongs to the I quadrant of the complex $k_z$ plane (i.e., ${\beta }_z{\alpha }_z>0$), whereas a physical backward one to the II quadrant (i.e., ${\beta }_z{\alpha }_z<0$). Moreover, to be physical, not only the wave needs to have a wavenumber with ${\alpha }_z>0$, but also its wavenumber has to respect the conditions stated in Table 1, which are also valid in the case of an observer along the negative z axis for which physical waves have ${\alpha }_z<0$. Further discussion about integration path and field representation is provided in Appendix A and references therein.

 

Fig. 2 Illustration of a physical complex wave for an observer along the positive z axis, according to its classification provided in Table 1. (A) and (B) represent forward and backward bound waves. (C) and (D) represent forward and backward leaky waves.

Download Full Size | PDF

Figure 2 shows how each physical wave outlined in Table 1 evolves in the analyzed linear chain, for the case in which the observer is along the positive z-direction. If the wave is proper, as defined in the previous section, it vanishes for $\rho \to \infty$ (as can be observed in Figs. 2(A), 2(B) and 2(D), where the wave decays along the radial direction). If the wave is improper, instead, it grows for $\rho \to \infty$, as shown in Fig. 2(C). Note that in this case, the leaky wave existence domain is limited in space as described at the end of Sec. 2.2, and therefore it does not violate any physical principle [41,42]. If the wave is forward, the phase travels towards the positive z-direction, as highlighted by the red arrows in Figs. 2(A) and 2(C), and attenuates while it travels. Whereas, if the wave is backward, the phase travels towards negative z-direction, as highlighted by the red arrows in Figs. 2(B) and 2(D), but still attenuates along positive z (for conservation of energy, a wave cannot grow along positive z). If the wave is bound, then no energy is radiated away (as in Figs. 2(A) and 2(B)); radiation phenomena are present when the wave is leaky, as shown in Figs. 2(C) and 2(D) by the black arrows.

2.4 Conditions for guidance and radiance

An important part of this work lies in the fact that we show which modes are physical (i.e., that are excitable by a localized source), and consequently employable in practical applications. For this reason, the conditions to discriminate physical bound and radiating (leaky) modes are here stated based on Table 1 and the discussion in Appendix A.

A bound mode is a mode that has all its Floquet wavenumbers outside the “radiating spectral region”, i.e., the spectral region for which$\left|{\beta }_{z,p}\right|<k$. Therefore, a bound mode is a mode with $\left|{\beta }_{z,p}\right|>k$ for any p. Taking into account the physical condition, considering the wave harmonic with ${\beta }_z$ in the fundamental BZ, an observer along the positive z-direction would be reached by a bound wave (excited by a localized source near the origin) either when (i) it is proper and $k<{\beta }_z$ and ${\alpha }_z$ is positive (forward wave), or when (ii) it is proper and ${\beta }_z<-k$ and ${\alpha }_z$ is positive (backward wave), i.e., the wavenumber is either in the first or in the second quadrant of the complex $k_z$ plane. Improper bound modes cannot be launched by a source in Fig. 2.

A leaky mode is a mode that has at least one Floquet wavenumber within the radiating spectral region, or $\left|{\beta }_{z,p}\right|<k$ for some p. Taking into account the physical condition, and considering the harmonic ${\beta }_z$ in the fundamental BZ, an observer along the positive z-direction would be reached by a radiating wave (excited by a localized source near the origin) either when (i) it is improper and $0<{\beta }_z<k$ and ${\alpha }_z$ is positive (forward wave), or when (ii) it is proper and $-k<{\beta }_z<0$ and ${\alpha }_z$ is positive (backward wave).

A guided mode (leaky or bound) travels a long distance when the imaginary part of its wavenumber is small with respect to the host wavenumber, or ${\alpha }_z<<k$ (i.e., low decay).

For completeness, we want to emphasize that the higher order Floquet harmonics outside the radiating spectral region (with $\left|{\beta }_{z,p}\right|>k$,$p\ne 0$) can be either proper or improper non radiating waves, and only the ones represented by poles captured from the path deformation described in Appendix A can be excited by a source in Fig. 2, and thus be physical. Also, as mentioned in Sec. 2.3, in the following sections we restrict our analysis to the case of chain period shorter than half wavelength, where only wavenumbers in the fundamental BZ may represent fast (leaky) waves.

3. Mode analysis

In this section, complex modes in a linear chain of silver nanospheres embedded in free space (i.e., $k=k_0$) are shown for both transverse and longitudinal polarization states (T- and L-pol, respectively), accounting also for metal losses. The radius of the nanospheres is r = 25 nm and the spatial period is d = 75 nm. The adopted Drude model parameters for silver are ${\epsilon }_{\infty }=5$, plasma frequency ${\omega }_p=1.37\times {10}^{16}\text{rad/s }$ and damping factor $\gamma =27.3\times {10}^{12}{\text{s}}^{-1}$ [46,47]. In the L-pol case, each nanosphere has an induced dipole moment along the z-direction (i.e., $p=p_z\widehat{z}$). In the T-pol case, each nanosphere has an induced dipole moment orthogonal to z.

For both polarizations, if $k_z$ is a complex wavenumber solution of $\det \underset{\_}{A}\left(k_z\right)=0$, then also $-k_z$ is a solution. In absence of losses (case not reported here for the sake of brevity) a four-quadrant symmetry holds, thus the correspondent complex conjugate solutions, $\pm k_z^{*}$, would be solutions as well. When losses are present, the complex conjugate solutions are modified from $\pm k_z^{*}$ correspondingly to the amount of losses introduced in the model. Because of space limitation, and without deterring information, we show here only the fundamental p = 0 FW harmonic of complex modes. The periodic condition $k_{z,p}=k_z+2\pi p/d$ would determine the behavior of FW harmonics with wavenumbers in other BZs. Therefore, the focus is on the wavenumber solutions belonging to the I and II quadrants of the complex$k_z$ plane in the dispersion diagrams in the fundamental BZ, i.e., those solutions with ${\alpha }_z\ge 0$, which are the only ones that can be captured in the path deformation described in Appendix A when observing along the positive z axis. Note that this subset has the same information as the complete set of modes (i.e., a complete basis for the field) in the linear chain because of the symmetry property stating that $k_z$ and $-k_z$ are both solutions.

3.1 Transverse polarization (T-pol)

Figure 3 shows the dispersion diagrams for the real and imaginary parts of the wavenumber of the FW harmonic in the fundamental BZ of the modes with T-pol and relatively small imaginary part. Figure 4 shows the wavenumber evolutions in the complex $k_z$ plane varying frequency, using two different normalizations. We want to clarify that the reported curves in Figs. 3 and 4 are relative to both physical and nonphysical solutions (because together they form a complete basis for the field); physical ones are then distinguished by using Table 1.

 

Fig. 3 Dispersion diagram for T-pol. (a) Real part and (b) imaginary part of the wavenumber of the fundamental FW $k_z={\beta }_z+i{\alpha }_z$.

Download Full Size | PDF

 

Fig. 4 Modal wavenumbers in the complex relative $k_z$ plane for T-pol, normalized with respect to (a) the spatial period d, and (b) the host wavenumber $k=k_0\sqrt{{\epsilon }_h}$ (here ${\epsilon }_h=1$). Note that in (b), crossing “1” in the abscissa means crossing the light line. The shown curves are parametric curves of frequency, and the arrows indicate direction of increasing frequency.

Download Full Size | PDF

Modes ‘Proper 1’, ‘Proper 2’, ‘Improper 1’ and ‘Improper 2’ in Figs. 3 and 4 are computed using the method described in [33], where only m = 0 and −1 are considered, and adapted to the present case as in Appendix B. When the wavenumber dispersion curves change complex quadrant, they can be continued with dispersion curves computed with different m-values (see appendix B). We also show the solution ‘Improper 3′, computed with m = +1, as a continuation of the m = 0 red curve into the II and IV quadrants (Fig. 4), which also represents the wavenumber close to the complex conjugate of the ‘Improper 1’, according to the symmetry properties of the modal solutions described in Sec. 3. Though these waves computed with m values different from 0 and −1 are correct mathematical solutions, their physical validity requires further investigation, as explained in Appendix B. For an observer along the positive z-direction, according to Table 1 and Sec. 2.4, and looking at the complex $k_z$ plane I and II quadrants (the upper half) in Fig. 4 (i.e., the only part in which physical solutions could be found), the mode ‘Proper 1’ is backward and physical in the entire II quadrant (specifically, it is a physical leaky mode when $-k<{\beta }_z<0$, and physical bound mode when ${\beta }_z<-k$); similarly, the mode ‘Improper 1’ is a physical forward leaky mode in the I quadrant when ${\beta }_z<k$, whereas it is a nonphysical bound mode when ${\beta }_z>k$. The mode ‘Proper 2’ is a physical forward bound mode only in the I quadrant when ${\beta }_z>k$, whereas the mode ‘Improper 2’ is always nonphysical, being in the II quadrant (remember that improper modes in the II quadrant cannot be captured, according to Appendix A).

From the point of view of modal wavenumber evolution when increasing frequency, the wavenumber of the blue mode ‘Proper 1’ goes from the II to the I quadrant of the complex $k_z$ plane (or from the IV to the III due to the symmetry of the modal wavenumber solutions with respect to the origin stated in Sec. 3, i.e., $k_z$ and $-k_z$are both solutions, even in lossy cases), thus it crosses the imaginary axis (exhibiting its lowest ${\alpha }_z$ before this transition) and so a branch cut positioned there (as shown in Appendix A). Therefore, this modal wavenumber passes from the top to the bottom Riemann sheet, i.e., it becomes the one of the blue mode ‘Improper 1’. This transition is well observed in Fig. 4. It is noteworthy that due to this change in the Riemann sheet of the modal wavenumber, the physical backward leaky mode ‘Proper 1’ (II quadrant), when frequency increases, transitions into the physical forward leaky mode ‘Improper 1’ (I quadrant). The wavenumber of the red mode ‘Proper 2’ is almost over the light line ${\beta }_z=k$ at low frequencies. Increasing frequency, it departs from the light line in the bound region until its trajectory is inverted (Fig. 4), and then it goes from the I to the II quadrant (or from the III to the IV due to symmetry), and thus, similarly to what described for the mode ‘Proper 1’, crosses a branch cut, transitioning into the wavenumber of the mode ‘Improper 3′. With smaller losses, the wavenumber of the mode ‘Proper 2’ would be able to depart farther away from the light line in the I quadrant, reaching large values of ${\beta }_z$, and thus causing short plasmonic wavelengths.

3.2 Longitudinal polarization (L-pol)

Figure 5 shows the dispersion diagrams for the real and imaginary parts of the wavenumber of the FW harmonic in the fundamental BZ of the modes with L-pol and relatively small imaginary part. Figure 6 shows the wavenumber evolutions in the complex $k_z$ plane varying frequency, using two different normalizations. We want to clarify again that the reported curves in Figs. 5 and 6 are relative to both physical and nonphysical solutions; physical ones are then distinguished by using Table 1.

 

Fig. 5 As in Fig. 3, for L-pol.

Download Full Size | PDF

 

Fig. 6 As in Fig. 4, for L-pol. Inset in (b) shows the modal wavenumbers zero crossing.

Download Full Size | PDF

Modes ‘Proper 1’ and ‘Improper 1’ are computed using the method described in [33], and adapted to the present case as in Appendix B, where only m = 0 and −1 are considered. When the wavenumber dispersion curves change complex quadrant, they can be continued with different m-values (see appendix B). We also show the curve related to the ‘Improper 2’ mode obtained with m = +1, which represents the continuation of the ‘Proper 1’ solution into the complex $k_z$ plane quadrants II and IV (Fig. 6). This wavenumber also represents the wavenumber close to the complex conjugate solution symmetric to the ‘Improper 1’ solution, according to the symmetry properties described in Sec. 3. The physical validity of the waves found with m-values different than 0 and −1, requires further investigation, as explained in Appendix B. For an observer along the positive z-direction, according to Table 1 and Sec. 2.4, and looking at quadrants I and II in the upper half part of the complex $k_z$ plane in Fig. 6, the mode ‘Proper 1’ is a physical forward bound mode in the I quadrant when ${\beta }_z>k$; similarly, the mode ‘Improper 1’ is a physical forward leaky mode in the I quadrant when $0<{\beta }_z<k$. The mode ‘Improper 1’ in the II quadrant is always nonphysical.

From the point of view of modal wavenumber evolution, at high frequencies, the wavenumber of the blue mode ‘Proper 1’ is in the bound region in the I quadrant, as in Fig. 6(b). Decreasing frequency, the real part of its wavenumber increases while its imaginary part decreases up to about $kd/\pi \approx 0.45$, then both the real and imaginary parts of this modal wavenumber decrease. Around $kd/\pi \approx 0.39$, this wavenumber enters in the radiating spectral region (i.e., $\left|{\beta }_z\right|<k$). There, its imaginary part transitions from positive to negative values (see the inset in Fig. 6). This transition translates into the wavenumber passing from the I to the IV quadrant of the complex $k_z$ plane (or from the III to the II due to symmetry), thus crossing the $\mathrm{Re}\left(k_z\right)$ axis and so a branch cut located parallel to the $\mathrm{Re}\left(k_z\right)$ axis (in Appendix A, this branch cut is slightly above the axis when small ambient losses are present; it would be exactly on the axis in absence of ambient losses), and it changes Riemann sheet, i.e., it becomes the wavenumber of the blue mode ‘Improper 2’. This transition is well observed in Fig. 6. The wavenumber of the forward green mode ‘Improper 1’ starts in the I quadrant of the complex plane in Fig. 6. Increasing frequency, its imaginary part goes to zero while its real part remains under the light line (i.e., $\left|{\beta }_z\right|>k$) as in the inset in Fig. 6, thus it does not encounter any branch cut, as can be seen in Appendix A, and remains improper with negative growing imaginary part (this modal wavenumber transitions from the I to the IV quadrant, or from the III to the II due to symmetry).

Notice for example that in [21] modes ‘Improper 1’ and ‘Proper 1’ have been described as the same mode, discontinuous for a small frequency range; whereas in [26] they have been described as one continuous longitudinal mode.

4. Analysis of modal guidance and radiance

4.1 Bound modes for transverse and longitudinal polarizations

Accounting for the presence of losses in the silver nanospheres, the only bound mode for T-pol that travels without large decay is the mode ‘Proper 2’ in Figs. 3(a) and 3(b). This mode is bound, forward and physical (as described in Sec. 3.1); it presents growing imaginary part for increasing frequency, with very small positive attenuation constant near the light line for $kd/\pi \le 0.4$ (this mode lies on top of the light line for $kd/\pi <0.37$). In the entire domain of its physical existence, this mode can be guided very well, in the sense that it can travel long distances from the excitation source. The other physical bound mode for T-pol is the ‘Proper 1’ that has the imaginary part of its wavenumber which is large in comparison to the host wavenumber k. This mode is bound, backward and physical; it presents first decreasing imaginary part while frequency increases, a minimum at $kd/\pi \approx 0.4$ at which ${\alpha }_{z}d/\pi \approx 0.17$, and then increasing imaginary part.

The only bound mode for L-pol that travels without large decay is the mode ‘Proper 1’ in Figs. 5(a) and 5(b). This mode is bound, forward and physical (as described in Sec. 3.2); it presents growing imaginary part for increasing frequency, and can be excited in all its domain of physical existence, although it would be slowly decaying while traveling only when it has very small attenuation constant (e.g., near the light line at about $kd/\pi \approx 0.39$).

4.2 Leaky modes for transverse and longitudinal polarizations

The radiating leaky modes for T-pol are described next. The physical modes found in the array of silver nanospheres here analyzed radiate, but they cannot provide very directive radiation since the imaginary part of their wavenumber is large with respect to the host wavenumber k [48]. In Figs. 3 and 4, the mode ‘Proper 1’ is a leaky backward mode, whereas the mode ‘Improper 1’ is a leaky forward mode. We will show in Sec. 5 the possibility of modifying the imaginary part of the wavenumber of these two modes by varying geometrical parameters of the array.

The only radiating leaky mode for L-pol is the mode ‘Improper 1’ in Figs. 5 and 6. This mode is leaky and forward. However, notice that this mode lies very close to the light line; since the dipoles approximating the nanospheres are longitudinally polarized (i.e., along the z-direction), their far-field radiation pattern would have the classical $\sin \theta$ shape, with a far-field null corresponding to the direction along which the dipoles are aligned. Indeed, in the case of the analyzed linear chain, when ${\beta }_z=k$ the radiation should take place along the z-direction and this is not possible due to the aforementioned null, whereas when ${\beta }_z<k$ for ${\beta }_z\approx k$, a very small radiation (still due to the dipole $\sin \theta$ radiation pattern) could take place.

A detailed analysis of how to adopt these leaky modes to obtain directive radiation at optical frequencies has been recently shown in [7].

5. Parametric analysis of bound and leaky modes

In this section, we analyze the evolution of bound and leaky modes by varying the radius r of each nanosphere and the spatial period d of the array.

For the T-pol case, the dependence with respect to the array geometrical parameters of the backward ‘Proper 1’ (always physical, bound for ${\beta }_z<-k$, leaky for $-k<{\beta }_z<0$) and the forward ‘Improper 1’ (physical leaky for $0<{\beta }_z<k$, nonphysical bound for ${\beta }_z>k$) is illustrated in Fig. 7 . The linear chain under analysis here is the same as the one in Sec. 3 where we vary one parameter at the time. In Fig. 7(a), d = 75 nm and r varies from 15 nm to 35 nm; whereas in Fig. 7(b), r = 25 nm and d varies from 55 nm to 100 nm. It is observed in Fig. 7(a) that increasing the radius of each nanosphere to 35 nm (blue lines), both the ‘Proper 1’ and ‘Improper 1’ modes present lower losses for a wider frequency range than the array in Sec. 3 (green lines). Also, notice that by increasing the nanospheres’ radius, the mode ‘Proper 1’ is “shifted” towards lower frequencies, whereas the mode ‘Improper 1’ is “shifted” towards higher frequencies in Fig. 7(a). It is noteworthy that for r = 35 nm the backward proper mode has a very small imaginary part in correspondence of $kd/\pi \approx 0.38$ (where this mode is bound); if excited at this particular frequency, this mode would travel long distances from the excitation source. Moreover, ${\beta }_z$ is quite large, thus causing short plasmonic wavelengths. In Fig. 7(b), we show the mode dependence on the array period: it can be noticed that the imaginary part ${\alpha }_z$ increases with increasing periodicity. Notice that this is in agreement with what it has been previously shown in [20].

 

Fig. 7 Dependence of the ‘Proper 1’ bound mode (solid line) and ‘Improper 1’ leaky mode (dashed line) for T-pol in Fig. 3 with respect to the geometrical parameters of the linear array. The analyzed array is the same as in Sec. 3. In (a) d = 75 nm and r varies from 15 nm to 35 nm. In (b) r = 25 nm and d varies from 55 nm to 100 nm.

Download Full Size | PDF

In the L-pol case, the dependence with respect to the array geometrical parameters of the forward ‘Improper 1’ (physical leaky only for $0<{\beta }_z<k$, nonphysical bound for ${\beta }_z>k$) is illustrated in Fig. 8 . The performed parameterization is the same as in Fig. 7. As for the T-pol case previously analyzed, it can be observed in Fig. 8(a) that increasing the radius of each nanosphere to 35 nm (blue lines), the ‘Improper 1’ mode presents lower losses for a wider frequency range than the array in Sec. 3 (green lines). Furthermore, notice that for large nanospheres’ radius the frequency range for which the mode ’Improper 1’ is physical is wider (this range is about 0.01 for the red curve and about 0.21 for the blue curve, using the scale $kd/\pi$ in Fig. 8(a)) and “shifted” towards lower frequencies. In Fig. 8(b), we modify the periodicity of the array with respect to the case in Fig. 5, and it is noticed that the imaginary part ${\alpha }_z$ decreases for decreasing d. Again, this result is in agreement with what it has been previously shown in [20]. Also, notice that, similar to the case in Fig. 8(a), for small periodicity d, the frequency range for which the mode ’Improper 2’ is physical is wider than the other cases with larger periodicity.

 

Fig. 8 Parameter analysis as in Fig. 7, but for the ‘Improper 1’ leaky mode for L-pol in Fig. 5.

Download Full Size | PDF

6. Conclusion

This work provides an extensive and comprehensive description of modal analysis in a representative case of a 1D-periodic array of metallic nanospheres in terms of bound/leaky, proper/improper, physical/nonphysical wave conditions. The subsets of physical and nonphysical waves have been identified and the general principle of their identification has been provided. A physical wave is one that can be excited by a realistic localized source, and it is classified in terms of path deformation in the complex wavenumber domain. Though the theory is involved, a simple Table summarizing the classification in terms of longitudinal and transverse wavenumbers, based on spectral theory, has been given. This is the first time such a complete physical characterization in a chain of plasmonic nanospheres has been presented, and to the authors‘ knowledge, this analysis for complex modes is the most comprehensive carried out so far. The formulation and classification are general. The high subwavelength regime and resonance of the array nanospheres are responsible for the interesting modal dispersion described in this paper; for example, it is very interesting to note that the T-pol physical backward leaky wave “Proper 1” in Fig. 4, increasing frequency, transitions into a physical forward leaky wave ‘Improper 1’. Knowledge of complex mode propagation is useful for the development of innovative sensors, field confinement and field enhancement structures, and highly directive antennas. Mode diagrams can be engineered by modifying the array parameters, thus allowing the design of waveguiding structures by using physical bound modes (such as ‘Proper 2’ in T-pol, and ‘Proper 1’ in L-pol) and radiating structure by using physical leaky modes (such as ‘Improper 1’ in L-pol).

Appendix A: on the physical existence of modes excited by a localized source

Assume that the linear array of nanospheres is excited by a single dipolar source, ${\mathbf{P}}^i(\mathbf{r}\text{'})={\widehat{\mathbf{p}}}^i\delta \left(\mathbf{r}\text{'}-{\mathbf{r}}_0\right)$ at ${\mathbf{r}}_0=x_0\widehat{\mathbf{x}}+y_0\widehat{\mathbf{y}}+z_0\widehat{\mathbf{z}}$, where ${\widehat{\mathbf{p}}}^i$, having units [Cm], is the impressed unit-amplitude electric dipole source. The total field produced by ${\mathbf{P}}^i$ and the chain, at r, can be represented as

(8)$$\mathbf{E}(\mathbf{r},{\mathbf{r}}_0)=\frac{d}{2\pi }{\displaystyle \underset{-\infty }{\overset{\infty }{\int }}{e}^{i{k}_{z}z}}{\mathbf{e}}_0^{\infty }\left(x,y,{\mathbf{r}}_0,k_z\right)d{k}_z,$$

where ${\mathbf{\text{e}}}_0^{\infty }$ is the 0-th Floquet expansion term of the total field ${\mathbf{E}}^{\infty }\left(\mathbf{r},{\mathbf{r}}_0,k_z\right)={\displaystyle {\sum }_p{e}^{i{k}_{z,p}z}}{\mathbf{e}}_p^{\infty }\left(x,y,{\mathbf{r}}_0,k_z\right)$, generated by the chain excited by a phased array of point sources ${\mathbf{p}}^{i,\infty }\left({\mathbf{r}}^{\prime },k_z\right)={\widehat{\mathbf{p}}}^i{\displaystyle {\sum }_{n=-\infty }^{\infty }\delta \left({\mathbf{r}}^{\prime }-{\mathbf{r}}_0-nd\widehat{\mathbf{z}}\right)}{e}^{i{k}_{z}nd}$ (phased replicas of the real source ${\mathbf{\text{P}}}^i$ at ${\mathbf{\text{r}}}_0$) according to the array scanning method (ASM), outlined in [41,44,45], where a reader is referred to also for a direct expression of ${\mathbf{\text{e}}}_0^{\infty }$. The integrand term ${\mathbf{\text{e}}}_0^{\infty }$ possesses spectral singularities: poles and branch points. Branch points are in correspondence of $k_{zb,p}=\pm k-2\pi p/d$ such that$k_{\rho ,p}=\sqrt{k^2-k_{zb,p}^2}=0$. The term with p = 0 is within the fundamental BZ.

Figure 9 shows the path deformation onto the steepest descent paths (SDPs) and the capturable regions for poles representing physical modes in the complex $k_z$ plane when an observer is along positive z. For an observer along negative z, the path on the real axis should be deformed onto the lower half of the complex $k_z$ plane. These SDPs, in correspondence of branch points at $k_{zb,p}=k-2\pi p/d$, represent the “space wave” or continuous spectrum terms, and a precise discussion of these terms is not within the scope of this paper (more details in [41,42,45]). The thick dark blue lines in Fig. 9 represent branch cuts, and when the integration path crosses a cut, a change in the Riemann sheet is performed (i.e., from proper to improper sheet and vice versa). The residues of the captured poles represent the physical waves (those that can be excited by a source ${\mathbf{\text{P}}}^i$ at ${\mathbf{\text{r}}}_0$). Note that each mode is represented by an infinite set of poles (Floquet spatial harmonics) periodic along the $\mathrm{Re}\left(k_z\right)$ axis. As can be noticed in Fig. 9, looking only at the spatial harmonics in the fundamental BZ, for an observer along the positive z-direction, only poles relative to backward proper (i.e., belonging to the top Riemann sheet) modes (leaky or bound) can be captured in the II quadrant (so called region R1) of the complex $k_z$ plane. In the I quadrant, instead, only poles relative to forward proper bound (region R3, bound region, below the light line in the $k\beta \text{ plane}$) and forward improper (i.e., belonging to the bottom Riemann sheet) leaky modes (region R2, radiating region, above the light line in the $k\beta \text{ plane}$) can be captured (see for example Appendix A in [41] or [42] for more details regarding spectral singularities and spectral path deformation).

 

Fig. 9 Complex $k_z$ wavenumber plane (assuming $k<\pi /d$). The field produced by a dipole source is represented by the integration path along the real axis. For an observer along the positive z axis, the path can be deformed in the upper half space. The poles encountered represent waves (each mode is represented by an infinite set of residues). The SDPs around the branch points at $k_z=k_{zb,p}=k-2\pi p/d$ represent the “space wave” or continuous spectrum terms [41,44,45,49]. Physical modes (only the harmonic in the fundamental BZ, i.e., p = 0, is indicated) are those relative to poles captured in the three indicated regions (R1, R2 and R3). In each region, it is indicated the type of physical pole that can be captured (physical wave harmonic). PB = proper bound, PL = proper leaky, and IL = improper leaky wave.

Download Full Size | PDF

The procedure outlined above leads to the following field representation for the total field at r excited by a source at $r_0$, as a sum of two wave species:

(9)$$\mathbf{\text{E}}(\mathbf{r},{\mathbf{r}}_0)={\displaystyle \sum _n}{\mathbf{\text{E}}}_n^{\mathrm{mode}}(\mathbf{r},{\mathbf{r}}_0,k_z^n)+{\mathbf{\text{E}}}_{tot}^{\mathrm{sp}}(\mathbf{r},{\mathbf{r}}_0).$$

Here, ${\mathbf{\text{E}}}_{tot}^{\mathrm{sp}}(\mathbf{r},{\mathbf{r}}_0)$denotes a “spatial wave” collecting the contributions from all the SDPs, and ${\mathbf{\text{E}}}_n^{\mathrm{mode}}(\mathbf{r},{\mathbf{r}}_0,k_z^n)$ corresponds to the n-th mode with wavenumber $k_z^n$ [41,44,45], which is expressed as in Eq. (7) and its characterization is the main purpose of this paper. Note that the physical property of the modes depends on the location of the corresponding pole and its capability of being captured. In summary physical waves are classified according to their wavenumbers as in Table 1.

Appendix B: Dyadic Ewald representation of the periodic GF for linear chains

The electric-field dyadic Green’s function ${\underset{\_}{\stackrel{⌣}{G}}}^{\infty }\left(k_z\right)$ for the phased periodic array of nanospheres in Eq. (4) is represented according to the Ewald method [32–34] as

(10)$${\underset{\_}{\stackrel{⌣}{G}}}^{\infty }\left(k_z\right)={\underset{\_}{G}}_{\text{spectral}}\left(k_z\right)+{\underset{\_}{\stackrel{⌣}{G}}}_{\text{spatial}}\left(k_z\right),$$

where the spectral and the regularized spatial parts are

(11)$${\underset{\_}{G}}_{\text{spectral}}\left(k_z\right)=\frac{1}{{\epsilon }_h{\epsilon }_0}\left[k^2{G}_{\text{spectral}}\left(k_z\right)\underset{\_}{I}+{\underset{\_}{H}}_{\text{spectral}}\left(k_z\right)\right],$$

(12)$${\underset{\_}{\stackrel{⌣}{G}}}_{\text{spatial}}\left(k_z\right)=\frac{1}{{\epsilon }_h{\epsilon }_0}\left[k^2{\stackrel{⌣}{G}}_{\text{spatial}}\left(k_z\right)\underset{\_}{I}+{\underset{\_}{\stackrel{⌣}{H}}}_{\text{spatial}}\left(k_z\right)\right].$$

The Ewald representation for the periodic scalar GF [32,33] (specialized for ρ = 0 and z = 0, and with regularized spatial term) is

(13)$$G_{{}^{\text{spectral}}}\left(k_z\right)=\frac{1}{4\pi d}{\displaystyle \sum _{p=-\infty }^{\infty }{E}_1^G\left(-\frac{k_{\rho ,p}^2}{4E^2}\right)},$$

(14)$${\stackrel{⌣}{G}}_{{}^{\text{spatial}}}\left(k_z\right)=\frac{1}{8\pi }{\displaystyle \sum _{\begin{array}{l}n=-\infty \\ \left(n\ne 0\right)\end{array}}^{\infty }\frac{e^{i{k}_{z}nd}}{r_n}f\left(r_n\right)}+\frac{f^{\prime }\left(0\right)-2ik}{8\pi },$$

with

(15)$$f\left(r_n\right)=e^{-ik{r}_n}\text{erfc}\left({\beta }_{-}\right)+e^{ik{r}_n}\text{erfc}\left({\beta }_{+}\right),{\beta }_{\pm }=r_{n}E\pm i\frac{k}{2E},$$

where $f^{\prime }$is the prime derivative of f with respect to its argument and E is the Ewald parameter as in [32,33], which is in general chosen as $E=\sqrt{\pi }/d$ (other choices are discussed in [32]). Furthermore, $E_1^G\left(\nu \right)=E_1\left(\nu \right)-i2\pi m$ is the 1st-order generalized exponential integral, where $E_1\left(\nu \right)={\displaystyle {\int }_{\nu }^{\infty }{t}^{-1}{e}^{-t}dt}$, $-\pi \le \arg \nu \le \pi$ is the standard exponential integral [33] and $m=0/\pm 1$ is an integer, associated to the proper/improper character of the cylindrical p-th harmonic. In this particular implementation of the GF, $r_n=\left|n\right|d$ is the distance between the point where we evaluate the GF (at $r=0$) and the n-th dipolar source at ${\mathbf{r}}_n^{\prime }=nd\widehat{z}$, and $\text{erfc}\left(\nu \right)=2{\pi }^{-1/2}{\displaystyle {\int }_{\nu }^{\infty }{e}^{-t^2}dt}$ is the complementary error function [32,33]. Moreover

(16)$${\underset{\_}{H}}_{\text{spectral}}\left(k_z\right)=\frac{1}{4\pi d}{\displaystyle \sum _{p=-\infty }^{\infty }{\underset{\_}{F}}_{{}_{\text{spectral,}p}}\left(k_z\right)},$$

(17)$${\underset{\_}{\stackrel{⌣}{H}}}_{\text{spatial}}\left(k_z\right)=\frac{1}{8\pi }{\displaystyle \sum _{\begin{array}{l}n=-\infty \\ \left(n\ne 0\right)\end{array}}^{\infty }{e}^{i{k}_{z}nd}{\underset{\_}{F}}_{{}_{\text{spatial,}n}}\left(k_z\right)+}\frac{1}{8\pi }\frac{f^{‴}\left(0\right)+2i{k}^3}{3}\underset{\_}{I},$$

where

(18)$${\underset{\_}{F}}_{{}_{\text{spectral,}p}}\left(k_z\right)=\left[-2E^2{E}_2^G\left(\frac{-k_{\rho ,p}^2}{4E^2}\right)\left(\underset{\_}{I}-\widehat{z}\widehat{z}\right)-E_1^G\left(\frac{-k_{\rho ,p}^2}{4E^2}\right)k_{z,p}^2\widehat{z}\widehat{z}\right],$$

(19)$${\underset{\_}{F}}_{{}_{\text{spatial,}n}}\left(k_z\right)=\left(\frac{f^{\prime }\left(r_n\right)}{r_n^2}-\frac{f\left(r_n\right)}{r_n^3}\right)\underset{\_}{I}+\left(\frac{f^{″}\left(r_n\right)}{r_n}-\frac{3f^{\prime }\left(r_n\right)}{r_n^2}+\frac{3f\left(r_n\right)}{r_n^3}\right)\widehat{z}\widehat{z},$$

and$f^{″}$ and $f^{‴}$ are the second and third derivatives of f in Eq. (15) with respect to its argument, respectively, and $E_2^G\left(\nu \right)=e^{-\nu }-\nu E_1^G\left(\nu \right)$.

The discussion that follows is related to the choice of the m value in the generalized exponential integral $E_1^G\left(\nu \right)$. We first recall that a purely spectral representation of the periodic Green’s function is in terms of outgoing cylindrical waves [32] described by the Hankel function of the first kind $H_0^{\text{(1)}}$; therefore, to discern between proper and improper Riemann sheets of the $k_z$ plane, we need to analyze the properties of the Hankel function with argument $k_{\rho ,p}\rho$, describing the p-th spectral harmonic, whose asymptotic expansion for large argument is given by

(20)$$H_0^{\text{(1)}}\left(k_{\rho ,p}\rho \right)\sim \sqrt{\frac{2}{\pi k_{\rho ,p}\rho }}{e}^{i\left(k_{\rho ,p}\rho -\frac{\pi }{4}\right)}\sim \frac{e^{i{k}_{\rho ,p}\rho }}{\sqrt{k_{\rho ,p}\rho }}\left(-\pi <\text{arg }{k}_{\rho ,p}<2\pi \right),$$

as described in [40], page 364, where also the argument range $\left(-\pi ,2\pi \right)$ is provided.

From Eq. (20), considering $k_{\rho ,p}={\beta }_{\rho ,p}+i{\alpha }_{\rho ,p}$ and assuming that the general argument of the Hankel function, including argument phases in other Riemann sheets, is dictated by $k_{\rho ,p}=\xi e^{m\pi i}$ with $0<\text{arg }\xi <\pi$, we identify the following mapping (the first two conditions have been previously given in [33])

(21)$$\{\begin{array}{c}-\pi <\arg \left(\xi e^{m\pi i}\right)<0,m=-1\text{ Improper}\\ 0<\arg \left(\xi e^{m\pi i}\right)<\pi ,m=0\text{ Proper }\\ \pi <\arg \left(\xi e^{m\pi i}\right)<2\pi ,m=+1\text{ Improper}\end{array}.$$

The indexes m = 0 and −1 refer to the principal Riemann sheet of the Hankel function [40] with argument $k_{\rho ,p}\rho =\xi e^{m\pi i}\rho$ decaying (proper) and growing (improper) radially. The choice with m = +1 in Eq. (21) is relative to cylindrical waves growing radially (improper). Notice that when analyzing the FW harmonic with modal longitudinal wavenumber $k_z$ in the fundamental BZ, i.e., with p = 0, improper modes are those with $-\pi <\arg k_{\rho ,0}<0$ or $\pi <\arg k_{\rho ,0}<2\pi$, and are thus computed by using m = −1 or m = +1 only for p = 0, as it has been described in [41,45] for planar 1D-periodic structures. Accordingly, when the wavenumber $k_z$ is in the fundamental BZ, all other radial wavenumbers $k_{\rho ,p}$ with p ≠ 0 in Eq. (13) should be always such that $0<\arg k_{\rho ,p}<\pi$ (proper choice), and thus evaluated with m = 0. Furthermore, for m values other than 0 and −1, according to the formula $\pi m<\arg k_{\rho ,p}<\pi (m+1)$, the $k_z$ solutions are always found to be improper according to the analytic continuation of the Hankel function and its asymptotic behavior for large argument [40]

(22)$$\begin{array}{l}{H}_0^{\text{(1)}}\left(k_{\rho ,p}\rho \right)=H_0^{\text{(1)}}\left(\xi \rho \right)-2m{J}_0\left(\xi \rho \right)\\ \sim \frac{e^{i\text{Re}\left(\xi \right)\rho }{e}^{-\text{Im}\left(\xi \right)\rho }}{\sqrt{\xi \rho }}\left(1-m\right)-m\frac{e^{-i\text{Re}\left(\xi \right)\rho }{e}^{\text{Im}\left(\xi \right)\rho }}{\sqrt{\xi \rho }},\end{array}$$

where each index $m\ne 0,-1$ refers to non-principal Riemann sheets of the Hankel function.

The physical validity of these last improper solutions related to non-principal Riemann sheets of the Hankel function involved in cylindrical problems need further analysis, therefore they are not shown in this paper and will be the focus of future work. Nevertheless, here we present a detailed physical description of the modes previously reported (for example in [21,26]), distinguishing between proper and improper ones, forward and backward, bound and leaky, in a well defined framework. We would like to point out that the study of longitudinal complex modal propagation in open or non-uniform cylindrical domains, like a dielectric rod or a waveguide with an inner dielectric rod (simpler cases than the present one because they are not periodic), is an old subject [50–53]. However, not much has been done in the study of modal solutions in cylindrical problems involving other Riemann sheets of the Hankel function. Several papers have been discussing the physical validity (individual excitability) of complex modes, even reaching contradicting conclusions [50–53], as discussed in [54]. Recently, the authors of [54] claimed they sorted out old misconceptions regarding the physical validity of certain complex modes. Here we have provided a complete description regarding the principal Riemann sheet of the Hankel function but as already said further studies are required to analyze the solutions obtained from non-principal Riemann sheets. These are necessary to establish the wavenumber symmetry in the $k_z$ complex plane, and indeed, according to [55], the complex conjugate symmetry of improper solutions can be obtained only considering solutions in adjacent Riemann sheets of the Hankel function (i.e., with different m values). What described in [55] has been confirmed here for the mode ‘Improper 1’ in T-pol in Figs. 3 and 4 (where we also found the wavenumber close to the complex conjugate solution correspondent to the mode ‘Improper 2’ by using m = +1, not shown), and for the mode ‘Improper 1’ in L-pol in the I and III quadrants of the complex $k_z$ plane in Figs. 5 and 6 (where we also found the wavenumber close to the complex conjugate solution correspondent to the mode ‘Improper 1’ in the II and IV quadrants of the complex $k_z$ plane by using m = +1, not shown). In this paper, in accordance with [55], symmetry of the proper solutions has been found in the proper plane (i.e., with m = 0).

Acknowledgments

Support by NSF Award #CMMI-1101074 is acknowledged. The authors would also like to thank Prof. D.R. Jackson, University of Houston, TX, Dr. R. Shore and Dr. A. Yaghjian, HANSCOM-AFRL, MA, for useful discussions.

References and links

1. P. Alitalo, C. Simovski, A. Viitanen, and S. Tretyakov, “Near-field enhancement and subwavelength imaging in the optical region using a pair of two-dimensional arrays of metal nanospheres,” Phys. Rev. B 74(23), 235425 (2006). [CrossRef]  

2. S. Steshenko, F. Capolino, S. A. Tretyakov, and C. R. Simovski, “Super-Resolution and Near-Field Enhancement with Layers of Resonant Arrays of Nanoparticles,” in Applications of Metamaterials, F. Capolino, ed. (CRC Press, 2009), p. 4.1.

3. C. R. Simovski, A. J. Viitanen, and S. A. Tretyakov, “Sub-wavelength resolution in linear arrays of plasmonic particles,” J. Appl. Phys. 101(12), 123102 (2007). [CrossRef]  

4. C. Simovski, S. Tretyakov, and A. Viitanen, “Subwavelength imaging in a superlens of plasmon nanospheres,” Tech. Phys. Lett. 33(3), 264–266 (2007). [CrossRef]  

5. S. Steshenko, F. Capolino, P. Alitalo, and S. A. Tretyakov, “Effective model and investigation of the near-field enhancement and subwavelength imaging properties of multilayer arrays of plasmonic nanospheres,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 84(1), 016607 (2011). [CrossRef]   [PubMed]  

6. A. F. Koenderink, “Plasmon nanoparticle array waveguides for single photon and single plasmon sources,” Nano Lett. 9(12), 4228–4233 (2009). [CrossRef]   [PubMed]  

7. X. Liu and A. Alu, “Subwavelength leaky-wave optical nanoantennas: directive radiation from linear arrays of plasmonic nanoparticles,” Phys. Rev. B 82(14), 144305 (2010). [CrossRef]  

8. J. Beermann, S. M. Novikov, K. Leosson, and S. I. Bozhevolnyi, “Surface enhanced Raman imaging: periodic arrays and individual metal nanoparticles,” Opt. Express 17(15), 12698–12705 (2009). [CrossRef]   [PubMed]  

9. F. Liu, Z. Cao, C. Tang, L. Chen, and Z. Wang, “Ultrathin diamond-like carbon film coated silver nanoparticles-based substrates for surface-enhanced Raman spectroscopy,” ACS Nano 4(5), 2643–2648 (2010). [CrossRef]   [PubMed]  

10. I. Firkowska, S. Giannona, J. A. Rojas-Chapana, K. Luecke, O. Brustle, and M. Giersig, “Biocompatible Nanomaterials and Nanodevices Promising for Biomedical Applications ” in Nanomaterials for Application in Medicine and Biology, M. Giersig, and G. B. Khomutov, eds. (Springer, Berlin, 2008), p. I.1.

11. A. J. Haes and R. P. Van Duyne, “A nanoscale optical biosensor: sensitivity and selectivity of an approach based on the localized surface plasmon resonance spectroscopy of triangular silver nanoparticles,” J. Am. Chem. Soc. 124(35), 10596–10604 (2002). [CrossRef]   [PubMed]  

12. I. E. Sendroiu and R. M. Corn, “Nanoparticle diffraction gratings for DNA detection on photopatterned glass substrates,” Biointerphases 3(3), FD23–FD29 (2008). [CrossRef]   [PubMed]  

13. A. L. Fructos, S. Campione, F. Capolino, and F. Mesa, “Characterization of complex plasmonic modes in two-dimensional periodic arrays of metal nanospheres,” J. Opt. Soc. Am. B 28(6), 1446–1458 (2011). [CrossRef]  

14. S. A. Maier, P. G. Kik, and H. A. Atwater, “Optical pulse propagation in metal nanoparticle chain waveguides,” Phys. Rev. B 67(20), 205402 (2003). [CrossRef]  

15. S. Y. Park and D. Stroud, “Surface-plasmon dispersion relations in chains of metallic nanoparticles: an exact quasistatic calculation,” Phys. Rev. B 69(12), 125418 (2004). [CrossRef]  

16. R. A. Shore and A. D. Yaghjian, “Travelling electromagnetic waves on linear periodic arrays of lossless spheres,” Electron. Lett. 41(10), 578–580 (2005). [CrossRef]  

17. R. A. Shore and A. D. Yaghjian, ““Traveling electromagnetic waves on linear periodic arrays of lossless penetrable spheres,” IEICE Trans. Commun. E88B(6), 2346–2352 (2005). [CrossRef]  

18. D. S. Citrin, “Plasmon-polariton transport in metal-nanoparticle chains embedded in a gain medium,” Opt. Lett. 31(1), 98–100 (2006). [CrossRef]   [PubMed]  

19. A. F. Koenderink and A. Polman, “Complex response and polariton-like dispersion splitting in periodic metal nanoparticle chains,” Phys. Rev. B 74(3), 033402 (2006). [CrossRef]  

20. A. Alù and N. Engheta, “Theory of linear chains of metamaterial/plasmonic particles as subdiffraction optical nanotransmission lines,” Phys. Rev. B 74(20), 205436 (2006). [CrossRef]  

21. R. A. Shore and A. D. Yaghjian, “Complex Waves on 1D, 2D, and 3D Periodic Arrays of Lossy and Lossless Magnetodielectric Spheres,” AFRL-RY-HS-TR-2010–0019, (Air Force Research Laboratory, Hanscom, MA 2010).

22. K. B. Crozier, E. Togan, E. Simsek, and T. Yang, “Experimental measurement of the dispersion relations of the surface plasmon modes of metal nanoparticle chains,” Opt. Express 15(26), 17482–17493 (2007). [CrossRef]   [PubMed]  

23. T. Yang and K. B. Crozier, “Dispersion and extinction of surface plasmons in an array of gold nanoparticle chains: influence of the air/glass interface,” Opt. Express 16(12), 8570–8580 (2008). [CrossRef]   [PubMed]  

24. A. Alù and N. Engheta, “Effect of small random disorders and imperfections on the performance of arrays of plasmonic nanoparticles,” N. J. Phys. 12(1), 013015 (2010). [CrossRef]  

25. S. Campione and F. Capolino, “Linear and Planar Periodic Arrays of Metallic Nanospheres: Fabrication, Optical Properties and Applications,” in Selected Topics in Metamaterials and Photonic Crystals, A. Andreone, A. Cusano, A. Cutolo, and V. Galdi, eds. (World Scientific Publishing, 2011), pp. 141–194.

26. M. Conforti and M. Guasoni, “Dispersive properties of linear chains of lossy metal nanoparticles,” J. Opt. Soc. Am. B 27(8), 1576–1582 (2010). [CrossRef]  

27. C. F. Bohren and D. R. Huffman, Absorption and Scattering of Light by Small Particles (Wiley, 1983).

28. S. Steshenko and F. Capolino, “Single Dipole Approximation for Modeling Collections of Nanoscatterers,” in Theory and Phenomena of Metamaterials, F. Capolino, ed. (CRC Press, 2009), p. 8.1.

29. D. E. Muller, “A Method for Solving Algebraic Equations Using an Automatic Computer,” Math. Tables Other Aids Comput. 10(56), 208–215 (1956). [CrossRef]  

30. W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery, Numerical recipes: the art of scientific computing (Cambridge University Press, 2007).

31. “IMSL Fortran Numerical Library,” (Visual Numerics Corporate Headquarters, 2500 Wilcrest Drive, Suite 200, Houston, TX), www.vni.com.

32. F. Capolino, D. R. Wilton, and W. A. Johnson, “Efficient computation of the 3D Green's function for the Helmholtz operator for a linear array of point sources using the Ewald method,” J. Comput. Phys. 223(1), 250–261 (2007). [CrossRef]  

33. V. R. Komanduri, F. Capolino, D. R. Jackson, and D. R. Wilton, “Computation of the one-dimensional free-space periodic green's function for leaky waves using the ewald method,” in Proc. URSI Gen. Ass.(Chicago, IL, 2008).

34. F. T. Celepcikay, D. R. Wilton, D. R. Jackson, and F. Capolino, “Choosing splitting parameters and summation limits in the numerical evaluation of 1-D and 2-D periodic Green's functions using the Ewald method,” Radio Sci. 43(6), RS6S01 (2008). [CrossRef]  

35. J. D. Jackson, Classical Electrodynamics (Wiley, 1998).

36. V. A. Markel, V. N. Pustovit, S. V. Karpov, A. V. Obuschenko, V. S. Gerasimov, and I. L. Isaev, “Electromagnetic density of states and absorption of radiation by aggregates of nanospheres with multipole interactions,” Phys. Rev. B 70(5), 054202 (2004). [CrossRef]  

37. D. W. Mackowski, “Calculation of total cross sections of multiple-sphere clusters,” J. Opt. Soc. Am. A 11(11), 2851–2861 (1994). [CrossRef]  

38. R. Ruppin, “Evaluation of extended Maxwell-Garnett theories,” Opt. Commun. 182(4–6), 273–279 (2000). [CrossRef]  

39. W. T. Doyle, “Optical properties of a suspension of metal spheres,” Phys. Rev. B Condens. Matter 39(14), 9852–9858 (1989). [CrossRef]   [PubMed]  

40. M. Abramowitz and I. A. Stegun, Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables (Dover Publications, 1965).

41. F. Capolino, D. R. Jackson, and D. R. Wilton, “Field Representations in Periodic Artificial Materials Excited by a Source,” in Theory and Phenomena of Metamaterials, F. Capolino, ed. (CRC Press, 2009), p. 12.1.

42. L. B. Felsen and N. Marcuvitz, Radiation and Scattering of Waves (IEEE Press, 1994).

43. P. Baccarelli, S. Paulotto, and C. Di Nallo, “Full-wave analysis of bound and leaky modes propagating along 2D periodic printed structures with arbitrary metallisation in the unit cell,” IET Proc. Microwaves, Antennas Propag. 1(1), 217–225 (2007). [CrossRef]  

44. F. Capolino, D. R. Jackson, D. R. Wilton, and L. B. Felsen, “Comparison of methods for calculating the field excited by a dipole near a 2-D periodic material,” IEEE Trans. Antenn. Propag. 55(6), 1644–1655 (2007). [CrossRef]  

45. F. Capolino, D. R. Jackson, and D. R. Wilton, “Fundamental properties of the field at the interface between air and a periodic artificial material excited by a line source,” IEEE Trans. Antenn. Propag. 53(1), 91–99 (2005). [CrossRef]  

46. A. Alù, A. Salandrino, and N. Engheta, “Negative effective permeability and left-handed materials at optical frequencies,” Opt. Express 14(4), 1557–1567 (2006). [CrossRef]   [PubMed]  

47. I. El-Kady, M. M. Sigalas, R. Biswas, K. M. Ho, and C. M. Soukoulis, “Metallic photonic crystals at optical wavelengths,” Phys. Rev. B 62(23), 15299–15302 (2000). [CrossRef]  

48. A. A. Oliner and D. R. Jackson, “Leaky-wave antennas,” in Antenna Engineering Handbook, J. Volakis, ed. (McGraw Hill, 2007), p. 11.1.

49. F. Capolino, D. R. Jackson, and D. R. Wilton, “Mode excitation from sources in two-dimensional EBG waveguides using the array scanning method,” IEEE Microw. Wirel. Compon. Lett. 15(2), 49–51 (2005). [CrossRef]  

50. P. J. B. Clarricoats and K. R. Slinn, “Complex modes of propagation in dielectric-loaded circular waveguide,” Electron. Lett. 1(5), 145–146 (1965). [CrossRef]  

51. J. D. Rhodes, “General constraints on propagation characteristics of electromagnetic waves in uniform inhomogeneous waveguides,” in Proc. Inst. Electr. Eng. 118, 849–856 (1971).

52. T. F. Jablonski, “Complex modes in open lossless dielectric waveguides,” J. Opt. Soc. Am. A 11(4), 1272–1282 (1994). [CrossRef]  

53. T. Rozzi, L. Pierantoni, and M. Farina, “General constraints on the propagation of complex waves in closed lossless isotropic waveguides,” IEEE Trans. Microw. Theory Tech. 46(5), 512–516 (1998). [CrossRef]  

54. R. Islam and G. V. Eleftheriades, “On the Independence of the Excitation of Complex Modes in Isotropic Structures,” IEEE Trans. Antenn. Propag. 58(5), 1567–1578 (2010). [CrossRef]  

55. J. R. James, “Leaky waves on a dielectric rod,” Electron. Lett. 5(11), 252–254 (1969). [CrossRef]