# Supplemental materials for “Topological Invariant for Multi-Band Non-
hermitian Systems with Chiral Symmetry”

Chun-Chi Liu  National Laboratory of Solid State Microstructures, Department
of Physics, Nanjing University, Nanjing 210093, China  Liu-Hao Li  National
Laboratory of Solid State Microstructures, Department of Physics, Nanjing
University, Nanjing 210093, China  Jin An<EMAIL_ADDRESS>National
Laboratory of Solid State Microstructures, Department of Physics, Nanjing
University, Nanjing 210093, China  Collaborative Innovation Center of Advanced
Microstructures, Nanjing University, Nanjing 210093, China

<EMAIL_ADDRESS>Boundary Conditions(OBCs) and Their Simplification

The eigenequation of $H(\beta)$ can be written as:

$\sum_{m=-M_{2}}^{M_{1}}T_{m}\psi(j-m)=E\psi(j),$ (1)

where $M_{1}<j<L-M_{2}$. When $1\leq j\leq M_{1}$, the eigenequation becomes
the left-boundary conditions, which are given by:

$\begin{split}&\sum_{m=-M_{2}}^{M_{1}-1}T_{m}\psi(M_{1}-m)=E\psi(M_{1}),\\\
&\sum_{m=-M_{2}}^{M_{1}-2}T_{m}\psi(M_{1}-1-m)=E\psi(M_{1}-1),\\\
&\quad\quad\vdots\\\ &\sum_{m=-M_{2}}^{0}T_{m}\psi(1-m)=E\psi(1).\\\
\end{split}$ (2)

When $L-M_{2}+1\leq j\leq L$, the eigenequation becomes the right-boundary
conditions, which are:

$\begin{split}&\sum_{m=-M_{2}+1}^{M_{1}}T_{m}\psi(L-M_{2}+1-m)=E\psi(L-M_{2}+1),\\\
&\sum_{m=-M_{2}+2}^{M_{1}}T_{m}\psi(L-M_{2}+2-m)=E\psi(L-M_{2}+2),\\\
&\quad\quad\vdots\\\ &\quad\sum_{m=0}^{M_{1}}\quad T_{m}\psi(L-m)=E\psi(L).\\\
\end{split}$ (3)

Figure 1: A one-dimensional open chain with $L$ unit cells represented by
solid circles, where each unit cell contains $2N$ internal degrees of freedom.
Here the hollow circles at both ends of the chain are the extra unit cells
added for the discussion of open-boundary conditions.

The two sets of equations contain $M_{1}$ and $M_{2}$ equations respectively.
In order to make the eigenequation Eq.(1) be obeyed for all unit cells,
including those near the left and right boundaries, one can introduce extra
unit cells as shown in Fig.1. Combination of these new equations and the above
two sets of equations lead to the following simplified left-boundary
conditions:

$\begin{split}&T_{M_{1}}\psi(0)=0,\\\
&T_{M_{1}}\psi(-1)+T_{M_{1}-1}\psi(0)=0,\\\ &\quad\vdots\\\
&T_{M_{1}}\psi(-M_{1}+1)+T_{M_{1}-1}\psi(-M_{1}+2)+...+T_{1}\psi(0)=0,\end{split}$
(4)

and right-boundary conditions:

$\begin{split}&T_{-M_{2}}\psi(L+1)=0,\\\
&T_{-M_{2}}\psi(L+2)+T_{-M_{2}+1}\psi(L+1)=0,\\\ &\quad\vdots\\\
&T_{-M_{2}}\psi(L+M_{2})+T_{-M_{2}+1}\psi(L+M_{2}-1)+...+T_{-1}\psi(L+1)=0.\end{split}$
(5)

If both $T_{M_{1}}$ and $T_{-M_{2}}$ are invertible, then according to the
first equations of the two sets of equations, we have: $\psi(0)=\psi(L+1)=0$.
From the second equations, we further have $\psi(-1)=\psi(L+2)=0$. Along the
similar argument, we would find the following equivalent boundary equations:

$\begin{split}&\psi(0)=\psi(-1)=...=\psi(-M_{1}+1)=0,\\\
&\psi(L+1)=\psi(L+2)=...=\psi(L+M_{2})=0.\\\ \end{split}$ (6)

These boundary conditions can be understood as that the wave function
$\psi(j)$ at the extra unit cells denoted by the hollow circles in Fig.1
should all be zero. Let’s consider the system described by $H(\rho\beta)$
whose wave function is denoted by $\widetilde{\psi}(j)$. The real-space
Hamiltonian $\mathcal{H}\rho$ of $H(\rho\beta)$ and the real-space Hamiltonian
$\mathcal{H}$ of $H(\beta)$ are related by a similarity transformation
$\mathcal{H}\rho=S^{-1}\mathcal{H}S$, with $S={\rm
diag}(\rho^{1},\rho^{2},...,\rho^{L})\otimes I_{2N}$. Accordingly, we have
$\widetilde{\psi}(j)=\rho^{-j}\psi(j)$ and the hopping matrices of
$H(\rho\beta)$ become $T_{m}\rho^{-m}$, so for $\widetilde{\psi}(j)$ we can
obtain the similar OBCs to Eq.(6).

<EMAIL_ADDRESS>Maximization of Matrices $\mathcal{M_{L}}$ and $\mathcal{M_{R}}$

Introduce a continuous parameter $\Lambda$, and $R_{\pm}$ can be viewed as
functions of $\Lambda$: $R_{\pm}(\beta,\Lambda)$ obeying
$R_{\pm}(\beta,0)=R_{\pm}(\beta)$. Here we shall prove in this section that
for a gapful open system with $|\beta_{p}(0)|\neq|\beta_{p+1}(0)|$, while
keeping the system gapful, one can always vary the parameter $\Lambda$ to
continuously modify the positions of the $\beta$ roots of
det$R_{\pm}(\beta,\Lambda)=0$ to make the rank of matrices of
$\mathcal{M_{L}}$ and $\mathcal{M_{R}}$ in the main text maximized, without
breaking the chiral symmetry. This kind of variation of the system is
obviously a topologically equivalent one, so the final state should share the
same topological invariant $\nu$ with the initial one.

By solving the eigenequation of $E=0$:

$H(\beta)\phi_{i}=0\rightarrow\begin{cases}R_{+}(\beta)v_{i}=0\\\
R_{-}(\beta)u_{i}=0,\end{cases}$ (7)

we get the $\beta$ roots and their corresponding eigenvectors:
$\beta^{+}_{i}$, $\phi_{i}^{+}=(0,v_{i})^{T}$($\beta^{-}_{i}$,
$\phi_{i}^{-}=(u_{i},0)^{T}$), $i=1$, $2$, …, $N(M_{1}+M_{2})$. From all the
eigenvectors $v_{i}$ of $R_{+}$, we construct a $N(M_{1}+M_{2})\times
N(M_{1}+M_{2})$ matrix $\mathcal{M}_{N}^{+}$:

$\mathcal{M}_{N}^{+}=(\Phi^{+}_{1},\Phi^{+}_{2},...,\Phi^{+}_{N(M_{1}+M_{2})}),$
(8)

with the column vector
$\Phi^{+}_{i}=(v_{i},(\beta_{i}^{+})^{-1}v_{i},...,(\beta_{i}^{+})^{-M_{1}-M_{2}+1}v_{i})^{T}=(1,(\beta_{i}^{+})^{-1},...,(\beta_{i}^{+})^{-M_{1}-M_{2}+1})^{T}\otimes
v_{i}$.

When $N=1$ or in the two-band case, the matrix $\mathcal{M}_{1}^{+}$ becomes:

$\begin{split}\mathcal{M}_{1}^{+}=\left(\begin{array}[]{cccc}1&1&...&1\\\
(\beta_{1}^{+})^{-1}&(\beta_{2}^{+})^{-1}&...&(\beta_{M_{1}+M_{2}}^{+})^{-1}\\\
\vdots&\vdots&\vdots&\vdots\\\
(\beta_{1}^{+})^{-M_{1}-M_{2}+1}&(\beta_{2}^{+})^{-M_{1}-M_{2}+1}&...&(\beta_{M_{1}+M_{2}}^{+})^{-M_{1}-M_{2}+1}\\\
\end{array}\right),\end{split}$ (9)

which is a Vandermonde matrix, whose determinant is a well-known polynomial
$\mathcal{P}_{V}$ of $(\beta_{i}^{+})^{-1}$:

$\text{det}\mathcal{M}_{1}^{+}=\prod_{1\leq i<j\leq
M_{1}+M_{2}}((\beta_{i}^{+})^{-1}-(\beta_{j}^{+})^{-1})\equiv\mathcal{P}_{V}((\beta_{1}^{+})^{-1},(\beta_{2}^{+})^{-1},...,(\beta^{+}_{M_{1}+M_{2}})^{-1}),$
(10)

If all $\beta_{i}^{+}$ are different from each other, then obviously
det$\mathcal{M}_{1}^{+}\neq 0$. Otherwise, there exist repeated $\beta$ roots,
then by fixing the other roots we can easily do a parameter change to make the
repeated roots split without closing the gap, since the two roots
$\beta_{p}(0)$ and $\beta_{p+1}(0)$ with unequal magnitude are still left
unchanged. So for any state of the system, either this state or its
topologically equivalent one would always obey: det$\mathcal{M}_{1}^{+}\neq
0$. For later convenience, we call vectors like
$(1,(\beta_{i}^{+})^{-1},...,(\beta_{i}^{+})^{-M_{1}-M_{2}+1})^{T}\equiv\Phi_{\beta^{+}_{i}}$
Vandermonde vectors. The space spanned by $n$ $M$-dimensional Vandermonde
vectors with different $\beta$ values is $n$-dimensional if $n<M$, and is
$M$-dimensional otherwise.

When $N\geq 2$, each column vector $v_{i}$ has $N$ components and
$R_{+}(\beta)$ is a $N\times N$ matrix. The $(kl)$ element of $R_{+}(\beta)$
is a polynomial of $\beta$:

$R_{kl}(\beta)=\sum_{m=-M_{2}}^{M_{1}}t^{kl}_{m}\beta^{-m},$ (11)

where $k$, $l$ =$1$, $2$,…, $N$. Design a set of continuous functions of
$\Lambda$: $\beta^{+}_{i}(\Lambda)$, satisfying
$\beta^{+}_{i}(0)=\beta^{+}_{i}$, and

$\begin{split}&\beta^{+}_{1}(1)=\beta^{+}_{2}(1)=...=\beta^{+}_{N}(1)=\widetilde{\beta}^{+}_{1},\\\
&\beta^{+}_{N+1}(1)=\beta^{+}_{N+2}(1)=...=\beta^{+}_{2N}(1)=\widetilde{\beta}^{+}_{2},\\\
&\quad\quad\vdots\\\
&\beta^{+}_{(M_{1}+M_{2}-1)N+1}(1)=\beta^{+}_{(M_{1}+M_{2}-1)N+2}(1)=...=\beta^{+}_{(M_{1}+M_{2})N}(1)=\widetilde{\beta}^{+}_{M_{1}+M_{2}}.\end{split}$
(12)

Here $\widetilde{\beta}^{+}_{i}$ are the target $\beta$ with different
absolute values satisfying:
$|\beta^{+}_{(i-1)N+1}|\leq|\widetilde{\beta}^{+}_{i}|\leq|\beta^{+}_{iN}|$.
So all $\beta^{+}_{i}(\Lambda)$ are expected to be divided into $M_{1}+M_{2}$
groups, and each group contains $N$ $\beta$ roots which would merge into one
single $N$-fold root as $\Lambda$ is varied to be $1$. This is shown
schematically in Fig.2. Now we treat all the coefficients $t^{kl}_{m}$ as
unknown variables which can be viewed as functions of $\Lambda$:
$t^{kl}_{m}(\Lambda)$. They can be found by solving the equations:
det$R_{+}(\beta^{+}_{i},\Lambda)=0$, with $i=1$, $2$, …, $N(M_{1}+M_{2})$ and
imposing the initial conditions: $t^{kl}_{m}(0)=t^{kl}_{m}$. The number of the
non-linear equations is $N(M_{1}+M_{2})$, while the number of the unknown
variables is $N^{2}(M_{1}+M_{2}+1)$. So for any fixed value of $\Lambda$, we
generally have multiple solutions of $t^{kl}_{m}(\Lambda)$, but for
$\Lambda=0$, we have only one solution: $t^{kl}_{m}(0)=t^{kl}_{m}$. Thus we
have found a solution of $R_{\pm}(\beta,\Lambda)$, which is a continuous
function of $\Lambda$. When $\Lambda=1$, each $\widetilde{\beta}_{i}$ is a
$N-$fold root of det$R_{+}(\beta^{+}_{i},\Lambda=1)=0$, indicating that its
eigen-space should be $N$-dimensional and so
$R_{+}(\beta^{+}_{i},\Lambda=1)\equiv 0$. Therefore, as varying $\Lambda$ from
$0$ to $1$, one would continuously transform $R_{+}(\beta,\Lambda)$ from
$R_{+}(\beta)$ to $R_{+}(\beta,\Lambda=1)\propto\beta^{-NM_{1}}\prod_{1\leq
i\leq M_{1}+M_{2}}(\beta-\widetilde{\beta}_{i}^{+})\times I_{N}$. When
$\Lambda=1$, eigenvectors $v_{i}$ can be chosen to be
$(v_{1},v_{2},...,v_{N(M_{1}+M_{2})})=(I_{N},I_{N},...,I_{N})$. It’s easy to
find that
$\mathcal{M}_{N}^{+}(\Lambda=1)=\mathcal{\widetilde{M}}_{1}^{+}\otimes I_{N}$,
where $\mathcal{\widetilde{M}}_{1}^{+}$ is the Vandermonde matrix of
$(\widetilde{\beta}^{+}_{i})^{-1}$, and so
det$\mathcal{M}_{N}^{+}(\Lambda=1)=\mathcal{P}^{N}_{V}((\widetilde{\beta}_{1}^{+})^{-1},(\widetilde{\beta}_{2}^{+})^{-1},...,(\widetilde{\beta}^{+}_{M_{1}+M_{2}})^{-1})\neq
0$

