{'fraction_non_alphanumeric': 0.047552751794648684, 'fraction_numerical': 0.01407439634544268, 'mean_word_length': 4.161801033011453, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Boundaries of boundaries: a systematic approach to lattice models with solvable boundary states of arbitrary codimension

We present a generic and systematic approach for constructing D-dimensional lattice models with exactly solvable d-dimensional boundary states localized to corners, edges, hinges and surfaces. These solvable models represent a class of"sweet spots"in the space of possible tight-binding models---the exact solutions remain valid for any tight-binding parameters as long as they obey simple locality conditions that are manifest in the underlying lattice structure. Consequently, our models capture the physics of both (higher-order) topological and non-topological phases as well as the transitions between them in a particularly illuminating and transparent manner.


I. INTRODUCTION
Topological phases of matter form one of the central topics in condensed matter physics due to their exotic properties, such as the presence of robust boundary states [1][2][3]. The most prominent manifestation of these phases comes in the form of topological insulators, which have a gapped D-dimensional bulk spectrum while supporting a dissipationless current on their (d = D − 1)-dimensional boundaries. Celebrated examples of these modes are the chiral, edge states of quantum Hall and Chern insulators [4][5][6][7][8][9], the helical, edge states of two-dimensional Z 2 insulators [10][11][12][13][14], and the Dirac surface states of threedimensional strong topological insulators [15]. Closely related to these models are topological semimetals whose bulk is instead gapless but nevertheless host robust states on their boundaries. A paradigmatic example of this phase is the Weyl semimetal, which features Weyl nodes in its low-energy description, and supports the existence of Fermi arcs on its surfaces [16][17][18][19][20][21][22].
Recently, an extended family of topological models were proposed in the form of higher-order topological insulators [23][24][25][26][27][28][29][30][31][32][33]. In this case, in-gap modes localize to boundaries with codimension (D − d) > 1, e.g., the zerodimensional corners of two-or three-dimensional latices, or the one-dimensional hinges of three-dimensional lattices. The existence of these (D − d)th-order modes is closely related to the presence of crystalline symmetries that leave the higher-order boundaries invariant but relate the adjacent edges or surfaces to one another [30].
Both in the case of first-and higher-order topological phases the appearance of robust states on ddimensional boundaries is related to the D-dimensional bulk physics via a topological invariant. This celebrated bulk-boundary correspondence is one of the guiding principles for topological phases of matter and underlines the importance of the boundary modes as an indicator for the presence of this phase. While solutions for first-order topological boundary states have been derived in several specific models [34][35][36][37][38][39][40][41][42][43][44] and several general approaches have been developed to retrieve them [45][46][47][48], there is a surprising lack of methods to find analytical solutions for these boundary wave functions that are straightforward and transparent, and can be used to describe modes of any codimension. Such a method is not only of theoretical relevance but also provides practical insight on how to engineer lattices that support these states.
In this paper, building on the brief list of specific examples given in our recent publication Ref. 49, we give a detailed account of a generic method on how to construct such lattice models, which fulfills the aforementioned criteria by allowing for the engineering of the boundary dispersion, giving access to exact solutions for the boundary-mode wave functions and providing the tools to choose the desired localization of these modes. We make explicit use of exact, destructive interference, which is naturally present on a large family of lattices and allows us to write down exact wave-function solutions whose details depend almost solely on a reoccurring structure of the total system, sublattice A, and whose eigenvalues stem from the "on-site" Hamiltonian on A. The localization of these modes is dictated by the indirect connection between different A sublattices via B (s) sublattices [see Figs. 1 and 2], which can be readily adjusted to obtain the desired localization. We thus have the power to choose the codimension, the localization and the dispersion of the exact wave-function solution.
This method has already been investigated and proven to be extremely useful to describe first-order topological phases on frustrated lattice models of the type shown in Fig. 1(a) to find the chiral, edge states of Chern insulators and Fermi arcs of Weyl semimetals [50][51][52]. In the current paper this method is extended to a more general family of lattices and to topological phases of any order. In Ref. 49 we briefly showed how to obtain exact solutions for corner modes of two-dimensional breathing kagome and three-dimensional breathing pyrochlore, as well as, hinge states in two different three-dimensional lattices, one formed of stacked honeycomb lattices and one consisting of stacked Rice-Mele chains [49]. Here, we discuss Breathing pyrochlore Hinge insulators Strong TI [54], Weyl semimetals [50] TABLE I: Examples of D-dimensional models that host ddimensional boundary modes. Solvable models corresponding to each of these can be constructed using the framework presented in this work. Details are given here or in the indicated references.
these results in great detail, and provide a generic and indepth description on how to construct these models. We also introduce an additional example in which we realize chiral, hinge states on the pyrochlore lattice by stacking two-dimensional Z 2 insulators and connecting them in such a way that the helical edge states gap out and only chiral modes on the hinges survive. This technique of creating higher-order topological phases has recently also been put forward in Ref. 31, where a D-dimensional model with second-order modes is created by stacking (D − 1)-dimensional systems with an alternating topological invariant. In contrast, we couple our topologically non-trivial systems to each other via a trivial layer consisting of B sublattices, which, due to the disappearing amplitude of the wave function, does not actively contribute to the boundary wave function. We have exemplified some of the different D-dimensional models with d-dimensional solvable boundary states in Table I. This paper is ordered as follows. In Sec. II, we introduce our method of constructing lattices that have exact wave-function solutions in its full generality. This is followed by a detailed description of the Hamiltonians on each of the individual lattice structures in Figs. 1(a), 1(b), 1(c), and 2, an explicit derivation of the exact wavefunction solutions, as well as concrete and relevant examples to each of the lattices in Secs. III-VI, respectively. We conclude in Sec. VII.

II. GENERAL METHOD
We start by constructing a general family of lattices whose primary building blocks consist of A, B, B , . By allowing no direct hoppings between the A sublattices, these lattices naturally support the exact disappearance of the wave-function amplitude of n wave functions on the B (s) sublattices with n the number of degrees of freedom on the A sublattices. The associated eigenvalues Ai of these wave functions are the eigenvalues of the Hamiltonian h A , which is the "on-site" Hamiltonian on each A sublattice. Wave functions of this form read where k || is the quasi momentum parallel to the direction of the open boundary conditions, i = 1, 2, . . . , n labels the solution, N i ( k || ) is the normalization factor, m (s) labels the A sublattices in the a s+1 direction with a total of creates an electron on sublattice A in unit cell {m (s) } with energy Ai , and r (s) i ( k || ) is obtained by using exact, destructive interference on the B (s) sublattices.
The localization of this wave function on each unit cell is determined by which only has non-zero weight on the A sublattices. We thus find that if |r  Solutions on lattices of the form of Fig. 1(a) are studied in great detail in the context of frustration in Ref. 50, and are generalized in Sec. III. The boundary modes are captured by Eq. (1) with s = 0, and localize to the A sublattice m = 1 (M ) when |r i ( k || )| < 1 (> 1). While the schematic lattice exactly corresponds to the famed Su-Schrieffer-Heeger (SSH) model [53], which we indeed use as an example in Sec. III to illustrate our findings, by dimensional extension it also describes two-and threedimensional systems with each of the sites now corresponding to one-dimensional, periodic chains and twodimensional, periodic planes, respectively [50]. We can thus use our results to find exact solutions for a plethora of different systems, such as for the chiral edge states of Chern insulators and Fermi arcs of Weyl semimetals [50] as well as for the helical edge states of two-dimensional Z 2 insulators and Dirac surface states of three-dimensional strong topological insulators [54] as summarized in Table I.
The schematic lattices in Figs. 1(b) and 1(c) are straightforward extensions of the lattice in Fig. 1(a) by stacking and connecting chains of a similar structure as the chain in Fig. 1(a). In Fig. 1(b), we stack chains consisting of A (red) and B (green) sublattices, and connect them via B (blue) sublattices. On the resulting structure we find that the wave function given by Eq. (1) with s = 0, 1 is an exact solution and corresponds to zero-dimensional corner modes or one-dimensional hinge states in which case each sublattice is associated with a one-dimensional periodic chain. While the lattice in Fig. 1(b) corresponds to the Lieb lattice, in Sec. IV we show an explicit example for the corner modes on the breathing kagome lattice, which only differs from the Lieb lattice by the inclusion of hopping terms between certain B and B sublattices. In that section, we also solve explicitly for the hinge modes on the pyrochlore lattice. The Lieb lattice is treated explicitly in a forthcoming publication [54], where we solve the complete eigensystem by making use of a symmetry that relates the energies E(k x , k y ) = E(−k x , k y ) = E(k x , −k y ). To construct the lattice in Fig. 1(c) we employ a similar construction as before, where we now stack chains consisting of A and B, A and B , and A and B (purple) sublattices in the a s+1 directions such that we can realize zero-dimensional corner states, which are given by Eq. (1) with s = 0, 1, 2. In Sec. V, we provide an explicit example for this schematic lattice in the form of the breathing pyrochlore model, which, similar to the breathing kagome lattice, is readily obtained by allowing for hoppings between certain B, B and B sublattices.
The lattice in Fig. 2 is also closely related to the chain in Fig. 1(a), and consists of A, B, and B sublattices. As we do not allow for any hoppings between the A sublattices, the wave function in Eq. (1) with s = 0, 1 is an exact solution on this lattice. The chains formed by the A and B sublattices (in red and green) and the B sublattices (in blue) are connected in such a way that the lattice has a mirror symmetry around m = (M + 1)/2, and consequently, r i ( k || ) = −1, ∀i. The modes captured by the exact wave-function solution in Eq. (1) thus localize to the A sublattice {m, m } = {1, (M + 1)/2} or {M, (M + 1)/2} depending on whether |r i ( k || )| < 1 or > 1, respectively. In Sec. VI we show how we can use this construction to find exact solutions for hinge modes both in space and time.
In the following, we describe each of these lattices in great detail by studying their concomitant tight-binding Hamiltonians, writing down the relevant form of the exact wave-function solution in Eq. (1), and deriving explicit expressions for r

III. EXACTLY SOLVABLE BOUNDARY STATES WITH CODIMENSION ONE
In this section we show how to find exact boundary-state wave functions on lattices with boundaries of codimen- sion one as shown schematically in Fig. 1(a). Solutions on lattices of this form are extensively described in Ref. 50 in the context of frustrated lattices, and are generalized here to a wider class of lattices. We supplement the general findings with an explicit example in the form of the SSH chain [53]. We start by describing generic, one-dimensional structures, which consist of alternating zero-dimensional A and B lattices, where the A lattices are only coupled to adjacent B lattices through the matrices h A,B;± . Within each A lattice there are n internal degrees of freedom, which may correspond to orbital, spin, or lattice-site degrees of freedom, and h A describes the physics within each A lattice. We denote the corresponding creation operators with c † Ai,m with m the unit-cell index and i = 1, . . . , n. Each B lattice hosts a single degree of freedom with energy h B , and the corresponding creation operator is given by b † m with m the unit cell index, which is the same as the unit-cell index of the A lattice to its left. The Fourier-transformed creation operators are defined as where L denotes the total number of unit cells in the periodic system, and k 1 = k · a 1 with a 1 the lattice vector shown in Fig. 1(a). We find that the Bloch Hamiltonian is given by and For simplicity we have assumed that there is no hopping among the different B sites, which does not restrict the generality of our conclusions as they are not altered by the inclusion of such hopping terms. Next, we consider the corresponding system with open boundary conditions in a 1 that starts and terminates with an A lattice as shown in Fig. 1(a). In this case the Hamiltonian with a total of M A sublattices reads where the row vector is given by . . , c † An , and the matrix H M is of the form which is shown schematically in Fig. 3. Given this structure, it is straightforward to construct exact wavefunction solutions that localize at the boundaries formed by the A sublattices at the ends. For simplicity, we assume that h A is diagonal, i.e. h A = diag( A1 , . . . , An ), which can be achieved by performing a unitary transformation. We continue by taking linear combinations of c † Ai,m |0 and choose their weights in such a way that they interfere destructively on all the B sublattices. Since the on-site energy is the same on all the A lattices, i.e., Ai , we indeed find an eigenstate of the Hamiltonian. In particular, for each i = 1, . . . , n, we make the following ansatz, which corresponds to an exponentially decaying wave function where N i is the normalization constant and |ψ i fulfils the energy eigenvalue equation The explicit equation for r i follows from destructive interference on the B lattices and is given by where the labels {1, i} refer to the matrix elements of the matrices h † A,B;± [68]. The function r i dictates the localization of the wave function, such that while the wave function is completely delocalized in the system when |r i | = 1, we find that it localizes to the right, m = 1, and left end, m = M , when |r i | < 1 and |r i | > 1, respectively. Note that we are able to solve this equation because the number of equations matches the number of variables n, which is why we restrict ourselves to the case where the B lattices host a single orbital, i.e., n B = 1. Moreover, the simple form of this equation for r i follows from writing the Hamiltonian h A in its diagonal basis. If h A is expressed in a different basis, we find that the ansatz in Eq. (4) and the equation for r i in Eq. (7) explicitly depend on the normalized eigenvectors φ i of h A with eigenvalues Ai , such that they read where φ i,j corresponds to the jth component of the ith eigenvector φ i . To elucidate these results, we work out an explicit example in the form of the SSH chain [53], which is also mentioned briefly in Ref. [49]. This model consists of alternating A and B sites, which are coupled through the nearest-neighbor hopping parameters t 1 and t 2 as shown in Fig. 4(a), such that the A lattices host a single degree of freedom, i.e. n = 1. The corresponding Bloch Hamiltonian reads For an open chain with A sites at its ends, we find that the Hamiltonian is given by Eq. (3) with the exact wavefunction solution in Eq. (4). The terms in the Hamiltonian are where we use the notation introduced above such that Eq. (7) now reads where we have dropped the subscript i = 1 of r to simplify notation, and we find r = −t 1 /t 2 . The wave function for the open chain in Eq. (4) with eigenvalue h A = 0 is depicted in Fig. 4(b), where the mechanism that results in the exact disappearance of the wave-function amplitude of the B sites is explicitly illustrated. The band spectrum is shown in Fig. 5(a), and we see that there is an exact zero-energy mode in red for all values of r. This result is in complete agreement with the well-known appearance of an exponentially localized end mode at m = 1 (M ) for |t 1 | < |t 2 | (|t 1 | > |t 2 |), and disappearance of this mode at |t 1 | = |t 2 | (|r| = 1) in the SSH chain with an odd number of sites (oddSSH).  We note that in the case with no broken unit cell, i.e., an even number of total sites (evenSSH), there are two zero-energy (up to algebraic finite-size corrections to the energy) end modes when |t 1 | < |t 2 | and no end modes when |t 1 | ≥ |t 2 | as shown in Fig. 5(b). In the case |t 1 | < |t 2 |, the two ends of the evenSSH chain can locally be mapped onto the end m = 1 of the oddSSH chain, which indeed features an end mode when |t 1 | < |t 2 |. This mapping breaks down, however, when |t 1 | ≥ |t 2 | as the putative solutions interfere strongly with each other due to their increasing amplitude away from their respective ends. In this case, the exact solutions are no longer good approximate eigenstates of the even length system. From the preceding discussion, we observe that while our solutions are only exact on lattices that terminate with A sublattices on both ends, the quantity r i accurately predicts the presence of boundary modes also in the case where there are no broken unit cells, i.e., a termination with an A sublattice on one end and a B sublattice on the other end, in which case |r i | < 1 (|r i | > 1) signals the presence (absence) of boundary modes on both boundaries simultaneously.
We emphasize that the method described in this section is not restricted to one-dimensional models and can be straightforwardly generalized to higher dimensions as long as the boundaries of the model have codimension one. Indeed, by allowing the A and B sublattices in Fig. 1(a) to be (d = 1)-or (d = 2)-dimensional, we form (D = 2)-and (D = 3)-dimensional lattices on which we can realize edge and surface modes, respectively. In these cases one performs a partial Fourier transformation to reduce the model to a family of one-dimensional chains parametrized by k || , which is the momentum parallel to the boundary. This leads to the following, generalized exact solutions of the form of Eq. (1) which differs from Eq. (4) only by the explicit inclusion of momentum dependence. In Ref. 50, this equation was retrieved for the chiral edge states of a Chern insulators, and the Fermi arcs of a Weyl semimetal on frustrated lattices with n = 2 and 3, and n B = 1 (see also Table I).

IV. EXACTLY SOLVABLE BOUNDARY STATES WITH CODIMENSION TWO
We now shift our attention to models exhibiting solvable eigenstates that localize to boundaries with codimension two; these include the so-called second-order topological states. In particular, we find exact wave-function solutions for the zero-dimensional corner states in twodimensional crystals of the type shown in Fig. 1(b), and use this to show in detail how to obtain corner states in the breathing kagome lattice [49]. In addition to the examples mentioned in our previous paper [49], we also include a new example of chiral, hinge states on the pyrochlore lattice.
In this discussion, we begin by considering twodimensional models with zero-dimensional corners for the sake of simplicity but, as before, the method can be easily applied to higher-dimensional systems with boundaries of codimension two of which an explicit example is shown at the end of this section. The idea is that we start with M one-dimensional chains of the form discussed in the previous section, which we then couple to each other in the a 2 -direction by inserting intermediate B sublattices, such that the corners are formed by A sublattices, as is shown in Fig. 1(b). As before, the individual A sublattices host n degrees of freedom and are described by h A , whereas both the B and B sublattices host a single orbital with energies h B and h B , respectively. Moreover, the A sublattices hybridize in the a 1 ( a 2 ) direction only with the B (B ) sublattices, which is described by the matrices h A,B;± (h A,B ;± ). The creation operators on the Fourier-transformed creation operators are given by where L and L denote the total number of unit cells in the periodic system in the directions a 1 and a 2 , respectively, and k i = k · a i . Using this notation we find that the full Bloch Hamiltonian is given by and As in the previous section, we do not consider hoppings amongst the B and B sublattices for reasons of simplicity while the inclusion of such terms does not alter any of our conclusions as is made explicit in the examples that follow below. Next, we consider the same system with open boundary conditions in a 1 and a 2 with A sublattices at its ends with M (M ) A blocks in the a 1 ( a 2 ) direction as shown in Fig. 1(b). For this model, the Hamiltonian reads and the matrix H M,M reads 7 The matrix H M , which is given by Eq. (3) and shown schematically in Fig. 3, accounts for the physics within each row m composed of M A sublattices and M − 1 B sublattices. The matrix H B describes each row m composed of M B sublattices, and can in principle assume any form. However, for simplicity we assume there are no interactions between B sublattices, such that its Hamiltonian is simply given by The matrices H AB,B ;± take care of the hybridization between neighboring rows, and read In a similar spirit to before, we now construct exact wave functions that localize to the corners, where we again assume that h A is written in a diagonal form, i.e., h A = diag( A1 , . . . , An ), such that the corner modes have energies Ai . From the previous section, we know how to construct the exact wave-function solutions for each single row m , which localize to the left or right ends of the row, i.e., m = 1 or M , respectively. These n wave-function solutions are given by which differ from the solutions in Eq. (4) by the introduction of the unit-cell index m and where N i and r i are given by Eqs. (5) and (7), respectively. Next, we go one step further by taking linear combinations of |Ψ i;m , and choosing their weights such that there is also destructive interference on the B sublattices. This leads to the following ansatz where the normalization constant N i is given by Eq. (5) with r i → r i and M → M . The value of r i is found by solving the following equation that corresponds to destructive interference on the B sites such that |Ψ i indeed is a solution to the eigenvalue equation These results correspond to what was found in Ref. 49, and we find that these wave functions remain exact as long as one does not include direct hoppings between the A sublattices, while one is free to include any type of disorder or hoppings within each of the sublattices B and B .
To illustrate these results, we consider the breathing kagome model shown in Fig. 7(a), which was also studied in Refs. 49, 55, 56. The arrangement of the sublattices A, B, and B shown in red, green, and blue, respectively, follows that of the lattice in Fig. 1(b) and the only difference is the coupling between certain B and B sublattices. The A sublattices only couple to the neighboring B (B ) sublattices with alternating, nearest-neighbor hoppings t 1 and t 2 (t 1 and t 3 ), while neighboring B and B sites are coupled via alternating hoppings t 1 and t 2 . Note that these hoppings between B and B can be arbitrarily chosen as they do not influence the form of the exact wavefunction solution. The corresponding Bloch Hamiltonian reads which closely resembles the Hamiltonian in Eq. (13), and h A,B and h A,B given in Eqs. (14) and (15), with Considering open boundary conditions in the a 1 and a 2 direction, we find that the Hamiltonian is given by Eq. (19) with H AB,B ;± reading which differ from the matrices in Eqs. (20) and (21) by the inclusion of extra hopping terms between certain B and B sublattices. The exact wave function is given in Eq. (22) and is unaltered by the extra hopping terms as it only depends on the details of the Hamiltonian on sublattice A and the hoppings between the A and B (B ) sites. The terms in the Hamiltonian are given in Eq. (27), such that Eqs. (7) and (23), which correspond to destructive interference on the B and B sites, respectively, read where we have dropped the subscript 1 on r and r to simplify notation, and which yields r = −t 1 /t 2 and r = −t 1 /t 3 . The exact wave function is shown explicitly in Fig. 7(b) and the energy spectrum for t 2 = t 3 is shown in Fig. 7 Note that these results were also summarized in Ref. 49. Next, we consider the pyrochlore lattice shown in Fig. 8(a). This lattice can be interpreted as a dimensional extension of the kagome lattice in Fig. 7(a) by considering the A, B, and B sites in the latter as onedimensional, periodic chains, which are connected to each other in such a way that they form the pyrochlore lattice. The general structure in Fig. 1(b) thus again applies. We implement two different tight-binding models on this lattice of which the first realizes an all-in-all-out spin configuration on the tetrahedra, which is relevant for the pyrochlore irridates Eu 2 Ir 2 O 7 [57] and Nd 2 Ir 2 O 7 [58]. The Bloch Hamiltonian is given by Eq. (25) supplemented by Eqs. (14), (15), and (26) with where t, t ⊥ , and t ⊥ are nearest-neighbor hopping parameters, and the model is periodic in k || = k · a || with a || the lattice vector shown in Fig. 8(a). This Hamiltonian renders a Chern-insulating phase on each kagome-layer cut in the pyrochlore lattice. Assuming open boundary conditions in the a 1 and a 2 -direction such that we obtain the lattice in Fig. 8(a), we obtain the Hamiltonian in Eqs. (19), (28), and (29) parametrized by k || with the Hamiltonian terms given in Eq. (31). Note that we have not written h A in its diagonal form, such that we cannot simply use Eqs. (7) and (23) to find r i and r i , with i = 1, 2. Instead, we first find the eigenvectors of h A , which are given by φ 1 (k || ) = (1, e ik || /2 ) T / √ 2 and φ 2 (k || ) = (1, −e ik || /2 ) T / √ 2 with eigenvalues t √ 2 cos(k || /2) and −t √ 2 cos(k || /2), respectively. In terms of these eigenvectors, we find that the equations corresponding to destructive interference on the B and B sites read with i = 1, 2, such that where we have substituted the index i for ± where + (−) corresponds to i = 1 (i = 2). The energy spectrum is shown in Fig. 8(b) with the bands highlighted in red corresponding to the exact solution given in Eq. (22) with r i and r i specified above. However, the bands are completely blurred by the bulk bands in blue. In order to isolate the bands corresponding to the exact wave-function solution, we continue by implementing a Kane-Mele-type Hamiltonian on each kagome layer [10,11], such that each of the two sites in the A lattices and the B sites host two degrees of freedom instead of one, while the B sites still host a single orbital per site. The corresponding Bloch Hamiltonian is given by Eq. (25) supplemented by Eqs. (14), (15), and (26) with where t, t ⊥ and t ⊥ are nearest-neighbor hopping parameters, t 2 is a next-nearest-neighbor hopping parameter, and V and V 1 are on-site potentials. Assuming open boundary conditions in a 1 and a 2 as before, we find the Hamiltonian given in Eqs.
where j = 1 (j = 2) if i = 1, 2 (i = 3, 4). From the explicit form of φ i (k || ) and h † A,B ;± , we find that r i (k || ) = 0 (= −∞) if i = 1, 2 (i = 3, 4), such that the wave-function solution in Eq. (22) reduces to and we find that two of the exact solutions are precisely localized to the layer m = 1 while the other two solutions are localized to the layer m = M . The band spectrum for this model is shown in Fig. 8(c), where the perpendicular hopping between the layers m is taken to be t

V. EXACTLY SOLVABLE BOUNDARY STATES WITH CODIMENSION THREE
Finally, we repeat the same steps as before to find exact solutions for the zero-dimensional corner modes of the three-dimensional model shown in Fig. 1(c), which is a dimensional extension of the general lattice in Fig. 1(b). We demonstrate the validity of our findings by applying it to retrieve zero-energy corner modes in the breathing pyrochlore lattice. In a similar spirit to before, we start with M copies of the two-dimensional crystal discussed in the previous section for which we can find exact boundary states. We couple these models by inserting intermediate B sublattices, which host a single orbital only at an energy h B as shown in Fig. 1(c). Again, it is crucial that the A sublattices only talk to neighboring B, B , and B sublattices.  a 1 , a 2 and a 3 , respectively, with a i the lattice vectors shown in Fig. 9. The corresponding Fourier-transformed creation operators are given by where L, L , and L denote the total number of unit cells in the direction a 1 , a 2 , and a 3 , respectively, and k i = k · a i . This leads to the Bloch Hamiltonian and with h A,B and h A,B given in Eqs. (14) and (15), respectively.
Next, we provide the Hamiltonian that describes the same system with open boundary conditions in all three directions as shown in Fig. 1(c).
with d † m ,m and e † m ,m defined in Eqs. (16) and (18) and schematically shown in Fig. 9. The matrices H M,M and H B on the diagonal correspond to the physics within each two-dimensional lattice labelled by m that are stacked in the a 3 direction, and are given by Eq. (19) and with Again, we have not included hoppings that hybridize B, B , and B sublattices for the sake of simplicity. We continue by finding exact corner-state solutions for this structure. Note that in the previous section we found exact wave-function solutions within each single plane m , which typically localize to one of the four corners.
For the three-dimensional system, these are given by with i = 1, . . . , n, and differ only from Eq. (22) by an extra unit-cell index m . As before, we take linear combinations of these wave functions such that they also interfere destructively on the B sublattices. This leads to the following wave-function solutions where N i is the normalization constant given by Eq. (5) with r i → r i and M → M , and r i solves the following equation ensuring As an illustrative example we now consider the breathing pyrochlore lattice shown in Fig. 10(a), which was also studied in Refs. 49, 55. This lattice can be thought of as consisting of a stack of breathing kagome lattices that are coupled via intermediate triangular lattices forming tetrahedral structures. The sites in the intermediate triangular lattice play the role of the B sites. The relevant hopping parameters are shown in Fig. 10(a), and we again see that certain B, B , and B sites are connected with each other. The Bloch Hamiltonian for this model reads which closely resembles Eq. (35) and h A,B , h A,B , h A,B , h B,B given by Eqs. (14), (15), (36) and (26), respectively, with which differ from Eqs. (38) and (39) by the inclusion of extra hopping terms between B and B (B and B ) sublattices.
The exact wave function is given in Eq. (40) and is unaltered by the extra hopping terms as it only depends on the details of the Hamiltonian on sublattice A and the hoppings between the A and B (B ) ((B )) sites. The terms in the Hamiltonian are given in Eq. (46), such that Eqs. (7), (23), and (41), which correspond to destructive interference on the B, B , and B sites, respectively, read t 1 + r t 2 = 0, t 1 + r t 3 = 0, t 1 + r t 4 = 0, (47) where the subscript 1 on r, r and r is dropped to simplify the notation, such that r = −t 1 /t 2 , r = −t 1 /t 3 , and r = −t 1 /t 4 . The exact wave function is shown explicitly in Fig. 10

VI. EXACTLY SOLVABLE HINGE STATES ON MIRROR SYMMETRIC LATTICES
Lastly, we focus on finding wave functions that localize at mirror-symmetric corners and hinges with codimension two of the type shown in Fig. 2. The structures we consider are composed of two different types of onedimensional chains, which are shown in red and green, and in blue corresponding to a chain made up of zerodimensional A and B lattices and zero-dimensional B lattices, respectively. Compared to the structure shown in Fig. 1(b) the main difference is that now the B chain consists of two sublattices (in Fig. 2 this is shown with two blue sites instead of only one). Earlier we remarked that it is crucial for the retrieval of exact wave-function solutions that both the B and B sublattices host a single degree of freedom. There are, however, some instances where this condition can be relaxed; if we couple the AB and B chains in a mirror-symmetric fashion we again retrieve wave-function solutions of the form of Eq. (22). Moreover, we again remark that this construction can be easily expanded to higher dimensions as is demonstrated in the explicit examples that follow.
We assume that the A sublattices host n degrees of freedom and B one. Moreover, we refer to the two inequivalent sublattices of the B chains as B I and B II sites, which are located between A and B sublattices, respectively, and carry n I and n II degrees of freedom, respectively. We require that the A lattices only hybridize with the B and B I sites and not with the B II sites. The corresponding creation operators are denoted by c † Ai,m,m , b † m,m , b † I,i 1 ,m,m , and b † II,i II ,m,m , respectively, with i = 1, . . . , n, i I = 1, . . . , n I , and i II = 1, . . . , n II . We assign to each B site the same unit-cell index as the A site to its left, and to each B I (B II ) site the same index as the A (B) site below. The corresponding Fouriertransformed creation operators are defined as before [see Eqs. (10)- (12)] and we find that the Bloch Hamiltonian is given by In contrast to the notation used so far, we have not included subindices ± in the n × n I matrix h A,B I and the 1 × n II matrix h B,B II , since the coupling in the a 2direction between the A and B I (B and B II ) sites is symmetric. For simplicity, we have assumed that there is no hopping among the B and B I lattice sites, but as before including these terms does not alter our findings. Next, we consider the system with open boundary conditions with A lattices at its corners. Here, we provide the generic Hamiltonian associated with the lat-tice shown in Fig. 2 α and Ψ † is a row-vector of creation operators, which reads The full Hamiltonian is shown schematically in Fig. 11.
To find the exact corner-state wave functions, we need to solve the equations that correspond to destructive interference on the B and B I sites. Since, the B I site may in principle host multiple orbitals, we have to solve 1+n I equations. Note that as we do not hybridize A lattices with the B II sites, we get destructive interference for free for the latter. Starting with the ansatz given in Eq. (22), we obtain the following equations ) l,i = 0, l = 1, . . . , n I . (57) Note that the second set of equations is trivially solved by r i = −1, which is due to the mirror symmetry. Therefore, we find that the mode described by the wave function in Eq. (22)  In Ref. 49, we showed that we can realize models on this lattice structure that are topologically protected by virtue of the mirror Chern number C m (k 2 ) = (C + (k 2 ) − C − (k 2 ))/2, where C +(−) (k 2 ) corresponds to the Chern number in the even (odd) sector of the mirror eigenvalues. If C m (0) + C m (π) = 1 mod 2, where k 2 = 0, π are the mirror-invariant slices in the Brillouin zone, there are an odd number of chiral, boundary modes at the boundaries that are symmetric under a 2 → − a 2 . If we thus require that on each AB sublattice we implement a model that has a non-zero Chern number opposite to that of the model on B sublattices, we find that the boundary modes of neighboring sublattices gap out and only those states on the layer m = (M + 1)/2 remain. In this case, while the "ordinary" Chern numbers C(0) and C(π) are zero, the mirror Chern number yields the desired situation, i.e., C m (0) + C m (π) = 1 mod 2, implying the presence of boundary modes exponentially localized to {m, m } = {1, (M + 1)/2} and {M, (M + 1)/2}. We supplement the above findings with two explicit examples, which were only briefly discussed in Ref. 49. As a first example we consider a stack of honeycomb lattices, where each of the A, B, B I and B II sublattices represents a one-dimensional periodic chain with one degree of freedom such that each of the individual m AB and B chains in Fig. 2 forms a honeycomb lattice with a zigzag edge on one side and a bearded edge on the other as shown in Fig. 12(a). On each honeycomb lattice, we realize a Haldane-like tight-binding model [5] as shown in Fig. 12(a) for which we find an exact hingestate solution. The corresponding Fourier-transformed Hamiltonian is given by Eqs. (48)- (52) with parameters as shown Fig. 12(a). The term ±2∆ sin(k || ) originates from purely imaginary next-nearest-neighbor hoppings ±i∆ in the direction a || , and realizes a Cherninsulating phase on each individual honeycomb lattice when ∆ = 0 while an ordinary semi-metallic phase is observed for ∆ = 0. The opposite sign of ∆ in the AB and B planes ensures that the two honeycomb lattices have opposite Chern number, such that the two different planes possess counter-propagating edge states. Even though this is a three-dimensional model, we can view it as a collection of two-dimensional models parametrized by k || . In particular, we may consider periodic boundary conditions in the a || direction, and open in the directions a 1 and a 2 with the termination shown in Fig. 11. Using the notation introduced above, we have n = n I = n II = 1, and the Hamiltonian is given in Eqs. (53)-(55) with the hopping terms in Eq. (58). Eqs. (56) and (57) thus where we dropped the index 1 in r and r for simplicity. This leads to r(k || ) = − 1 + e −ik || and r = −1. The band spectrum is shown in Fig. 12(b), and the chiral hinge state described by Eq. (22) has eigenenergy h A = −2∆ sin(k || ) shown in red.
As a second example, we consider a stack of Rice-Mele chains [59]. Each A, B, B I , and B II sublattice is formed by a single site with one degree of freedom such that we literally consider the lattice as given in Fig. 2. We implement the Rice-Mele model on each individual AB and B chain as shown in Fig. 12(c). The Bloch Hamiltonian is given by Eqs. (48)- (52) with where s = −1, t and s are nearest-neighbor hopping parameters,δ(t) = δ(t)/t = cos(2πt/T )/t and T the period of driving as shown in Fig. 12(c). Considering open boundary conditions, we find that the Hamiltonian is given by Eqs. (53)-(55) with the hopping terms in Eq. (58) and with n = n I = n II = 1. The band spectrum is shown in Fig. 12(d). This leads to the following for Eqs. (56) and (57) (t + δ(t)) + r(t) (t − δ(t)) = 0, s + r s = 0, where we dropped the index 1 for simplicity. Hence, at each instantt we find precisely one corner mode with energy s sin(2πt/T ), or a hinge mode in time, and its associated wave-function solution given in Eq. (22) with r(t) = −(t + δ(t))/(t − δ(t)) and r = −1.

