Poloidal and toroidal plasmons and fields of multilayer nanorings

TL;DR

This paper analyzes the surface plasmon modes of multilayer toroidal nanorings, deriving dispersion relations, resonance frequencies, and Green's functions, with applications to controlling electromagnetic fields and quantum emitters.

Contribution

It introduces a comprehensive theoretical framework for poloidal and toroidal plasmons in multilayer nanorings, including new Janus nanoring configurations and generalized Green's functions.

Findings

Derived canonical dispersion relations for layered nanorings.

Identified distinct forward and backward mode coupling in tori.

Calculated resonance frequencies for gold, silver, aluminum, and silicon nanorings.

Abstract

Composite and janus type metallo-dielectric nanoparticles are increasingly considered as a means to control the spatial and temporal behavior of electromagnetic fields in diverse applications such as coupling to quantum emitters, achieve invisibility cloaks, and obtain quantum correlations between qubits. We investigate the surface modes of a toroidal nano-structure and obtain the canonical plasmon dispersion relations and resonance modes for arbitrarily layered nanorings. Unlike particle plasmon eigenmodes in other geometries, the amplitudes of the eigenmodes of tori exhibit a distinct forward and backward coupling. We present the plasmon dispersion relations for several relevant toroidal configurations in the quasistatic limit and obtain the dominant retarded dispersion relations of a single ring for comparison, discuss mode complementarity and hybridization, and introduce two new types of toroidal particles in the form of janus nanorings. The resonance frequencies for the first few dominant modes of a ring composed of plasmon supporting materials such as gold, silver, and aluminum are provided and compared to those for a silicon ring. A generalized Green's function is obtained for multilayer tori allowing for calculation of the scattering response to interacting fields. Employing the Green's function, the scalar electric potential distribution corresponding to individual poloidal and toroidal modes in response to an arbitrarily polarized external field and the field of electrons is obtained. The results are applied to obtain the local density of states and decay rate of a dipole near the center of the torus.

Figures

Equations

$$0={\int }_{\phi =0}^{2\pi }{\int }_{\eta =0}^{2\pi }{\int }_{\mu =0}^{\infty }{h}_{\mu }{h}_{\eta }{h}_{\phi }\rho (\mu ,\eta ,\phi )d\mu d\eta d\phi ,$$

$$\rho (\mu ,\eta ,\phi )=\frac{{ϵ}_0\delta (\mu -{\mu }_1)}{h_{\mu }^2}\sum _{m,n}{T}_n(\eta )V_m(\phi )W_{mn}({\mu }_1),$$

$$L=\frac{nm}{2}\int \eta \dot{\sigma }ds-\frac{1}{2}\int \frac{\sigma (\mathbf{r})\sigma ({\mathbf{r}}^{\prime })}{|\mathbf{r}-{\mathbf{r}}^{\prime }|}dsd{s}^{\prime },$$

$$R=r_1\cosh {\mu }_1=r_i\cosh {\mu }_i,$$

$${\mu }_i={\cosh }^{-1}\left(\frac{d_1}{d_1+\cdots +d_i}\cosh {\mu }_1\right),i=2,\dots ,k.$$

$$\Phi (\mu ,\eta ,\phi ,t)$$

$$+\sum _{j=1}^{k-1}\Theta ({\mu }_j-\mu )\Theta (\mu -{\mu }_{j+1}){\Phi }_{j+1}(\mu ,\eta ,\phi ,t)$$

$$+\Theta ({\mu }_k-\mu ){\Phi }_{k+1}(\mu ,\eta ,\phi ,t),$$

$${\Phi }_j$$

$$\displaystyle=f(\mu,\eta)\sum_{m=-\infty}^{\infty}\sum_{n=-\infty}^{\infty}\bigl{[}(1-\delta_{j,1})C_{mn}^{\text{j}}(t)P_{n-\frac{1}{2}}^{m}(\cosh\mu)$$

$$\displaystyle+(1-\delta_{j,k+1})D_{mn}^{\text{j}}(t)Q_{n-\frac{1}{2}}^{m}(\cosh\mu)\bigr{]}e^{in\eta}e^{im\varphi},$$

$$f(\mu ,\eta )=\sqrt{\cosh \mu }\left[1-\frac{1}{2}\frac{\cos \eta }{\cosh \mu }+\mathcal{O}\left({\left(\frac{\cos \eta }{\cosh \mu }\right)}^2\right)\right],$$

$${\Phi }_j\approx \left(\sqrt{\cosh \mu }-\frac{\cos \eta }{2\sqrt{\cosh \mu }}\right)\sum _{n,m}{L}_{mn}^j{e}^{in\eta }{e}^{im\phi },$$

$$\Phi_{j}\approx\sum_{n,m}\Bigl{(}I_{m,n-1}^{j}-4\cosh\mu\,I_{m,n}^{j}+I_{m,n+1}^{j}\Bigr{)}e^{in\eta}e^{im\varphi}.$$

$${{\Phi }_i|}_{\mu ={\mu }_i}$$

$${\epsilon }_i{\frac{\partial {\Phi }_i}{\partial \cosh \mu }|}_{\mu ={\mu }_i}$$

$${\mathcal{C}}_n=\left[\begin{array}{c}{C}_n^2\\ C_n^3\\ ⋮\\ C_n^{k+1}\end{array}\right]\text{ and }{\mathcal{D}}_n=\left[\begin{array}{c}{D}_n^1\\ D_n^2\\ ⋮\\ D_n^k\end{array}\right],$$

$${\mathbb{P}}_n{\mathcal{C}}_n={\mathbb{Q}}_n{\mathcal{D}}_n,$$

$${\mathcal{C}}_n={\mathbb{P}}_n^{-1}{\mathbb{Q}}_n{\mathcal{D}}_n.$$

$${\mathcal{W}}_{n+1}-{\text{R}}_n{\mathcal{W}}_n+{\mathcal{W}}_{n-1}=0,$$

$$\begin{array}{c}{\mathcal{W}}_n={\mathrm{J}}_n{\mathcal{D}}_n,\\ {\mathrm{J}}_n=({\mathbb{P}}_n^{\prime }{\mathrm{E}}_2{\mathbb{P}}_n^{-1}-{\mathbb{Q}}_n^{\prime }{\mathrm{E}}_1{\mathbb{Q}}_n^{-1}){\mathbb{Q}}_n,\end{array}$$

$${\mathrm{R}}_n={\mathrm{D}}_{\mu }+{\mathbb{A}}_n{\mathbb{B}}_n^{-1},$$

$$\begin{array}{c}{\mathbb{A}}_n={\mathbb{P}}_n{\mathrm{E}}_2{\mathbb{P}}_n^{-1}-{\mathbb{Q}}_n{\mathrm{E}}_1{\mathbb{Q}}_n^{-1},\\ {\mathbb{B}}_n={\mathbb{P}}_n^{\prime }{\mathrm{E}}_2{\mathbb{P}}_n^{-1}-{\mathbb{Q}}_n^{\prime }{\mathrm{E}}_1{\mathbb{Q}}_n^{-1}.\end{array}$$

$${\mathcal{X}}_n={\mathcal{W}}_n+{\mathcal{W}}_{-n}\text{ and }{\mathcal{Y}}_n={\mathcal{W}}_n-{\mathcal{W}}_{-n},$$

$${\mathcal{X}}_2-({\text{R}}_1-2{\text{R}}_0^{-1}){\mathcal{X}}_1$$

$${\mathcal{X}}_{n+1}-{\text{R}}_n{\mathcal{X}}_n+{\mathcal{X}}_{n-1}$$

$${\mathcal{Y}}_2-{\text{R}}_1{\mathcal{Y}}_1$$

$${\mathcal{Y}}_{n+1}-{\text{R}}_n{\mathcal{Y}}_n+{\mathcal{Y}}_{n-1}$$

$$X_{n+1}-{\text{R}}_n{X}_n+X_{n-1}=0,n=1,2,\cdots ,$$

$$\frac{X_2}{X_1}$$

$$V_{n-1}=\frac{I}{{\text{R}}_n-V_n},n=2,3,\cdots .$$

$$V_1=\frac{I}{{\text{R}}_2}-\frac{I}{{\text{R}}_3}-\cdots -\frac{I}{{\text{R}}_n-V_n}\cdot$$

$$\frac{I}{{\text{R}}_2}-\frac{I}{{\text{R}}_3}-\frac{I}{{\text{R}}_4}-\cdots .$$

Full text

I Introduction

Collective electronic effects in solid tori and toroidal shells by virtue of their topology exhibit unique electromagnetic response. Consequently, electronic excitation in nanorings has been of increasing importance in nano-optics aizp:rev , metamaterial mari , trapping of cold molecules salhi and are of great potential for achieving enhanced spectroscopies mich ; xu ; nikoob ; hayne ; grand , where giant enhancement in the near field of the surface modes of similar nano-structures associated with the collective electronic excitation is being utilized. Ordered or structured embedded nanoparticles in metamaterial as in photonic bandgap material, and material with negative index of refraction rely strongly on the shape of the structures vesel ; pendry1 . Increasingly, the local curvature and geometry of the nanoparticles are being capitalized upon in tip and surface-enhanced spectroscopies dier , in spectral tuning of optical response halas , in optical antenna bharad and chemical and biological sensing lee . Moreover, during the last decade, the implications of the emerging field of transformation optics leon ; pendry have driven the exploration of novel anisotropic materials and their realizations in spherical, cylindrical, and toroidal geometries within the optical regime. In these studies, to achieve the required anisotropy, mainly the homogenization principle of laminates milton has been resorted to, and the associated proof of principle of the modeled concepts has been attempted in the quasistatic regime kadic ; schit . In particular, in the case of toroidal structures, significant difficulties arise when employing such approaches due to the vectorial nature of the resulting three-dimensional electromagnetic problem. Specifically, the degree of the underlying challenges can be appreciated in the approximation employed by Kadic et al in their work on plasmonic analogs of electromagnetic wormholes kadica , where the availability of effective properties of metallo-dielectric coatings would have been beneficial.

Employing advanced fabrication, both metallic and nonmetallic nanorings are being considered. Metallic nanorings may be fabricated by various lithographic techniques such as vacuum evaporation and self-assembly. Aizpurua et al. fabricated gold nanorings on glass substrate using colloidal lithography and investigated their light scattering properties and found that the excitation wavelengths of metallic nanoparticles depend strongly on the particle geometry allowing one to choose and tailor the shape of a particle to achieve excitation spectra on demand aizp:quant . Similarly, while Kong et al. demonstrated the fabrication and application of zinc oxide nanorings in investigating fundamental physical phenomena, such as the Aharonov-Bohm oscillations in exciton luminescence besides possible use as nanoscale sensors, transducers, and resonators kong , Kuyucak et al. used the toroidal domain as a realistic model important for the conductance of ions in biological channels kuk . Many forms of toroidal structures are emerging, such as toroidal micelles of polystyrene-block-Poly(acrylic acid) c:liu , and toroidal carbon nanotubes l:liu . Similarly, various fabrication methods are being reported including a strategy for generating monodisperse toroidal particles through solidification of droplets of polymer solution in a microchannel b:wang , and very recently, vortex ring-derived particles were used toward mass production of larger rings an . Other more elaborate configurations continue to be envisioned and studied as reported by Chen et al. in their theoretical work on the mechanical properties of connected carbon nanorings n:chen . Also in recent development yi , Janus particles have been employed as a new tool for imaging and sensing in biological systems.

