1. Introduction

Electromagnetic resonances enable precise control and enhancement of spontaneous emission and other light–matter interactions. While it is well known that resonances can enhance spontaneous emission rates via the celebrated Purcell effect [1–3] by confining light to small volumes for long times, recent work [4–6] suggests that giant enhancements can occur via the less familiar Petermann effect [7–11]. The Petermann enhancement factor is a measure of non-orthogonality of the modes in non-Hermitian systems and it appears to diverge when two modes coalesce at an exceptional point (EP)—an exotic degeneracy in which two modes share the same frequency and mode profile [12, 13]. In recent years, there has been an explosion of interest in EPs due to the many interesting and counter-intuitive phenomena associated with them, e.g., unidirectional reflection and transmission [14–16], topological mode switching [17–19], intrinsic single-mode lasing [20, 21], and lasers with unconventional pump dependence [22–24]. An understanding of spontaneous emission at EPs is essential for their implementation in optical devices, but the existing theory is limited to one-dimensional [25, 26] or discrete oscillator systems [27].

In this paper, we present a general theory of spontaneous emission near EPs. Our theory extends beyond spontaneous emission to any light–matter interaction described by the local density of states (LDOS) or, more precisely, any situation in which one analyzes the contribution of a given resonance to the emission of a source, such as narrowband thermal emission [28–31], absorption [32], perfect coherent absorption [33–35], and nonlinear harmonic generation [36]. Whereas traditional theories of spontaneous emission imply infinite enhancement factors at EPs (since the Petermann factor diverges), we use a modified Jordan-form-based perturbation theory to derive (finite) bounds on the enhancement at second- and higher-order EP systems. We show that line narrowing leads to a maximum enhancement of 4 in passive systems with second-order EPs and significantly larger enhancements (exceeding 400 ×) in gain-aided and higher-order EP systems. Our analytical results are presented in Sec. 2, where we express the emission rate at an EP in terms of the degenerate mode and its corresponding Jordan vector. This derivation assumes negligible dispersion, but we show in appendix B that the effect of dispersion amounts to merely modifying the normalization of the resonant modes, changing the results quantitatively but not qualitatively. Then, in Sec. 3, we demonstrate the implications of our theory via a concrete numerical example of coupled plasmonic resonators. Motivated by the fact that an EP is associated with a double pole in the Green’s function, we find specific locations where the emission lineshape becomes a squared Lorentzian, with peak amplitude scaling as Q², where Q is the resonance quality factor (a dimensionless measure of the resonance lifetime). We show that the enhancement at the EP is thus, potentially, much larger than the Purcell factor, which scales linearly with Q. Then, in Sec. 4, we derive bounds on the maximal enhancement at an EP, and we explore these bounds using a periodic system, which allows us to tune gain, loss, and degeneracy independently. Our theory provides a quantitative prescription for achieving large enhancements in practical optical systems, which is applicable to arbitrary geometries and materials and can be implemented with the recent experimental realizations of EPs [16–18, 20–22, 37, 38].

Traditional enhancement formulas fail at EPs since they are based on non-degenerate perturbation theory, which is invalid at EPs. Standard perturbation theory relies on Taylor expansions of differentiable functions while, near EPs, eigenvalues change non-analytically in response to small matrix perturbations. Instead, one needs to use a Jordan-form-based perturbative expansion [39, 40]. Although Jordan-vector perturbation theory is well known in linear algebra, along with related results on resolvent operators, these algebraic facts have not previously been applied to analyze Purcell/Petermann enhancement or LDOS in a general EP setting. By using such a modified expansion, we obtain a quantitative formula for the LDOS, which is a measure of how much power a dipole source can radiate [41] or, equivalently, a measure of its spontaneous emission rate. Note that similar expansion methods were previously used to evaluate the Green’s functions at EPs [42–44]; however, these works were limited to one-dimensional and paraxial systems, and were not applied to study spontaneous emission. An alternative semi-analytic approach was presented in [25]. In this work, the authors applied a scattering matrix formulation to model a simple one-dimensional system and analyzed its emission properties under $\mathcal{P}\mathcal{T}$-symmetry conditions. With proper modifications, such an analysis could be generalized to handle more complicated one-dimensional structures (e.g., continuously varying media or complex layered media). However, our general formulas [Eqs. (9) and (10)] can be directly applied to any system with an EP (e.g., three-dimensional photonic or plasmonic structures). Our results demonstrate that the unique spectral properties at EPs are general and do not rely on certain symmetry or dimensionality. Moreover, our theory enables modeling complex experimental apparatuses, performing numerical optimization and design, and deriving bounds on the enhancement, thereby clarifying the usefulness and limitations of EPs for enhancing light matter interactions.

Formally, the Petermann factor is inversely proportional to the “unconjugated norm” of the resonant mode ${\displaystyle \int dx\epsilon {\mathbf{E}}_n^2}$, where E_n is the mode profile and ε is the dielectric permittivity (with modifications to this “norm” when treating dispersive media [45].) At an EP, the unconjugated norm vanishes, ${\displaystyle \int dx\epsilon {\mathbf{E}}_n^2}=0$ [13] (a property also called “self-orthogonality” [12]), and the Petermann factor diverges. In fact, the Petermann factor can only diverge at an EP. This is because the Petermann factor is proportional to the sensitivity of an eigenvalue to perturbations [46] (its “condition number”), and that sensitivity can only diverge when two eigenvectors coalesce (i.e., at an EP) by the Bauer–Fike theorem [46]. This implies that our theory is applicable to any system exhibiting a giant Petermann factor. Specifically, in any laser and optical parametric oscillator (OPO) system where a giant Petermann factor was identified [47–49], there must have been a nearby “hidden” EP.

2. Local density of states and Green’s function expansions

In the following section, we give some background on LDOS calculations in non-degenerate systems, i.e., systems without EPs (Sec. 2.A), and then we review perturbation theory for systems with EPs (Sec. 2.B). Finally, in Sec. 2.C, we present a condensed derivation of our key analytical result—a formula for the LDOS at an EP [Eqs. (9) and (10)].

2.A. LDOS formula for non-degenerate resonances

The spontaneous emission rate of a dipolar emitter, oriented along the direction ê_µ, is proportional to the local density of states (LDOS) [50–52], which can be related to the dyadic Green’s function G via [41, 50, 53]

Computationally, one can directly invert Maxwell’s equations to find G and evaluate Eq. (1), but this provides little intuitive understanding. A modal expansion of the Green’s function, when applied properly, can be more insightful. Away from an EP, one can use the standard modal expansion formula for non-dispersive media [54]:

In many cases of interest, one can get a fairly accurate approximation for the LDOS by including only low-loss resonances in the Green’s function expansion [Eq. (2)] since only those contribute substantially to the emission spectrum. Under this approximation (i.e., considering only resonances ω_n = Ω_n − iγ_n which lie close to the real axis in the complex plane, with γ_n ≪ Ω_n), the spectral lineshape of the LDOS reduces to a sum of Lorentzian functions, weighted by the local field intensity:

