1. Introduction

Plasmonic metasurfaces (PMs) are a preferred avenue as sensors for the label-free detection of biological or chemical analytes [1–5]. PMs are made of two-dimensional arrays of plasmonic scatterers and support global resonances that arise from the coupling of local resonances, which occur on each scatterer due to the presence of localized surface plasmons. In recent years, this governing principle has been used to design several metasensors [5–8].

Both the sensing detection limit and accuracy are greatly increased by using narrow width spectral features [9]. Most often, it is obtained by combining a super-radiating resonator, bearing an electric dipolar momentum and acting as a radiative or bright mode, with a sub-radiating resonator bearing a quadrupolar electric or magnetic dipolar momentum and acting as a trapped or dark mode (DM) [9–12].

Entering the territory of coupled resonators with generic radiative coupling, several concepts are of order: (i) Associated with dark mode is the concept of electromagnetically induced transparency (EIT), which grants narrow-band resonances in the middle of broad absorption peaks; (ii) When coupling and losses are involved in two resonators, exceptional points (EP) may arise [13–20], as has been popularized in Parity-Time (PT) symmetric system [21–23]. Here the “passive-PT” depiction applies, whereby a loss difference in the two resonators combined with adequate coupling strength gives rise to eigenmodes that either break or preserve the symmetry, and coalesce at the very EP.

In this work, we discuss the role of EPs but only briefly address their role in EIT in Supplementary Material (SM), for the sake of completeness. Rather, we mostly retain the role of EPs as structuring points of the parameter space that we shall depict, namely the coupling–detuning plane, called CDP in short. Based on a temporal coupled-mode theory (t-CMT), we show that this space provides a general (“universal”) map that greatly helps navigating the resulting Fano-resonance configurations (control of spectral response asymmetry and its relation to DM and EP).

This result provides a global vision on the problem of sensing using metasurfaces based on bright and dark mode resonators coupling, and thus makes it possible to synthesize the strategies of a large body of the experimental work carried out in this field. The universality of the scattering matrix solutions evidenced in our work provides a general recipe for an unlimited diversity of optimal designs, answering bespoke sensors demands with tailored spectral responses, especially in the wavelength range from microwave to 10 THz.

Special attention is given to explicitly take into account the open nature of a system due to the antenna radiation properties of plasmonic resonators. In particular, it is pointed out that the loss engineering is an important factor in the control of the contrast of Fano resonance.

2. Temporal coupled mode theory model

The ingredients of the t-CMT developed below are (i) a set of two resonators that are coupled radiatively to outer space (ii) a set of coefficients (matrices), partly related, that describes the coupling of waves to resonators and also the non resonant background (essential to Fano-type operation). The system of coupled resonators is schematically represented in Fig. 1.

 

Fig. 1. Sketch of single-layer metasurface. The coupling to a specific incoming and outgoing waves (scattering matrix) are not represented.

Download Full Size | PDF

Based on [24–28], a set of coupled ordinary differential equations then embodies the t-CMT as follows:

Here, lowercase letters designate vectors whereas uppercase bold letters designate matrices. a describes resonance amplitude, ${{\text s}_ + }$ incoming waves, and ${{\text s}_ - }$ outgoing waves. In the parenthesis of Eq. (1) right hand side, there are three matrices detailed later: Ω describes the Hermitian part of the resonator dynamics (not coupled to outside), ${{\text X}}_\textrm{A}$ corresponds to partial decay rates due to absorption, ${{ \Gamma }_\textrm{R}}$ corresponds to partial decay rates due to radiation, and ${{\text K}_ \pm }$ denotes the coupling between resonances and incoming or outgoing waves. ${\text C}$ describes the non-resonant contribution to the scattering matrix.

Due to reciprocity and time-reversal symmetry of the radiative coupling [24–29], these matrices obey additional relations and constraints as follows:

It is convenient to express the general result through the effective Hamiltonian:

In the frequency domain, the scattering matrix that relates incoming and outgoing waves is then simply given by:

For the effective non-Hermitian Hamiltonian H_eff, we consider a generic model of two coupled resonators that lie on a single metasurface layer:

Here ω_j is the resonance frequency of a single and uncoupled resonator, κ is the near-field coupling coefficient between both resonators, γ₁ and γ₂ are the partial decay rates due to radiation towards port 1 and 2. The off-diagonal elements $\sqrt {{\mathrm{\gamma }_\textrm{1}}{\mathrm{\gamma }_\textrm{2}}} $ describe the indirect coupling, or in other terms the radiative coupling, between the resonances induced by the interaction between each resonator and the free space. This is an important constraint that is key to several simplifications. We assume no absorption for the moment (X_A= 0).