The expanding use of toroidal structures and fields is diverse. For example, virtually free from dissipation, all-dielectric metamaterials capable of supporting toroidal dipolar excitations has also been proposed basharin . Similarly, the excitation of toroidal dipolar moment in a novel metal-dielectric-metal nanostructure under the radially polarized light was investigated bao . A toroid-like nanostructure metamaterial comprising asymmetric double-bar magnetic resonators was designed to demonstrate the toroidal dipolar response in the optical regime dong1 . Also the toroidal dipole response was numerically investigated on a torus-like structure obtained by rotating the planar double-ring structure into a multifold double-ring torus-like metamaterial dong2 . Similar investigations have been performed in metamaterials with their unit cell including three magnetic resonators, metal-dielectric-metal, consisting of two Ag rods and a SiO2 spacer tang . Experiments and numerical simulations were conducted to investigate the feasibility of tuning the optical response of a dimer type nanoantenna using plasmonic nanorings panaretos . Treating the charged ionic channels as deformable torus with a circular or elliptical cross section, analytical solution of the equilibrium ion distribution for a toroidal model of a ionic channel was discussed enriquez . Other examples include obtaining an exact solution for the electrostatic potential of a family of conducting charged toroids with no hole in the toroid lekner ; obtaining the surface states of a toroidal topological insulator supporting both bound-states and charged zero-modes fonseca ; and analysis of a uniformly magnetized torus and its derivative self-intersecting shapes beleggia .

Interestingly, while a plethora of applications are emerging, early treatment of the toroidal problem has been reported in the context of plasma modes and stability. Specifically, Tokamak confinement of fusion plasma typically entails the generation and magnetic containment of a ring of plasma hinton . Toroidal plasma stability and control have been a major focus in the field of hot fusion research mukhov . These early works that appear to have been unnoticed by recent investigators, most likely due to a length-scale difference of at least $10^{6}$, provide much of the necessary mathematical developments. However, to date in the context of collective electronic effects at the toroidal surfaces, that is, surface plasmon excitation or particle plasmons, an accurate and comprehensive account on the toroidal response is not available. Furthermore, there are discrepancies in what has been reported, obscuring the application of the reported results to materials of current interest in nanoscience. The interior and exterior fields of an isolated dielectric ring subject to a uniform axial electrostatic field was first obtained by Love love:jmp . In a follow on article, Love also obtained the modes of a cold plasma love:plasma . Very recently, Di Biasio et al. calculated the polarizibility of various dielectric materials including a ring-shaped dielectric biasio . The quasi-static modes of a cold plasma obtained by Love were claimed by Avramov et al. to account for only half of the existing modes when taking into consideration the parity of the modes under coordinate inversion avramov . A more recent attempt to obtain the plasmon dispersion relations for a gold nanoring was made by Mary et al. However, the reported solutions were obtained by making an assumption that does not conform to a diagonalization or perturbation scheme to solve the resulting transcendental dispersion equation mary . Interestingly, employing a theorem established by Lord Kelvin, F. H. Safford, obtained the solution to the toroidal problem by using a geometric inversion and transformation from conical into toroidal surfaces saff . Similarly, an alternative separation scheme, where the coordinate dependence of the boundary conditions is altered, was discussed by M. Andrews andrew . While Love considered the Greens function to obtain the potential distribution of a torus in a static field, Scharstein and Wilson employed the static thin-wire kernel approximation in an integral equation to calculate the excitation of a conducting ring schars . Unlike the Laplace equation, the Helmholtz equation is not separable in the toroidal systems. Furthermore, the case of a multiply coated toroids have not previously been treated. While not straightforward, numerical techniques, in principle, can be employed to obtain the dispersion relations in the fully retarded case, both analytical and semi-analytical approaches in the quasi-electrostatic limit have proved to be useful. Tsuchimoto et al. employed the boundary element method to solve an eigenvalue problem derived from the Helmholtz equation to obtain the dispersion relations of a cold plasma for both circular and oval cross sections tsuch . Using the Hertz vector to obtain a scalar Helmholtz equation, Janaki and Dasgupta attempted a partial separation for the azimuthal solutions followed by expansion to obtain approximate eigenmodes of a toroidal cavity janaki .

Motivated by the above introduction, we proceed to present the organization of this article. In section II, we begin by formulating the general case of a $k$-layered toroidal structure ($k=1,2,\cdots$) and obtain the corresponding plasmon dispersion relations. These results, being reported for the first time, warrant an in-depth investigation pap1 due to a three-term difference equation that arises and its relation to continued fractions and infinite determinants. We will simply utilize the results of this investigation without giving a full account, which is out of the scope of this article. We will then specialize, in subsections II.A-II.B, the obtained results for a series of cases of current importance beginning with a $k$-layered toroidal structure dispersion relations, an isolated solid torus and its complementary toroidal void in II.B.1. Here, we will clarify the discrepancies noted above among the results reported by Love love:plasma , and Avaramov et al. avramov Significant simplification can be achieved by use of perturbation theory to derive the dispersion relations allowing to further calculate the effects of a planar substrate. Thus, using this approach, we investigate the more realistic case of a nanoring residing on a dielectric substrate, which is generally known to cause a red-shift in the plasmon energies. Furthermore, we will compute the dispersion relations for the specific materials gold, aluminum, silver, and silicon typically used in experiments. Proceeding with the treatment of the specialized cases, we will then obtain, in II.B.2, the plasmon dispersion relations for a bimetallic (i.e., singly coated) torus, an isolated solid shell, and a dielectric core-metal shell. Again, here as in the previous subsection, we will consider application of perturbation theory and use of specific materials. In subsection II.B.3, we will study two multilayer cases – a three layered ring in vacuum, and a four layered system. As a final example, we introduce two types of janus nanorings in II.B.4 and discuss how one may formulate the resulting problem. Section III is dedicated to the calculation and visualization of the various contributing modes to the response of the nanoring to external fields, where both uniform and nonuniform external fields are treated. Here, employing the Greens function approach, the three-term difference equation is solved for the two source terms corresponding to the external fields. We then extend the Greens function approach for the general multilayered tori. The results, which in principle can be used to construct the dyadic Green’s function and consequently also the local density of states (LDOS), are verified computationally using the finite element (FEM) and finite difference time domain (FDTD) methods. Here, to illuminate the obtained modes, we present, for the dominating modes, the computationally determined resonances of the single solid ring, and discuss two cases of janus nanorings. The scattering and absorption cross sections as well as the vibrational eigenmodes that may be excited following the photon scattering and plasmon decay for the same cases are also computed and discussed. Concluding remarks are provided in section IV.

II Plasmon Dispersion Relations

The multilayer ring depicted in Fig. 1 possess finite surfaces and therefore the spectrum of the normal modes representing the oscillations of any induced polarization charge will necessarily be discrete. A domain corresponding to a multilayer metal-dielectric nanoring may be modeled using the toroidal coordinate system $(\mu,\eta,\varphi)$. Here, $\mu\in[0,\infty)$, $\eta\in[0,2\pi]$, and $\varphi\in[0,2\pi]$, and the surface of a torus is defined by $\mu=\mu_{1}$, wherein the circle associated with the toroidal cross section grows from zero radius to infinite radius (with center at infinity) as $\mu$ decreases progressively from $\infty$ to [math]. Constant $\mu_{1}$ follows a circle of the minor radius $r=a/\sinh\mu_{1}$ centered at the major radius $R=a\coth\mu_{1}$ as the coordinate $\eta$ changes from [math] to $2\pi$ meaning that for a given $\mu=\mu_{1}$, $\eta$ traces a circle of fixed radius $r$ as it changes from [math] to $2\pi$. The coordinate $\varphi$ is the azimuthal angle about the $z$ symmetry axis. As $\varphi$ changes from [math] to $2\pi$, the circle is revolved around the $z$ axis by $360^{\circ}$ counterclockwise, thus generating a toroidal shape.

Prior to formulating the normal modes corresponding to charge density oscillations in the multi-domain depicted in Fig. 1, we note that the charge density $\rho=-\epsilon_{0}\Delta\Phi$ must satisfy

[TABLE]

where $h_{\mu}$, $h_{\eta}$, $h_{\varphi}$ are scale factors given by Eq. (21), and $\Phi$ denotes the electrostatic potential. For a single ring:

[TABLE]

where $T$ and $V$ are separated solutions and $W_{mn}(\mu_{1})$ is the Wronskian of the interior and exterior solutions, in conjunction with Eq. (1). In many cases the equation of motion for $\rho$, or rather the surface charge density $\sigma$, since here $\rho(\mu,\eta,\varphi)=\delta(\mu-\mu_{0})\sigma(\eta,\varphi)/h_{\mu_{0}}$, can be obtained and solved to calculate the eigenfrequencies board . While Eq. 1 is shown here to generate a condition for the surface modes (as described in SM VI) that is consistent with the structure of the dispersion relations, the eigenmode coupling in this case prevents an analytical expression for the eigenfrequencies. Similarly, viewing the plasmons as incompressible irrotational deformations of the conduction electron gas of density $n$, one may study the dynamics of small deformations of the electron fluid from the following Lagrangian  prodan ; board :

[TABLE]

where $m$ is the electron mass and the velocity potential satisfies $\Delta\eta=0$. However, as discussed in the next section, the amplitude coupling of the eigenmodes also prevents this formulation to generate analytical expressions for plasmon frequencies of the ring system.

II.1 General derivation for a $k$-layered ring

Nonretarded plasmon dispersion relations describing the resonant surface modes of a material domain may be obtained from a transcendental equation that is generated when imposing the quasi-static boundary conditions at the bounding surfaces of the domain within which the scalar electric field satisfies the Laplace equation board . As shown in Fig. 1 and Fig. 2, a $k$ layered torus can be described as a solid torus with toroidal surface $\mu=\mu_{1},$ frequency dependent dielectric function $\varepsilon_{1}$ and minor radius $r_{1}=d_{1},$ together with $k-1$ sublayers of toroidal shells (representing various material domains), each with corresponding thickness $d_{i}$ and dielectric function $\varepsilon_{i}$, where $i=2,\dots,k$. The first region corresponds to interior of the solid torus and is given by $\mu\geq\mu_{1}.$ The remaining $k-1$ concentric toroidal shells share the same major-axis $R$ with the solid torus. Therefore, one can compute the boundary surfaces between the layers via the relation

[TABLE]

where $r_{1}=d_{1}=R/\cosh\mu_{1},$ $r_{i}=d_{1}+\dots+d_{i},$ and $i=2,\dots,k$ (see Fig. 2). Since $r_{1}<r_{2}<\dots<r_{k},$ the relation (3) gives $\mu_{1}>\mu_{2}>\dots>\mu_{k},$ where

[TABLE]

We can thus describe the regions enumerated from 2 to $k$ by $r_{i-1}\leq r\leq r_{i}$ ($i=2,\dots,k$), respectively. Relation (3) transforms the inequalities given in $r$ into the corresponding reversed inequalities in $\mu$ given by $\mu_{i}\leq\mu\leq\mu_{i-1},$ where each $\mu_{i}$ satisfies (4). Finally, the $(k+1)^{th}$ region, which lies outside the $k$-layered torus, is described by $0\leq\mu\leq\mu_{k}$ (see Fig. 2 for the case $k=3$). The general form of the quasi electrostatic potential satisfying $\Delta\Phi=0$ is given by

[TABLE]

and

[TABLE]

where $f(\mu,\eta)=\sqrt{\cosh\mu-\cos\eta},$ $j=1,\dots,k+1$, and $\Theta$ is the Heaviside function. Utilization of the delta functions assures consistency with the toroidal harmonics $P_{n-\frac{1}{2}}^{m}(\cosh\mu)$ and $Q_{n-\frac{1}{2}}^{m}(\cosh\mu)$ being unbounded for $\mu\to\infty$ and $\mu\to 0$, respectively, and $C_{mn}^{\text{j}}(t)$ and $D_{mn}^{\text{j}}(t)$ are the non-retarded mode amplitudes at time $t$.

Prior to invoking the boundary conditions, it is noteworthy to observe that the consequence of the nonlinearity in the weight function $f$ for the interrelation of the surface mode amplitudes may be regarded as a manifestation of the surface charges affecting each other across the toroidal opening. This conjecture may be supported by noting that the dominant modes of the $mn^{\text{th}}$ harmonics may be studied within the first-order approximation:

[TABLE]

in (II.1) to obtain

[TABLE]

