Topological edge states in two-gap unitary systems: A transfer matrix approach

We construct and investigate a family of two-band unitary systems living on a cylinder geometry and presenting localized edge states. Using the transfer matrix formalism, we solve and investigate in details such states in the thermodynamic limit. Analitycity considerations then suggest the construction of a family of Riemman surfaces associated to the band structure of the system. In this picture, the corresponding edge states naturally wind around non contractile loops, defining by the way a topological invariant associated to each gap of the system.


Introduction
First discovered in the context of the quantum Hall effect [1,2,3], boundary states turn out to be the hallmark of topological properties that can emerge in any dimensions and for various symmetry classes [4,5]. It was also realized that these topological properties can be achieved in miscellaneous physical systems beyond solids [6,7] leading for instance to the discovery of chiral edge states in cold-atom [8], electromagnetic [9,10,11] and acoustic [12,13,14,15] lattices.
These topological properties characterize the bulk bands and are encoded by a topological invariant whose integer value cannot change unless bands touch, or equivalently, unless the gap closes. For instance, in two dimensions and in absence of particular symmetry, this topological invariant is the first Chern number [16,17]. As the edge states live in the gaps spectrum, it follows that they cannot be removed or added unless a topological transition of the bulk bands happens when the gap closes. This gives a topological robustness to the edge states.
Recently, similar topological properties have been found in unitary systems, namely physical systems whose behavior is described by a unitary operator rather than a Hermitian operator. It follows that their spectrum is periodic, unlike an energy spectrum which is bounded. Among them are the Floquet systems that are periodic in time [18,19]. These dynamical systems, such like periodically driven solids [20,21], shaken cold-atom gases [22], photonic lattices [23,24] or discrete-time quantum walks [25] are fully described by their (unitary) time evolution operator. Importantly, it was shown that beyond the analogy with equilibrium (Hermitian) systems, they can exhibit edge states whereas the usual bulk topological invariants vanish. This is understood as an exotic topological property characterized by a novel invariant assigned to a gap (rather than a band) and that accounts for the time evolution over a period [26,27].
However, there are unitary systems which are not Floquet systems but still exhibit edge states. This is the case of photonic and microwave networks which can instead be described by (unitary) scattering matrices [28,29,30]. How to describe the topological origin of these edge states? Can it still be related to a bulk property?
To answer these questions we construct, in section 2, a generic model of a two-gap unitary system on a cylinder geometry. This model achieves the specific situation for which the topological invariants of the bands (namely the Chern number in our case) vanish. Then, in section 3, we apply the transfer matrix method to investigate the appearance of edge states and discuss their simultaneous existence in the two gaps (thus guarantying a vanishing Chern number of each band). In its seminal paper, Hatsugai showed that the transfer matrix approach provides a deep understanding of the topological nature of the edge states in the quantum Hall phase [31]. In particular, the transfer matrix allows one to focus directly on the gaps where the edge states live rather than the bands only. We follow the same strategy and adapt this method to the unitary case in section 4. Our analysis reveals the underlying topologically non-trivial structure of the edge states and justifies the definition of a topological invariant assigned to a gap instead of a band. Finally, several examples of physical systems ruled by the present model are discussed in section 5.
2 Two-gap unitary models with topological edge states

