Non-Hermitian systems and topology: A transfer-matrix perspective

TL;DR

This paper introduces a comprehensive transfer-matrix framework to analyze non-Hermitian topological systems, elucidating phenomena like the skin effect and boundary modes through analytical methods and topological invariants.

Contribution

It provides a unifying analytical approach using transfer matrices to understand non-Hermitian topological phases and their boundary phenomena, connecting spectral properties with topological invariants.

Findings

Transfer matrix determinant indicates non-Hermitian skin effect.

Real-space exceptional points scale with system size.

Topological invariants are derived for non-Hermitian systems.

Abstract

Topological phases of Hermitian systems are known to exhibit intriguing properties such as the presence of robust boundary states and the famed bulk-boundary correspondence. These features can change drastically for their non-Hermitian generalizations, as exemplified by a general breakdown of bulk-boundary correspondence and a localization of all states at the boundary, termed the non-Hermitian skin effect. In this article, we present a completely analytical unifying framework for studying these systems using generalized transfer matrices -- a real-space approach suitable for systems with periodic as well as open boundary conditions. We show that various qualitative properties of these systems can be easily deduced from the transfer matrix. For instance, the connection between the breakdown of the conventional bulk-boundary correspondence and the existence of a non-Hermitian skin effect, previously observed numerically, is traced back to the transfer matrix having a determinant not equal to unity. The vanishing of this determinant signals real-space exceptional points, whose order scales with the system size. We also derive previously proposed topological invariants such as the biorthogonal polarization and the Chern number computed on a complexified Brillouin zone. Finally, we define an invariant for and thereby clarify the meaning of topologically protected boundary modes for non-Hermitian systems.

Figures

Equations

$$\mathcal{H}$$

$$+c_{n+\ell ,\alpha }^{†}[{\mathbf{t}}_{\mathrm{R},\ell }^{†}{]}_{\alpha \beta }{c}_{n,\beta }^{})+c_{n,\alpha }^{†}[{\mathbf{t}}_0^{}{]}_{\alpha \beta }{c}_{n,\beta }^{}].$$

$$\mathcal{H}({\mathbf{k}}_{\perp })$$

$$|\Psi ⟩=\sum _{n=0}^N{\Psi }_n{\mathbf{c}}_n^{†}|\Omega ⟩,$$

$$J_{\mathrm{L}}^{}{\Psi }_{n+1}+M{\Psi }_n+J_{\mathrm{R}}^{†}{\Psi }_{n-1}=\epsilon {\Psi }_n.$$

$$J_{\mathrm{R}}=J_{\mathrm{L}}=J,J^2=0,$$

$$J{\Psi }_{n+1}+M{\Psi }_n+J^{†}{\Psi }_{n-1}=\epsilon {\Psi }_n,$$

$${\mathcal{H}}_{\mathrm{B}}(\mathbf{k})=J({\mathbf{k}}_{\perp }){\mathrm{e}}^{i{k}_x}+M({\mathbf{k}}_{\perp })+J^{†}({\mathbf{k}}_{\perp }){\mathrm{e}}^{-i{k}_x}.$$

$${\Psi }_n=\mathcal{G}J{\Psi }_{n+1}+\mathcal{G}{J}^{†}{\Psi }_{n-1},$$

$$J=V\Xi W^{†},$$

$$V^{†}V=W^{†}W={\text{openone}}_r,V^{†}W=0,$$

$${\mathbf{\alpha }}_n=V^{†}{\Psi }_n,{\mathbf{\beta }}_n=W^{†}{\Psi }_n,$$

$${\Psi }_n=\mathcal{G}V\Xi {\mathbf{\beta }}_{n+1}+\mathcal{G}W\Xi {\mathbf{\alpha }}_{n-1}.$$

$${\mathbf{\alpha }}_n$$

$${\mathbf{\beta }}_n$$

$$\Phi_{n+1}=T\Phi_{n},\qquad\Phi_{n}\equiv\left(\begin{array}[]{c}{\boldsymbol{\beta}}_{n}\\ {\boldsymbol{\alpha}}_{n-1}\end{array}\right),$$

$$T=\left(\begin{array}{cc}{\Xi }^{-1}\cdot {\mathcal{G}}_{vw}^{-1}& -{\Xi }^{-1}\cdot {\mathcal{G}}_{vw}^{-1}\cdot {\mathcal{G}}_{ww}\cdot \Xi \\ {\mathcal{G}}_{vv}\cdot {\mathcal{G}}_{vw}^{-1}& \left({\mathcal{G}}_{wv}-{\mathcal{G}}_{vv}\cdot {\mathcal{G}}_{vw}^{-1}\cdot {\mathcal{G}}_{ww}\right)\cdot \Xi \end{array}\right).$$

$${\Phi }_n=T^n{\Phi }_0,\forall n\in \mathbb{Z},$$

$$\det T=\det ({\mathcal{G}}_{vw}^{-1}{\mathcal{G}}_{wv})=\frac{\det {\mathcal{G}}_{wv}}{\det {\mathcal{G}}_{vw}}.$$

$$\tilde{\mathcal{H}}={\mathcal{H}}^{†}⟹\tilde{\mathcal{G}}(\epsilon )={\mathcal{G}}^{†}({\epsilon }^{\ast })$$

$${\mathcal{G}}^{†}(\epsilon )={\left[({\epsilon }^{\ast }\text{openone}-M{)}^{-1}\right]}^{†}={\left(\epsilon \text{openone}-M^{†}\right)}^{-1}=\mathcal{G}(\epsilon ),$$

$$\det T=\frac{\det {\mathcal{G}}_{wv}^{}(\epsilon )}{\det {\mathcal{G}}_{wv}^{†}({\epsilon }^{\ast })}=\frac{\det {\mathcal{G}}_{wv}^{}(\epsilon )}{{\left[\det {\mathcal{G}}_{wv}({\epsilon }^{\ast })\right]}^{\ast }}.$$

$$J=\mathcal{U}{J}^T{\mathcal{U}}^{†},M=\mathcal{U}{M}^{\ast }{\mathcal{U}}^{†}.$$

$${\mathcal{G}}^{\ast }(\epsilon )={\left(\epsilon \text{openone}-{\mathcal{U}}^{†}M\mathcal{U}\right)}^{-1}={\mathcal{U}}^{†}\mathcal{G}(\epsilon )\mathcal{U}.$$

$$V=\mathcal{U}{W}^{\ast },W=\mathcal{U}{V}^{\ast },$$

$$J=V\Xi W^{†}=\mathcal{U}{W}^{\ast }\Xi V^T{\mathcal{U}}^{†}=\mathcal{U}{J}^T{\mathcal{U}}^{†}$$

$${\mathcal{G}}_{vw}^{\ast }(\epsilon )$$

$$=V^{†}{\mathcal{U}}^T{\mathcal{U}}^{†}\mathcal{G}(\epsilon )\mathcal{U}{\mathcal{U}}^{\ast }W={\mathcal{G}}_{wv}(\epsilon ),$$

$$\det T=\frac{\det {\mathcal{G}}_{wv}^{}(\epsilon )}{\det {\mathcal{G}}_{wv}^{\ast }(\epsilon )}=\frac{\det {\mathcal{G}}_{wv}^{}(\epsilon )}{{\left[\det {\mathcal{G}}_{wv}({\epsilon }^{\ast })\right]}^{\ast }},$$

$$\Xi =\text{diag}\left\{{\xi }_1{\text{openone}}_2,{\xi }_2{\text{openone}}_2,\dots ,{\xi }_{r/2}{\text{openone}}_2\right\}.$$

$$\Sigma \equiv \text{diag}\left\{\mathcal{J},\dots ,\mathcal{J}\right\},\mathcal{J}=\left(\begin{array}{cc}0& 1\\ -1& 0\end{array}\right).$$

$$V=\mathcal{U}{W}^{\ast }\Sigma ,W=\mathcal{U}{V}^{\ast }\Sigma ,$$

Full text

I Introduction

A fundamental tenant of quantum mechanics is the reality of the spectra of operators that describe observables, which is typically achieved by demanding these operators to be Hermitian. Abandoning Hermiticity, however, has proved useful in constructing effective descriptions of dissipative systems Rotter (2009); Rotter and Bird (2015), where non-Hermitian operators encode the interactions with the environment, so that the imaginary part of their spectra can be assigned physical meaning. For instance, the imaginary part of the “energy” can be interpreted as the inverse lifetime of a (quasi-)particle Gamow (1928); Landau and Lifshitz (1981).

The study of non-Hermitian systems has primarily been driven by experiments in photonics Rüter et al. (2010); Peng et al. (2016, 2014); Lin et al. (2011); Regensburger et al. (2012); Feng et al. (2014); Hodaei et al. (2014); Chen et al. (2017); Hodaei et al. (2017); Gao et al. (2015); Zhen et al. (2015); Weimann et al. (2017); Zeuner et al. (2015); Schomerus (2013); Poli et al. (2015), where non-Hermiticity can be realized by judiciously incorporating gain and loss El-Ganainy et al. (2018); Lu et al. (2014); Malzard et al. (2015). These setups thus provide concrete realizations of non-Hermitian lattice models, such as a photonic analog of the Su-Schrieffer-Heeger model with topologically protected mid-gap states Weimann et al. (2017); Zeuner et al. (2015); Schomerus (2013). Furthermore, non-Hermitian photonic systems can be engineered to operate at “exceptional points” at which they exhibit intriguing phenomena such as unidirectional transmission Peng et al. (2016, 2014), one-sided invisibility Lin et al. (2011); Regensburger et al. (2012), single-mode lasing Feng et al. (2014); Hodaei et al. (2014), and enhanced sensitivity to perturbations Chen et al. (2017); Hodaei et al. (2017). Similar realizations of non-Hermitian models are also possible in other experimental setups, such as mechanical Bender et al. (2013), acoustic Fleury et al. (2015); Shi et al. (2016), electronic Schindler et al. (2011), and ultracold atomic Kreibich et al. (2014) systems.

Theoretically, non-Hermitian Hamiltonians have been used to describe condensed matter systems such as Majorana fermions in topological superconductors San-Jose et al. (2016); Avila et al. (2018), finite lifetime quasiparticles in heavy-fermion systems Kozii and Fu (2017); Yoshida et al. (2018), and bosonic superconductors Lieu (2018a), as well as to simulate the out-of-equilibrium systems described by a Lindblad master equation Daley (2014); Plenio and Knight (1998). In addition, certain symmetries of non-Hermitian Hamiltonians, such as a parity-time ($\mathrm{PT}$) symmetry Bender (2005); *bender_pt_symm_rev; Bender and Boettcher (1998); *bender_pt_symm; Jones-Smith and Mathur (2010) or a more general pseudo-Hermiticity Mostafazadeh (2002); *mostafazadeh_pseudo_hermitian_qm, ensure the reality of its spectrum. These Hamiltonians have garnered significant interest in mathematical physics as an analytic continuation of quantum mechanics to the complex plane Bender et al. (2002). Recently, a new direction of research has been established by investigating these systems from the perspective of topological phases Alvarez et al. (2018); Hu and Hughes (2011); Yuce (2015a); Harter et al. (2016); Yuce (2015b); Xiong (2018); Kunst et al. (2018a); Budich et al. (2019); Carlström and Bergholtz (2018); Xu et al. (2017); Yao and Wang (2018); Gong et al. (2018); Shen et al. (2018); Esaki et al. (2011); Lee (2016); Yin et al. (2018); Kawabata et al. (2018a); Yao and Wang (2018); Yao et al. (2018); Herviou et al. (2018); Kawabata et al. (2018b).

Noninteracting topological phases of matter have been of much theoretical Bernevig (2013); Shen (2013); Hasan and Kane (2010); Qi and Zhang (2011); Bernevig et al. (2006); Kane and Mele (2005); Qi et al. (2008) and experimental König et al. (2007); Hsieh et al. (2008); Xu et al. (2014) interest over the last decade. Lacking a local order parameter, these phases are characterized by features that stay unchanged under continuous deformations, such as a quantized bulk topological invariant and the appearance of robust states on their boundaries. A particularly profound result is the bulk-boundary correspondence, which establishes a direct link between the bulk invariants and the presence of robust boundary modes. Mathematically, this correspondence can be thought of as a relation between the continuous spectrum and the point spectrum of the system with open boundary conditions.

Non-Hermitian analogs of topological phases often exhibit drastically different physics from their Hermitian counterparts. An example is the existence of the aforementioned exceptional points (EPs) or more general exceptional structures Carlström and Bergholtz (2018); Budich et al. (2019); Xu et al. (2017), where a spectral degeneracy is accompanied by a coalescence of the corresponding eigenstates Kato (1995); Heiss (2012). Another remarkable feature is the possibility of a marked difference between the spectra of systems for periodic and open boundary conditions (hereafter PBC and OBC, respectively), in stark contrast to the Hermitian systems. This is accompanied by a piling up of “bulk” states at the boundaries for a finite system, a phenomenon termed the non-Hermitian skin effect Xiong (2018); Kunst et al. (2018a); Yao and Wang (2018). For gapped systems with robust boundary modes, this difference in spectra as well as nature of states for PBC and OBC signify a breakdown of the bulk-boundary correspondence.