where $L_{mn}^{j}$ stands for the expression given in the bracket of (II.1). Setting $L_{mn}^{j}=-4\sqrt{\cosh\mu}\,I_{mn}^{j}$, the above may be rewritten as:

[TABLE]

This result is also consistent with the lack of distinct $mn$ modes (as described in SM II), where complete modes are rather seen to be of the form $m\sum_{n}$ type. In fact, the electron gas occupies a non simply-connected region of space with genus one, where the potential varies across the central hole. The electron motion in a region of the torus may therefore couple to the electron motions at a different location, leading to the described mode coupling. Such mode coupling is absent for electron motion in simply-connected regions with genus zero such as spheres, ellipsoids, cylinders, etc.

We now proceed and impose the ($2k$) boundary conditions:

[TABLE]

where $i=1,\dots,k.$ Since $\{e^{im\varphi}\}$ is a complete orthogonal system, we can treat Eqs. (8) and (9) separately for each integer $m$. With $m$ maintained fixed, one may therefore suppress its notation in upcoming equations. Moreover, we adopt the simplified notations $P_{n}^{i},$ $Q_{n}^{i},$ $f_{i},$ $C^{i}_{n}$, $D^{i}_{n}$ for $P_{n-\frac{1}{2}}^{m}(\cosh\mu_{i})$, $Q_{n-\frac{1}{2}}^{m}(\cosh\mu_{i})$, $f(\mu_{i},\eta)$, $C^{i}_{mn}(t)$, $D^{i}_{mn}(t)$ respectively, and let $P_{n}^{{}^{\prime}i}$, $Q_{n}^{{}^{\prime}i}$, $f^{\prime}_{i}$ denote the derivatives of $P_{n}^{i}$, $Q_{n}^{i}$, $f_{i}$ with respect to $\cosh\mu$ evaluated at $\mu=\mu_{i}$. Introducing the $k\times 1$ vectors

[TABLE]

and applying the boundary condition (8), one obtains the relation

[TABLE]

where the $k\times k$ matrices $\mathbb{P}_{n}$ and $\mathbb{Q}_{n}$ are defined by Eqs. (2) and (3). Since both bidiagonal matrices $\mathbb{P}_{n}$ and $\mathbb{Q}_{n}$ have non-zero diagonal entries, they are invertible. Thus, one can solve, say, $\mathcal{C}_{n}$ in terms of $\mathcal{D}_{n}$ to obtain

[TABLE]

Using the substitution (12) into the second boundary condition (9) and after some involved algebra pap1 , one can show that $\mathcal{D}_{n}$ satisfies the vector three-term recurrence relation:

[TABLE]

for $n=0,\pm 1,\pm 2,\cdots$, where

[TABLE]

and

[TABLE]

where

[TABLE]

For the definition of the $k\times k$ matrices $\mathbb{P}^{\prime}_{n},$ $\mathbb{Q}^{\prime}_{n},$ $\mathrm{D}_{\mu},$ $\mathrm{E}_{1},$ and $\mathrm{E}_{2},$ see Eqs. (4)–(7). Note that Eqs. (12)–(16) contain all that is needed to calculate the quasi electrostatic potential $\Phi$ (Eq. (II.1)) and normal modes for a $k$-layered ring in the absence of any applied field.

II.2 Surface modes of $k$-layered ring configurations

II.2.1 Multilayered dispersion relations, an isolated solid torus, and a toroidal void

Plasmon dispersion relations for a $k$-layered ring

Noting love:plasma that the toroidal harmonics $P_{n-\frac{1}{2}}^{m}(z)$ and $Q_{n-\frac{1}{2}}^{m}(z),$ and their derivatives with respect to $z$ are symmetric in $n$, it follows from Eqs. (15) and (16) that $\text{R}_{n}=\text{R}_{-n}$ for all $n=\pm 1,\pm 2,\cdots$. Using this together with substitutions

[TABLE]

one can transform the bi-infinite vector three-term recurrence (13) into two semi-infinite vector three-term recurrences

[TABLE]

and

[TABLE]

It is easily seen that the above system displays a single vector three-term recurrence relation together with two initial conditions given by Eqs. (17) and (19). This observation leads us to consider the more general format of the matrix three-term recurrence:

[TABLE]

with $\mathcal{X}_{n}$ and $\mathcal{Y}_{n}$ being treated as a single column of the $k\times k$ matrix $X_{n}.$ It follows that (17) and (19) give the two initial conditions

[TABLE]

where $I$ is the $k\times k$ identity matrix and the quotient of two $k\times k$ matrices $A$ and $B$ with $B$ non-singular is defined by $A/B:=AB^{-1},$ a convention that will be used throughout the the paper.

Let us assume that the sequence of $k\times k$ nonsingular matrices $\{X_{n}\}$ is a solution to the matrix three-term recurrence(21). Multiplying Eq. (21) from the right by $X_{n}^{-1}$ and defining $V_{n}$ by $V_{n}=\dfrac{X_{n+1}}{X_{n}},$ we arrive at the first order non-linear matrix recurrence:

[TABLE]

Successive iteration of Eq. (24) yields the finite matrix continued fraction:

[TABLE]

It turns out that the limit of the left-hand side of Eq. (25) exists as $n\to\infty;$ i.e., one can show Perron_Mat ; lev the convergence of infinite matrix continued fraction (MCF):

[TABLE]

The answer as to which unique value the MCF (26) converges is addressed in the following result:

It follows that the MCF (26) converges if and only if the matrix three-term recurrence (21) has a minimal solution. Moreover, the minimal solution $X_{n}$ satisfies

[TABLE]

The proof of the above fact is quite involved and is beyond the scope of this paper but is reported elsewhere pap1 .

In the scalar case of (21), a solution $X_{n}$ is called minimal Gau if $\lim_{n\to\infty}\frac{X_{n}}{Y_{n}}=0$ holds for any other linearly linearly independent solution $Y_{n}.$ This notion can be generalized Ah1 ; Ah2 to the matrix case of Eq. (21). As the above result indicates, a matrix three-term recurrence may very well have no minimal solution. These type of non-minimal solutions are called the dominant solutions in accordance with the terminology for the scalar case.

We turn our attention back to the Eq. (25). Since the MCF (26) converges, the matrix three-term recurrence (21) has the minimal solution, say, $X_{n}$ and one can let $n\to\infty$ in (25). Assuming that $\mathcal{X}_{n},$ defined by Eqns. (17) and (18), constitute a single column of $X_{n},$ it follows from (27) and (22) that

[TABLE]

A similar argument applied to $\mathcal{Y}_{n},$ defined by Eqns. (19) and (20), together with (23) implies

[TABLE]

Equations (28) and (29) are the sought dispersion relations for a $k$-layered ring configurations.

We end this section by making some important remarks regarding the obtained dispersion relations. First of all, the existence of the minimal solution does not guarantee its derivation. In fact, almost all initial conditions generate a dominant solution to the matrix three-term recurrence (21). Simply put, to find the minimal solution of (21) is equivalent to finding the exact analytic expression for the limit of MCF (26), which is not possible. So the follow up question is how we could assume the minimality of either $\mathcal{X}_{n}$ or $\mathcal{Y}_{n}$ in obtaining the dispersion relations, where the minimality of $\mathcal{X}_{n}$ or $\mathcal{Y}_{n}$ is referred to the case in which one of these two vectors constitutes a single column of the minimal solution, say, $X_{n}$ of (21). The answer to this question lies in the fact that the convergence of the MCF (26) is independent of the particular choice of the dielectric values $\varepsilon_{1},\varepsilon_{2},\dots,\varepsilon_{k+1}.$ These dielectric values are implicitly embedded as diagonal matrices $\text{E}_{1}$ and $\text{E}_{2}$ in the definition of $\text{R}_{n}$ defined by the Eq. (15). Thus, by assuming the minimality of either $\mathcal{X}_{n}$ or $\mathcal{Y}_{n},$ one can seek those sets of dielectric values $\varepsilon_{1},\varepsilon_{2},\dots,\varepsilon_{k+1}$ which satisfy (28) and (29), respectively. In doing so, one finds the surface normal modes of a $k$-layered ring. Finally, recall that for two minimal solutions, say, $Y_{n}$ and $Z_{n}$ of (21), we have $\frac{Y_{2}}{Y_{1}}=\frac{Z_{2}}{Z_{1}}.$ By considering the initial conditions (22) and (23), it follows that both $\mathcal{X}_{n},$ and $\mathcal{Y}_{n}$ can not be minimal at the same time; that is, no unique set of dielectric values $\varepsilon_{1},\varepsilon_{2},\dots,\varepsilon_{k+1}$ may satisfy both dispersion relations. In conclusion, we have two independent sets of dispersion relations (28) and (29) corresponding to the assumptions of $\mathcal{X}_{n},$ and $\mathcal{Y}_{n}$ being the minimal solution, respectively.

Clarification of discrepancies in existing models

As noted in the introduction, within the context of fusion and Tokamak physics, Love love:plasma and Avramov et al. avramov explored the surface modes of a toroidal plasma. More recently, within the context of condensed matter and nearfield optics, Mary et al. mary attempted to obtain the surface plasmon dispersion relations for an isolated solid nanoring. In this section, by considering the $k=1$ case, that is, a solid torus without layers, we will delineate the derivation of the associated modes and address the differences among the reported results.

Thus, let us consider a solid ring with the dielectric function $\varepsilon_{1}$ immersed in a medium with a dielectric function $\varepsilon_{2}.$ Before we proceed, it should be emphasized that all the involved vectors and matrices introduced in the derivation of the dispersion relations are now simply scalars. Invoking Eqs. (12)–(14), it follows from (8) that the inside and outside potentials $\Phi_{1}$ and $\Phi_{2}$ are given by

[TABLE]

and

[TABLE]

respectively, where $u=\cosh\mu$ and the summation index $m$ has been suppressed in $\mathbb{P}_{n},$ $\mathbb{Q}_{n},$ $\text{J}_{n},$ and $\mathcal{W}_{n}.$

Love love:plasma derives the three-term recurrence (13) for $\mathcal{W}_{n},$ but only considers the transformation $\mathcal{X}_{n}=\mathcal{W}_{n}+\mathcal{W}_{-n}.$ As a result, only the dispersion relation (28) is obtained. We now know that $\mathcal{W}_{n}=\frac{1}{2}\mathcal{X}_{n}+\frac{1}{2}\mathcal{Y}_{n},$ which clearly indicates that one should also consider $\mathcal{Y}_{n}$ as a minimal solution. This has led us to the second dispersion relation (29) discussed earlier. Love, considering the more involved transformation $V_{n}=\dfrac{\mathcal{X}_{n+1}}{\mathcal{X}_{n}}$ at once, might have missed the second dispersion relation.

On the other hand, Avramov et al. avramov , aiming to study the high-frequency surface waves on a toroidal isotropic plasma, correctly obtains both dispersion relations (28) and (29) for a single toroidal plasma.

It is easily seen that Eqs. (17)–(18) and (19)–(20) can be equivalently expressed as two infinite homogeneous systems of linear equations:

[TABLE]

[TABLE]

Avramov et al. avramov describes that the dispersion relations (28) and (29) are obtained using the statement that the infinite systems (32) and (33) possess nontrivial solutions if the determinants of the corresponding infinite matrices vanish. Contrary to the finite case, the mentioned statement is incorrect in the infinite case. In general, the determinant of an infinite matrix may not exist. In fact, there are infinite matrices with non-zero or non-existing determinant for which the corresponding homogenous systems still possess non-trivial solutions. A simple example of such a case is provided in SM V. Moreover, any formal calculation of the determinants of (32) and (33) does not result in a continued fraction resembling either (28) or (29). A similar incorrect interpretation of the determinant of an infinite matrix and non-triviality of the solutions also occurs in Love love:plasma .