Near the resonant frequencies, ω ≈ Ω_n, the peak of the LDOS scales linearly with the resonance quality factor $Q_n\equiv \frac{{\mathrm{\Omega }}_n}{2{\gamma }_n}$, leading to the celebrated Purcell enhancement factor [1]. On the other hand, the “unconjugated norm,” $({\mathbf{E}}_n^L,{\mathbf{E}}_n^R)$, which appears in the denominator of Eq. (3) leads to the Petermann enhancement factor, defined as $({\mathbf{E}}_n^{R\ast },{\mathbf{E}}_n^R)({\mathbf{E}}_n^{L\ast },{\mathbf{E}}_n^L)/{\left|({\mathbf{E}}_n^L,{\mathbf{E}}_n^R)\right|}^2$ [5]. (In non-Hermitian systems, the mode profiles (E_n) are complex and the Petermann factor is, generally, greater than one. At the extreme case of an EP, the unconjugated norm in the denominator vanishes and the enhancement factor diverges. However, this divergence does not properly describe LDOS or spontaneous emission at EPs since Eq. (3) is invalid at the EP. That is because the derivation of Eq. (2) assumes that the set of eigenvectors of the operator $\widehat{\mathcal{A}}$ spans the Hilbert space, but this assumption breaks down at the EP. In order to complete the set of eigenvectors of $\widehat{\mathcal{A}}$ into a basis and obtain a valid expansion for the Green’s function and the LDOS at the EP, we introduce in the following section additional Jordan vectors [39, 46, 60].

2.B. Jordan vectors and perturbation theory near EPs

At a (second order) EP, the operator ${\widehat{\mathcal{A}}}_0$ is defective—it does not have a complete basis of eigenvectors and is, therefore, not diagonalizable. However, one can find an eigenvector $({\mathbf{E}}_0^R)$ and an associated Jordan vector $({\mathbf{J}}_0^R)$, which satisfy the chain relations [39]:

Near the EP, on can find a pair of nearly degenerate eigenvectors and eigenvalues that satisfy

2.C. LDOS formula at exceptional points

Near the EP (i.e., for small but non-zero p), one can use the non-degenerate expansion formula Eq. (2) to compute G. In order to compute G at the EP, we substitute the perturbative expansions for λ_± and E_± [Eq. (7)] into Eq. (2) and take the limit of p going to zero, namely:

Crucially, Eqs. (9) and (10) yield finite results at the EP resonance frequency Ω_EP. Moreover, a squared Lorentzian prefactor, $\frac{{\mathrm{\Omega }}_n}{2\pi }{\left(\frac{{\gamma }_n}{{\left(\omega -{\mathrm{\Omega }}_n\right)}^2+{\gamma }_n^2}\right)}^2$, replaces the traditional Lorentzian prefactors near non-degenerate resonances, $\frac{1}{\pi }\frac{{\gamma }_n}{{\left(\omega -{\mathrm{\Omega }}_n\right)}^2+{\gamma }_n^2}$ [compare with Eq. (3)]. This unique spectral feature follows directly from the existence of a double pole in the modified expansion formula for G [the first term in Eq. (9)]. As shown in the following section, a squared Lorentzian lineshape implies a narrower emission peak and greater resonant enhancement in comparison with a non-degenerate resonance at the same complex frequency. The strength of our formulation [Eqs. (9) and (10)] is that it is applicable to arbitrary structures and can therefore be used to design and optimize complex three-dimensional geometries with EPs. Moreover, it clarifies the conditions under which LDOS enhancement is maximal; essentially, determined by the relative magnitude of the two terms in Eq. (9), depending on the location of the emitter.

3. Properties of the LDOS at EPs

In this section, we explore the spectral emission properties at EPs [following Eqs. (9) and (10)] via a numerical model system of coupled plasmonic resonators.

3.A. Example: EPs in plasmonic systems

A convenient approach for tailoring the LDOS in practice is by using plasmonic resonances in metallic structures, which can enable ultra-high LDOS enhancements and are widely used in various applications [61]. In this subsection, we numerically explore a system of two plasmonic resonators with approximate parity-time $(\mathcal{P}\mathcal{T})$ symmetry. $\mathcal{P}\mathcal{T}$ symmetric systems in optics are characterized by having balanced distributions of gain and loss, and are known to possess EPs at critical gain/loss values at which the mode profiles undergo spontaneous symmetry breaking [62]. Our numerical setup is shown on the left-hand side of Fig. 1(a). It includes two rods with metallic cores (ε = −2.3 + 0.0001i) and a silica coating (ε = 2.1316) surrounded by air (ε = 1). The dimensions of the structure (b = 0.643a, t₁ = 0.536a, t₂ = 0.16a, and t₃ = 0.268a) were tuned in order to have two nearly degenerate resonances that are spectrally separated from the neighboring resonances of the structure. We implement outgoing boundary conditions using perfectly matched layers (PML) [41]. In order to preserve the approximate $\mathcal{P}\mathcal{T}$ symmetry of the system, gain and loss are added symmetrically to the outer parts of the coating. By brute-force optimization, we find that an EP occurs when the gain/loss value is |Im ε| ≈ 0.0002551 and the background permittivity of the upper half of the silica coating is slightly shifted to ε ≈ 2.129. The right-hand side of Fig. 1(a) depicts the trajectories of the two eigenvalues (red and blue curves) that merge at the EP (orange dot) upon varying the gain/loss parameter |Im ε|. More details on the numerical optimization procedure are given in appendix D.

 

Fig. 1 Plasmonic system with exceptional points. (a) Left: Two plasmonic resonators with silver cores and a silica coating, with gain and loss (red/blue) added to the outer sides of the coating. Right: Eigenvalue trajectories in the complex plane upon increasing gain/loss. The trajectories merge at the EP (orange dot). (b) The Petermann factor, $({\mathbf{E}}_n^{R\ast },{\mathbf{E}}_n^R)({\mathbf{E}}_n^{L\ast },{\mathbf{E}}_n^L)/{\left|({\mathbf{E}}_n^L,{\mathbf{E}}_n^R)\right|}^2$, diverges at the EP (left) while the LDOS_x [evaluated using Eq. (2)] remains finite (right), since the giant contributions of the terms in Eq. (2) have opposite signs (blue/red curves, scaled by 10⁻³).

Download Full Size | PDF

We discretize Maxwell’s equations using a finite-difference frequency-domain (FDFD) approach [63] and compute the LDOS by taking imaginary part of the Green’s function. The Green’s function is found via three methods: (i) directly inverting Maxwell’s partial differential equation, (ii) using the non-degenerate expansion formula, Eq. (2), which is valid away from the EP, and (iii) using the degenerate formula, Eq. (9), at the EP. In principle, one can use the non-degenerate formula, Eq. (2), very close to the EP to compute the LDOS, relying on cancellation between the terms to yield a finite result. However, non-degenerate perturbation theory obfuscates the finite enhancement at the EP, confusing previous estimates of LDOS, whereas our formulation, based on degenerate perturbation theory, naturally captures the finite enhancement.

Our numerical results are shown in Fig. 1(b). The plot on the left shows the Petermann factor, which diverges at the EP, while the right-hand side plot shows the finite LDOS. We excite TM modes and, therefore, compute the LDOS for x-polarized modes (denoted by LDOS_x). The red and blue curves show the two terms that contribute to the sum in Eq. (2). Upon approaching the EP, the individual contributions to the sum diverge with opposite signs, while their sum (black curve) remains finite. The red and blue curves are scaled by 10⁻³ for ease of presentation.

3.B. Simplified model for the LDOS

Although our general formula for the LDOS at the EP can be directly applied to the plasmonic system of Fig. 1, it is useful to introduce a simplified model to interpret the results. First, we project Maxwell’s operator, $\widehat{\mathcal{A}}={\epsilon }^{-1}\nabla \times \nabla \times$, onto the two-dimensional subspace spanned by the nearly degenerate eigenvectors near the EP, thus producing a reduced 2 × 2 matrix, $\widehat{\mathcal{A}}$ (which gives a meaningful description of the system as long as the emission spectrum is dominated by the two coalescing resonances). Second, we employ an approach similar in spirit to coupled-mode theory (CMT) [64, 65], which involves expressing the modes of the two-rod system (the coupled system) in terms of modes of two reduced systems (the uncoupled systems), containing only one or the other rod. Such an approach is valid whenever the rod–rod separation (t₁) in the coupled system is large enough so that the frequency splitting induced by the coupling [κ depicted in Fig. 1(a)] is smaller than the uncoupled resonance frequencies. Denoting by U and V the matrices whose columns are the right and left eigenvectors of the uncoupled system ( ${\mathbf{E}}_{1,2}^R$, and ${\mathbf{E}}_{1,2}^L$ corresponding to rod 1 and 2 respectively), we find in appendix E that the reduced operator is

Next, we obtain a simplified formula for the LDOS at the EP. The reduced matrix $\widehat{\mathcal{A}}$ has an EP at the critical gain/loss value: η = κ. Denoting the defective matrix by $\widehat{\mathcal{A}}$₀ and the identity opertor by 11, we can relate the full Green’s function at the EP to $\widehat{\mathcal{A}}$₀ via:

The advantage of this formulation becomes apparent when evaluating the LDOS in close proximity to one of the resonators. For instance, near the gain region [e.g., at r₀ in Fig. 1(a)], the lossy mode intensity, $\propto \left|{\mathbf{E}}_2^R{\mathbf{E}}_2^L\right|$, is negligible, so it follows from Eq. (12) that the LDOS is proportional to the first diagonal entry of the resolvent

 

Fig. 2 Spectral properties of LDOS at EPs. (a) LDOS peak vs. quality factor Q for five coupling values κ, showing quadratic (linear) scaling for large (small) κ values. Q is varied by adding gain to the coating, while κ is varied by changing the rod–rod separation. Inset: Red and blue points from the main plot, on a log-log scale. (b) Normalized linewidth (FWHM) at the EP, Γ/Γ₀, vs. normalized coupling, κ/γ, computed using: 1. FDFD-discretization of Maxwell’s equations (dots), and 2. the CMT-based linewidth formula, Eq. (16) (black line). The limit as κ/γ → ∞ is shown in purple (dashed line). (Γ₀ ≡ 2γ is the FWHM of a Lorentzian curve, see text.) Inset: LDOS_x at the EP, computed via FDFD (black) and via CMT [Eq. (14), red dashed line].

Download Full Size | PDF

3.C. Quadratic scaling and linewidth narrowing

In this subsection, we apply our CMT-based simplified formulas [Eqs. (12)–(14)] to evaluate the LDOS_x near the upper rod [see Fig. 1(a)]. As can be seen from Eq. (14), the LDOS peak value scales as

Another consequence of a squared Lorentzian lineshape is a narrower emission peak, compared to that of a standard Lorentzian spectrum. The full-width half maximum (FWHM) of a standard Lorentzian curve, $\frac{1}{\pi }\frac{{\gamma }_n}{{\left(\omega -{\mathrm{\Omega }}_n\right)}^2+{\gamma }_n^2}$, is Γ₀ = 2γ (where Ω₀ ± Γ₀/2 is the frequency at which the Lorentzian drops to half of its maximal value). On the other hand, following Eq. (14), the FWHM of the LDOS near an EP with Im[ω_EP] = γ is

4. LDOS enhancement at EPs

In the previous section, we showed that a squared Lorentzian lineshape can lead to enhanced emission rates and reduced linewidth. In this section, we quantify the enhancement at the EP and study its bounds. To demonstrate our results, we employ another numerical example: a periodic waveguide (Sec. 4.A). This structure allows us to independently tune gain/loss and degeneracy and, therefore, demonstrate the impacts of the two effects separately. We treat both passive (Sec. 4.B) and active (Sec. 4.D) systems (i.e., systems without and with gain respectively), and then generalize our results (in Sec. 4.C) to higher-order EPs, which form at the coalescence of multiple resonance.

4.A. Example: EPs in periodic structures

To demonstrate our results in this section, we numerically explore the periodic structure shown in Fig. 3(a). The modes of a periodic system are Bloch wavefunctions of the form E(r) = E_k (r)e^ik·r, where E_k (r) is a periodic function and k is the wavevector [66]. At each k, the mode E_k (r) solves an eigenvalue problem of the form [64, 65]: $\widehat{\mathcal{A}}(\mathbf{k}){\mathbf{E}}_{nk}={\omega }_{nk}^2{\mathbf{E}}_{nk}$, where $\widehat{\mathcal{A}}(\mathbf{k})\equiv {\epsilon }^{-1}(\nabla +i\mathbf{k})\times (\nabla +i\mathbf{k})\times$ and n labels the band. The figure of merit for spontaneous emission in periodic structures is the LDOS per wavevector k and field component µ, which is a measure of the power expended by a Bloch-periodic dipole source with a particular k-vector, in the presence of an electromagnetic field polarized along direction µ. We abbreviate it as LDOS_k (also called the mutual density of states [67]). The LDOS_k can be integrated over k to obtain the LDOS from an isolated (non-periodic) point source in the periodic structure [68]. However, the effect of the EP is much clearer in the integrand (LDOS_k) than in the integral (LDOS), and so we focus here on the former for illustration purposes, exploiting the fact that k allows us to control how close we are to an EP without altering losses.

 

Fig. 3 Passive periodic waveguides with EPs. (a) Top left: Periodic waveguide with outgoing boundary condition in the transverse direction. Bottom left: Field patterns of TM modes at k = 0 (only the real part is shown). Top right: Degenerate mode and Jordan vector at k_EP. Bottom right: Resonance [ω_n (k)] vs. k-vector, showing an EP (orange dot) at k_EPa/2π ≈ 1.46 × 10⁻⁴. (b) Top plot, positive y-axis: LDOS_k at r₀ for four k-values (k_EP, k₁, k₂, k_∞, marked on the lower plot). Negative y-axis: Resonances [ω_n (k)] in the complex plane. Bottom: Normalized LDOS_k peak (M_k/M_∞) vs. deviation from the EP (Δk) for four structures with different Q values, showing 4-fold enhancement at the EP. (M_k is the LDOS peak at k.)

Download Full Size | PDF

Our example system consists of a waveguide with periodic index modulation along its central axis, $\widehat{x}$ [Fig. 3(a)]. PML are used to truncate the transverse (ŷ) direction. The design parameters are: ε₁ = 12, ε₂ = 13.137, d = 0.51a, and t = 0.25a. These parameters were chosen so that the corresponding one-dimensional system (with the same ε₁, ε₂, d and infinite thickness, t → ∞) has nearly degenerate second and third frequency bands near k_x = 0 [guaranteed by choosing parameters close to the quarter-wave plate condition: $\sqrt{{\epsilon }_1}d=\sqrt{{\epsilon }_2}(a-d)$]. We force an EP for TM-polarized modes (E_z, H_x, H_y) by fine-tuning the wavevector k_x and the permittivity contrast δε ≡ ε₂ − ε₁ (see appendix D, Fig. 6). Figure 3(a) also depicts the real and imaginary parts of the coalescing eigenvalues (blue and red curves respectively). At k = 0 (also called the Γ point), Maxwell’s eigenvalue problem is complex-symmetric and, consequently, the eigenvectors ${\mathbf{E}}_1^R$ and ${\mathbf{E}}_2^R$ are orthogonal (see mode profiles in the lower-left corner, having even/odd symmetry when x is flipped to −x). At the EP, the eigenmodes ${\mathbf{E}}_{1,2}^R$ merge into a single degenerate mode ${\mathbf{E}}_0^R$ (upper-right panel). We also compute the Jordan vector $({\mathbf{J}}_0^R)$ with our novel iterative method [69], and use it to verify our formula for the LDOS at the EP [Eq. (10)].

4.B. Passive structures with EPs

Figure 3(b) depicts the LDOS_k near the periodic waveguide at several wavevectors (k_EP, k₁, k₂, k_∞ marked on the lower plot). Far from the EP (at k_∞), the LDOS_k is a sum of two non-overlapping Lorentzian curves. As k approaches the EP (e.g., at k₁ and k₂), the resonance peaks begin to overlap and the LDOS_k peaks increase due to the growing Petermann factor. Most importantly, for k values near but not equal to k_EP, the traditional modal expansion formula [Eq. (2)] approaches the limiting Jordan-form-based formula [Eq. (9)]. Physically, this means that structures with nearby EPs can be approximated by truly defective structures, making this analysis useful for experimental systems, which can only be close to but not exactly at an EP due to fabrication and calibration imperfections. Computationally, this implies that as long as Maxwell’s operator A(k) is not exactly defective, one can use Eq. (2) to evaluate the LDOS_k, when properly canceling the divergent terms in the sum.

The lower plot in Fig. 3(b) compares the enhancement, M_k/M_∞, for structures with varying quality factors, where we introduce the definitions M_k ≡ γ · max_ω[LDOS_k (ω)] and M_∞ ≡ γ · max_ω[LDOS_∞(ω)]. We change Q by modifying the permittivity contrast δε: Radiation losses decrease with decreasing index contrast of the periodic modulation [70], with the limit of zero contrast corresponding to infinite Q. By plotting the normalized LDOS_k peak (M_k/M_∞) vs. deviation from the EP (Δk ≡ k − k_EP), we find that the enhancement at the EP (M_EP/M_∞) is always four-fold, regardless of Q. This follows from a sum rule which states that the spectrally integrated LDOS (and, therefore, also the LDOS_k) is a constant [71]. It implies that the maximum LDOS at an EP [i.e., the peak of a square Lorentzian $\frac{M_{\text{EP}}{\gamma }^3}{{\left[{\left(\omega -{\mathrm{\Omega }}_{\text{EP}}\right)}^2+{\gamma }^2\right]}^2}$] is four times larger than the maximum LDOS at a non-degenerate resonance [i.e., the peak of a simple Lorentzian $\frac{M_{\infty }\gamma }{{\left(\omega -{\mathrm{\Omega }}_{\text{EP}}\right)}^2+{\gamma }^2}$].

4.C. Passive structures with higher order EPs

Motivated by recent interest in higher-order EPs [72–75], we generalize our results from the previous section to nth order degeneracies (i.e., EPs that form when n degenerate eigenvectors merge). In this case, we define the enhancement factor as the ratio of the LDOS peak at the EP and at a reference point with n non-degenerate resonances (generalizing our earlier definition for M_EP/M_∞). Following the arguments from Sec. 3.B, we expect to find a squared Lorentzian emission curve at second-order EPs, a cubic Lorentzian curve at third-order EPs, and a Lorentzian to the nth power, $L_n(\omega )=\frac{M_n{\gamma }^{2n-1}}{{\left[{\left(\omega -{\mathrm{\Omega }}_{\text{EP}}\right)}^2+{\gamma }^2\right]}^n}$, at nth order EPs. (This is essentially equivalent to a known result on the rate of divergence of the norm of the resolvent matrix as an nth order EP is approached [46].) From the sum rule mentioned above [71], the spectrally integrated LDOS at an nth order EP, $S_n={\displaystyle \int d\omega L_n(\omega )=\frac{M_n\sqrt{\pi }\mathrm{\Gamma }[n-\frac{1}{2}]}{\mathrm{\Gamma }[n]}}$, is equal to the integrated LDOS before the n resonances merge (here, Γ[n] is the gamma function). Realizing that the enhancement at the EP is maximal when merging n identical resonances, a bound can be computed by solving S_n = nS₁. We find that the enhancement at the EP is at most $M_n/M_1=\frac{\sqrt{\pi }\mathrm{\Gamma }[n+1]}{\mathrm{\Gamma }[n-\frac{1}{2}]}$. For example, at third order passive EPs (n = 3), the enhancement is at most 8-fold (as some of us recently confirmed in [76]).

These results imply that higher order EPs could potentially provide a new route for achieving order-of-magnitude enhancement of the LDOS and order-of-magnitude narrowing of the emission linewidth (in contrast to traditional methods, which typically aim to maximize the traditional Purcell factor by increasing the quality factor and reducing the mode volume [1]). However, this result does not necessarily mean that higher-order EPs will yield a larger LDOS than the best lower-order EP or single resonance. The reason is that the degrees of freedom that one would use to bring n resonances together might otherwise be employed to enhance the Q of an individual resonance. In practice, there will be a tradeoff between the quality of individual resonances and the order of the EP.

4.D. Active structures

We show in this section that much greater enhancements can be achieved in active systems, i.e., by introducing gain. Figure 4 compares four periodic waveguides with different index contrasts (δε), corresponding to different passive quality factors, Q_p. Gain and loss are added to each of the waveguides in order to force EPs, all of which share the same active quality factor Q_a. The empty (filled) markers in Fig. 4(a) are the EP resonances (ω_EP) in the complex plane before (after) adding gain/loss. As shown in Fig. 4(b), the structure with the smallest passive quality factor Q_p has the largest LDOS_k enhancement at the EP since it requires more gain in order to attain the same Q_a. Structures with smaller initial Q_p values can lead to even greater enhancements, and we present such a case in appendix F. In principle, the relative enhancement at the EP, M_EP/M_∞, is not bounded in this computational model. However, in practice, giant relative enhancements do not necessarily imply giant absolute LDOS values (since the value of the LDOS at the reference point may be very small). Moreover, when adding enough gain to bring the system to the lasing threshold, stimulated emission eventually limits line narrowing and LDOS enhancement at the EP.

 

Fig. 4 LDOS_k enhancement in active structures. (a) Resonances with different passive loss rates γ_p (hollow symbols) and equal active loss rates γ_a (full symbols). Inset: Structure from Fig. 3(a) with gain/loss added to the waveguide. (b) Normalized LDOS_k peak (M_k =M_∞) vs. Δk, evaluated at r₀ [marked in Fig. 3(a)], after adding gain/loss. The enhancement at the EP increases with the overall gain.

Download Full Size | PDF

Similar to the analysis of the plasmonic example in Sec. 3, we can introduce a simplified 2 × 2 model to interpret the results [analogous to Eq. (11)]. We project the full Maxwell’s operator $\widehat{\mathcal{A}}$ onto the subspace spanned by the two modes at k_x = 0 [shown in the left-lower corner of Fig. 3(a)] and we obtain

5. Discussion

The theory presented in this work provides a quantitative prescription for achieving large spontaneous emission rates using EPs, potentially exceeding by orders of magnitude those attained with standard non-degenerate resonances. Such enhancements could be useful for various applications, including fluorescent and Raman sensing [77], high-power low-coherence light sources [78], and sources with tunable coherence [79]. Moreover, by extending the current linear theory to account for EPs in nonlinear systems, our approach could be applied to study the properties of lasers at EPs—a topic that has recently drawn great attention in the optics community [20–22, 80].

Our formulation for the LDOS at the EP [Eqs. (9) and (10)] establishes that the enhancement generally consists of two terms, one that scales linearly and one that scales quadratically with the quality factor. In Secs. 3,4, we verified this scaling argument via two numerical examples, and we employed simplified 2 × 2 models to estimate the coefficients of the quadratic terms in the LDOS formula. (In the plasmonic system, we found that the coefficient was the coupling constant κ, while in the periodic example, the coefficient was the passive decay rate γ_p.) More generally, we show in appendix G that for arbitrary low-loss systems, the quadratic coefficient is bounded by

Finally, our theory extends beyond spontaneous emission, and can be applied to a broader class of phenomena described by the LDOS. For example, we anticipate similar enhancements in higher-harmonic generation rates in nonlinear media (e.g., Kerr media). In that case, the input lower-harmonic field (multiplied by the nonlinear susceptibility) acts as a source to the higher-harmonic field. To lowest order in the nonlinearity, the converted power is found by convolving the Green’s function with the input signal. This is analogous to our formulation of spontaneous emission (which involved convolving a dipole source with the Green’s function), and would therefore result in similar enhancements [81]. More generally, a similar treatment could produce enhancements in related quantities at EPs in other areas of physics (e.g., exciton-polariton [38] mechanical systems [18], and also leak-wave antennas both at radio-frequency and visible frequencies [82, 83]). Finally, our theory can be generalized to study scattering and extinction problems. We find, however, that even though EPs can produce special spectral features in scattering cross-sections, they do not give rise to giant scattering enhancements, since the scattered intensity is bounded by the incoming intensity (and the bound of perfect scattering can be easily achieved with a single non-degenerate resonance [84]).

Appendix

A. Non-degenerate Green’s function modal expansion formula

In this section we review the derivation of the Eq. (2) in the main text. Our derivation is similar to standard methods for Hermitian eigenvalue problems [54], with the necessary modifications for treating general non-Hermitian systems (most importantly, using the unconjugated “inner product” [12] between left and right eigenmodes).

Given a non-magnetic medium with dielectric permittivity ε, the fields E and currents J in the medium are related via Maxwell’s frequency-domain partial differential equation:

(19)$$({\widehat{\mathcal{A}}}_{\mathbf{r}}-{\omega }^2\mathbf{1})\mathbf{E}(\mathbf{r})=i\omega \mathbf{J}(\mathbf{r}).$$

where ${\widehat{\mathcal{A}}}_{\mathbf{r}}\equiv \frac{1}{\epsilon }{\nabla }_{\mathbf{r}}\times {\nabla }_{\mathbf{r}}\times$. The response to arbitrary currents can be found by convolving the dyadic Green’s function with the current sources: $\mathbf{E}(\mathbf{r})={\displaystyle \int d{\mathbf{r}}^{\prime }\mathbb{G}(\mathbf{r},{\mathbf{r}}^{\prime })}\mathbf{J}({\mathbf{r}}^{\prime })$, where $\mathbb{G}$ is defined via:

(20)$$({\widehat{\mathcal{A}}}_{\mathbf{r}}-{\omega }^2\mathbf{1})\mathbb{G}(\mathbf{r},{\mathbf{r}}^{\prime })=-\delta (\mathbf{r}-{\mathbf{r}}^{\prime })\mathbf{1}.$$

In this section, we expand G in terms of right and left resonant modes ( ${\mathbf{E}}_n^R$ and ${\mathbf{E}}_n^L$), which are outgoing solutions of the partial differential equations

(21)$$\begin{array}{l}({\widehat{\mathcal{A}}}_{\mathbf{r}}-{\omega }_n^2\mathbf{1}){\mathbf{E}}_n^R(\mathbf{r})=0\\ ({\widehat{\mathcal{A}}}_{\mathbf{r}}^T-{\omega }_n^2\mathbf{1}){\mathbf{E}}_n^L(\mathbf{r})=0.\end{array}$$

When the set of eigenvectors of Eq. (21) forms a complete basis of the Hilbert space, one can introduce the completeness relation, which consists of expanding the Dirac delta function as

(22)$$\delta (\mathbf{r}-{\mathbf{r}}^{\prime })\mathbf{1}={\displaystyle \sum _n{\mathbf{E}}_n^R(\mathbf{r})\otimes {\mathbf{E}}_n^L}({\mathbf{r}}^{\prime }),$$

where ⊗ is the outer/tensor product u ⊗ v = uv^T. The question of completeness of eigenmodes in non-Hermitian open systems has not been proven in general, for arbitrary three-dimensional systems. However, since in this work we always evaluate the Green’s function in close proximity to the resonators and at frequencies close to the resonance frequencies, the eigenmodes which overlap spectrally and spatially with the emitter give a good approximation for the LDOS, justifying the use of Eq. (22).

Similarly, we wish to find an expansion formula for $\mathbb{G}$ or, more explicitly, find the coefficients a_n(r′) in the series

(23)$$\mathbb{G}(\mathbf{r},{\mathbf{r}}^{\prime })={\displaystyle \sum _n{\mathbf{E}}_n^R}(\mathbf{r})\otimes {\mathbf{a}}_n({\mathbf{r}}^{\prime }).$$

To this end, we substitute Eq. (22) and Eq. (23) into Eq. (20) and obtain

(24)$${\displaystyle \sum _n\left({\mathcal{A}}_{\mathbf{r}}-{\omega }^2\right){\mathbf{E}}_n^R(\mathbf{r})\otimes {\mathbf{a}}_n({\mathbf{r}}^{\prime })=-}{\displaystyle \sum _n{\mathbf{E}}_n^R(\mathbf{r})\otimes {\mathbf{E}}_n^L({\mathbf{r}}^{\prime })}.$$

Next, we multiply both sides of the equation by ${\left[{\mathbf{E}}_m^L(\mathbf{r})\right]}^T$. Using the relation: ${\left[{\mathbf{E}}_m^L(\mathbf{r})\right]}^T{\mathcal{A}}_{\mathbf{r}}={\omega }_m^2{\left[{\mathbf{E}}_m^L(\mathbf{r})\right]}^T$, integrating over r, and invoking the bi-orthogonality relation for non-Hermitian systems [85]:

(25)$${\displaystyle \int dr{\mathbf{E}}_m^L(\mathbf{r})\cdot {\mathbf{E}}_n^R(\mathbf{r})}={\delta }_{m,n},$$

we obtain ${\mathbf{a}}_m({\mathbf{r}}^{\prime })=\frac{{\mathbf{E}}_m^L({\mathbf{r}}^{\prime })}{{\omega }^2-{\omega }_m^2}$. Finally, substituting this result in Eq. (23) we obtain an eigenmode expansion of the dyadic Green’s function [which reduces to Eq. (2) in the text]:

(26)$$\mathbb{G}(\mathbf{r},{\mathbf{r}}^{\prime },\omega )={\displaystyle \sum _n\frac{1}{{\omega }^2-{\omega }_n^2}}\frac{{\mathbf{E}}_n^R(\mathbf{r})\otimes {\mathbf{E}}_n^L({\mathbf{r}}^{\prime })}{{\displaystyle \int dr{\mathbf{E}}_n^R(\mathbf{r})\cdot {\mathbf{E}}_n^L(\mathbf{r})}}.$$

The integral in the denominator of Eq. (26) is one, but we keep it here for comparison with Eq. (2).

Last, we note that in reciprocal media ε = ε^T, and there exists a simple relation between left and right eigenvectors: ${\mathbf{E}}_n^L=\epsilon {\mathbf{E}}_n^R$. More generally, the left and right eigenvectors of a symmetric generalized eigenvalue problem (EVP): $\mathbb{A}{\mathbf{E}}_n^R={\lambda }_n\mathbb{B}{\mathbf{E}}_n^R$, ${\mathbb{A}}^T{\mathbf{E}}_n^L={\lambda }_n{\mathbb{B}}^T{\mathbf{E}}_n^L$, with A = A^T and B = B^T are related via ${\mathbf{E}}_n^L=\mathbb{B}{\mathbf{E}}_n^R$. To see this, rewrite the EVP for the left vectors as ${\left({\mathbb{B}}^{-1}\mathbb{A}\right)}^T{\mathbf{E}}_n^L=\mathbb{A}({\mathbb{B}}^{-1}{\mathbf{E}}_n^L)={\lambda }_n\mathbb{B}({\mathbb{B}}^{-1}{\mathbf{E}}_n^L)$ which shows that ${\mathbb{B}}^{-1}{\mathbf{E}}_n^L$ is a right eigenvector. Note that this relation holds also with Bloch-periodic boundary conditions (and the surface-term correction found in [86] does not appear in our formulation). Although the matrix $\mathbb{A}$ is no longer symmetric in the k-periodic problem, it satisfies the relation $\mathbb{A}{\left(k\right)}^T=\mathbb{A}(-k)$, and since we are essentially relating k-right eigenvectors to (−k)-left eigenvectors, the relation above remains unchanged.

B. Effect of dispersion on the LDOS formula

In this appendix, we consider the effects of dispersion on the Green’s function near and at the EP. In accordance with previous work on quasi-normal modes in dispersive media [87, 88], we find that the Green’s function has non-diagonal contributions $\propto {\mathbf{E}}_{\pm }^L\otimes {\mathbf{E}}_{\mp }^R$ near the EP. However, exactly at the EP, the Green’s function has exactly the same form as the non-dispersive degenerate formula [Eq. (9)], with dispersion affecting only the normalization of the degenerate mode (E₀) and Jordan vector (J₀).

In the same spirit as our derivation of the non-dispersive formula, we expand the Green’s function $\mathbb{G}$ in eigenmodes:

(27)$$\mathbb{G}(\mathbf{r},{\mathbf{r}}_0,\omega )={\displaystyle \sum _m{\mathbf{E}}_m^R}(\mathbf{r})\otimes {\mathbf{\alpha }}_m(\omega ,{\mathbf{r}}_0),$$

and use simple algebraic manipulations to express the coefficients α_m in terms of the modes. Recall that $\mathbb{G}$ is defined as the solution to the partial differential equation:

(28)$$\left[\nabla \times \nabla \times -{\omega }^2\epsilon (\mathbf{r},\omega )\right]\mathbb{G}(\mathbf{r},{\mathbf{r}}_0,\omega )=-\delta (\mathbf{r}-{\mathbf{r}}_0)\mathbf{1}$$

By multiplying both sides of Eq. (28) from the left by ${\mathbf{E}}_n^L(\mathbf{r})$ and integrating over dr, we obtain

(29)$${\mathbf{E}}_n^L({\mathbf{r}}_0)={\displaystyle \int dr{\mathbf{E}}_n^L}(\mathbf{r})\left[{\omega }^2\epsilon (\mathbf{r},\omega )-{\omega }_n^2\epsilon (\mathbf{r},{\omega }_n)\right]\mathbb{G}(\mathbf{r},{\mathbf{r}}_0,\omega )$$

Then, by substituting Eq. (27) into Eq. (29), we find

(30)$${\mathbf{E}}_n^L({\mathbf{r}}_0)={\displaystyle \sum _m{\displaystyle \int dr{\mathbf{E}}_n^L}(\mathbf{r})\left[{\omega }^2\epsilon (\omega )-{\omega }_n^2\epsilon ({\omega }_n)\right]{\mathbf{E}}_m^R}(\mathbf{r})\otimes {\mathbf{\alpha }}_m(\omega ,{\mathbf{r}}_0)\equiv {\displaystyle \sum _m({\mathbf{E}}_m,{\mathbf{E}}_n){\mathbf{\alpha }}_m(\omega ,{\mathbf{r}}_0)}.$$

(Note that this result was also derived in [88] by invoking Lorentz reciprocity.) Since we are interested in emission from emitters in close proximity to the resonators and in frequencies near the resonant frequencies, we may assume that a finite set of N eigenvectors adequately describes the system’s response [similar to our assumption in the non-dispersive derivation in Eq. (22)]. We introduce the vector s(r) = {E₁(r),…, E_N (r)} and the matrix ${\mathbb{O}}_{m,n}\equiv ({\mathbf{E}}_m,{\mathbf{E}}_n)$, and rewrite Eq. (30) as: $\mathbf{s}({\mathbf{r}}_0)=\mathbb{O}(\omega )\alpha (\omega ,{\mathbf{r}}_0)$ or equivalently $\alpha (\omega ,{\mathbf{r}}_0)={\mathbb{O}}^{-1}(\omega )\mathbf{s}({\mathbf{r}}_0)$. With this notation, Eq. (27) can be rewritten as: $\mathbb{G}=\alpha {(\omega ,{\mathbf{r}}_0)}^T\mathbf{s}(\mathbf{r})=\mathbf{s}{({\mathbf{r}}_0)}^T{({\mathbb{O}}^{-1})}^T\mathbf{s}({\mathbf{r}}_0)$ or

(31)$$\mathbb{G}(\mathbf{r},{\mathbf{r}}_0,\omega )={\displaystyle \sum _{nm}{({\mathbb{O}}^{-1})}_{nm}{\mathbf{E}}_n^R(\mathbf{r})}\otimes {\mathbf{E}}_m^L({\mathbf{r}}_0)$$

where

(32)$${\mathbb{O}}_{ij}={\displaystyle \int dr}{\mathbf{E}}_i^L[{\omega }^2\epsilon (\mathbf{r},\omega )-{\omega }_i^2\epsilon (\mathbf{r},{\omega }_i)]{\mathbf{E}}_j^R.$$

Equation 31 and Eq. (32) are the main result of this appendix—an eigenmode expansion of the Green’s function in the presence of dispersion, generalizing Eq. (2) from the main text. In the non-dispersive limit, this expression reproduces our previous result since ${\mathbb{O}}_{mn}=({\mathbf{E}}_m,{\mathbf{E}}_n)=({\omega }^2-{\omega }_n^2){\displaystyle \int dr\epsilon {\mathbf{E}}_m^L\cdot }{\mathbf{E}}_n^R=({\omega }^2-{\omega }_n^2){\delta }_{mn}$ and $\mathbb{G}={\displaystyle {\sum }_n{({\omega }^2-{\omega }_n^2)}^{-1}}{\mathbf{E}}_n^R(\mathbf{r})\otimes {\mathbf{E}}_n^L({\mathbf{r}}_0)$.

Next, let us calculate the limit of Eq. (31) when two modes E_± coalesce at an EP. Keeping only the terms corresponding to E_± in Eq. (31), expanding both ω and ω_± in Taylor series around ω₀ and introducing the notation:

(33)$${\omega }^2\epsilon (\omega )-{\omega }_{\pm }^2\epsilon ({\omega }_{\pm })=(\omega -{\omega }_{\pm })\frac{d({\omega }^2\epsilon )}{d\omega }|{ }_{\omega 0}\equiv (\omega -{\omega }_{\pm }){({\omega }^2\epsilon )}_0^{\prime },$$

we find that ${\mathbb{O}}_{ij}=(\omega -{\omega }_j)={\displaystyle \int dr{({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_i^L\cdot }{\mathbf{E}}_j^R$, and the matrix inverse is

(34)$${\mathbb{O}}^{-1}=\left(\begin{array}{cc}\frac{1}{\omega -{\omega }_{+}}\frac{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{-}^L\cdot {\mathbf{E}}_{-}^R}}{N}& -\frac{1}{\omega -{\omega }_{+}}\frac{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{+}^L\cdot {\mathbf{E}}_{-}^R}}{N}\\ -\frac{1}{\omega -{\omega }_{-}}\frac{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{-}^L\cdot {\mathbf{E}}_{+}^R}}{N}& \frac{1}{\omega -{\omega }_{-}}\frac{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{+}^L\cdot {\mathbf{E}}_{+}^R}}{N}\end{array}\right),$$