The problem of restoring the bulk-boundary correspondence by defining a bulk invariant that can predict the existence of topologically protected boundary modes is rather subtle Gong et al. (2018); Shen et al. (2018). Various topological invariants have been defined using a generalization of the conventional Berry connection by replacing the standard inner product with a biorthogonal one Lee (2016); Esaki et al. (2011); Yin et al. (2018); Shen et al. (2018); Kawabata et al. (2018a), but they often fail to predict the existence of robust boundary modes Kawabata et al. (2018a); Herviou et al. (2018); Yao and Wang (2018); Yao et al. (2018). This failure can be traced back to the computation of the topological invariant using the continuous spectrum for a periodic system as opposed to that for a system with open boundaries, a distinction that does not exist for Hermitian systems. Indeed, bulk invariants computed taking this into account for certain specific models Yao and Wang (2018); Yao et al. (2018); Yao and Wang (2018); Lee and Thomale (2019) have been shown to correctly predict the existence of the topological boundary modes. A correct understanding of the bulk-boundary correspondence for non-Hermitian systems thus necessitates an understanding of the eigenstates of a finite system with OBC.

In this paper, we present just such an approach by constructing generalized transfer matrices Hatsugai (1993); Lee and Joannopoulos (1981); Tauber and Delplace (2015); Dwivedi and Chua (2016) for quasi one-dimensional non-Hermitian tight-binding models. We show that various qualitative features of these models can be readily gleaned off from the determinant of the transfer matrix without any numerical exact diagonalizations. For a given tight-binding model, we can thus directly answer questions such as:

• •

Is there a difference between the PBC and OBC spectra? If yes, is it always accompanied by the non-Hermitian skin effect?

• •

Where do the EPs occur in a finite system? Are they at the same parameter values as those in a periodic system?

• •

How does one define a bulk topological invariant that predicts the existence of robust boundary modes?

The transfer matrix approach also facilitates analytic computation of the full spectrum and wave functions for arbitrary finite non-Hermitian systems with OBC, quantities which so far have only been accessible by numerical computations. The implementation of the symmetries of the tight-binding model on the transfer matrix provides a new lens to view systems with additional symmetry, which can be used to explain similarities between, for instance, Hermitian and $\mathrm{PT}$-symmetric systems.

More concretely, for the transfer matrix $T$, we show that a necessary and sufficient condition for the equality of the bulk spectra for PBC and OBC is $\left|\det T\right|=1$. We prove that this condition is always satisfied for Hermitian Hamiltonians as well as $\mathrm{PT}$-symmetric Hamiltonians in the $\mathrm{PT}$-unbroken phase, thereby explaining the observed qualitative similarities in their behavior. The “bulk” states for OBC are shown to vary as $\left|\Psi_{n}\right|\sim\left|\det T\right|^{n/2}$, so that they are localized at the left/right boundary of the system for $\left|\det T\right|\lessgtr 1$. Thus, the difference between the PBC and OBC spectra is always accompanied by the non-Hermitian skin effect. Finally, for $\det T\to 0,\infty$, the propagation using the transfer matrix becomes unidirectional, which corresponds to a “real space EP” whose order scales with the system size, as previously observed in numerical computations Martinez Alvarez et al. (2018).

The correct bulk topological invariant can then be computed by using these decaying states for OBC, which corresponds to deforming the Brillouin zone by adding a complex part to the momentum, as derived for particular cases from ad hoc methods in Refs. Yao and Wang (2018); Yao et al. (2018); Yao and Wang (2018). A more geometric picture follows from the algebraic nature of the construction, which is used to construct a Riemann surface associated with the complex energy Hatsugai (1993). We show that the deformed Brillouin zone used to compute the bulk invariant is then associated with one set of noncontractible loops, while the other set of noncontractible loops are associated with the boundary modes. We thus get a winding number associated with the boundary modes and thereby clarify the meaning of a “topologically protected” boundary mode for a non-Hermitian system, where the notion of a gap may be ill defined. Finally, our formalism extends the real-space biorthogonal polarization Kunst et al. (2018a) to more general lattice topologies than those considered in Refs. Kunst et al. (2017, 2018a, 2018b); *flore_hoti_two.

Interestingly, the distinction between the PBC and OBC bulk spectra turns out to be quite transparent for cases where the transfer matrix $T(\varepsilon)$ is $2\times 2$ with eigenvalues $\rho_{\pm}(\varepsilon)$. We show that the PBC bulk bands contain energies $\varepsilon\in\mathbb{C}$ for which $\left|\rho_{a}(\varepsilon)\right|=1$ for one of $a=\pm$, the other naturally satisfying $\left|\rho_{a}(\varepsilon)\right|=\left|\det T\right|$. On the other hand, the OBC bulk bands contain energies where $\left|\rho_{a}(\varepsilon)\right|=\sqrt{\left|\det T\right|}$ for both $a=\pm$. The curves in the complex $\varepsilon$ plane for these two conditions are generically different for $\left|\det T\right|\neq 1$, as we depict schematically in Fig. 1. This explains the difference between the bulk spectra for PBC and OBC for generic non-Hermitian Hamiltonians. On the other hand, we have $\left|\det T\right|=1$ for Hermitian as well as $\mathrm{PT}$-symmetric systems, so that the two curves merge in this case, leading to identical bulk spectra for PBC and OBC.

The rest of this paper is organized as follows: In Sec. II, we construct the generalized transfer matrix for non-Hermitian tight-binding models and obtain several general results for the spectra. These results are further specialized to a particularly analytically tractable case of $2\times 2$ transfer matrices in Sec. III, and an associated energy Riemann surface is constructed. A set of explicit examples illustrating the previously derived general results are presented in Sec. IV. We finally conclude and place this work in a broader context in Sec. V. Various nonessential details of the calculations are relegated to the appendices.

Notation: We denote the set of $n\times n$ real or complex matrices as $\mathrm{Mat}(n,\mathbb{R})$ and $\mathrm{Mat}(n,\mathbb{C})$, respectively. We denote the spectrum of a matrix $M$ by $\text{Spec}\!\left[M\right]$.

II Transfer matrices

Transfer matrices arise naturally in discrete calculus as a representation of recursion relations. Since tight-binding models are essentially composed of hopping, i.e., shift, operators acting on the wave functions, the Schrödinger equation for a one-dimensional system can be reduced to a set of recursion relations, which can then be recast into a transfer-matrix equation Hatsugai (1993); Dwivedi and Chua (2016). Thus, starting with a $d$-dimensional system, we impose PBC along $(d-1)$ directions and OBC along the remaining direction, along which the transfer matrix is computed. Choosing the direction of OBC, we can analytically explore the boundary states for various possible boundaries.

II.1 General setup

Consider a system in $d$ spatial dimensions with OBC along $x$ and PBC along the remaining $(d-1)$ directions, which are parametrized by the transverse quasi-momentum $\mathbf{k}_{\perp}\in\mathbb{T}^{d-1}$. This system can also be interpreted as a family of one-dimensional chains parametrized by $\mathbf{k}_{\perp}$. Explicitly, we consider a system described by a general tight-binding non-Hermitian Hamiltonian

[TABLE]

Here, the matrices $\mathbf{t}_{\mathrm{L},\ell}$ ($\mathbf{t}_{\mathrm{R},\ell}$) denote the hopping to the left (right) and $\mathbf{t}_{0}$ is the intra-unit-cell term. For Hermitian systems, these matrices satisfy $\mathbf{t}_{\mathrm{L},\ell}=\mathbf{t}_{\mathrm{R},\ell}$ and $\mathbf{t}_{0}^{\dagger}=\mathbf{t}_{0}$. The hopping depends only on the distance owing to translation invariance and $R<\infty$ is the range of hopping. We have $q$ internal degrees of freedom, e.g., spin, orbital, or sublattice, per unit cell. The explicit dependence on $\mathbf{k}_{\perp}$ is suppressed to avoid notational clutter, however, all parameters should be assumed to depend on $\mathbf{k}_{\perp}$, unless stated otherwise.

We reduce this Hamiltonian to a nearest-neighbor form Lee and Joannopoulos (1981) by bundling together $\mathfrak{n}\geq qR$ degrees of freedoms into a supercell, whose creation (annihilation) operators are denoted by $\mathbf{c}^{\dagger}$ ($\mathbf{c}^{\phantom{\dagger}}\!\!$). This definition is not unique, and one may indeed choose arbitrarily large supercells with nearest-neighbor hopping. The Hamiltonian reduces to

[TABLE]

with the hopping matrices $J_{\mathrm{L},\mathrm{R}}$ and the on-site matrix $M$, where the latter encodes the hopping between degrees of freedom within the supercell as well as the on-site energies. An arbitrary single-particle state is given by

[TABLE]

with $\ket{\Omega}$ the fermionic vacuum state and $\Psi_{n}\in\mathbb{C}^{\mathfrak{n}}$ the wave function for each supercell. The single-particle Schrödinger equation $\mathcal{H}\ket{\Psi}=\varepsilon\ket{\Psi}$ thus reduces to the recursion relation

[TABLE]

We seek to express this as a transfer-matrix equation for cases where $J_{\mathrm{L},\mathrm{R}}$ may be singular.

In this paper, we take $M$ to be arbitrary, possibly non-Hermitian, while we demand that the hopping matrices satisfy

[TABLE]

For a Hermitian system, $M^{\dagger}=M$ and $J_{\mathrm{R}}=J_{\mathrm{L}}$, so that we have lifted the Hermiticity condition on the on-site matrix $M$ but not the hopping matrix $J$. The nilpotence of $J$ implies that no sublattice site within the supercell has hoppings to both left and right adjacent supercells. This can always be ensured by choosing a large enough supercell (also see Ref. (Dwivedi and Chua, 2016, Appendix B)). Mathematically, we need this condition to ensure that the singular vectors of $J$ [Eq. (10)] form an orthonormal set. Under these assumptions, the recursion relation becomes

[TABLE]

which corresponds to the (full) Bloch Hamiltonian

[TABLE]

In practice, we simply use this equation to identify $M$ and $J$ as the coefficients of $\mathrm{e}^{ik_{x}}$ and 1, respectively, to compute the transfer matrix for propagation along $x$.

II.2 Constructing the transfer matrix

We construct the generalized transfer-matrix representation of the recursion relation in Eq. (6) following Ref. (Dwivedi and Chua, 2016, Sec. II), which we briefly describe here. The recursion relation can be rewritten as

[TABLE]

where $\mathcal{G}=(\varepsilon\openone-M)^{-1}$ is the on-site Green’s function, which is nonsingular except when $\varepsilon$ is an eigenvalue of $M$. Next, we compute a reduced singular value decomposition (SVD) (Strang, 2009, Sec. 6.3)

[TABLE]

where $\Xi=\text{diag}\{\xi_{1},\dots,\xi_{r}\}$ with $r=\text{rank}\,{J}$ and the singular values $\xi_{i}$ are real and positive. Physically, this signifies that a suitable unitary transform of the Hamiltonian reduces it to a form where the consecutive supercells have exactly $r$ hoppings, with the magnitude of the corresponding hopping strengths given by the singular values $\xi_{i}$’s. The $r$ corresponding left and right singular vectors are assembled in the $\mathfrak{n}\times r$ matrices $V$ and $W$, which satisfy

[TABLE]

where the orthogonality of $V$ and $W$ follows from $J^{2}=0$, which also ensures that $r\leq\mathfrak{n}/2$.

As the vectors in $V$ and $W$ form an orthonormal set, they can be extended 111 We refer to Ref. (Dwivedi and Chua, 2016, Sec. II.B) for details. Also note that this breaks down if $J_{\mathrm{L}}\neq J_{\mathrm{R}}$.

to an orthonormal basis of $\mathbb{C}^{\mathfrak{n}}\ni\Psi_{n}$. We then define the coefficients of $\Psi_{n}$ in this basis as

[TABLE]

in terms of which Eq. (8) becomes

[TABLE]

Multiplying to the left by $V^{\dagger}$ and $W^{\dagger}$, we find

[TABLE]

where we have defined $\mathcal{G}_{AB}=B^{\dagger}\,\mathcal{G}\,A\in\mathrm{Mat}(r,\mathbb{C})$ with $A,B\in\{V,W\}$. This system of equations can be rewritten as

[TABLE]

where the $2r$-dimensional transfer matrix is given by

[TABLE]

