Surface plasmon polariton resonant sensors [1] have found wide-ranging applications, particularly for nanoscale biosensing [2], where they allow label-free optical detection of binding events between molecules (e.g.,  antibodies and antigens [3]), protein interactions [4] or exosomes [5]. The original prism-based geometries such as the Kretschmann [6] and Otto [7] configurations are versatile and precise, but rely on free-space optics and are thus relatively bulky. Surface plasmons lend themselves to extreme confinement of light, which could be implemented in very small footprint devices, motivating extensive research investigating strategies for integrating surface plasmon sensors with chip-based [8,9] and fiber-based [10,11] circuitry.

Experimental and theoretical work on waveguide-based plasmonic sensors goes through several decades [12,13]: we refer the reader to Refs. [14–16] for a selection of recent reviews. Such devices are commonly composed of rectangular [9,17], cylindrical [18–21], or microstructured [22,23] dielectric waveguides (e.g.,  composed of silica [21], PMMA [22], or silicon nitride [9]), adjacent to one or several metallic nanostructures (e.g.,  nanofilms [21,23], nanowires [24], and nanoantennas [9]), which are in contact with a region to be sensed (e.g.,  a liquid [21] or a gas [25,26]). Because dielectric waveguides and plasmonic films typically differ by orders of magnitude in lateral dimensions, it is generally challenging to couple light efficiently between them. One approach is to tailor the geometry of each waveguide such that individual uncoupled propagation constants are equal at a particular wavelength [27,28], i.e., they are phase matched, analogous to what occurs in sensors that rely on dielectrics alone [29–32]. Figure 1 (left) shows a concept schematic of an example plasmonic waveguide sensor. After coupling light into the dielectric core, a wavelength-dependent excitation/propagation of the modes in the sensing region occurs, leading to characteristic transmission spectra (Fig. 1, right). In typical sensing schemes, the phase-matching (PM) wavelength is associated with a local transmission dip due to directional coupling; because this condition is sensitive to the refractive index, shifts in the transmitted spectrum thus contain information on changes in the analyte. Compared to all-dielectric sensors, plasmonic sensors can exploit plasmonic modes with extremely short evanescent tails, which allow for exquisite sensitivity to small refractive index changes within a few tens of nanometers [2]. However, because plasmonic systems are inherently lossy (i.e., non-Hermitian [33,34]), they feature subtle and counterintuitive eigenmode (EM) topologies [35], whose excitation and propagation are far from trivial [36], and thus demand careful consideration. Figure 1 (right) shows different achievable transmission spectra that can result by varying wavelength and interaction length. Depending on the analyte refractive index, such sensors can behave as “effective” parity-time ($\mathcal{P}\mathcal{T}$) symmetric (EPTS), or effective $\mathcal{P}\mathcal{T}$ broken (EPTB) systems [18,37], whose resonant behavior is markedly different when the analyte length $L$ is changed [18]. Furthermore, non-Hermitian systems can support exceptional points (EPs) [35], which in some experiments can be used for enhanced sensing [38]; what role EPs play in the specific context of plasmonic waveguide sensors has, to the best of our knowledge, yet to be discussed in detail. Most commonly, experimental reports of plasmonic waveguide sensors are accompanied by mode simulations to explain the overall measured features by judiciously selecting the mode that dominates the loss spectrum [20,23]; conversely, several numerical reports are not supported by detailed comparisons with experiments [39,40] (often because the designed devices, although realistic, are challenging to fabricate). In such cases, sensing performance—e.g.,  how much a plasmonic resonance shifts as the analyte index changes—is deduced from the properties of individual 2D mode calculations, rather than considering the 3D excitation, propagation, and interference of all participating modes. Simulations of plasmonic sensor implementations [41] using commercially available solvers (which typically rely on finite element or finite difference time domain methods) often require fine meshing and large devices, making them computationally demanding and time-consuming. While full device transmission spectra have been discussed for selected configurations [18,41,42], a comprehensive study of how the inferred sensing properties fare against experimentally measurable quantities as a function of the key parameters is still missing.

Here we show, with 1D modes undergoing 2D propagation, that full transmission characteristics of non-Hermitian plasmonic waveguide sensors can be reproduced with a straightforward EM model [28] that relies on modal calculations. This model crucially takes into account both coupling in and out of the device, as well as propagation over a finite length. Our results are validated by full vector finite-element method (FEM) calculations (COMSOL). These fast computations enable us to calculate the full transmission characteristics, such as the extinction ratio, resonance width, and sensitivity, as a function of all key parameters, including wavelength, interaction length, and analyte index. This approach allows a direct comparison of full device performance (i.e., which considers mode excitation and propagation) with commonly used approaches that use modal calculations alone (i.e., which do not consider excitation and propagation). We show that many proxies extracted from mode calculations or coupled mode theory alone, such as PM wavelength, are not representative of a device’s performance and can provide misleading sensitivity and detection limit (DL) values. Our study reveals a number of additional key features: (i) the DL is a property of a specific device, and crucially depends on the physical sensor length; (ii) the lowest DLs occur in regions where directional coupling occurs, as a result of beating between hybrid modes that produce narrow spectral widths; (iii) single-mode filtering at the multimode sensor output is fundamental for achieving such a narrow spectral feature; (iv) the high sensitivity at the EP can be seen in the modal beating of the hybrid EMs in the plasmonic region. These results can be immediately adapted to more realistic systems formed by 2D modes undergoing 3D propagation [27].