Let us note for convenience δ=(ω₂-ω₁)/2, ${\bar{\mathrm {\gamma} }\ =\ }{{({{\mathrm{\gamma }_\textrm{1}}\textrm{ + }{\mathrm{\gamma }_\textrm{2}}} )} / \textrm{2}}$, η=(γ₂-γ₁)/2 and ${\mathrm{\gamma }_\textrm{0}}\textrm{ = }\sqrt {{\mathrm{\gamma }_\textrm{1}}{\mathrm{\gamma }_\textrm{2}}} $. They respectively represent the detuning, the (arithmetic) average radiative loss, the radiative loss contrast and the “far-field” radiative coupling (the geometric average loss).

The eigenmodes of H_eff are poles of the transmission matrix, the inner term of Eq. (7), and correspond to sharp features in the spectral response. They can be classified as either symmetric (even, ω₊), or antisymmetric (odd, ω_-). As follows directly from Eqs. (1,8–10) and using above notations we get:

3. Engineering of Fano resonances

Let us have some first insight on the control of the Fano resonance asymmetry. We shall demonstrate that a symmetric spectral response of a Fano resonance corresponds to the condition of operation along the PT-symmetry phase transition path. In turn, the largest Fano resonance asymmetry is obtained when operating away from the PT-symmetry path and close to the DM condition also identified as bound state in continuum (BIC) [30–34].

We see from Eq. (11) that the loss difference η is combined with the detuning δ and that the “radiative coupling” γ₀ is combined with the near-field coupling κ. Furthermore, it is generically possible for the two squared complex numbers in the square root to cancel out (say, δ+jη and κ+jγ₀ just differ by π/2 in argument, while having same module), thus causing eigenvalue degeneracy that corresponds to the exceptional point (EP). The real and imaginary parts of ω₊ and ω_- eigenvalues coalesce at the EP when the following conditions are simultaneously fulfilled:

For didactic purposes, we formally extend these solutions to the negative values of the coupling coefficient κ. The crossing of the EP corresponds to the PT-symmetry phase transition when the operating point of the system in the two-dimensional κ-δ space (the CDP) evolves along the path delimited by the relation:

As depicted in Fig. 2(a) it corresponds to a dotted line passing through the point of origin and forming an angle φ=-atan(γ_0/η) with the CDP abscissa axis. Substituting Eq. (13) in Eq. (11) the characteristic expression for PT-symmetry phase transition is easily obtained:

 

Fig. 2. Logarithmic-scale colormap representation of the imaginary part of the antisymmetric mode eigenfrequency Im(ω_-) in the κ-δ CDP. a) Parameters of the coupled resonators system used in numerical modeling: ω₁= 10, γ₁= 0.8, γ₂= 0.2. All values are given in normalized frequency units. The contour curves corresponding to Im(ω_-) = 0.001; 0.01; 0.1 are plotted as thin dashed lines. The dotted line at angle φ corresponds to the evolution along PT-symmetric path. The two black bold dots indicate exceptional points EP _±. The three bold white dots correspond to CDP points whose spectral responses are shown in Fig. 4; b) Parameters of the coupled resonators system used in numerical modeling: ω₁= 10, γ₁= 0.5, γ₂= 0.5 in normalized frequency units. The two black black bold dots indicate exceptional points EP _±. The three bold white dots correspond to the points in CDP whose spectral responses are shown in Fig. 3.

Download Full Size | PDF

By introducing a generalized “radial” evolution parameter ρ such that ρ cos(φ)= κ corresponding to the distance between the CDP origin and a given point on the PT-symmetry phase transition path, Eq. (14) comes into the more convenient form:

Obviously, in the CDP, solutions in the range |ρ|<|$\bar{\mathrm{\gamma}}$| correspond to the broken PT-symmetric phase, while those where |ρ|>|$\bar{\mathrm{\gamma}}$| correspond to the PT-symmetric phase. The coordinates {κ,δ} of exceptional points EP_- and EP₊ corresponding to |ρ|=|$\bar{\mathrm{\gamma}}$| are then {η,-γ₀} and {-η,γ₀}, respectively. The colormap of real and imaginary components of eigenfrequencies as well as the evolution of these components along different paths in the CDP are represented in Figs. S1 and S2 of SM, respectively.