Heuristic construction on a finite size lattice
We construct heuristically a two-gap unitary model for a strip geometry that exhibits topologically protected edge states. For simplicity, we shall treat the cases of 0 or 1 edge state, but the generalization to several edge states is straightforward. To this end, we consider a system with two degrees of freedom -that we refer to as A and B -in the cylinder geometry such that the (dimensionless) quasi-momentum k ∈ U (1) is a well-defined continuous parameter in the periodic direction, whereas the lattice remains finite in the other one, as sketched in figure 1 (a). A state |Ψ(k) of the system is then given by a 2N -component vector (A 1 (k), B 1 (k), . . . , A N (k), B N (k)) T , which, for the scope of this study, is an eigenvector of a unitary matrixŨ We would like the phases (k) ∈ S 1 of the eigenvalues ofŨ (2N ) to display two-gaps and two edge states lying in these gaps. Formally, the simplest way to obtain such a phase spectrum is to impose a quasi-diagonal form forŨ (k) such as where the bulk part iσ y yields two flat bands = ±π/2 (in the range ∈] − π, π]), whereas e ik and e −ik guaranty the existence of propagating modes at the boundaries A 1 and B N respectively, (see figure 1 (b)). The phase spectrum (k) ofŨ 0 (k) can Figure 1: (a) Sketch of the cylindrical geometry for the strip of width N . The two arrows represent the chiral edge states at the boundary n = 1 (purple) and N (red). Phase spectrum (k) obtained by a direct diagonalization ofŨ (k) with N = 10 for θ = 0 (b) and θ = π/8 (c).
then be obtained as a particular case of a more general unitary matrix where U µ ∈ U (2) and 1 is a scalar. One gets explicitlỹ where, for concreteness, we choose In Eq. (4), Λ is a 2 × 4 matrix coupling components A and B at sites n and n + 1 for 1 ≤ n < N (hence describing the bulk properties), whereas first and last lines are constraints for A and B at sites 1 and N (hence yielding the boundary conditions). This model allows us to study the fate of the edge states when varying θ. First, note thatŨ (k) corresponds toŨ 0 (k) for θ = 0, up to a unitary transformation that does not change the spectrum. Then, for θ = 0, the bulk bands acquire a dispersion (figure 1 (c)) and eventually touch at (k = 0, = 0) and (k = π, = π) for θ = π/4, thus implying a transition toward a gapped phase with no edge state (π/4 < θ < π/2). A similar matrixŨ (k) was derived to investigate the propagation of electromagnetic modes in arrays of optical resonators [28,29], and can also me adapted to describe discrete-time quantum walks for a spin-1/2 particle or either a time-dependent tight-binding model, as we discuss in section 5.
In the following, we investigate the topological properties of this model with the transfer matrix formalism.

Transfer matrix formalism
The transfer matrix formalism is a standard method to tackle a large variety of problems. The starting point consists in relating the wave function amplitudes on adjacent sites of a lattice by a matrix T . In our case, this translates as From the general form of the matrixŨ (k) in Eq. (4), one can infer the relation Note that the relation (7), valid for 1 ≤ n < N describes the bulk; it does not take into account the boundary conditions (n = 1, N ) that will be given in section 3.2.1. We can then deduce the expression of T which, for the specific parameters (5), reads By construction, the matrix T l describes how the amplitude of a state evolves from the edge n = 1 to a site n = l of the bulk. This information is then encoded into the eigenvalues λ ± ( ) of the transfer matrix T . In particular, either |λ ± (k, )| = 1, which corresponds to an oscillatory (or delocalized) bulk state at (k, ); or one of the two eigenvalues satisfies |λ(k, )| < 1 -note that λ + = λ −1 − since det T = 1and no delocalized state can appear i.e. there is a band gap. The eigenvalues of T being solutions of the characteristic polynomial λ 2 − tr T λ + det T = 0, they read where In other words, for the model (5), the eigenvalues of the transfer matrix are fully determined by its trace. Finally, the two eigenvectors v − and v + , associated respectively to the eigenvalues λ − and λ + of the matrix T = {T ij }, can be written as In the following, we show how the projected bulk bands, the edge states and their topological properties can be inferred from these eigenvalues and eigenvectors.

Band structure from the transfer matrix
The transfer matrix formalism allows one to solve explicitly the initial problem (1), namely to reconstruct the band structure and the corresponding states. We focus on the thermodynamic limit N → ∞ of the infinite cylinder: in that case the bulk bands become a continuum, whereas the edge states persist.

Bands and gaps
The bulk bands (k) correspond to regions of the diagram (k, ) for which the eigenvalues of the transfer matrix satisfy |λ ± (k, )| = 1. For the model we consider, this region is determined by the solutions of f 2 (k, ) 1 according to Eq. (10). This can be found explicitly by using the expression of f (k, ) in Eq. (11). The result is shown in figure 2 and compared to a direct diagonalization ofŨ . Consistently, the gaps correspond to regions of the diagram (k, ) for which the eigenvalues of the transfer matrix satisfy |λ ± (k, )| = 1, that is when f 2 (k, ) > 1. This inequality has two solutions f (k, ) < −1 and f (k, ) > 1 yielding two different domains in the diagram (k, ). Since these two domains are separated by a band region (f 2 1), they thus correspond to the two gaped regions that we will refer to as gap g − when f < −1 and g + when f > 1. The gaps g − and g + are respectively delimited by f = −1 or f = +1. We therefore end up with a criteria that distinguishes the bands and the two gaps directly from the trace of the transfer matrix The explicit dispersion relation of the two gap edges are then inferred from Eq. (11) G −,t (k) = arccos (sin θ cos θ(cos k + 1)) , G −,b (k) = − arccos (sin θ cos θ(cos k + 1)) G +,b (k) = arccos (sin θ cos θ(cos k − 1)) , G +,t (k) = − arccos (sin θ cos θ(cos k − 1)).
and plotted in Fig. 2.