where $N\equiv \left[{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{+}^L\cdot {\mathbf{E}}_{+}^R}\right]\left[{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{-}^L\cdot {\mathbf{E}}_{-}^R}\right]-\left[{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{+}^L\cdot {\mathbf{E}}_{-}^R}\right]\left[{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{-}^L\cdot {\mathbf{E}}_{+}^R}\right]$. Substituting Eq. (34) into Eq. (31), we find that the Green’s function near the EP is

(35)$$\begin{array}{l}\mathbb{G}(\mathbf{r},{\mathbf{r}}_0,\omega )=\frac{{\mathbf{E}}_{+}^R\otimes {\mathbf{E}}_{+}^L}{\omega -{\omega }_{+}}\frac{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{-}^L\cdot {\mathbf{E}}_{-}^R}}{N}+\frac{{\mathbf{E}}_{-}^R\otimes {\mathbf{E}}_{-}^L}{\omega -{\omega }_{-}}\frac{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{+}^L\cdot {\mathbf{E}}_{+}^R}}{N}-\\ \frac{{\mathbf{E}}_{-}^R\otimes {\mathbf{E}}_{+}^L}{\omega -{\omega }_{+}}\frac{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{+}^L\cdot {\mathbf{E}}_{-}^R}}{N}-\frac{{\mathbf{E}}_{+}^R\otimes {\mathbf{E}}_{-}^L}{\omega -{\omega }_{-}}\frac{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{-}^L\cdot {\mathbf{E}}_{+}^R}}{N}.\end{array}$$

