---
abstract: 'We investigate the supersymmetry (SUSY) structures for inductor-capacitor circuit networks on a simple regular graph and its line graph. We show that their eigenspectra must coincide (except, possibly, for the highest eigenfrequency) due to SUSY, which is derived from the topological nature of the circuits. To observe this spectra correspondence in the high frequency range, we study spoof plasmons on metallic hexagonal and lattices. The band correspondence between them is predicted by a simulation. Using terahertz time-domain spectroscopy, we demonstrate the band correspondence of fabricated metallic hexagonal and lattices.'
author:
- Yosuke Nakata
- Yoshiro Urade
- Toshihiro Nakanishi
- Fumiaki Miyamaru
- Mitsuo Wada Takeda
- Masao Kitano
nocite: '[@*]'
title: ' Supersymmetric correspondence in spectra on a graph and its line graph: From circuit theory to spoof plasmons on metallic lattices '
---

\#1[[~~\#1~~]{}]{} \#1[[\#1]{}]{}

\#1[[Ø u ‘=‘Ø=]{}]{} \#1[\#1|]{} \#1[\#1|]{} \#1[|]{} \#1[|]{} \#1 \#1\#2[|]{} \#1\#2\#3[||]{} \#1[[ ]{}]{} \#1 \#1\#2 \#1[\_]{} \#1[\^]{}

Introduction
============

Supersymmetry (SUSY) is a conjectured symmetry between fermions and bosons. Although the concept of SUSY was introduced in high-energy physics and remains to be experimentally confirmed, the underlying algebra is also found in quantum mechanics. When the SUSY algebra is applied to the field of quantum mechanics it is called supersymmetric quantum mechanics (SUSYQM) [@Cooper1994]. The algebraic relations of SUSY link two systems that at first glance might seem to be very different. The linkage through SUSY can be utilized to construct exact solutions for various systems in quantum mechanics. Recently, SUSYQM has been applied to construct quantum systems enabling exotic quantum wave propagations: reflectionless or invisible defects in tight-binding models [@Longhi2010] and complex crystals [@Longhi2013a], transparent interface between two isospectral one-dimensional crystals [@Longhi2013], reflectionless bent waveguides for matter-waves [@Campo2014], and disordered systems with Bloch-like eigenstates and band gaps [@Yu2015].

The SUSY structure was also found in other physics fields besides quantum mechanics, e.g., statistical physics through the Fokker-Planck equations [@Bernstein1984]. Through the similarity between quantum-mechanical probability waves and electromagnetic waves, the SUSY structure can be formulated for electromagnetic systems. Electromagnetic SUSY structures have been found in one-dimensional refractive index distributions [@Chumakov1994; @Miri2013], coupled discrete waveguides [@Longhi2010; @Miri2013], weakly guiding optical fibers with cylindrical symmetry [@Miri2013], planar waveguides with varying permittivity and permeability [@Laba2014], and non-uniform grating structures [@Longhi2015]. Even a quantum optical deformed oscillator with $\mathrm{SU}(1,1)$ group symmetry and its SUSY partner were constructed as a classical electromagnetic system [@Zuniga-Segundo2014].

The SUSY transformation generates new optical systems whose spectra coincide with those of the original system (except possibly for the highest eigenvalue of the fundamental mode of original or generated systems). The SUSY transformations have been utilized to synthesize mode filters [@Miri2013] and distributed-feedback filters with any desired number of resonances at the target frequencies [@Longhi2015]. The scattering properties of the optical systems paired by the SUSY transformation are related to each other [@Longhi2010; @Miri2013]. It is possible to design an optical system family with identical reflection and transmission characteristics by using the SUSY transformations [@Miri2014]. A reflectionless potential derived from the trivial system by SUSY transformation was applied to design transparent optical intersections [@Longhi2015a]. Moreover, SUSY has also been intensively investigated in non-Hermitian optical systems. If a system is invariant under the simultaneous operations of the space and time inversions, it is called $\mathcal{PT}$-symmetric. The SUSY transformation for the $\mathcal{PT}$-symmetric system allows for arbitrarily removing bound states from the spectrum [@Miri2013a]. In addition, non-Hermitian optical couplers can be designed [@Principe2015]. By using double SUSY transformations, the bound states in the continuum were also formulated in tight-binding lattices [@Longhi2014; @Longhi2014a] and continuous systems [@Correa2015]. The SUSY transformation in the $\mathcal{PT}$-symmetric system can also reduce the undesired reflection of one-way-invisible optical crystals [@Midya2014].

From an experimental perspective, it is still challenging to extract the full potential of electromagnetic SUSY because of fabrication difficulties. However, using dielectric coupled waveguides, researchers have realized a reflectionless potential [@Szameit2011], interpreted as a transformed potential derived from the trivial one by a SUSY transformation [@Longhi2010], and SUSY mode converters [@Heinrich2014]. The SUSY scattering properties of dielectric coupled waveguides have also been observed [@Heinrich2014a].

As we have described so far, many studies have been done for the electromagnetic SUSY, but their focusing point is mainly limited to dielectric structures. Recent progress of plasmonics [@AlexanderMaier2007] and metamaterials [@Solymar2009] using metals in optics demands further studies of SUSY for metallic systems. To design and analyze the characteristics of metallic structures, intuitive electrical circuit models are very useful, because they extract the nature of the phenomena despite reducing the degree of freedom for the problem [@Nakata2012a]. Actually, a circuit-theoretical design strategy called [*metactronics*]{} has been proposed even in the optical region [@Engheta2007] and the circuit theory for plasmons has also been developed [@Staffaroni2012]. If we could design circuit models enabling exotic phenomena, they open up new possibilities for application to higher frequency ranges due to the scale invariance of Maxwell equations. Thus, in this paper we develop how SUSY appears in inductor-capacitor circuit networks and demonstrate the SUSY correspondence in the high frequency region. In particular, we focus on the SUSY structure for inductor-capacitor circuit networks on a graph and its line graph.

This article is organized as follows. In Sec. \[sec:2\], we start by introducing the graph-theoretical concepts and formulate a general class of inductor-capacitor circuit network pairs related through SUSY, derived from the topological nature of the graphs representing the circuits. In Sec. \[sec:3\], we theoretically and experimentally demonstrate the SUSY eigenfrequency correspondence for paired metallic lattices in the terahertz frequency range. In Sec. \[sec:4\], we summarize and conclude the paper.

Theory \[sec:2\]
================