We first consider the plasmonic waveguide sensor shown in Fig. 2(a), where modes propagate in $z$, and different materials are distributed in $x$. Here we will limit ourselves to the analysis of a specific pure 2D geometry, which was designed to exemplify the different topological situations of EMs for analytes with the refractive index of aqueous solutions, as well as the various sensor performance and design considerations one can encounter. All waveguides are assumed to be infinite in $y$, providing physical insight into coupling mechanisms over a large parameter space and relatively fast calculation time. The system is here formed by a 1D slab of ${\mathrm{SiO}}_2$ (waveguide width: $d=1.95\mathrm{\mu m}$), coated with a finite length gold nanofilm on one side (thickness $t=30\mathrm{nm}$), and surrounded by an analyte (RI: $n_a=1.32–1.42$) elsewhere. This dielectric waveguide is multimode, supporting three to six modes between 400 and 800 nm. While we will only excite the fundamental mode at the input, having multiple modes at the output is of importance for our later analysis. The relevant plasmonic mode here is the long-range surface plasmon (LR-SPP) [43], whose real part can be designed to match with the fundamental mode of the dielectric waveguide. We choose this geometry as a starting point of the discussion for two reasons: first, because it most closely resembles commonly used fiber-based structures [15]; second, because long-range plasmons have a cutoff wavelength [44] that is highly sensitive to the environment and that can also be harnessed for enhancing sensors. For this gold film thickness, the (cutoff-free) short-range surface plasmon’s (SR-SPP’s) effective index and loss are too high to couple to the dielectric mode and are considered separately for a suitably modified geometry in Section 4. The interaction length $L$ corresponds to the length where the gold film and the analyte overlap. The modes in the sensing section are excited by the fundamental mode of the silica slab at $z=0$. At $z=L$, the superposition of modes of the sensor in turn excites modes of the silica slab. Experimentally, measurements consider the resulting transmitted wavelength-dependent intensity, shown schematically in Fig. 2(b). These measurements are characterized by a transmission minimum $T_{\min }$, a resonant wavelength ${\lambda }_R$, and a spectral width $\delta \lambda$, all of which depend on $n_a$.

There are two dominant interpretations of the cause of the resonant dip in waveguide plasmonic sensors. The first is that the dominant mode inside the sensor is lossy, with a loss peak that depends on the analyte’s index. In that interpretation, the wavelength-dependent loss of this dominant mode is measured. The second is that two modes are excited in the sensing section, and the output results from their interference—in short, that the plasmonic sensor acts as a directional coupler, with the PM wavelength being dependent on the analyte index. In this paper we will show in which circumstances each of these interpretations is correct, their limitations, and how to use their understanding to optimize sensing performance.

Figure 2(b) shows the main features of typical transmission spectra as $n_a$ increases. The position of the resonance wavelength ${\lambda }_R$ is a function of $n_a$, and the sensor’s sensitivity $S$ is defined by the shift of resonant wavelength per change in analyte refractive index (RI), $$S=\frac{\mathrm{d}{\lambda }_R}{\mathrm{d}{n}_a}.$$(1)

A refractive index sensor’s overall performance is best characterized by its DL, which is the smallest detectable change in RI $\delta n$, which generally depends on a specific user’s experimental configuration. Wu et al. [31] extended White and Fan’s [45] heuristic formula for the DL in the context of measurements that rely on resonant transmission spectra, most relevant for the present case, $$\delta n\approx \frac{1}{S}\frac{\delta \lambda }{1.5{(\mathrm{SNR})}^{0.25}},$$(2)where SNR is the signal-to-noise ratio in linear units near the transmission minimum $T_{\min }$, and $\delta \lambda$ is the spectral width at twice the minimum transmission, i.e., the $T_{\min }+3\mathrm{dB}$ limit, as shown in Fig. 2(b). The smallest detectable $\delta n$ would stem from a combination of narrow resonance width, high sensitivity, and low instrument noise. Since the SNR is dictated by the instrument, for the remainder of this work we shall consider only $S$ and $\delta n$. In general, a high $S$ and a small $\delta \lambda$ do not occur in the same region of parameter space. While $S$ could be inferred from modal calculations (in lossless dielectric directional couplers, for example, it can be extracted from the shifts in the PM wavelength [29–31,46]), $\delta \lambda$ also depends on the interaction length $L$, and thus requires that the modes’ excitation and propagation through the sensor be considered in detail. In lossy systems, even calculating $S$ can be difficult, as losses and interference compete, so that the minimum transmission does not necessarily occur at the PM wavelength or at the wavelength of maximum loss for any one mode.

We begin by considering all relevant bounded EMs shown in the Fig. 2(a) schematic. The propagation constants, as well as the electric and magnetic fields of each mode are obtained by numerically solving a complex transcendental dispersion equation resulting from enforcing boundary conditions between the layers [44]. The material dispersions for silica [47] and gold (Drude model in Ref. [48]) are taken into account, and the analyte index is taken as a wavelength-independent constant as labeled. Figure 3(a) shows a schematic of the modes considered. The two hybrid modes in the sensing region of Fig. 2(a) (propagation constants: ${\beta }_{i=\mathrm{1,2}}$) are shown as solid lines. For comparison, the isolated (uncoupled) EMs—supported by an equivalent dielectric waveguide without a gold film (propagation constant: ${\beta }_d$), or by an equivalent gold nanofilm sandwiched between silica/analyte on each side—are shown as dashed lines.

Figures 3(b)–3(d) show the associated EMs’ dispersion curves ($n_{\mathrm{eff},i}={\beta }_i/k_0$, $k_0=2\pi /\lambda$) for three values of $n_a$. To aid physical intuition, we show each mode’s real part $\Re e(n_{\mathrm{eff},i})$ [top of Figs. 3(b)–3(d)] and loss ${\alpha }_i$ [in $\mathrm{dB}/\mathrm{\mu m}$, bottom of Figs. 3(b)–3(d)], which is related to the imaginary part of $\Im \mathit{m}(n_{\text{eff},\mathit{i}})$ via $${\alpha }_i[\mathrm{dB}/\mathrm{\mu m}]=10{\log }_{10}\{\exp [2\Im \mathit{m}({\beta }_i)\times 1\mathrm{\mu m}]\}.$$(3)