Our final remark concerns with the parity of the dispersion modes avramov . Using the standard decomposition of $\mathcal{W}_{n}=\frac{1}{2}\mathcal{X}_{n}+\frac{1}{2}\mathcal{Y}_{n}$ into symmetric ($\mathcal{X}_{n}$) and antisymmetric ($\mathcal{Y}_{n}$) terms, together with the symmetry of $\mathbb{P}_{n},$ $\mathbb{Q}_{n},$ and $\text{J}_{n},$ one can split each of the sums indexed by $n$ in (30) and (31) into two separate terms with respect to $\sin n\eta$ and $\cos n\eta$ as:

[TABLE]

where

[TABLE]

and $f=f(\mu,\eta).$ The expressions for $\Phi_{2,m}^{c}$ and $\Phi_{2,m}^{s}$ can be found in a similar fashion with $Q_{n-\frac{1}{2}}^{m}(u)$ and $\frac{1}{\text{J}_{n}}$ replaced by $P_{n-\frac{1}{2}}^{m}(u)$ and $\frac{\mathbb{Q}_{n}}{\mathbb{P}_{n}\text{J}_{n}},$ respectively. The reason behind this decomposition is to identify the even and odd parities of $\Phi_{i,m}^{c}$ and $\Phi_{i,m}^{s}$ with respect to the toroidal inversion or point reflection:

[TABLE]

Clearly, for $i=1,2$:

[TABLE]

Under this transformation, $\Phi_{i,m}^{c}$ ($i=1,2$) is an even function for $m=\pm 2,\pm 4,\cdots$ and an odd function for $m=\pm 1,\pm 3,\cdots.$ A similar statement holds for $\Phi_{i,m}^{s}$ ($i=1,2$) with reversed parity for being even or odd.

In conclusion, one can use the parity argument in the following sense. By assuming the minimality of $\mathcal{X}_{n},$ the normal modes are calculated from the dispersion relation  (28) while $\mathcal{Y}_{n}$ becomes a dominant solution. Now for a fixed $m,$ if $m$ is even (odd), then the even parity corresponds to the solution $\mathcal{X}_{n}$ ($\mathcal{Y}_{n}$) and the odd parity corresponds to the solution $\mathcal{Y}_{n}$ ($\mathcal{X}_{n}$). An exact similar argument applied to the case when $\mathcal{Y}_{n}$ is assumed as the minimal solution yields the dispersion relation  (29) with $\mathcal{X}_{n}$ as a dominant solution, where the odd and even parity can be recognized. As shown in the next subsection, however, this method works only for small aspect ratio where there is a significant difference between the two separate sets of the normal modes. For larger aspect ratios, this difference is almost unrecognizable.

Exact numerical solutions

The results of the numerical solution of the dispersion relations (28) and (29) for a single solid torus are shown in Fig. 3.

To facilitate comparison love:plasma , the plots are obtained for the mode values of $m=1,2,3,4$ with $n=1,2,3$. The values of the associated Legendre functions $P_{n}^{i},\hskip 2.84526ptQ_{n}^{i}$ in Eqs. (28) and (29) are computed using the algorithm described in Gil and Segura segura:cpc ; gil:jcp ; gil:cpc . Numerical analysis may be performed on Eq. (13) with $k=1$ in the spherical and cylindrical limits wherein the aspect ratio tends to 1 and $\infty$, respectively. As $Q_{n}^{1}$ diverges as $\mu\rightarrow 0$ and $P_{n}^{1}$ diverges as $\mu\rightarrow\infty$, one can numerically observe that the coefficient terms associated with $\mathcal{D}_{k}$ for $k=n-1$, $n$, $n+1$ in Eq. (13) diverge in the spherical limit and converge to zero in the cylindrical limit. This observation is consistent with the earlier note on the surface charge influence across the toroidal opening ($n$, $n-1$, and $n+1$ decouple). Also in Fig. 1 in SM II we present various limiting cases of a ring.

The presented numerical scheme may be employed directly to also obtain the eigenfrequencies of a single toroidal cavity or alternatively we may consider using the geometric complementarity and the sum rule introduced by S. Apell et al. sumrule . The surface plasmon sum rule may therefore be used to check the validity of numerical calculations performed on the ring configurations studied. Here, to describe the metallic media we have employed the frequency $\omega$ dependent dielectric function:

[TABLE]

where $\omega_{p}$ is the plasma frequency, $\Gamma$ is the electron collision rate (damping), and $\varepsilon_{\infty}$ corrects for the high frequency response due to interband transitions of bound electrons (numerically, e.g. for gold, $\hbar\Gamma=0.069$ eV, $\omega_{p}=8.95$ eV, and $\varepsilon_{\infty}=9.5$). For the sake of discussions, we occasionally use: $\varepsilon_{\infty}=1$ and $\Gamma=0$, that is, free electron gas. Denoting the frequencies of the surface modes of the toroidal void with $\omega_{v}$ and those of the solid ring with $\omega_{r}$, the inset in Fig. 4 depicts the complementarity of the modes, which when using Eq. (38) without damping reads:

[TABLE]

The satisfaction of the plasmon sum rule for the ring-cavity pair is shown in Fig. 4, where Eq. (28) was solved for $\varepsilon_{1}=\varepsilon_{r},\varepsilon_{2}=1$, that is, a vacuum bounded solid ring, and subsequently for $\varepsilon_{1}=1,\varepsilon_{2}=\varepsilon_{v}$, that is, a toroidal vacuum cavity. In Fig. 4, the dispersion relations Eq. (28) were obtained for $m=4$ and $n=0$.

Given that the wave equation is not separable in the toroidal system, the results in Fig. 3 are particularly useful. Further information on the eigenmodes and photonic and plasmonic response may be obtained computationally. In addition to computing the eigenfrequencies by solving an eigenvalue equation over an appropriately defined computational domain, we may compute the dominant branches of the retarded dispersion relations from the scattering properties of the nanorings. To begin with, employing the FEM method, in Fig. 2 in SM VII we present the computationally determined transient response of a gold nanoring to linearly polarized incident photons, and for comparison in Fig. 5, the quasi-static response to different polarizations.

Parametric studies based on such results allow the determination of the retarded dispersion relations. In Fig. 6 we present a comparison of the nonretarded dispersion relations to some of the dominant retarded dispersion relations of a gold nanoring obtained computationally using FDTD. The FDTD converts the computational domain occupied by the ring and the surrounding medium into a lattice with a basic cell (the Yee cell) over which the field components are staggered in space and time. The fields are then evolved (time marching via leapfrog scheme) while satisfying Maxwell’s equations. This procedure obtains retarded fields and quantities. As can be seen, the agreement is quite good for larger aspect ratios, which is reasonable in light of Fig. 3. We also observe that, for the same toroidal mode number, the agreement improves with higher poloidal mode numbers.

Perturbation technique

Since the dispersion relations cannot be obtained analytically from the continued fractions, approximate analytical solutions may be sought employing perturbation theory. Here we show that a decoupling of the poloidal surface modes as presented by the three-term recursion can be achieved using the perturbation approach. To demonstrate the mode decoupling, after specializing Eq. (II.1) for $k=1$, and decomposing the parity of the poloidal solutions, we consider the even surface modes i.e., $\cos n\eta$ modes with fixed $m$ (the odd modes $\sin n\eta$ may be treated similarly). Setting $\varepsilon=\varepsilon_{1}/\varepsilon_{2}$, $u=\cosh\mu$, $u_{1}=\cosh\mu_{1}$, for convenience, the conditions Eq. (8) and Eq. (9), yield the Fourier modes:

[TABLE]

[TABLE]

where we defined $b_{-1}=0$. This equation can be solved using perturbation theory. We shall work in the scaling limit $u_{1},m\to\infty$, with $m/u_{1}$ held fixed. In this limit,

[TABLE]

[TABLE]

and the right-hand side of Eq. (II.2.1) vanishes. Therefore, Eq. (II.2.1) reduces to

[TABLE]

at $u=u_{1}$, and the modes decouple. For the $n=n_{0}$ mode, we obtain

[TABLE]

where $I$ and $K$ are the modified Bessel functions of the first and second kind, respectively. For first order in perturbation theory, set $b_{n_{0}}=1$, and expand

[TABLE]

where

[TABLE]

Notice that $C_{n_{0}}^{0}=0$. At $\mathcal{O}(1/u_{1})$, Eq. (II.2.1) becomes

[TABLE]

therefore

[TABLE]

where $(n\neq n_{0}\pm 1)$. At $\mathcal{O}(1/u_{1}^{2})$, Eq. (II.2.1) becomes

[TABLE]

For $n=n_{0}$, this gives

[TABLE]

In the special case $n_{0}=0$, the above reduces to

[TABLE]

Having obtained $C_{n_{0}}^{2}$, we then extract $\varepsilon_{mn}^{2}$ from Eq. (II.2.1).

A comparison between zeroth order and next-to-leading order is shown in Fig. 7. When corrections are included, the analytic formulas agree well with numerical results down to small aspect ratios. They disagree near the minimum aspect ratio $u_{1}=1$, indicating the effect of mode amplitude coupling in this regime, and thus need for inclusion of higher orders.

Plasmon redshift due to substrate effect

From an experimental point of view, most measurement efforts will consider use of a supporting substrate. It is known that, as in the case of thin solid films, generally surface plasmon energy redshifts when a dielectric substrate is considered. Similar modifications have been reported for a hyperboloidal nanostructures near a dielectric surface pass . Employing the perturbation approach, here we calculate the redshift of the plasmon energies by obtaining the first order correction of the vacuum bounded plasmon dispersion due to the substrate.

Suppose that the ring is placed at $z=0$ plane on a dielectric substrate (occupying the $z<0$ region) characterized by $\varepsilon^{\prime}$. Let $\Phi_{0}(\vec{r})$, $\vec{r}=(x,y,z)$ in Cartesian coordinates be the potential without the substrate. Then the potential with the substrate is found using the method of images to be

[TABLE]

We define the toroidal coordinates $u_{+}=\cosh\mu_{+}$, $\eta_{+}$ related to the Cartesian coordinates $\rho=\sqrt{x^{2}+y^{2}}$, $z$ by

[TABLE]

For the image ring, we similarly define toroidal coordinates $u_{-}=\cosh\mu_{-}$, $\eta_{-}$ by

[TABLE]

Concentrating on the $\cos n\eta$ modes and fixed $m$, we have for $z>0$,

[TABLE]

respectively for inside and outside the ring. Now for the image we have:

[TABLE]

Continuity at $u_{+}=u_{1}$ implies

[TABLE]

as before, because the potential of the image is already continuous. Continuity of the $u_{+}$-component implies

[TABLE]

at $u_{+}=u_{1}$. Notice that both $u_{-}$ and $\eta_{-}$ depend on $u_{+}$. To solve this, we shall treat the right-hand side as a perturbation. At zeroth order, the modes decouple. After fixing $n=n_{0}$, we obtain the zeroth order estimate

[TABLE]

Thus, in the zeroth order, the substrate effects on plasmon energies are absent. To calculate the redshift, we proceed to the first-order calculation. Introducing the expansion (bookkeeping) parameter $\xi$, we re-write Eq. (II.2.1) as

[TABLE]

We fix $b_{n_{0}}=1$, and expand

[TABLE]

where $(n\neq n_{0})$ and $\overline{\varepsilon}_{mn_{0}}^{1}$ is the first-order correction due to the substrate.

At $\mathcal{O}(\xi)$, Eq. (II.2.1) becomes

[TABLE]

Multiplying by $\cos n_{0}\eta_{+}$ and integrating over $[0,\pi]$, we deduce the sought correction to the $n=n_{0}$ mode due to the substrate,

[TABLE]

which may be simulated to obtain the amount of the shift induced in the plasmon energies. While the correction vanishes for $\varepsilon^{\prime}\to 1$, that is, for a vacuum bounded isolated ring, it diverges for $\varepsilon^{\prime}\to-1$. Interestingly, for the surface plasmons on a plane-bounded semi-infinite metal, the quasi-static resonance occurs when the real part of the dielectric function is -1.

