1. Introduction

The history of plasmonics can be traced to the early 1900s with the description of radio wave propagation along the surface of a conductor [1–4]. Early work by pioneers in the field have paved the way for the current, rapid growth of research into plasmonics [5–11]. While the field has continued to grow, it has been developed and refined with respect to a widening field of potential applications including microscopy, computing, and biomedicine [3,4,12–15].

Symmetric thin-film waveguides have been extensively studied for decades and reviews may be found in [3,4]. The field of symmetric thin-film waveguides is presently experiencing a rapid increase in the development of practical applications such as optical circuits, biomedical sensing, and photovoltaics [14,16–22]. Symmetric thin-film waveguides are noted for high power confinement while asymmetric thin-film designs have a nonuniform confinement. On the other hand, asymmetric structures have smaller attenuation as films become thinner and may be practical within the optical frequency range [23–25].

Analytic results, constructed in an ad-hoc manner and therefore less generally applicable, may be found in earlier works on thin-film structures [26–32]. To date, most research into thin-film waveguides has been focused on the development of numerical methods of solution to the dispersion relations [33–38].

In this work, a generic methodology for the construction of closed-form approximations is developed and applied to three- and multi-layer symmetric designs. Closed-form results are found that describe the physics of the propagation and confinement of plasmons in the thin-film regime. The closed-form results are referenced, for convenience, to that of Dionne et al. ([39–42]). The results expose the rich dependence of the propagation length and attenuation on material properties, wavelength, and thin-film physical characteristics. The parametric basis for improved three-layer, MIM (metal-insulator-metal) propagation and attenuation properties over its complement IMI (insulator-metal-insulator) is fully elucidated. Numerical results in [39–42] and elsewhere [41,43,44] are reconstructed using the closed-form analysis. For multilayer designs, in direct contrast to three-layer structures, it is found that the propagation and attenuation in the space of the cladding thickness and wavelength are governed by proximity to singularities associated with the zero limit of cladding dielectric constant imaginary part. Enhanced propagation and attenuation at desired ranges of wavelength are uncovered at specific material and thin-film structural properties. The overall results provide a generic approach to the design of symmetric thin-film waveguides and are aimed at providing parametric estimates that yield desirable properties over given frequency ranges. These results may also be applied as general design and material guidelines and extended to consideration of device construction in applications such as waveguides, and hyperlenses [45].

2. Symmetric three-layer dispersion equations

Methods for the closed-form analysis of the symmetric thin-film, three-layer dispersion equations ([39] and [40]) are developed for the IMI (insulator-metal-insulator) configuration. The approach is then applied to the analysis of MIM and comparison is made between IMI and MIM device performance. Finally, symmetric thin-film multi-layer devices such as MIMIMIM are analyzed using the method and the sources of fundamental differences in behaviour of three-layer and multi-layer devices are examined. To fix ideas, the dielectric is assumed to be $\textrm {SiO}_{\textrm {2}}$ and the cladding is $\textrm {Ag}$ which follows ([39] and [40]) and this is without loss of generality regarding this class of materials.

2.1 Dimensional dispersion equations

The IMI dispersion equations, following [39], are

2.2 Dimensionless dispersion equations

In this section the dimensionless form of the dispersion equations is found and reasoning for the choice of dimensionless variables and their relevance to the analysis is provided.

The dispersion relations, (1) and (2), are non-dimensionalized and reduced to a minimal dependence on two, dimensionless groups and a single dependent variable.

The waveguide vector component Eq. (3) may be recast as

The thickness is rescaled as

A dimensionless axial waveguide component variable is defined as

The difference between the dielectric constants $\tilde \epsilon _1$ and $\tilde \epsilon _2$ mediates the relationship between the metal and insulator axial waveguide components $k'_{z_1}$ and $k'_{z_2}$. From (4), with (6), (8), and (9),

The waveguide axial variable relation in (16) is further simplified by removing $k_{z_1}'$ in favor of the surrogate $z$ (12) and substituting from (13) to (16) yields

The waveguide component equation in (18) is substituted to the dispersion relations in (1) and (2) along with the definitions in (8),  (12), and (13), and reduces to

The form in (20) is further reduced to a minimum dependence on the dependent variable $z$ and two dimensionless groups via

For the IMI configuration, the scaled dielectric difference $\alpha$ is small compared with unity near peaks in propagation length for the IMI device and materials considered here (not shown) and it is useful to introduce the order one variable

Applying (21), (22), and (23) to (19) and (20) yields the final form