Last, we want to take the limit of Eq. (35) as the two modes E_± approach the EP. We expand ω_± and E_± around the EP using degenerate perturbation theory (using the notation of Sec. 2.B):

(36)$$\begin{array}{l}{\omega }_{\pm }\approx {\omega }_0\pm \sqrt{p}\mathrm{\Delta },\\ {\mathbf{E}}_{\pm }\approx {\mathbf{E}}_0\pm \sqrt{p}\mathrm{\Delta }{\mathbf{J}}_0.\end{array}$$

Generalizing our derivation of the non-dispersive formula, we choose the normalization conditions ${\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_0^L\cdot {\mathbf{J}}_0^R}={\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{J}}_0^L\cdot {\mathbf{E}}_0^R}$ =1 and ${\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{J}}_0^L\cdot {\mathbf{J}}_0^R}=0$, while the condition ${\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_0^L\cdot {\mathbf{E}}_0^R}=0$ is automatically satisfied due to self-orhotognality. With this normalization, we can calculate integrals in Eq. (35):

$$\begin{array}{l}{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{\pm }^L{\mathbf{E}}_{\pm }^R}=\cancel{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_0^L{\mathbf{E}}_0^R}}\pm 2\sqrt{p}\mathrm{\Delta }{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_0^L{\mathbf{J}}_0^R+p{\mathrm{\Delta }}^2\cancel{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{J}}_0^L{\mathbf{J}}_0^R}}+\mathcal{O}(p^{3/2})}\\ {\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_{\pm }^L{\mathbf{E}}_{\pm }^R}=\cancel{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_0^L{\mathbf{E}}_0^R}}-p{\mathrm{\Delta }}^2\cancel{{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{J}}_0^L{\mathbf{J}}_0^R}}+\mathcal{O}(p^{3/2})\\ N=-4p{\mathrm{\Delta }}^2{\displaystyle \int {({\omega }^2\epsilon )}_0^{\prime }{\mathbf{E}}_0^L{\mathbf{J}}_0^R}+\mathcal{O}(p^{3/2})\end{array}$$