In the exactly similar way, from all the eigenvectors $u_{i}$ of $R_{-}$ one
can construct the $N(M_{1}+M_{2})\times N(M_{1}+M_{2})$ matrix
$\mathcal{M}_{N}^{-}$:

$\mathcal{M}_{N}^{-}=(\Phi^{-}_{1},\Phi^{-}_{2},...,\Phi^{-}_{N(M_{1}+M_{2})}),$
(13)

with the column vectors
$\Phi^{-}_{i}=(u_{i},(\beta_{i}^{-})^{-1}u_{i},...,(\beta_{i}^{-})^{-M_{1}-M_{2}+1}u_{i})^{T}=(1,(\beta_{i}^{-})^{-1},...,(\beta_{i}^{-})^{-M_{1}-M_{2}+1})^{T}\otimes
u_{i}=\Phi_{\beta^{-}_{i}}\otimes u_{i}$. We also have:
det$\mathcal{M}_{N}^{-}=\mathcal{P}^{N}_{V}((\widetilde{\beta}_{1}^{-})^{-1},(\widetilde{\beta}_{2}^{-})^{-1},...,(\widetilde{\beta}^{-}_{M_{1}+M_{2}})^{-1})\neq
0$, where $\widetilde{\beta}_{i}^{-}$ are the target $\beta$ at $\Lambda=1$ in
this case.

Figure 2: Continuous Variation of $R_{+}(\beta)$, where the $\beta$ roots of
det$R_{+}(\beta)=0$ are varied continuously. The $NM$ roots are equally
divided into $M$ groups with the $N$ roots in each group being finally merged
into one repeated root, where $M=M_{1}+M_{2}$.

Under this situation, we discuss why the ranks of the matrices
$\mathcal{M_{L}}$ and $\mathcal{M_{R}}$ in the main text have been maximized.
We take $\mathcal{M_{L}}$ as an example to illustrate this fact. If $w=0$,
then $p_{1}=p_{2}=p/2=NM_{1}$, $q_{1}=q_{2}=q/2=NM_{2}$. From
$\mathcal{M}_{N}^{+}(\Lambda=1)$, choose the first $NM_{1}$ column vectors
$\Phi^{+}_{i}=\Phi_{\beta^{+}_{i}}\otimes v_{i}$. For each vector we remove
the last $NM_{2}$ components and replace $v_{i}$ by $\phi_{i}=(0,v_{i})^{T}$
to obtain $\Phi^{+}_{i,\mathcal{L}}$. These new vectors are comprised of $N$
sets of $M_{1}$-dimensional Vandermonde vectors with each set spanning an
independent $M_{1}$-dimensional space. From $\mathcal{M}_{N}^{-}(\Lambda=1)$,
also choose its first $NM_{1}$ column vectors
$\Phi^{-}_{i}=\Phi_{\beta^{-}_{i}}\otimes u_{i}$. For each vector we also
remove the last $NM_{2}$ components but replace $u_{i}$ by
$\phi_{i}=(u_{i},0)^{T}$ to obtain $\Phi^{-}_{i,\mathcal{L}}$. These new
vectors are also comprised of $N$ sets of $M_{1}$-dimensional Vandermonde
vectors with each set spanning an independent $M_{1}$-dimensional space.
Obviously, two spaces spanned by $\Phi^{+}_{i,\mathcal{L}}$ or
$\Phi^{-}_{i,\mathcal{L}}$ are orthogonal. Since the matrix
$\mathcal{M_{L}}=(\Phi_{1,\mathcal{L}}^{+},...,\Phi_{NM_{1},\mathcal{L}}^{+},\Phi_{1,\mathcal{L}}^{-},...,\Phi_{NM_{1},\mathcal{L}}^{-})$
is comprised of $\Phi^{+}_{i,\mathcal{L}}$ and $\Phi^{-}_{i,\mathcal{L}}$, we
have
Rank$\mathcal{M_{L}}=$Rank$(\Phi_{1,\mathcal{L}}^{+},...,\Phi_{NM_{1},\mathcal{L}}^{+})$+Rank$(\Phi_{1,\mathcal{L}}^{-},...,\Phi_{NM_{1},\mathcal{L}}^{-})=NM_{1}+NM_{1}=p$.
Thus, when $w=0$, Rank$\mathcal{M_{L}}$ is maximized to be $p$. Similarly, it
can be shown that Rank$\mathcal{M_{R}}$ is maximized to be $q$.

When $w\neq 0$, $p_{1}=NM_{1}+w$, $p_{2}=NM_{1}-w$, $q_{1}=NM_{2}-w$ and
$q_{2}=NM_{2}+w$. We only consider the case of $0<w<N$ for simplicity. Let’s
make a slight modification on the continuous variation of $\Lambda$. We change
the target $\beta$ roots slightly: only the $(M_{1}+1)$th equation of Eq.(12)
is modified to be:

$\begin{split}&\beta^{+}_{M_{1}N+1}(1)=\beta^{+}_{M_{1}N+2}(1)=...=\beta^{+}_{M_{1}N+w}(1)=\widetilde{\beta}^{+}_{M_{1}+1},\\\
&\beta^{+}_{M_{1}N+w+1}(1)=\beta^{+}_{M_{1}N+w+2}(1)=...=\beta^{+}_{(M_{1}+1)N}(1)=\widetilde{\beta^{\prime}}^{+}_{M_{1}+1},\end{split}$
(14)