Figure 3(e) also shows a detailed 3D plot of the dispersion of the hybrid EMs as a function of $n_a$. Note in particular the transition between a regime where $\Re e(n_{\mathrm{eff}})$ cross and $\Im \mathit{m}(n_{\mathrm{eff}})$ anticross [Fig. 3(b)], and a regime where $\Re e(n_{\mathrm{eff}})$ anticross and $\Im \mathit{m}(n_{\mathrm{eff}})$ cross [Fig. 3(d)], separated by an EP [Fig. 3(c)] where the propagation constants coalesce, with equal real and imaginary parts for both hybrid modes [35]. This kind of transition is characteristic of coupled non-Hermitian systems and is frequently encountered when designing plasmonic waveguide sensors [18,40]. The transition through the EP is similar to those in $\mathcal{P}\mathcal{T}$ symmetric systems [37], for example in situations including perfectly balanced optical gain and losses [37]. In such $\mathcal{P}\mathcal{T}$ symmetric situations, the eigenvalues can either be real and follow PT symmetry, or form complex conjugate pairs and thus break $\mathcal{P}\mathcal{T}$ symmetry. These two regimes are separated by the EP, where the eigenvalues coalesce. Structures with a global net loss, e.g.,  plasmonic [18,38] and leaky [51] waveguides, share many of the same features as non-Hermitian systems with no net loss: the main difference is that eigenvalues are shifted along the positive imaginary axis with respect to the perfectly loss-balanced case. Coupled lossy systems can thus be more rigorously classified as having eigenvalues that are either EPTS or EPTB, as labeled in Fig. 3(a). To simplify the discussion, we point out that in the present configuration and near regions where plasmonic resonances are measured, eigenvalues coalesce at the EP when $n_a=n_{\mathrm{EP}}=1.33$, $n_a<1.33$ supports EPTB modes, and $n_a>1.33$ supports EPTS modes, as highlighted by Figs. 3(b)–3(e). In the EPTB regime near the resonance, both hybrid EMs have field distributions with sizeable overlap with the dielectric core region and have comparable loss: they are both excited by the incoming dielectric mode, and their interference leads to directional coupling and energy exchange between the dielectric waveguide and metal film [27,28]. In the EPTS regime, the incoming dielectric mode predominantly excites the low-loss hybrid mode, effectively leading to monotonic exponential decay resulting from the transmission of the excited lossy mode [18].

One important question that arises when designing plasmonic refractive index sensors is the following. How can the location and shift of the resonant wavelength ${\lambda }_R$ be inferred from the modal calculations of Fig. 3? Various approaches may be found in the literature, as a result of the many possible choices available, which result from the transition between EPTS and EPTB regions near resonance. Inspecting the loss of hybrid mode 1 in the EPTB region [e.g.,  Fig. 3(b), $n_a=1.32$] might suggest a simple interpretation: the fundamental mode of the dielectric waveguide couples only to the lowest loss mode, and the resonant transmission $T(\lambda )$ (in dB) can be computed as $$T(\lambda )=\min [{\alpha }_1,{\alpha }_2]\times L,$$(4)where ${\alpha }_i$ is the loss of hybrid mode $i$ (in $\mathrm{dB}/\mathrm{\mu m}$) and $L$ is the interaction length (in $\mathrm{\mu m}$). In other words, the loss is computed from the minimum of the two hybrid mode loss curves of Fig. 3(b), and both ${\lambda }_R$ and $\delta \lambda$ follow immediately from ${\alpha }_1$ and are independent of $L$. $T_{\min }$ thus occurs where the loss of this curve is maximum, as highlighted by the label “max. loss” in Fig. 3(b). We will show that this approach is adequate in the EPTB region, but many reports extend this reasoning to the EP and EPTS regions [i.e., by computing the loss from the minimum of the two hybrid mode loss curves in Fig. 3(c) and Fig. 3(d) at each $n_a$] without accounting for mode excitation, propagation, and interference. This reasoning would imply that the resonance minimum $T_{\min }$ occurs at the wavelength where the two loss curves of the hybrid EMs in Fig. 3(d) intersect, i.e., the “loss matching” wavelength ${\lambda }_{\mathrm{LM}}$ [49], as indicated in Fig. 3(d). One other choice that can be borrowed from the dielectric waveguide literature to predict ${\lambda }_R$ is the “PM” wavelength ${\lambda }_{\mathrm{PM}}$, i.e., where the real parts of the uncoupled effective indices cross, as indicated by “PM” in Fig. 3(b). The PM wavelength would be a reasonable estimate of the resonant wavelength of the directional coupler consisting of the plasmonic and dielectric waveguides if these two waveguides were weakly coupled. However, plasmonic refractive index sensors use waveguides that are often in close proximity, invalidating this approximation. Furthermore, because the plasmonic guide has a cutoff, the PM point ceases to exist at higher analyte indices. In the EPTS region, where directional coupling occurs, one alternative choice is the wavelength ${\lambda }_{\mathrm{\Delta }{n}_{\mathrm{eff}}^{\min }}$ at the minimum effective index difference, i.e., where $\Re e(\mathrm{\Delta }{n}_{\mathrm{eff}}^{\min })=\min [\Re e(n_{\text{eff},2}-n_{\text{eff},1})]$, as highlighted in Fig. 3(d). This wavelength informs where coupling between adjacent waveguides is strongest [50] and resonant energy transfer is most efficient. Finally, another proxy for the resonant wavelength one can use is the cutoff wavelength of the plasmonic mode ${\lambda }_{\text{cutoff}}$ as indicated by “cutoff” in Fig. 3(b). This is perhaps less used in the plasmonic literature, but has been used to successfully analyze the high sensitivity of all-dielectric photonic crystal fiber sensors with liquid analyte satellites [29]. Note that the dependence of these parameters on $n_a$, should it correlate with ${\lambda }_R$, would only provide information on $S$ via Eq. (1), but that $\delta \lambda$ and $\delta n$ require knowledge of the full transmission spectrum, which generally cannot be obtained simply from modal dispersion curves alone.