The first two terms in Eq. (35) have the same form as the non-dispersive expansion [Eq. (8)]: $\frac{1}{2\mathrm{\Delta }\sqrt{p}}\frac{{\mathbf{E}}_{+}^R\otimes {\mathbf{E}}_{+}^L}{\omega -{\omega }_{+}}-\frac{1}{2\mathrm{\Delta }\sqrt{p}}\frac{{\mathbf{E}}_{-}^R\otimes {\mathbf{E}}_{-}^L}{\omega -{\omega }_{-}}+\mathcal{O}(p^{3/2})$. Upon approaching the EP, each contributions diverges with an opposite sign, and a finite contribution remains. The last two terms scale as $\sqrt{p}$ near the EP, and vanish at the EP. Finally, we obtain

(37)$${\mathbb{G}}_{\text{EP}}(\mathbf{r},{\mathbf{r}}_0,\omega )=\frac{{\mathbf{E}}_0^R(\mathbf{r})\otimes {\mathbf{E}}_0^L({\mathbf{r}}_0)}{{\left(\omega -{\omega }_0\right)}^2}+\frac{{\mathbf{J}}_0^R(\mathbf{r})\otimes {\mathbf{E}}_0^L({\mathbf{r}}_0)+{\mathbf{E}}_0^R(\mathbf{r})\otimes {\mathbf{J}}_0^L({\mathbf{r}}_0)}{\omega -{\omega }_0}.$$

Equation (37) implies that the green’s function at the EP in dispersive media has the same form as the non-dispersive formula [Eq. (9) in the text], but the normalization of the modes changes.

C. Convergence of the “unconjugated norm”

In this appendix, we show that when perfectly matched layers (PML) are used to implement outgoing boundary conditions, the “unconjugated norm” of a scattering eigenmode, ${\displaystyle \int dx\epsilon (\mathbf{x}){\mathbf{E}}_n}{(\mathbf{x})}^2$, converges to a finite result as the PML thickness tends to infinity. Scattering eigenmodes are solutions to Maxwell’s eigenvalue problem, ${\epsilon }^{-1}\nabla \times \nabla \times {\mathbf{E}}_n={\omega }_n^2{\mathbf{E}}_n$ [55], with outgoing radiation conditions. These solutions (also called “leaky modes” [89]) have complex frequencies ω_n which lie in the lower-half of the complex plane (Im[ω_n] < 0) [90] and, consequently, the modal amplitudes (|E_n|) grow unboundedly at large distances from the structure. Even though the modal amplitude diverges, we show that the unconjugated inner product is finite and independent of the PML parameters, as long as the PML works effectively (i.e., designed so that outgoing waves are attenuated exponentially inside the PML and are not reflected at the air-PML interface.) We provide a simple proof for a one dimensional geometry. An abstract proof for the independence of the unconjugated norm on the PML parameters appeared in [88] (using a generalized definition of the norm, which includes dispersion). An alternative one-dimensional proof was given in [91], where the authors used analytic continuation of the coordinates, similar to the coordinate stretching method that we use here.

We consider the one-dimensional geometry depicted in Fig. 5(a), which consists of a slab of thickness a and refractive index n, embedded in air. (Generalizations to more complex one-dimensional structures straightforwardly follow.) We truncate the computational cell by placing PML at a finite distance from the slab (|x| = L/2), and we impose metallic boundary conditions at the cell boundary $\left(x=\pm \frac{N}{2}\right)$. The PML can be viewed as an analytic continuation of the spatial coordinate into the complex plane [92], $\tilde{x}=x+if(x)$, where the derivative of f satisfies $\frac{df}{dx}=\frac{\sigma (x)}{\omega }$. (The latter condition guaranties that fields oscillating at different frequencies ω will be attenuated inside the PML at rate that is independent of ω.) For concreteness, we choose: $\sigma (x)=\frac{{\sigma }_{0}x}{d}$, where $d=\frac{N-L}{2}$ is the PML thickness. Consequently, the coordinate stretching transformation inside the PML is

(38)$$\tilde{x}=x+i\frac{{\sigma }_0{\left(x-\frac{L}{2}\right)}^2}{2\omega d}.$$

 

Fig. 5 Leaky modes in one dimension. (a) A slab in air with PML at the cell boundary. (b) The intensity of the leaky mode increases exponentially in the air region and is attenuated in the PML.

Download Full Size | PDF

 

Fig. 6 Forcing EPs in the periodic waveguides. (a) Real and (b) imaginary parts of the eigenvalues of Maxwell’s operator, as a function of the wavevector k_x and the permittivity contrast δε ≡ ε₂ − ε₁, for the geometry from Fig. 3(a). The EP is found by numerically minimizing the distance between the two eigenvalue sheets (red and blue surfaces). Cyan dots: Eigenvalues for fixed k_x = k_EP and varying δε. Black dots: Eigenvalues for fixed δε = δε_EP and varying k_x.

Download Full Size | PDF

Scattering solution for this geometry can be written explicitly as