Eigenequation for inductor-capacitor circuit networks 
------------------------------------------------------

We consider an inductor-capacitor circuit network on a simple directed graph $G=(\mathcal{V},\mathcal{E})$, where $\mathcal{V}$ and $\mathcal{E}$ are the sets of vertices and directed edges, respectively. The modifier [*simple*]{} means that there are no multiple edges between any vertex pair and no edge (loop) that connects a vertex to itself. The number of the edges connected to a vertex $v$ of a graph is called the degree of $v$. A [*regular*]{} graph is a graph whose every vertex has the same degree. We assume that $G$ is an $m$-[*regular*]{} graph with all vertices having degree $m$. The capacitors, all with the same capacitance $C$, are connected between each vertex $v\in \mathcal{V}$ and the ground. Coils, all with the same inductance $L$, are loaded along all $e \in \mathcal{E}$. An example of $G$ and the inductor-capacitor circuit network on it are shown in Fig. \[fig:lc\_ladder\](a) and (b).

![\[fig:lc\_ladder\] (a) Example of simple $3$-regular directed graph. (b) Inductor-capacitor circuit network on the graph. (c) Line graph of the graph shown in (a). (d) Inductor-capacitor circuit network on the line graph (c). ](lc_ladder.eps){width="86mm"}

For $v\in \mathcal{V}$ and $e\in \mathcal{E}$, the incidence matrix $\mathsf{X}=[X_{ve}]$ of a directed graph $G$ is defined as follows: $X_{ve}=-1$ ($e$ enters $v$), $X_{ve}=1$ ($e$ leaves $v$), otherwise $X_{ve}=0$.

Using vector notation, we represent the current distribution $J_e$ flowing along $e\in \mathcal{E}$ as a column vector $\vct{J}=[J_e]^\mathrm{T}$. The charge distribution is denoted by $\vct{q}=[q_v]^\mathrm{T}$ with a stored charge $q_v$ at $v\in \mathcal{V}$. The charge conservation law is given by $$\dot{\vct{q}}=-\mathsf{X} \vct{J},  \label{eq:1}$$ where the time derivative is represented by the dot. The scalar potential $\Phi_v$ at $v \in {\mathcal V}$ must satisfy Faraday’s law of induction, so we have $$\dot{\vct{J}}=\frac{1}{L} \mathsf{X}^\mathrm{T}\vct{\Phi} ,  \label{eq:2}$$ with $\vct{\Phi}=[\Phi_v]^\mathrm{T}$. The scalar potential $\vct{\Phi}$ is written as $$\vct{\Phi}=\mathsf{P} \vct{q},  \label{eq:3}$$ with a potential matrix $\mathsf{P}$. In our case, $\mathsf{P}$ is given by $$\mathsf{P}=C^{-1}\mathsf{I},  \label{eq:4}$$ where we use the identity matrix $\mathsf{I}$.

From Eqs. (\[eq:1\])–(\[eq:4\]), we obtain $$\ddot{\vct{q}}=-{\omega_0}^2\mathsf{X}\mathsf{X}^\mathrm{T} \vct{q}, \nonumber$$ with $\omega_0=1/\sqrt{LC}$. Assuming $\vct{q}=\tilde{\vct{q}}\exp(-\ii \omega t)+\cc$, we have an eigenequation $$\mathsf{L}\tilde{\vct{q}}= \left(\frac{\omega}{\omega_0}\right)^2\tilde{\vct{q}}  \label{eq:5}$$ with the Laplacian $\mathsf{L}=\mathsf{X}\mathsf{X}^\mathrm{T}$. We introduce an adjacency matrix $\mathsf{A}=[A_{vw}]$, where $A_{vw}$ is 1 if $v, w\in \mathcal{V}$ are connected by an edge, otherwise 0. From a direct calculation, we can write $\mathsf{L}$ by $\mathsf{A}$ as $$\mathsf{L}=\mathsf{X}\mathsf{X}^\mathrm{T}=-\mathsf{A}+m\mathsf{I},  \label{eq:6}$$ where $m$ is the degree of the vertex of $G$. Note that $\mathsf{L}$ is independent of the direction of the edges in $G$ because $\mathsf{L}$ is expressed in terms of $\mathsf{A}$ and $\mathsf{I}$.

The directed graph $G$ can be also regarded as an undirected graph. For $e\in \mathcal{E}$, we can make an undirected edge $\bar{e}$, where the bar operator ignores the direction of the edge. Then, we have $\bar{G}=(\mathcal{V},\bar{\mathcal{E}})$ with $\bar{\mathcal{E}}=\{\bar{e}|e\in \mathcal{E}\}$. We can also define the undirected incidence matrix $\bar{\mathsf{X}}=[\bar{X}_{ve}]$ as $\bar{X}_{ve}=1$ ($e$ and $v$ are connected), otherwise $\bar{X}_{ve}=0$. Using $\bar{\mathsf{X}}$, $\mathsf{A}$ is written as follows [@Biggs1994]: $$\bar{\mathsf{X}}\bar{\mathsf{X}}^\mathrm{T}=\mathsf{A}+m \mathsf{I}.  \label{eq:7}$$ From Eqs. (\[eq:6\]) and (\[eq:7\]), we obtain $$\mathsf{L}=- \bar{\mathsf{X}}\bar{\mathsf{X}}^\mathrm{T}+2m\mathsf{I}.  \label{eq:8}$$

SUSY correspondence in spectra on a simple regular graph and its line graph
---------------------------------------------------------------------------

Next, we introduce the line graph concept [@Biggs1994]. The line graph $L(G)=(\mathcal{V}\sub{L}, \mathcal{E}\sub{L})$ of a directed graph $G$ is constructed as follows. Each edge in $G$ is considered to be a vertex of $L(G)$. Two vertices of $L(G)$ are connected if the corresponding edges in $G$ have a vertex in common. There are two possible choices for the direction of each edge in $L(G)$ and we adapt one of them. From here on, we only consider $L(G)$ of a simple $m$-regular graph $G$. In this case, the line graph $L(G)$ is a simple $m\sub{L}$-regular graph. The degree $m\sub{L}$ can be represented by $m$. For a vertex $v \in \mathcal{V}$ included in $e\in \mathcal{E}$, there are $m-1$ edges $e'\in \mathcal{E}\setminus\{e \}$ connected to $v$. Then, we obtain $$m\sub{L}=2(m-1).  \label{eq:9}$$ Note that $(m/2)\#\mathcal{V}= \#\mathcal{E}=\#\mathcal{V}\sub{L}=(2/m\sub{L})\#\mathcal{E}\sub{L}$ is satisfied for a finite graph $G$, where $\#\mathcal{S}$ represents the numbers of the elements of the set $\mathcal{S}$. Figure \[fig:lc\_ladder\](c) is an example of the line graph of the graph $G$ shown in Fig. \[fig:lc\_ladder\](a). Figure \[fig:lc\_ladder\](d) is the inductor-capacitor circuit network on $L(G)$.