We now evaluate the transmission spectra associated with the waveguide configuration shown in Fig. 2, accounting for the modes’ excitation and propagation along the device length $L$ as $n_a$ is varied. Here, we always consider the case where the input mode is the fundamental mode of the isolated dielectric waveguide (i.e., a silica core and an analyte cladding), corresponding to typical experimental conditions [15]. The fundamental mode couples to the two hybrid EMs of the plasmonic sensor, which then propagate through the sensing region before exciting the modes of the output dielectric waveguide. Although this simple model does not account for reflections and scattering at the boundary between the dielectric and hybrid waveguide, these have a negligible impact on the overall transmission spectra at the configurations considered, as we discuss below. For 1D modes propagating along $z$, the total electric and magnetic fields in the sensing region can be written as a superposition of its supported EMs, $$\mathbf{E}(x,z)=a_1{\mathbf{E}}_1(x)\exp (i{\beta }_{1}z)+a_2{\mathbf{E}}_2(x)\exp (i{\beta }_{2}z),$$(5)$$\mathbf{H}(x,z)=a_1{\mathbf{H}}_1(x)\exp (i{\beta }_{1}z)+a_2{\mathbf{H}}_2(x)\exp (i{\beta }_{2}z),$$(6)where the subscript $i=\mathrm{1,2}$ labels the two hybrid EMs in the sensing region. Here we use an implicit $\exp (-i\omega t)$ time dependence, and the convention that actual fields are the real parts of the complex fields. $\mathbf{H}$ and $\mathbf{E}$ are the transverse total magnetic and electric field vectors at any point in the plasmonic waveguide, respectively, ${\mathbf{H}}_i$ and ${\mathbf{E}}_i$ are the magnetic and electric field distributions of the hybrid EMs, ${\beta }_i={\beta }_i^R+i{\beta }_i^I$ are their propagation constants, and $a_i$ are the modal amplitudes that determine the contribution of the respective EMs to the total field. Strictly speaking, the above expansions are approximate, since contributions from higher-order modes and radiation modes have been neglected, but we shall be verifying their validity through comparisons with full-field finite-element calculations. Following the convention used in Ref. [52], modal fields are normalized so that $$\frac{1}{2}\int [{\mathbf{E}}_i(x)\times {\mathbf{H}}_j(x)]\mathrm{d}x={\delta }_{i,j}.$$(7)

This normalization does not use the complex conjugate of $\mathbf{H}$, since the latter is problematic in waveguides including material losses [52]. Note, however, that this normalization is equivalent to the complex conjugated version for purely lossless waveguides.

Modes can always be normalized to satisfy Eq. (7), which can then be used to calculate the complex modal amplitudes at input, leading to $$a_i=\frac{1}{2}\int [{\mathbf{E}}_i(x)\times {\mathbf{H}}_{\mathrm{in}}(x)]\mathrm{d}x,$$(8)where in this case ${\mathbf{H}}_{\mathrm{in}}(x)={\mathbf{H}}_1^d(x)$ is the magnetic field of the fundamental mode of the dielectric waveguide, corresponding to the input ($z=0$), as shown in Fig. 2(a). Note that because Eq. (7) is unconjugated and the waveguide includes losses, the sum of ${|a_i|}^2$ is not equal to the total power carried in the waveguide. At the output boundary between the plasmonic sensor and the dielectric waveguide ($z=L$), the field is a superposition of the EMs of the plasmonic waveguide, i.e., ${\mathbf{H}}_{\mathrm{out}}(x)=\mathbf{H}(x,L)$ via Eq. (6). The transmitted amplitude is then obtained by projecting this field onto the modes of the silica waveguide via $$t_i=\frac{1}{2}\int [{\mathbf{E}}_i^d(x)\times {\mathbf{H}}_{\mathrm{out}}(x)]\mathrm{d}x,$$(9)where ${\mathbf{E}}_i^d$ is the electric field distribution of each dielectric waveguide EM. In the purely lossless dielectric waveguide, the normalization Eq. (7) is identical to power normalization, and thus the total transmission is given by $T(\lambda )={\sum }_i|t_i(\lambda {)|}^2$; in cases where the output is a single mode waveguide, $T={|t_1|}^2$.