Boundary conditions
So far, the bulk part ofŨ have been encoded within the transfer matrix. To characterize the existence of edge states, one also need to specify the boundary conditions. Whereas in Hermitian systems, a standard boundary condition on a lattice consists in imposing the vanishing of the wavefunction, such a procedure does not apply here since it is not compatible with unitarity. As for the bulk transfer matrix, the boundary conditions are instead inferred from the original eigenvalues problem onŨ . From Eq. (4), we relate for each edge n = 1, N the amplitudes A n and B n with the phase and the coefficients of U 1 as A geometric way to reformulate these boundary conditions is to introduce two boundary vectors b 1 and b N that are parallel to where we have defined These two boundary vectors are related by the transfer matrix as This equation constrains the allowed values of couple (k, ) for finite size N , thus yielding the corresponding solution of the initial problem 1 (1). Next, we impose that in the thermodynamic limit (N → ∞), b N does not depend on N anymore. Thus, in that limit, we can write T N b 1 = T N −1 b 1 so that b 1 becomes an eigenvector of the transfer matrix, given by (12). Then, by decomposing b 1 as b 1 = αv + + βv − , we end up with (where we have set λ + = λ for commodity). Let us assume in addition that |λ| < 1.
Then |λ| N −1 vanishes and |λ −(N −1) | tends to infinity. It follows that β = 0 and thus b 1 is proportional to v + . We can show in the same way that b N is proportional to v − when N → ∞. This result will allow us to fully characterize the edge states in the next section. Besides, a 1 and a N beeing also eigenstate of T in that case, the corresponding state |Ψ of the initial problem is exponentially decreasing or increasing with typical length ln λ, namely it is localized at one edge of the system.

Existence and dispersion relation
Concretely, in order to know which one of the λ ± 's satisfies |λ| < 1, one has to fix f , which amounts to fix the gap. Let us for instance focus on the gap g − (i.e. f < −1). Then, it is clear that |λ + | > 1 and thus |λ − | < 1. Therefore, in the gap g − , an edge state at n = 1 is given by the proportionality relation between the boundary vector b 1 and the eigenvector v − , whereas the edge state at n = N is given by the proportionality relation between the boundary vector b N and the eigenvector v + . The opposite situation occurs in the gap g + , for which f > 1, and thus |λ − | > 1 and |λ + | < 1, meaning that an edge state localized at the edge n = 1 (n = N ) is now associated with the eigenvalue λ + (λ − ) and the eigenvector v + (v − ). The existence of an each edge state and its dispersion relation is therefore given independently for each gap and for each boundary n = 1, N by the zeros of the function as summarized in table 1. Let us treat the case of the edge state n = 1. Each of the two equations I 1,± (k, ) = 0 can be split into two equations corresponding to the real part and the imaginary part of I 1,± (k, ). These two equations can be treated simultaneously, and after some algebra, one gets Whereas the explicit dispersion relation of the edge state localized at edge n = 1 is straightforwardly obtained by inverting Eq. (20), its domain of validity in k is constrained by the equation (21) which depends on the gap g ± through the contribution ± of the square root. Concretely, the edge state only exists in the region (k, ) satisfying by Q(k, ) > 0 for the gap g + and in the region Q(k, ) < 0 for the gap g − . A similar procedure can be performed for the other edge so that the spectrum of the strip in the thermodynamic limit is recovered. As displayed in figure 3 this procedure provides the dispersion relation of the edge states (see also figure 6 (a) for the reconstruction of the full spectrum). Moreover, it shows that edge states cannot exist if the sign of Q does not change in the gap. It reveals the importance of the role of the branch ± of the square root in Eq. (21). This is a crucial point to understand the topological robustness of the edge states that will be developed in section 4.