Plasmon redshift for a ring on a substrate with finite thickness

Thin solid films are ubiquitous within plasmonics and nano-optics. We therefore suppose that the substrate has finite thickness $d$, that is, it possesses two planar interfaces and therefore two dispersion branches, which in the metallic case correspond to symmetric and anti-symmetric electron density oscillations pass . For thick films (large $d$), the two branches degenerate into a single branch, which in the nonretarded limit approaches -1. It is perhaps best to solve this problem iteratively. We start with the solution for the semi-infinite slab, i.e., for $z>0$

[TABLE]

for $-d<z<0$

[TABLE]

By design, this potential satisfies the correct boundary conditions at the interface $z=0$. For correct boundary conditions at $z=-d$, we take the image of $\Phi_{II}^{(0)}$ with respect to $z=-d$ and obtain the modified potential

[TABLE]

where $\Phi_{III}$ is defined in the region $z<-d$, satisfying $\Phi_{II}^{(1)}=\Phi_{III}^{(1)}$, and $\varepsilon^{\prime}\frac{\partial\Phi_{II}^{(1)}}{\partial z}=\frac{\partial\Phi_{III}^{(1)}}{\partial z}$ at $z=-d$. Next, we correct $\Phi_{II}$ and $\Phi_{I}$ so that they satisfy the correct boundary conditions at $z=0$. $\Phi_{I}^{(0)}$ and $\Phi_{II}^{(0)}$ already do, so we only need to consider the second image and its image with respect to $z=0$. We obtain

[TABLE]

satisfying $\Phi_{I}^{(2)}=\Phi_{II}^{(2)}$, and $\frac{\partial\Phi_{I}^{(2)}}{\partial z}=\varepsilon^{\prime}\frac{\partial\Phi_{II}^{(2)}}{\partial z}$ at $z=0$. Alternating between $z=0$ and $z=-d$, we thus build the potential as a series of images, which can be subsequently used to obtain a similar correction as in Eq. (66).

Plasmon dispersion in metallic and semiconductor nanorings

Experiments involving surface enhanced spectroscopies, nanophotonics and plasmonics typically employ metallic nanorings made of aluminum, gold, and silver as well as semiconductors such as silicon and germanium. For the case of vacuum bounded modes (i.e., $\varepsilon_{2}=1$), the plasmon resonance frequencies as a function of aspect ratio can provide useful information regarding the characteristics of the excitation photons impinging on the nanorings. To obtain the dispersion relations for specific materials, we employ an interpolations scheme between the solutions of Eq. (28) with $k=1$, and the experimental data composed by Palik palik .

This results in the dispersion relations shown in Fig. 8 for the modes $m=1$, $n=0$ and $m=4$, $n=1$. The usefulness of these results can be appreciated when seeking the response of the system. For example, for an aspect ratio $\cosh\mu=2$, following the nonretarded resonance values for the dominant mode $m=1$, $n=0$, presented in Fig. 8, the corresponding absorption cross sections are displayed in the insets of Fig. 8. The fully retarded cross sections in the insets of Fig. 8 were here computed by numerically solving Maxwell’s equations over a spatial domain containing the solid gold ring in vacuum using the FDTD method fdtd . As can be seen the agreement is quite good between the two cases. Also the FEG dispersion curves presented in Fig. 8 exhibit qualitatively similar patterns but quantitatively are different due to the effect of $\omega_{p}$ values.

It is further useful to study the behavior of the modes with respect to poloidal variations. Therefore, we fix $\cosh\mu=2$ and consider a solid ring in vacuum ($\varepsilon_{2}=1$) to study the corresponding resonance values of $\varepsilon_{1}$ satisfying the dispersion relation Eqs. (28), (29) with $k=1$ or, for comparison, the cylindrical limit solution love:plasma . The results are shown as $1-\varepsilon$ in Fig. 9(a) for various $n$ at the toroidal mode $m=4$. The obtained solutions are then matched with the experimental data palik to determine the resonance frequency spectrum of specific materials. This is as shown in Figs. 9(b), 9(c) for aluminum and silver. Working in the cylindrical limit, it is straightforward to obtain solutions for comparison with the exact solutions using the continued fraction Eq. (28). Thus, for the sake of sensing the asymptotic behavior $n=0,1,\cdots,8$ in the cylindrical limit suffice. Comparison with the exact solutions for $n=0,1,\cdots,6$ clearly reveals the convergence to the asymptotic limit. Clearly, while $1-\varepsilon$, for higher values of $n$, decreases asymptotically to 2 (and the decrease is much more rapid as $\cosh\mu\to\infty$), the metallic resonance frequencies grow larger for both aluminum and silver nanorings, as one may expect (implying higher energy photons for excitation of higher $n$ modes).

II.2.2 A singly coated toroid and its complementary void

Surface modes of coated rings in vacuum

The specific case of a coated ring is of particular interest due to its potential use in modeling an oxide layer or a thin film or a monolayer of a nonlinear material in strong exciton-plasmon coupling studies, or for generating thin walled (shells) representing carbon based molecules or a hollow cage of atoms such as fullerenes, etc.

Plasmon excitation in a singly coated torus embedded in an external medium or vacuum may be studied through analysis of the dispersion relations for the surface modes in a layered torus with fixed $\mu_{1}$, $\mu_{2}$. We begin by considering variations in $\varepsilon_{1}$, $\varepsilon_{2}$ for the system in vaccum, i.e., $\varepsilon_{3}=1$. The dispersion relations associated with this model are given by setting $k=2$ in Eqs. (28) and (29). To check the correctness of the numerical simulation in case of a multilayered torus, we first perform a limiting study of dielectric function to compare the dispersion relations for the $k=2$ case to those of the single solid torus (with $k=1$ in Eqs. (28) and (29), respectively). The results are shown in Fig. 10 for the modes $m=2$ and $n=1,2,3$ considered for fixed aspect ratios $\cosh\mu_{1},\cosh\mu_{2}$. These solutions in the $\varepsilon_{1}\varepsilon_{2}$ plane are obtained by using the corresponding dispersion relations Eqs. (28) and (29) with $k=2$. We then consider the values when $\varepsilon_{1}=\varepsilon_{2}=\varepsilon$. This corresponds to the case of a single solid torus since the similar $\varepsilon_{1},\varepsilon_{2}$ values in a $2$-layered torus corresponds to a single solid torus with dielectric constant $\varepsilon$ and fixed aspect ratio $\cosh\mu_{2}$. The common epsilon values $1-\varepsilon_{1}=1-\varepsilon_{2}=1-\varepsilon$ are marked on the $\varepsilon_{1}=\varepsilon_{2}$ line.

The values marked on the $\varepsilon_{1}=\varepsilon_{2}$ line in Fig. 10 agree well with the values generated from the solutions for a single solid ring as shown in Fig. 3. For example in Fig. 10 for the mode $m=2$ and $n=1$, the common $\varepsilon$ values are marked at $2.34$, $2.433$ for Eqs. (28), (29), respectively. The fixed outer aspect ratio in this case is $\cosh\mu_{2}=1.6$. The agreement is seen in Fig. 3, at the fixed aspect ratio $\cosh\mu_{2}=1.6$ for the mode $m=2$ and $n=1$. It can be seen that these marked values are in agreement with the values corresponding to the single solid torus solutions selected based on the ($m,n$) at the desired aspect ratio. This is true for all other values marked in Fig. 10.

Another interesting numerical experiment regards a singly coated torus with a fixed $\mu_{1}$ but varying $\varepsilon_{1}$, $\varepsilon_{2}$, and $\mu_{2}$. As a specific example of a physical arrangement modeled by this case, one may note a closed form of coated carbon nanotube, for which plasmon excitation may be studied in a similar fashion as that for a nanotube stokli . Thus, considering a coated ring in vacuum ($\varepsilon_{3}=1$), Fig. 11 shows the dispersion relation (Eq. (28) with $k=2$) for mode $m=4$ and $n=0$ in the $\varepsilon_{1}\varepsilon_{2}$ plane for a fixed aspect ratio $\cosh\mu_{1}=5$. As can be seen a redistribution of the dielectric values occur as the outer aspect ratio $\cosh\mu_{2}$ approaches the inner $\cosh\mu_{1}$, that is, the $\varepsilon_{1}$ range narrows down while the $\varepsilon_{2}$ range broadens up within the considered dielectric constant value ranges.

Surface modes of coated rings by perturbation technique

As in the case of a single ring, significant simplification can be obtained when treating the eigenmodes of the coated ring using the perturbation technique. Thus, writing the function $f$ in Eq. (II.1) as $f=\sqrt{u-\xi\cos\eta}$, and specializing Eq. (II.1) for $k=2$, and decomposing the parity of the poloidal solutions, we consider the odd surface modes, i.e., $\sin n\eta$ modes with fixed $m$ (the even modes $\sin n\eta$ may be treated similarly). At zeroth order, $\xi=0$ and the modes decouple. Let us concentrate on the $n$th mode. Using the conditions in Eq. (8) and Eq. (9) for $i=1,2$, and assuming $\varepsilon_{v}=\varepsilon_{1}/\varepsilon_{3}$ and $\varepsilon_{c}=\varepsilon_{2}/\varepsilon_{3}$, we have:

[TABLE]

and

[TABLE]

[TABLE]

This system has a solution if the determinant vanishes,

[TABLE]

where

[TABLE]

The solutions can thus be obtained following these results, as discussed for the special case of a thin toroidal membrane in the next section.

Plasmon hybridization and sum rules for single toroidal shell

The utility of the calculated expressions describing the generalized dispersion relations can be appreciated when considering more elaborate toroidal configurations. Here, we consider a single shell cavity structure, for which we can obtain the eigenvalues corresponding to a variety of parameter configurations.

These modes describe the out-of-phase and in-phase linear combinations of the induced screening charges when the shell is subjected to an electromagnetic disturbance such as an incident photon, electron, or ion. Therefore, polarization effects prevent the direct application of the previously considered sum rule, that is, the sum of the square frequencies of the complementary cavity, the shell, and the solid ring does not add up to those of the square bulk $\omega^{2}_{p}$. Instead, the frequencies of the lower energy symmetric $\omega_{s-}$ and the higher energy antisymmetric $\omega_{s+}$ modes form the new sum rule:

[TABLE]

where $\omega_{r}$ and $\omega_{v}$ now correspond to the plasmon frequencies of the original surfaces making up the shell, that is, the frequencies of a solid ring with a radius $r$ equal to the radius of the outer shell surface, and a cavity with a radius equal to the radius of the shell’s inner surface. The continued fractions in Eqs. (28) and (29) for $k=2$, generate all the needed frequencies to verify the sum rule. Furthermore, within the zeroth order, one readily extracts, from Eq. (77), a second order polynomial behavior for $\varepsilon$, from which we obtain the frequencies of the symmetric and antisymmetric charge density oscillations:

[TABLE]

where the energy splitting for a vacuum bounded shell is given by:

[TABLE]

where $\mathcal{A}$ and $\mathcal{B}$ are as given in the SMVIII.

Since the shell possesses two distinct surfaces, any excited modes will couple unless for large thicknesses of the shell so that any engendered screen charges on the inner and outer surfaces would decouple. In the latter case, the structure of the modes would be expected to approach those of either $\omega_{r}$ or $\omega_{v}$. These energy considerations can be captured within the plasmon hybridization picture prodan , as shown in Fig. 12. Characterizing the plasmon mode of the bulk with the energy state $|\omega_{p}\rangle$, and forming a single interface, the plasmon sum rule was shown to be satisfied by the energy of the new surface modes $|\omega^{mn}_{r}\rangle$ and $|\omega^{mn}_{v}\rangle$. Forming a second interface to create the shell, the modes $|\omega^{mn}_{r,v}\rangle$ hybridize to form the states $|\omega^{mn}_{s\pm}\rangle$, as described in Fig. 12. For $(m,n)=(4,1)$ as an example, the numerical annotation in Fig. 12 shows, the higher mode energy of the void $\omega_{v}/\omega_{p}=0.754$, the lower energy of the solid ring $\omega_{r}/\omega_{p}=0.646$, and an energy splitting of $(\omega_{s+}-\omega_{s-})/\omega_{p}=0.34$ between the hybridized states. Specializing the exact Eqs. (28) and (29) for $k=2$, $m=3$, and $n=1,2,3$, as well as the 0-order perturbative solutions Eq. (79), Fig. 13 demonstrates that the plasmon sum rule is satisfied in all cases.