The rank of $J$, and hence the size of the transfer matrix is independent of the choice of a supercell (Dwivedi and Chua, 2016, Appendix B).

Given $\Phi_{0}$, we can propagate it with the transfer matrix $T$ as

[TABLE]

provided $T$ is invertible, i.e., $\det T\neq 0$. We explicitly compute

[TABLE]

A distinct possibility for non-Hermitian systems is $\left|\det T\right|\to 0,\infty$ when $\left|\det\mathcal{G}_{wv}\right|\to 0$ and $\left|\det\mathcal{G}_{vw}\right|\to 0$, respectively. Note that these two cases are dual to each other, since if $\left|\det T\right|\to\infty$ for some parameters, we can compute the transfer matrix for translation in the opposite direction, whose determinant would then tend to zero. Physically, this corresponds to unidirectionality in the system, since it means that the states can be propagated only in one direction.

The construction above computes the transfer matrix for a right eigenstate. We can perform a similar construction of a transfer matrix for the left eigenstates by considering the action of $\mathcal{H}$ on an arbitrary single particle bra state $\bra{\Psi}$, instead of the ket in Eq. (3). Alternatively, we note that the left eigenvectors of $\mathcal{H}$ are related to the right eigenvectors of $\mathcal{H}^{\dagger}$ by a conjugate transpose. Thus, we can repeat the computation above with a new Bloch Hamiltonian

[TABLE]

to get the transfer matrix for the left eigenstates of $\mathcal{H}$.

II.3 Special cases

The transfer matrix possesses additional structure if the original Hamiltonian is Hermitian or $\mathrm{PT}$-symmetric, as we now show.

II.3.1 Hermitian systems

For Hermitian systems, the Bloch Hamiltonian satisfies $\mathcal{H}_{\mathrm{B}}^{\dagger}(\mathbf{k})=\mathcal{H}_{\mathrm{B}}^{\phantom{\dagger}}(\mathbf{k})$. For the Bloch Hamiltonian defined in Eq. (7), this implies that $M^{\dagger}=M$ with no additional condition on $J$. We compute $\mathcal{G}^{\dagger}(\varepsilon)\equiv\left[\mathcal{G}(\varepsilon^{\ast})\right]^{\dagger}$ as

[TABLE]

so that $\mathcal{G}_{AB}^{\dagger}(\varepsilon^{\ast})=\mathcal{G}_{BA}^{\phantom{\dagger}}(\varepsilon)$ and Eq. (17) reduces to

[TABLE]

Thus, for $\varepsilon\in\mathbb{R}$, i.e., the regime of physically relevant energies for Hermitian systems, $\det T=\exp[2i\arg\mathcal{G}_{wv}(\varepsilon)]$ lies on the unit circle. As expected, this reproduces Ref. (Dwivedi and Chua, 2016, Eq. (26)).

II.3.2 $\mathrm{PT}$-symmetric systems

$\mathrm{PT}$-symmetry is implemented as $\mathcal{PT}=\mathcal{U}K$ with $\;\mathcal{U}\in\mathrm{U}(\mathfrak{n})$ and $K$ the complex conjugation, so that a $\mathrm{PT}$-symmetric system satisfies $\mathcal{U}\,\mathcal{H}_{\mathrm{B}}^{\ast}(\mathbf{k})\,\mathcal{U}^{\dagger}=\mathcal{H}_{\mathrm{B}}^{\phantom{\dagger}}(\mathbf{k})$. Imposing this on the Bloch Hamiltonian in Eq. (7), we find

[TABLE]

Using the condition on the on-site matrix, we can compute $\mathcal{G}^{\ast}(\varepsilon)\equiv\left[\mathcal{G}(\varepsilon^{\ast})\right]^{\ast}$ as

[TABLE]

We next derive a condition on the singular vectors $V$ and $W$ that satisfy the condition on $J$. We here need to distinguish the two cases corresponding to $\left(\mathrm{PT}\right)^{2}=\pm 1$, which are discussed in Appendix A.

$\left(\mathrm{PT}\right)^{2}=+1$:

In this case, $\mathcal{U}=\mathcal{U}^{T}$ and in Appendix B, we show that $V,W$ must satisfy

[TABLE]

which is consistent, since $\mathcal{U}\mathcal{U}^{\ast}=\mathcal{U}\mathcal{U}^{\dagger}=\openone$. Furthermore,

[TABLE]

as desired. We can now compute

[TABLE]

so that Eq. (17) reduces to

[TABLE]

which, as in the Hermitian case, lies on the unit circle for $\varepsilon\in\mathbb{R}$, i.e., in the $\mathrm{PT}$-unbroken phase.

$\left(\mathrm{PT}\right)^{2}=-1$:

In this case, $\mathcal{U}^{T}=-\mathcal{U}$ is even dimensional, as shown in Appendix A. Alternatively, this must be the case since $\mathcal{U}\in\mathcal{U}(\mathfrak{n})\implies\left|\det\mathcal{U}\right|=1$, while the determinant vanishes for any odd-dimensional antisymmetric matrix. As we show in Appendix B, the singular values of $J$ must also come in doubly degenerate pairs, so that $\text{rank}\,J$, i.e., the number of nonzero singular values of $J$, is even, and we can write

[TABLE]

We now define

[TABLE]

Here, $\Sigma$ is antisymmetric and satisfies $\Sigma^{2}=-\openone$ and $[\Sigma,\Xi]=0$, the latter being the case because $\Xi$ is proportional to the identity matrix in each $2\times 2$ block. In Appendix B, we show that $V,W$ satisfy

[TABLE]

which is consistent since

[TABLE]

Finally, we can compute

[TABLE]

Thus,

[TABLE]

since $\det\left[-\Sigma^{2}\right]=\det\openone=1$. As in Hermitian case, $\det T$ lies on the unit circle for $\varepsilon\in\mathbb{R}$.

In conclusion, the presence of either Hermiticity or a $\mathrm{PT}$-symmetry implies the unimodularity of the transfer matrix. This is the precise sense in which the two systems behave in a similar fashion. Other symmetries of non-Hermitian Hamiltonians may also lead to this similarity with Hermitian systems, e.g., for parity-particle-hole (CP) symmetry which takes $\mathcal{H}_{\mathrm{B}}(\mathbf{k})\to-\mathcal{U}\,\mathcal{H}_{\mathrm{B}}^{\ast}(\mathbf{k})\,\mathcal{U}^{\dagger}$, we find

[TABLE]

so that $T$ is unimodular if $\varepsilon\in i\mathbb{R}$.

II.4 Spectra and states

The spectrum of the transfer matrix for a given $(\varepsilon,\mathbf{k}_{\perp})$ contains information about the possible states for that specific energy $\varepsilon$. This can also be thought of as a discrete scattering problem, where for an incoming “plane wave” of a given energy, the spectrum of the transfer matrix contains information about the fate of that plane wave as it propagates through the system. The eigenstates of the systems can then be thought of as the stationary or standing-wave solutions. Given a boundary condition, the task then is to find the values $(\varepsilon,\mathbf{k}_{\perp})$ that are compatible with such stationary solutions.

For condensed-matter systems, the most common boundary conditions to consider are periodic (PBC) and open (OBC) ones. In the following, we start with a ring with $N$ supercells realizing PBC and consider an interpolation between these two cases by tuning the strength of one of the bonds continuously to zero.

II.4.1 Periodic boundary condition

For a periodic system with $N$ supercells, $\Psi_{n}=\Psi_{n+N}$, so that using Eq. (14), we must have

[TABLE]

Thus, the system with PBC has a state for a given $(\varepsilon,\mathbf{k}_{\perp})$ iff $1\in\text{Spec}\!\left[T^{N}(\varepsilon,\mathbf{k}_{\perp})\right]$, which reduces to

[TABLE]

for some $\ell\in\{0,\dots,N-1\}$. As $N\to\infty$, the set of these points is dense on the unit circle. Thus, the bulk band for a given $\mathbf{k}_{\perp}$ is the closed, compact set of $\mathbb{C}\ni\varepsilon$ for which at least one eigenvalue $\rho$ of $T(\varepsilon,\mathbf{k}_{\perp})$ lies on the unit circle. Setting $\rho=\mathrm{e}^{ik_{x}}$ and $\Phi_{0}=\boldsymbol{\varphi}$ as the corresponding eigenvector (or one of the eigenvectors, if the corresponding eigenspace is degenerate), we write

[TABLE]

which is simply Bloch’s theorem for periodic systems.

We next set the hopping matrix connecting $\Psi_{1}$ and $\Psi_{N}\equiv\Psi_{0}$ as $\kappa J$ for some $\kappa\in\mathbb{R}$. Then, we may interpolate continuously between PBC and OBC by tuning $\kappa$ from one to zero. Following the approach of Ref. (Dwivedi, 2018, Sec. III.C.2), we write the modified recursion relation in Eq. (6) for $n=0,1$ as

[TABLE]

Multiplying to the left with $V^{\dagger}$ and $W^{\dagger}$ as earlier, these reduce for $\kappa\neq 0$ to

[TABLE]

respectively, where $K_{\mathrm{L}}=\text{diag}\left\{\openone_{r},\kappa\openone_{r}\right\}$ and $K_{\mathrm{R}}=\text{diag}\left\{\frac{1}{\kappa}\openone_{r},\openone_{r}\right\}$. Using $\Phi_{N}=T^{N-2}\Phi_{2}$, we get

[TABLE]

We finally set $\boldsymbol{\varphi}=K_{\mathrm{L}}\Phi_{1}$ to obtain

[TABLE]

Thus, we have a state iff $1\in\text{Spec}\!\left[KT^{N}(\varepsilon,\mathbf{k}_{\perp})\right]$. For $\kappa=1$, i.e., $K=\openone_{2r}$, we recover Eq. (26), which can be reduced to a condition on the spectrum of $T$ as opposed to that of $T^{N}$, and can thus be readily generalized to $N\to\infty$. This is convenient, since $T^{N}$ is generally difficult to compute analytically. For arbitrary $\kappa$, we have been able to obtain the $N\to\infty$ limit only when $r=1$ using an explicit form of $T^{N}\!\!$, as described in Sec. III.1.

II.4.2 Open boundary condition

For OBC, we need to take the limit $\kappa\to 0$, for which Eq. (30) is singular. To remedy this, we multiply to the left by $K_{\mathrm{R}}^{-1}$ to get

[TABLE]

which is well-behaved as $\kappa\to 0$. Setting $\kappa=0$, we find

[TABLE]

where ${\boldsymbol{\alpha}}_{N}$ and ${\boldsymbol{\beta}}_{1}$ are arbitrary. This is equivalent to the Dirichlet boundary condition used in Ref. (Dwivedi and Chua, 2016, Sec. II.D.3), where one starts with an infinite chain and sets $\Psi_{0}=\Psi_{N+1}=0$.

To solve this condition for $(\varepsilon,\mathbf{k}_{\perp})$, the general strategy is to find solutions to the eigenvalue problem

[TABLE]

and to then expand $\Phi_{1}$ and $\Phi_{N+1}$ in terms of these eigenvectors. We first consider the case where $T$ is diagonalizable, so that $\boldsymbol{\varphi}_{\ell}$ form a (generically non-orthogonal) basis of $\mathbb{C}^{2r}$. The condition in Eq. (33) then becomes

[TABLE]

This can be further reduced by projecting down to the sectors where the left hand side of these equations vanishes. Explicitly,

[TABLE]

where the projectors $\mathcal{P}_{{\boldsymbol{\alpha}},{\boldsymbol{\beta}}}\colon\mathbb{C}^{2r}\to\mathbb{C}^{r}$ are defined as $\mathcal{P}_{\boldsymbol{\alpha}}=\left(0,\openone_{r}\right)$ and $\mathcal{P}_{\boldsymbol{\beta}}=\left(\openone_{r},0\right)$. This is a set of $2r$ complex homogeneous linear equations in $2r$ variables $\mathbf{a}=\{a_{1},a_{2},\dots,a_{2r}\}$, which can be recast into a matrix equation of the form $\mathcal{R}\cdot\mathbf{a}=\mathbf{0}$, which, by Cramer’s rule, has a nontrivial solution iff

[TABLE]

where we have defined

[TABLE]

Since $\mathcal{R}$ is defined only in terms of the eigenvalues and eigenvectors of $T$, we obtain a condition for states that satisfy OBC purely in terms of $(\varepsilon,\mathbf{k}_{\perp})$, which can be solved to get the set of energies for which the system with OBC has an eigenstate.

On the other hand, if $T$ is non-diagonalizable or defective, we need to augment the set of eigenvectors with the generalized eigenvectors to form a basis of $\mathbb{C}^{2r}$, which can then be used to expand $\Phi_{0}$. However, the action of the transfer matrix on these eigenvalues is more complicated than in the previous case, so that the associated conditions take the general form

[TABLE]