Gaps correspondence
Before discussing the topological nature of the edge states, we emphasize that they always appear simultaneously in the two gaps for any generic unitary matrixŨ of the form of Eq. (3) (provided the spectrum ofŨ is gapped). This is a particularly interesting property of unitary systems since the Chern number of the Bloch bands of the periodic system in both x and y directions is guaranteed to vanish in this case [26]. However, the conditions to achieve such a remarkable topological property are not clear and remain an open question. Here, we point out that the transfer matrix formalism is particularly useful to reveal a correspondence between the gaps that constrains the existence of an edge state simultaneously in each gap.
To make explicit such a property, we define a map υ acting in parameters space υ : (k, ) → (υk, υ ) that relates an edge state at the boundary n = 1 that exists in the gap g + to an edge state localized at the same boundary but existing in the gap g − . We aim at determining under which conditions I n,+ (k, ) and I n,− (υk, υ ) vanish simultaneously. By swapping the gaps, the function f changes sign as f (υk, υ ) = −f (k, ) so that the eigenvalues become It is clear that a map υ acting on the eigenvalues of the transfer matrix as in Eq. (22) actually changes the sign of the trace of the transfer matrix. This allows us to express such a constraint on the transfer matrix itself as with σ z the standard Pauli matrix and where Υ is allowed to be either the complex conjugation operator κ or the identity. It is then easy to check that, for the model (5), the map υ : (k, ) → (π − k, π − ) (defined modulo 2π) satisfies Eq. (23) with Υ = κ, so that Eq. (22) is satisfied as well. It follows that an eigenvector v ± (k, ) of the transfer matrix T (k, ) (see Eq. (12)) is transformed as since T 12 (υk, υ ) = T * 12 (k, ) and λ ± (k, ) = λ * ± (k, ) for the model (5). A relation of proportionality between I 1,± (k, ) and I 1,∓ (υk, υ ) can finally be inferred provided that the boundary vector b 1 (k, ) -being itself an eigenvector of the transfer matrix in the thermodynamic limit -is also transformed as (as it can be checked explicitly from Eqs. (8) and (12)). Indeed, according to Eqs. (24) and (25) one gets and thus I 1,± (k, ) and I 1,∓ (υk, υ ) vanish simultaneously. This shows that if an edge state localized at edge n = 1 exists in the gap g + (defined by the zeros of I 1,+ ), then so does an edge state localized at the same edge in the gap g − (defined by the zeros of I 1,− ).
4 Topological property of the edge states from the transfer matrix

Motivations
We would like to specify the topological nature of the edge states regardless the topological properties of the bands which are anyway always characterized by a vanishing Chern number as discussed above. A natural approach then consists in assigning a winding number to the edge states themselves [32,33,34]. When θ = 0, thenŨ (k) reduces to a block-diagonal matrix where the contributions of the bulk and of each edge can be distinguished. In that case, we notice that each edge term h ± (k) ≡ e ±ik is a U (1)−valued scalar so that its winding is well defined [26] and we get 1 How to construct such a quantity for arbitrary θ wheñ U (k) is not block-diagonal anymore? How to assign a specific gap to this winding? It is clear that the dispersion relations of the edge states obtained before will be a key ingredient. However they are not defined for every k, not even periodic when considering one edge, and multivalued when considering both edges in one relation. Hence the winding number of such quantity is ill defined.
A possibility to overcome this problem is to consider instead the functions I n,± (k, ), whose zeros fully define the edge states dispersion relations. These functions depend explicitly on the eigenvalues λ ± of the transfer matrix that contain the square root f 2 − 1. Therefore, their first derivative diverge around f = ±1, that is when (k, ) touch the bands edges. This divergence of the first derivative appears in the wave function too, which -remembering that the boundary vectors are also eigenvectors of the transfer matrix -is proportional to (A n , B n ) T ∝ λ n−1 ± b 1,N . Thus, to properly define a topological invariant, one first needs to find a way to obtain smooth solutions for the wave-functions when changing the branch of the square root ± f 2 − 1 → ∓ f 2 − 1.
The natural tool to deal with this issue is the theory of Riemann surfaces, as used by Hatsugai in 1993 to correctly define the winding number of the edge states in the (Hermitian) Harper model [31]. To construct the Riemann surfaces associated to our problem, we adapt the approach by Hatsugai to the unitary case where the eigenenergies are replaced by exp(−i ) which are U (1)-valued. Substituting f by its expression, this approach suggests to consider quantity µ 2 k = cos k − cos sin θ cos θ 2 − 1 for → z ∈ C rather than the square root itself which is not analytic in the complex plan. Such equation can be rewritten where the expression of the four C-valued functions φ j (k) can be obtained explicitly.
Here we just need to know that φ j (k) ∈ U (1), φ j (k) = φ j (k) for j = j (the gaps are not closing, except at θ = π/4 which has to be avoided) and φ j (k) is smooth in k as soon as 0 < θ < π/2 (the bands are not degenerated to a line). Note that when e −i = φ j (k) one has µ k = 0 that is f 2 = 1: this corresponds to the gap/band edge. The construction of the Riemann surface consists in two steps. One is the standard construction of an elliptic curve and its geometric interpretation as a compact surface obtained by gluing two Riemann sheets corresponding to the two branches ± of the square root. The second step, specific to unitary systems, accounts for the periodic structure of the phases through the dependence in exp ±i (and its powers) appearing in µ 2 k .