where
$|\beta^{+}_{M_{1}N+1}|\leq|\widetilde{\beta}^{+}_{M_{1}+1}|\leq|\beta^{+}_{M_{1}N+w}|<|\beta^{+}_{M_{1}N+w+1}|\leq|\widetilde{\beta^{\prime}}^{+}_{M_{1}+1}|\leq|\beta^{+}_{(M_{1}+1)N}|$.
So here the $(M_{1}+1)$th group of $\beta^{+}_{i}(\Lambda)$ roots would merge
into two repeated roots $\widetilde{\beta}^{+}_{M_{1}+1}$ and
$\widetilde{\beta^{\prime}}^{+}_{M_{1}+1}$, where the former is a $w$-fold
root while the latter is a $(N-w)$-fold root. The eigenspace of the two roots
formed by their corresponding eigenvectors $v_{i}$ is $w$\- or $(N-w)$\-
dimensional, respectively. From $\mathcal{M}_{N}^{+}(\Lambda=1)$, choose the
first $p_{1}$ column vectors $\Phi^{+}_{i}=\Phi_{\beta^{+}_{i}}\otimes v_{i}$.
For each vector we remove the last $NM_{2}$ components and replace $v_{i}$ by
$\phi_{i}=(0,v_{i})^{T}$, so we obtain $\Phi^{+}_{i,\mathcal{L}}$. These new
vectors have already contained $N$ sets of $M_{1}$-dimensional Vandermonde
vectors with each set spanning an independent $M_{1}$-dimensional space. So
Rank$(\Phi_{1,\mathcal{L}}^{+},...,\Phi_{p_{1},\mathcal{L}}^{+})$ has been
maximized to be $NM_{1}=p/2$. From $\mathcal{M}_{N}^{-}(\Lambda=1)$, also
choose its first $p_{2}$ column vectors
$\Phi^{-}_{i}=\Phi_{\beta^{-}_{i}}\otimes u_{i}$. For each vector we also
remove the last $NM_{2}$ components but replace $u_{i}$ by
$\phi_{i}=(u_{i},0)^{T}$, so we obtain $\Phi^{-}_{i,\mathcal{L}}$. If proper
representation of $v_{i}$ is chosen, these new vectors are found to be
comprised of $N-w$ sets of $M_{1}$-dimensional Vandermonde vectors with each
set spanning an independent $M_{1}$-dimensional space, together with $w$ sets
of $(M_{1}-1)$-dimensional Vandermonde vectors with each set spanning an
independent $(M_{1}-1)$-dimensional space. So
Rank$(\Phi_{1,\mathcal{L}}^{-},...,\Phi_{p_{2},\mathcal{L}}^{-})=NM_{1}-w=p/2-w$.
Therefore, we have
Rank$\mathcal{M_{L}}=$Rank$(\Phi_{1,\mathcal{L}}^{+},...,\Phi_{p_{1},\mathcal{L}}^{+})$+Rank$(\Phi_{1,\mathcal{L}}^{-},...,\Phi_{p_{2},\mathcal{L}}^{-})=p/2+(p/2-w)=p-w$,
indicating that it has been maximized. Thus, when $w>0$, Rank$\mathcal{M_{L}}$
is maximized to be $p-w$. Similarly, it can be shown that
Rank$\mathcal{M_{R}}$ is simultaneously maximized to be $q-w$.

So we have proved that for a state of the open system denoted by $w$ in the
main text, by continuously varying $R_{\pm}(\beta)$ without closing the gap,
the ranks of $\mathcal{M}_{L}$ and $\mathcal{M}_{R}$ in the main text can be
maximized to be $p-|w|$ and $q-|w|$ respectively.

<EMAIL_ADDRESS>Theorem on $\pm 1$ change at topological transitions of the Winding
Number Defined on One of the SubGBZs

In this section we prove that when the non-hermitian system is sufficiently
near the topological transition point, if one replaces the closed loop
$\mathcal{L_{\beta}}$ in the integral of the topological invariant $\nu$
defined in the main text, by one of the subGBZs, the change of the new winding
number $\nu^{{}^{\prime}}$ at the transition point would be
$\Delta\nu^{{}^{\prime}}=\pm 1$. The subGBZ chosen corresponds to the branch
of the OBS closest to the $E=0$ point.

When the system under OBCs is near the topological phase transition, the
system must be nearly gapless and then one part of its OBS must be approaching
the $E=0$ point. Under parameter change, the transition process can be
summarized into two scenarios, as shown in Fig.3, where one is the innerpoint
connection and the other is the endpoint connection. At the critical point
where the system is gapless, as shown in Fig.3 (a2) and (b2), the state of the
isolated point $E=0$ belongs to one pair of the open arcs of the OBS and its
$\beta_{p}(0)$ and $\beta_{p+1}(0)$ would simultaneously fall on the
corresponding subGBZ. Now consider the situation when the system is
sufficiently close to the transition point.

Figure 3: Two typical topological transition processes under parameter change
in a chiral symmetric non-hermitian system under OBCs, where gap closes and
then opens. Each pair of arcs in each figure represents the corresponding pair
of the OBSs on $E$ plane closest to $E=0$. The isolated solid circle at $E=0$
in (a3) and (b3) denotes the topological protected zero modes with $\nu$
denoting the corresponding topological invariant. (a1)-(a3): innerpoint
touching process. (b1)-(b3): endpoint touching process.

The subGBZ as a closed loop, is composed of the $\beta_{p}$ and $\beta_{p+1}$
points of the energy branch of the OBS closest to $E=0$. Choose a special
state $E_{0}$ on the branch which is nearest to $E=0$. If the system is
sufficiently close to the transition point, then $E_{0}$ would be approaching
$E=0$ so that its $\beta_{j}(E_{0})$ is nearly equal to $\beta_{j}(0)$.

Let’s first consider the case where the innerpoints $\pm E_{0}$ are going to
be touching with each other, as shown in Fig.3(a1)-(a2). Expand the polynomial
$f(E,\beta)\equiv$det$(H(\beta)-E)=0$ near $E_{0}$ and its $\beta_{j}(E_{0})$.
We have

$\frac{\partial f}{\partial E}\bigg{|}_{(E_{0},\beta_{j}(E_{0}))}\Delta
E+\frac{\partial
f}{\partial\beta}\bigg{|}_{(E_{0},\beta_{j}(E_{0}))}\Delta\beta=0.$ (15)

So the small deviations $\Delta E$ and $\Delta\beta$ are connected with each
other. For any two kinds of deviations ($\Delta E,\Delta\beta$), and ($\Delta
E^{{}^{\prime}},\Delta\beta^{{}^{\prime}}$), one must have

${\Delta E^{{}^{\prime}}}/{\Delta
E}={\Delta\beta^{{}^{\prime}}}/{\Delta\beta}.$ (16)

This is actually the conformal property of analytic functions in complex
analysis. Now take $\Delta E$ as the deviation along the OBS, and $\Delta
E^{{}^{\prime}}$ as the deviation along the line segment, as shown in Fig.4
(a1). If $j$ is chosen to be $p$ or $p+1$, we would have two deviations
$\Delta\beta$, starting from $\beta_{p}$ or $\beta_{p+1}$ respectively, as
shown in Fig.4 (a2). What’s important is that $\Delta\beta$ at $\beta_{p}$ and
$\Delta\beta$ at $\beta_{p+1}$ are pointing oppositely along the
subGBZPhysRevB.105.045422 , namely, one is clockwise while the other is
anticlockwise. Therefore, according to the conformal property Eq.(16), the
deviation $\Delta\beta^{{}^{\prime}}$ at $\beta_{p}$($\beta_{p+1}$) would be
pointing towards the interior(exterior) of the subGBZ. So among two points
$\beta_{p}(0)$ and $\beta_{p+1}(0)$, the subGBZ only encloses $\beta_{p}(0)$.