It is easy to observe that the distance $\bar{\mathrm{\gamma}}$ between the origin and the EP does not depend on the ratio between γ₁ and γ₂ as long as their sum, i.e. the trace of the matrix Γ_R, is constant. As shown in Appendix A, the variation of the ratio γ₁/γ₂ corresponds indeed to the CDP rotation around origin. Since the distance between arbitrary eigenvalues remains constant under the rotation operation, this characteristic reflects the isometric properties of the Hamiltonian solutions given by Eq. (11). On the other hand, a variation of the trace of the matrix Γ_R is equivalent to a uniform contraction or expansion in the CDP. In other words, the solution of Eq. (11) are homothetic with respect to the variation $\bar{\mathrm{\gamma}}$. Therefore, the eigenfrequency solutions obtained for a particular case of Γ_R matrix are general and valid for arbitrary values of γ₁ and γ₂ by means of elementary appropriate transformations of the CDP by linear extension and rotation.

The universality of the solutions obtained greatly simplifies the problem, as it can be reduced to the analysis of the system behavior along the PT-symmetry phase transition path on the base of a single selected example. To calculate the transmission and reflection coefficients of the PMs we use the t-CMT formalism [24–28]. The scattering matrix is then:

Note that this particular case corresponds to the situation where the angle φ=π/2, i.e. the PT-symmetry phase-transition path coincides with the δ-axis in CDP, as represented in Fig. 2(b). Taking into account that ω₂= ω₁+ 2δ and considering the resonance frequency of the first resonator fixed, the evolution of the spectral response as a function of the detuning of the resonance frequency of the second resonator is illustrated in Fig. 3. Since ρ=δ, it corresponds to operation at different points along the PT-symmetry phase transition path. For the selected examples the operation points corresponds to: (a) and (b) broken PT-symmetric phase [(a) ρ=0.01 $\bar{\mathrm{\gamma}}$, (b) ρ=0.5$\bar{\mathrm{\gamma}}$]; (c) EP (ρ=$\bar{\mathrm{\gamma}}$); (d) PT-symmetric phase (ρ=1.5$\bar{\mathrm{\gamma}}$).

 

Fig. 3. Scattering matrix coefficients spectral response at different points along PT-symmetric path. a) ρ=0.01$\bar{\mathrm{\gamma}}$ − broken PT and also near vicinity of the dark mode; b) ρ=0.5$\bar{\mathrm{\gamma}}$ − broken PT; c) ρ=$\bar{\mathrm{\gamma}}$ − EP; d) ρ=1.5$\bar{\mathrm{\gamma}}$ − PT-symmetric. Parameters of the coupled resonators system used in numerical modeling: ω₁= 10, γ₁= γ₂=$\bar{\mathrm{\gamma}}$=0.5. All values are given in normalized frequency units.

Download Full Size | PDF

As can be seen, Fano resonance as well as the accompanying EIT effect are observed in all cases, the spectral responses being qualitatively similar. Their particularity is to be perfectly symmetric with respect to the frequency of the EIT maximum. According to Eq. (17) the maximum of EIT, i.e. S₁₁= 0, occurs when ${\mathrm{\omega}}=\hat{{\mathrm{\omega}}}=\left({\mathrm{\omega}}_1+{\mathrm{\omega}}_2\right) / 2$ while the minima, i.e. S₁₂= 0, occur at ${\mathrm{\omega}}=\breve{{\mathrm{\omega}}}_{-}={\mathrm{\omega}}_1$ and ${\mathrm{\omega}}=\breve{{\mathrm{\omega}}}_{+}={\mathrm{\omega}}_2={\mathrm{\omega}}_1+2 {\mathrm{\delta}}$ (using a hat ^ for maxima and an anti-hat ˘ for minima). It is evident that there is no abrupt change that one would naively expect in the vicinity of the EP. The half-width δ of the EIT resonance varies linearly with path evolution parameter ρ, and becomes infinitesimally narrow in the close vicinity of the origin point. This reflects the fact that, as shown in Fig. 2(b), at the origin point one obeys Im(ω_-) = 0, i.e. the DM (BIC) condition. The presence of the transition from broken to PT-symmetric path can be however evidenced in absorption spectra, as detailed in Fig. S6 of SM.

When the operating point evolves on a trajectory different from the PT-symmetric path, then the Fano resonance becomes more classically asymmetrical, especially when located in the near vicinity of the BIC path condition given by the expression [34]:

It corresponds to the blue diagonal wedges in Fig. 2(a) forming an angle ψ = π/2 - φ with the CDP ordinate axis (additional discussion on BIC behavior is presented in Sec. 2 of SM). To illustrate this, in virtue of the isometric properties of the eigenfrequencies of the effective Hamiltonian given by Eq. (3), we consider the situation illustrated in Fig. 4, where the operating point evolves from EP to BIC along the CDP path κ=const for the system considered in Fig. 2(a).

 