Elliptic curves
From now on, we consider a fixed momentum k and forget about it in the notations. We start with the standard construction of elliptic curves by considering the Riemann surface with the φ j satisfying the previous assumptions. This is an almost standard construction [35]. First, this is a compact Riemann surface since it can be covered by open subset that are homeomorphic to open subset in C, with holomorphic transition function on the intersections. The singularities that appear when Z → 0 and ∞ are not essential because the corresponding neighborhood can also be identified with open subsets of C. There is a nice interpretation of this curve as a double covering of C, or its compactified version [35], as we will see in the next section.
The complex square root cannot be defined analytically over the whole complex plane C, nor on its compact version, namely the Riemann sphere S = C ∪ ∞ [35]. A natural branch cut on C (or S) is defined such that W (Z) < 0. We know from Eq. (29) that W (Z) is real for Z ∈ U (1) and that the sign of W changes at φ i where µ = 0. By construction, the φ i are located on the equator of the Riemann sphere. We chose the order of the φ i 's such that W (Z) < 0 for Z ∈ [φ 1 , φ 2 ] and Z ∈ [φ 3 , φ 4 ] (this notation is ambiguous in U (1) but that does not matter in the following). As a result W (Z) is real and positive for Z on the equator of each Riemann sphere as well, but outside the branch cuts. To get an analytic structure for R 0 , consider two copies of such cut Riemann spheres S + and S − and set the convention that µ = + √ W > 0 (respectively − √ W < 0) on such a region of the Riemann sphere S + (S − ). Hence, by construction of these two copies, one travels from one sphere to the other by crossing a branch cut, such that µ goes smoothly from positive to negative values [35]. The square root is now analytic over R 0 (instead of S) which is a torus of genus 1 (instead of a sphere). The Riemann surface R 0 is locally homeomorphic to the complex plane, but not globally since the topology is different. This is the price to pay to have smooth functions. But this also gives a direct geometrical interpretation of the winding of the edge states [31], as we will see.