where $f_{\ell\ell^{\prime}}(N)$ are products of polynomials and exponentials in $N$. In the case of $T$ diagonalizable, these reduce to $f_{\ell\ell^{\prime}}=\rho_{\ell}^{N}\delta_{\ell\ell^{\prime}}$.

In the following, we elucidate this idea for a simple case. Recall that if $\rho\in\text{Spec}\!\left[T\right]$ is a doubly degenerate eigenvalue with a single eigenvector $\boldsymbol{\varphi}_{1}$, then the corresponding generalized eigenvector $\boldsymbol{\varphi}_{2}$ is defined by the relations Strang (2009)

[TABLE]

Given $\Phi_{1}=a_{1}\boldsymbol{\varphi}_{1}+a_{2}\boldsymbol{\varphi}_{2}$ for some $a_{1,2}\in\mathbb{C}$, the transfer matrix acts as

[TABLE]

Thus, we identify

[TABLE]

so that $\mathbf{f}$ is the $N^{\text{th}}$ power of the Jordan normal form of $T$ in the eigenspace of $\rho$. Note that $f_{21}$ has picked up an additional term linear in $N$. In general, we may get terms that grow or decay as $N^{k}\rho^{N-k}$, where $k$ is the difference between the algebraic and geometric multiplicity of an eigenvalue of $T$. Thus, for OBC, the nondiagonalizability of $T$ gives rise to a family of states whose localization is not purely exponential, but has a polynomial decay. This would clearly be most apparent if the repeated eigenvalue lies on the unit circle. Another interesting case is when $\rho=0$, where we get a state that decays to zero within a finite number of steps, independent of the system size.

III Restricting to rank 1

In this section, we restrict the formal discussions of Sec. II to systems with $r=1$, which encompasses many tight-binding models of interest and has the advantage that the relevant computations are analytically tractable. Note that this condition is not directly related to either the range of hoppings or the number of degrees of freedom in a supercell; instead $r=1$ signifies that there exists a local unitary transform under which the Hamiltonian reduces to a form where the consecutive supercells are connected by a single bond.

For $r=1$, the transfer matrix $T\in\mathrm{Mat}(2,\mathbb{C})$ can be written as

[TABLE]

where $\mathcal{G}_{ab}\in\mathbb{C}$ and $\xi\in\mathbb{R}^{+}$ is the (only) singular value of $J$. The eigenvalues of $T$ are

[TABLE]

where

[TABLE]

In Appendix C, we show that $\mathcal{G}_{ab}$’s are rational functions of $\varepsilon$, with the numerator a polynomial in $\varepsilon$ of order $\mathfrak{n}$ for $\mathcal{G}_{vv},\mathcal{G}_{ww}$ and order $\mathfrak{n}-1$ for $\mathcal{G}_{vw},\mathcal{G}_{wv}$. We next specialize the results of Sec. II.4 to the present case and use them to explain various interesting aspects of non-Hermitian systems such as the skin effect and real-space EPs. We also construct a Riemann surface associated with $\varepsilon$, which can be used to define topological invariants for the boundary states.

III.1 Boundary conditions and spectra

We split this discussion between bulk and boundary spectra. For a given transverse momentum, the bulk spectrum is generically a set of closed curves in the complex plane, which can generically be written as $\varepsilon=F(\phi)$ with $\phi\in[0,2\pi]$ and $F$ periodic in $\phi$. The boundary spectrum, on the other hand, is a discrete set of points on the complex plane.

III.1.1 Bulk spectra

For a system of $N$ supercells and PBC, we use Eq. (26) to write the condition for the existence of a Bloch state as

[TABLE]

where $\phi=2\pi\ell/N,\,\ell\in\{0,\dots N-1\}$, and the $N\to\infty$ limit is taken by setting $\phi\in[0,2\pi)$. Since the numerator of $\Delta$ is a polynomial in $\varepsilon$ of order $\mathfrak{n}$, we obtain $\mathfrak{n}$ complex solutions for $\varepsilon$ for each $\phi$ and $\mathbf{k}_{\perp}$. Scanning over $\phi$, we get the PBC bulk spectrum. We reiterate that if $\boldsymbol{\varphi}$ is the eigenvector of $T$ associated with $\mathrm{e}^{i\phi}$, then the corresponding bulk states are given by $\Phi_{n}=\mathrm{e}^{in\phi}\boldsymbol{\varphi}$, which are precisely the Bloch states.

We next turn to the condition for OBC [cf. Eq. (33)], which can be rewritten in the present case as

[TABLE]

for some $\mathfrak{r}\in\mathbb{C}$. This can be used to derive a Cramer’s condition, as was done in Sec. II.4. However, for $r=1$, we can explicitly compute $T^{N}$ Sánchez-Soto et al. (2012) and use it to derive conditions involving only the transfer matrix in the $N\to\infty$ limit. As shown in Appendix D, for $\Gamma\neq 0$,

[TABLE]

where $z=\Delta/2\sqrt{\Gamma}$ and $U_{n}(z)$’s are the Chebyshev polynomials of the second kind, explicitly defined in Eq. (98). Combining this with Eq. (42), we derive the condition for OBC as

[TABLE]

The behavior of the right hand side as $N\to\infty$ strongly depends on $z$. If $z$ is real and $z\in[-1,1]$, we set $z=\cos\phi$ for some $\phi\in[0,\pi]$ and use Eq. (98) to rewrite Eq. (44) as

[TABLE]

The right hand side has poles at $\phi=\ell\pi/(N-1)$ and zeros at $\phi=\ell\pi/N$ with $\ell=0,1,\dots N-1$. Thus, Eq. (45), solved for $\phi$, has $N$ solutions, one lying in each open set $(\ell\pi/N,\ell\pi/(N-1))$. As $N\to\infty$, we get a dense set of solutions, which is our bulk band for OBC. Thus, the condition for the bulk band in terms of $(\varepsilon,\mathbf{k}_{\perp})$ can be written as

[TABLE]

for some $\phi\in[0,\pi]$. Since $\Delta$ is a polynomial in $\varepsilon$ of order $\mathfrak{n}$, this equation has $\mathfrak{n}$ complex solutions for $\varepsilon$ for each $\phi$. Scanning over $\phi$, we thus get the OBC bulk bands generically as a set of $\mathfrak{n}$ closed curves on $\mathbb{C}$. Furthermore, the corresponding eigenvectors can be computed analytically as

[TABLE]

which yield the correct boundary vectors for $n=0,N$.

To obtain further insight into the meaning of this condition, we substitute Eq. (46) in Eq. (39) to conclude that the eigenvalues of the transfer matrix are $\rho_{\pm}=\sqrt{\Gamma}\mathrm{e}^{\pm i\phi}$. The equality of magnitude of the eigenvalues can be alternatively inferred as follows: if $\left|\rho_{+}\right|\neq\left|\rho_{-}\right|$, then for any initial vector $\Phi_{0}$ which is a linear combination of the corresponding eigenvectors, the eigenvector corresponding to the larger magnitude of the eigenvalue will dominate as $N\to\infty$. Thus, a matching condition like Eq. (42) can be satisfied only if the two eigenvalues are equal in magnitude, thereby leading to Eq. (46). This argument is essentially identical to that of Ref. Yao and Wang (2018) with their decay exponents $\beta$ being equal to the eigenvalues of the transfer matrix.

For transfer matrices with $\Gamma=0$, excluded in the above derivation, we find $T^{n}=\Delta^{n-1}T$. Substituting this in Eq. (42) results in $\Delta=0$, which is equivalent to $\Gamma\to 0$ limit of Eq. (46). Thus, we get a bulk state for OBC iff

[TABLE]

Since this is independent of $\phi$ unlike Eq. (45), we get a discrete set of $\mathfrak{n}$ points instead of $\mathfrak{n}$ bands, i.e., each bulk band collapses to a single energy eigenvalue. Furthermore,

[TABLE]

and $\Phi_{n}=\mathbf{0}\;\forall n>2$. Thus, we have a single state for each band, which is localized at the left boundary and has a finite support.

The condition for the bulk states for PBC and OBC can be written concisely as

[TABLE]

with the cosine of a complex angle with $\phi\in[0,2\pi)$ and $\zeta\in\mathbb{R}$. For $\kappa=0$ (OBC) and $\kappa=1$ (PBC), we get $\zeta=0$ and $\frac{1}{2}\log\Gamma$, respectively. We can extend these further by continuously tuning between these two values of $\zeta$ as discussed in Sec. II.4. In this setup, for $0<\kappa<1$, we find some intermediate $\zeta=\zeta_{\kappa}$ that interpolates between [math] and $\frac{1}{2}\log\left|\Gamma\right|$. We derive an approximate expression for $\zeta_{\kappa}$ in Appendix D.

III.1.2 Boundary spectra

The boundary states are obtained as additional discrete solutions to Eq. (44). For $z\notin[-1,1]$, the $N\to\infty$ limit of the right hand side of Eq. (44) is finite, so that Eq. (44) reduces to the condition

[TABLE]

The solutions to these equations give us the boundary spectrum, but more care is needed to physically interpret them. The problem stems from the fact that for $N\to\infty$, we have essentially ignored the boundary condition at the other end, thereby effectively treating the system as semi-infinite. We need to additionally ensure that the state so obtained decays into the bulk. Thus, only those solutions of Eq. (51) describe a physical left boundary mode for which the corresponding eigenvalue of the transfer matrix satisfies $\left|\rho_{\mathrm{L}}\right|<1$, and a similar condition for the right boundary mode.

The boundary states can alternatively be obtained in a more straightforward manner by starting with a semi-infinite system and demanding that the boundary vector is an eigenvector of the transfer matrix, as in Ref. (Dwivedi and Chua, 2016, Sec. III.A). More explicitly, a left boundary state is obtained when

[TABLE]

A similar calculation for the right boundary results in $\mathcal{G}_{ww}=0$ and $\rho_{\mathrm{R}}=\xi\mathcal{G}_{wv}$, in agreement with Eq. (51). We can alternatively write the expressions for boundary spectra as a special case of the equation

[TABLE]

since $\boldsymbol{\varphi}^{T}\mathcal{J}\boldsymbol{\varphi}=0,\,\forall\boldsymbol{\varphi}\in\mathbb{C}^{2}$. Setting $\boldsymbol{\varphi}=(1,0)^{T}$ or $(0,1)^{T}$, we recover the boundary state conditions computed above. In writing this equation, we have ignored the decay condition, so that we obtain physical states (in ED, for instance) only for a subset of the solutions of Eq. (52). On the other hand, a solution to this equation exists for all $\mathbf{k}_{\perp}$. For two-dimensional systems where $\mathbf{k}_{\perp}\in S^{1}$, this fact can be used to define closed curves corresponding to the boundary states on a Riemann surface, as we show in Sec III.3.

III.2 Aspects of non-Hermiticity

We now discuss several intriguing aspects of non-Hermitian systems that can be readily deduced from the knowledge of its transfer matrix.

III.2.1 PBC vs OBC bulk spectra

The difference between the PBC and OBC bulk spectra can be easily visualized by plotting the magnitudes of the eigenvalues of the transfer matrix, as shown in Fig. 2. The bulk bands for PBC are then given by the intersection of the plane $\left|\rho\right|=1$ with the eigenvalues (black lines in Fig. 2), while those for OBC are given by the intersection with the plane where the two eigenvalues are equal in magnitude, i.e., $\left|\rho\right|=\sqrt{\left|\Gamma\right|}$ (blue lines). For $\left|\Gamma\right|\neq 1$, these two planes do not coincide, so that their intercepts, i.e., the curves on the complex $\varepsilon$ plane corresponding to the bulk bands, can be different in the two cases.

On the other hand, these two curves become identical if $\left|\Gamma\right|=1$. We can also see this analytically, since setting $\Gamma=\mathrm{e}^{-2i\chi}$ for some $\chi\in[0,2\pi)$, the conditions for PBC and OBC bulk modes in Eqs. (41) and (46) both reduce to 222 Here we have replaced $(\phi-\chi)$ with $\phi$ in the condition for PBC [cf. Eq. (41)] since we are scanning over $\phi$.

[TABLE]

with $\phi\in[0,2\pi)$. Thus, in the large system limit, the bulk spectra for a system with PBC and OBC are identical iff the transfer matrix is unimodular. In Sec. II.3, we showed that Hermiticity or $\mathrm{PT}$-symmetry implies unimodularity of the transfer matrices for physically relevant energies. This may, however, also be true in more general settings.

This difference between the PBC and OBC bulk spectra can lead to an interesting situation, where as one tunes a parameter (or $k_{y}$), the PBC bulk gap closes while the OBC bulk bands remain gapped. This scenario is depicted in Fig. 2. In this case, any topologically nontrivial boundary states, if present, will also remain qualitatively unchanged, since they cannot be removed without closing the gap between the two bands connected by them, i.e., the OBC bulk bands. Thus, this presents an instance of a dramatic breakdown of the conventional bulk-boundary correspondence.