In the context of mathematics, it is known that the spectra of the Laplacians for the graph and its line graph are related to each other [@Shirai2000]. For the convenience of the readers, we rederive this property in a simple manner and apply it to the inductor-capacitor circuit networks. The Laplacian of $L(G)$ is written as $$\mathsf{L}\sub{L}=\mathsf{X}\sub{L}{\mathsf{X}\sub{L}}^\mathrm{T}=-\mathsf{A}\sub{L}+m\sub{L}\mathsf{I}\sub{L},  \label{eq:10}$$ with the identity matrix $\mathsf{I}\sub{L}$, the incidence matrix $\mathsf{X}\sub{L}$, and the adjacency matrix $\mathsf{A}\sub{L}$ for $L(G)$. The adjacency matrix of $L(G)$ is represented as follows [@Biggs1994]: $$\mathsf{A}\sub{L}=\bar{\mathsf{X}}^\mathrm{T} \bar{\mathsf{X}}-2\mathsf{I}\sub{L}.  \label{eq:11}$$ From Eqs. (\[eq:10\]) and (\[eq:11\]), we have $$\mathsf{L}\sub{L}=-\bar{\mathsf{X}}^\mathrm{T} \bar{\mathsf{X}}+(m\sub{L}+2)\mathsf{I}\sub{L}.
  \label{eq:12}$$ Now, we consider the composite system of $L(G)$ and $G$. Then the composite Laplacian $\mathcal{L}\sub{c}$ is given by $$\mathcal{L}\sub{c}=-\mathcal{K}\sub{c} +2m\mathcal{I}\sub{c},  \label{eq:13}$$ with $$\mathcal{K}\sub{c}=
    \begin{bmatrix}
   \bar{\mathsf{X}}^\mathrm{T}\bar{\mathsf{X}}&0\\
   0&\bar{\mathsf{X}}\bar{\mathsf{X}}^\mathrm{T}
  \end{bmatrix},\
  \mathcal{I}\sub{c}=
    \begin{bmatrix}
   \mathsf{I}\sub{L}&0\\
   0&\mathsf{I}
  \end{bmatrix},
\nonumber$$ where we have used Eqs. (\[eq:8\]), (\[eq:9\]), and (\[eq:12\]). The composite operator $\mathcal{K}\sub{c}$ is written as $$\mathcal{K}\sub{c}=\mathcal{Q}\mathcal{Q}^\dagger+\mathcal{Q}^\dagger\mathcal{Q},  \label{eq:14}$$ where the symbol $\dagger$ represents the Hermitian conjugate, and we define the supercharge as $$\mathcal{Q} =
\begin{bmatrix}
 0 & 0\\
 \bar{\mathsf{X}} & 0
\end{bmatrix}.$$ These operators satisfy the superalgebra [@Cooper1994]: $$[\mathcal{K}\sub{c},\mathcal{Q}]=[\mathcal{K}\sub{c},\mathcal{Q}^\dagger]=0,$$ $$\{\mathcal{Q},\mathcal{Q}^\dagger\}=\mathcal{K}\sub{c},\ \{\mathcal{Q},\mathcal{Q}\}=\{\mathcal{Q}^\dagger,\mathcal{Q}^\dagger\}=0,$$ where $\{\mathcal{A},\mathcal{B}\}$ and $[\mathcal{A},\mathcal{B}]$ are the anticommutator and the commutator, respectively. Therefore, the eigenspectra of the inductor-capacitor circuit networks on the simple regular graph and its line graph must coincide except, possibly, for the highest eigenfrequency. Actually, if we have eigenvector $\vct{x}$ satisfying $\mathsf{L}\vct{x}=E\vct{x}$ with the eigenvalue $E$, we obtain $\bar{\mathsf{X}}^\mathrm{T}\vct{x}$, satisfying $\mathsf{L}\sub{L}(\bar{\mathsf{X}}^\mathrm{T}\vct{x})=E(\bar{\mathsf{X}}^\mathrm{T}\vct{x})$. Then, we have an eigenvector $\bar{\mathsf{X}}^\mathrm{T}\vct{x}$ for $\mathsf{L}\sub{L}$ when $\bar{\mathsf{X}}^\mathrm{T}\vct{x}\ne \vct{0}$ ($E\ne 2m$). The eigenvalue $E$ of $\mathsf{L}$ and $\mathsf{L}\sub{L}$ must satisfy $E\leq 2m$, because $\bar{\mathsf{X}}\bar{\mathsf{X}}^\mathrm{T}$ and $\bar{\mathsf{X}}^\mathrm{T}\bar{\mathsf{X}}$ are positive-semidefinite. For an eigenvector $\vct{y}_{2m}$ of $\mathsf{L}$ ($\mathsf{L}\sub{L}$) with eigenvalue $E=2m$, the partner mode cannot be obtained by multiplying $\bar{\mathsf{X}}^\mathrm{T}$ ($\bar{\mathsf{X}}$), because of $\bar{\mathsf{X}}^\mathrm{T}\vct{y}_{2m}=0$ ($\bar{\mathsf{X}}\vct{y}_{2m}=0$). The condition for complete spectral coincidence of the eigenvalues of $\mathsf{L}$ and $\mathsf{L}\sub{L}$ is discussed in Appendix \[sec:appA\]. Note that quantum tight-binding models represented by Eqs. (\[eq:8\]) and (\[eq:12\]) are isospectral except, possibly, for the highest eigenenergy, but an accurate tuning of the on-site potential satisfying Eqs. (\[eq:8\]) and (\[eq:12\]) is usually difficult to achieve. The significant point of SUSY for the inductor-capacitor circuit networks is that the on-site potential tuning is accomplished naturally.