We calculate the wavelength-dependent modal amplitudes $a_i$, the propagation constants ${\beta }_i$, and all relevant fields, which in combination are used to calculate the amplitudes $t_i$ at the output interface of the plasmonic sensor region. Figure 4(a) shows the intensity transmitted by the fundamental mode [i.e., $T(\lambda )=|t_1(\lambda {)|}^2$] for three examples of interaction lengths $L=25\mathrm{\mu m}$ (solid blue curves), $L=37.5\mathrm{\mu m}$ (solid orange curves), and $L=50\mathrm{\mu m}$ (solid yellow curves), for increasing analyte index as labeled, associated with the mode calculations of Fig. 3. Before proceeding with a more detailed analysis, we can already note a number of important features. In the EPTB region ($n_a=1.32$), the transmission spectrum monotonically and exponentially decreases with increasing length, and the resonant wavelength remains nominally unchanged. The transmission spectrum near the EP ($n_a=1.33$) exhibits similar overall features, but with a larger overall loss. Note that the EP wavelength (${\lambda }_{\mathrm{EP}}=564\mathrm{nm}$) does not correspond to the resonant wavelength at $T_{\min }$ (${\lambda }_R=570\mathrm{nm}$), as highlighted by the dashed circle in Fig. 4(a). Further increasing the analyte index to $n_a=1.36$ leads to a transmission spectrum that oscillates with propagation length due to interference between the excited EMs, noting that here $T_{\min }$ and $\delta \lambda$ are smallest at the intermediate length ($L=37.5\mathrm{\mu m}$). One immediate consequence is that the lowest DLs $\delta n$ as per Eq. (2) will occur in EPTS regions where directional coupling is supported, and for specific analyte lengths. Finally, at larger analyte indices ($n_a=1.38$ and $n_a=1.40$), the resonant transmission spectrum becomes increasingly complicated: ${\lambda }_R$ and $\delta \lambda$ depend on length, and the EMs’ excitation, propagation, and loss contribute to the overall transmission in ways that are challenging to predict by inspecting the modal dispersion alone. This is particularly important from an applications standpoint, since higher analyte indices fundamentally yield larger sensitivities [29]. Finally, we note that in contrast to the EM theory presented here, a perturbative coupled mode analysis [53], while providing insight into overall non-Hermitian behavior, cannot rigorously be applied to the present case, because (i) the isolated plasmonic mode cuts off before longer wavelengths regions where hybrid modes are supported and where the highest sensitivities can be reached, and (ii) the two waveguides are physically connected, and thus no evanescent field is present. Analogous to previous analyses of similar systems, the EM method presented here is the most straightforward method to compute the mode properties and transmitted spectra.

Figure 4 also includes a comparison of the transmission spectra obtained via the EM method described above with full vector FEM simulations (COMSOL), shown as dashed lines. A port boundary condition at the input ensures that only the fundamental TM mode of the waveguide is excited. A port boundary condition at the output provides the option of considering either the amount of power in the fundamental mode [54] [Fig. 4(a)], or the total power transmitted by the waveguide [Fig. 4(b)]. Perfectly matched layers at every external boundary suppress any reflections in the simulation volume, absorbing all power carried by the dielectric waveguide at the output. We find an excellent agreement between our EM method and the FEM calculations in all cases, provided meshing of the gold region is fine enough $(<10\mathrm{nm})$. The EM method accurately predicts the transmission spectra with computations that are orders of magnitude faster than FEM, making large parameter sweeps more practical. Our FEM calculations indicate that the reflected power back into the input waveguide is $<0.03\%$ for $L>10\mathrm{\mu m}$ for the analyte indices and wavelengths discussed, resulting in negligible Fabry–Perot effects. Furthermore, the scattered power into radiation modes is $<2\%$ for the configurations considered. Such effects only lead to slight offsets in the total transmitted power, which a posteriori validates our EM model assumptions in the sensor configurations considered. Such small reflection and scattering effects are likely a result of the gold film being 2 orders of magnitude thinner than the dielectric waveguide, and located near the evanescent tail of the incoming dielectric mode. These properties are typical for plasmonic waveguide sensors but might need revisiting when approaching other designs.

For completeness, Fig. 4(a) shows the computed transmission spectrum from Eq. (4), in the example case of $L=37.5\mathrm{\mu m}$, as a dashed-dotted line. While Eq. (4) predicts the transmission in regions where the isolated mode is in the topological vicinity of the coupled modes (i.e., off-resonance), and may be considered adequate in finding an approximate location for ${\lambda }_R$ for EPTB regions, this model fails to be accurate at resonance for all cases, and fails entirely in the EPTS region.

Figure 4(b) shows the same calculations when considering the total power transmitted at output, i.e., $T(\lambda )={\mathrm{\Sigma }}_i|t_i(\lambda {)|}^2$ for two representative values of $n_a$ in the EPTB and EPTS regions. While this only marginally impacts the EPTB region ($n_a=1.32$), this is more significant in the EBTS region, where directional coupling is supported, because the output power is distributed among the available dielectric waveguide modes, increasing $T_{\min }$ and decreasing $\delta \lambda$. Figures 4(c) and 4(d), respectively, show the calculated axial component of the Poynting vector at resonance in the EPTB regime ($n_a=1.32$, ${\lambda }_R=560\mathrm{nm}$) and EPTS regime ($n_a=1.36$, ${\lambda }_R=610\mathrm{nm}$), using the FEM. The EPTB regime is dominated by absorption upon propagation as a result of metallic losses in the gold film; in the EPTS regime [Fig. 4(d)], a higher fraction of power is transferred to the gold surface upon propagation due to directional coupling. Note in particular that the output fields of the Fig. 4(d) configuration can couple to higher-order modes of the dielectric waveguide, which leads to the aforementioned reduction $\delta \lambda$ in the absence of a mode filter, as per Fig. 4(b).

From the above analysis, it immediately follows that, in order to minimize $\delta \lambda$ and thus $\delta n$, it is important to filter out higher-order modes at output, for example, by splicing the sensor with a single-mode fiber [10]. Henceforth, we shall therefore consider such filtering to be implemented and only consider power in the fundamental mode.