III.2.2 The non-Hermitian skin effect

To study the skin effect, we need to look at the asymptotic behavior of the states for systems with PBC and OBC. For systems with PBC,

[TABLE]

independent of $n$, as one would expect for Bloch waves. On the other hand, for OBC, we have

[TABLE]

If $\left|\Gamma\right|\neq 1$, the “bulk states”, or more precisely, the states associated with the continuum spectrum for the system with PBC, decay into the bulk. These states are localized on the left boundary for $\left|\Gamma\right|<1$ and on the right boundary for $\left|\Gamma\right|>1$. Thus, the existence of the non-Hermitian skin effect can be deduced simply from the value of $\left|\det T\right|$.

Combining this with the result from the previous subsection, we note that the phenomena of the skin effect and the difference between the PBC and OBC bulk spectra are intimately linked, since they are both governed by the same condition. More explicitly, a non-Hermitian system exhibits the skin effect iff the PBC and OBC bulk spectra are different.

III.2.3 Exceptional points

The Bloch and real-space EPs can also be understood in the transfer matrix formalism. We have a Bloch EP if the condition for the bulk states, i.e., Eq. (41), solved for $\varepsilon$, has a repeated root. The multiplicity of the roots sets the order of the EP. On the other hand, the real-space EPs are obtained when $\left|\Gamma\right|\to 0,\infty$. We remark that the order of the real-space EP is $(N-1)$, where $N$ is the system size, so that we can get EPs with arbitrarily high order for a given Bloch Hamiltonian. On the other hand, the order of the Bloch EPs is limited by the dimensionality of the Bloch Hamiltonian. Thus, if $\left|\Gamma\right|\neq 0,\infty$, then the maximum order of an EP in the real-space spectrum is restricted by the dimensionality of the Bloch Hamiltonian, where we make use of the fact that when $\Gamma\neq 0,\infty$ the a nonunitary similarity transform of the original Hamiltonian yields a Hamiltonian for which $\Gamma=1$ Yao and Wang (2018).

III.2.4 Bulk topological invariants

A bulk-boundary correspondence for non-Hermitian gapped systems can be defined if one computes the “bulk invariant” using the eigenvectors for a system with boundaries as opposed to that for a periodic system. In the present setup, the eigenvectors for PBC vary with position as $\mathrm{e}^{in\phi}$ where $\phi$ can be identified with $k_{x}$ while those for OBC go as $\sqrt{\Gamma}\mathrm{e}^{in\phi}$. This suggests the topological invariants should be computed using a modified Bloch Hamiltonian with the replacement

[TABLE]

For instance, for gapped two-dimensional systems, a modified Chern number can be computed by integrating the biorthogonal Berry curvature on the modified “Brillouin zone” with coordinates $(\phi,k_{y})$.

This approach was indeed shown to predict the existence of edge modes correctly in Refs. Yao and Wang (2018); Yao et al. (2018), albeit only for specific tight-binding models. Our setup thus provides a direct way of analytically generalizing the topological invariants for Hermitian Hamiltonians to non-Hermitian Hamiltonians to arbitrary lattice models without resorting to continuum limits or numerical computations.

III.2.5 Biorthogonality condition

The case of rank-1 systems subsumes the non-Hermitian models discussed in Ref. Kunst et al. (2018a), whose boundary modes can be obtained analytically by construction. As an indicator of the existence of boundary modes, a biorthogonal polarization was proposed, defined in terms of $\mathfrak{p}\equiv\left|\widetilde{\rho}_{\mathrm{L}}^{\ast}\rho_{\mathrm{L}}^{\phantom{\ast}}\right|$, i.e., the product of decay exponents of the left and right eigenstates of the Hamiltonian, localized at the left boundary. It was shown that the boundary states merge into the bulk band when $\left|\mathfrak{p}\right|=1$.

We now derive this quantity using the transfer matrix formalism by considering a semi-infinite non-Hermitian system on $\mathbb{Z}^{+}$. Let $\Psi$ be the right eigenstate of the Hamiltonian for a left boundary mode, with the decay exponent $\rho_{\mathrm{L}}=-[\xi\mathcal{G}_{vw}(\varepsilon_{\mathrm{L}})]^{-1}$, where $\varepsilon_{\mathrm{L}}$ satisfies $\mathcal{G}_{vv}(\varepsilon_{\mathrm{L}})=0$. For the corresponding left eigenstate, we need the transfer matrix $\widetilde{T}$ for $\widetilde{\mathcal{H}}=\mathcal{H}^{\dagger}$ in terms of which the decay exponent is given by $\widetilde{\rho}_{\mathrm{L}}=-[\xi\widetilde{\mathcal{G}}_{vw}(\varepsilon_{\mathrm{L}})]^{-1}$. Using Eq. (18) to relate $\mathcal{G}$ to $\widetilde{\mathcal{G}}$, we find

[TABLE]

Next, we note that the bulk and boundary bands merge for a given $(\varepsilon,\mathbf{k}_{\perp})$ if the conditions for both bulk and boundary states for OBC are simultaneously satisfied. Thus, we seek to solve $\mathcal{G}_{vv}=0=\Delta-2\sqrt{\Gamma}\,\cos\phi$ for some $\phi$. We combine these to get

[TABLE]

which can be rearranged as

[TABLE]

This is solved by $\sqrt{\mathfrak{p}}=\mathrm{e}^{\pm i\phi}$, which is equivalent to demanding that $\left|\mathfrak{p}\right|=1$, precisely what was obtained in Ref. Kunst et al. (2018a). Note that the exact same condition is obtained by alternatively considering $\left|\widetilde{\rho}_{\mathrm{R}}^{\ast}\rho_{\mathrm{R}}^{\phantom{\ast}}\right|$ for the right boundary.

III.3 $\varepsilon$-Riemann surface

The algebraic structure of the transfer matrix naturally lends itself to the construction of a Riemann surface associated with the complex energy. Explicitly, the map $\varepsilon\mapsto\rho$ of in Eq. (39) is not analytic for $\varepsilon\in\mathbb{C}$, since there are square-root singularities at the zeros of $Q(\varepsilon)\equiv\Delta^{2}(\varepsilon)-4\Gamma(\varepsilon)$. Since $\Delta$ and $\Gamma$ are rational functions in $\varepsilon$, so is $Q(\varepsilon)$, with the numerator being a polynomial of order $2\mathfrak{n}$ (see Appendix C for details). Thus, $Q(\varepsilon)$ has exactly $2\mathfrak{n}$ complex roots, which must be connected in pairs by $\mathfrak{n}$ branch cuts. Since these zeros are points where the two eigenvalues of the transfer matrix are degenerate, i.e., where

[TABLE]

we define the branch cuts to lie along the bulk spectrum for OBC. More explicitly, we define the branch cuts as the curves in the $\varepsilon$ plane for which $\rho_{\pm}=\sqrt{\Gamma}\mathrm{e}^{\pm i\phi}$. For example, in Fig. 2, we have $\mathfrak{n}=2$, and the four zeros of $Q(\varepsilon)$ are denoted by dark blue dots, with the branch cuts lying along the blue solid lines.

Two copies of $\mathbb{C}$ are glued along these branch cuts and a compact Riemann surface $\mathfrak{R}$ is then obtained by one-point compactifying these sheets into Riemann spheres and gluing them. By the Riemann-Hurwitz lemma, we deduce that $\mathfrak{R}$ has genus $(\mathfrak{n}-1)$, i.e., one less than the number of Bloch bands. Thus, in the case of Fig. 2, the Riemann surface is a 2-torus. An explicit mapping from the $\varepsilon$-plane with two branch cuts to a torus can be implemented via the elliptic integrals, as shown in Ref. (Dwivedi and Chua, 2016, Sec. IV.A).

This construction particularly caters to a system with OBC. For each $\mathbf{k}_{\perp}$, the continuum states run precisely along the branch cuts. Furthermore, for a two-dimensional system where $\mathbf{k}_{\perp}=k_{y}\in S^{1}$, the boundary modes are a map $S^{1}\to\mathfrak{R}$, which can be classified by a winding number. This is the topological invariant associated with the boundary states. This approach can be used to define a “topologically protected boundary mode” for non-Hermitian systems as modes with nonzero winding numbers, in absence of the conventional definition for Hermitian systems in terms of directed zero crossings.

III.4 A generic two-band model

We now illustrate the ideas discussed in this section by explicit computations on a generic two-band model. We consider a $d$-dimensional system described by an arbitrary Bloch Hamiltonian of the form

[TABLE]

where ${\boldsymbol{\eta}}\colon\mathbb{T}^{d-1}\to\mathbb{C}^{3}$ depends on $\mathbf{k}_{\perp}$, $\boldsymbol{\sigma}=\left(\sigma^{x},\sigma^{y},\sigma^{z}\right)$ is the vector of Pauli matrices, and

[TABLE]

The eigenvalues of the Bloch Hamiltonian are

[TABLE]

where

[TABLE]

with ${\boldsymbol{\eta}}_{R}$ and ${\boldsymbol{\eta}}_{I}$ the real and imaginary parts of ${\boldsymbol{\eta}}$, respectively. Note that $\eta$ here is not the usual norm of ${\boldsymbol{\eta}}\in\mathbb{C}^{3}$, i.e., $\eta^{2}\neq{\boldsymbol{\eta}}\cdot{\boldsymbol{\eta}}^{\ast}$.

To compute the transfer matrix, we identify the hopping and on-site matrices as coefficients of $\mathrm{e}^{ik_{x}}$ and $1$ in the Bloch Hamiltonian, so that

[TABLE]

The SVD results in $J=\xi\,\mathbf{v}\cdot\mathbf{w}^{\dagger}$, with

[TABLE]

The on-site Green’s function is

[TABLE]

Writing $\mathcal{G}$ as a matrix for the given definitions of $\mathbf{v}$ and $\mathbf{w}$, we identify

[TABLE]

Using Eq. (38), the transfer matrix is

[TABLE]

We compute

[TABLE]

in terms of which the eigenvalues of $T$ are given by Eq. (39).

For PBC, the energies of Bloch states are given by Eq. (41), which can be simplified to get

[TABLE]

We note that this expression is identical to Eq. (59) with the identification $\phi\to k_{x}$, as expected. For OBC, the bulk states are given by Eq. (42), which simplifies to

[TABLE]

These states are localized near the left boundary if $\left|\Gamma\right|<1$, i.e, if

[TABLE]

Similarly, they are localized near the right boundary if $\left|\Gamma\right|>1$, i.e, if $\text{Im}[\eta_{x}^{\ast}\eta_{y}]>0$.

When $\Gamma=0,\infty$, i.e., $\eta_{x}=\pm i\eta_{y}$, we get the real-space EP, where the bulk bands collapse to two points with energies $\varepsilon=\pm\sqrt{1+\eta_{z}^{2}}$. The corresponding states are all localized on the leftmost/rightmost supercell for $\eta_{x}=\pm i\eta_{y}$. Since $\eta_{x}\pm i\eta_{y}$ is the intracell hopping, these real-space EPs occur when the system has a completely unidirectional hopping, so that all states are piled up at one end of the system. Finally, we note that Eqns. (66) and (67) become identical if the transfer matrix is unimodular, as follows from Eq. (53).

The boundary states are given by Eq. (52), so that the boundary spectra and the corresponding decay exponents become

[TABLE]

The left boundary state exists for $\mathbf{k}_{\perp}$ where $\left|\eta_{x}+i\eta_{y}\right|>1$, while the right one exists if $\left|\eta_{x}-i\eta_{y}\right|>1$. Using this boundary spectrum, we can also compute

[TABLE]

which signals that the boundary states merge into the bulk bands for $\left|\mathfrak{p}\right|=1$, i.e., for $\left|\eta_{x}^{2}+\eta_{y}^{2}\right|=1$. This is identical to the result obtained from the decay exponents.

To illustrate these transfer-matrix calculations, we briefly discuss a well-known example of the one-dimensional Su-Schrieffer-Heeger (SSH) model Su et al. (1979), whose various non-Hermitian variations have been studied in the literature Malzard et al. (2015); Xiong (2018); Yao and Wang (2018); Lee (2016); Kunst et al. (2018a); Yin et al. (2018); Herviou et al. (2018); Lee and Thomale (2019). We consider a non-Hermitian SSH model with an asymmetry between the left and right intra-cell hoppings Kunst et al. (2018a); Yin et al. (2018). Thus, explicitly we consider the Bloch Hamiltonian

[TABLE]

As shown in Ref. Kunst et al. (2018a), conventional bulk-boundary correspondence is broken for this model, and the bulk states pile up at the ends. We now show how the same result is derived using the transfer matrix.