$$E_n=\{\begin{array}{cc}{A}_n{e}^{-i{\omega }_{n}x-({\sigma }_0/2d){(x+L/2)}^2}+B_n{e}^{-i{\omega }_{n}x+({\sigma }_0/2d){(x+L/2)}^2}& -N/2<x<-L/2\\ A_n{e}^{-i{\omega }_{n}x}+B_n{e}^{i{\omega }_{n}x}& -L/2<x<-a/2\\ e^{-i{\omega }_{n}nx}+e^{i{\omega }_{n}nx}& -a/2<x<a/2\\ A_n{e}^{i{\omega }_{n}x}+B_n{e}^{-i{\omega }_{n}x}& a/2<x<L/2\\ A_n{e}^{-i{\omega }_{n}x-({\sigma }_0/2d){(x-L/2)}^2}+B_n{e}^{-i{\omega }_{n}x+({\sigma }_0/2d){(x-L/2)}^2}& L/2<x<N/2\end{array},$$

where the coefficients A_n, B_n and resonant frequency ω_n are found by requiring continuity of the field (E_n) and its derivative (dE_n/dx) at the interfaces, while imposing the boundary condition: E_n = 0 at ±N/2. We obtain

$$\begin{array}{l}{A}_n=\frac{e^{inka/2}+e^{-inka/2}}{e^{ika/2}-C_n{e}^{-ika/2}},\\ B_n=-A_n{C}_n,\end{array}$$

where we introduced $C_n\equiv e^{i{\omega }_{n}N-{\sigma }_{0}d}$.

Next, we compute the unconjugated norm, ${\displaystyle \int \epsilon E_n^2}$, and study its convergence. Introducing the antiderivative function

(39)$$J(x)\equiv {\displaystyle {\int }_{}^{x}d{x}^{\prime }\epsilon {\left(A_n{e}^{i{\omega }_n{x}^{\prime }}+B_n{e}^{-i{\omega }_n{x}^{\prime }}\right)}^2},$$

the unconjugated norm can be written as

(40)$$ℐ(d)\equiv {\displaystyle {\int }_{\text{all}}dx\epsilon E^2}={\displaystyle {\int }_{\text{Cav}}dx\epsilon E^2+2\left[J\left(\frac{L}{2}\right)-J\left(\frac{a}{2}\right)\right]}+2\left\{J\left[\tilde{x}(\frac{N}{2})\right]-J\left[\tilde{x}(\frac{L}{2})\right]\right\}.$$

We now show that whenever the condition Eq. (41) (below) holds, $ℐ$ is finite and independent of the PML parameters in the limit of d → ∞. Since the coordinate stretching factor is zero at the air-PML boundary $\left[\tilde{x}\left(\frac{L}{2}\right)=\frac{L}{2}\right]$, the antiderivative terms at L/2 cancel $(J\left(\frac{L}{2}\right)=J\left[\tilde{x}\left(\frac{L}{2}\right)\right])$ Next, consider the boundary term, $J\left[\tilde{x}\left(\frac{N}{2}\right)\right]$. Straightforward integration of Eq. (39) yields three terms, all of which decay exponentially as d → ∞, provided that |C_n| decays exponentially. [The first is $\frac{A_n^2{e}^{2i{\omega }_n\tilde{x}(N/2)}}{2i{\omega }_n}=\frac{A_n^2{C}_n}{2i{\omega }_n}$, the second term is $\frac{B_n^2{e}^{2i{\omega }_n\tilde{x}(N/2)}}{-2i{\omega }_n}=\frac{A_n^2{C}_n}{-2i{\omega }_n}$, and the last term is $2A_n{B}_n\tilde{x}(N/2)=-2A_n^2{C}_n\cdot \mathcal{O}(d)$.] Introducing ${\omega }_n={\omega }_n^{\prime }-i{\omega }_n^{″}$ and α = d/L, one finds from the definition of C_n that it decays exponentially whenever

(41)$$[{\omega }_n^{″}(\alpha +2)+{\sigma }_0]d\gg 1.$$

Evaluating the remaining terms in Eq. (40), ${\displaystyle {\int }_{\text{Cav}}dx\epsilon E^2}$ and $J\left(\frac{a}{2}\right)$ we obtain:

$$\underset{d\to \infty }{\lim }ℐ(d)=\frac{2}{n{\omega }_n}[\sin (n{\omega }_{n}a)+n{\omega }_{n}a]+\frac{4i}{{\omega }_n}\cos \left(\frac{n{\omega }_{n}a}{2}\right).$$

The result is finite, and independent of the location of the air-PML interface or the cell boundary, thus completing the proof.

D. Computational details on forcing EPs in periodic waveguides

In order to force an EP (i.e., an accidental degeneracy of two nearby resonances), we need to satisfy a single complex equation (ω_m = ω_n), which can be done by searching two real parameters. In the plasmonic example in Sec. 3, we search the two-dimensional parameter space spanned by the gain/loss parameter (|Im ε|) and the (real part) of the refraction index in the upper-half of the silica coating. In the periodic example from Sec. 4, we search the two-dimensional space spanned by the wavevector k_x and the permittivity contrast δε. Figure 6 presents our numerical results for the periodic example. As shown, we find a degenerate resonance at ω_EP = 0.3924377 − 0.00004119303i when k_EPa/2π ≈ 1.468 × 10⁻⁴ and δε_EP ≈ 1.137.

E. Derivation of the simplified model for the LDOS

In this section, we derive Eq. (11) from the main text, which provides a simplified formula for the LDOS in Sec. 3. Our approach is similar to coupled-mode theory, originally developed for photonic waveguides [64, 65]. Our derivation consists of projecting Maxwell’s operator [describing the coupled-rod system from Fig. 1(a)] onto the subspace spanned by the modes of the corresponding uncoupled-system (i.e., a system in which the rod–rod separation is infinite).

The field of the (original) coupled system satisfies Maxwell’s equation:

(42)$$\widehat{\mathcal{A}}\mathbf{E}\equiv \frac{1}{\epsilon }\nabla \times \nabla \times \mathbf{E}={\omega }^2\mathbf{E}.$$

Denoting by ε_ℓ the permittivity of a system that contains only a single rod [where ℓ = 1 (2) indicates the upper (lower) rod in Fig. 1(a)], the right eigenvectors of the uncoupled system satisfy

(43)$$\frac{1}{{\epsilon }_{\ell }}\nabla \times \nabla \times {\mathbf{E}}_{\ell }^R={\omega }_{\ell }^2{\mathbf{E}}_{\ell }^R,$$

with a similar definition for left eigenvectors.

Next, let us compute the projection of $\mathcal{A}$ onto the subspace spanned by E₁ and E₂, denoted by $\widehat{\mathcal{A}}$. The diagonal terms of $\widehat{\mathcal{A}}$ are

(44)$${\displaystyle \int {\mathbf{E}}_j^L\frac{1}{\epsilon }\nabla \times \nabla \times {\mathbf{E}}_j^R\approx }{\displaystyle \int {\mathbf{E}}_j^L\frac{1}{{\epsilon }_j}\nabla \times \nabla \times {\mathbf{E}}_j^R\approx {\omega }_j^2.}$$

The first approximation follows from the fact that the mode profile ${\mathbf{E}}_j^R$ is significant only in close proximity to resonator j, and in that region ε = ε_j. The second equality follows from Eq. (43). Now, let us define Δε_j ≡ ε − ε_j. The off-diagonal terms of $\widehat{\mathcal{A}}$ are

(45)$$\begin{array}{l}{\displaystyle \int {\mathbf{E}}_i^L\frac{1}{\epsilon }\nabla \times \nabla \times {\mathbf{E}}_j^R\approx }{\displaystyle \int {\mathbf{E}}_i^L\frac{{\epsilon }_j}{\epsilon }\frac{1}{{\epsilon }_j}\nabla \times \nabla \times {\mathbf{E}}_j^R=}\\ {\omega }_j^2{\displaystyle \int {\mathbf{E}}_i^L\frac{\epsilon -\mathrm{\Delta }{\epsilon }_j}{\epsilon }{\mathbf{E}}_j^R=2\omega j\left(-\frac{{\omega }_j}{2}{\displaystyle \int {\mathbf{E}}_i^L\frac{\mathrm{\Delta }{\epsilon }_j}{\epsilon }{\mathbf{E}}_j^R}\right)\equiv 2{\omega }_j{\kappa }_{ij}}\end{array}$$

Assuming that the structure is approximately symmetric in the ŷ direction (i.e., under exchanging i ↔ j), we have κ_ij ≈ κ_ji ≈ κ. In the limit of low losses/gain, γ ≪ Ω_EP, we obtain Eq. (11):

(46)$$\widehat{A}=V^T\widehat{\mathcal{A}}U=\left(\begin{array}{cc}{\left({\omega }_{\text{EP}}-i\eta \right)}^2& 2{\mathrm{\Omega }}_{\text{EP}}\kappa \\ 2{\mathrm{\Omega }}_{\text{EP}}\kappa & {\left({\omega }_{\text{EP}}+i\eta \right)}^2\end{array}\right).$$

where U is the matrix whose columns are ${\mathbf{E}}_1^R$ and ${\mathbf{E}}_2^R$ and V is the matrix whose columns are ${\mathbf{E}}_1^L$ and ${\mathbf{E}}_2^L$

F. Giant LDOS enhancement in active periodic systems

In this section, we demonstrate substantial LDOS enhancements by forcing EPs in periodic waveguides with significant amounts of gain. Similar to Sec. 3 in the text, we consider a waveguide with periodic index modulation along $\widehat{x}$ and outgoing boundary conditions in the transverse y-direction [Fig. 3(a)]. The parameters of the system, prior to adding gain, are ε₁ = 12, ε₂ ≈ 14.43, and d = 0.4807a.

Figure 7 shows our numerical results. At each wavevector, we compute the LDOS_k peak value (M_k) vs. deviation from the EP (Δk ≡ k − k_EP) when adding increasing amounts of gain to the waveguide (i.e., fixing Q_p while increasing Q_a). The LDOS is evaluated at r₀ [Fig. 3(a)] by direct inversion of Maxwell’s equations. We show here enhancements of ≈ 400, while higher values can easily be obtained; the enhancement is essentially not bounded in this computational model. However, in reality, it is bounded by quantum noise near threshold [93]. The side panels in the figure demonstrate that the LDOS peak scales quadratically with Q_a near the EP [Fig. 7(a)] and linearly with Q_a away from the EP [Fig. 7(c)]. Finally, we note that when evaluating the LDOS near the center of the computational cell (i.e., at x ≈ 0), the lineshape changes dramatically, and it actually has a minimum at the resonance frequency and two side peaks, whose amplitude scales as $Q_a^2$ (not shown).

 

Fig. 7 LDOS enhancement in active periodic waveguides. Middle: LDOS peak (M_k) vs. Δk for increasing amounts of gain in the waveguide from Fig. 3(a). Left: log₁₀(M_EP) vs. log₁₀(Q_a) for the data from the middle panel near k_EP, showing quadratic scaling with Q_a. Right: log₁₀(M_∞) vs. log₁₀(Q_a) for the data from the middle panel at k_∞, showing linear scaling with Q_a. Inset: Eigenvalues move vertically in the complex plane upon adding gain.

Download Full Size | PDF

G. Theoretical limit of LDOS enhancement at an EP

In this section, we derive Eq. (18) in the text. Let us first define an effective mode amplitude 〈E, E〉 ≡ ∫_Cdx |E|², where C denotes a finite region containing the geometry (e.g., the last scattering surface). In order to obtain an upper bound on the LDOS enhancement at the EP, we need to estimate the quantities $〈{\mathbf{E}}_0^R,{\mathbf{E}}_0^R〉/〈{\mathbf{J}}_0^R,{\mathbf{J}}_0^R〉$ and $〈{\mathbf{E}}_0^L,{\mathbf{E}}_0^L〉/〈{\mathbf{J}}_0^L,{\mathbf{J}}_0^L〉$, which determine the relative magnitude of the two terms in the expansion formula for G at the EP [Eq. (9)]. Let us decompose the complex-symmetric Maxwell’s operator into: $\widehat{\mathcal{A}}={\widehat{\mathcal{A}}}^{\prime }+i{\widehat{\mathcal{A}}}^{″}$, where ${\widehat{\mathcal{A}}}^{\prime }\equiv \frac{\widehat{\mathcal{A}}+\widehat{\mathcal{A}}*}{2}$ and ${\widehat{\mathcal{A}}}^{″}\equiv \frac{\widehat{\mathcal{A}}-\widehat{\mathcal{A}}*}{2i}$, and the asterisk denotes complex conjugation. In many cases of interest, one can assume $\Vert {\widehat{\mathcal{A}}}^{″}\Vert \ll \Vert {\widehat{\mathcal{A}}}^{\prime }\Vert$ (under an appropriate matrix norm). In such cases, one can use defective perturbation theory (Sec. 2.B) to expand the eigenmodes ${\mathbf{E}}_{\pm }^R$ and eigenvalues ω_± of ${\widehat{\mathcal{A}}}^{\prime }$ in terms ${\mathbf{E}}_0^R$, ${\mathbf{J}}_0^R$ and λ_EP. Using Eq. (6), we obtain

(47)$${\mathbf{E}}_0^R\approx ({\mathbf{E}}_{+}^R+{\mathbf{E}}_{-}^R)/2$$

(48)$${\mathbf{J}}_0^R\approx ({\mathbf{E}}_{+}^R-{\mathbf{E}}_{-}^R)/(2{\lambda }_1{p}^{\frac{1}{2}}).$$

Since ${\widehat{\mathcal{A}}}^{\prime }$ is a real Hermitian operator, it has real and orthogonal eigenvectors E_±. Using Eqs. (47) and (48) and assuming 〈E₊, E₊〉 ≈ 〈E₋, E₋〉 and 〈E₋, E₊〉 = 0, one obtains

(49)$$\begin{array}{l}〈{\mathbf{E}}_0^R,{\mathbf{E}}_0^R〉\approx 〈{\mathbf{E}}_{+},{\mathbf{E}}_{+}〉/2\\ 〈{\mathbf{J}}_0^R,{\mathbf{J}}_0^R〉\approx 〈{\mathbf{E}}_{+},{\mathbf{E}}_{+}〉/2{\left|{\lambda }_1\right|}^{2}p.\end{array}$$

Substituting the explicit perturbative expansion Eq. (7), we find

(50)$$\frac{〈{\mathbf{E}}_0^R,{\mathbf{E}}_0^R〉}{〈{\mathbf{J}}_0^R,{\mathbf{J}}_0^R〉}\approx {\left|{\lambda }_1\right|}^{2}p\approx \left|\frac{({\mathbf{E}}_0^L,{\widehat{\mathcal{A}}}^{″},{\mathbf{E}}_0^R)}{({\mathbf{J}}_0^L,{\mathbf{E}}_0^R)},\right|$$

(with an equivalent expression for the left eigenvector). Next, recall the definition of $\widehat{\mathcal{A}}\equiv \frac{1}{\epsilon }\nabla \times \nabla \times$. In the limit of low losses, ∥Im ε/ε∥ ≪ 1, we can approximate ${\widehat{\mathcal{A}}}^{″}\approx -\frac{\mathrm{Im}\epsilon }{{\epsilon }^2}\nabla \times \nabla \times$ and, consequently

(51)$$\begin{array}{l}|({\mathbf{E}}_0^L,{\widehat{\mathcal{A}}}^{″},{\mathbf{E}}_0^R)|\approx \left|{\omega }_{\text{EP}}^2{\displaystyle \int {\mathbf{E}}_0^L(\mathrm{Im}\epsilon /\epsilon ){\mathbf{E}}_0^R}\right|=\\ \left|{\omega }_{\text{EP}}^2{\displaystyle \int \mathrm{Im}\epsilon {\mathbf{E}}_0^{R^2}}\right|\approx |{\omega }_{\text{EP}}^2{\displaystyle \int \mathrm{Im}\epsilon {\left|{\mathbf{E}}_0^R\right|}^2|\le {\left|{\omega }_{\text{EP}}\right|}^2\max |\mathrm{Im}\epsilon |},\end{array}$$

where, in going from the first to the second line, we used the relation ${\mathbf{E}}_0^R={\epsilon }^{-1}{\mathbf{E}}_0^L$ (which holds for Maxwell’s eigenvalue problem) and, in the following approximations, we used the property ${\mathbf{E}}_0^R={\left|{\mathbf{E}}_0\right|}^2$ (valid for low-loss systems) and the normalization condition ∫|E₀|² = 1. This completes the proof of Eq. (18).

Funding

Army Research Office(ARO) (W911NF-13-D-0001); U.S. Department of Energy (DOE) (DE-SC0001299); National Science Foundation (NSF) (DMR-1307632, DMR-1454836); United States–Israel Binational Science Foundation (BSF) (2013508).

Acknowledgments

AP, FH, and SGJ were partially supported by the Army Research Office through the Institute for Soldier Nanotechnologies under Contract No. W911NF-13-D-0001. MS (reading and analysis of the manuscript) was supported by MIT S3TEC EFRC of the U.S. DOE under Grant No. DE-SC0001299. CWH was supported by the NSF grant No. DMR-1307632. AWR was supported by the NSF grant No. DMR-1454836. BZ was partially supported by the United States–Israel Binational Science Foundation (BSF) under award No. 2013508. We would like to thank Yi Yang, Aristeidis Karalis, Wonseok Shin, Nick Rivera, and David Liu for helpful discussions.

References and links

1. E. M. Purcell, “Spontaneous emission probabilities at radio frequencies,” Phys. Rev. 69, 681 (1946).

2. H. Yokoyama and K. Ujihara, Spontaneous Emission and Laser Oscillation in Microcavities (CRC Press, 1995), Vol. X.

3. S. V. Gaponenko, Introduction to Nanophotonics (Cambridge University, 2010). [CrossRef]  

4. W. J. Firth and A. M. Yao, “Giant excess noise and transient gain in misaligned laser cavities,” Phys. Rev. Lett. 95, 073903 (2005). [CrossRef]   [PubMed]  

5. M. V. Berry, “Mode degeneracies and the Petermann excess-noise factor for unstable lasers,” J. Mod. Opt. 50, 63–81 (2003). [CrossRef]  

6. S.-Y. Lee, J.-W. Ryu, J.-B. Shim, S.-B. Lee, S. W. Kim, and K. An, “Divergent Petermann factor of interacting resonances in a stadium-shaped microcavity,” Phys. Rev. A 78, 015805 (2008). [CrossRef]  

7. K. Petermann, “Calculated spontaneous emission factor for double-heterostructure injection lasers with gain-induced waveguiding,” IEEE J. Quant. Elect. 15, 566–570 (1979). [CrossRef]  

8. M. A. van Eijkelenborg, Å. M. Lindberg, M. S. Thijssen, and J. P. Woerdman, “Resonance of quantum noise in an unstable cavity laser,” Phys. Rev. Lett. 77, 4314 (1996). [CrossRef]   [PubMed]  

9. A. M. van der Lee, N. J. van Druten, A. L. Mieremet, M. A. van Eijkelenborg, Å. M. Lindberg, M. P. van Exter, and J. P. Woerdman, “Excess quantum noise due to nonorthogonal polarization modes,” Phys. Rev. Lett. 79, 4357 (1997). [CrossRef]  

10. A. E. Siegman, “Excess spontaneous emission in non-Hermitian optical systems. I. Laser amplifiers,” Phys. Rev. A 39, 1253–1263 (1989). [CrossRef]  

11. A. E. Siegman, “Excess spontaneous emission in non-Hermitian optical systems. II. Laser oscillators,” Phys. Rev. A 39, 1264–1268 (1989). [CrossRef]  

12. N. Moiseyev, Non-Hermitian Quantum Mechanics (Cambridge University, 2011).

13. W. D. Heiss, “The physics of exceptional points,” J. Phys. A - Math. Theo. 45, 444016 (2012). [CrossRef]  

14. Z. Lin, H. Ramezani, T. Eichelkraut, T. Kottos, H. Cao, and D. N. Christodoulides, “Unidirectional invisibility induced by PT -symmetric periodic structures,” Phys. Rev. Lett. 106, 213901 (2011). [CrossRef]  

15. B. Peng, S. K. Özdemir, F. Lei, F. Monifi, M. Gianfreda, G. L. Long, S. Fan, F. Nori, C. M. Bender, and L. Yang, “Parity-time-symmetric whispering-gallery microcavities,” Nat. Phys. 10, 394–398 (2014). [CrossRef]  

16. L. Feng, Y.-L. Xu, W. S. Fegadolli, M.-H. Lu, J. E. Oliveira, V. R. Almeida, Y.-F. Chen, and A. Scherer, “Experimental demonstration of a unidirectional reflectionless parity-time metamaterial at optical frequencies,” Nat. Mater. 12, 108–113 (2013). [CrossRef]  

17. J. Doppler, A. A. Mailybaev, J. Böhm, A. G. U. Kuhl, F. Libisch, T. J. Milburn, P. Rabl, N. Moiseyev, and S. Rotter, “Dynamically encircling an exceptional point for asymmetric mode switching,” Nature 537, 76–79 (2016). [CrossRef]   [PubMed]  

18. H. Xu, D. Mason, L. Jiang, and J. G. E. Harris, “Topological energy transfer in an optomechanical system with exceptional points,” Nature 537, 80–83 (2016). [CrossRef]   [PubMed]  

19. D. Heiss, “Mathematical physics: Circling exceptional points,” Nat. Phys. 12, 823–824 (2016). [CrossRef]  

20. L. Feng, Z. J. Wong, R.-M. Ma, Y. Wang, and X. Zhang, “Single-mode laser by parity-time symmetry breaking,” Science 346, 972–975 (2014). [CrossRef]   [PubMed]  

21. H. Hodaei, M.-A. Miri, M. Heinrich, D. N. Christodoulides, and M. Khajavikhan, “Parity-time–symmetric microring lasers,” Science 346, 975–978 (2014). [CrossRef]   [PubMed]  

22. M. Brandstetter, M. Liertzer, C. Deutsch, P. Klang, J. Schöberl, H. E. Türeci, G. Strasser, K. Unterrainer, and S. Rotter, “Reversing the pump dependence of a laser at an exceptional point,” Nat. comm. 5, 1 (2014). [CrossRef]  

23. M. Liertzer, L. Ge, A. Cerjan, A. D. Stone, H. E. Türeci, and S. Rotter, “Pump-induced exceptional points in lasers,” Phys. Rev. Lett. 108, 173901 (2012). [CrossRef]   [PubMed]  

24. B. Peng, S. K. Özdemir, S. Rotter, H. Yilmaz, M. Liertzer, F. Monifi, C. M. Bender, F. Nori, and L. Yang, “Loss-induced suppression and revival of lasing,” Science 346, 328–332 (2014). [CrossRef]   [PubMed]  

25. G. Yoo, H.-S. Sim, and H. Schomerus, “Quantum noise and mode nonorthogonality in non-Hermitian $\mathcal{P}\mathcal{T}$-symmetric optical resonators,” Phys. Rev. A 84, 063833 (2011). [CrossRef]  

26. W. D. Heiss and G. Wunner, “Fano-Feshbach resonances in two-channel scattering around exceptional points,” Eur. Phys. J. D 68, 1–6 (2014). [CrossRef]  

27. A. M. van den Brink, K. Young, and M. H. Yung, “Eigenvector expansion and Petermann factor for ohmically damped oscillators,” J. Phys A: Math. Gen. 39, 3725 (2006). [CrossRef]  

28. M. D. Zoysa, T. Asano, K. Mochizuki, A. Oskooi, T. Inoue, and S. Noda, “Conversion of broadband to narrowband thermal emission through energy recycling,” Nat. Photon. 6, 535–539 (2012). [CrossRef]  

29. X. Liu, T. Tyler, T. Starr, A. F. Starr, N. M. Jokerst, and W. J. Padilla, “Taming the blackbody with infrared metamaterials as selective thermal emitters,” Phys. Rev. Lett. 107, 045901 (2011). [CrossRef]   [PubMed]  

30. M. Laroche, R. Carminati, and J.-J. Greffet, “Coherent thermal antenna using a photonic crystal slab,” Phys. Rev. Lett. 96, 123903 (2006). [CrossRef]   [PubMed]  

31. C. Khandekar, W. Jin, O. D. Miller, A. Pick, and A. W. Rodriguez, “Giant frequency-selective near-field heat transfer in PT-symmetric structures,” Phys. Rev. B 94, 115402 (2016). [CrossRef]  

32. T. Søndergaard, S. M. Novikov, T. Holmgaard, R. L. Eriksen, J. Beermann, Z. Han, K. Pedersen, and S. I. Bozhevolnyi, “Plasmonic black gold by adiabatic nanofocusing and absorption of light in ultra-sharp convex grooves,” Nat. Comm. 3, 969 (2012). [CrossRef]  

33. W. Wan, Y. Chong, L. Ge, H. Noh, A. D. Stone, and H. Cao, “Time-reversed lasing and interferometric control of absorption,” Science 331, 889–892 (2011). [CrossRef]   [PubMed]  

34. Y. D. Chong, L. Ge, H. Cao, and A. D. Stone, “Coherent perfect absorbers: time-reversed lasers,” Phys. Rev. Lett. 105, 053901 (2010). [CrossRef]   [PubMed]  

35. Y. Sun, W. Tan, H.-Q. Li, J. Li, and H. Chen, “Experimental demonstration of a coherent perfect absorber with pt phase transition,” Phys. Rev. Lett. 112, 143903 (2014). [CrossRef]   [PubMed]  

36. H. Linnenbank, Y. Grynko, J. Förstner, and S. Linden, “Second harmonic generation spectroscopy on hybrid plasmonic/dielectric nanoantennas,” Light Sci. Appl. 5, e16013 (2016). [CrossRef]  

37. B. Zhen, C. W. Hsu, Y. Igarashi, L. Lu, I. Kaminer, A. Pick, S.-L. Chua, J. D. Joannopoulos, and M. Soljačić, “Spawning rings of exceptional points out of Dirac cones,” Nature 525, 354–358 (2015). [CrossRef]   [PubMed]  

38. T. Gao, E. Estrecho, K. Y. Bliokh, T. C. H. Liew, M. D. Fraser, S. Brodbeck, M. Kamp, C. Schneider, S. Hofling, Y. Yamamoto, F. Nori, Y. S. Kivshar, A. G. Truscott, R. G. Dall, and E. A. Ostrovskaya, “Observation of non-Hermitian degeneracies in a chaotic exciton-polariton billiard,” Nature 526, 554–558 (2015). [CrossRef]   [PubMed]  

39. A. P. Seyranian and A. A. Mailybaev, Multiparameter Stability Theory With Mechanical Applications (World Scientific Publishing, 2003), Vol. XIII. [CrossRef]  

40. D. C. Brody and E.-M. Graefe, “Information geometry of complex hamiltonians and exceptional points,” Entropy 15, 3361–3378 (2013). [CrossRef]  

41. A. Taflove, A. Oskooi, p., and S. G. Johnson, Advances in FDTD Computational Electrodynamics: Photonics and Nanotechnology (Artech House, 2013), http://arxiv.org/ftp/arxiv/papers/1301/1301.5366.pdf.

42. G. W. Hanson, A. I. Nosich, and E. M. Kartchevski, “Green’s function expansions in dyadic root functions for shielded layered waveguides,” Prog. Electromagn. Res. 39, 61–91 (2003). [CrossRef]  

43. E. Hernàndez, A. Jàuregui, and A. Mondragòn, “Jordan blocks and Gamow-Jordan eigenfunctions associated with a degeneracy of unbound states,” Phys. Rev. A 67, 022721 (2003). [CrossRef]  

44. E. Hernandez, A. Jauregui, and A. Mondragon, “Degeneracy of resonances in a double barrier potential,” J. Phys. A 33, 4507–4523 (2000). [CrossRef]  

45. J. D. Jackson, Classical Electrodynamics, Edition III (John Wiley and Sons, 1999).

46. L. N. Trefethen and M. Embree, Spectra and Pseudospectra: The Behavior of Nonnormal Matrices and Operators (Princeton University, 2005).

47. S. Longhi, “Enhanced excess noise in laser cavities with tilted mirrors,” Opt. Lett. 25, 811–813 (2000). [CrossRef]  

48. G. D’Alessandro and C. B. Laforet, “Giant noise amplification in synchronously pumped optical parametric oscillators,” Opt. Lett. 34, 614–616 (2009). [CrossRef]  

49. G. D’Alessandro and F. Papoff, “Giant subthreshold amplification in synchronously pumped optical parametric oscillators,” Phys. Rev. A 80, 023804 (2009). [CrossRef]  

50. F. Wijnands, J. B. Pendry, F. J. Garcia-Vidal, P. M. Bell, L. M. N. Moreno, and P. J. Roberts, “Green’s functions for Maxwell’s equations: Application to spontaneous emission,” Opt. Quant. Elect. 29, 199–216 (1997). [CrossRef]  

51. Y. Xu, R. K. Lee, and A. Yariv, “Quantum analysis and the classical analysis of spontaneous emission in a microcavity,” Phys. Rev. A 61, 033807 (2000). [CrossRef]  

52. K. Sacked, Optical Properties of Photonic Crystals (Springer Science & Business Media, 2004), Vol. LXXX.

53. A. Lagendijk and B. A. V. Tiggelen, “Resonant multiple scattering of light,” Phys. Rep. 270, 143–215 (1996). [CrossRef]  

54. G. B. Arfken and H. J. Weber, Mathematical Methods for Physicists (Elsevier Academic Press, 2006).

55. J. D. Joannopoulos, S. G. Johnson, J. N. Winn, and R. D. Meade, Photonic Crystals: Molding the Flow of Light (Princeton University, 2011).

56. P. M. Morse and H. Feshbach, Methods of Theoretical Physics (McGraw-Hill, 1953), Vol. I.

57. A. E. Siegman, Lasers (University Science Books, 1986).

58. H. E. Türeci, A. D. Stone, and B. Collier, “Self-consistent multimode lasing theory for complex or random lasing media,” Phys. Rev. A 74, 043822 (2006). [CrossRef]  

59. R. G. Newton, Scattering Theory of Waves and Particles (Springer Science & Business Media, 2013).

60. S. H. Weintraub, “Jordan canonical form: Application to differential equations,” Synth. Lect. Math. Statist. 1, 1–85 (2008). [CrossRef]  

61. S. A. Maier, Plasmonics: Fundamentals and Applications (Springer Science & Business Media, 2007).

62. C. E. Rüter, K. G. Makris, R. El-Ganainy, D. N. Christodoulides, M. Segev, and D. Kip, “Observation of parity–time symmetry in optics,” Nat. Phys. 6, 192–195 (2010). [CrossRef]  

63. A. Christ and H. L. Hartnagel, “Three-dimensional finite-difference method for the analysis of microwave-device embedding,” IEEE Trans. Microw. Theory Techn. 35, 688–696 (1987). [CrossRef]  

64. W. P. Huang, “Coupled-mode theory for optical waveguides: an overview,” J. Opt. Soc. Am. A 11, 963–983 (1994). [CrossRef]  

65. K. Okamoto, Fundamentals of Optical Waveguides (Academic Press, 2010).

66. N. W. Ashcroft and N. D. Mermin, Solid State Physics (Holt, Rinehart and Winston, 1976).

67. R. C. McPhedran, L. C. Botten, J. McOrist, A. A. Asatryan, C. M. de Sterke, and N. A. Nicorovici, “Density of states functions for photonic crystals,” Phys. Rev. E 69, 016609 (2004). [CrossRef]  

68. F. Capolino, D. R. Jackson, and a. B. F. D. R. Wilton, “Comparison of methods for calculating the field excited by a dipole near a 2-d periodic material,” IEEE Trans. Antennas Propag. 55, 1644–1655 (2007). [CrossRef]  

69. A. P. F. Hernàndez and S. G. Johnson, “Scalable computation of jordan chains,” arXiv:1704.05837 (2017).

70. H. Benisty, D. Labilloy, C. Weisbuch, C. J. M. Smith, T. F. Krauss, D. Cassagne, A. Beraud, and C. Jouanin, “Radiation losses of waveguide-based two-dimensional photonic crystals: Positive role of the substrate,” App. Phys. Lett. 76, 532–534 (2000). [CrossRef]  

71. S. M. Barnett and R. Loudon, “Sum rule for modified spontaneous emission rates,” Phys. Rev. Lett. 77, 2444 (1996). [CrossRef]   [PubMed]  

72. E. M. Graefe, H. J. Korsch, and A. E. Niederle, “A non-hermitian symmetric bose? hubbard model: eigenvalue rings from unfolding higher-order exceptional points,” J. Phys. A. Math. Theor. 41, 255206 (2008). [CrossRef]  

73. J.-W. Ryu, S.-Y. Lee, and S. W. Kim, “Analysis of multiple exceptional points related to three interacting eigenmodes in a non-hermitian hamiltonian,” Phys. Rev. A 85, 042101 (2012). [CrossRef]  

74. W. D. Heiss and G. Wunner, “Resonance scattering at third-order exceptional points,” J. Phys. A. Math. Theor. 48, 345203 (2015). [CrossRef]  

75. K. Ding, G. Ma, M. Xiao, Z. Q. Zhang, and C. T. Chan, “Emergence, coalescence, and topological properties of multiple exceptional points and their experimental realization,” Phys. Rev. X 6, 021007 (2016).

76. Z. Lin, A. Pick, M. Lončar, and A. W. Rodriguez, “Enhanced spontaneous emission at third-order dirac exceptional points in inverse-designed photonic crystals,” Phys. Rev. Lett. 117, 107402 (2016). [CrossRef]   [PubMed]  

77. N. Colthup, Introduction to Infrared and Raman Spectroscopy (Elsevier, 2012).

78. J. G. Fujimoto, C. Pitris, S. A. Boppart, and M. E. Brezinski, “Optical coherence tomography: An emerging technology for biomedical imaging and optical biopsy,” Neoplasia 2, 9–25 (2000). [CrossRef]   [PubMed]  

79. J. D. McKinney, M. A. Webster, K. J. Webb, and A. M. Weiner, “Characterization and imaging in optically scattering media by use of laser speckle and a variable-coherence source,” Opt. Let. 25, 4–6 (2000). [CrossRef]  

80. P. Miao, Z. Zhang, J. Sun, W. Walasik, S. Longhi, N. M. Litchinitser, and L. Feng, “Orbital angular momentum microlaser,” Science 353, 464–467 (2016). [CrossRef]   [PubMed]  

81. A. Pick, Z. Lin, W. Jin, and A. W. Rodriguez, “Enhanced frequency conversion and nonlinear purcell effect at exceptional points,” (In preparation).

82. D. R. Jackson and A. A. Oliner, Leaky-Wave Antennas (Wiley Online Library, 2008).

83. F. Monticone and A. Alù, “Leaky-wave theory, techniques, and applications: From microwaves to visible frequencies,” Proc. IEEE 103, 793–821 (2015). [CrossRef]  

84. O. D. Miller, A. G. Polimeridis, M. T. H. Reid, C. H. H. Chia, B. G. DeLacy, J. D. Joannopoulos, M. Soljačić, and S. G. Johnson, “Fundamental limits to optical response in absorptive systems,” Opt. express 24, 3329–3364 (2016). [CrossRef]   [PubMed]  

85. A. E. Siegman, “Excess quantum noise in nonnormal oscillators,” Frontiers of Laser Physics and Quantum Optics (Springer, 2000). [CrossRef]  

86. B. Vial, F. Zolla, A. Nicolet, and M. Commandré, “Quasimodal expansion of electromagnetic fields in open two-dimensional structures,” Phys. Rev. A 89, 023829 (2014). [CrossRef]  

87. M. Perrin, “Eigen-energy effects and non-orthogonality in the quasi-normal mode expansion of maxwell equations,” Opt. Express 24, 27137–27151 (2016). [CrossRef]   [PubMed]  

88. C. Sauvan, J.-P. Hugonin, I. S. Maksymov, and P. Lalanne, “Theory of the spontaneous optical emission of nanosize photonic and plasmon resonators,” Phys. Rev. Lett. 110, 237401 (2013). [CrossRef]   [PubMed]  

89. J. Hu and C. R. Menyuk, “Understanding leaky modes: slab waveguide revisited,” Adv. Opt. Photonics 1, 58–106 (2009). [CrossRef]  

90. R. E. Collin, Field theory of guided waves (McGraw-Hill, 1960).

91. P. T. Leung, S. Y. Liu, and K. Young, “Completeness and orthogonality of quasinormal modes in leaky optical cavities,” Phys. Rev. A 49, 3057 (1994). [CrossRef]   [PubMed]  

92. W. C. Chew and W. H. Weedon, “A 3d perfectly matched medium from modified maxwell’s equations with stretched coordinates,” Microw. Opt. Technol. Lett. 7, 599–604 (1994). [CrossRef]  

93. A. Pick, A. Cerjan, D. Liu, A. W. Rodriguez, A. D. Stone, Y. D. Chong, and S. G. Johnson, “Ab initio multimode linewidth theory for arbitrary inhomogeneous laser cavities,” Phys. Rev. A 91, 063806 (2015). [CrossRef]