VII. DISCUSSION
In this work, we have presented a systematic and coherent method on how to construct lattices on which generic local tight-binding models feature exactly solvable states exponentially localized to boundary motifs of arbitrary codimension and with a dispersion that can be tuned at will by tuning the local hopping amplitudes. By adjusting the hopping parameters that determine the connection of the A sublattices to the B (s) sublattice(s), we have full access to all single particle properties, such as the (D − d) localization lengths, of these modes. The validity and exact form of our solutions are completely independent of the topological properties of the model as well as the size of the system, and we have demonstrated that while they are no longer exact solutions, they also provide useful insights in the phase diagram of models whose boundaries preserve unit cells, i.e., boundaries made up of A sublattices on one end and B (s) sublattices on the other. While we have restricted ourselves to the dimensions relevant to our three-dimensional physical world in the investigation of different concrete examples, this method can be straightforwardly extend to describe systems of any dimension D with boundaries of any dimension d < D. The results presented in this paper were summarized in a shorter previous work in Ref. 49, and have been extended and explored here in more detail. We point out that our bottom-up approach to engineer lattices with the wished for properties is extremely versatile and transparent, and provides useful insights both in the experimental and theoretical context. Especially, in the latter context our method shines a rare but clear light on the mechanisms underlying the appearance of boundary modes of any codimension. While several other methods have been introduced to find exact wavefunction solutions in systems with open boundaries with codimension one [45-48, 65, 66], none of these methods have been generalized to describe higher-order phases so far, and it is not a priori obvious that it is at all possible. Indeed, while all of these methods have different strengths, the clear advantage of our approach is that its generalization to boundaries of any codimension is simple and straightforward. This is also advantageous in the experimental setting in the sense that our method allows for the proposal of simple lattice structures that realize interesting phenomena. Indeed, the type of lattices that we propose can be straightforwardly realized in artificial lattices such as cold atom systems, optical experiments, microwave cavities and acoustic metamaterials [9,[60][61][62][63][64].
While this work does not explicitly focus on establishing the topological properties of the boundary modes we realize beyond looking at whether the state corresponds to an in-gap state and changes localization upon a bandgap closing, we point out that the topological robustness of several of the modes we have realized can be straightforwardly established. Indeed, for models with boundaries of codimension one, topological invariants have been studied in numerous examples, e.g., the end modes of the SSH chain are protected by chiral symmetry, which allows for the definition of a winding number, whereas the chiral modes of a Chern insulator are protected by a non-zero Chern number. We make use of the latter when realizing chiral hinge states in the pyrochlore lattice as well as constructing mirror-symmetric lattices that feature chiral hinge modes, which are protected by the mirror Chern number. To establish the topological robustness of higher-order models in general, one may make explicit use of the crystalline symmetries that protect these modes to define invariants [32], whereas a more generic approach is to use the Wilson-loop formalism to show topological protection [28]. A quick analysis, however, results in the conclusion that the lack of a spectral symmetry, e.g., E k = −E k , in the spectrum of the breathing kagome and pyrochlore models prevents the odd number of corner modes we find at half filling from being topologically protected. Indeed, when introducing a defect on the A sublattice, the zero-energy mode may be pushed away into the bulk, while the localization and exact form of the wave function is unchanged when the defect is introduced on any of the B (s) sublattices. This simple consequence of our construction clearly demonstrates the stability of boundary modes even when they are not explicitly topologically protected; while altering the bulk spectrum, a defect on one or more B (s) sublattice site does not alter the properties of the boundary mode. This finding is especially relevant in the case of large lattices, i.e., M (s) → ∞, where the explicit role of the A and B (s) sublattices is washed out away from the boundaries. As a last note, we point out that upon turning of certain couplings such that the breathing kagome and pyrochlore lattices take the form of the Lieb lat-tice [cf. Fig. 1(b)] and its three-dimensional cousin [cf. Fig. 1(c)], respectively, the spectral symmetry is restored by virtue of which the corner modes in these lattices captured by Eq. (1) are protected.
The advantage of the transparency of our method has in fact already proven to bear fruit. Notably, it was very recently demonstrated that topological phases in the context of non-Hermitian systems, which set themselves apart from topological phases in Hermitian models by the general breaking down of conventional bulk-boundary correspondence, can be directly described in a system with open boundary conditions using biorthogonal quantum mechanics [67]. This important insight was strongly inspired by the use of our exact solutions which readily generalize to the non-Hermitian setting. This work thus serves as a supplement to our previous works and can be viewed as a comprehensive guidebook on how to engineer lattices from the bottom up with the desired topological properties up to any order.