The variable $Z$ is $O(1)$ (’order one’) for the IMI configuration and the balance between the scale of the first and second terms is explicitly controlled by the parameters $\alpha$ and $\tilde \epsilon$. This form will also be seen to be convenient starting point for the analysis of the three-layer, MIM configuration.

The focus here will be on the analysis of the asymmetric mode $L_+$ (24) which provides the main results [39] and [40] and methods applied here to $L_+$ are easily transferred to $L_-$ (25).

2.3 Analysis of IMI $L_{+}$ (Eq. (24))

The goal, as stated earlier with respect to the surface plasmon propagation length and its skin depth, is to to build closed-form approximations to these in regions of greatest practical interest.

The numerical method is briefly summarized for clarity and follows a stable iterative scheme with arbitrary solution accuracy like that in [36]. Equation (24) is rewritten as

The numerical solution of (26) is found by constructing a coupled pair of equations whose iterative solution $Z$ is a stable fixed point. Given the function

The pair (27) and (29) is an exact representation of the full Eq. (24) or its equivalent (26) and provide an iterative set

Note that the coupled pair developed above are equivalent to the iteration

A qualitative estimate $Z_0$ begins by noting that the scaled dielectric difference $\alpha$ is small compared with unity and the dielectric ratio is treated as $\tilde \epsilon = O(1)$. Thus, a balance between the first and second terms in (26) or (28) requires $Z \approx 1$ in order to reduce the magnitude of the $\sqrt {Z^2-1}$ in the first term and thus $Z_0=1$. This qualitative estimate of $Z_0$ is formalized as the zeroth term of a regular perturbation [46] of $Z$ in the parameter $\alpha$ used below in the closed-form approximation.

2.4 Closed-form approximation

The required balance between the two terms in (26) for small $\alpha$ implies that $Z$ may be expanded in (26) as a regular perturbation in $\alpha$; $Z = Z_0 + \alpha Z_1 + \alpha ^2 Z_2 + \cdots$ and the two-term expansion is

The propagation length $L = \lambda /2|\Im (k_x')|$, given $k_x'$ (14), $k_{z_1}'$ (13), and $Z$ (23) leads to

From the definition of $\Delta \tilde \epsilon$ (17), and the scaled dielectric difference $\alpha$ (22), and noting that $\tilde \epsilon _1 + \alpha ^2/d'^2 = \tilde \epsilon _2$

Further expanding for small $\alpha$ and replacing $\alpha ^2 = d'^2(\tilde \epsilon _2 - \tilde \epsilon _1)$ from (17) and (22)

Since the dielectric constant $\tilde \epsilon _2$ is real, the imaginary part of the wave component $k_x'$ resides in the second term. After some algebra, the propagation length $L = \lambda /2|\Im (k_x')|$ may be reduced to the dimensional form

2.5 Discussion

It is desirable to achieve both large propagation length and small attenuation. A primary difficulty with achieving both these aims, in the IMI configuration, is immediately apparent from the propagation length $L \propto \lambda ^3/d^2$ (38) and the skin depth $P \propto \lambda ^2/d$ (40). These demonstrate both quantities have correlated behaviour in the wavelength $\lambda$ and insulator thickness $d$: increasing or decreasing the wavelength or the depth always improves one feature at the expense of the other. Note that the quasi-bound results (not shown) [39] near $\lambda =300nm$ are recovered from (40).

At larger wavelengths, the metal dielectric constant becomes large with a real part that is an order of magnitude larger than its imaginary part. On the other hand, the insulator dielectric constant $\tilde \epsilon _2$ is order one. Thus, at larger $\lambda$ the propagation length is

Note that although the size of $\Re (\tilde \epsilon _1)$ may be generally available to limit the skin depth, the ratio $\lambda / \Re (\tilde \epsilon _1) \gg 1$ for the metal [39] considered here.

A complementary set of possibilities unfold below for the MIM device. This is due to a swapping of roles in $\tilde \epsilon _1$ and $\tilde \epsilon _2$ with a dependence on a re-defined dielectric ratio $\tilde \epsilon = \tilde \epsilon _2/\tilde \epsilon _1$. This ratio is now small in magnitude while the scaled dielectric difference magnitude is unchanged at large wavelengths.

3. Analysis of MIM $L_+$ (Eq. (24))