Since the Bloch Hamiltonian takes the same form as Eq. (57), we readily identify ${\boldsymbol{\eta}}=(t,i\gamma,0)$. Using Eq. (64), we get the transfer matrix

[TABLE]

so that

[TABLE]

Thus, the bulk spectra are given by

[TABLE]

The two bulk spectra are different, i.e., $\varepsilon_{\textrm{PBC}}\neq\varepsilon_{\textrm{OBC}}$ for $\gamma\neq 0$. The system exhibits the non-Hermitian skin effect, since the bulk states for OBC are localized at the left/right end when $t\gamma\lessgtr 0$, as follows from Eq. (68). The transition between these two cases occurs at $\gamma=0$, where the system is Hermitian and hence bulk states extend across the system. Finally, we get real space EPs for $t=\pm\gamma$, where the bulk bands collapse to $\varepsilon=\pm 1$.

Moreover, from Eq. (69), we find that the system exhibits boundary modes in the gap with energies $\varepsilon_{\mathrm{L}}=\varepsilon_{\mathrm{R}}=0$ and decay exponents $\rho_{\mathrm{L}}=-1/(t-\gamma)$ and $\rho_{\mathrm{R}}=-(t+\gamma)$. The biorthogonal polarization in $\mathfrak{p}$ in Eq. (70) then reads as $\mathfrak{p}=t^{2}-\gamma^{2}$ such that the boundary states attach to the bulk bands when $t=\pm\sqrt{\gamma^{2}+1},\pm\sqrt{\gamma^{2}-1}$, in accordance with the results in Ref. (Kunst et al., 2018a, Eq. (9)). Finally, we note that the associated bulk invariant can be computed as the Berry phase computed around a modified Brillouin zone, which is obtained by the replacement

[TABLE]

This expression corresponds exactly to the result derived in Ref. Yao and Wang (2018).

IV Two-Dimensional Examples

We now apply the ideas discussed above to a variety of lattice models for two-dimensional topological phases. The PBC and OBC spectra are computed analytically and compared to those computed using numerical exact diagonalization (ED). We also discuss the topology associated with the boundary states in terms of the energy Riemann surface.

IV.1 A non-Hermitian Chern insulator

We consider a non-Hermitian generalization of the two-dimensional Chern insulator Yao and Wang (2018); Yao et al. (2018); Kawabata et al. (2018a), for which we take

[TABLE]

where $\mathbf{h}=(h_{x},\,h_{y},\,h_{z})\in\mathbb{R}^{3}$. Physically, $h_{x}$ and $h_{y}$ represent an anisotropy in the phase and amplitude of the intracell left and right hoppings, respectively, while $h_{z}$ represents an on-site gain and loss on alternative sublattices.

For $\mathbf{h}=\mathbf{0}$, i.e., the Hermitian limit, the system is gapless for $m=0,\pm 2$, a trivial insulator for $\left|m\right|>2$, and a topological insulator with Chern number $\pm 1$ for $\left|m\right|<2$. For OBC along $x$, the topological phase exhibits modes localized on the edge. Using the transfer-matrix method, we can compute the edge spectra as $\varepsilon_{\mathrm{L},\mathrm{R}}=\mp\sin k_{y}$, with the corresponding decay exponents being $\rho_{\mathrm{L}}=\cos k_{y}-m$ and $\rho_{\mathrm{R}}=1/(\cos k_{y}-m)$, respectively (Dwivedi and Chua, 2016, Sec. III.D.3). Demanding that the edge modes decay into the bulk, we deduce that they exist near $k_{y}=0$ for $0<m<2$ and near $k_{y}=\pi$ for $-2<m<0$. In the following analysis, we set $m=1.4$, for which we get the celebrated Chern insulator spectrum, as shown in Fig. 3.

For the non-Hermitian generalization, we find

[TABLE]

Thus, the transfer matrix is unimodular if $h_{y}=0$ and non-unimodular otherwise. Since these two cases exhibit qualitatively different behaviors, we shall distinguish between them in the following analysis. with the OBC bulk states localized on the left/right edge when

[TABLE]

as follows from Eq. (68). Thus, the system exhibits the non-Hermitian skin effect only when the non-Hermiticity is along $\sigma^{y}$, corresponding to asymmetric hopping within the unit cell. Furthermore, for a fixed $h_{y}$ and $\left|m\right|<1$, there are OBC bulk states localized at both ends of the system, corresponding to different ranges of $k_{y}$.

The bulk states can be computed explicitly from Eqs. (41) and (46). The edge states occur at energies $\varepsilon_{\mathrm{L},\mathrm{R}}=\mp(\sin k_{y}+ih_{z})$, with the associated decay exponents $\rho_{\mathrm{L}}=\eta_{x}+i\eta_{y}$ and $\rho_{\mathrm{R}}=1/(\eta_{x}-i\eta_{y})$, respectively. We now set the terms in $\mathbf{h}$ to $\gamma\in\mathbb{R}^{+}$ one by one and apply the results of Sec. III.4 to deduce the behavior of the OBC spectrum.

IV.1.1 Non-unimodular transfer matrix

We begin with the most interesting case, viz., that with a non-unimodular transfer matrix, by setting $h_{y}=\gamma$. In this case, the system exhibits the non-Hermitian skin effect as well as a difference in the PBC and OBC spectra. Interestingly, in certain parameter ranges, the PBC spectrum is actually gapless, while the OBC spectrum remains gapped with a robust edge mode in the gap (as also pointed out in Ref. Kawabata et al. (2018a) from ED.). The robustness of the edge mode is clear from its winding on the energy Riemann surface, which remains unchanged throughout this transition, as shown in Fig. 4. We also visualize the transition in the PBC spectrum by plotting the complex PBC bulk spectrum as a function of $k_{y}$, which forms a surface whose topology changes from two cylinders to a “pair-of-pants”.

The bulk spectra for PBC and OBC are given by

[TABLE]

where $A=2+m^{2}-\gamma^{2}-2m\cos k_{y}$. We note that $\varepsilon_{\text{OBC}}$’s are either real or come in complex-conjugate pairs, which can also be traced back to the pseudo-Hermiticity of the real-space Hamiltonian Yao and Wang (2018); Yao et al. (2018). The edge spectra given by $\varepsilon_{\mathrm{L},\mathrm{R}}=\mp\sin k_{y}$ with purely real decay exponents. Next, the Bloch Hamiltonian [cf. Eq. (59)] exhibits second-order Bloch EPs at

[TABLE]

Finally, this system also exhibits a pair of real-space EPs at $k_{y,\text{EP}}=\cos^{-1}(m\mp\gamma)$ by setting $\left|\Gamma\right|$ to $0,\infty$, whose order is $N-1$, where $N$ is the system size. At these points, each bulk band collapses to a single point with energy

[TABLE]

As expected, aside from the qualitative difference, the Bloch and real-space EPs occur at different values of $k_{y}$ for a given $m$ and $\gamma$.

We can now analytically deduce the behavior of this system as the non-Hermitian term is turned on. In the following, take $1<m<2$, so that we start in a topological phase for $\gamma=0$. Tuning $\gamma$ up, we nucleate a real-space EP at $k_{y}=0$ when $\gamma=m-1$, for which all the states are localized at the left edge. Further increase in $\gamma$ splits this EP into two real-space EPs at $\pm\cos^{-1}(m-\gamma)$, which move out and merge again at $k_{y}=\pi$ when $\gamma=m+1$. On the other hand, we nucleate a Bloch EP at $k_{y}=0$ for $\gamma=2-m$, which splits into two EPs that merge at $k_{y}=\pi$ when $\gamma=m$. For a full phase diagram obtained numerically, see Ref. Kawabata et al. (2018a).

In Fig. 5, we plot the PBC and OBC spectrum for the non-Hermitian Chern insulator for a fixed $m$, and we choose two values of $\gamma$ in two different phases: one with only real-space EPs and one with both real-space and Bloch EPs. The spectrum was computed numerically using ED for a finite system with PBC/OBC. We also plot the curves obtained by solving the equations for the bulk spectra for $\phi=0,\pi$ (blue solid lines) and $\phi=\pm\pi/2$ (blue dashed lines), which follow various contours of the numerically computed bulk bands. We also plot the analytically computed edge spectrum $\varepsilon_{\mathrm{L},\mathrm{R}}(k_{y})$, only a part of which (corresponding to the decay condition on the eigenvalues) are seen for the particular termination used for the OBC calculation. The spectra in Fig. 5 for PBC and OBC show vastly different qualitative behavior. For PBC, the system goes from gapped to gapless. This effect can be seen more clearly in a 3D plot of the complex bulk band energies as a function of $k_{y}$, as shown in the left and middle panel in Fig. 4. The bulk band topology clearly changes from two cylinders, which can be “flattened” into two bands, to the pair of pants, which cannot be flattened. The OBC spectrum is qualitatively unaffected by this transition. Indeed, we note that the edge states run along a noncontractible loop on the $\varepsilon$-Riemann surface in both cases, as shown in the right panel of Fig. 4. We point out that this behavior was previously observed in Ref. Kawabata et al. (2018a).

The difference between PBC and OBC spectra here can be intuitively understood as the manifestation of a preferred hopping direction, which leads to a pileup of the continuum states at the edges and thereby to an extreme sensitivity to boundary conditions Kunst et al. (2018a). Moreover, when $\gamma$ is chosen such that the hopping in one direction is completely turned off, we get real-space EPs, where the bulk bands indeed collapse to single points, as shown in Fig. 5. The corresponding eigenstates have a finite support, independent of the system size. These EPs are thus associated with an extreme form of unidirectionality.

IV.1.2 Unimodular transfer matrix

We recall that unimodularity of the transfer matrix implies identical qualitative behavior for the spectrum for PBC and OBC (cf. Sec. III.2.1), and thus we only plot band spectra for the latter in this section without loss of information. Moreover, as in this case no real-space EPs may appear, we may make use of the eigenvalues of the Bloch Hamiltonian in Eq. (59) to determine the location of EPs in the spectrum of OBC.

We first set $h_{z}=\gamma$. With a rotation of $\boldsymbol{\sigma}$, the corresponding Bloch Hamiltonian is equivalent to the Bloch Hamiltonian in the case of $h_{y}=\gamma$ with $k_{x}$ and $k_{y}$ interchanged. These two models are thus equivalent up to a $\pi/2$-rotation from a Bloch Hamiltonian perspective, and the EPs are given by expressions identical to the case of $h_{y}=\gamma$ [cf. Eqs. (75) and (76)] with the roles of $k_{x}$ and $k_{y}$ interchanged. In particular, we find the same behavior for the EPs as we tune $\gamma$, while the systems look completely different from a real-space perspective. The edge spectra are $\varepsilon_{\mathrm{L},\mathrm{R}}=\mp\left(\sin k_{y}+i\gamma\right)$, and have now picked up an imaginary part, so that the edge modes now have a finite lifetime. The opposite sign of the imaginary part in the energy of these states is explained by the fact that they are primarily localized on alternate sublattices. Their decay exponents, however, stay real. In Fig. 7, we plot the spectra for OBC with the same parameters as for the previous case, i.e., $m=1.4$ and $h_{z}=0.75$. We note that the EPs appear at $k_{y}=0$, which is indeed suggested by Eqs. (75) and (76), once the roles of $k_{x}$ and $k_{y}$ are interchanged.

We finally set $h_{x}=\gamma$, so that our model is the usual lattice Dirac equation with a complex mass. The Bloch Hamiltonian, and hence the bulk spectra for both PBC and OBC, exhibit second-order EPs when

[TABLE]

The edge spectra $\varepsilon_{\mathrm{L},\mathrm{R}}=\mp\sin k_{y}$ is real, but the corresponding decay exponents now pick up an imaginary part. We find Bloch EPs when the imaginary part of the above equations disappears and the real part is confined to $[-1,1]$. In Fig. 7, we plot the spectra for OBC with appropriate parameter values and we indeed find Bloch EPs at these values of $k_{y}$.

IV.2 A non-Hermitian 2D Dirac semimetal

In this section, we consider a non-Hermitian lattice model with $\mathrm{PT}$ symmetry, viz., a two-dimensional Dirac semimetal. This is essentially a two-dimensional stacking of the $\mathrm{PT}$-symmetric Su-Schrieffer-Heeger chains studied in Refs. Schomerus (2013); Lieu (2018b). Explicitly, we consider the model of Sec. III.4 with

[TABLE]

with $m,\gamma\in\mathbb{R}$. The $\mathrm{PT}$ operation is implemented by $\mathcal{PT}=\sigma^{x}K$. Physically, the non-Hermitian term $i\gamma\sigma_{z}$ in the Bloch Hamiltonian can be understood as a gain on one of the site types and a loss on the other type.