For shells made of specific materials, the plasmon dispersion relations may be obtained following the same numerical procedure as before. For example, for aluminum and silver shells in vacuum, the corresponding resonance frequencies are shown for $\cosh\mu_{1}=5$ in Fig. 14. Employing the presented approach, other useful shell configurations may be studied, as discussed further in SM IX. Notably, the case of two shells ($k=5$) is amenable to the hybridization picture, in which the plasmon modes $|\omega^{mn}_{s\pm}\rangle$ are seen to further hybridize and form states $|\omega^{mn}_{s\pm\pm}\rangle$ creating highly blue-shifted antibonding-antibonding $|\omega^{mn}_{s++}\rangle$ and red-shifted bonding-bonding $|\omega^{mn}_{s--}\rangle$ plasmon modes.

II.2.3 Metal-dielectric stratified toridal medium

Similar to applications capitalizing on the photon scattering and propagation properties in layered media (e.g., interference filters and field enhancement in dielectric or metalo-dielectric stacks lereu ), insight into the mode structures and resonance spectra of multilayered tori can be useful. Using the general approach presented here, preliminary information can indeed be obtained for specific cases that are sufficiently practical to be of relevance in applications, such as, a three layered and a four layered torus. For example, with a fixed $\varepsilon_{2}=1$, fixed aspect ratios $\cosh\mu_{1},\cosh\mu_{3}$, and variable aspect ratio $\cosh\mu_{2}$, the dispersion relations (Eq. (28) with $k=3$) describing the toroidal mode $m=4$ and $n=0$ may be obtained for the three layered torus in vacuum ($\varepsilon_{4}=1$). Fig. 15 displays the results in the $\varepsilon_{1}\varepsilon_{3}$ plane for fixed aspect ratios of $\cosh\mu_{1}=5,\cosh\mu_{3}=3$.

It is seen that as $\cosh\mu_{2}\to\cosh\mu_{3}$, the eigenvalues redistribute, that is, the $\varepsilon_{1}$ range narrows down while that of $\varepsilon_{3}$ broadens up. In Fig. 15 one can make an additional observation that there seems to be a crossing point or small region in the neighborhood of $\varepsilon_{1}=-4$, $\varepsilon_{3}=-2$ where there is always a solution irrespective of (or weakly dependent on) the aspect ratio. This can be potentially interesting for fabrication and design consideration.

Similarly, one may consider a four layered torus in which two concentric shells are separated by vacuum so that $\mu_{1}$, $\mu_{4}$ are fixed and $\varepsilon_{1}=\varepsilon_{3}=1$. The results are shown in the $\varepsilon_{2}\varepsilon_{4}$ plane for various $\cosh\mu_{2}$ and $\cosh\mu_{3}$ values as shown in the Fig. 16. The dispersion relations were obtained from Eq. (28) with $k=4$ for the toroidal mode $m=4$ and $n=0$. Other effects such as the dependence on the distance or layer thickness, as shown in Fig. 4(a) in SM IX, or obtaining positive values satisfying Eq. (28), as shown in Fig. 4(b) in SM IX may be studied. As can be seen, these results suggest that the geometric parameters $\cosh\mu_{i}$, as well as the material parameters $\varepsilon_{i}$ may be used to engineer the eigenvalue distribution that best serves a given application.

II.2.4 Janus rings

Janus nanoparticles honegger ; zentgraf possess unique optical and electronic properties that are finding important applications. We here introduce and analytically formulate two new janus systems with different asymmetric properties, which have not appeared previously. The introduced janus particles may be considered special cases of the inhomogeneous toroidal domain shown in Fig. 17. For the sake of illuminating their electromagnetic response, we here suffice with a computational determination of their fields as a full treatment of their surface modes are beyond the scope of the present article, for reasons noted below.

The $z$-type janus ring is formed when $\eta_{r}-\eta_{l}=\pi$, that is, composition of two half rings described with dielectric functions $\varepsilon_{1}(\omega)$ and $\varepsilon_{2}(\omega)$, respectively, as depicted in Fig. 17 and Fig. 18. To formally define the problem, we consider the potentials inside and outside such a janus torus:

[TABLE]

The boundary conditions which $\Phi_{\text{in}}^{1}$, $\Phi_{\text{out}}^{2}$ must satisfy on $\mu_{1}\leq\mu<\infty$ reduce to

[TABLE]

While on $\mu=\mu_{1}$ the potential equations $\Phi_{\text{in}}^{1}$, $\Phi_{\text{in}}^{2}$ and $\Phi_{\text{out}}$ must satisfy

[TABLE]

The $y$ and $z$ are as given in Eq. (16) in SM III. In case of $z$-type janus rings we choose $\eta=0$, $z=0$. To obtain the normal modes and resonance frequencies associated with the charge density oscillations, it would be prudent to first establish whether for each $m=0,1,2,\cdots,$ the functions \big{\{}Q_{n-\frac{1}{2}}^{m}\big{\}}_{n=0}^{\infty} are pairwise orthogonal in $L^{2}(\sqrt{\cosh\mu-1}d\mu)$, the space of square integrable functions with respect to the measure $\sqrt{\cosh\mu-1}d\mu$. To the best of our knowledge, such an orthogonality relation is not known. Also care must be taken with singularities that may arise from the boundary condition Eq. (83), see also Eq. (IV)-(IV) in SM IV. However, bypassing the analytical complexity, we may obtain the eigenmodes and field distributions computationally employing FEM or the FDTD.

For example, for janus rings made up of gold palik and silicon palik , employing FEM to solve the frequency domain form of Maxwell’s wave equation, we computationally visualize the induced field and current distributions for an incident various linearly polarized fields of wavelength 652.6 nm, as shown in Fig. 18. The field and current density confinement and localization properties of the janus effect are clearly seen from the comparison to the uniform Au ring for two polarization.

The dissimilarity in the material domains in the $\varphi$-type janus ring takes the form of two half rings with transition regions satisfying $\varphi_{r}-\varphi_{l}=\pi$, and described by dielectric functions $\varepsilon_{1}(\omega)$ and $\varepsilon_{2}(\omega)$, respectively, as depicted in Fig. 17, and Fig. 18. In Eq. (II.2.4) we now have $0\leq\varphi<\pi$ and $\pi\leq\varphi<2\pi$. To satisfy the boundary conditions, we consider $\varphi=0$, $y=0$. The same arguments as in the case of $z$-type janus rings would be true for $\varphi=0$ case as well, but now one needs to establish the orthogonality condition for functions \big{\{}Q_{n-\frac{1}{2}}^{m}\big{\}}_{m=0}^{\infty} with $n$ fixed. Altogether Eqs. (82)–(83) and (84)–(85) shall result in three difference equations which can be transformed into a matrix vector product form as in Eq. (14). One may now follow the same procedure described in section II to obtain proper dispersion relations with needed symmetry checks performed. For the $\varphi$-type janus rings, field modifications under the same conditions as the $z$-type, are shown in Fig. 18. Noting that the contrasts have been adjusted individually in Fig. 18 for visualization purposes, the silicon domain of the Au-Si janus rings (both types) is marked with minimum field and current when compared with the uniform material (Au-Au) cases. Following the above formulation, the plasmon dispersion relations may be obtained.

III Eigenmodes and field distribution

Nanofocusing and field gradient engineering in multilayer rings are of particular interest due to toroidal and poloidal geometric features of the metallic, dielectric, and metallo-dielectric transitions. In the quasi-static limit, photon scattering by the multilayered torus can be formulated as an inhomogeneous Neumann problem. Here, considering a multilayer nanoring, we obtain a three-term difference amplitude equation with source terms characterizing the external radiation. We will then employ a Green function approach to solve for the amplitudes and thus the potential and fields. We will now proceed to obtain the appropriate Green’s function for the general case of a multilayer ring.

III.1 Green function for the response of a stratified toroidal medium to nonuniform fields

In order to study the scattering properties of the studied system, we begin by treating a $k$-layered ring in presence of electric monopoles and dipoles of charge $q_{s}$ at $(\mu_{0,s},\eta_{0,s},\varphi_{0,s})$ outside the toroidal boundary $\mu=\mu_{k}>\mu_{0,s}$. Expanding the electrostatic Green function for free space and superimpose the potential of $q_{s}$ to the exterior potential of the source-free ring, we can write:

[TABLE]

where $s$ counts the number of charges, $\delta$ denotes Kronecker delta, $j=1,\cdots,k+1$ and $\eta^{j}_{{}_{mn,s}}$, $\eta^{\prime j}_{{}_{mn,s}}$ are phase factors that are not necessarily coherent with $\eta_{0,s}$ and

[TABLE]

and where $\varepsilon_{0}$ is the permittivity of free space. As in the case of multiple charges kuk:cpc , Eq. (III.1) together with the Green’s function discussed in this section, allows us to compute inside and outside potentials in case of a multilayer nanoring structures, when multiple charges are placed outside the toroidal boundary $\mu=\mu_{k}$.

One can evaluate the potential represented in Eq. (III.1) by computing $C^{\text{j}}_{mn}$, $D^{\text{j}}_{mn}$ for individual charges. Thus in the following equations we adopt the notation $C^{\text{j}}_{mn}=C^{{}^{\prime}\text{j}}_{mn,s}e^{-in\eta^{j}_{{}_{mn,s}}},D^{\text{j}}_{mn}=D^{{}^{\prime}\text{j}}_{mn,s}e^{-in\eta^{\prime j}_{{}_{mn,s}}},K_{mn}=K^{{}^{\prime}}_{mn,s}e^{-in\eta_{0,s}}$. We suppress index $m$ in the same way as described in Section II.A. Introducing the $k\times k$ matrix

[TABLE]

and the $k\times 1$ vector

[TABLE]

and applying the boundary condition Eq. (8), one obtains the relation

[TABLE]

Using the substitution (90) into the second boundary condition (9) and after some algebra pap1 , one can show that $\mathcal{D}_{n}$ satisfies the vector three-term recurrence relation:

[TABLE]

for $n=0,\pm 1,\pm 2,\dots,$ where

[TABLE]

and all other used notations, except for $\mathcal{C}_{n}$ given by Eq. (90), are in accordance with the ones already introduced in Section II.A. We note that the three-term recurrence (91) holds for each fixed $m=0,\pm 1,\pm 2,\cdots$.

Generalizing Green’s function method love:jmp ; kuk , we obtain $\mathcal{W}_{n}$ in Eq. (91) by noting that

[TABLE]

forms a solution to Eq. (91) provided that for each $N$, the $k\times k$ matrices $G_{n,N}$ satisfy

[TABLE]

for $n=0,\pm 1,\pm 2,\dots,$ where $I$ denotes the $k\times k$ identity matrix. Using Kronecker delta’s definition at $n=N-1,\ N,\ N+1$, Eq. (III.1) can be written as:

[TABLE]

For other $n$ values, we have

[TABLE]

Suppressing the index $N,$ the three term matrix recurrences (98) and (99) display a single recurrence

[TABLE]

where either $n\geq N+2$ or $n\leq N-2$. Since the MCF

[TABLE]

is the tail of the converging MCF (26), it follows from the italic statement preceding Eq. (27) that the three-term recurrence (98) has a minimal solution. Letting $G_{n,N}$ denote this minimal solution, Eq. (27) yields

[TABLE]

Denoting the right-hand side of (102) by $\boldsymbol{\gamma}_{n+1},$ we get

[TABLE]