Figure 4: Two deviation paths in $E$ plane and the corresponding paths of
their $\beta_{p}$ and $\beta_{p+1}$ in $\beta$ plane for the two cases when
the system is sufficiently near the critical point, as shown in Fig.3. Here,
$E_{0}$ is the point on the OBS which is chosen to be the nearest one to
$E=0$. One path on $E$ plane is the straight line segment connecting $E_{0}$
to $E=0$, the other path starting from $E_{0}$ is a small arc along the same
OBS.

For the case of the endpoints $\pm E_{0}$ touching with each other, as shown
in Fig.4 (b1)-(b2), since $\beta_{p}(E_{0})=\beta_{p+1}(E_{0})$, one has
$\partial f/\partial\beta\mid_{(E_{0},\beta_{p}(E_{0}))}=\partial
f/\partial\beta\mid_{(E_{0},\beta_{p+1}(E_{0}))}=0$. So the expansion of the
polynomial $f(E,\beta)$ near $E_{0}$ and its $\beta_{p}(E_{0})$ must be made
in the second order of $\beta$ and we obtain:

$\begin{split}&\frac{\partial f}{\partial
E}\bigg{|}_{(E_{0},\beta_{p}(E_{0}))}\Delta
E+\frac{1}{2}\frac{\partial^{2}f}{\partial\beta^{2}}\bigg{|}_{(E_{0},\beta_{p}(E_{0}))}\Delta\beta^{2}=0.\end{split}$
(17)

So similarly for two kinds of deviations, we must have:

${\Delta E^{{}^{\prime}}}/{\Delta
E}=({\Delta\beta^{{}^{\prime}}}/{\Delta\beta})^{2},$ (18)

where the two kinds of deviations $\Delta E$ and $\Delta E^{{}^{\prime}}$ are
taken similarly as before, as shown in Fig.4 (b1). The two corresponding
deviations of $\Delta\beta^{{}^{\prime}}$ is found to be pointing oppositely.
Thus for $\beta_{p}(0)$ and $\beta_{p+1}(0)$ points, the subGBZ only encircles
$\beta_{p}(0)$.

Therefore for both cases, we have proved that among $\beta_{p}(0)$ and
$\beta_{p+1}(0)$, the subGBZ as a closed loop only encircles $\beta_{p}(0)$,
with $\beta_{p+1}(0)$ being kept outside. Consider a topological transition
process with $\Delta\nu=\pm 1$. Before or after the transition, the system
would always be gapful, but is characterized by different topological
invariant $\nu$. Let’s denote $\beta^{i}_{p}(0)$($\beta^{f}_{p}(0)$) and
$\beta^{i}_{p+1}(0)$($\beta^{f}_{p+1}(0)$) for the system just before(after)
the phase transition. According to our theorem, at the transition point, one
of the $\beta_{p}(0)$ and $\beta_{p+1}(0)$ points is from det$R_{+}=0$, while
the other is from det$R_{-}=0$, and they would exchange to induce the
topological transition. So during the phase transition, $\beta^{i}_{p}(0)$
would continuously change to be $\beta^{f}_{p+1}(0)$, while
$\beta^{i}_{p+1}(0)$ would continuously change to be $\beta^{f}_{p}(0)$. Then
the change of the phase winding along the subGBZ for det$R_{+}$ would be minus
that for det$R_{-}$. This indicates that the change of the new winding number
$\nu^{{}^{\prime}}$ defined on the subGBZ would be $\pm 1$. So we complete the
proof of the statement at the beginning of this section. We remark that for a
two-band system, its GBZ naturally becomes the subGBZ discussed here.

<EMAIL_ADDRESS>Model Example Exhibiting Relationship between the Invariant $\nu$
and the Edge States

In this section, by studying an explicit model in detail, we will follow the
way in the theorem proof in the main text, step by step, to find out the edge
states for the model system under OBCs. On the other hand, by direct
diagonalization of the Hamiltionian for the corresponding open system, we
obtain the numerical results on the edge states of $E=0$. In the following, we
will compare these numerical results with our analytical calculations, to
illustrate our conclusions on the topologically protected edge states(TPESs)
given in the main text.

When a $2N$-band non-hermitian system with chiral symmetry takes the form
Eq.(3) in the main text, it is convenient to introduce the two sublattices A
and B, where the first(last) $N$ internal degrees of freedom belong to
sublattice A(B), as shown in Fig.5. We choose the model in the main text,
which is a four-band model with only the nearest-neighbor hoppings, implying
$N=2$, $M_{1}=M_{2}=1$ and $p=q=4$.
First, the parameter $\lambda$ is chosen to be $\lambda=0.1$, corresponding to
a topologically nontrivial state with $\nu=-1$. By solving det$H(\beta)=0$, we
get eight $\beta$ roots of $E=0$, which are given by:

$\begin{split}\begin{array}[]{cccc}\beta_{1}=\beta_{1}^{-}=-0.195,&\quad\beta_{2}=\beta_{1}^{+}=-0.401,&\quad\beta_{3}=\beta_{2}^{-}=-0.495,&\beta_{4}=\beta_{3}^{-}=0.811,\\\
\beta_{5}=\beta_{2}^{+}=-1.149,&\beta_{6}=\beta_{3}^{+}=2.014,&\quad\beta_{7}=\beta_{4}^{-}=-2.442,&\quad\beta_{8}=\beta_{4}^{+}=-5.634,\end{array}\end{split}$
(19)

and accordingly their corresponding eigenmodes are:

$\begin{split}\begin{array}[]{cccc}\phi_{1}=\phi_{1}^{-}=\left(\begin{array}[]{c}-0.569\\\
-0.822\\\ 0\\\ 0\\\
\end{array}\right),&\phi_{2}=\phi_{1}^{+}=\left(\begin{array}[]{c}0\\\ 0\\\
-0.664\\\ -0.748\\\
\end{array}\right),&\phi_{3}=\phi_{2}^{-}=\left(\begin{array}[]{c}-0.084\\\
-0.997\\\ 0\\\ 0\\\
\end{array}\right),&\phi_{4}=\phi_{3}^{-}=\left(\begin{array}[]{c}-0.464\\\
0.886\\\ 0\\\ 0\\\ \end{array}\right),\\\
\phi_{5}=\phi_{2}^{+}=\left(\begin{array}[]{c}0\\\ 0\\\ 0\\\ 1.000\\\
\end{array}\right),&\quad\phi_{6}=\phi_{3}^{+}=\left(\begin{array}[]{c}0\\\
0\\\ -0.409\\\ 0.912\\\
\end{array}\right),&\phi_{7}=\phi_{4}^{-}=\left(\begin{array}[]{c}-0.714\\\
-0.700\\\ 0\\\ 0\\\
\end{array}\right),&\phi_{8}=\phi_{4}^{+}=\left(\begin{array}[]{c}0\\\ 0\\\
-0.648\\\ -0.761\\\ \end{array}\right).\end{array}\end{split}$ (20)