If $G$ is a periodic graph with lattice vectors $\{\thrvct{a}_i\}$, we can show that the spectral coincidence except possibly for the highest eigenfrequency holds for [*each wave vector.* ]{} We define parallel translation with $\thrvct{a}_i$ as $\mathsf{P}^\mathcal{E}_i$ and $\mathsf{P}^\mathcal{V}_i$ for edges and vertices, respectively. From the translational symmetry, we have $\bar{\mathsf{X}}\mathsf{P}^\mathcal{E}_i=\mathsf{P}^\mathcal{V}_i\bar{\mathsf{X}}$. For a Bloch vector $\vct{x}_{\thrvct{k}}$ satisfying $\mathsf{P}^\mathcal{E}_i \vct{x}_{\thrvct{k}}=\exp({\ii\thrvct{k}\cdot\thrvct{a}_i}) \vct{x}_{\thrvct{k}}$, we have $\mathsf{P}^\mathcal{V}_i(\bar{\mathsf{X}}\vct{x}_{\thrvct{k}})=\exp(\ii\thrvct{k}\cdot\thrvct{a}_i)(\bar{\mathsf{X}}\vct{x}_{\thrvct{k}})$. This means that $\bar{\mathsf{X}}\vct{x}_{\thrvct{k}}$ is also a Bloch vector. Similar discussion can be applied to $\bar{\mathsf{X}}^\mathrm{T}$. Then, $\bar{\mathsf{X}}$ and $\bar{\mathsf{X}}^\mathrm{T}$ map Bloch vectors to Bloch vectors without changing $\thrvct{k}$. Therefore, the decomposition shown in Eq. (\[eq:13\]) is valid in the subspace of Bloch vectors with a wave vector $\thrvct{k}$. This means that the spectral coincidence except, possibly, for the highest eigenfrequency holds for each wave vector.

Examples \[subsec:c\]
---------------------

### finite case\[subsub:1\]

For the graph shown in Fig. \[fig:lc\_ladder\](a), we have $$\bar{\mathsf{X}}=\begin{bmatrix}
  0 &1 &1 &1 &0 &0\\
  1 &0 &1 &0 &1 &0\\
  1 &1 &0 &0 &0 &1\\  
  0 &0 &0 &1 &1 &1
         \end{bmatrix}\nonumber.$$ Then, we get $\omega=2\omega_0,\ 2\omega_0,\ 2\omega_0,\ 0$ for the inductor-capacitor circuit network on $G$. On the other hand, $\omega=\sqrt{6}\omega_0,\ \sqrt{6}\omega_0,\ 2\omega_0,\ 2\omega_0,\ 2\omega_0,\ 0$ are obtained for the inductor-capacitor circuit network on $L(G)$. We can see that all angular eigenfrequencies for $G$ are included in those for $L(G)$.

### Infinite case

![\[fig:kagome\_hex\] (color online) (a) Hexagonal lattice. (b) lattice. (c) Dispersion relation for a hexagonal inductor-capacitor circuit network. (d) Dispersion relation for a inductor-capacitor circuit network. (e) Eigenmodes of the higher band at the $\Gamma$ and $\mathrm{M}$ points for the hexagonal lattice. (f) Eigenmodes of the middle band at the $\Gamma$ and $\mathrm{M}$ points for the lattice. ](kagome_hex.eps){width="86mm"}

Here, we consider a hexagonal lattice as $G$ \[Fig. \[fig:kagome\_hex\](a)\]. The line graph $L(G)$ is a lattice \[Fig. \[fig:kagome\_hex\](b)\]. To see the spectra coincidence directly, we calculate the angular eigenfrequencies $\omega$. At first, we calculate the eigenvalue $\alpha$ and $\alpha\sub{L}$ for $\mathsf{A}$ and $\mathsf{A}\sub{L}$, respectively. Due to the Bloch theorem, it is enough to calculate them in the restricted space $\mathcal{W}_{\thrvct{k}_\parallel}$ of waves with wave vector $\thrvct{k}_\parallel$. For the hexagonal lattice, we have two vertices $v_p\in \mathcal{V}$ $ (p=1,2)$ in a unit cell. The vertex displaced from $v_p$, with $i\thrvct{a}_1+j\thrvct{a}_2$, is denoted by $v_p^{(i,j)}$ for $(i,j)\in \mathbb{Z}^2$, where $\thrvct{a}_1$ and $\thrvct{a}_2$ are lattice vectors for $G$. Now, we define $\vct{\Psi}_1 (\thrvct{k}_\parallel)=\sum_{(i,j)\in\mathbb{Z}^2} \exp(\ii \thrvct{k}_\parallel\cdot(i\thrvct{a}_1+j\thrvct{a}_2) ){\vct{v}_1^{(i,j)}}$ and $\vct{\Psi}_2 (\thrvct{k}_\parallel)=\sum_{(i,j)\in\mathbb{Z}^2} \exp(\ii \thrvct{k}_\parallel\cdot(i\thrvct{a}_1+j\thrvct{a}_2) ){\vct{v}_2^{(i,j)}}$, where $\{\vct{v}_{p}^{(i,j)}|(i,j)\in\mathbb{Z}^2, p=1,2\}\subset \mathcal{H}$ is a complete orthogonal basis of the Hilbert space $\mathcal{H}$. The vector subspace $\mathcal{W}_{\thrvct{k}_\parallel}$ is spanned by $\vct{\Psi}_1 (\thrvct{k}_\parallel)$ and $\vct{\Psi}_2 (\thrvct{k}_\parallel)$. The action of $\mathsf{A}$ in the restricted space $\mathcal{W}_{\thrvct{k}_\parallel}$ is represented by a $2\times 2$ matrix $\mathsf{A}({\thrvct{k}_\parallel})=[A_{ij}({\thrvct{k}_\parallel})]$, satisfying $\mathsf{A}\vct{\Psi}_i({\thrvct{k}_\parallel})= \sum_{j=1}^{2}\vct{\Psi}_j({\thrvct{k}_\parallel}) A_{ji}({\thrvct{k}_\parallel})$. Diagonalizing $\mathsf{A}({\thrvct{k}_\parallel})$ we have $$\alpha(\thrvct{k}_\parallel)=\pm\sqrt{3+2F(\thrvct{a}_1, \thrvct{a}_2;\thrvct{k}_\parallel)},  \label{eq:15}$$ with $F(\thrvct{u}_1, \thrvct{u}_2; \thrvct{k}_\parallel)=\cos (\thrvct{k}_\parallel\cdot\thrvct{u}_1)+\cos( \thrvct{k}_\parallel\cdot\thrvct{u}_2)+\cos \thrvct{k}_\parallel\cdot(\thrvct{u}_1-\thrvct{u}_2)$. By applying a similar calculation to the lattice, we obtain $$\alpha\sub{L}(\thrvct{k}_\parallel) =-2, 1\pm\sqrt{3+2F(\thrvct{a}_1\sur{L}, \thrvct{a}_2\sur{L};\thrvct{k}_\parallel)},  \label{eq:16}$$ with lattice vectors $\thrvct{a}_1\sur{L}$ and $\thrvct{a}_2\sur{L}$ for $L(G)$. Using Eqs. (\[eq:5\]), (\[eq:6\]), and (\[eq:15\]), we obtain $$\frac{\omega}{\omega_0}=\sqrt{3\pm\sqrt{3+2F(\thrvct{a}_1, \thrvct{a}_2;\thrvct{k}_\parallel)}}  \label{eq:17}$$ for the hexagonal lattice. From Eqs. (\[eq:5\]), (\[eq:10\]), and (\[eq:16\]), the lattice also has the dispersion relation $$\frac{\omega}{\omega_0}=\sqrt{6},  \sqrt{3\pm\sqrt{3+2F(\thrvct{a}_1\sur{L}, \thrvct{a}_2\sur{L};\thrvct{k}_\parallel)}}.  \label{eq:18}$$ The obtained dispersion relations are shown in Figs. \[fig:kagome\_hex\](c) and (d). The lower two bands are identical as we expected. Note that these bands are determined only by the product $LC$ and are independent of the ratio $L/C$. We also show examples of the eigenmodes for the hexagonal and lattices in Figs. \[fig:kagome\_hex\](e) and (f).