Punctured torus for the unitary problem
We would like to apply the general theory discussed above for Z = exp (−i ). To do so we first need to extend the real variable to the complex plane. The phases being defined modulo 2π, a natural extension of their domain of definition is the complex cylinder C = {z = + iη | ( , η) ∈ S 1 × R}, as depicted in figure 4. Then we define the map ϕ : that sends this cylinder to the complex plane by preserving the circles (see Fig. 4). However, note that this map is not analytic since it has two essential singularities at η → +∞ and η → −∞, mapped respectively to ∞ and 0 in C. This is specific to unitary models where the phase is U (1)-valued whereas such singularities do not appear when doing E → z for real energy E of Hermitian systems [31]. These singularities will stay all along the construction and will be actually necessary. Indeed, the image of C by ϕ is the Riemann sphere S * punctured of two singular points 0 and ∞. This surface is not compact anymore, but outside this two points the function is still analytic such that around the image of U (1) (the equator in this picture), phase terms exp(−i ) can be extended in an analytic way.
One can now apply the general method of the previous section to construct the Riemann surface by replacing the two cut Riemann spheres S ± by two punctured cut Riemann spheres (S * ) ± . The Riemann surface R, defined as some pull-back of R 0 by ϕ, is obtained by gluing (S * ) + and (S * ) − together by the branch cuts W (e −iz ) < 0 which correspond to the regions of the two bands (at fixed k) as shown in Fig. 5. These bands, together with the gaps g ± , constitute the equator of each punctured sphere as they originate from the real part of the complex variable z (blue circle in figure 4). All the quantities used, besides the square root, only involve polynomials in e −iz and e iz that are perfectly smooth on R since they are so on each copy of S * . The square root of µ 2 k is an analytic function on the Riemann surface R and so are the eigenvalues of the transfer matrix. In particular, one can write λ = f + µ that corresponds to λ + when e −iz is in a gap and the corresponding µ is positive, and similarly for λ − with µ negative. In this picture λ contains both square roots and one can go analytically from one to the other using the continuation in the complex numbers. Consequently the edge states of the system are now described in an analytic way on R, and the full square root is given by one single formula instead of two. Geometrically, the Riemann surface R is a torus with four punctured points (see Fig. 5). Each gap is present in two copies that generates loops along the torus. Indeed, the two copies of g + generate the outer loop + (in green) and the two copies of g − generate the inner loop − (in blue). Importantly, due to the four singular points 0 + , ∞ + , 0 − and ∞ − that cannot be crossed, these two loops are not equivalent (or homotopic) as in the standard torus since they cannot be deformed one to each other without crossing a singular point. Hence the Riemann surface R keeps track of the two distinct gaps g − and g + of the unitary problem. Finally note that such relative position (interior/exterior) is completely arbitrary and just depends on the way we draw the construction of the torus. Indeed the two configurations of torus are topologically equivalent by "twisting" the full torus.
The construction of the Riemann surface was performed at fixed quasi-momentum k. When varying k, the size of the gaps changes so that one obtains a family of Riemann surfaces R k . Still, the topology remains the same as these surfaces are all homotopic one to another (as long as the gaps do not close): one can deform all this family to the same torus R for all k.