Figure 5: Open chain and its unit cells which can be divided into two parts:
sites belonging to sublattice A and those belonging to sublattice B.

Here the superscript $\pm$ indicates that the corresponding $\beta$ roots come
from $R_{\pm}$. So the matrices
$\mathcal{M_{L}}=(\phi_{1}^{+},\phi_{1}^{-},\phi_{2}^{-},\phi_{3}^{-})$ and
$\mathcal{M_{R}}=(\beta_{2}^{+}\phi_{2}^{+},\beta_{3}^{+}\phi_{3}^{+},\beta_{4}^{+}\phi_{4}^{+},\beta_{4}^{-}\phi_{4}^{-})$
can be explicitly expressed as

$\begin{split}\mathcal{M_{L}}=\left(\begin{array}[]{cccc}0&-0.569&-0.084&-0.464\\\
0&-0.822&0.996&0.886\\\ -0.664&0&0&0\\\ -0.748&0&0&0\\\
\end{array}\right),\quad\mathcal{M_{R}}=\left(\begin{array}[]{cccc}0&0&0&1.744\\\
0&0&0&1.710\\\ 0&-0.824&3.649&0\\\ -1.149&1.840&4.292&0\\\
\end{array}\right).\end{split}$ (21)

It’s easy to check that Rank$\mathcal{M_{L}}=$Rank$\mathcal{M_{R}}=3$, namely,
they have been maximized so either $\mathcal{M_{L}}|\mathcal{L}\rangle=0$ or
$\mathcal{M_{R}}|\mathcal{R}\rangle=0$ has one independent solution, which can
be expressed as
$|\mathcal{L}\rangle=(c_{1}^{+},c_{1}^{-},c_{2}^{-},c_{3}^{-})^{T}$ and
$|\mathcal{R}\rangle=(c_{2}^{+},c_{3}^{+},c_{4}^{+},c_{4}^{-})^{T}$. The
solutions are found to be: $|\mathcal{L}\rangle=(0,0.335,0.765,-0.550)^{T}$
and $|\mathcal{R}\rangle=(0.922,0.377,0.0852,0)^{T}$. Notice that in these
expressions the decimal numbers shown here are all retained to a few decimal
places, but in the actual calculations they have all been kept in high
precision to ensure the accuracy of the calculating results.

From the two solutions one can see that $|\mathcal{L}\rangle$ only contains
the coefficients of $\phi^{-}$, while $|\mathcal{R}\rangle$ only contains
those of $\phi^{+}$. This indicates that $|\mathcal{L}\rangle$ is a left-
localized state residing on sublattice A, while $|\mathcal{R}\rangle$ is a
right-localized one residing on sublattice B. By using these coefficients, the
corresponding real-space normalized wave function $\psi_{\mathcal{L}}$ and
$\psi_{\mathcal{R}}$ can be written as:

$\begin{split}&\psi_{\mathcal{L}}(j)\propto\sum_{j=1}^{L}c_{1}^{-}(\beta_{1}^{-})^{j}\phi_{1}^{-}+c_{2}^{-}(\beta_{2}^{-})^{j}\phi_{2}^{-}+c_{3}^{-}(\beta_{3}^{-})^{j}\phi_{3}^{-},\\\
&\psi_{\mathcal{R}}(j)\varpropto\sum_{j=1}^{L}c_{2}^{+}(\beta_{2}^{+})^{j-L}\phi_{2}^{+}+c_{3}^{+}(\beta_{3}^{+})^{j-L}\phi_{3}^{+}+c_{4}^{+}(\beta_{4}^{+})^{j-L}\phi_{4}^{+}.\end{split}$
(22)

The space distributions of these two wave functions are shown in Fig.6
(a)-(b). For a finite system under OBCs, we have two final real-space wave
functions $\psi_{\pm}$ which are the linear combinations of
$\psi_{\mathcal{L}}$ and $\psi_{\mathcal{R}}$:

$\psi_{\pm}=\mathcal{A}\psi_{\mathcal{L}}\pm\mathcal{B}\psi_{\mathcal{R}},$
(23)

where $\psi_{+}$ and $\psi_{-}$ are related by chiral symmetry:
$\mathcal{S}\psi_{\pm}=-\mathcal{A}\psi_{\mathcal{L}}\pm\mathcal{B}\psi_{\mathcal{R}}=-\psi_{\mp}$,
with $\mathcal{S}$ the chiral operator. The coefficients $\mathcal{A}$ and
$\mathcal{B}$ are determined by the matrices $\mathcal{M_{RL}}$ and
$\mathcal{M_{LR}}$, which take the following forms:

$\begin{split}&\mathcal{M_{RL}}=((\beta_{1}^{+})^{L+1}\phi_{1}^{+},(\beta_{1}^{-})^{L+1}\phi_{1}^{-},(\beta_{2}^{-})^{L+1}\phi_{2}^{-},(\beta_{3}^{-})^{L+1}\phi_{3}^{-}),\\\
&\mathcal{M_{LR}}=((\beta_{2}^{+})^{-L}\phi_{2}^{+},(\beta_{3}^{+})^{-L}\phi_{3}^{+},(\beta_{4}^{+})^{-L}\phi_{4}^{+},(\beta_{4}^{-})^{-L}\phi_{4}^{-}).\end{split}$
(24)

Substitute $|\mathcal{L}\rangle=(0,0.335,0.765,-0.550)^{T}$ and
$|\mathcal{R}\rangle=(0.922,0.377,0.0852,0)^{T}$ into the equation

$\frac{\mathcal{B}}{\mathcal{A}}=\frac{\sqrt{\langle\mathcal{R}|\mathcal{M_{RL}|L\rangle}}}{\sqrt{\langle\mathcal{L}|\mathcal{M_{LR}|R\rangle}}},$
(25)