The analysis of the MIM configuration only requires the interchanging of $\tilde \epsilon _1$ and $\tilde \epsilon _2$ [40] in the dispersion and waveguide equations ((1)–(3)) while the subscripts ’1’ and ’2’ respectively continue to refer to the metal and the insulator. The analysis is focused here on the anti-symmetric mode at large wavelength as was the case in [40]. Notational simplicity of the dimensionless problems across the IMI, MIM, and multi-layer configurations is maintained by re-stating dimensionless variables for each configuration.

The dimensionless MIM wave component Eqs. (14) and (15) are

The waveguide axial variable relation in (45) is further simplified, as before, the waveguide axial component variable from (13)

The MIM dielectric ratio and scaled dielectric difference respectively are

The $L_+$ form of the MIM dispersion equations

The primary focus in the analysis of the MIM configuration will be at large $\lambda$ following [40]. At large $\lambda$, the dielectric ratio $\tilde \epsilon \ll 1$ and striking a balance between the two terms in (50) requires $Z=o(1)$. Inspection of (50) shows that $Z=O(\sqrt {\tilde \epsilon /\alpha })$ where $Z=o(1)$ is also consistent with $\tilde \epsilon \ll \alpha$ with $\alpha$ treated as an $O(1)$ constant. Setting $U=Z/(\sqrt {\tilde \epsilon /\alpha })$ and substituting to (50)

The numerical method follows as in Section 2.3 and referring to (51) if

The zeroth approximation to $U$ in (52) is $U_0=O(1)$ with $F=O(\sqrt {\alpha \tilde \epsilon })$.

The solution $U$ may be written as a regular perturbation expansion in $\tilde \epsilon$: $U = U_0 + \sqrt {\tilde \epsilon } U_1 + \tilde \epsilon U_2 + \cdots$, with $\tilde \epsilon \ll 1$ and $\tilde \epsilon \ll \alpha$. It follows by inspection from (51) that $U_0$ may be chosen in the first quadrant $U_0 = (1+\mathbf {i})/\sqrt {2}$ without loss of generality.

3.1 Closed-form approximation

The MIM propagation length and skin depth are examined at large $\lambda$ for $\tilde \epsilon \ll 1$ and $\alpha =O(1)$, as explained above: $U = U_0 + \sqrt {\tilde \epsilon } U_1 + \tilde \epsilon U_2 + \cdots$. Substituting this regular perturbation to (51)

Using $U$ from (54) with $k_x'$ from (43) and approximating for $\tilde \epsilon \ll 1$, the wave component $k_x'$ simplifies to