Winding number
When k spans S 1 , an edge state may cross a gap by moving from one band to the other one, whereas the other edge state, located at the other boundary, crosses the gap in the opposite direction. On the Riemann torus R, these two edge states span one of the two loops ± and can thus be qualified as topological : the winding of their pair cannot change value unless the gap closes. This defines a topological invariant Figure 5: Construction of the Riemann surface R by gluing two punctured Riemann spheres (S * ) + and (S * ) − . The branch cuts, in red, correspond to the bands. They and are delimited by the φ i according to Eq. (29) and separate the two gap g − (blue) from the gap g + (green). For convenience, we have flipped one sphere in a way such that when arriving at a branch cut we are going to the other surface, stay to the same hemisphere (north or south), as illustrated by the arrows. The gaps form two non-contractile loops − and + on R which are non-homotopic one to each other because of the essential singularities ∞ ± and 0 ± . for each gap that counts algebraically the number of chiral edge states. As noticed by Hatsugai [31], this winding number is nothing but the intersection number of the curve spanned by the edge state dispersion relation (k) (after having identified all the Riemann surfaces R k ) and some "vertical line" of the torus, e.g. one of the branch cuts.
To compute explicitly this number W , it is particularly convenient to deal with one continuous function˜ (k) for the dispersion relation of the pair of edge states rather than a multi-valued function as it is the case when the two edge states cross in the gap (see figure 6 (a)). Such a single-valued function can always be obtained by shifting one of the two edge state's dispersion relation, as shown in figure 6 (b). This is performed by adding a phase-shift to the boundary vector (e.g. b 1 ) coefficients as τ 1 → τ 1 e iϕ and ρ 1 → ρ 1 e iϕ so that the unitarity of U 1 (and thenŨ ) is preserved 2 . Importantly, the transfer matrix is not affected by this transformation, hence both the bulk bands and the construction of the Riemann surface remain the same.
It is enough to focus on the gap g − only. The single-valued dispersion relatioñ − (k) of the pair of edge states is smooth and periodic, so that the winding number of exp(−i˜ − (k)) always vanishes. However in the Riemann surface picture it can be used to wind around nontrivial loops. In this gap,˜ − (k) touches the bands in k * , for G −,t/b given in (14). That way, one gets µ k > 0 for k ∈]k * 1 , k * 2 [ and µ k < 0 for k ∈]0, k * Figure 6: (a) Bulk bands and edge states obtained in the thermodynamic limit with the transfer matrix method for θ = π/8. (b) The boundary state at edge n = 1 (blue) has been shifted (ϕ = π) so that the two edge states touch at k * 1 and k * 2 . The dispersion relation of a pair of edge states (in each gap) is single-valued. (c) Appearance of trivial edge states for θ = 3π/8 when complex phases are added to the model (5).
where D − (k) measures the relative distance between the dispersion relation of the edge state and a bound of the gap g − . Let us choose the lower bound as a conven- The function F − is well defined as soon as the gap does not close (θ = π/4). It is 2π-periodic in k, continuous in k * 1 (since it is 0) but not continuous at k * 2 . However the function is continuous, and has non trivial winding number Since F − is periodic, then, if F − were also continuous, its winding number would be 0. However F − is only piecewise continuous, such that where we have dropped the gap index because of the gap correspondence discussed in section 3.2.3. From an effective point of view, the computation of W [w] amounts to count the number of times ∂ k µ k changes sign at points k * i such that F (k * i ) is discontinuous, that is D(k * i ) = 1. On can thus write W [w] as an intersection number 3 Note that the other convention would reverse the sign of the winding number.
where the derivative ∂ k µ k is well defined on the Riemann surface. In practice, when going from a lower edge state (µ k < 0) to an upper one (µ k > 0) while increasing k, the derivative of µ k is positive and the winding number increases by 1, and conversely decreases by 1 when going from the upper to the lower Riemann sheet. It is clear that this analysis is still valid for several edge states and without the gap correspondence; the invariant W thus counting the number of topological edge states in a given gap. Besides, it also distinguishes topological edge states from nontopological "accidental" edge states that may appear for certain set of parameters, as illustrated in figure 6 (c). Indeed these edge states are associated to contractile loops on R and thus do not benefit of any topological robustness. Their winding (36) is clearly zero.
Finally note that the loops around the essential singularities of R might be related to previous invariant W . For example, a loop around 0 − on Figure 5 is homeotopic to the one passing by φ 1 → φ 2 → φ 3 → φ 4 → φ 1 and traveling along left solid red curve, dashed blue part of − , right solid red curve and solid green part of + , respectively. This loop corresponds to a global path in the spectrum, crossing both bands and both gaps but with opposite localization on the edges. In particular, the sum of loops around facing singularities from distinct Riemann sheet (e.g. 0 − and ∞ + ) are homeotopic to + ∪ − and hence related to gap invariant W .