Fig. 4. Scattering matrix coefficients spectral response at different points along the κ=0.3 path. a) δ=0.4 - EP; b) δ=0.1; c) δ=−0.2 close to BIC ; d) δ=−0.5. Parameters of the coupled resonators system used in numerical modeling: ω₁= 10, γ₁= 0.8, γ₂= 0.2. All values are given in normalized frequency units.

Download Full Size | PDF

As can be seen, at EP, the spectral response shown in Fig. 4(a) is symmetric despite the strong asymmetry of the radiative losses, as in the case shown in Fig. 3 where γ₁= γ₂. This result goes beyond current wisdom on this topic and demonstrates the universality of eigenfrequency solutions due to isometric and homothetic properties of scattering matrix. As one moves away from PT-symmetric phase transition path Fano resonance becomes strongly asymmetrical [Fig. 4(b)] and extremely sharp in close proximity to BIC [Fig. 4(c)]. Beyond the BIC condition, the sharpness of the Fano edge is gradually reduced [Fig. 4(d)]. The PT-symmetric phase transition path and BIC path corresponding to two radically different types of Fano resonance behavior are therefore natural orthogonal axis in the CDP for scattering matrix eigenvalues.

At this stage, we could wonder what happens if we restore nonzero absorption losses, e.g. χ=χ₁= χ₂> 0. Intuitively, the balance that governs the situation is the one between Im(ω_-) and χ. This means that a substantial impact occurs only in a “ribbon” around the blue wedge of Fig. 2(a). Logically, no line sharper than the value dictated by χ can be observed. However, having in mind the generality of Fig. 2 versus rotations (but no more versus homothetical changes), let us stress that the sharpness depends on where the BIC wedge (even smeared out by absorption losses) is intersected. Hence, it is convenient to use as a reference the ψ=0 case corresponding to Fig. 2(b), because in this case, varying the detuning does correspond to the most direct crossing of the DM wedge (vertical trajectory in CDP) and EIT type spectral response. We report in Appendix B many analytical expressions that can be developed in this case.

In terms of the resonance itself, the contrast is the most affected property. For instance, it vanishes for the exact BIC condition (a sharp line becomes nil). More elaborate modeling results taking into account the presence of absorption losses are presented in Fig. S6 of SM. In summary, the EIT resonance contrast in the presence of absorption losses gradually decreases when approaching the BIC condition.

We can introduce the ratio Q=$\bar{\mathrm{\gamma}}$/ρ as the measure of the spectral sharpness of the EIT resonance. Within the limits of low absorption losses (χ<<$\bar{\mathrm{\gamma}}$), the following approximate expressions can be used:

As it follows from Eq. (19), to get an observable Fano resonance of width ρ, absorption losses as low as 1/Q² are needed. This explains the inherent difficulty in achieving highly resonant behavior using plasmonic resonators in the optical and even in the THz domain, where metal-related ohmic losses are sizable. More detailed discussion on the impact of absorption losses is presented in Sec. 3 of SM.

4. Comparison of t-CMT model predictions with experimental results

To confirm the validity of our approach, we compared the results of the t-CMT model with those reporting an ultra-high-Q Fano resonance for a metasurface with a unit cell composed of an array of asymmetrically-split ring resonators (ASRs) [8], each cut in two pieces, and operating in the THz domain, schematically shown as insets in Figs. 5 and 6. Asymmetry in the ASRs controlling the detuning δ of the two resonators (the half-rings) is introduced by displacing by distance d the lower gap gradually from the central vertical axis. To carry out the comparison, the graphs represented in Fig. 2 of Ref. 8 were digitized, for accurate enough extraction of the experimental data as well as the full 3D numerical simulations data.

 

Fig. 5. Transmission amplitude for experimental and 3D numerical simulations data extracted from Ref. 8 as well as t-CMT model as function of asymmetry parameter variation. a) d = 1µm; b) d = 3µm; c) d = 5µm; d) d = 10µm. The inset shows the direction of electric field polarization.

Download Full Size | PDF

 

Fig. 6. Transmission amplitude for experimental and 3D numerical simulations data extracted from Ref. 8 as well as t-CMT model as function of asymmetry parameter variation. a) d = 1µm; b) d = 3µm; c) d = 5µm; d) d = 10µm. The inset show the direction of electric field polarization.

Download Full Size | PDF