For the backward recurrence (100) with $n\leq N-2,$ a similar argument presented in section II.B.1.a together with the substitution $V_{n}=X_{n-1}X_{n}^{-1}$ gives the following first order non-linear matrix recurrence:

[TABLE]

Using the fact $\text{R}_{n}=\text{R}_{-n},$ the backward MCF:

[TABLE]

equals

[TABLE]

The last MCF converges using the same argument discussed above for the MCF (101). Thus, the backward recurrence has a minimal solution. Denoting this minimal solution by $G_{n,N},$ it follows from Eq. (104) that

[TABLE]

Similarly, denoting the right-hand side of (105) by $\boldsymbol{\zeta}_{n-1},$ we get

[TABLE]

Setting $n=N+1$ and $n=N-1$ in (103) and (106), respectively, one may now find $G_{N+1,N},$ $G_{N,N},$ and $G_{N-1,N}$ via Eqs. (95)–(97) as kuk :

[TABLE]

Finally, for each fixed $n,$ one can use Eqs. (103) and (106) together with Eqs. (107)–(109) to obtain $\mathcal{W}_{n}$ via (93). Since $\mathrm{J}_{n}$ in Eq. (14) is invertible, one can rewrite Eq. (14) in terms of $\mathcal{D}_{n}$ as

[TABLE]

Thus one can obtain the $k\times 1$ vector $\mathcal{D}_{n}$ as considered in Eq. (1). On the same lines, since $\mathbb{P}_{n}$ is invertible, one can use Eq. (90) to obtain another $k\times 1$ vector $\mathcal{C}_{n}$ as presented in Eq. (1). The obtained vector values can now be substituted into the potential Eq. (III.1).

III.2 Green function for the difference equation of the multilayer nanoring in a uniform field

III.2.1 Polarization parallel to the symmetry axis of the ring

The general form of the potential for a uniform field polarized parallel to $z$ axis outside the toroidal boundary $\mu=\mu_{k}$ is given by

[TABLE]

where $j=1,\cdots,k+1$ and $\delta$ denotes Kronecker delta and

[TABLE]

similar to the uniform electrostatic field defined in Love love:jmp . Because of the axial symmetry in case of a uniform field along $z$ axis, the potential is independent of the coordinate $\varphi$. Hence we do not have $m$-sums in Eq. (III.2.1). Also there are no phase differences in the $\eta$ solutions and hence they are just given by $e^{in\eta}$ everywhere. One could now follow the same procedure defined through Eqs. (88)–(109) to obtain Eqs. (107)–(109) but with $K_{n}$ as defined in Eq. (112). The rest of the procedure to obtain $\mathcal{D}_{n}$ and $\mathcal{C}_{n}$ is discussed following Eq. (109).

III.2.2 Polarization not parallel to the symmetry axis of the ring

The general form of the potential corresponding to a field with a $x-z$ polarization axis outside the toroidal boundary $\mu=\mu_{k}$ is given by

[TABLE]

where $j=1,\cdots,k+1$ and $\delta_{j,k+1}$ denotes Kronecker delta and

[TABLE]

The derivation of the applied potential

[TABLE]

in Eq. (III.2.2) is discussed in SM III. The method to obtain Green’s function and coefficients in Eq. (III.2.2) follows the same lines to the procedure discussed in previous subsections.

III.3 Poloidal and toroidal modes in response to specific excitations

III.3.1 Response of a vacuum bounded solid ring to nonuniform fields

The torus develops induced surface charges in response to the field of an external charge such as an electron at $(\mu_{0},\eta_{0},\varphi_{0})$ outside the toroidal boundary $\mu=\mu_{1}>\mu_{0}$, which can be calculated after solving the second-order difference equation (91) using the Green function with $k=1$. The combined potentials now read:

[TABLE]

where

[TABLE]

and $K_{mn}=K^{{}^{\prime}}_{mn}e^{-in\eta_{0}}$ with $K^{{}^{\prime}}_{mn}$ as defined in Eq. (III.1) and $\theta(x)=1$ if $x>0$, and 0 otherwise. The total potential obtained by summation of $m$ and $n$ modes as in Eqs. (III.3.1) and (III.3.1) is shown in Fig. 5 in SM X. The poloidal and toroidal mode potential distribution induced by a nearby dipole using Eq. (III.1) by considering only relevant modes is shown in Fig. 19. The potential (encompassing contributions from toroidal and poloidal modes in Eq. (III.1)) is shown in Fig. 20 and in the inset of Fig. 21. The plotting procedure for these figures is briefly discussed in SM X.

The utility of the results may be appreciated if we consider forming a dipole of moment $p$ and from its particular placements at ${\bf r}_{0}=(\mu_{0},\eta_{0},\varphi_{0})$ and orientation build the dyadic Green’s function:

[TABLE]

where $c$ is the speed of light and $\omega$ is the frequency at which $\varepsilon$ is specified.

For example, setting ${\bf r}_{0}$ to form a dipole parallel to the symmetry axis of the ring, then the partial ($z$) local density of states (LDOS), which is directly proportional to the spontaneous decay rate $\Gamma$, can be obtained from

[TABLE]

where $\alpha=3\hbar/2\omega p^{2}$. Fig. 21 shows the modification of the radiative decay rate due to the presence of the gold ring. The higher local density of states corresponding to the resonant response of the ring to the two polarization states of the exciting field (parallel and perpendicular to the symmetry axis of the ring) is clearly observable. Employing Eq. (III.1) to obtain the response of the ring to an electric dipole with an arbitrary moment $\bar{d}$ we may thus obtain the partial local density of states.

III.3.2 Response of a vacuum bounded solid ring to uniform fields

In case of potential for a uniform field polarized along the $z$-axis, the corresponding potentials $\Phi_{1}$ and $\Phi_{2}$ for the regions inside and outside $\mu=\mu_{1}$ respectively, after solving second-order difference equation (91) using the Green function with $k=1$ are given as

[TABLE]

where $d_{{}_{n,N}}$ is as shown in Eq. (III.3.1) but with $m=0$ considered and $K_{n}$ is as defined in Eq. (112).

The toroidal and poloidal mode potential distribution obtained by considering relevant modes in Eqs. (118)–(III.3.2) is further discussed in SM X.

III.4 Scattering cross section

To calculate the scattering cross section for a single ring, here we shall consider the equations of a uniform field of an arbitrary polarization (subject to the conditions Eqs. (8), (9) with $k=1$). Similar procedure discussed in section III.A is followed and finally after solving the second-order difference equation Eq. (91) using the Green function with $k=1$, the corresponding potentials $\Phi_{1}$ and $\Phi_{2}$ for the regions inside and outside $\mu=\mu_{1}$ respectively, are obtained:

[TABLE]

where $d{{}_{m,n,N}}$ is as defined in Eq. (III.3.1) with

[TABLE]

The derivation of the applied potential

[TABLE]

in Eq. (III.4) is discussed in SM III.

The induced surface charge density on the ring in the limit of dipolar contribution ($m=0,1$) is therefore

[TABLE]

which is connected to the volume charge density via

[TABLE]

yielding the induced dipole moment $\vec{p}$

[TABLE]

where $x$, $z$ are Eqs. (14), (16) in SM III and $h_{\mu}$, $h_{\eta}$, $h_{\varphi}$ are presented in Eq. (21) in SM IV. The induced dipole gives rise to a scattered field such that the total scattering cross section can be calculated. Using this result, the simulated spectrum are shown in Fig. 22 for a vacuum bounded gold ring. For comparison, the computationally obtained cross sections for both types of janus nanorings and a metal nanosphere are presented.

IV Conclusion

In conclusion, the introduction of a composite system of ring toroids to model multilayered metal, dielectric, metal-dielectric, and janus nanorings successfully produced the sought surface modes that may be exhibited by such a system. In particular, the acquisition of the plasmon dispersion relations either analytically, up to a three-term difference equation including its associated continued fraction and infinite tridiagonal determinants, or numerically, visualized in the form of specific solution branches, was shown to produce the resonance spectra of the multilayered ring system. Specifically, using the simulated resonance values corresponding to the dominant (dipole) modes of a simple ring isolated in vacuum, reasonable agreement with spectral peak positions from the computationally determined absorption cross sections of the ring was demonstrated using FDTD. The agreement was demonstrated for several materials (gold, silver, aluminum and silicon) using experimental values of their dielectric functions. Furthermore, the analysis of the initial conditions for the obtained three-term difference equation and the consideration of its convergence properties, clarified the ambiguity of the reported surface modes in the literature. Therefore, exact quasi-static surface plasmon dispersion relations were calculated for a composite ring and can be simulated following Eq. (28) and Eq. (29) presented here for various configurations of current interest in plasmonics, nanophotonics, and metamaterial research. The agreement with the computationally obtained retarded dispersion relations in the case of a single ring, further implies that the presented results do offer working estimates of the resonance spectra of solid nanorings. Additionally, for multiple cases simulating the analytical results from a first order perturbation theory and visualizing and comparing the leading dispersion branches demonstrated the viability and accuracy of obtaining the exact dispersion relations via the presented results. Additionally, the perturbation solutions were instrumental for calculating the substrate induced plasmon shift, as well as for analytical observation of the energy splitting due to hybridization between the induced bonding and anti-bonding screening charges at the two interfaces of a toroidal shell. The presented formulation and modeling of the multilayered rings also proved useful in obtaining the Green’s function for rings with an arbitrary number of layers. Therefore, field distribution and scattering properties of composite nanorings may be studied when considering interaction with arbitrarily polarized photons, which was shown to amount to modifying a source term appearing in the three-term difference equation. The correctness of the obtained Green’s function was demonstrated for a single ring interacting with a uniform field with polarization parallel to the ring symmetry axis and perpendicular to it, as well as for a nonuniform field, resulting for example from electrons and dipoles. Therefore, following the results presented in Sec. III.3, one may obtain the LDOS for applications requiring the determination of the spontaneous decay rate of an emitter, such as a two level or four level quantum system in presence of multilayer nanorings and toroidal nanoshells. We demonstrated this viability by determining the partial LDOS and the corresponding spontaneous decay rate for the simpler system of a transition dipole within the nearfield of a gold ring. In addition to providing the resonance spectra, plasmon dispersion relations and fields, the presented results are of importance for use in validating computational models, in particular since retarded quantities are analytically unavailable due to the inseparability of the wave equation in the toroidal geometry. Finally, the introduction of the two types of ring based janus particles may lead to new opportunities in nanofabrication, which could trigger novel applications in sensing, quantum applications salhi , and metamaterials. In the presented view, in the case of $z$-type janus rings, for any non-negative integer $m$, one will have to establish a new orthogonality relation of the following kind

[TABLE]

for all non-negative integers $n$ and $n^{\prime}$, where $C_{mn}$ only depends on $m$ and $n$. This orthogonality, which to the best of our knowledge is being reported for the first time, may need to hold in case of $\varphi$-type janus ring but now with fixed $n$, varying $m$ and weight function $\sqrt{\cosh\mu-\cos\eta}$ considered. The proof of this orthogonality will be the subject of a forthcoming article.

V Acknowledgements

This work was supported by the laboratory directed research and development (LDRD) fund at Oak Ridge National Laboratory (ORNL). ORNL is managed by UT-Battelle, LLC, for the US DOE under contract DE-AC05-00OR22725.

Supplemental Materials (SM)

I

Formulas in SM I hold for $n=0,\pm 1,\pm 2,\cdots$.

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

II

It follows from Eq. (II.1) and Eq. (II.1) that the general solution is a linear combinations of the modes

[TABLE]

where $L_{n-\frac{1}{2}}^{m}(\mu)$ denotes either of the associated Legendre functions $Q_{n-\frac{1}{2}}^{m}(\cosh\mu)$ or $P_{n-\frac{1}{2}}^{m}(\cosh\mu)$. Letting ∗ stand for the complex conjugation, we have

[TABLE]

where

[TABLE]