Application to physical models
As mentioned in section 2, the unitary matrixŨ (k) defined in Eq. (3) maps on several physical two-dimensional systems in a cylinder geometry. First, let us notice that from the factorized form (3), it turns out thatŨ actually describes an oriented square lattice similar to the Ho-Chalker model [36] as depicted in Fig. 7. In that case, U 1 and U 2 can be interpreted as scattering matrices that describe the coherent reflection and transmission processes at the nodes of the network. It was shown by Chong and collaborators that such an oriented network actually also describes the propagation of electromagnetic modes in arrays of optical [28,29] or micro-wave [30] resonators beyond the tight-binding model 4 .
Interestingly, this model of a (static) network, also maps on other dynamical Floquet systems as we now show. To see it explicitly, let us factorizeŨ (k) and replace the scattering parameters by their value (Eq. (5)). One gets Clearly,Ũ (k) reveals a quantum protocol that consists in six steps acting on a twolevel system. By repeating periodically this protocol,Ũ (k) can be interpreted as the Floquet operator (evolution operator after one period of time) of a discrete-time quantum walk (see Ref. [37] for a pedagogical introduction) of a quantum system consisting of spin-1/2 particles located at the nodes of a lattice. The first five steps are given byŨ 1 (k) which is block-diagonal in the basis of the position across the (a) Width π/2 θ π/2 θ − π/2 cylinder's width. Thus, these operations are local in position as sketched in Fig. 7 (b) and correspond successively to various shifts and spin-rotations. Finally, the non block-diagonal operatorŨ 2 is applied, so that the corresponding operation (a spin-rotation by an angle 2θ − π around the y axis), is applied on a two-level quantum state which is delocalized on sites n for spin down and n + 1 for spin up. During this last step, one spin at each edge of the strip is left unchanged. In contrast with the oriented network model which is static, this describes a (Floquet) dynamical process.
This model can finally be mapped onto a time-dependent tight-binding model on a square lattice, where the hopping amplitudes are successfully switched-on and off, in the spirit of previous Floquet toy models [18,26,27] as illustrated in Fig. 7 (c). Unlike the two previous models, one can even write down explicitly the corre-sponding Bloch Hamiltonian where g ± (k) = (cos k/2, ± sin k/2, 0). By giving a step profile to the time evolution of the couplings as depicted in Fig. 7 (d), then U 1 (k) and U 2 can be seen as the evolution operators The Hamiltonian (39) has a simple interpretation sketched in Fig. 7 (c) and (d).
First, a coupling of amplitude J 1 is switched-on between second nearest neighbors (A and B sites). Next, a purely imaginary intra-dimer coupling is switched-on with an amplitude J 2 . Then another coupling of amplitude J 3 is switched-on between second nearest neighbors. And finally, a purely imaginary inter-dimer coupling is switched-on with an amplitude J 4 . During this last step, sites at each side of the strip are left uncoupled 5 . Considering that these four steps are repeated periodically in time, then the unitary matrixŨ (k) is nothing but the Floquet operator in a cylinder geometry and the phase can thus be interpreted as the quasi-energy of the periodic dynamics, as for the previous model. These three models, described by the same unitary matrixŨ (k) and thus the same transfer matrix, share the same topological properties independently of whether a Hamiltonian or a periodic dynamics can be associated to them. This illustrates the generality of the framework we have used along the paper.

Discussion
The transfer matrix approach allows the definition of a gap topological invariant, as opposed to a bulk topological invariant defined for the bands of the periodic system, e.g. the Chern number. This is particularly useful for unitary systems for which the Chern numbers can vanish while the system still exhibits topological edge states, as illustrated all along this paper.
The topological nature of the edge states is revealed by the Riemann surface when taking into account the boundary conditions. Note that the construction of the Riemann surface only requires the transfer matrix that describes the (projected) bands in the thermodynamic limit. It thus only contains information about the bulk (one-dimensional) system. This interplay between the bulk information and the boundary conditions is encoded into the quantity Q(k, ) (defined in Eq. (21)) whose change of sign for every (k 0 ) in a gap at fixed k 0 guaranties the existence of an edge state, as shown in figure 3. It follows that a topological transition occurs when the gap closes at (k 0 , 0 ) if and only if Q(k 0 , 0 ) changes sign.
Besides, the topological invariant W appears explicitly as an obstruction in defining a continuous and periodic argument for the phase eigenvalues on the underlying Riemann surface R. This invariant takes a geometrical meaning by counting the number of times a pair of edge states (located at two opposite edges) winds around one of the two loops ± of the punctured torus R. The fact that these two loops are non-homotopic results from the existence of essential singularities that are absent in Hermitian systems for which a similar construction was first performed [31].
Other gap invariants have been proposed in two dimensions, especially in the context of Floquet systems which are described by a time ordered evolution operator [26,27]. In that case, the additional time parameter allows for the definition of a gap topological invariant for the bulk (two dimensional) system, irrespective of the boundary conditions. In contrast, the approach developed in this paper applies for any unitary system in a finite size (cylindrical) geometry, irrespective of the existence of a time dependent dynamics. Furthermore, other topological indexes associated to edge states of Hermitian systems have been already proposed [32,33,34]. Their generalization to unitary systems and the relation with the present index is a natural direction for future investigations. Besides, it would be interesting to bridge the transfer matrix approach with the scattering matrix approach that also provides a topological characterization of edge states in both Hermitian [38,39] and unitary [40] systems. Unlike the transfer matrix which somehow probes a bulk property of the reduced (one-dimensional) system, the scattering matrix requires to connect the system to a lead and thus only probes the edges properties. Moreover, the winding number that is defined is assigned to one edge whereas the two edges are required in the Riemann torus picture.
Finally it would be very interesting to adapt this method to other exotic topological phases, in particular when time-reversal symmetry [27,41] or dissipation [42,43] plays a crucial role.