Band correspondence between metallic hexagonal and lattices \[sec:3\]
=====================================================================

In the previous section, we developed inductor-capacitor circuit networks that are related through SUSY. As an example, we saw that the bands of hexagonal and inductor-capacitor circuit networks are isospectral by SUSY, except for the highest band of the lattice. In this section, we examine this correspondence for realistic system. It is known that bar-disk resonators composed of metallic disks connected by metallic bars can be qualitatively modeled by the inductor-capacitor circuit networks discussed in Sec. \[sec:2\], because charges on the disks are coupled dominantly by the current flowing along the bars [@Nakata2012; @Kajiwara2016]. The modes on the bar-disk resonators are called spoof plasmons. Here, we study the spoof plasmons of metallic hexagonal and lattices whose designs are shown in Fig. \[fig:kagome\_hex\_lattice\_design\].

![\[fig:kagome\_hex\_lattice\_design\] Designs for (a) metallic hexagonal lattice and (b) metallic lattice. The following parameters are used: $d=10\,\U{\mu m}$, $r=150\,\U{\mu m}$, $b=800/\sqrt{3}\,\U{\mu m}$, and thickness $h=30\,\U{\mu m}$.](kagome_hex_lattice_design.eps){width="86mm"}

Simulation\[sec:3-1\]
---------------------

![\[fig:hex\_kagome\_bands\] (color online) Dispersion relations obtained by simulation for a (a) metallic hexagonal lattice and (b) metallic lattice. The fitting parameters for the theoretical models discussed in Appendix \[sec:appB\] are given by $\omega_0\sur{hex}=2\pi \times 0.107\,\U{THz}$, $\eta\sur{hex}_0=0.0916$ for the hexagonal lattice and $\omega_0\sur{kag}=2\pi\times 0.101\,\U{THz}$ and $\eta\sur{kag}_0=0.142$ for the lattice. ](hex_kagome_band.eps){width="86mm"}

![image](eigenmodes.eps){width="172mm"}

We perform an eigenfrequency analysis for the metallic hexagonal and lattices by the finite element method solver (<span style="font-variant:small-caps;">Comsol Multiphysics</span>). The parameters of the structures for Fig. \[fig:kagome\_hex\_lattice\_design\] are as follows: bar width $d=10\,\U{\mu m}$, radius of disks $r=150\,\U{\mu m}$, distance between nearest disks $b=800/\sqrt{3}\,\U{\mu m}$, and thickness $h=30\,\U{\mu m}$. In each simulation, the finite thickness metallic lattice parallel to $z=0$ is located in $z\in [-h/2,h/2]$. The unit cell in the $xy$ plane is the rhombus spanned by the lattice vectors and denoted by $U$. To reduce the degrees of freedom, we use the mirror symmetry with respect to $z=0$. A simulation domain with the material parameters of a vacuum is set in $U\times [0,6l]$ with $l=\sqrt{3} b=800\,\U{\mu m}$, and a perfect magnetic conductor condition imposed on the surface $z=0$. Half of the structure in $z\in [0,h/2]$ is engraved in the simulation domain and a perfect electric conductor (PEC) boundary condition is imposed on the structure surface. A perfect matched layer (PML) in $U\times [5l, 6l]$ with a PEC boundary at $z=6l$ is used to truncate the infinite effect. The periodic boundary condition with a phase shift (Floquet boundary conditions with a wave vector $\thrvct{k}_\parallel$) is applied to $\partial U \times [0,6l]$. Changing $\thrvct{k}_\parallel$ along the Brillouin zone boundary, we calculate the eigenfrequencies. To remove the modes which are not localized near the metallic surface [@Parisi2012], we select the modes with $$\xi=\frac{\int_{U\times [9l/2, 6l]} |\tilde{\thrvct{E}}|^2 \dd V}{\int_{U\times [h/2, l/2]} |\tilde{\thrvct{E}}|^2 \dd V}<1,\nonumber$$ where the complex amplitude of the electric field of the mode is denoted by $\tilde{\thrvct{E}}$.

The calculated eigenfrequencies for the metallic hexagonal and lattices are shown in Fig. \[fig:hex\_kagome\_bands\] as circles. Note that some points are missing because unphysical modes located near PML accidentally exist or couple with the modes. As we explained earlier, we eliminated such modes with $\xi\geq 1$. In Fig. \[fig:hex\_kagome\_bands\], we observe the lower two band correspondence between the metallic hexagonal and lattices. The bands of the metallic hexagonal lattice are about 5% higher than those of the lattice. However, we can say that the band correspondence is qualitatively established. The detailed theoretical models for fitting curves are discussed later in Appendix. \[sec:appB\]. Figure \[fig:eigenmodes\] shows the electric flux density $D_z$ on $z=h/2$ of the specific modes, where $D_z$ corresponds to the surface charge on the metal. These mode profiles agree with the theoretically calculated eigenmodes shown in Figs. \[fig:kagome\_hex\](e) and (f).