one can find
$|\mathcal{B}/\mathcal{A}|\approx|\beta^{-}_{3}\beta^{+}_{2}|^{L/2}=|\beta_{4}\beta_{5}|^{L/2}$.
Here since $|\beta_{4}\beta_{5}|<1$, when $L$ is sufficiently large, it will
be $|\mathcal{B}|\ll|\mathcal{A}|$, indicating the contribution of
$\psi_{\mathcal{L}}$ to $\psi_{\pm}$ is much greater than that of
$\psi_{\mathcal{R}}$. So for this parameter, it looks like that the two TPESs
are all localized at the left boundary, and when $L\rightarrow\infty$,
$|\mathcal{B}|$ will be approaching zero so that $\psi_{+}$ and $\psi_{-}$
linked by chiral symmetry would coalesce into one single state
$\psi_{\mathcal{L}}$. This is actually the origin of the defectiveness of the
TPESs. In Fig.6 (c) we show the numerical matrix diagonalization result of the
TPESs. The two edge states are nearly overlapping, and are found to be
localized at left boundary while at right boundary the wave function only has
a very small weight. This is consistent with our analysis. If we perform a
similarity transformation on this model to make
$|\beta_{4}\beta_{5}/\rho^{2}|=1$, i.e., $\rho=|\beta_{4}\beta_{5}|^{1/2}$,
then this model described by $H(\beta)$ would be transformed to one described
by $H(\rho\beta)$ whose wave function $\widetilde{\psi}(j)=\rho^{-j}\psi(j)$
and so its ratio $|\mathcal{B}/\mathcal{A}|\approx 1^{L}=1$. In our analysis
the two TPESs in this situation would be localized at both ends of the open
chain, consistent with the numerical results as shown in Fig.6 (d).

Figure 6: Space distributions of the edge states with $\nu=-1$ and
$\lambda=0.1$. (a) Edge state $\psi_{\mathcal{L}}$. (b) Edge state
$\psi_{\mathcal{R}}$. (c) Numerical matrix diagonalization result of the zero-
energy modes, and (d) the corresponding result after a similarity
transformation with $\rho=|\beta_{4}\beta_{5}|^{1/2}$, where data on
sublattice A(B) are colored in blue(red). We have 800 lattice sites since here
$L=200$ and $N=2$.

Meanwhile, notice that a negative $\nu=-1$ means $|\mathcal{L}\rangle$ only
has $\phi^{-}$ components, while $|\mathcal{R}\rangle$ only has $\phi^{+}$
ones, indicating that $|\mathcal{L}\rangle$($|\mathcal{R}\rangle$) is a
left(right)-localized A(B)-sublattice edge state. However, for a positive
$\nu=1$, the $|\mathcal{L}\rangle$($|\mathcal{R}\rangle$) state would be a
left(right)-localized B(A)-sublattice edge state, as mentioned in the main
text. This is because that when $\nu=+1$, all nonzero components in
$|\mathcal{L}\rangle$($|\mathcal{R}\rangle$) are $\phi^{+}$($\phi^{-}$)
instead of $\phi^{-}$($\phi^{+}$).

Figure 7: Space distributions of the zero-energy modes obtained by the
numerical matrix diagonalization, with $\nu=1$, where data on sublattice A(B)
are colored in blue(red). (a) $\lambda=-0.6$. (b) $\lambda=-0.9$. Here
$L=200$.

In Fig.7, we show the case of $\nu=1$ with $\lambda=-0.6$ and $\lambda=-0.9$
respectively. The eight $\beta$ values for the two parameters are ordered as
usual as follows:

$\begin{split}&|\beta_{1}|=|\beta^{+}_{1}|<|\beta_{2}|=|\beta^{+}_{2}|<|\beta_{3}|=|\beta^{-}_{1}|<|\beta_{4}|=|\beta^{+}_{3}|<1<|\beta_{5}|=|\beta^{-}_{2}|<|\beta_{6}|=|\beta^{+}_{4}|<|\beta_{7}|=|\beta^{-}_{3}|<|\beta_{8}|=|\beta^{-}_{4}|,\\\
&|\beta_{1}|=|\beta^{+}_{1}|<|\beta_{2}|=|\beta^{+}_{2}|<|\beta_{3}|=|\beta^{+}_{3}|<|\beta_{4}|=|\beta^{-}_{1}|<1<|\beta_{5}|=|\beta^{+}_{4}|<|\beta_{6}|=|\beta^{-}_{2}|<|\beta_{7}|=|\beta^{-}_{3}|<|\beta_{8}|=|\beta^{-}_{4}|,\end{split}$
(26)

where the former corresponds to $\lambda=-0.6$ while the latter corresponds to
$\lambda=-0.9$. When $\lambda=-0.6$,
$|\mathcal{B}/\mathcal{A}|\approx|\beta_{4}\beta_{5}|^{L/2}=|\beta^{+}_{3}\beta^{-}_{2}|^{L/2}$
and when $\lambda=-0.9$,
$|\mathcal{B}/\mathcal{A}|\approx|\beta_{3}\beta_{6}|^{L/2}=|\beta^{+}_{3}\beta^{-}_{2}|^{L/2}$.
In both cases, $|\mathcal{B}/\mathcal{A}|\approx 0.997^{L}$. In our numerical
calculation, since $L=200$, $|\mathcal{B}/\mathcal{A}|\approx 0.74$, which
makes $|\mathcal{A}|$ slightly larger than $|\mathcal{B}|$. So each TPES is
localized at both boundaries of the open chain. But when $L\rightarrow\infty$,
only the left-localized components would be remaining. By setting $L=400$ we
redo the numerical calculation of the edge states in Fig.6 (d) and Fig.7 (a),
and show the corresponding results in Fig.8. In Fig.8 (a), the space
distribution is nearly unchanged since $|\mathcal{B}/\mathcal{A}|\propto 1$ is
also kept unchanged while in Fig.8 (b), the left-localized components are even
more dominant than the right-localized ones since the ratio
$|\mathcal{B}/\mathcal{A}|$ in this case becomes smaller. These results show
that the predicted ratio between the coefficients $\mathcal{A}$ and
$\mathcal{B}$ is quite reasonable.

Figure 8: (a) The same as Fig.6 (d), and (b) the same as Fig.7 (a), except
that $L=400$.

Now we choose $\lambda=1$ to study the case of $\nu=-2$. The eight $\beta$
values are ordered as follows:

$|\beta_{1}|=|\beta^{-}_{1}|<|\beta_{2}|=|\beta^{-}_{2}|<|\beta_{3}|=|\beta^{-}_{3}|<|\beta_{4}|=|\beta^{-}_{4}|<1<|\beta_{5}|=|\beta^{+}_{1}|<|\beta_{6}|=|\beta^{+}_{2}|<|\beta_{7}|=|\beta^{+}_{3}|<|\beta_{8}|=|\beta^{+}_{4}|.$
(27)