The power of the EM model is that it allows us to rapidly obtain resonance wavelength, spectral width, and extinction ratio from rapid calculations of the realistic transmission spectrum $T(\lambda )$ of the full device. Having noted the salient features of plasmonic sensors at different non-Hermitian regimes, and having validated our EM model, we now quantify the key parameters $S$, $\delta \lambda$, and $T_{\min }$ in detail, particularly with regard to their dependence on $L$ and $n_a$.

We begin by considering the performance of the sensor as a function of $n_a$ for the representative lengths $L$ considered so far. Henceforth, we will only consider the EPTS regime, which can yield the smallest $\delta \lambda$ due to directional coupling. We repeat the calculations for the three analyte lengths shown in Fig. 4, but with a much finer resolution on $n_a$. The resulting transmission spectra, as a function of $n_a$ and $\lambda$, are shown in Figs. 5(a)–5(c) for $L=25\mathrm{\mu m}$, $L=37.5\mathrm{\mu m}$, and $L=50\mathrm{\mu m}$, respectively. We immediately note that, although in all cases the resonances qualitatively redshift for increasing analyte index, their shape, sharpness, and location can change significantly. To quantify this further, Figs. 5(d)–5(f) show the associated resonant wavelength ${\lambda }_R$, 3-dB width $\delta \lambda$, and transmission minimum $T_{\min }$ as a function of $n_a$, for the same lengths as in Figs. 5(a)–5(c): $L=25\mathrm{\mu m}$ (blue), $L=37.5\mathrm{\mu m}$ (orange), and $L=50\mathrm{\mu m}$ (yellow). Figure 5(b) shows that the resonant wavelength is dependent on the length of the chosen device, particularly at longer lengths, as anticipated in our preliminary analysis of Fig. 4(a) at $n_a=1.38$ and $n_a=1.4$. Furthermore, high analyte indices show a fluctuation in the resonant transmission wavelength as the length of the sensor increases, and even multiple transmission minima—see also the spectra in Fig. 4(a), $n_a=1.40$. This behavior is challenging to interpret, as a result of the wavelength dependence of the complex dispersion profiles, where mode excitation, propagation, and losses all contribute to the total transmission spectrum in a nontrivial manner. The EPTB regime is dominated by the excitation of the lowest-loss hybrid EM over the entire wavelength range (i.e., ${|a_1|}^2>0.9$), with minimal contributions from the highest loss mode (i.e., ${|a_2|}^2<0.1)$. In contrast, the EPTS regime is characterized by broad wavelength regions (bandwidth: 50–100 nm) where both modes are excited (i.e., ${|a_{\mathrm{1,2}}|}^2>0.4$), analogous to earlier reports of broadband plasmonic directional couplers [28]. As a result, fixing the length as per Fig. 5 shows complete coupling only for certain combinations of $n_a$ and $\lambda$; elsewhere, incomplete coupling occurs, and local transmission minima occur over the associated wavelength region, with much larger $\delta \lambda$ values than the optimal case, and fluctuating resonant wavelength ${\lambda }_R$, as shown in Fig. 5(c) and the yellow line in Fig. 5(d).

To quantify this effect, Fig. 5(e) shows $\delta \lambda$ as a function of $n_a$: in the EPTS regime considered, it is possible to achieve a subnanometer $\delta \lambda$ for specific combinations of $L$ and $n_a$. This spectral width minimum is associated with a transmission minimum $T_{\min }$, as shown in Fig. 5(f) and would require a high enough SNR to be measured (here: $>60\mathrm{dB}$). However, small deviations in $n_a$ can change $\delta \lambda$ significantly for fixed $L$.

Both the sensitivity $S$ and resonant width $\delta \lambda$ thus depend on device length $L$ and analyte index $n_a$. Therefore, an optimized plasmonic waveguide sensor design (small $\delta n$) crucially requires a judicious combination of $n_a$ (chosen by the user) and $L$ (determined by the choice of target $n_a$ range) to ensure that $\delta \lambda$ is minimized. Most importantly, engineering and analyzing such designs require implementing the model presented here to appropriately account for hybrid mode excitation, propagation, and loss.

Recall that while the evaluation of $\delta \lambda$ demands full transmission spectra, recommended design criteria often infer the sensitivity $S$ from modal calculations, as a valuable proxy with fast calculation time [55], informing the first steps of a plasmonic sensor design. This approach relies on calculating the wavelength ${\lambda }_i$ at which a particular condition is met, as a function of $n_a$, and computing the associated sensitivity via $S_i=\mathrm{d}{\lambda }_i/\mathrm{d}{n}_a$. Because our fast numerical tool allows us to obtain ${\lambda }_R$ from $T(\lambda )$, we can now also compare how the sensitivity inferred from different conditions, obtained from mode calculations alone, compares with the $S=\mathrm{d}{\lambda }_R/\mathrm{d}{n}_a$, obtained from propagation calculations. Continuing our earlier analysis, we only consider the EPTS region ($n_a>1.33$), which can achieve the lowest $\delta n$ by complementing high sensitivities with subnanometer $\delta \lambda$ at specific $L$.