![\[fig:hdbr\_kdbr\] (color online) Microphotographs of the (a) metallic hexagonal and (b) metallic lattices.](hdbr_kdbr.eps){width="86mm"}

Experiment
----------

![image](transmission.eps){width="172mm"}

To investigate the dispersion relation experimentally, we fabricated the metallic hexagonal and lattices by etching and performed transmission measurement on them using the terahertz time-domain spectroscopy technique. The samples made of stainless steel (SUS304) are shown in Fig. \[fig:hdbr\_kdbr\]. The geometrical parameters of these samples are the same as those for the simulation model in Sec. \[sec:3-1\]. The area where structures are patterned is $4\,\U{cm}\times 4\,\U{cm}$. The terahertz beam is generated by a spiral antenna and collimated by a combination of a hyper-hemispherical silicon lens and a Tsurupica$^\text{\textregistered}$ lens. The beam diameter is set to $13\, \U{mm}$ by an aperture. Wire-grid polarizers are located near the emitter and detector, which are adjusted so that the emitted and detected fields have the same linear polarization. The transmission spectrum $T(\omega)$ in the frequency domain is obtained from $T(\omega)=|\tilde{E}(\omega)/\tilde{E}\sub{ref}(\omega)|^2$, where $\tilde{E}(\omega)$ and $\tilde{E}\sub{ref}(\omega)$ are Fourier transformed electric fields with and without the sample. To scan the Brillouin zone, power transmission spectra are measured with changing incident angle $\theta$ from $\theta=0^\circ$ to $60^\circ$ with step $2.5^\circ$. Here, the magnitude of the in-plane wave vector $\thrvct{k}_\parallel$ is given by $k_\parallel = (\omega/c) \sin \theta$, where $c$ is the speed of light.

To observe the higher band of the hexagonal lattice and the middle band of the lattice, the incident waves are set as follows: (i) transverse electric (TE) modes in $\Gamma$–$\mathrm{K}$ scan and (ii) transverse magnetic (TM) modes in $\Gamma$–$\mathrm{M}$ scan. Figures \[fig:transmission\](a) and (b-1) show the power transmission spectra for these incident waves entering into the metallic hexagonal and lattices, respectively. The calculated eigenfrequencies are shown simultaneously as circles in Fig. \[fig:transmission\]. We can see that the transmission dips form a band from $0.15$ to $0.3\,\U{THz}$. The calculated eigenfrequencies are located around the experimental transmission dips. Thus, the SUSY band correspondence for the second band is experimentally demonstrated.

The highest band for the metallic lattice can be observed for differently polarized incident waves. Figure \[fig:transmission\](b-2) shows the transmission spectra for the metallic lattice, where the incident waves are set as (i) TM modes in the $\Gamma$–$\mathrm{K}$ scan and (ii) TE modes in the $\Gamma$–$\mathrm{M}$ scan. In Fig. \[fig:transmission\](b-2), we can see the flat band reported in Ref. . Note that the frequencies of the lowest band modes are under the light line, so it is impossible to excite them by free-space plane waves. To excite the lowest band modes, another method, e.g., attenuated total reflection measurement, is needed [@Kajiwara2016].

Discussion \[sec:3.3\]
----------------------

![\[fig:tuned\_kagome\] (color online) Comparison between eigenfrequencies for a tuned metallic lattice with $d=28\,\U{\mu m}$ and metallic hexagonal lattice with $d=10\,\U{\mu m}$. The other parameters are $r=150\,\U{\mu m}$, $b=800/\sqrt{3}\,\U{\mu m}$, and thickness $h=30\,\U{\mu m}$. ](tuned_kagome.eps){width="86mm"}

In the previous subsections, the band correspondence between the metallic hexagonal and lattices was demonstrated, but a $5\%$ discrepancy between the bands was also observed. Here, we investigate the possibility to compensate empirically for the discrepancy. We calculated the eigenfrequencies for the metallic lattice with the bar width $d=28\,\U{\mu m}$. The other parameters are the same as the previous one. The calculated results are shown in Fig. \[fig:tuned\_kagome\] as circles. To compare with the previous result, the eigenfrequencies for the metallic hexagonal lattice with $d=10\,\U{\mu m}$ are also plotted in Fig. \[fig:tuned\_kagome\]. We can see the improvement of band correspondence between the metallic hexagonal lattice and tuned metallic lattice.

Conclusion \[sec:4\]
====================

In this paper, we showed that the inductor-capacitor circuit networks on a simple regular graph and its line graph are related through SUSY, and their spectra must coincide (except possibly for the highest eigenfrequency). The SUSY structure for the circuits was derived from the topological nature of the graphs. To observe SUSY correspondence of the bands in the high frequency range, we investigated the metallic hexagonal and lattices. The band correspondence between them was predicted by a simulation. We performed terahertz time-domain spectroscopy for these metallic lattices and observed the band correspondence. Finally, we proposed an empirical tuning method to reduce the discrepancy of the corresponding bands of the metallic hexagonal and lattices.

The theoretical results are formulated for the inductor-capacitor circuit networks and independent of the implementations. Therefore, our results is also applicable to transmission-line systems such as microstrip. The SUSY correspondence in the spectra of inductor-capacitor circuit networks has the potential to extend mode filters [@Miri2013] and mode converters [@Heinrich2014] in two-dimensional (spoof) plasmonic systems.

The present research was supported by the JSPS KAKENHI Grant No. 25790065 and Grant-in-Aid for JSPS Fellows No. 13J04927. Two of the authors (Y. N. and Y. U.) were supported by JSPS Research Fellowships for Young Scientists.

Condition for complete spectral coincidence of eigenvalues of $\mathsf{L}$ and $\mathsf{L}\sub{L}$\[sec:appA\]
==============================================================================================================

In this Appendix, we consider condition for complete spectral coincidence of the eigenvalues of $\mathsf{L}$ and $\mathsf{L}\sub{L}$. Here, we mainly consider the finite graph cases. At first, we introduce a bipartite graph. A bipartite graph is a graph whose vertex set can be separated into two disjoint sets $\mathcal{V}_1$ and $\mathcal{V}_2$ such that every edge is connected between a vertex in $\mathcal{V}_1$ and that in $\mathcal{V}_2$. Using the concept of a bipartite graph, we obtain the following lemma:

Consider a connected graph $G=(\mathcal{V}, \mathcal{E})$ satisfying $\#\mathcal{V}>0$ and $\#\mathcal{E}>0$. Let $\bar{\mathsf X}$ be an undirected incidence matrix of $\bar{G}$. In this case, $\rank \bar{\mathsf X}< \#\mathcal{V}$ is satisfied if and only if $G$ is bipartite.

$\Leftarrow$: If the graph is bipartite, $\mathcal{V}$ is represented by the disjoint union of $\mathcal{V}_1$ and $\mathcal{V}_2$ as $\mathcal{V}=\mathcal{V}_1 \sqcup \mathcal{V}_2$. The $v_i$-component row vector of $\bar{\mathsf X}$ is denoted by $\bar{\vct{X}}_{i}$, where $v_i \in \mathcal{V}$ $(i=1,2,\cdots, n)$, and $n=\#\mathcal{V}$. Without loss of generality, we can assume $v_i \in \mathcal{V}_1$ $(i=1,2,\cdots, l)$ and $v_i\in \mathcal{V}_2$ $(i=l+1, l+2,\cdots, n)$. Because the graph is bipartite, $\sum_{i=1}^{l}\bar{\vct{X}}_{i}=\sum_{i=l+1}^{n}\bar{\vct{X}}_{i}$ is satisfied.\
$\Rightarrow$: If we assume $\rank \bar{\mathsf X}\ne \#\mathcal{V}$, $$\sum_{i=1}^{n}c_i \bar{\vct{X}}_{i}=0 \label{eq:19}$$ is satisfied for $c_i \in \mathbb{R}$ and not every $c_i=0$. At first we show $c_i\ne 0$ for all $i\in\{1,2,\cdots,n\}$. We assume $c_{q}=0$ for $q \in \{1,2,\cdots,n\}$. For arbitrary $r \in \{1,2,\cdots,n\}$, we consider a path $(v_{f(1)}, e_{1}, v_{f(2)}, e_2, \cdots, e_{p}, v_{f(p+1)})$ from $v_{q}$ to $v_{r}$, where $p$ is the path length, $f$ is a function from $\{1,2,\cdots, p+1\}$ to $\{1,2,\cdots, n\}$ satisfying $f(1)=q$ and $f(p+1)=r$, and the edge $e_{s}$ connects vertices $v_{f(s)}$ and $v_{f(s+1)}$ $(s=1,2,\cdots, p)$. Considering $e_{s}$-column of Eq. (\[eq:19\]), we have $c_{f(s)}=-c_{f(s+1)}$. Then, we obtain $c_i=0$ for all $i$. This leads a contradiction because we assumed not every $c_i=0$. Therefore, we have $c_i\ne0$ for all $i$. We assume $c_i>0$ $(i=1,2,\cdots, l)$, and $c_i<0$ $(i=l+1,l+2,\cdots, n)$, without loss of generality. For a given arbitrary $e\in \mathcal{E}$, we consider the column component about $e$ of Eq. (\[eq:19\]). Then, we find $e$ is connected between a vertex in $\mathcal{V}_1=\{v_1,v_2,\cdots,v_l\}$ and that in $\mathcal{V}_2=\{v_{l+1},v_{l+2},\cdots,v_{n}\}$. This shows $G$ is bipartite. This proof is based on Ref. .

From this lemma, we can prove the following theorem:

Consider a simple $m$-regular connected graph $G=(\mathcal{V}, \mathcal{E})$ with $\#\mathcal{E}>0$. There exists at least one mode with $\omega/\omega_0=\sqrt{2m}$ for the inductor-capacitor circuit network on $G$ if and only if $G$ is bipartite.

Let $\bar{\mathsf X}$ be an undirected incidence matrix of $\bar{G}$. $G$ is bipartite $\Leftrightarrow$ $\rank{\bar{\mathsf X}}<\#\mathcal{V}$ $\Leftrightarrow$ $\dim \ker{\bar{\mathsf{X}}^\mathrm{T}}=\#\mathcal{V}-\rank{\bar{\mathsf X}}>0$ $\Leftrightarrow$ $^\exists \vct{x}$ satisfying $\bar{\mathsf{X}}^\mathrm{T}\vct{x}=0$ $\Leftrightarrow$ $^\exists \vct{x}$ satisfying $\mathsf{L}\vct{x}=2m\vct{x}$.

On the other hand, we can formulate the condition for presence of eigenvalue $E=2m$ of $\mathsf{L}\sub{L}$ as follows:

Consider a simple $m$-regular graph $G=(\mathcal{V}, \mathcal{E})$ with $\#\mathcal{V}>0$. There exists at least one mode with $\omega/\omega_0=\sqrt{2m}$ for the inductor-capacitor circuit network on $L(G)$ if $m> 2$.

Let $\bar{\mathsf X}$ be an undirected incidence matrix of $\bar{G}$. We have $\dim \ker{\bar{\mathsf X}}=\#\mathcal{E}-\rank{\bar{\mathsf{X}}}\geq \#\mathcal{E}-\#\mathcal{V}=\frac{m-2}{2}\#\mathcal{V}>0$ for $m>2$. Then, $^\exists \vct{x}$ satisfying $\bar{\mathsf{X}}\vct{x}=0$. Finally, $^\exists \vct{x}$ satisfying $\mathsf{L}\sub{L}\vct{x}=2m\vct{x}$.

From these theorems, all spectra of the inductor-capacitor circuit networks on a simple $m$-regular ($m>2$) connected bipartite graph and its line graph completely coincide. On the other hand, we find that $\mathcal{L}\sub{c}$ have a non-degenerated eigenvector with eigenvalue $E=2m$ for a simple $m$-regular ($m>2$) connected non-bipartite graph $G$ (e.g. Sec. \[subsub:1\]). For $m=2$, all spectra of the inductor-capacitor circuit networks on a simple connected $m$-regular graph $G$ and its line graph $L(G)$ coincide because of $G=L(G)$. In this case, there are $\omega/\omega_0=\sqrt{2m}$ modes if and only if $\#\mathcal{V}$ is even.

To analyze an infinite periodic graph with lattice vectors $\{\thrvct{a}_i\}$, we take a supercell spanned by $\{N_i \thrvct{a}_i\}$, $N_i>1$. We impose a periodic boundary condition (called Born–von Karman boundary condition) on the sides of the supercell. Note that this boundary condition just leads to discretization of wave vectors in the Brillouin zone. Now, the infinite graph is reduced to a finite graph and we can use the theorems. For example, the hexagonal lattice is a simple 3-regular connected bipartite graph. Then, the spectra of the inductor-capacitor circuit networks on hexagonal and lattices include $\omega/\omega_0=\sqrt{6}$, simultaneously. Note that the spectra of the inductor-capacitor circuit networks on hexagonal and lattices completely coincide, but their dispersion relations do not.

