Topological superconductivity enhanced by exceptional points

Majorana zero modes (MZMs) emerge as edge states in topological superconductors and are promising for topological quantum computation, but their detection has so far been elusive. Here we show that non-Hermiticity can be used to obtain dramatically more robust MZMs. The enhanced properties appear as a result of an extreme instability of exceptional points to superconductivity, such that even a vanishingly small superconducting order parameter already opens a large energy gap, produces well-localized MZMs, and leads to strong superconducting pair correlations. Our work thus illustrates the large potential of enhancing topological superconductivity using non-Hermitian exceptional points.

In the ideal situation, MZMs emerge at zero energy and are well separated from the quasicontinuum of other states by a topological gap, which creates an isolated subspace for the topological qubit.However, in practice, the topological gap strongly depends on material properties and is usually smaller than the superconducting order parameter in the parent superconductor.This considerably limits a realistic Majorana-based computation [2,57], which urgently calls for platforms with much larger values of the topological gap.
With the advent of non-Hermitian quantum mechanics, it has been predicted that system properties can be drastically enhanced due to the strong sensitivity of non-Hermitian systems [58][59][60][61][62].This sensitivity is due to the presence of non-Hermitian degeneracies, also known as exceptional points (EPs), which are points where both eigenvalues and eigenstates coalesce [59,.The sensitivity gets even larger when increasing the number N of coalescing energy levels, where N is also referred to as the order of the EP.In fact, for a dimensionless perturbation x, the spectrum changes as x 1/N [58][59][60][61][62], revealing that even an infinitesimally small x can induce a sizable effect for sufficiently large N .A tantalizing question arises here if it may be possible to exploit this strong sensitivity of non-Hermitian systems to enhance the topological gap protecting the MZMs.
In this work, we show that non-Hermiticity can dramatically enhance topological superconductivity, producing much more robust MZMs.In particular, we consider the non-reciprocal Kitaev chain and show that the presence of a high-order EP in the normal state makes the Lattice representation of the HNK Hamiltonian with lattice sites (red dots) and terms connecting different sites (dashed/dotted lines).
system highly unstable to superconductivity, such that even a vanishingly small superconducting order parameter is sufficient to open a sizable gap protecting the MZMs, resulting in exceptionally enhanced topological superconductivity.With sizable superconductivity being hard to achieve experimentally, our results open a compelling avenue for designing systems with much more robust MZMs.
Model and its normal state.-Toinvestigate and most clearly elucidate exceptionally enhanced topological superconductivity, we consider a non-Hermitian extension of the simplest known topological superconductor: the Kitaev chain.Further, as we want to investigate systems with strong non-Hermitian sensitivity, we seek a very simple model where the normal state presents a high-order EP.An ideal candidate for the latter is the Hatano-Nelson (HN) model [82,92,93], where simple non-reciprocal hoppings give rise to high order EP [82,92,93].Combining the two models we arrive at the Hatano-Nelson-Kitaev (HNK) chain, which is illustrated in Fig. 1 and given by where c † r (c r ) creates (annihilates) a spinless fermion at site r, t is the nearest neighbor hopping, µ the chemical potential, ∆ the p-wave superconducting order parameter, and Γ encodes the non-reciprocity.
The Hermitian, or Kitaev, limit, found at Γ = 0, hosts a topological phase for |µ/t| < 2, characterized by a bulk topological invariant and MZMs at the endpoints of any finite chain [3][4][5][6].Inclusion of the non-reciprocal Γ into the Kitaev chain has already been considered [28][29][30].Even though non-Hermitian systems can present new topological phases or extend the topological phases in Hermitian systems [82,94], the topological phase diagram for the HNK chain was found to remain, being characterized by the non-Hermitian version of the same invariant [28] and showing localized MZMs.It has additionally been shown that the non-Abelian braiding is preserved in non-Hermitian systems [69,95].These earlier results seem to indicate that one does not gain much from making the Kitaev chain non-Hermitian.However, an entirely overlooked aspect is how the normal-state EPs can dramatically enhance the topological superconducting phase, albeit the phase diagram boundaries do not change.
Before showing that the EP in the normal state spectrum enhances the superconducting phase, we need to describe this EP.The normal state of the HNK chain in just the HN model [82,92,93] or Eq. ( 1) with ∆ = 0. Due to the non-reciprocal hopping, the HN model presents remarkably different properties depending on the boundary conditions.In the presence of periodic boundary conditions (PBC), it hosts non-Hermitian topological phases but no EPs.In contrast, the open boundary conditions (OBC) system does not have a topological phase but has an EP of the order of the system size L [82,92,93].As we are interested in the EP effects and MZMs, we focus on OBC in the main text but discuss the PBC system in the Supplementary Material (SM) [96] (see also Refs.[97][98][99][100][101][102] therein).In general, for non-Hermitian systems, the bulk spectrum for OBC can be captured by the non-Bloch Hamiltonian h, which takes values in a complex generalization of the Brillouin zone [103][104][105].In many cases, including the HN model, we can then obtain an analytical expression for the eigenvalues of h and, consequently, of the OBC spectrum [103][104][105].However, for the HNK chain, to the best of our knowledge, there exists no analytical solution, see SM [96].Therefore, we are here forced to primarily rely on a numerical solution of the Boguliubov-de Gennes (BdG) Hamiltonian [106,107], where , for an open chain with L sites.The real/imaginary parts of the spectrum of Eq. (1) in the normal state are shown in Fig. 2(a)/(c) as a function of Γ/t.For |Γ/t| < 1, the energies are completely real, while for |Γ/t| > 1, the energies are purely imaginary.The presence of regions with completely real energies in the system is due to pseudo-Hermitian (PH) symmetry, while the developing of a spectrum with complex conjugated pairs marks the spontaneous breaking of this symmetry [108][109][110][111][112].These PH-preserved and PHbroken phases are separated by an EP of the order of the system size, where all energy levels coalesce [82,92,93]  [113].Exceptional gap.-Having understood the HNK system without superconductivity and its EP, we next add a small superconducting order parameter ∆.Doing so, we arrive at a surprising conclusion: already adding a vanishingly small ∆ dramatically changes the energy levels.In Fig. 2(b)/(d), we show the real/imaginary part of the spectrum for the extremely small ∆ = 10 −10 t (see SM [96] for complementary data at different L).For such small values of ∆, most of the spectrum barely changes from the non-superconducting case.The exception is for Γ-values around the normal-state EP, which disappears and is replaced by a continuum of bands, which joins pairwise in multiple EPs.In addition, there are also two zero energy states, two MZMs, clearly isolated from the rest of the other modes, easily seen when plotting the complex energies (red filled circles) in the inset of Fig. 2(d).We further notice that although these MZMs are clearly isolated modes, there still exist bulk modes with, separately, Re(ϵ) = 0 or Im(ϵ) = 0. Thus, one may naively think the system is gapless, just as in the Hermitian limit (blue circles), but thanks to the complex energy spectrum, there exists a finite region in the complex plane around the MZMs where no bulk states are found.This spectrum configuration corresponds to a point gap, which provides a non-Hermitian topological protection of the MZMs [82].As a consequence, even for a vanishingly small ∆, we find a finite minigap δ = |ϵ QP − ϵ MZM |, with ϵ QP being the lowest quasi-particle mode, protecting the MZMs.The finite minigap exists for a finite range of Γ, marked by the vertical green dashed lines in Fig. 2(b,d).
Although we find no general analytical expression for the minigap protecting the MZM, we can gain understanding by perturbatively adding superconducting pairing at the normal-state EP.The details are given in the SM [96] and here we focus on the results.At Γ/t = 1, the (bulk) energy levels ϵ m are given, to the lowest order in ∆/t, by with m = 0, 1, . . ., 2L − 2. From this result, we see that for sufficiently large L, L−1 ∆/t ∼ 1 for any ∆, which directly explains the exceptional sensitivity of the energy levels to any perturbation ∆.This result is in agreement with the recently reported sensitivity of non-Hermitian systems [66,114], but they were not considering superconductivity or other electronic ordering.This approximate energy level expression also clarifies why the energies are complex since each of them is a point on a circle in the complex plane.The radius of this circle also automatically gives the minigap, δ/t ≈ (∆/t) since the MZMs are at zero energy.In Fig. 3(a,b), we compare these perturbative and thus approximate analytical results for the energy levels (a) and minigap (b) with the earlier obtained exact numerical results using Γ/t = 1 and ∆/t = 10 −10 for many different values of L. We find that analytical (lines) and numerical (points) results agree so well with each other that they are not distinguishable in the figure.We also directly see that as L increases, the radii of the energy level circles increase, leaving a larger minigap.In particular, we notice that the (L − 1)-root expression in Eq. ( 2) means that there is a power-law behavior of the minigap δ, as we show in the SM [96].This is very different from the Hermitian limit (black points, lines) where the minigap is always δ ∼ ∆.We thus establish, both numerically and analytically, that the normal-state EP makes the system exceptionally unstable in the thermodynamic limit to superconductivity, generating a finite minigap even for vanishingly small ∆.It is the coalescence of all levels at the normal-state EP that generates this extreme instability as a collective phenomenon.
The above analytical analysis is only valid at the normal-state EP at Γ/t = 1, but the instability towards superconductivity extends for a finite range of Γ around the EP.We illustrate this in Figs.3(c,d), where we plot the minigap δ as a function of both Γ and ∆ for both L = 30 (c) and L = 100 (d).The strong system size dependence is also visible here, as a sizable minigap is present in a much larger region for the larger L system, clearly emphasizing the collective aspect of the effect.The spectrum is also affected by other values of the parameters, but MZMs are still present in a large region of the parameter space.For instance, although we here only report data for µ = 0, for ∆ = 10 −10 t the enhancement appears in the large regime |µ/t| ≲ 1.To not deviate from the main result, we discuss these details in the SM [96] and in Supplementary Videos (SVs) [115].
Taken together, we find that already a vanishingly small ∆ around the normal-state EP changes the features of the spectrum in the complex plane, opening a sizable point gap that produces isolated MZMs.We coin this remarkable sensitivity exceptionally enhanced topological superconductivity.The sensitivity is a collective effect, arising due to the accumulation of all system modes at the EP and thus growing with system size.This can be viewed as a non-Hermitian equivalent of how the accumulation of density of states in flat band systems creates superconductivity, as well as other electronic ordered states, as recently illustrated in e.g.twisted bilayer graphene at the magic angle [116][117][118][119] or bi-and trilayer graphene in electric fields [120,121].
Robust MZMs and superconductivity.-Exceptionaltopological superconductivity is not restricted to influencing the topological gap, but next we show how it also manifest in the spectral weight and superconducting pair correlations.These are encoded in the Nambu Green's function G [88,101,122,123] [124] where G e and G h are the normal electron and hole Green's functions, respectively, while F eh and F he are the anomalous contributions containing the superconducting pair correlations or amplitudes.From G, we obtain all two-particle observables of the system, see SM [96] for details.First, we calculate the spectral weight , which represents the density of states at the frequency (energy) ω (ℏω).
Note that Tr here denotes the sum over both electron and hole components and spatial coordinates.In Figs.4(a,d), we plot A for the HNK system.At ∆ = 0 (a), A is only finite when PH-symmetry is preserved, and we find mainly non-zero real energy values.The only visible peak occurs at the normal-state EP (|Γ/t| = 1) due to its accumulation of modes.For a vanishingly small ∆ = 10 −10 t (d), the same behavior, with A reflecting the non-zero real part of the spectrum, is now only found in a narrower region of Γ.Instead, there is a huge peak at ω = 0 that extends until the real energies disappear.This is the same region, marked by green dashed lines in Fig. 2, where we find a point gap separating the MZMs from the remaining modes.Thus, this peak in A is nothing else than a zero bias peak (ZBP) produced by the MZMs [5,6].Moreover, since the sizable point gap is associated with finite imaginary energies, which broadens the poles for all modes but the MZMs, we find a very clear ZBP with no other low-energy spectral weight, again signaling the exceptional enhancement of the topological phase.The exceptional enhancement is not only seen in the spectrum and spectral weight but also in the spatial properties of the system.To illustrate this, we consider the local density of states (LDOS), ρ(ω; r) = −Im [G e (ω; r, r)] /π.In Fig. 4(b), we plot ρ as a function of position along the chain for ∆ = 10 −10 t in the Hermitian limit (Γ = 0) and find, as expected, essentially uniform occupation across energies and position.In contrast, tuning to the normal-state EP at Γ = t in Fig. 4(c), we find very pronounced and well-localized peaks at zero energy at each end of the chain.This demonstrates that non-Hermiticity also strongly promotes MZM localization, which, together with the large minigap, establishes the existence of robust MZMs.Farther away from the EP, we still find MZMs, but they then decay further into the bulk, see SV [115], since the coherence length of the MZMs is proportional to δ −1 [5,6].
Finally, we also show that the superconducting pair correlations, expressed through F eh , are also exceptionally enhanced due to the normal-state EP.Keeping a vanishingly small ∆ = 10 −10 t, we find in the Hermitian limit in Fig. 4(c) that the local pair correlations are of the same order as ∆ and delocalized in the lattice, as expected.In contrast, for Γ = t in Fig. 4(f), we find an excessively strong F , up to 20 orders of magnitude larger than ∆.Such a strong response shows that, although the superconducting order parameter ∆ is vanishingly small, there exists strong superconducting pairing, further cementing the notion of exceptionally enhanced topological superconductivity.Another remarkable feature is the localization of F eh at the right edge due to the Non-Hermitian skin effect [73,[103][104][105][125][126][127][128][129][130][131][132][133][134][135][136].If instead plotting F he , we find a similar localization at the left edge, illustrating that electron and hole components localize on opposite edges.
Concluding remarks.-Weshow that even a vanishingly small superconducting order parameter ∆ in a nonreciprocal Kitaev chain can create exceptional topological superconductivity with robust MZMs, due to a normalstate EP of the order of the system size.The robustness of the MZMs stems from the opening of an exceptionally enhanced point gap in the system, which both generates the minigap protecting the MZMs and strongly localized MZMs.Moreover, the superconducting pair correlations are also exceptionally enhanced.While the most dramatic behavior occurs for vanishingly small ∆ close to the normal-state EP, any non-Hermiticity still enhances MZM robustness and superconductivity.Our results are also stable to changes in the chemical potential.We expect that this relative insensitivity to parameter values also makes the results stable to weak disorder.In terms of a direct experimental realization, cold atom systems are promising as they can host both p-wave pairing [137,138] and non-reciprocity [80,93].Our results are also transferable to other superconductors as they only require the opening of a point gap close to an EP in the normal state.Additionally, we speculate that the same reasoning can be applied to other correlated phases that are amplified by state degeneracies, simple examples being the Stoner mechanism for ferromagnetism [139] and the Peierls instability [140].Finally, an interesting outlook is to understand the role of particle-hole symmetry in exceptionally enhanced superconductivity, since the enhancement of the gap with a power 1/(L − 1), while preserving MZMs, seems to be related to the analysis of perturbations of EPs in the presence of chiral and sublattice symmetry in Ref. [141].
We are grateful for discussions with E. J.In the main text, we show that for open boundary conditions (OBC), the presence of a high order exceptional point (EP) in the normal state Hatano-Nelson-Kitaev (HNK) chain (i.e.just the Hatano-Nelson (HN) model) leads to an extreme, or exceptional, sensitivity to superconductivity.Here we show that due to the absence of such an EP for periodic boundary conditions (PBC), there exists no such enhancement of superconductivity in the periodic version of the HNK chain.With PBC, the Hamiltonian is block-diagonal in momentum space and reads such that we can write the BdG Hamiltonian as From here, we obtain the energy spectrum directly The real part of the spectrum is equal to the energy spectrum of the (Hermitian) Kitaev chain, and the imaginary part is just linear in Γ.The Green's function, in this case, is also easy to obtain and is given by −2i∆ sin(ka) 2i∆ sin(ka) ℏω + 2iΓ sin(ka) + µ + 2t cos(ka) (7) where B(ω, k) = ℏ 2 ω 2 + 2iℏωΓ sin(ka) − 4(∆ 2 + Γ 2 ) sin 2 (ka) − [µ + 2t cos(ka)] 2 .We note that in this case, for ∆ small, F eh/he ∝ ∆, directly illustrating that the superconducting amplification seen for OBC in the main text is always absent for PBC.
The fact that the periodic system is barely modified by a small ∆ is also clearly seen when we consider how the spectrum and the spectral weight change with ∆, as illustrated in Fig. 5.This figure is to be compared with Fig. 2 of the main text, which shows the same quantities with similar parameters, but with OBC.For ∆ = 0, we seem to have a continuum of bands in both the real, Fig. 5(a), and imaginary, Fig. 5(d), part of the spectrum.However, the analysis of the spectrum in the complex plane, inset of Fig. 5(d), shows that the system has a point gap, but now without any Majorana zero modes (MZMs) in the middle, as expected for a periodic system.This point gap does not reflect any interesting properties in the spectral weight, see Fig. 5(e), aside from the van Hove singularities at ℏω = ±2t, related to the borders of the Brillouin zone, and the spectrum of the PBC Kitaev chain at Γ = 0.The inclusion of a vanishingly small ∆ = 10 −10 t does not alter the spectrum, see Figs. 5 (b,e) or the spectral weight, see Fig. 5(h).For a sufficiently large ∆ = 0.1t, we find an opening of a small gap in the real part of the energy in Fig. 5(c), but with an unchanged imaginary energy part, see Fig. 5(f), which leads to line gap of a similar size as ∆ for all Γ, see inset of Fig. 5(f).
To summarize, the abrupt change of the electronic properties, even for a vanishingly small value of the superconducting order parameter ∆ observed in the main text (see Figs. 2, 3, and 4 of the main text) for OBC is clearly not present for PBC.This is manifest both in the superconducting pair correlations, obtained analytically, and the numerical analysis of the spectrum and the spectral weight.The fact that the HNK system with PBC is not exceptionally sensitive to superconducting pairing can be traced back to the lack of EPs in the spectrum for the normal state (HN model).

NON-BLOCH FORMALISM
In this section of the SM we discuss the non-Bloch formalism and motivate why it is difficult to obtain analytical results for the HNK chain.Some general properties of this model, and of non-Hermitian superconductors in general, are explained in detail in Ref. [28].Here, we reproduce some of these results for pedagogical reasons.In non-Hermitian systems, there is a breaking of the traditional bulk-boundary correspondence, such that it cannot be used.Nevertheless, there exist different approaches to obtaining the spectrum and the states, such as the transfer matrix [97,108], and the non-Bloch formalism [103][104][105].In particular, the non-Bloch formalism makes that the exp(ik) of the Fourier transform is substituted by a complex number z.Using it, the non-Bloch Hamiltonian of the HNK chain is of the form [28] The secular equation det ϵ 1 − h(z) = 0 then becomes a polynomial equation in ϵ and z.In our case Therefore, the energies are given by while the secular equation is a fourth-order polynomial equation for z (11) The four solutions for z can be sorted according to their absolute value.The physical solutions are given by the solutions z 2 and z 3 , which satisfy |z 2 | = |z 3 |.In some cases, like the HN model, |z| does not depend on ϵ.In this case, we can write z = exp [ik(ϵ)]|z| and obtain an effective Bloch Hamiltonian called a surrogate Hamiltonian.However, for the HNK chain, |z| is a function of ϵ, |z|(ϵ), such that the OBC modes can only be obtained numerically.Even though there exist no analytical solutions for the bands, z(ϵ) = z −1 (−ϵ) due to particle-hole symmetry [28], leading to electron and holes presenting a non-Hermitian skin effect (NHSE) in different edges of the system.The fact that |z|(ϵ) explains why we in the main text have to resort primarily to a numerical approach in our description of exceptionally enhanced superconductivity, although we can still derive key properties of the energy levels and minigap analytically.

FINITE SIZE EFFECTS
In this section of the SM we extend the analysis presented in Figs.2(b,d) of the main text to show how the spectrum at vanishingly small ∆ = 10 −10 t changes for different sizes.In particular, Fig. 6 shows the spectrum and spectral weight for the OBC HNK system with ∆ = 10 −10 t for increasing values of L from left to right using L = 10, 20, 30, 40, 50.We clearly see how the exceptional sensitivity of the non-superconducting HN system to the introduction of superconductivity heavily depends on system size as the normal-state EP becomes more and more degenerate with an increasing number of lattice sites.As such it is a collective effect due to the accumulation of all the system modes at the EP in the normal state.In Supplementary Video (SV) 1, we continuously vary L to even more clearly illustrate this effect.The exceptional sensitivity of the non-superconducting HN system to the introduction of superconductivity is a collective effect due to the accumulation of all the system modes at an EP in the normal state.Therefore, it is a phenomenon that depends heavily on system size as the EP becomes more and more degenerate with an increasing number of lattice sites.

Lattice with an odd number of sites
In both the main text and in the rest of this SM, we consider lattices with an even number of lattice sites because systems with an odd number of lattice sites can have spurious zero modes that are not MZMs.This is due to the combination of PH-symmetry, which constraints the energies to appear in (ϵ, ϵ * ) pairs, and particle-hole symmetry, which causes the energies to appear in (ϵ, −ϵ) pairs, which make lattices with an odd number of sites to always host zero energy modes.These modes can be trivial, if degenerate with spatially overlapping eigenstates, or isolated MZMs.In Fig. 7, we show a version of Fig. 6, i.e. we display the spectrum and spectral weight of the system for increasing system sizes, but for lattices with an odd number of sites, still keeping ∆ = 10 −10 t.Overall, Figs. 6 and 7 are quite similar.But for the odd-numbered lattices, we find an additional ZBP for all Γ/t, which clearly illustrates the spurious, non-MZM, zero-energy modes in odd-numbered systems.We also note that there are no new EPs points induced by superconductivity in contrast to the even-numbered systems.However, the normal-state EP of the HN model is still present, which is the root of the exceptionally enhanced superconductivity, as it is independent of evenor odd-numbered lattices.In Fig. 8 we additionally show the minigap δ, local density of states ρ, and superconducting pair correlations F for a lattice with L = 31 sites and ∆ = 10 −10 t.This figure should be directly compared with Figs. 3 and 4 of the main text that displays the same quantities and for the same parameters but for a lattice with L = 30 sites.The behavior of δ as a function of ∆ and Γ, Fig. 8(a), is the same as the one shown in Fig. 3(c), as is the behavior of ρ, Fig. 8(c) and F , Fig. 8(e) in the PH-broken phase.This illustrates that exceptionally enhanced topological superconductivity is present in all systems, independent of them having an even or odd number of sites.However, in the Hermitian limit, Figs.8(b,d), the same quantities now show additional peaks at ω = 0 due to the spurious zero modes.We note, however, that these modes are not localized at the edges but are fully delocalized across the chain, and thus per definition not spatially isolated MZMs.Interestingly, they are modes with a wavelength of two lattice parameters, which causes the dotted patterns seen in Figs.8(b,d).The analysis present in this SM shows how the system size dramatically changes the superconducting gap and spectral properties of the system due to the large degeneracy at the normal-state EP, which further supports the analysis present in the main text.

APPROXIMATE ANALYTICAL EXPRESSION FOR SPECTRUM AT Γ = t
In this section of the SM, we demonstrate how we obtain the analytical results in Eq. ( 2) of the main text and discuss this result in more detail.Since the HNK Hamiltonian is not a Toeplitz matrix, a generic analytic expression for the eigenvalues and eigenstates is difficult to obtain in general.However, it turns out that we can still perform an approximate analytical analysis at Γ = t.
For ∆ = 0, H BdG is a direct sum of the electron and hole Hamiltonians H BdG = (H e ⊕ H h ) /2.For Γ = t, both Hamiltonians assume a Jordan form typical of EPs.As a matter of fact, the secular equation shows that all the 2L energies are equal to zero, explicitly demonstrating the very high degeneracy of the EP.The degeneracy is proportional to the number of sites L and has been named an 'infernal point' in Ref. [98].
The inclusion of even a vanishingly small value of a superconducting order parameter ∆ combines the electron and hole components of the Hamiltonian.Explicitly, H BdG that was block diagonal above, now has the block form where ∆ is the matrix of p-wave pairing in real space Besides connecting the electron and hole components, ∆ makes such that H BdG is no longer a Jordan block.Although the exact expression of the determinant cannot, to the best of our knowledge, be obtained analytically, the expression of the secular equation can be expressed in terms of a series for small ∆, such that we can obtain a scaling relationship for the energies.We describe this procedure in detail later, but here we first comment on the result and compare it with our numerical data.
The secular equation can be approximated by which has as solutions either the MZMs with ϵ = 0 or the bulk states where m = 0, 1, . . ., 2L − 3 labels a different 2L − 2-th root of (−1) L .These energies then describe points on a circle of radius r = L−1 ∆/t in the complex energy plane.This is the approximate analytical result reported in the main text in Eq. (2).Although this result is similar to the response of an EP to perturbations [66,114], it here depends on L − 1 instead of L since the modes related to the MZMs remain at zero energy.Due to the spectrum changing as L−1 ∆/t, for L → ∞, any ∆ creates a sizable value of a gap, which is the essence of exceptionally enhanced superconductivity.In Fig. 9 we compare the exact numerical results with these approximate analytical expressions for the bulk energies and also the extracted minigap to complement the results reported in the main text in Fig. 3(a,b).In Fig. 9(a), we compare Eq. ( 16) with the numerical data for ∆ = 10 −10 t, for many different system sizes L. We notice that for all system sizes, there are ϵ = 0 solutions corresponding to the MZMs, while the bulk modes form circles in the complex plane with a radius very well described by the analytical expression in Eq. ( 16).There is a very minor eccentricity of these circles, but it is reduced for the larger system size, indicating that the subleading corrections to Eq. ( 16) are only having any effect at small system sizes.The radius of the circle directly gives the minigap δ as it is the energy difference between the MZM (at zero energy) and the lowest-lying bulk modes.In Fig. 9(b) we plot the minigap δ as a function of ∆ on a log-log plot for the same systems sizes.The dots are the values obtained by the numerical diagonalization, while dashed lines are the function L−1 ∆/t.We find that the analytical expression reproduces extremely well the power-law behavior of the numerical data for not only small ∆, but also larger ∆ for large system sizes.Again, this means subleading corrections are not relevant for larger system sizes.

Form of the determinant for different system sizes
The expression of the secular determinant is complicated and, as far as we know, does not have a simple solution.Still, we can gain insight by keeping only terms of order ∆ 2 .The resulting expression is different for even and odd numbers of sites, involving different powers of ϵ, but, nevertheless, we can understand why the spectrum can be approximated by Eq. ( 16).
The determinant of a matrix can generally be obtained using the method of cofactors.Although this method does not generate an equation relating the determinants for different systems sizes, in contrast to what happens to Toeplitz matrices [99], we can still use this method to extract the determinant of the system for all small system sizes and from there clearly motivate the expression for any L. We will first perform a change of basis in H BdG equivalent to changing the Nambu vector to T , such that every two rows/columns describe the electron and holes states in each site.This change of basis is made to make the calculation of the determinant more convenient since the matrix becomes more sparse.For convenience, we also introduce the notation Below we explicitly extract D L for L = 2, 3 followed by the expressions for large even and odd L.

sites
We start with a system of just 2 sites.In this case, the determinant is a 4 × 4 matrix that can be split into 3 × 3 determinants that are easily computable: Thus the energies are given by ϵ 2 = (∆/2) 2 , such that ϵ = ±∆/2.Notice that for this small system, the modifications of the energy levels are only linear in ∆.

sites
For three sites, the determinant is needed of a 6 × 6 matrix, which makes its calculation cumbersome.Using the method of cofactors again, it can be split into two 5 × 5 determinants Again, using cofactors to split these into determinants of 4 × 4 matrices, we find the first term given by while the second is given by Notice that D 3 is not simply connected to D 2 , which indeed appears in the first term of Eq.( 20).This is a general feature of D L , which illustrates the increasing complexity.Nevertheless, for three sites, we can collect all terms and re-express the determinant as from which we obtain the energy levels ϵ = 0 and ϵ = ± ∆ 2 /2 ± it∆.
Notice that the second solution changes as √ ∆ for small ∆.This solution can be obtained by approximating The same rationale can be used to obtain approximate solutions for any L, which show the scaling of the bulk modes with the system size.For L > 3, there is no simple analytical solution, but this type of expansion to the lowest order in ∆ motivates the scaling observed in the numerical data.
Thermodynamic limit, even L Moving beyond the analytically exact expressions for the smallest system sizes, we find that for an even L, the determinant D L takes to the lowest order in ∆ the form (checked explicitly in Mathematica up to 100 sites) The term multiplying ∆ 2 in parenthesis can be re-expressed as (26) In the thermodynamic limit (L → ∞), we can use the expression for the series to approximate the secular equation by This expression is valid in the thermodynamic limit L → ∞ up to order O (∆/t) 2 .
Thermodynamic limit, odd L For an odd L, we instead get a different expression for the determinant in lowest order in ∆ (checked explicitly in Mathematica up to 99 sites) Again, we can re-express Notice that this expression does not include the term proportional to ϵ 2L−6 .In the thermodynamic limit, we can use the expression of the series to finally arrive at the equation This expression is valid in the thermodynamic limit L → ∞ up to order O (∆/t) 2 .

Scaling and MZMs
From the results above, we see that to the lowest order of ∆ we arrive at complicated polynomial equations, which in the thermodynamic limit can be approximated by the equations, ( 32) and ( 28), for even and odd L, respectively.We notice that these two equations have two important similarities.First, both present two solutions with ϵ = 0.This is exactly the MZMs.Second, the scaling of ϵ with system size can be obtained by taking the log of Eqs. ( 28) and (32).In both cases, assuming ϵ to be of form exp(iθ)∆ x , we obtain where y = (2L − 10)x or y = (2L − 8)x for odd or even L, respectively, and m describes a branch of the log of (−1) L .We thus have approximately the solutions with m restricted to be between 0 and 2L − 3 to give the different 2L − 2 bulk solutions of the Hamiltonian.We here remark that the values of the exponents and the phases obtained in this analysis can also be obtained by just disregarding the terms in the right-hand side of Eq. ( 33) in the secular equation, which leads to Eq. ( 16).
To summarize, the analytical approximation explicitly shows how the enhancement of superconductivity occurs due to the large degeneracy at the normal-state EP, such that, for a large enough system, even a vanishingly small value of ∆ is sufficient to open a gap of the order of t, as discussed in section Exceptional gap of the main text.

FINITE CHEMICAL POTENTIAL
Most results in the main text are reported at a chemical potential µ/t = 0.The inclusion of a finite chemical potential µ/t eventually destroys the topological phase in the Hermitian Kitaev chain.In particular, for |µ/t| > 2, the system is a topologically trivial superconductor, regardless of the values of ∆ and Γ, as shown by topological invariants in Ref. [28].For finite µ, but still within the topological phase, there are also some changes.We illustrate this behavior in Fig. 10 by plotting the spectrum and spectral weight for ∆ = 10 −10 t with finite µ, increasing from the left to right.We see that for µ = 0.5t, Figs.10(a,f,k), the real and imaginary part of the spectrum changes in comparison with the µ = 0 case in Fig. 2(b, d) and Fig. 4(d) of the main text, but the system still present a point gap with isolated MZMs, which are also clearly visible in the spectral weight as a ZBP.This illustrates explicitly that all our results in the main text are stable in a substantial parameter regime centred around µ = 0.Even for µ = t, Figs.10(b,g,l), we still have MZMs for a small Γ region around the normal-state EP, but the structure of the energies in the complex plane now starts to deviate considerably in the PH-broken phase.There is now also an apparent breaking of particle-hole symmetry, seen in the different structures of modes with negative and positive real parts of the energy.It is not until µ = 1.5t,Figs.10(c,h,m), that there no longer exist values of Γ where a ZBP occurs in the spectral weight.Thus, we find evidence of exceptionally enhanced superconductivity for values of the chemical potential all the way up to |µ/t| ≈ t, although the regime of Γ-values for which this enhancement occurs decreases for increasing µ.At µ = 2t, the point of the Hermitian topological phase transition, Figs.10(d,i,n) shows how the system is now gapless only in the Hermitian limit (Γ = 0).The spectral weight also only contains signatures of the spectrum of the system in the PH-preserved phase.Finally, at µ = 2.5t, Figs.10(e,j,o), there is a full energy gap with no midgap states.In SV4, we continuously vary µ for ∆ = 10 −10 t to more clearly show the effect of µ in the spectrum and spectral weight.Finally, for completeness, we consider the effect of µ on the spectrum and spectral weight for a more sizable ∆.
In Fig. 11, we show the spectrum and spectral weight for increasing µ for the large value of ∆ = t, where also the Hermitian Kitaev chain has clearly defined MZMs.We see then that for all |µ/t| ̸ = 2, the system presents a line gap, see insets of Figs.11(f-j) and hosts MZMs, while for |µ/t| > 2 there are no midgap states in the spectral function.This is precisely the phase diagram of the HNK system, also obtained in Ref. [28].In SV5, we continuously change µ for ∆ = t to show its effect more clearly.In summary, the exceptionally enhanced superconductivity due to the normal-state EP is slowly diminishing with increasing µ, but still persists all the way up to |µ/t| ≈ 1.This shows that our results in the main text are not finetuned to µ = 0.With increasing ∆, we find MZMs also for larger µ/t, but this is then probing the straightforward non-Hermitian extension of the Hermitian Kitaev chain.
where ρ e/h is the electron/hole local density of states.Notice that in the main text, we used ρ = ρ e due to the convention that the local density of states is related to the occupation of electrons only.We see that the local density of states shows the occupation of the electronic levels for a given frequency.Using the above equations and noting that for non-Hermitian systems the resolution of the identity is written using the right |m⟩ where we use again that we have a complex conjugated spectrum.We note that this is similar to the definition of the local density of states in a Hermitian system, but with the difference that ψ R ̸ = ψ L in general.
For the superconducting pair amplitudes, we simply use the anomalous components of G, resulting in Note that this expression does not become a δ function, but it is still highly peaked on states with real energies.As a final remark, we notice that the MZMs present distinguished peaks in all the above quantities since they are fully real.We also note that F presents an asymmetry between left and right edges due to the fact that the electron component of ψ R and the hole component of ψ L are localized on the same edge due to the non-Hermitian skin effect.

FIG. 1 .
FIG. 1.Lattice representation of the HNK Hamiltonian with lattice sites (red dots) and terms connecting different sites (dashed/dotted lines).

FIG. 2 .
FIG. 2. Spectrum of the OBC HNK system as a function of Γ/t with ∆ = 0 (a,c) and ∆ = 10 −10 t (b,d).First (second) row shows the real (imaginary) part of the energy spectrum.Insets in (c,d) show the energy levels in the complex plane in the Hermitian limit (Γ = 0, blue dashed lines and open circles) and in the PH-broken phase (Γ = 1.2t, red dashed lines and filled circles).Green dashed lines mark the region with MZMs.Other parameters: µ = 0, L = 30.

FIG. 5 .
FIG.5.Spectrum and spectral weight of the PBC HNK system as a function of Γ/t (same as Fig.2in the main text but for PBC).(a,d,g) have ∆ = 0, (b, e, f )∆ = 10 −10 t, and (c,f,j) ∆ = 0.1t, where the first (second) row shows the real (imaginary) part of the spectrum, while the last row shows the log of the spectral weight |A|.Insets of (d,e,f) show the energy levels in the complex plane for the system in the Hermitian limit (Γ = 0, blue dashed lines and open circles) and in the PH-broken phase (Γ = 1.2t, red dashed lines and filled circles).Other parameters: 100 k-points (corresponding to a lattice with L = 100 sites) and µ = 0.

FIG. 6 .
FIG. 6. Spectrum and spectral weight of the OBC HNK system as a function of Γ/t for L = 10 (first column, (a,f,k)), L = 20 (second column, (b,g,l)), L = 30 (third column, (c,h,m)), L = 40 (fourth column, (d,i,n)), and L = 50 (fifth column, (e,j,o)).The first (second) row shows the real (imaginary) part of the spectrum, while the third row shows the log of the spectral weight |A|.Insets in (f-j) show the energy levels in the complex plane for the system in the Hermitian limit (Γ = 0, blue dashed lines and open circles) and in the PH-broken phase (Γ = 1.2t, red dashed lines and filled circles).Other parameters: ∆ = 10 −10 t and µ = 0.

FIG. 7 .
FIG. 7. Spectrum and spectral weight of the OBC HNK system as a function of Γ/t for L = 11 (first column, (a,f,k)), L = 21 (second column, (b,g,l)), L = 31 (third column, (c,h,m)), L = 41 (fourth column, (d, i,n)), and L = 51 (fifth column, (e,j,o)).The first (second) row shows the real (imaginary) part of the spectrum, while the third row shows the log of the spectral weight |A|.Insets in (f-j) show the energy levels in the complex plane for the system in the Hermitian limit (Γ = 0, blue dashed lines and open circles) and in the PH-broken phase (Γ = 1.2t, red dashed lines and filled circles).Other parameters: ∆ = 10 −10 t and µ = 0.

FIG. 9 .
FIG. 9. Comparison between the approximate analytical expressions for (a) bulk energies ϵ and (b) minigap δ and the exact numerical data for many different system sizes L. (a) shows ϵ in the complex plane using ∆ = 10 −10 t.(b) shows the scaling of δ with ∆.Symbols indicate numerical data, dashed lines represent the analytical expressions.Other parameters: Γ = t, µ = 0.

FIG. 11 .
FIG.11.Spectrum and spectral weight of the OBC HNK system as a function of Γ/t for µ = 0 (first column, (a,f,k)), µ = t (second column, (b,g,l)), µ = 1.5t (third column, (c,h,m)), µ = 2.0t (fourth column, (d,i,n)), and µ = 2.5t (fifth column, (e,j,o)).First (second) row shows the real (imaginary) part of the spectrum, while third row (k-o) shows the log of the spectral weight |A|.Insets in (f-j) show the energy levels in the complex plane for the system in the Hermitian limit (Γ = 0, blue dashed lines and open circles) and in the PH-broken phase (Γ = 1.2t, red dashed lines and filled circles).Parameters: L = 30 and ∆ = t.
Bergholtz, J. C. Budich, and T. Yoshida.We acknowledge financial support from the Knut and Alice Wallenberg Foundation and the European Research Council (ERC) under the European Union Horizon 2020 research and innovation programme (ERC-2017-StG-757553).J. C. acknowledges financial support from the Swedish Research Council (Vetenskapsrådet Grant No. 2021-04121).
R. Arouca, Jorge Cayao, and Annica M. Black-SchafferDepartment of Physics and Astronomy, Uppsala University, Uppsala, SwedenPERIODIC BOUNDARY CONDITIONS AND LACK OF EXCEPTIONAL POINTS