When referring to our analytical model, the case of the polarization of the electric field along the split arms, shown as insets in Fig. 5, corresponds approximately to an operation along the PT-symmetry phase transition path. The spectral response reported in Ref. 8 and reproduced in Fig. 5 indeed presents an essentially symmetric EIT resonance. The spectral bandwidth and the contrast of the EIT resonance both increase with the growth of detuning parameter d, which is in full agreement with the prediction of our model [Fig. S7 of SM]. The t-CMT model was then used to fit the experimental data extracted from Ref. 8. The agreement between t-CMT modeling and experimental results is similar to that given by full 3D numerical simulations. Note that we deliberately do not introduce any correction factors liable to improve the quality of fit in our model, such as the frequency dependence of absorption losses. Our goal is to demonstrate that even a basic model with fixed parameters does capture the essential physics governing the behavior of two coupled resonators. The fitted parameters of coupled resonators system used in t-CMT modeling are provided in Table 1.

Table 1. Parameters of coupled resonators system used in numerical modeling shown in Fig. 5

View Table | View all tables in this article

The other interesting situation for our study occurs when the polarization of the electric field, shown as inset in Fig. 6, is perpendicular the split arms. It corresponds to operation in the vicinity of the DM condition with asymmetric cliff-like Fano resonance. The variation of the detuning parameter d is then equivalent to the tuning of the decay-rate for the anti-symmetric dark mode appearing from the asymmetry of the structure, as detailed in Ref. 35. The t-CMT fit results presented in Fig. 6 are in good agreement with those of experimental and numerical modeling in Ref. 8, thus confirming the validity of our analytical approach.

The parameters of coupled resonators system used in numerical modeling are provided in Table 2. Note how different they are from those of Table 1, reflecting the distinct ways the half-rings eigenmodes are different for the corresponding mode symmetry.

Table 2. Parameters of coupled resonators system used in numerical modeling shown in Fig. 6

View Table | View all tables in this article

With respect to the variety of more specific approaches reported in the literature, ours is a synthetic view on the engineering of the Fano resonance spectral shape that was lacking to the best of our knowledge. Empirically, a majority of the work targeting DM regime has focused on the intentional introduction of asymmetry in the radiative-decay using subradiant and superradiant resonators [4,6,7,10,12,29] while a minority proceeded through detuning their resonant frequencies [8,9,11], but the fact that there is a unifying common picture had not emerged from these works.

5. Summary and conclusions

The aim of the article is to provide guidelines for the engineering of Fano resonances hinging on BIC and EP conditions. It has been shown that the eigenfrequency solutions obtained for a particular case are general and valid for arbitrary values of scattering matrix by means of an appropriate transformation of the κ-δ coordinates by linear expansion and rotation. The power of our analytical approach lies in its ability to provide, on the basis of a mathematical formalism, the prediction of the spectral behavior of the system without going through more complex and time-consuming complete 3D numerical simulations which are still necessary, but rather at next stages, when translating the mathematical model into a practical design.

Specifically, we demonstrated the control of asymmetry and sharpness of the Fano resonance obtained through navigation between BIC and exceptional point corresponding to a PT-symmetric phase transition. In particular, we demonstrate the possibility of achieving a fully symmetric Fano resonance through appropriate engineering of near-field coupling and resonator detuning despite the use of bright and dark mode resonators with significantly different radiative rates. This non-trivial result, unanticipated in previously reported studies, illustrates the generality of the scattering matrix solutions highlighted in our study. To prove the validity of our approach, we compared it with experimental and full 3D numerical simulations results published in the literature and found them to be in good agreement with those of our analytical model. Furthermore, the same formalism can also be applied to describe the behavior of complementary metasurfaces by simply setting for the substrate r = 1 and t = 0.

Special care in the analysis was given to explicitly take into account the open nature of a system due to the antenna radiation properties of plasmonic resonators. It has been shown that radiative loss engineering is an important factor that determines the spectral sharpness (or quality factor Q) of the EIT resonance. In order to obtain even a very modest Q≈10, the level of absorption losses must be two orders of magnitude lower than the radiative losses. This explains the inherent difficulty in achieving highly resonant behavior using plasmonic resonators in the optical and even in the upper THz domain, where metal-related ohmic losses are important.

These results provides a global vision on the problem of sensing using metasurfaces based on bright and dark mode resonators coupling, and thus makes it possible to synthesize the strategies of a large body of the experimental work carried out in this field.

The physics and design methods described here have a broad validity that can be applied to a variety of systems, appreciably differing from metasurfaces. We notably think of systems composed of resonant nanostructures of various types interacting with waveguides (e.g. biosensors, plasmonic components, optical antennas…).