Detailed theoretical model for the metallic hexagonal and lattices \[sec:appB\]
===============================================================================

In this appendix, we derive detailed theoretical models for fitting the eigenfrequencies of the metallic hexagonal and lattices. For the metallic hexagonal and lattices, the circuit models treated in Sec. \[sec:2\] are approximately valid. To improve the model accuracy, we have to take into account capacitive couplings between disks. Considering only nearest neighboring couplings, we modify Eq. (\[eq:4\]) as $\mathsf{P}=C^{-1}( \mathsf{I} + \eta \mathsf{A})$. Here, $C^{-1} \eta \mathsf{A}$ represents the capacitive coupling between the adjacent disks [@Kajiwara2016]. Then, we obtain $$\frac{\omega}{\omega_0}=\sqrt{ \big[m-\alpha(\thrvct{k}_\parallel)\big]\big[1+\eta \alpha(\thrvct{k}_\parallel)\big] }$$ for the metallic hexagonal lattice and $$\frac{\omega}{\omega_0}=\sqrt{ \big[m\sub{L}-\alpha\sub{L}(\thrvct{k}_\parallel)\big]\big[1+\eta \alpha\sub{L}(\thrvct{k}_\parallel)\big] }$$ for the metallic lattice. Generally, the coupling constant $\eta$ depends on $\omega$ as $\eta=\eta_0 \exp[\ii (\omega/c) b]$, where $b$ is the distance between the nearest disks [@Yeung2011]. The imaginary part of $\eta$ represents the resistive component. If we only focus on the real part of the eigenfrequencies, we may ignore the resistive term (note that we have already ignored the imaginary part of $L$ and $C$). Then, we assume $\eta=\eta_0 \cos[(\omega/c) b]$. Using the real part of the eigenvalues calculated by simulation, we numerically minimize the error $$\mathrm{err}\sur{hex}(\omega_0, \eta_0)=\sum_{(\omega,\thrvct{k}_\parallel) \in \text{\{data points\}} }g\big(\omega_0, \eta_0;\omega, m, \alpha(\thrvct{k}_\parallel)\big)^2$$ for the hexagonal lattice and $$\mathrm{err}\sur{kag}(\omega_0, \eta_0)=\sum_{(\omega,\thrvct{k}_\parallel) \in \text{\{data points\}} }g\big(\omega_0, \eta_0;\omega, m\sub{L}, \alpha\sub{L}(\thrvct{k}_\parallel)\big)^2$$ for the lattice, where we define $$g(\omega_0, \eta_0;\omega, m, \alpha)=\omega- \omega_0\sqrt{\left[m-\alpha\right]\left[1+\alpha \eta_0 \cos\left(\frac{\omega b}{c}\right)\right]}.$$ The obtained fitting parameters are as follows: $\omega_0=\omega_0\sur{hex}=2\pi \times 0.107\,\U{THz}$, $\eta_0=\eta\sur{hex}_0=0.0916$ for the hexagonal lattice, and $\omega_0=\omega_0\sur{kag}=2\pi\times 0.101\,\U{THz}$ and $\eta_0=\eta\sur{kag}_0=0.142$ for the lattice. Because the magnetic coupling and higher order effects (beyond the nearest capacitive coupling) are included in these parameters, the parameters for the hexagonal and lattice can be different. The dispersion curves with these fitting parameters are shown in Fig. \[fig:hex\_kagome\_bands\]. These curves agree with the simulated data despite the simplicity of the model.

[35]{}ifxundefined \[1\][ ifx[\#1]{} ]{}ifnum \[1\][ \#1firstoftwo secondoftwo ]{}ifx \[1\][ \#1firstoftwo secondoftwo ]{}““\#1””@noop \[0\][secondoftwo]{}sanitize@url \[0\][‘\
12‘\$12 ‘&12‘\#12‘12‘\_12‘%12]{}@startlink\[1\]@endlink\[0\]@bib@innerbibempty [****,  ()](\doibase 10.1016/0370-1573(94)00080-M) [****,  ()](\doibase 10.1103/PhysRevA.82.032111) [****, ()](\doibase 10.1016/j.aop.2013.03.011) [****,  ()](\doibase 10.1209/0295-5075/102/40008) [****,  ()](\doibase 10.1038/srep05274) [****,  ()](\doibase
  10.1038/ncomms9269) [****,  ()](\doibase 10.1103/PhysRevLett.52.1933) [****,  ()](\doibase 10.1016/0375-9601(94)00616-4) [****, ()](\doibase 10.1103/PhysRevLett.110.233902) [****,  ()](\doibase 10.1103/PhysRevA.89.033826) [****, ()](\doibase 10.1088/2040-8978/17/4/045803) [****,  ()](\doibase
  10.1364/OE.22.000987) [****, ()](\doibase 10.1364/OPTICA.1.000089) [****, ()](\doibase 10.1364/OL.40.000463) [****,  ()](\doibase 10.1103/PhysRevA.87.043819) [****,  ()](\doibase
  10.1038/srep08568) [****, ()](\doibase 10.1364/OL.39.001697) [****, ()](\doibase 10.1103/PhysRevA.89.052132) [****,  ()](\doibase 10.1103/PhysRevA.92.023839) [****, ()](\doibase 10.1103/PhysRevA.89.032116) [****,  ()](\doibase
  10.1103/PhysRevLett.106.193903) [****,  ()](\doibase
  10.1038/ncomms4698) [****,  ()](\doibase
  10.1364/OL.39.006130) @noop [**]{} (, , ) @noop  () [****,  ()](\doibase 10.1002/pssb.201248154) [****,  ()](\doibase 10.1126/science.1133268) [****,  ()](\doibase
  10.1016/j.photonics.2011.12.002) @noop  () [****,  ()](\doibase 10.1090/S0002-9947-99-02497-6) [****,  ()](\doibase 10.1103/PhysRevB.85.205128) @noop [****,  ()]{} [****,  ()](\doibase 10.1364/OE.20.016690) [****,  ()](\doibase 10.1109/TCS.1976.1084251) [****,  ()](\doibase 10.1109/TMTT.2011.2163803)