Figure 6(a) shows the loss-matching wavelength ${\lambda }_{\mathrm{LM}}$ (purple), the PM wavelength ${\lambda }_{\mathrm{PM}}$ (light blue), the wavelength where real part of the effective index difference is minimum ${\lambda }_{\mathrm{\Delta }{n}_{\mathrm{eff}}^{\min }}$ (green), and the plasmonic cutoff wavelength ${\lambda }_{\text{cutoff}}$ (black), as defined in Fig. 3 and associated text. Note in particular that ${\lambda }_{\mathrm{PM}}$ does not exist for $n_a>1.35$ (i.e., the real parts of the isolated EMs do not cross). Analogously, the imaginary parts of the hybrid EMs do not cross for $n_a>1.385$, and ${\lambda }_{\mathrm{LM}}$ does not exist in that region. Because ${\lambda }_R$ depends on $L$, we show ${\lambda }_R$ for both $L=25\mathrm{\mu m}$ (dark blue) and $L=37.5\mathrm{\mu m}$ (orange), obtained from the transmission calculations of Fig. 5. All curves show superlinear wavelength shifts with increasing $n_a$. The associated sensitivities are shown in Fig. 6(b) on a logarithmic scale. For the two example lengths shown, we find that the cutoff wavelength of the uncoupled plasmonic mode (black line) appears to be the best proxy for the sensitivity of the full (coupled) dielectric-plasmonic waveguide sensor (blue and orange lines), predicting the devices’ sensitivity within a factor of $\sim 3$ over the entire $n_a$ range. However, this is not the case if a wider range of lengths $L$ is considered. The shaded gray region in Fig. 6(a) and Fig. 6(b), respectively, shows the ranges of possible ${\lambda }_R$ and associated sensitivities for $L=10–2000\mathrm{\mu m}$, the latter of which spans more than an order of magnitude. This example serves to reinforce the message that no single criterion obtained from mode dispersion calculations can be used to accurately infer device sensitivity, which is crucially dependent on the device length.

We now calculate the transmission spectra for representative values of $n_a$, but with a much finer resolution on $L$. The resulting transmission spectra, as a function of $L$ and $\lambda$, are shown in Fig. 7(a). Figures 7(b)–7(d), respectively, show the associated resonant wavelength ${\lambda }_R$, 3-dB width $\delta \lambda$, and transmission minimum $T_{\min }$ as functions of $L$ for each analyte index $n_a$ of Fig. 7(a) as labeled. Figure 7(b) confirms that the resonant wavelength is nominally constant for all $L$ and $n_a\le 1.33$ (in the EPTB regime). For $n_a>1.33$ (in the EPTS regime), however, ${\lambda }_R$ is dependent on length and can fluctuate significantly, as per our earlier analysis.

Figure 7(c) shows that the spectral width $\delta \lambda$ decreases monotonically with increasing length for resonances occurring in the EPTB regime (i.e., $n_a=1.30$ and $n_a=1.32$) because only coupling to the lowest loss mode dominates, and resonant interference does not occur. In contrast, the EPTS regime supports resonant coupling effects from the dielectric core to the gold surface—or equivalently, resonant interference between the two hybrid modes—leading to subnanometer $\delta \lambda$ at specific interaction lengths $L$. It is striking to note that in the EPTS region, the wavelength at which complete coupling occurs [i.e., blue color map regions in Fig. 7(a)] can also depend on $L$. This effect becomes more prominent as $n_a$ increases; see, for example, the case of $n_a=1.4$, which is due to the differential wavelength-dependent losses between the two beating modes, leading to a change in the wavelength with the most complete destructive interference. Once again, $\delta \lambda$ minima are associated with local $T_{\min }$ minima, as shown in Fig. 7(c). In the EPTB regime, a monotonic decrease in $T_{\min }$ is observed, consistent with a behavior dominated by the lowest loss mode alone. In the EPTS regime, however, $T_{\min }$ oscillates as a result of directional coupling.

Let us now analyze how directional coupling is affected by $n_a$ and its relation to the EP. To simplify the discussion, we consider $L_R$ as the location of the first transmission minimum, highlighted by a dotted line in Fig. 7(d) for the case $n_a=1.34$, which roughly (but not exactly) corresponds to half a beat length (the two differ because of modal losses). Inspecting Fig. 7(d) already suggests that $L_R$ depends on $n_a$, and Fig. 8(a) quantitatively plots this dependence. As per our earlier analysis, the inset of Fig. 8(a) shows that the wavelength at which $L_R$ occurs is itself a function of $n_a$. Figure 8(a) clearly shows that $L_R$ increases dramatically close to the EP at $n_a=1.33$ and gradually reduces as the analyte index increases. We now show that this property is a direct consequence of the dispersion topology of Fig. 3(e) and indeed manifests key properties of perturbations near the EP.

Sensors that operate close to the EP have attracted much attention in recent years because the splitting of the coalesced eigenvalues is proportional to the square root of the change in the environment’s relative permittivity $\mathrm{\Delta }\epsilon$, with a much higher slope than the linear splitting for conventional Hermitian systems [56] and can thus be exploited for enhanced sensing [56–58]. The difference in effective index can be measured using an interferometric setup as a resonant dip with wavelength exquisitely dependent on analyte index, which one would expect is exactly what happens in our sensor through the beating between the two modes over a finite length. However, near the EP, the small splitting corresponds to beat lengths much larger than the loss length of the modes, so that dips are dominated by loss, rather than interference, and the benefits of the EP are lost. Intriguingly, to restore the benefits one could consider compensating losses by adding gain, which would get the sensor geometry closer to a true $\mathcal{PT}$-symmetric system.

Nonetheless, the characteristic square root dependence near the EP is readily found from the hybrid EM calculations of Fig. 3. To illustrate this, Figure 8(b) shows a plot of the associated minimum eigenvalue splitting, i.e., $\Re e(\mathrm{\Delta }{n}_{\mathrm{eff}}^{\min })=\min [\Re e(n_{\text{eff},2}-n_{\text{eff},1})]$, for each $\mathrm{\Delta }\epsilon =n_a^2-n_{\mathrm{EP}}^2$. For small values of $\mathrm{\Delta }\epsilon$, the eigenvalue splitting follows a square root dependence (black dashed line). For large values of $\mathrm{\Delta }\epsilon$, this is no longer the case because perturbations to the EP are large. Note that $\Re e(\mathrm{\Delta }{n}_{\mathrm{eff}}^{\min })$ occurs at wavelengths longer than the cutoff of the isolated plasmonic mode—a further indication of how removed this geometry becomes from the perturbative treatment.