For $\gamma=0$, we recover the Hermitian limit. In this case, the model is gapped and trivial if $\left|m\right|>2$, while for $\left|m\right|<2$ we get two Dirac points in the 2D Brillouin zone at $\mathbf{k}=\left(0,\pm\cos^{-1}(m-1)\right)$ for $0<m<2$ and $\mathbf{k}=\left(\pi,\pm\cos^{-1}(m+1)\right)$ for $-2<m<0$, which is indeed shown in the top left panel of Fig. 8. Turning on the non-Hermitian term $\gamma\neq 0$, these Dirac points broaden into curves of EPs (or exceptional lines). Using the eigenvalues of the Bloch Hamiltonian in Eq. (59), we compute that these lie along

[TABLE]

These lines separate the $\mathrm{PT}$-unbroken and $\mathrm{PT}$-broken phases as is shown in the phase diagram in the top right panel of Fig. 8 for $m=1.4$ and $\gamma=0.5$. Explicitly, we have a $\mathrm{PT}$-unbroken phase, i.e., real energies, if the left-hand side is positive, and $\mathrm{PT}$ is broken otherwise.

From the transfer-matrix perspective, we find

[TABLE]

and $\Gamma=1$. The latter implies that the bulk spectra for PBC and OBC are identical in both $\mathrm{PT}$-unbroken and $\mathrm{PT}$-broken phases 333 Note that in general, this is ensured by the $\mathrm{PT}$-symmetry only in the $\mathrm{PT}$-unbroken phase. . The bulk spectrum for both PBC and OBC is given by $\Delta=2\cos\phi$, i.e., the Bloch spectrum. The edge states satisfy $\varepsilon_{\mathrm{L},\mathrm{R}}=\mp i\gamma$, so that we get a gain for the left edge state and loss for the right one. This is expected, since each of the edge states obtained above is primarily localized on one of these two types of sites. We plot the spectra for OBC in the bottom row of Fig. 8 with $m=1.4$ and $\gamma=0.5$. We find Bloch EPs for four different values of $k_{y}$ as predicted by the phase diagram (cf. top right panel of Fig. 8]).

IV.3 A non-Hermitian Hofstadter model

We finally consider a non-Hermitian generalization of the Hofstadter model Hatsugai (1993), variations of which have also been studied in Refs. Matveenko and Ouvry (2014); Chernodub and Ouvry (2015). Recall that the Hofstadter model is essentially a hopping model on a square lattice with hopping of equal magnitude across each link and phases corresponding to a rational flux of $2\pi\phi$ with $\phi=p/q$ threading each plaquette. We introduce non-Hermiticity in this model either by adding on-site terms $i\gamma_{n}$ corresponding to absorption/decay and by staggering the magnitude of left hopping vs right hopping by $\delta_{n}$. Explicitly, assuming translation invariance and PBC along $y$, we consider the Hamiltonian

[TABLE]

where $\gamma_{n},\delta_{n}\in\mathbb{R}$. The original Hofstadter model is periodic with period $q$. To recover this periodicity as well as $J_{\mathrm{L}}=J_{\mathrm{R}}$ required for our transfer matrix construction, we choose $\gamma_{n}=\gamma_{n\,(\text{mod}\,q)}$, $\delta_{n}=\delta_{n\,(\text{mod}\,q)}$ and $\delta_{q}=0$, while the remaining $(2q-1)$ parameters are arbitrary. We can now write the hopping and on-site matrices explicitly. The hopping matrix $J$ has all entries equal to zero except $J_{1,q}=1$ and satisfies $J^{2}=0$ for all $q>1$. On the other hand, $M$ has $2t_{n}\cos(k_{y}-2\pi n\phi)+i\gamma_{n}$ as its diagonal entries and $(1\pm\delta_{n})$’s on the first diagonal, with $\gamma_{n},\delta_{n}\in\mathbb{R}$. Explicitly, for the simplest nontrivial case of $\phi=1/3$, we set

[TABLE]

and we find

[TABLE]

We plot the spectrum with OBC and PBC in Fig. 9. When we choose $\delta_{1}=\pm 1$ (or $\delta_{2}=\pm 1$), we find $\Gamma=0$ or $\infty$ and the continuum bands in the spectrum for the OBC shrink to exceptional lines of order $(N-1)$.

V Discussion

In this paper, we construct a generalized transfer matrix for non-Hermitian noninteracting tight-binding models and show that various peculiarities of non-Hermitian models are related to simple and readily computable features of the transfer matrix. For instance, the unimodularity of the transfer matrix, a property of Hermitian systems as well as $\mathrm{PT}$-symmetric systems in the $\mathrm{PT}$-unbroken phase, is shown to be related to a bulk-boundary correspondence, while a departure from unimodularity is related to a difference between the PBC and OBC spectra as well as the non-Hermitian skin effect, thereby establishing a formal connection between these two phenomena. These results are illustrated through various examples, which are analytically tractable and highlight the power of this method. For a particular class of systems where the transfer matrix is two-dimensional, we find that the singularity of the transfer matrix is accompanied by the appearance of real-space EPs in the OBC spectrum at which all states are confined to the boundary.

We further find that the topological invariants proposed in Refs. Yao and Wang (2018, 2018); Yao et al. (2018), which make use of a deformation of the Brillouin zone to complex quasimomentum, can be naturally understood and generalized in the language of transfer matrices. Explicitly, for $2\times 2$ transfer matrices, they are obtained straightforwardly by replacing $\mathrm{e}^{ik}$ in the Bloch Hamiltonian with $\sqrt{\det T}\mathrm{e}^{ik}$, i.e., by replacing the eigenvalue of the transfer matrix for PBC with that for OBC. Moreover, the real-space invariant proposed in Ref. Kunst et al. (2018a) in the form of the biorthogonal polarization can also be readily obtained making use of transfer matrices. Additionally, we show that, at least for two-dimensional systems, we can assign a topological invariant to the edge states by identifying the edge spectra $\varepsilon_{\mathrm{L},\mathrm{R}}(k_{y})$ as closed loops on the energy Riemann surface. If these loops are noncontractible, the edge modes can only be removed if the bulk gap for OBC collapses, which may be independent of the bulk gap for PBC for non-Hermitian systems. This is indeed the case for the non-Hermitian Chern insulator model studied in Sec. IV.1.1, where a gap closing in the PBC spectrum leaves the edge state unaffected. Therefore, the transfer matrices, which give access to the eigensystem of the model with open boundary conditions, provide a crucial insight to establish the topology of these models.

Interestingly, the extension to non-Hermitian Hamiltonians of many of the results previously obtained in Ref. Dwivedi and Chua (2016) for Hermitian models involves several aspects that were not needed or not quite visible in the Hermitian case. For instance, the bulk spectra for PBC and OBC correspond to very different mathematical conditions, which happen to coincide when the transfer matrix is unimodular. A further attraction of this generalization is the possibility of complex energies, which lends a physical significance to the construction of an energy Riemann surface, which was introduced for Hermitian systems purely for mathematical convenience.

The explicit computations for the rank-1 case discussed in Sec. III do not readily generalize to higher ranks. In certain cases, however, additional structure in the model can be leveraged. For instance, if the transfer matrix happens to be symplectic, then the spectra can be obtained by using the restriction on the eigenvalue problem that it imposes (Dwivedi and Chua, 2016, Appendix C). In absence of such additional structure, the spectra can be computed directly using the results of Sec. II.4, in particular, Eqs. (26) and (37). Both of these require the diagonalization of the transfer matrix and hence should be numerically, if not analytically, tractable. A more algebraic approach to these cases, especially the construction of associated Riemann surfaces and definition of boundary invariants, requires further study.

We emphasize that the transfer-matrix approach is also useful for systems not described by tight-binding models. For instance, transfer matrices have been extensively used to study localization in the phenomenological Chalker-Coddington network models Evers and Mirlin (2008). A non-Hermitian version thereof was also studied in Ref. Brouwer et al. (1997) for a one-dimensional periodic chain with an imaginary vector potential. We believe that the insights gleaned from our extension of transfer-matrix formalism to non-Hermitian systems would prove useful in diverse contexts.

One particularly interesting direction for further investigation is the implications of symmetries on the transfer matrix. Indeed, one of the central parts of the study of Hermitian topological phases has focused on their classification, e.g, gapped non-interacting fermionic Hamiltonians belong to one of the 10 equivalence (Altland-Zirnbauer) classes based on their antiunitary symmetries Altland and Zirnbauer (1997); Schnyder et al. (2008). The non-Hermitian analogs of this “ten-fold way” are described by Bernard and LeClair Bernard and LeClair (2001); *bernard-leclair_book, and allow for many more symmetries. Indeed, recent investigations of non-Hermitian Bloch Hamiltonians have resulted in their classification in terms of genuinely non-Hermitian symmetry classes Gong et al. (2018); Kawabata et al. (2018c), while in previous studies explicit non-Hermitian topological phases with a trivial Hermitian limit have already been constructed Budich et al. (2019), and the role of time-reversal symmetry was also investigated in Ref. Sato et al. (2012). The transfer-matrix approach can shed further light on the general classification of non-Hermitian Hamiltonians, since as we show in this paper for $\mathrm{PT}$ symmetry, the symmetries of the Hamiltonian may be implemented in a nontrivial manner on the transfer matrix. Moreover, the possible difference between spectra obtained for PBC and OBC makes a classification of systems based on a real-space approach highly relevant.

The discrepancy between periodic and open spectra necessitates the access to exact solutions of the real-space Hamiltonian in order to probe the topological aspects of a system. The transfer-matrix approach, being a purely real-space construction, is the ideal platform for such an endeavor. We thus believe that transfer matrices provide a natural framework for a general understanding of non-Hermitian systems.

Acknowledgements.

We thank Emil J. Bergholtz, Victor Chua, Henry Legg, Max Geier, Sebastian Diehl, and Thors Hans Hansson for useful conversations. FKK was funded by the Swedish Research Council (VR) and the Knut and Alice Wallenberg Foundation. VD was funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) – Projektnummer 277101999 – TRR 183 (project B03).

Appendix A Non-Hermitian Hamiltonians and $\mathrm{PT}$ symmetry

Non-Hermitian systems are described by Hamiltonians with $\mathcal{H}\neq\mathcal{H}^{\dagger}$, with the adjoint taken under the usual inner product on the Hilbert space. Consequently, the left and right eigenstates, defined by

[TABLE]

are not related by conjugation and $\braket{\psi_{m}}{\psi_{n}}\neq\delta_{mn}$, i.e., the eigenbasis is no longer orthonormal. However, $\braket{\phi_{m}}{\psi_{n}}=\delta_{mn}$, so that one defines a biorthogonal basis and uses it to compute the so-called biorthogonal expectation value Brody (2013) of an observable $\mathcal{O}$ as $\braket{\phi}{\mathcal{O}}{\psi}$. The non-orthogonality of eigenvectors may also lead to two or more eigenvectors becoming linearly dependent, so that the eigenstates do not span the Hilbert space. Such Hamiltonians are termed defective. For the Hamiltonian dependent on a set of parameters, their value for which the Hamiltonian is defective is termed an exceptional point (EP) Kato (1995); Heiss (2012), whose order is defined as the number of eigenvectors that coalesce at that EP.

Systems with a parity-time-reversal ($\mathrm{PT}$) symmetry Jones-Smith and Mathur (2010) form a particularly well-studied subset of non-Hermitian systems. $\mathrm{PT}$ is implemented as an antilinear and antiunitary operator $\mathcal{PT}=\mathcal{U}K$ with $\mathcal{U}$ a unitary matrix and $K$ denoting complex conjugation. Furthermore, $\left(\mathcal{PT}\right)^{2}=\pm\openone$, which corresponds to $\mathcal{U}^{T}=\pm\,\mathcal{U}$. Owing to its antilinearity, $\mathcal{PT}$ cannot have eigenvectors. Explicitly, this is because if $\psi$ were an eigenvector, so would $a\psi$ for any $a\in\mathbb{C}$, but

[TABLE]

However, a state may be left invariant up to a phase under $\mathrm{PT}$. A further constraint is imposed by

[TABLE]

so that $\left(\mathcal{PT}\right)^{2}=-\openone$, there are no states left invariant under $\mathcal{PT}$, while for $\left(\mathcal{PT}\right)^{2}=\openone$, one might have such states.

Under the $\mathrm{PT}$ operation, the Bloch Hamiltonian transforms as

[TABLE]

If $\mathcal{H}_{\mathrm{B}}$ is $\mathrm{PT}$-symmetric, the states must satisfy

[TABLE]