The dielectric constant $\tilde \epsilon _2$ is real and thus the propagation length $L=\lambda /|2\Im (k_x')|$, depends on the imaginary part of the second term in (55). After returning to dimensional variables, $L$ takes the form

The skin depth, or attenuation, $P=\lambda /|k_{z_2}'|$, depends on $k_{z_2}'$ (47) with (46) , which expanded for $\tilde \epsilon \ll 1$ is $|k_{z_2}'| = O(1/\tilde \epsilon ) \gg 1$ and to $O(1)$, and after returning to dimensional variables,

The comparison of the MIM propagation length $L$ (56) to the numerical solution is provided in Fig. 1 and the results approximate each other near peak propagation.

 

Fig. 1. The MIM structure is Ag/SiO$_2$/Ag with $\epsilon$$SiO2$ and $\epsilon$$Ag$ computed using formulas from the appendix of [39]. The solid and dashed lines are the respective numerical solution and closed-form approximation (56) for the noted propagation distance $L$. The thick and thin lines respectively are the propagation length for metal cladding thicknesses $d_1=100$nm and $d_1=50 nm$.

Download Full Size | PDF

3.2 Discussion

In the MIM configuration there is a complementary dependency on the dielectric constants compared to the IMI form: the IMI and MIM forms respectively depend upon (21) and (48). This reversal has a profound effect on the device properties since the metal dielectric constant has a real part with a magnitude that dominates the imaginary part at large wavelength. The end result is a MIM device with a simple form for the propagation length and skin depth: (i) wavelength dependence implicitly comes from the dielectric constants, (ii) there is a complementary dependence on dielectric constants allowing for simultaneous improvement of both propagation length and skin depth which contrasts with the IMI device, and (iii) a linear dependence on metal thickness $d$ and relatively small dependence on the wavelength at large wavelength due to (i) giving a desirable broadband response.

The applicability of the results to other MIM devices is considered for several cases taken from the literature and these appear in Table 1. The devices are constructed from a variety of materials and all cases have scaled dielectric differences $\alpha = O(1)$. The closed-form results from Section 3.1 closely approximate those from the literature.

Table 1. Three examples of the MIM configuration involving different MIM structures are presented. The closed-form results approximate those determined numerically.

View Table

4. Multilayer symmetric structures

Given the desirable features found in the symmetric MIM device, it is of interest to discover: (i) how the MIM device dispersion equations generalize in the asymptotic sense to the multi-layered case, and (ii) to what extent, and by what means might improved characteristics become available in a multi-layered device.

The methodologies developed for the three-layer configurations, as applied to the MIM device, are extended to repeating symmetric multilayer structures. This is illustrated here for a seven-layer MIMIMIM structure which represents a generalization of the MIM configuration.

4.1 Anti-symmetric $L_+$ mode dispersion equation

The dimensional form of the dispersion relation for the seven-layer MIMIMIM device is considered for the anti-symmetric $L_+$ mode [42].

The anti-symmetric dispersion relation is written in the form $H_1+H_2=0$ where $H_1$ and $H_2$ respectively refer to the ’uncoupled’ and ’coupled’ terms. For three-layer structures there are only the uncoupled terms and in this seven-layer device these correspond to

The coupled term $H_2$ may be written as

The coupling $H_2$ is not present in the three-layer structures where only two terms are must be balanced. The objective here is to compare the seven-layer waveguide to the previous analysis of the three-layer MIM structure and, therefore, the focus will at large wavelengths $\lambda$ where the magnitude of the metal complex dielectric real part dominates its imaginary part.

4.2 Closed-form approximation

The anti-symmetric dispersion equation $H_1+H_2=0$ is first rescaled following Section 2.1 and, fully restating to fix ideas, is

A closed-form approximation is found to the solution of $H_1+H_2=0$ as before. Applying (62) to the uncoupled and coupled terms, $H_1$ and $H_2$ respectively, and setting $\tilde \epsilon = \tilde \epsilon _2/\tilde \epsilon _1$ (48), with $|\tilde \epsilon | \ll 1$, shows that all three terms are small compared with unity.

The dominant balance indicates that the second term of the uncoupled component $H_1$ in (60) and coupled terms balance to give $z = O(\sqrt {\tilde \epsilon }) \ll 1$. The variable $z$ is rescaled to an $O(1)$ variable $W=z/\delta$ where $\delta = \sqrt {\tilde \epsilon }$ is for convenience.

Proceeding to $O(\delta ^3)$ in the coupled term $H_2$, treating $\alpha$ as $O(1)$, and working exactly with the uncoupled term $H_1$, yields

The first term in (63) is the $O(\delta ^3)$ approximation to the coupling $H_2$ while the last two terms are the exact restatement of the uncoupled form $H_1$ (60). The balance among terms in (63) is now clear: the first and third terms are $O(\delta )$ and form the dominant balance between the uncoupled and coupled terms.

The three terms in (63) respectively are $O(\delta )$, $O(\delta ^3)$, and $O(1)$ where $W$, $\alpha$, and $G(\beta )$ are assumed $O(1)$. The solution of (63) is expanded in a regular perturbation series as

The propagation length $L=\lambda /|2\Im (k_x')|$ where, from (67), ${k_x'}^2 = \tilde \epsilon _2 - {k_{z_2}'}^2$. Substituting $z$ from (67), expanding $k_x'$ for $\delta = \sqrt {\tilde \epsilon } \ll 1$

Since the insulator dielectric $\tilde \epsilon _2$ is real, the imaginary part of $k_x'$ resides in the second term of (68) and thus the propagation length $L=\lambda /|2 \Im (k_x')|$ written in dimensional form is

Proceeding similarly for the skin depth, or attenuation, $P = \lambda /|k_{z_2}'|$ and returning to dimensional variables leads to

4.3 Discussion

The primary difference between the three-layer MIM and seven-layer MIMIMIM devices is the existence of the coupling term $H_2$ in the latter that greatly enhances the propagation length for a judicious choice of cladding thickness. The coupling introduces a propagation length amplification factor $G(\beta )$ (64, 65) that enables increased propagation length over a broad band of higher frequencies.

The metal dielectric is considered here at large wavelength [42] and has $\Im (\tilde \epsilon _1) \ll |\Re (\tilde \epsilon _1)|$ in this region. Any reduction in $\Im (\tilde \epsilon _1)$ tends to increase the propagation length through its effect on $\Im (\sqrt {G})$ and therefore $\Im (\tilde \epsilon _1)$ is referred to here as a ’damping’. In the limit of zero damping, the parametric dependence of the propagation length $L$ (69) collapses to $L_0$ which depends upon the zero damping form of the scaled dielectric constant (71) and

The material considered here has two distinguishing features driving the large wavelength propagation length. Firstly, the thickness $d_{1n}$(74) has a weak dependence on $\lambda$ since $\Im (\tilde \epsilon _2) \ll \Re (\tilde \epsilon _2) \ll |\Re (\tilde \epsilon _1)|$ and, in particular, $|\Re (\tilde \epsilon _1)| \propto \lambda ^2$. Secondly, the damping is reduced as $\lambda$ is increased: in (65), the ratio $\Im (\tilde \epsilon _1)/(\tilde \epsilon _2 - \Re (\tilde \epsilon _1))$, $\Re (\tilde \epsilon _1)<0$, is monotonic decreasing for increasing wavelength. Taken together, these imply that increased wavelength is associated with increased proximity to the near constant $d_1=d_{1n}$ lines (74) where the propagation length is singular.

In Fig. 2, the propagation length $L$ (69) is shown as a function of cladding thickness $d_1$ and $\lambda$ near the line of singularity $d_{10}$ ((74), $n=0$, depicted as constant thickness, white line). The line of singularity (constant thickness white line), associated with zero damping, is approached as the wavelength is increased and this leads to the observed thickening of the saturated region ($L=60$nm) at increased wavelength. These singular choices of the metal thickness $d_1=d_{1n}$, $n=0,1,2,\ldots$ are approximate zones where a guaranteed minimum propagation level can be achieved for wavelengths that exceed a threshold. Propagation lengths at the scale of millimeters are available for relatively broad bands of wavelength near ($\lambda \approx 1000$nm, $d_1 \approx 33.2$nm) or ($\lambda \approx 1700$nm, $d_1 \approx 32.5$nm). Evaluation of the skin depth  (70) (not shown) in Fig. 2 leads to $P \ll L$ at large wavelength.

 

Fig. 2. The 7-layer device being examined is composed of alternating Ag/GaP layers, where $\epsilon$$Ag$ is computed using formulas from the appendix of [39] and $\epsilon$$GaP$ = 3.31 taken from [42]. The GaP dielectric layer thickness is 50nm. The constant width, solid white line is the undamped case(69) and the damped case (74) is the remaining greyscale (’shining’) white line. The propagation length $L$ (69) is shown in terms of cladding thickness $d_1$ and wavelength $\lambda$ near the line of singularity $d_{10}$ ((74), $n=0$, depicted as solid, constant thickness white line). Propagation lengths approaching millimeters are available for relatively broad bands of wavelength near ($\lambda \approx 1000$nm, $d_1 \approx 33.25$nm) or ($\lambda \approx 1700$nm, $d_1 \approx 32.75$nm). The nearly uniform difference between the two lines reflects the small damping, $\Im (\tilde \epsilon _1)$, over the range of frequencies in the figure.

Download Full Size | PDF

5. Conclusions

Closed-form approximations were developed via a novel asymptotic method constructed here for the symmetric multi-layer plasmonic dispersion equations. The original closed-form analysis recovers numerical results and illuminates fundamental differences in the physics of differences, found in propagation length and attenuation, between thin-film MIM and IMI devices. Tabulated results comparing the analytical results here to those found numerically in other work indicate close agreement to within 1%. Improved performance of the MIM versus IMI devices is discovered to primarily depend on a complementary dependence of the thin-film, MIM and IMI dispersion equations on the ratio of the material dielectric properties.

The analysis of multi-layer devices, demonstrated here for a seven-layer design, shows a remarkable departure from IMI and MIM devices to a propagation length and attenuation in the multi-layer design that is uniquely governed by a proximity to zero-damping lines of singularity. In these designs, the desired wavelength dependence may be tuned through judicious choices of the physical and material characteristics. The new results indicate that exploitation in multi-layer devices of the existence of zero-damping singularities remains and, as such, may be important for a range of applications.

The asymptotic method developed here has been demonstrated on thin-film, three- and multi-layer devices. The method provides insight into fundamental differences in the physics distinguishing three- and multi-layer thin-film designs. It has also been found that multi-layer hyperbolic metamaterial devices possess properties approximated by much simpler MIM devices [45]. Numerical studies such as these demonstrate that further analytical exploration is needed to explain the physics of their similar behaviour. Further application of asymptotic approaches, like those shown here, should allow for better understanding of the physics governing plasmonic devices and aid in development of improved designs with wider application.