One important subtlety, in the present context, is that the wavelength at which the minimum EM splitting occurs, shown in the inset of Fig. 8(b), is also a function of $\mathrm{\Delta }\epsilon$. Recalling that the beat length and EM splitting are related by $L_b={\lambda }_R/\mathrm{\Delta }{n}_{\mathrm{eff}}$—at least in the lossless case—Fig. 8(c) shows $\frac{1}{2}{\lambda }_R/L_R$ as a function of $\mathrm{\Delta }\epsilon$, as obtained from the transmission spectra of Fig. 7. We find a square root dependence (black dashed line) for small $\mathrm{\Delta }\epsilon$ perturbations, and deviations for larger values of $n_a$: the signature square root dependence on perturbation, close to the EP, is thus obtained from the device beat length. Interestingly, the deviation from the square root dependence is opposite to that predicted from mode calculations alone at a high-index region, once again as a result of the subtle interplay of mode excitation, interference, and loss—a further indication that mode dispersion calculations, in isolation, do not adequately predict the behavior of plasmonic sensors.

Having considered a plasmonic sensor that relies on the long-range surface plasmon, we now briefly discuss how our conclusions carry over to the case of short-range surface plasmons (SR-SPPs) [42,59]. SR-SPPs are ideally suited for chip-scale miniaturization due to their smaller effective modal area, lower group velocities, larger effective index, and no cutoff. In the present context, phase matching to the SR-SPP thus requires a different geometry than in previous sections, namely a higher-index core waveguide adjacent to the plasmonic film. For example, let us consider a dielectric waveguide formed by a silicon nitride core ($d=400\mathrm{nm}$, refractive index dispersion, Ref. [60]), surrounded by an infinite silica layer on one side, and a gold nanofilm of length $L$ on the other ($t=30\mathrm{nm}$), separated by a silica spacer (thickness, $s=111\mathrm{nm}$). This geometry possesses an EP at $\lambda =670\mathrm{nm}$ when the surrounding analyte index is $n_a=1.33$. Preliminary calculations indicate that many of the features and properties discussed so far carry over, with micrometer-scale interaction lengths due to the associated beat lengths [28]. As per our analysis so far, the spectral width $\delta \lambda$, resonant wavelength ${\lambda }_R$, and sensitivity $S$ all depend on the interaction length. Furthermore, regions of small $\delta \lambda$ occur for different combinations of $n_a$ and $L$, since they are linked to the characteristic beat length for a certain configuration. In all cases, an a priori quantitative analysis on the basis of mode dispersion calculations provides only an order-of-magnitude estimate of the sensitivity, and optimization of device performance requires a case-by-case analysis that considers propagation through the sensor to obtain the full DL. Note that sensitivities using SR-SPPs in this geometry are about 1 order of magnitude lower than in the previously discussed geometry, which used LR-SPPs. Recent chip-based waveguide sensors harnessing plasmonic Mach–Zehnder interferometry have shown that the resonance due to directional coupling can be fine-tuned by appropriately phase shifting a reference waveguide arm [61], which, in the present context, could be adapted to target a minimum $\delta \lambda$ in regions of high $S$.

In conclusion, we have comprehensively evaluated the sensing properties of plasmonic waveguide sensors by calculating their resonant transmission spectra in different regions of the non-Hermitian EM space. Our study highlights the limits of using modal dispersion calculations alone to predict plasmonic sensor performance and transmission spectra. These limits are easily addressed by using the same modal calculations in the framework of a model that accounts for both excitation and propagation of the EMs supported by the sensor. The resulting transmission calculations faithfully reproduce the transmission spectra, verified via a comparison with full-vector finite-element calculations, with the added benefit of allowing for a rapid sweep over three important parameters (wavelength and analyte index, but most importantly, device length), in turn revealing many important aspects that have so far eluded discussion in the context of practical devices. By increasing the resolution on the full wavelength and analyte parameter space, we showed that no single mode dispersion criterion can be used as a proxy for sensitivity. Indeed, the highest DLs occur where directional coupling is supported (via subnanometer spectral linewidths) and close to plasmonic cutoffs. The latter suggests revisiting sensor performance in cylindrical fibers/wires, close to the cutoff of high-order long-range cylindrical plasmonic modes. Near the EP, the hybrid plasmonic modes yield a characteristic square root dependence of the EM splitting with respect to the permittivity perturbation of the sensor, which in this context is identified through the sensor beat length. The square root dependence theoretically leads to high slopes of the coupling length versus refractive index, and thus high sensitivity. However, the small difference in effective index between beating modes in this supralinear region near the EP corresponds to beat lengths that are much larger than the modes’ loss length. In this case, the transmission dip is dominated by loss rather than interference. This makes it difficult to exploit the supralinear behavior near the EP in lossy systems. Theoretically, this issue could be avoided if one could compensate losses by introducing gain, thereby bringing the sensor closer to a true $\mathcal{P}\mathcal{T}$-symmetric system. Note that recent experiments on nanofluidic-core fiber platforms [62] showed that hybrid mode excitation and propagation can be directly visualized via sideways-detected fluorescence, which could be adapted in the present context for novel, single-wavelength sensing avenues that rely on beat-length measurements. Our analysis will find widespread applications in a variety of waveguide-based refractive index sensors, whose theoretical performance, in some instances, might require revisiting.

Acknowledgment. A. T. and B. T. K. thank Pranav A. Alavandi and Zachary J. R. Mann for fruitful discussions.