It is found that Rank$\mathcal{M_{L}}=$Rank$\mathcal{M_{R}}=2$, namely, they
have been maximized so either $\mathcal{M_{L}}|\mathcal{L}\rangle=0$ or
$\mathcal{M_{R}}|\mathcal{R}\rangle=0$ has two independent solutions, which
can be labeled as $|\mathcal{L}_{1}\rangle$, $|\mathcal{L}_{2}\rangle$ and
$|\mathcal{R}_{1}\rangle$, $|\mathcal{R}_{2}\rangle$, respectively. According
to degenerate perturbation theory, let’s solve
$\left(\begin{array}[]{cc}\mathcal{M_{L}}&\mathcal{M_{LR}}\\\
\mathcal{M_{RL}}&\mathcal{M_{R}}\end{array}\right)|\psi\rangle=0$ which is
equivalent to $M|\varphi\rangle=0$, with $M$ being the $2|\nu|\times 2|\nu|$
matrix given by:

$\begin{split}M=&\left(\begin{array}[]{cccc}0&M_{LR}\\\ M_{RL}&0\\\
\end{array}\right)\\\
=&\left(\begin{array}[]{cccc}0&0&\langle\mathcal{L}_{1}|\mathcal{M_{LR}}|\mathcal{R}_{1}\rangle&\langle\mathcal{L}_{1}|\mathcal{M_{LR}}|\mathcal{R}_{2}\rangle\\\
0&0&\langle\mathcal{L}_{2}|\mathcal{M_{LR}}|\mathcal{R}_{1}\rangle&\langle\mathcal{L}_{2}|\mathcal{M_{LR}}|\mathcal{R}_{2}\rangle\\\
\langle\mathcal{R}_{1}|\mathcal{M_{RL}}|\mathcal{L}_{1}\rangle&\langle\mathcal{R}_{1}|\mathcal{M_{RL}}|\mathcal{L}_{2}\rangle&0&0\\\
\langle\mathcal{R}_{2}|\mathcal{M_{RL}}|\mathcal{L}_{1}\rangle&\langle\mathcal{R}_{2}|\mathcal{M_{RL}}|\mathcal{L}_{2}\rangle&0&0\\\
\end{array}\right),\end{split}$ (28)

where $|\varphi\rangle=(l_{1},l_{2},r_{1},r_{2})^{T}$. So
$|\psi\rangle=(l_{1}|\mathcal{L}_{1}\rangle+l_{2}|\mathcal{L}_{2}\rangle,r_{1}|\mathcal{R}_{1}\rangle+r_{2}|\mathcal{R}_{2}\rangle)^{T}$.
We have four solutions:
$|\psi^{\pm}_{m}\rangle=(\mathcal{A}_{m}|\mathcal{L}_{m}^{\prime}\rangle,\pm\mathcal{B}_{m}|\mathcal{R}_{m}^{\prime}\rangle)^{T}$,
where $m=1$, $2$. Here each normalized $|\mathcal{L^{\prime}}_{m}\rangle$
($|\mathcal{R^{\prime}}_{m}\rangle$) is the linear superposition of
$|\mathcal{L}_{1}\rangle$ and
$|\mathcal{L}_{2}\rangle$($|\mathcal{R}_{1}\rangle$ and
$|\mathcal{R}_{2}\rangle$). $|\psi^{+}_{m}\rangle$ is also connected with
$|\psi^{-}_{m}\rangle$ by chiral symmetry:
$\mathcal{S}|\psi^{\pm}_{m}\rangle=-|\psi^{\mp}_{m}\rangle$. From the ratio
equation

$\frac{\mathcal{B}_{m}}{\mathcal{A}_{m}}=\frac{\sqrt{\langle\mathcal{R^{\prime}}_{m}|\mathcal{M_{RL}}|\mathcal{L^{\prime}}_{m}\rangle}}{\sqrt{\langle\mathcal{L^{\prime}}_{m}|\mathcal{M_{LR}}|\mathcal{R^{\prime}}_{m}\rangle}},\quad
m=1,2$ (29)

we get
$|\mathcal{B}_{1}/\mathcal{A}_{1}|\approx|\beta_{4}\beta_{5}|^{L/2}=|\beta^{-}_{4}\beta^{+}_{1}|^{L/2}=0.11567$
and
$|\mathcal{B}_{2}/\mathcal{A}_{2}|\approx|\beta_{3}\beta_{6}|^{L/2}=|\beta^{-}_{3}\beta^{+}_{2}|^{L/2}=0.00084$,
with $L=200$. The numerical result for one pair of the TPESs is shown in Fig.9
(a), while the other pair(data not shown) shows analogous behavior in this
case. We found that the predicted ratio is consistent very well with the
numerical matrix diagonalization result. As before, if we apply the similarity
transformation to this system, i.e., $H(\beta)\rightarrow H(\rho\beta)$, the
ratios become
$|\mathcal{B}_{1}/\mathcal{A}_{1}|\approx|\beta_{4}\beta_{5}/\rho^{2}|^{L/2}$
and
$|\mathcal{B}_{2}/\mathcal{A}_{2}|\approx|\beta_{3}\beta_{6}/\rho^{2}|^{L/2}$.
In Fig.9 (b), we set a suitable $\rho$ to make one of the ratios close to $1$,
so the corresponding pair of the TPESs would be localized at both ends of the
open chain. In this situation, the total number of TPESs would be three. In
Fig.9 (c)-(d), we set $\rho=2$ and $\rho=1/2$ respectively. In either case,
when $L$ is large, among $\mathcal{B}_{m}$ and $\mathcal{A}_{m}$, one would
always be much larger than the other, and the smaller one can be viewed as be
approaching zero, implying that the number of the TPESs would be equal to
$|\nu|=2$.

Figure 9: Space distributions of the edge states with $\nu=-2$ and
$\lambda=1$. (a) One pair of the nearly overlapping zero-energy modes obtained
by the numerical matrix diagonalization, and the corresponding result after a
similarity transformation with (b) $\rho=|\beta_{4}\beta_{5}|^{1/2}=0.989$,
(c) $\rho=2$ or (d) $\rho=1/2$, respectively. Here $L=200$.

Here we see that we can always choose a sufficiently small or large $\rho$ to
perform the similarity transformation so that the number of the TPESs can be
reduced to the absolute value of the topological invariant $\nu$ but is
definitely not less than $|\nu|$. If all ratios
$|\mathcal{B}_{m}/\mathcal{A}_{m}|$ share the same relationship, one can
naturally choose a common $\rho$ to make all edge states be localized at both
ends, so the number of the TPESs reaches its maximum value, $2|\nu|$. However,
in most cases of non-hermitian systems, this is usually unachievable, so for a
general non-hermitian open chain with chiral symmetry, the number of the TPESs
would be varying from $|\nu|$ to $2|\nu|$.

## References

  * [1] Deguang Wu, Jiao Xie, Yao Zhou, and Jin An.  Connections between the open-boundary spectrum and the generalized brillouin zone in non-hermitian systems.  Phys. Rev. B, 105:045422, Jan 2022.