where $\varepsilon\in\mathbb{C}$ in general. For $\left(\mathcal{PT}\right)^{2}=\openone$, a state may be invariant under $\mathrm{PT}$. In this case, the two eigenvalue equations in Eq. (80) contain the same eigenvector, so that $\varepsilon=\varepsilon^{\ast}$. When this is true for all eigenstates, the system is termed to be in a $\mathrm{PT}$-unbroken phase. On the other hand, if there are eigenstates that are not invariant under $\mathrm{PT}$, i.e., the $\mathrm{PT}$ symmetry is spontaneously broken, then we are in a $\mathrm{PT}$-broken phase. For $\left(\mathcal{PT}\right)^{2}=-\openone$, no state is left invariant by $\mathrm{PT}$, so that the system is always in a $\mathrm{PT}$-broken phase. The states satisfy an analog of Kramers’ theorem Jones-Smith and Mathur (2010), viz., all eigenstates come in orthogonal pairs whose energies are related by complex conjugation. A trivial consequence of this case is that the Bloch Hamiltonian must be even dimensional.

Appendix B $\mathrm{PT}$-symmetry and SVD

We derive the constraints imposed by the $\mathrm{PT}$-symmetry on the singular values and vectors of the hopping matrix. Recall that the reduced SVD of a matrix $J\in\mathrm{Mat}(\mathfrak{n},\mathbb{C})$ is defined as (Strang, 2009, Sec. 6.3)

[TABLE]

where $r=\text{rank}\,J$, $\xi_{n}>0$ are the singular values and $\mathbf{v}_{i},\mathbf{w}_{i}$ the corresponding left/right singular vectors. A pair of singular vectors $\mathbf{v},\mathbf{w}$ is defined by the relations

[TABLE]

Note that these expressions are manifestly invariant under a simultaneous phase rotation $\mathbf{v}\to\mathrm{e}^{i\theta}\mathbf{v},\;\mathbf{w}\to\mathrm{e}^{i\theta}\mathbf{w}$.

In presence of a $\mathrm{PT}$ symmetry, the hopping matrix satisfies $J^{T}=\mathcal{U}^{\dagger}J\mathcal{U}$ with $\mathcal{U}\in\mathrm{U}(\mathfrak{n})$, so that Eq. (82) becomes

[TABLE]

A complex conjugation leads to

[TABLE]

where $\widetilde{\mathbf{w}}=\mathcal{U}^{T}\mathbf{v}^{\ast}$ and $\widetilde{\mathbf{v}}=\mathcal{U}^{T}\mathbf{w}^{\ast}$. We thus find two sets of vectors satisfying the equation for a singular value $\xi$, so that either the two vectors are proportional, i.e.,

[TABLE]

or $\xi$ is degenerate as a singular value with two sets of left and right singular vectors. When the two vectors are proportional, we find

[TABLE]

Here, we consider the two possible cases: If $\mathcal{U}^{T}=\mathcal{U}$, then Eq. (83) holds iff $\left|\lambda\right|=1$. Setting $\lambda=\mathrm{e}^{2i\chi}$, we get $\mathbf{v}=\mathrm{e}^{2i\chi}\mathcal{U}^{T}\mathbf{w}^{\ast}$, and we can use the invariance of Eq. (82) under the phase rotation $\mathbf{v}\to\mathbf{v}\mathrm{e}^{i\chi},\,\mathbf{w}\to\mathbf{w}\mathrm{e}^{i\chi}$ to fix the phase of $\mathbf{v},\mathbf{w}$ such that $\mathbf{v}=\mathcal{U}\mathbf{w}^{\ast}$ and $\mathbf{w}=\mathcal{U}\mathbf{v}^{\ast}$. Continuing this for all singular vectors of $J$, we find

[TABLE]

which are the requisite conditions on the singular vectors of $J$.

On the other hand, if $\mathcal{U}^{T}=-\mathcal{U}$, then $\nexists\,\lambda\in\mathbb{C}$ for which Eq. (83) holds. Therefore, the singular value $\xi$ must be degenerate with the corresponding right and left singular vectors reading $\mathbf{w},\widetilde{\mathbf{w}}$ and $\mathbf{v},\widetilde{\mathbf{v}}$, respectively. In this degenerate sector, we set $\mathfrak{v}=(\mathbf{v},\widetilde{\mathbf{v}})$ and $\mathfrak{w}=(\mathbf{w},\widetilde{\mathbf{w}})$ to write the SVD of $J$ in this subspace as $\mathfrak{v}\,\openone_{2}\,\mathfrak{w}$. Using the definition of $\widetilde{\mathbf{v}}$ and $\widetilde{\mathbf{w}}$, this leads to

[TABLE]

$J$ thus falls apart into these $2\times 2$ $\mathfrak{v}\,\openone_{2}\,\mathfrak{w}$ blocks with degenerate singular values, so that $r=\text{rank}\,J$ must be even. Defining $\Sigma=\mathscr{J}\otimes\openone_{r/2}$, we find

[TABLE]

which are the requisite conditions on the singular vectors of $J$.

Appendix C Schur complement and inversion

Here, we explore some of the algebraic properties of the rank-1 transfer matrix. Given the on-site matrix $M\in\mathrm{Mat}(\mathfrak{n},\mathbb{C})$ and energy $\varepsilon\in\mathbb{C}$, the on-site Green’s function can be written as

[TABLE]

where $Q(\varepsilon)\equiv\det(\varepsilon\openone-M)$ is by definition a polynomial in $\varepsilon$ of order $\mathfrak{n}$, and $G(\varepsilon)\equiv\text{adj}(\varepsilon\openone-M)\in\mathrm{Mat}(\mathfrak{n},\mathbb{C})$ is the adjugate (i.e., the matrix of minors) (Strang, 2009, Sec. 4.4.1) of $M$, whose elements are polynomials in $\varepsilon$ of order $\leq\mathfrak{n}-1$. The transfer matrix can be written as

[TABLE]

where $G_{ab}(\varepsilon)\equiv Q(\varepsilon)\,\mathcal{G}_{ab}(\varepsilon),\;a,b\in\{\mathbf{v},\mathbf{w}\}$ are polynomials (instead of rational functions) in $\varepsilon$. The discriminant becomes

[TABLE]

and we are interested in its zeros. Naively, owing to the $G_{vv}^{2}G_{ww}^{2}$ term, the numerator is a polynomial of order $\leq\varepsilon^{4(\mathfrak{n}-1)}$ in $\varepsilon$. However, in the following we show that

[TABLE]

is a polynomial in $\varepsilon$ of order $\leq\mathfrak{n}-1$, so that the leading-order term in the numerator of $\Delta^{2}-4\Gamma$ arises from $Q^{2}$, rendering it a polynomial in $\varepsilon$ of order $2\mathfrak{n}$.

We begin by using the basis independence of the transfer-matrix computation to choose a basis of $\mathbb{C}^{\mathfrak{n}}$ in which $\mathbf{v}=(1,0,0,\dots)$ and $\mathbf{w}=(0,1,0,\dots)$, so that

[TABLE]

where $A\in\mathrm{Mat}(2,\mathbb{C})$ and the numerator of $f(\varepsilon)$ is simply $\det A$. Using the fact that for block matrices

[TABLE]

we can rewrite $f(\varepsilon)$ as

[TABLE]

where we have used the fact that $\det G=Q^{\mathfrak{n}-1}$. Using the inversion formula for block matrices (Dwivedi and Chua, 2016, Eq. (A8)), Eq. (88) becomes

[TABLE]

But we also have

[TABLE]

We thus identify $S^{-1}$ as the $(\mathfrak{n}-2)\times(\mathfrak{n}-2)$ lower right term in the block structure of $(\varepsilon\openone-M)$. Finally,

[TABLE]

where $\left(.\right)_{X}$ denotes the restriction to the subspace spanned by the orthogonal complement of $\mathbf{v}$ and $\mathbf{w}$. Thus, $f(\varepsilon)$ is a polynomial in $\varepsilon$ of order $\mathfrak{n}-2$, which proves our desired result.

Appendix D Explicit computations for $r=1$

We compute $T^{n}(\varepsilon)$ explicitly for $T\in\mathrm{Mat}(2,\mathbb{C})$ and arbitrary $n\in\mathbb{Z}$ and use it to derive explicit conditions on $\varepsilon$ for obtaining an eigenstate of a system with OBC.

D.1 Computing $T^{n}$

We start off with Cayley’s theorem, which states that a matrix satisfies its characteristic equation. Thus, $T\in\mathrm{Mat}(2,\mathbb{C})$ satisfies

[TABLE]

where $\Delta=\text{tr}\,T$ and $\Gamma=\det T$. For $\Gamma=0$, we simply get $T^{n}=\Delta^{n-1}T$. On the other hand, for $\Gamma\neq 0$, using Eq. (94) repeatedly, one can reduce $T^{n}=A_{n}T+B_{n}\openone$. Using $T^{n+1}=T\,T^{n}$, i.e.,

[TABLE]

we obtain a recursion relation for the coefficients

[TABLE]

These reduce to a three-term recursion for $A_{n}$ as

[TABLE]

with the initial condition $A_{1}=1$ and $A_{2}=\Delta$. For $\Gamma\neq 0$, setting $A_{n}=\Gamma^{(n-1)/2}a_{n}$, this reduces to

[TABLE]

with the initial conditions $a_{1}=1$ and $a_{2}=2z$. This is the defining relation for the Chebyshev polynomials of the second kind $U_{n}(z)$ (Bateman, 1953, Sec. 10.11), so that we identify $a_{n}=U_{n-1}(z)$. This leads to our final result

[TABLE]

which can be easily evaluated using the closed form expressions for the Chebyshev polynomials (Bateman, 1953, Sec. 10.11, Eqn. (2))

[TABLE]

where

[TABLE]

The former expression for $U_{n}(z)$ in Eq. (98) is useful for arbitrary $z\in\mathbb{C}$, while the latter is naturally more useful when $z\in\mathbb{R}$ and $\left|z\right|<1$.

D.2 Open Boundary conditions

For $\Gamma=0$, using $T^{n}=\Delta^{n-1}T$, we trivially get

[TABLE]

Thus, the bulk spectrum collapses to a single point given by $\Delta=\Gamma=0$. For $\Gamma\neq 0$, using the explicit form of $T$ [cf. Eq. (38)], the condition for OBC in Eq. (42) becomes

[TABLE]

where $\mathfrak{r}$ is arbitrary and we have defined

[TABLE]

The condition on $\varepsilon$ can now be written as

[TABLE]

This can be recast in another useful form by substituting Eq. (100) in Eq. (D.2). We get

[TABLE]

where we have used the recursion relation for the Chebyshev polynomial in Eq. (96) as well as the product formula

[TABLE]

In the last step, we have used the fact that $U_{0}(z)=1$.

The conditions for OBC can be further reduced in the large-$N$ limit. Using the first definition of Chebyshev polynomials from Eq. (98), we have

[TABLE]

For $N\to\infty$, we need to consider three cases. For $\left|\lambda\right|>1$, we can compute

[TABLE]

while for $\left|\lambda\right|<1$, we get

[TABLE]

Finally, for $\left|\lambda\right|=1$, setting $\lambda=\mathrm{e}^{i\phi}$, we get

[TABLE]

which does not have a well-defined limit as $N\to\infty$; instead, the right hand side oscillates wildly, since it has zeros at $\phi=k\pi/N$ and poles at $k\pi/(N-1)$. Thus, for any $q(\phi)$, we get $N$ solutions in $\phi\in[0,\pi)$, which become dense in $[0,2\pi)$ as $N\to\infty$. This is our bulk band for OBC. In terms of $z$, this also corresponds to setting $z=\cos\phi$.

We can also derive a condition for a boundary condition interpolating between PBC and OBC [cf. Sec II.4], for we demand that $1\in\text{Spec}\!\left[KT^{N}\right]$. Since $KT^{N}$ is a $2\times 2$ matrix, its two eigenvalues must be $1$ and $\det(KT^{N})=\det K\left(\det T\right)^{N}=\Gamma^{N}$, so that

[TABLE]

where $\zeta_{0}=\frac{1}{2}\log\Gamma$. We next compute the left hand side using Eq. (97) as

[TABLE]

where we have used

[TABLE]

Using the recursion relation for the Chebyshev polynomials to replace

[TABLE]

Eq. (103) can be reduced to

[TABLE]

Setting $z=\cos\chi$ for some $\chi=\phi+i\zeta;\,\phi,\zeta\in\mathbb{R}$ and using the definition of $U_{n}(\cos\chi)$, this can be rearranged to get

[TABLE]

The left hand side is now independent of $N$. For the bulk states, we shall require that the right hand side does not have a limit as $N\to\infty$ (as in the OBC case above). Thus, the condition for an eigenstate becomes $\Delta=2\sqrt{\Gamma}\cos(\phi+i\zeta)$ for some $\phi\in[0,\pi]$, with

[TABLE]

for $\kappa$ close to 1. Therefore, the spectrum for $\kappa\neq 0,1$, unlike for the PBC and OBC case, is in general quite sensitive to the system size.