Using the orthogonality of $\{e^{in\eta}\}$ in Eq. (II)

[TABLE]

For the sake of the discussions, we note that for $\cosh\mu\rightarrow 1$, the behavior of the toroidal harmonics resembles that of spherical harmonics as determined by using formulas $8.751-1$ and $8.751-2$ in Gradshteyn et al. Grad On the other hand, when $\cosh\mu\rightarrow\infty$, the toroidal harmonics resemble cylindrical harmonics and are represented using modified Bessel functions of the first and second kind of order $n$ which are the solutions to Laplace’s equation in cylindrical coordinates. Such cylindrical limit solution is discussed in Love love:plasma . Various limiting cases of a ring is shown in Fig. 1.

III

The potential for a uniform field polarized along $x-z$ axis in toroidal coordinates is given by:

[TABLE]

For a torus in the toroidal coordinate system $(\mu,\eta,\varphi),$ by Morse & Feshbach morse on page 1304,

[TABLE]

Using the absolute convergence of the above series, differentiation w.r.t $\cosh\mu$ yields

[TABLE]

Using 8.732, 8.734, and 8.736.4 from Gradshteyn & Ryzhik Grad

[TABLE]

Letting $z=\cosh\mu$ in the last series Eq. (13) and using the above relation together with the fact that

[TABLE]

we obtain the desired expansion

[TABLE]

From the known expansion love:jmp for $z$, we have

[TABLE]

Thus,

[TABLE]

where $K_{mn}$ is as given in Eq. (III.2.2).

Furthermore we have

[TABLE]

Similar to Eq. (III), we now have

[TABLE]

The potential for a uniform field along $x-y-z$ axis in toroidal coordinates is given by:

[TABLE]

where $K_{mn}$ is as given in Eq. (III.4) with the assumption that $\mathcal{E}_{3}=\beta\mathcal{E}_{2}$. Using the above obtained expression for the applied potential $\Phi_{ap}$, the superposed potential for a uniform field of arbitrary polarization can be written as

[TABLE]

for $0\leq\mu\leq\mu_{1}$.

IV

The scale factors for the toroidal coordinates are

[TABLE]

Using the chain rule on the partial derivative terms involved in Eq. (83) gives

[TABLE]

Furthermore, differentiating $x$, $y$, $z$ terms in Eqs. (14), (16) gives

[TABLE]

V

Consider the system

[TABLE]

The infinite determinant of the above matrix does not exist, but yet the system has non-trivial solutions: $f_{0}=1,f_{1}=-1,f_{2}=0,f_{3}=1,\cdots$.

VI

It is instructive to consider some general points regarding the character of the charge density in the absence of applied fields. Here, we shall consider a solid ring without external sources, and therefore the charge density (confined to the surfaces) should integrate to zero. Thus, the charge density $\rho$ must satisfy Eq. (1). We will now discuss the potential satisfying the Poisson equation

[TABLE]

The discussion related to the Poisson equation where distributed charges of density $q$ may be present in a region $\mu>\mu_{1}$ can be found elsewhere zotter . In case of a solid ring, the electrostatic potential given in Eq. (II.1) becomes

[TABLE]

and considering Eq. (II.1) in general terms, Eq. (VI) can be arranged as a linear combination of the product solutions

[TABLE]

where $\mathcal{A}_{mn}$ is the corresponding amplitudes, $U_{mn}(\mu)$ is an associated Legendre function $Q_{n-\frac{1}{2}}^{m}(\cosh\mu)$ or $P_{n-\frac{1}{2}}^{m}(\cosh\mu)$, $T_{n}(\eta)$ is $\sin{n\eta}$ or $\cos{n\eta}$, and $V_{m}(\varphi)$ is $\sin{m\varphi}$ or $\cos{m\varphi}$. One makes the potential continuous across the boundary $\mu=\mu_{1}$ by writing

[TABLE]

and

[TABLE]

where $\delta(\mu)$ is the Dirac delta function. The Laplacian of $\Phi$ in toroidal coordinate system is

[TABLE]

We shall insert the basis in Eq. (VI) which is

[TABLE]

into Eq. (VI) and obtain weis the corresponding Laplace equation

[TABLE]

An alternative separation andrew of Laplace’s equation in toroidal coordinates can be obtained by transforming Eq. (VI) into

[TABLE]

where $r$ is the distance from the $z$-axis given by

[TABLE]

and $x$ and $y$ are as given in Eqs. (14) and (16) in SM III. Further substituting bye $U=\Phi\sqrt{r}$ into Eq. (VI), one would obtain the Laplace equation

[TABLE]

The solutions of the above Laplace equation form a complete basis andrew and are given by

[TABLE]

where $U_{pq}(\mu)$ is $Q_{p-\frac{1}{2}}^{q}(\coth\mu)$ or $P_{p-\frac{1}{2}}^{q}(\coth\mu)$, $T_{q}(\eta)$ is $\sin{q\eta}$ or $\cos{q\eta}$, and $V_{p}(\varphi)$ is $\sin{p\varphi}$ or $\cos{p\varphi}$. This basis is more appropriate for the situations where the boundary conditions do not depend on $\eta$. In such cases $q=0$ is considered and the solution involves the Legendre functions $Q_{p-\frac{1}{2}}(\coth\mu)$ or $P_{p-\frac{1}{2}}(\coth\mu)$. The basis in Eq. (31) is particularly convenient for axially symmetric situations, in which case we set $m=0$ and the solution involves the Legendre functions $Q_{n-\frac{1}{2}}(\cosh\mu)$ or $P_{n-\frac{1}{2}}(\cosh\mu)$. However the equivalence of the two approaches to separation of the variables is discussed in SM B in Andrews andrew .

The charge density $\rho(\mu,\eta,\varphi)$ corresponding to the potential of Eq. (VI) is given by Eq. (2) and $W_{mn}(\mu_{1})$ for each fixed $m$ is

[TABLE]

Substituting $\rho$ given by Eq. (2) into Eq. (1) and evaluating the integral with respect to $\mu$, one obtains

[TABLE]

where

[TABLE]

Evaluating the integral in Eq. (38) with respect to $\varphi$, one obtains

[TABLE]

where

[TABLE]

That is, to infer that a solid ring without external sources is axially symmetric. Simplifying integral expression in Eq. (39), we obtain

[TABLE]

One can now substitute Eq. (12) in SM III into Eq. (VI) and together with the orthogonality property of $T_{n}(\eta)$ on $0\leq\eta\leq 2\pi$, the following sufficient condition

[TABLE]

needs to hold in order to satisfy Eq. (1) for $n=0,1,\cdots$. As a final remark, we note that Eq. (41) further signifies the fact that for a given aspect ratio $\cosh\mu_{1}$, there are no roots for the dispersion relation Eqs. (28) and (29) when $m=0$ is considered. This can be observed from the fact that the left hand side and the right hand side of the equation obtained by applying boundary condition Eq. (9) for $k=1$ cancel each other out with Eq. (41) considered at $m=0$ mode. Notice that the superscripts in Eq. (41) indicate the $k$ layered index, in this case, $k=1$ and not to be confused with $m$ modes.

VII

Employing the FEM method, in Fig. 2 we present the computationally determined transient response of a gold nanoring to linearly polarized incident photons.

VIII

In case of a toroidal shell we assume $\varepsilon_{v}=1$ and $\varepsilon_{c}=\varepsilon$ in Eq. (77). The determinant in Eq. (77) generates a second order polynomial in $\varepsilon$:

[TABLE]

with

[TABLE]

IX

We briefly discuss plasmon dispersion in metallic shells with dielectric core. Instructive modes include the resonance values $\varepsilon_{2}$ of a single shell (in vacuum) as a function of the outer aspect ratio, assuming $\varepsilon_{1}=1$ as shown in Fig. 3(a), where the inner aspect ratio is fixed at $\cosh\mu_{1}=5$. The $x$-axis represents the shell thickness $d_{2}$ which is $d_{2}=r_{2}-r_{1}$ while the $y$-axis denote the roots of the dispersion relation Eq. (28) in the form of a dielectric function $\varepsilon_{2}$ for a toroidal surface mode of $m=4$ and $n=0$. Similarly, in Fig. 3(b) the outer aspect ratio is fixed at $\cosh\mu_{2}=1.05$ while the inner varied. The $x$-axis represents the inner ring radius $d_{1}=r_{1}$.

We allow the inner radius $r_{1}$ to approach 0 in order to establish its relation with $\varepsilon_{2}$ values satisfying the dispersion relation Eq. (28) for a toroidal surface mode of $m=4$ and $n=0$ with $k=2$. As noticed in Fig. 3(b), smaller ranges of $\varepsilon_{2}$ values with negative magnitudes tend to satisfy the dispersion relation Eq. (28) when the distance $d_{1}$ approaches towards zero. But when $d_{1}$ approaches $R=$1\text{,}\mathrm{nm}$$ as shown in the plot, the $\varepsilon_{2}$ values satisfying the dispersion relation Eq. (28) tend to gradually increase. This phenomenon is opposite to the one described in Fig. 3(a) wherein the $\varepsilon_{2}$ values gradually increase for smaller shell size thickness $d_{2}$.

Using the dispersion relation Eq. (28) with $k=4$ for the toroidal mode $m=4$ and $n=0$, other effects such as the dependence on the distance or layer thickness, as shown in Fig. 4(a), or obtaining positive values satisfying Eq. (28), as shown in Fig. 4(b) may be studied.

X

We briefly discuss the procedure used in plotting the non-uniform and uniform fields. To plot the non-uniform field, one may convert Eq. (III.3.1) into real Fourier series. A real Fourier series function in general can be expressed as

[TABLE]

Expressing $\cos{nx}$ and $\sin{nx}$ in exponential form, we may rewrite Eq. (44) as

[TABLE]

where $c_{n}=\frac{1}{2}(a_{n}-ib_{n})$ and $c_{-n}=\frac{1}{2}(a_{n}+ib_{n})$ for $n>0$ and $c_{0}=\frac{1}{2}a_{0}$. Hence $a_{n}$ and $b_{n}$ in Eq. (44) are given by $a_{n}=c_{n}+c_{-n}$ and $b_{n}=i(c_{n}-c_{-n})$. We now convert Eq. (III.3.1) into real Fourier series with respect to $\eta$. Thus, considering $c_{n}$ to be

[TABLE]

together with the following relations

[TABLE]

we obtain

[TABLE]

[TABLE]

In the same fashion we can convert Eq. (III.3.1) into real Fourier series by considering $c_{n}$ to be

[TABLE]

and using the fact that

[TABLE]

one can compute $a_{n}$, $b_{n}$ at the respective toroidal, poloidal modes (considering other modes to be irrelevant) and then use Eq. (44) to obtain the potential distribution $\Phi_{1}$, $\Phi_{2}$ of the poloidal and toroidal modes induced by a nearby electron. Similar procedure as above but using Eq. (III.1) accordingly one can obtain the potential distribution at the respective toroidal, poloidal modes induced by a dipole. Considering all other modes via summation of $m$ and $n$, results in the total potential distribution as shown in Fig. 5.

To plot the uniform field, We convert Eq. (118) into real Fourier series with respect to $\eta$ by considering $c_{n}$ to be

[TABLE]

together with the following relation

[TABLE]

we obtain

[TABLE]

In similar fashion we convert Eq. (III.3.2) into real Fourier series by considering $c_{n}$ to be

[TABLE]

and using the fact that $K_{n}=-K_{-n}$, one can use Eq. (44) to plot $\Phi_{1}$ and $\Phi_{2}$ at the respective toroidal and poloidal modes as shown in Fig. 6. The coloring scheme across the two cross-sections in case of uniform field in Fig. 6 depends on the evenness/oddness of $n$-mode. When odd, the color scheme differs at $\varphi=0$ or $\pi$ noticing the opposite loops on the ring cross-section and vice-versa. The total potential is then obtained by summation of $n$ modes as in Eqs. (118)–(III.3.2) is shown in Fig. 7.