Fermionic tensor networks for higher order topological insulators from charge pumping

We apply the charge pumping argument to fermionic tensor network representations of d-dimensional topological insulators (TIs) to obtain tensor network states for (d+1)-dimensional TIs. We exemplify the method by constructing a two-dimensional projected entangled pair states (PEPS)for a Chern insulator starting from a matrix product state (MPS) in d=1 describing pumping in the Su-Schrieffer-Heeger (SSH) model. In extending the argument to second-order TIs, we build a three-dimensional TNS for a chiral hinge TI from a PEPS in d=2 for the obstructed atomic insulator (OAI) of the quadrupole model. The (d+1)-dimensional TNSs obtained in this way have a constant bond dimension inherited from the d-dimensional TNSs in all but one spatial direction, making them candidates for numerical applications. From the d-dimensional models, we identify gapped next-nearest neighbour Hamiltonians interpolating between the trivial and OAI phases of the fully dimerized SSH and quadrupole models, whose ground states are given by an MPS and a PEPS with a constant bond dimension equal to 2, respectively.


I. INTRODUCTION
Higher order TIs [1][2][3][4] have recently been introduced as a new class of symmetry protected topological systems generalizing the framework of TIs with surface states [5]. According to one definition, a TI of order n in d spatial dimensions has topological boundary modes at the (d − n)dimensional intersection of n crystal faces. In this terminology, conventional TIs such as Chern insulators [6] are of order n = 1 with protected boundary modes at (d − 1)-dimensional edges. On the other hand, second order TIs protected by mirror or rotation symmetries have zero-dimensional corner states in d = 2 dimensions (in which case they may exemplify obstructed atomic limits [7] with local Wannier states) and one-dimensional chiral or helical hinge states in d = 3 dimensions (and bulk bands without a local Wannier description). Prototypical examples of second order topological phases include the two-dimensional quadrupole model of Ref. [1], a natural extension [8] of the Su-Schrieffer-Heeger (SSH) model [9], and the three-dimensional chiral hinge insulator of Ref. [2], both of which have been experimentally observed in either materials [10] or mechanical [11], acoustic [12,13], photonic [14][15][16][17] and electrical [18][19][20] systems.
The bulk-boundary correspondence states that the topological properties of a system are reflected in the excitations at its physical boundary. For instance, a twodimensional Chern insulator is characterised by an integer number of gapless chiral edge modes [21]. Similarly, second order TIs in two dimensions possess gapless corner modes at the intersection of two edges compatible with the crystal symmetry [1]. In three dimensions, one type of second order TI is characterised by the presence or absence of a chiral hinge mode [2]. For strong TIs [22] as well as chiral topological phases, the universal features in the boundary energy spectrum are encoded in its entanglement spectrum (ES) [23] characterizing the bulk entanglement properties. The convenience of this bulk characterization makes the ES an important tool for the numerical analysis of topological phases. Recently it was also observed, as expected, that the ES of higher order TIs displays characteristic (d − n)-dimensional boundary states as long as the entanglement cut preserves the crystal symmetries of the phase [24][25][26].
A natural platform for the bulk-boundary correspondence is provided by tensor network states (TNSs) [27] in which the entanglement between physical particles is mediated through virtual particles hosting the lowerdimensional boundary theory [28]. Physical and virtual particles are related by a local tensor whose bond dimension determines the maximal amount of entanglement in the state [29,30]. The structure of the local tensor encodes the topological properties of the quantum state, and an analysis of this relation has led to valuable insight for the systematic understanding of one-dimensional symmetry protected topological (SPT) phases [31] and intrinsically ordered topological phases [32]. Free fermion systems can be described using Gaussian fermionic TNSs (GfTNSs) that are defined purely in terms of two-point correlation functions [33]. Starting from such Gaussian TNSs, interacting states can be constructed, whose topological properties derive from the features of the initial Gaussian model [34].
Many non-chiral topologically ordered phases in two dimensions [35,36] have a simple representation in terms of projected entangled pair states (PEPSs) [37] with a finite and constant bond dimension [38,39]. However, Gaussian PEPSs with chiral topological properties do not adequately describe the bulk of a gapped phase, since their correlation functions decay algebraically [40,41]. Indeed, the existence of gapped PEPSs with finite bond dimension for SPT phases of free fermions is forbidden in dimensions d > 1 [41,42]. The no-go theorem states that only Wannierizable phases corresponding to the product of one-and zero-dimensional systems can be exactly represented as a TNS.
This no-go theorem does not prevent TNSs from being a useful numerical description at finite system size for chiral topological phases like the fractional quantum Hall effect [43]. Hence, numerically efficient TNSs for such phases are valuable, especially in two and three dimensions. One method for the construction of TIs in d + 1 dimensions from d-dimensional TIs is given by charge pumping [44]. For instance, a two-dimensional Chern insulator is obtained from a charge pumping interpolation in the one-dimensional SSH model when the periodic time direction is identified with the momentum in the second spatial direction [45]. For second order TIs, dipole pumping in the quadrupole model defines a three-dimensional chiral hinge insulator [2,3]. In both cases, the zero-dimensional boundary modes of the ddimensional model give rise to one-dimensional chiral boundary modes in the (d + 1)-dimensional system.
In this article, we exemplify how the charge pumping argument applied to d-dimensional GfTNSs with constant bond dimension yields gapped GfTNSs for TIs in d+1 dimensions. By construction, the (d+1)-dimensional TNS has a constant finite bond dimension in a hybrid coordinate system with d real-space axes and one momentum axis in the additional dimension. In order to obtain a real-space tensor network for the (d + 1)-dimensional state, we apply to the hybrid TNS an inverse Fourier transform (FT) in the direction d + 1. As a result, the bond dimension of the real-space local tensor in this direction generically grows with the system size due to the non-locality of the FT, whereas it is identical to the finite bond dimension of the d-dimensional TNS in the other d directions. We apply this construction both to a matrix product state (MPS) [46] for the SSH model in order to obtain a PEPS for a Chern insulator, and to a Gaussian PEPS with finite bond dimension for the topological quadrupole model in order to obtain a threedimensional GfTNS for the second order chiral hinge state TI of Refs. [2,3]. Therefore, our approach provides us with a gapped TNS with one-dimensional chiral boundary states and a constant finite bond dimension in all but one of the spatial directions. This representation is therefore potentially useful for tensor network algorithms.
This article is structured as follows: In Sec. II we begin with a brief overview of charge pumping in the SSH model and fermionic TNSs, and continue by studying an MPS for the SSH ground state along the charge pump-ing interpolation. In the following Sec. III, we employ this MPS to construct a two-dimensional PEPS describing a Chern insulator. In Sec. IV we extend our method to the two-dimensional second order quadrupole insulator and the construction of a three-dimensional PEPS for the chiral hinge state higher order TI. Finally, we summarize our results and discuss remaining open questions in Sec. V.

II. FERMIONIC MPS FOR CHARGE PUMPING IN THE SSH MODEL
In this section we introduce a fermionic MPS that describes a charge pumping cycle in the SSH model corresponding to a Chern insulator with Chern number C = 1 in two spatial dimensions. We begin by briefly reviewing the SSH model, the charge pumping argument and its relation to Chern insulators in Sec. II A. We continue with a short introduction to fermionic MPSs in Sec. II B. In Sec. II C, we construct the MPS for the SSH model and study its parent Hamiltonian.
A. Chern insulator from charge pumping in the SSH model The SSH model describes a one-dimensional chain of spinless fermions with two orbitals denoted A and B per unit cell [9]. We consider the model at half filling where the number of particles is equal to the number N x of unit cells. For open boundary conditions, the Hamiltonian reads where we use the notationâ A,x andâ B,x for the fermionic annihilation operators of the orbitals A and B in unit cell x with x = 0, . . . , N x − 1, respectively. Here, t (0) denotes the hopping amplitude between A and B orbitals within the same unit cell, and t (1) the hopping amplitude between sites on neighboring unit cells (see Fig. 1(a)). For |t (1) | > |t (0) |, the hopping between unit cells dominates the hopping within unit cells and the SSH model is in a phase topologically different from the case |t (0) | > |t (1) |. This phase is protected by spatial inversion symmetry and characterised by fermionic edge modes [45] and a "filling anomaly" [8]. It is called obstructed atomic insulator (OAI) [7]. When the intra-cell hopping t (0) vanishes in the deep OAI phase, the system is dimerized since its bulk splits into decoupled two-site blocks. In this case, the edge excitations are created by the mode operatorŝ a † A,0 andâ † B,Nx−1 and have exactly zero energy in this specific model, but can generally be moved in energy.
On the other hand, for |t (0) | > |t (1) | the SSH model is trivial, with a dimerized point at t (1) = 0 and no state at zero energy.
(3) At t = −π (equivalent to t = π), the system is in an atomic phase with only the A orbitals occupied. For −π < t < 0, the coupling between A and B sites in the same unit cell is non-zero and the charge is transferred from left to right by the changing chemical potential until only the B orbitals are occupied at t = 0. At t = −π/2, the staggered chemical potential vanishes and the system is in the trivial dimerized phase of the SSH model. For 0 < t < π, the intra-unit cell coupling vanishes whereas the hopping between different unit cells is non-zero and the charge is transferred from left to right such that at t = π the system returns to the state with all A orbitals occupied. At t = π/2, the chemical potential vanishes and the system is in the OAI dimerized phase of the SSH model. The ground state of the pumping Hamiltonian H pump (t) is continuous as a function of time t if the chain has periodic and anti-periodic boundary conditions for N x odd and even, respectively.
The charge pumping interpolation of Eq. (3) can be used to generate a lattice model in one dimension higher with the topology of a Chern insulator [3]. Indeed, H pump (t) corresponds to the time-dependent Bloch Hamiltonian where σ 1 = ( 0 1 1 0 ), σ 2 = ( 0 −i i 0 ) and σ 3 = ( 1 0 0 −1 ) are the Pauli matrices and k x = π Nx (2j − N x ) ∈ [−π, π] for 0 ≤ j ≤ N x −1 the lattice momentum. Since the interpolation is cyclic, the time t ∈ [−π, π] may be interpreted as the lattice momentum k y of a second spatial direction y and Eq. (4) as the Bloch Hamiltonian of a two-dimensional system with closed boundaries in both directions. This Hamiltonian has a Chern number C = 1 due to the charge transport between unit cells induced by the interpolation for t ∈ [0, π] [3].
We note that the cycle of Eq. (3) can be deformed to a smooth charge pumping interpolation, without changing the topology or breaking the dimerization. In this cycle, the couplings µ(t), t (1) (t) and t (0) (t) and hence the Bloch Hamiltonian are smooth functions of the time t, i.e. they are infinitely often continuously differentiable with respect to t. As explained in the previous paragraph, a two-dimensional Chern insulator model is obtained by identifying the time t and the momentum k y . Due to the smoothness of the cycle, the real-space representation of the Chern insulator has couplings which decay faster than any polynomial [47]. The smooth cycle is obtained by smoothing out the non-analyticities of Eq. (3) at t = 0, π using smooth functions that interpolate between 0 and 1. This can be done in such a way that, like in Eq. (3), at each time t the smooth couplings satisfy [48].

B. Fermionic MPSs for a half-filled lattice
Similarly to bosonic MPSs for spin chains, fermionic MPSs describing chains of electrons are obtained by associating virtual particles to each physical particle which mediate the entanglement between different physical constituents [33]. In the case of fermionic tensor networks, the virtual particles obey fermionic statistics. A fermionic MPS with f physical complex fermionic modes per site and ξ virtual complex fermionic modes per lattice site and nearest-neighbour bond has physical dimension d p = 2 f and bond dimension D = 2 ξ . The state is fully characterized by the set of local maps which relate the virtual and physical particles associated to one lattice site. With respect to orthonormal bases {|i } with i = 0, . . . , d p − 1 for the physical Hilbert space and {|l }, {|r } with l, r = 0, . . . , D − 1 for the left and right virtual spaces on one site, respectively, the local maps are represented by d p local matrices A i lr of dimension D × D [49]. In order to fix the global parity of the state, we restrict ourselves to parity-even local tensors which preserve the fermion parity between the physical and virtual layers [49]. We choose the orthonormal basis for the physical Hilbert space such that all basis states |i have either even fermion parity |i| = 0 or odd fermion parity |i| = 1, and similarly for the virtual Hilbert spaces. A local tensor is parity-even if for all non-vanishing elements A i lr = 0, the total fermion parity |i| + |l| + |r| = 0 is even, such that for all configurations of i, l and r (where no summation is implied). If the MPS is constructed from tensors A[j] at each site j with this property, the state on a closed chain with N x sites is given by Here, |l| is the parity of the virtual basis state on the link between the first and last site, and = 1 or = 0 corresponding to periodic and anti-periodic boundary conditions for the many-body state, respectively [49]. Below, we construct a fermionic MPS with parity-even local tensors for the ground state of the SSH model along the dimerized charge pumping interpolation. Due to the parity symmetry of Eq. (5), this MPS necessarily has an even number of physical particles on a closed chain where all virtual bonds are contracted. On the other hand, the ground state of the SSH model is half-filled such that the number of particles is odd if the number of unit cells is odd. Therefore, to construct the parityeven SSH MPS we use a many-body basis built by acting with creation operators on the state |Ω that contains N x physical particles. To that end, we define new mode operators a A,x and a B,x for the physical fermions by performing a particle-hole transformation on all B orbitals of the SSH chain while leaving the A orbitals unaltered, The state |Ω is the vacuum of the new operators, given in terms of the original vacuum state |Ω withâ A,x |Ω = a B,x |Ω = 0 as It satisfies a A,x |Ω = a B,x |Ω = 0. Therefore, the new vacuum contains N x of the original fermions and thus is half-filled. The MPS of the SSH model along the interpolation is then defined with respect to the Fock states constructed from the modes a † A,x and a † B,x acting on the new vacuum |Ω . This particle-hole transformation is motived by our desire to write the MPS in terms of parity-even local tensors, which will in turn enable us to express the state as a GfTNS and thereby compute a Bloch parent Hamiltonian analytically.

C. MPS for the SSH model
In this Section, we study a fermionic MPS for the SSH model. In Sec. II C 1, we define the state in terms of its local tensors and identify the U(1) symmetry leading to a conserved number of particles. In Sec. II C 2, we derive a parent Hamiltonian for the MPS, allowing us to conclude in Sec. II C 3 that the MPS describes the SSH model along the charge pumping cycle of Eq. (3).

Local tensors and U(1) symmetry
In order to describe the ground state of the SSH model along the dimerized charge pumping interpolation, we consider a fermionic MPS with physical dimension 1. (a) Extended SSH model with a staggered chemical potential µ, nearest-neighbour hoppings t (0) and t (1) within and between unit cells, and next-nearest neighbour hoppings t corresponding to one fermionic mode per site and bond dimension D SSH = 2. This is the minimal bond dimension for an MPS along the charge pumping cycle, since the ES of an open SSH chain in the OAI and trivial dimerized phases has two degenerate levels w.r.t. all cuts between and within unit cells, respectively.
The fermionic MPS is translation invariant with a unit cell of two sites. It is therefore fully specified by the local matrices A i lr and B i lr with i, l, r ∈ {0, 1} for sites on the A and B sublattices, respectively. In terms of the local tensors the physical state on a closed chain is given as where = 1 or = 0 for periodic and anti-periodic boundary conditions. The physical many-body basis state is with the vacuum |Ω from Eq. (8). Guided by the dimerized limits to be discussed below, for the local MPS matrices we make the ansatz depending on parameters α, β, γ ∈ R. Note that the normalised quantum state defined by these local MPS matrices depends only on the two quotients α/γ and β/γ. Nonetheless, we choose to work with the parametrisation of Eq. (11) since the case γ = 0 can be treated more conveniently without divergences. In order to motivate the ansatz of Eq. (11) for the local matrices, let us see how the MPS from Eq. (9) can describe the ground state of the SSH model both in the trivial and the OAI dimerized phases with appropriate choices for the parameters α, β and γ. In the trivial dimerized phase of the SSH model, each unit cell decouples from the rest of the system and is in the state (a † A a † B −1)|Ω . In order to ensure the absence of entanglement on the bonds between unit cells, the blocked MPS matrices (A i A B i B ) lr for one unit cell should be non-zero only if l = r = 0. Then, the restriction of the MPS to the unit cell, Therefore, the MPS represents the trivial dimerized phase of the SSH model if the parameters are chosen as α = 0 and |β| = |γ| = 0. Similarly, in the OAI dimerized phase the chain splits into decoupled plaquettes composed of two adjacent sites from different unit cells. The state on each of the plaquettes is given by the superposition (a † B a † A + 1)|Ω such that the blocked MPS matrices for the plaquette should have only two non-zero entries (B 0 A 0 ) 00 = (B 1 A 1 ) 00 . This is achieved if β = 0 and |α| = |γ| = 0. The ansatz of Eq. (11) is chosen such as to allow an interpolation between these two cases.
The MPS of Eq. (9) with the parametrisation of Eq. (11) has a U(1) symmetry which ensures that the physical state lies at half filling of the chain. The U(1) rotation acts on a single complex fermion with the matrix U (ϕ) = 1 0 0 e iϕ . The local tensors for the two sublattices are invariant under a combination of U(1) rotations of the physical and virtual legs as sketched in Fig. 1(b) and (c) (see for instance Ref. [50] for MPSs with physical symmetries), If we choose ϕ = π, these identities correspond to the parity symmetry of Eq. (5), showing that the local MPS tensors are parity-even. For each nearest-neighbour bond, one virtual leg transforms with U (ϕ) and the other with U (ϕ) † in Eq. (12), such that the two U(1) rotations cancel if the bond is contracted. Therefore, the state on a chain with closed boundaries after the contraction of all virtual bonds is invariant under the physical part of the U(1) rotations of Eq. (12), given by staggered transformations U (ϕ) and U (ϕ) † on A and B orbitals, respectively. In-variance under this global symmetry implies that forcing the number of particles measured in terms of the original physical modes to be equal to the number of unit cells as required for the SSH ground state. In Eq. (12), the physical legs on the A and B sublattices transform as particles and holes, respectively, as expected due to the particle-hole transformation of Eq. (7). In addition, Eq. (12) implies that the virtual legs on the A and B sublattices transform as hole-like and particle-like degrees of freedom (DOFs), respectively.
In order to gain a better understanding of the parameters α, β and γ as well as of the systems described by the MPS from Eq. (9), it is helpful to consider a parent Hamiltonian for which it is the exact ground state. This Hamiltonian can be computed directly in terms of its Bloch representation once the MPS of Eq. (9) is expressed as a Gaussian fermionic TNS.

Parent Hamiltonian
Since the charge pumping Hamiltonian of Eq. (3) describes non-interacting fermions, its ground state is a fermionic Gaussian state. It is thus fully characterised by its covariance matrix (CM) (see Appendix A for a summary of our conventions) [51]. As we review in Appendix B, the tensor network formalism may be used to construct Gaussian fermionic TNSs which are ground states of free fermion Hamiltonians [33]. The CM in Fourier space of a translationally invariant Gaussian TNS is given by a simple expression which can often be evaluated analytically. Then, any positive function (k) > 0 on the Brillouin zone gives rise to a parent Hamiltonian with dispersion relation (k), whose Bloch representation is obtained by multiplying the CM by (k).
For all values of the parameters α, β and γ, the MPS of Eq. (11) corresponds to a Gaussian fermionic TNS whose Fourier CM is computed analytically in Appendix E. We show that this MPS with bond dimension D = 2 is the unique ground state of a longer-range SSH-like model with a staggered chemical potential µ and next-nearest neighbour hoppings t (2) A and t (2) B on the A and B sublattices, respectively, as sketched in Fig. 1(a). The coupling constants of the parent Hamiltonian H MPS depend on the parameters of the MPS as where t (0) and t (1) denote the hopping amplitudes between nearest-neighbour sites on the same and adjacent unit cells, respectively. Depending on the parameter values, H MPS describes different phases of the fermionic chain. If any two out of the three parameters α, β, γ vanish, the system is in an atomic state without entanglement between different sites. Indeed, if only γ is non-zero, all hopping constants vanish and the chemical potential is µ = +1 such that we obtain the state with all B sites occupied. On the other hand, if either only α or only β is non-zero, all hopping constants vanish and the chemical potential is µ = −1 such that the MPS is equal to the state with all A sites occupied. This implies that the atomic state with occupied A orbitals is obtained from two distinct parameter configurations (α, β, γ) = (1, 0, 0) and (α, β, γ) = (0, 1, 0) whose corresponding local MPS matrices are related by a virtual unitary gauge transformation with representation matrix σ 1 .
On the other hand, if only γ and β are non-zero, the system is in a dimerized phase where the next-nearest neighbour hopping as well as the nearest-neighbour hopping between unit cells vanishes, whereas the nearestneighbour hopping within unit cells is finite. Unless |γ| = |β|, inversion symmetry is broken by the nonzero chemical potential. For |γ| = |β|, we obtain the trivial dimerized phase of the SSH model. Similarly, if only γ and α are non-zero, all couplings vanish except for the nearest-neighbour hopping between unit cells and the chemical potential. For |γ| = |α|, inversion symmetry is restored and we obtain the OAI dimerized phase of the SSH model. Finally, if all three parameters are nonzero, both nearest-neighbour hopping constants t (0) and t (1) are non-zero and there is an additional next-nearest neighbour hopping t B which is odd under inversion.

Charge pumping interpolation
We now use the MPS of Eq. (11) to describe charge pumping by considering the evolving state obtained from time-dependent parameters α(t), β(t) and γ(t) with time t ∈ (−π, π]. In Eq. (14) we saw that the MPS can represent charge pumping with a dimerized nearest-neighbour as well as a longer-range non-dimerized Hamiltonian. For  simplicity, we focus on the dimerized case. For the MPS given by the parametrisation The ground state of the charge pumping cycle H pump (t) and its correlation functions are continuous along the entire interpolation. However, the local MPS tensors parametrised by φ pump are discontinuous at the point t = ±π along the cycle: α → 0 and β → 1 as t → −π from above, whereas α → 1 and β → 0 as t → π from below. In fact, at t = ±π the system is in the atomic state with all A orbitals occupied. As discussed in the paragraph beneath Eq. (14), there are two distinct configurations for the MPS parameters corresponding to this state which are attained for t = ±π.
For the MPS we are considering, the discontinuity in the parametrisation φ pump is related to the chiral edge mode of the Chern insulator defined by the interpolation. Indeed, charge pumping in the SSH model with open boundaries corresponds to a Chern insulator on a cylinder with periodic and open boundaries in the verti-cal and horizontal directions, respectively. In the topological phase with Chern number C = 1, each edge of the cylinder hosts a one-dimensional chiral mode [21]. These modes are reflected in the ES of the SSH chain along the charge pumping cycle. For instance, the single-particle ES [52] of half the chain with periodic boundaries computed from the representation of the state as a Gaussian fermionic TNS along φ pump is shown in Fig. 2(c). It has both a left-and a right-moving mode which are localized at the two virtual edges and which are degenerate in the SSH topological phase at t = π/2.
Using the representation of the MPS in terms of local tensors, we can also compute the many-body ES of half of an infinite chain with open boundaries. Indeed, the many-body ES is isometric to the spectrum of the logarithm of the normalised and symmetrised left-and right fixed point of the MPS transfer matrix [28]. In the case of the MPS from Eq. (11), the fixed point is a matrix of dimension 2 × 2 with eigenvalues Therefore, the many-body ES has two non-trivial levels that are related by normalisation and which describe the spinless fermion at the single virtual boundary of the half-infinite chain. If and only if α 4 = β 4 + γ 4 , the two eigenvalues from Eq. (16) are degenerate. This corresponds to a crossing of the Fermi level by the edge fermion. If the interpolation describes a Chern insulator with Chern number C = 1, a degeneracy of the eigenvalues from Eq. (16) should therefore occur exactly once along the cycle. For the parametrisation φ pump this happens in the SSH topological phase at t = π/2 with α = γ and β = 0. Indeed, α 4 < β 4 + γ 4 for t < π/2 and α 4 > β 4 + γ 4 for t > π/2. Hence, the non-trivial Chern number implies that for MPS of the form of Eq. (11) the parametrisation has to be discontinuous along the cycle in order to combine the two parts of the interpolation where α 4 ≷ β 4 + γ 4 with only a single degenerate point. We emphasize that this discontinuity in the MPS interpolation is required by the topology of the charge pumping cycle, even if the pumping cycle itself is infinitely often continuously differentiable. For example, the smooth deformation of the topological pumping cycle, whose existence is discussed in Sec. II A, has an MPS ground state of the form of Eq. (11). Due to the non-trivial topology of the cycle, this MPS has a discontinuous interpolation despite the smoothness of its Bloch Hamiltonian.
The MPS of Eq. (11) can also describe cyclic interpolations corresponding to topologically trivial twodimensional models.
For instance, the MPS with parametrisation for t ∈ (−π, π] corresponds to the ground state of a model with finite chemical potential and nearest-neighbour hop-ping t (1) between unit cells, whereas the hopping within unit cells and the next-nearest neighbour hopping vanish (see Fig. 2(b)). From the single-particle ES in Fig. 2(d) we see that the conduction and valence bands in this system are not connected by the edge modes such that the Chern number is zero. Indeed, both the MPS parametrisation φ triv and the corresponding Bloch parent Hamiltonian are smooth along the cycle.

III. CHERN INSULATOR PEPS FROM SSH CHARGE PUMPING
In this section we discuss how charge pumping can be used to construct tensor networks in d + 1 dimensions starting from d-dimensional TNSs. Specifically, we show how the MPS from Sec. II for the SSH model along a charge pumping cycle leads to a two-dimensional PEPS for a Chern insulator. In Sec. III B we define a hybrid real-momentum space PEPS with finite bond dimension for the Chern insulator. In Sec. III C, we perform an inverse FT in the direction d + 1 required to transform the state to a fully real-space representation. In Sec. III D we study how the resulting state can be expressed as a real-space PEPS, whose bond dimension is investigated in Sec. III E.

A. (d + 1)-dimensional TIs from charge pumping
Charge pumping provides a systematic method to obtain a TI in d+1 spatial dimensions from a d-dimensional TI. Indeed, if the Bloch Hamiltonian of the d-dimensional model is smooth along the charge pumping interpolation as a function of the cyclic time t ∈ (−π, π], t can be identified with the momentum k d+1 in the (d + 1) st direction of the (d + 1)-dimensional system. The time-dependent Hamiltonian of the d-dimensional model then gives the Bloch Hamiltonian of the (d + 1)-dimensional system as a function of k d+1 . For instance, charge pumping in the SSH model with d = 1 defines a two-dimensional Chern insulator with Chern number C = 1 (cf. Sec. II A).
We discretize the (d + 1) st dimension with a finite number N d+1 of lattice sites. The discrete momentum values in (−π, π] are k Let us express the identification of time t and momentum k d+1 as a relation between the mode operators of the hybrid and real-space systems. The d-dimensional system has annihilation operatorsâ τ,x (t (j) ) for the physical fermion on the orbital τ = 1, . . . , f on the unit cell x ∈ Z d , where f is the number of orbitals per unit cell. They depend on the discretized time value t (j) along the pumping interpolation at which the ddimensional model is evaluated. For example, in the SSH model we have mode operatorsâ A,x (t (j) ),â B,x (t (j) ) with x = 0, . . . , N x − 1.
The (d + 1)-dimensional model obtained from charge pumping has the same number of orbitals as the ddimensional model it is derived from. For instance, the Chern insulator constructed from the SSH model also has sublattices A and B. Therefore, the physical creation operators of the (d+1)-dimensional system in real space arê a τ,(x,x d+1 ) , where x d+1 = 0, . . . , N d+1 −1 is the real-space coordinate in the (d+1) st direction and τ = 1, . . . , f . The charge pumping construction of the (d + 1)-dimensional model requires periodic boundary conditions in the direction d + 1. We may therefore consider the FT of the mode operators in the direction d + 1, while keeping the real-space coordinate x in the other d dimensions, given byâ In terms of these mode operators, the identification of time t and momentum k d+1 is the expressed by the identityâ for all d-dimensional unit cells x and sublattices τ = 1, . . . , f . For example, in the case of the SSH charge pumping we haveâ τ,x (t (j) ) =â τ,(x,k (j) y ) for τ = A, B and x = 0, . . . , N x − 1.
The Bloch Hamiltonian of the (d + 1)-dimensional model is the time-dependent Hamiltonian of the ddimensional system. Hence, with the identification of Eq. (19), the (d + 1)-dimensional ground state is given by the direct product of the many-body ground states |ψ d (t) of the ddimensional model evaluated at the N d+1 discrete times along the interpolation. From Eq. (19) it is clear that this defines the ground state w.r.t. hybrid (d+1)-dimensional real-momentum space coordinates (x, k (j) d+1 ). In order to obtain the state in terms of (d + 1)dimensional real-space coordinates (x, x d+1 ), we need to apply the inverse of the FT of Eq. (18) to the hybrid state of Eq. (20). If the d-dimensional pumping Bloch Hamiltonian is smooth as a function of time, the real-space correlation functions are guaranteed to decay faster than any polynomial [47].

B. Hybrid real-momentum space Chern PEPS
We now specialize this construction, which we explained above in terms of generic free-fermionic TIs, to d-dimensional TIs which are described by a GfTNS at all times along their charge pumping cycle. We will thereby obtain a GfTNS for the (d + 1)-dimensional TI.
For the remainder of this section and for pedagogical purposes, we will focus on the one-dimensional case (i.e. d = 1) and the corresponding notation. We will mostly rely on charge pumping in the SSH model which is described by the MPS with bond dimension D SSH = 2 from Eq. (9). Hence, we will obtain a Gaussian fermionic PEPS for the Chern insulator with Chern number C = 1. Extensions to other models and higher dimensions are straightforward.
In case the d-dimensional ground state along the pumping cycle is given as a GfTNS, Eq. (20) allows to define the (d + 1)-dimensional ground state as a hybrid TNS, where the first d dimensions correspond to real space and the (d + 1) st dimension is expressed in momentum space. Indeed, the local tensor of the hybrid TNS at the position (x, k Thus, the virtual DOFs of the hybrid TNS in the first d directions are identical to those of the d-dimensional model. In other words, the identification of time and momentum from Eq. (19) holds not just for the physical mode operators, but also for the virtual mode operators in the first d dimensions.
On the other hand, due to the direct product in the direction of k d+1 in Eq. (20), the hybrid TNS does not need virtual legs in the direction d+1, and we say that the bond dimension in this direction is equal to one (implying that the contraction of this direction corresponds trivially to a product).
The hybrid Chern PEPS obtained from the SSH charge pumping MPS of Eq. (9) is sketched in Fig. 3(a) in the hybrid coordinate system where the horizontal axis describes x and the vertical axis corresponds to k y ≡ k d+1 . The local tensor for a site on the A sublattice at the position (x, k (j) y ) is given by the SSH local tensor A i lr (t (j) ) at the time t (j) along the cycle, and similarly for the B sublattice. In the horizontal direction, the two-dimensional hybrid state hence inherits both the translation invariance and the constant finite bond dimension D SSH of the MPS. As discussed in the previous paragraph, virtual legs in the vertical direction are not needed for the hybrid Chern PEPS and are hence not shown in Fig. 3(a).

C. Inverse Fourier transform
As discussed below Eq. (20), the hybrid state obtained from charge pumping is mapped to a (d + 1)-dimensional real space coordinate system by applying an inverse FT in the direction d + 1 to the fermionic mode operators. For GfTNSs -which have virtual in addition to physical DOFs -the inverse FT should be applied to both the physical mode operators and the virtual fermionic mode operators in the first d directions. For the hybrid Chern PEPS, we therefore apply the inverse FT in the vertical direction to the physical and the horizontal virtual legs. The extension of the FT to the horizontal virtual modes amounts to a virtual basis change and does not alter the physical state, but ensures that the real-space PEPS is invariant under vertical translations y → y + 1 of its physical and virtual legs (see also Appendix F). For the hybrid Chern PEPS, the correct definition of the inverse FT in the vertical direction y entails a subtlety. Indeed, recall that the SSH charge pumping MPS from Eq. (9) is defined w.r.t. mode operators a τ,x (t (j) ) related to theâ τ,x (t (j) ) used in Eq. (19) by the particlehole transformation of Eq. (7). For the Chern PEPS, we perform an analogous particle-hole transformation in two-dimensional real space and define new mode operators The operators a A,(x,y) , a † A,(x,y) , a B,(x,y) , a † B,(x,y) span the basis in which we want to express the real-space Chern PEPS.
The SSH particle-hole transformation of Eq. (7) acts within the set of operators at time t (j) . Since we identify time and momentum, this is equivalent to a particlehole transformation acting on the modes of the Chern PEPS in momentum space rather than real space like in Eq. (21). Due to the anti-unitarity of the particlehole transformation, the inverse FT relating the modes a B,x (t (j) ) and a B,(x,y) on the B sublattice therefore contains an additional complex conjugation. Hence, the physical mode operators a τ,x (t (j) ) of the SSH models, providing the basis for the hybrid Chern PEPS, and the mode operators a τ,(x,y) for the real-space Chern PEPS are related as Here, η A = 1 and η B = −1 for the physical modes which are particle-like on the A sublattice and hole-like on the B sublattice. The inverse FT of the horizontal virtual legs of the hybrid PEPS is analogous to Eq. (22) for the physical modes. Here, the particle-or hole-like character of the virtual modes is determined by their transformation under the U(1) symmetry of Eq. (12). Specifically, the left and right virtual legs on the A sublattice transform as holes, i.e. η L,A = η R,A = −1, whereas the left and right virtual legs on the B sublattice transform as particles such that η L,B = η R,B = 1.
Let us now study how the inverse vertical FTF acts on the hybrid Chern PEPS. Due to the translation invariance in the horizontal direction, it is enough to consider one column of the hybrid state given by the sites (x, k (j) y ) for a fixed horizontal position x and j = 0, . . . , N y − 1 (see Fig. 3(a)). UnderF, the hybrid column is mapped to one column of the two-dimensional real-space state, given by the sites (x, y) with y = 0, . . . , N y − 1.
Since the inverse FTF is non-local, for a generic interpolation there are long-range correlations in the realspace column {(x, y)} 0≤y≤Ny−1 of the PEPS. In the tensor network language, the real-space column is therefore described by a single tensor obtained from the application ofF to one column of the hybrid PEPS (in contrast to the hybrid column {(x, k (j) y )} 0≤j≤Ny−1 which is described by a product of individual tensors for each k (j) y , signifying the absence of correlations in the vertical direction). For a real-space column of A sites, we denote this tensor by A col as shown in Fig. 3(b). Similarly, we define a tensor B col for a real-space column of B sites. The tensors A col and B col have N y physical legs of dimension d p = 2 and N y left and right virtual legs of dimension D SSH = 2, which describe the physical and horizontal virtual DOFs of the real-space PEPS.
The SSH pumping MPS is a Gaussian MPS for free fermions. Thus, the column tensors A col and B col are also Gaussian states, which can be described by their covariance matrices (CMs). These CMs are computed in App. F for the general d-dimensional case. There we show that the CMs of the real-space columns A col and B col , defined by the application of the inverse FT of Eq. (22) to the physical and horizontal virtual legs of the hybrid columns, are given by the inverse FT of the timedependent CMs describing the local SSH tensors A i lr (t (j) ) and B i lr (t (j) ) along the pumping cycle. This result fully characterizes the column tensors A col and B col .

D. Translation invariant real-space PEPS
From the previous subsections, we are now ready to discuss how the two-dimensional real-space state can be expressed as a TNS with local tensors for each site. For that purpose, the column tensors A col and B col have to be decomposed into a column of local PEPS tensors A 2D and B 2D with both horizontal virtual legs of dimension D SSH and vertical virtual legs with dimension D y,A and D y,B , respectively. We require the local PEPS tensors to be identical on all sites of the same sublattice as shown in Fig. 3(c) in order to make the invariance under real-space translations explicit.
The decomposition of the column into PEPS tensors can be achieved using tools developed for onedimensional MPSs [53,54]. To that end, we define one-dimensional pure states |ψ 1D (A col ) and |ψ 1D (B col ) , whose many-body wave functions are the tensor elements of A col and B col , respectively. Therefore, |ψ 1D (A col ) effectively describes a chain of length N y , where the local Hilbert space at each site of the fictitious chain is the tensor product of the physical and horizontal virtual Hilbert spaces of A col at the corresponding physical site (and similarly for |ψ 1D (B col ) and B col ). Pictorially, |ψ 1D (A col ) and |ψ 1D (B col ) are obtained by merging at each site the physical and horizontal virtual legs of A col and B col into an effective physical index of dimensiond p = d p D 2 SSH per site.
Our goal is now to express |ψ 1D (A col ) and |ψ 1D (B col ) as translation invariant MPSs with periodic boundary conditions. Then, after unfolding the effective physical index into the physical and horizontal virtual indices of the PEPS, the local tensors of these MPSs will define the PEPS tensors A 2D and B 2D , respectively. Similarly, the MPS bond dimensions D y,A and D y,B are equal to the vertical bond dimensions of the PEPS local tensors A 2D and B 2D , respectively.
Numerically, the translation invariant MPSs for |ψ 1D (A col ) and |ψ 1D (B col ) can be obtained as follows. Therefore, the bond dimension grows at least linearly with the system size [54].

E. Vertical PEPS bond dimension
Above, we showed that the inverse FT of the hybrid PEPS can be decomposed into a (d+1)-dimensional realspace TNS, where the bond dimension in the first d directions is equal to that of the original d-dimensional state. The core question is the scaling of the bond dimension in the (d + 1) st real space dimension w.r.t. the system size in this direction.

Lower bound from ES
The bond dimension of a tensor network is intimately related to the entanglement between its physical particles: For a TNS describing a physical system with a sub- Schmidt decomposition of the one-dimensional effective state |ψ1D for the PEPS column with physical dimensioñ dp into the subsystem AL of the first L sites and its complement. When |ψ1D is written as a translation invariant MPS with bond dimension Dy, the virtual boundary of the subsystem AL, marked in blue, crosses two bonds of the MPS.
system A, the total dimension of all virtual legs at the boundary ∂A can be no smaller than the rank of the state's Schmidt decomposition into the DOFs of A and its complement. This constraint allows us to infer a lower bound for the vertical bond dimension of the real-space PEPS from the ES of the one-dimensional column states.
We consider a generic one-dimensional effective column state denoted |ψ 1D with an effective physical dimensioñ d p . Below, we will choose either |ψ 1D = |ψ 1D (A col ) or SSH . We assume that |ψ 1D is expressed as a translation invariant MPS with a bond dimension D y corresponding to the vertical PEPS bond dimension (see Fig. 4). Such a representation can be obtained for instance using the procedure described in the previous subsection. We now perform a Schmidt decomposition of the pure state |ψ 1D w.r.t. the subsystem A L of the first L sites 0 ≤ y ≤ L − 1, and denote the rank of the decomposition by r L (|ψ 1D ). As sketched in Fig. 4, the cut between A L and its complement crosses exactly two virtual bonds of the MPS, namely those between sites (L − 1, L) and (N y − 1, 0), each of dimension D y . Since the virtual dimension can be no smaller than the Schmidt rank for any subsystem size L, we obtain the bound The Schmidt rank r L (|ψ 1D ) may be obtained by counting the non-trivial levels in the single-particle ES of the Gaussian column state |ψ 1D . This ES can be computed numerically based on the explicit result for the CM of |ψ 1D given in Appendix F. Since the state |ψ 1D with physical dimensiond p has log 2 (d p ) free fermionic modes per lattice site, its single-particle ES w.r.t. the subsystem A L consists of levels 0 ≤ |λ i (L)| ≤ 1 with 1 ≤ i ≤ L log 2 (d p ) (see Eq. (A5) for a definition of the single particle ES and its levels |λ j |). For an improved numerical stability, we actually compute µ i (L) = 1 − λ i (L) 2 that may be conveniently extracted from the CM (see Eq. (A6)). Here, values µ i (L) = 1 and µ i (L) = 0 correspond to maximally entangled and non-entangled modes, respectively. The many-body Schmidt rank of the state |ψ 1D is therefore given by the exponential of the number of entangled modes in the single-particle ES with µ i (L) bigger than zero.

Exponential growth
The number of entangled modes in the single-particle ES of a column A col of unit cells of the real-space PEPS is displayed in Fig. 5(a) for two different cyclic interpolations of the MPS from Eq. (11): on one hand, the parametrisation φ pump of Eq. (15) giving rise to a twodimensional Chern insulator, and on the other hand the parametrisation φ triv of Eq. (17) corresponding to a topologically trivial two-dimensional state. The corresponding data of the column B col is identical.
For φ pump the number of entangled modes is equal to when L ≤ N y /2 (the spectra for L and N y − L are identical). In Appendix G we show that this is the maximal number of entangled modes which is compatible with the global U(1) symmetry of A col inherited from the SSH model MPS. Here, the factor 3 in Eq. (25) is a consequence of the column tensor having three fermions (one physical fermion and two virtual fermions) per site. The validity of Eq. (25) can be seen in Fig. 5(a) for the smallest system size N y = 28. For larger system sizes, the number of entangled modes shown in Fig. 5(a) is lower than Eq. (25) because of the finite numerical resolution, as shown by the ES in Fig. 5(b). We have checked that for the Chern PEPS obtained from the MPS corresponding to the smooth charge pumping cycle described in Sec. II A, the number of entangled modes is also given by Eq. (25). For φ triv , the number of entangled modes is given by 2L when L ≤ N y /2. This is lower than Eq. (25) for φ pump except at L = N y /2, where both numbers agree. As discussed in Appendix G, this reduction of entangled modes is due to the decoupling of all right virtual legs and is not related to topology. Indeed, for the trivial cycle, the parameter β(t) = 0 vanishes for the entire duration of the interpolation φ triv . Any small perturbation of φ triv by a non-zero β increases the number of entangled modes to Eq. (25). Moreover, these additional decoupled modes disappear when blocking the two columns A col and B col , i.e. by contracting the bonds connecting them. Hence, the column for AB unit cells has the same Schmidt rank for both φ pump and φ triv . The number of entangled modes from Eq. (25) corresponds to a maximal Schmidt rank max 1≤L≤Ny−1 for the state |ψ 1D (A col ) , and similarly for |ψ 1D (B col ) . The bound of Eq. (23) therefore implies that the vertical bond dimensions D y,A and D y,B grow exponentially as a function of N y . This is a generic feature of our construction and not related to the topology of the state: It originates in the non-locality of the inverse FT which couples all states in the exponentially growing Hilbert space of one column. Indeed, a generic quantum state requires an exponentially growing bond dimension to be represented exactly as a TNS.
Despite the faster decay of the single-particle entanglement energies for the topologically trivial state than for the Chern insulator in Fig. 5(b), this does not allow us to make any statement about differences in the growth of the vertical bond dimension required for an approximative PEPS for the two systems.

IV. FERMIONIC PEPS FOR TWO-DIMENSIONAL HIGHER ORDER TI
In this section, we study a Gaussian fermionic PEPS for the topological quadrupole model from Ref. [1], which is reviewed in Sec. IV A. The PEPS is defined in Sec. IV B, where we also provide its parent Hamiltonian. In Sec. IV C we discuss a three-dimensional PEPS with chiral hinge states obtained from a dipole pumping interpolation of the quadrupole model [2,3].

A. Second order quadrupole insulator
A two-dimensional second-order topological phase with a quantized bulk quadrupole moment was recently proposed theoretically [1] and has subsequently been observed experimentally in mechanical [11], photonic [14][15][16][17] and electrical [18,19] systems. This phase is described by a microscopic free-fermionic model with one spinless fermionic mode per site and a unit cell of 2 × 2 sites depicted in Fig. 6(a). We consider the system at half-filling where only the lowest two bands are occupied. With open boundary conditions, the nearest-neighbour Hamiltonian is given by wherex andŷ denote the unit vectors in the horizontal and vertical direction, respectively. The positions of the unit cells are x = xx + yŷ with 0 ≤ x ≤ N x − 1 and 0 ≤ y ≤ N y − 1 on a lattice with N x and N y unit cells in the horizontal and vertical direction. For τ = 1, . . . , 4, a † τ,x denotes the creation operator for a fermion in orbital x and t (1) y , where the subscripts x and y refer to hoppings in the horizontal and vertical direction, and the superscripts (0) and (1) indicate hopping between sites on the same and adjacent unit cells, respectively. The signs of the hopping amplitudes ensure that there is a flux π through every plaquette of the square lattice.
The quadrupole model of Eq. (27) is in a second-order topological OAI phase protected by the horizontal and vertical mirror symmetries M x and M y if t (0) . In this phase, corners host gapless protected states at the intersection of two gapped edges. Moreover, the system possesses a quantized bulk quadrupole moment, edge dipole moment and corner charge [3]. When the hopping amplitudes t is in a trivial phase with gapped edges and corners. If the hoppings t (1) x and t (1) y between adjacent unit cells vanish, the system is in a trivial dimerized phase where each unit cell decouples from the rest of the system.
Similarly to the charge pumping cycle discussed in Sec. II A, a dipole pumping cycle interpolating between the trivial and OAI dimerized phases of the quadrupole model can be defined by adding the chemical potential µ to all sites on the sublattices 1 and 2 and −µ to all sites on the sublattices 3 and 4 [3]. Therefore, the staggering pattern breaks the symmetries protecting the topological phase but preserves the C 2 rotation symmetry (see Fig. 6(a)). The dipole pumping cycle is obtained from C 4 symmetric hopping amplitudes t (0) ≡ t (0) y and a chemical potential µ evolving according to the same cyclic interpolation of Eq. (3) as for the SSH charge pump.
The properties of the dipole pumping cycle for the quadrupole model are analogous to those of the charge pumping cycle for the SSH model. In particular, the system remains in a dimerized state throughout the interpolation since at all times t ∈ [−π, π] only one of the two hopping amplitudes t (0) (t) and t (1) (t) is non-zero. At t = ±π and t = 0, the system is in an atomic state with only the orbitals 1 and 2 occupied at t = ±π and only the orbitals 3 and 4 occupied at t = 0. For −π < t < 0, the hopping within unit cells is non-zero, whereas the hopping between unit cells is non-zero for 0 < t < π. In both cases, charge gets transferred by the changing chemical potential. At t = −π/2 and t = π/2, the chemical potential vanishes such that the mirror symmetries are restored and the Hamiltonian corresponds to the trivial and OAI dimerized phase of the quadrupole model, respectively.
In the same manner as charge pumping relates the SSH model to a Chern insulator, dipole pumping induces a model with chiral hinge states from the quadrupole model [2,3]. The chiral hinge model is a threedimensional second order topological insulator whose one-dimensional protected boundary modes occur at the intersection of a pair of two-dimensional faces. The topology of the hinge model obeys a Z 2 classification protected by the product M x M y T of the horizontal and vertical mirror symmetries and time reversal T [2].

B. PEPS for the quadrupole model
In this subsection, we provide a Gaussian fermionic PEPS for the ground state of the quadrupole model. After giving the details of the construction in Sec. IV B 1, we compute its parent Hamiltonian in Sec. IV B 2 and discuss its ES in Sec. IV B 3.

Construction
In order to construct a Gaussian PEPS for the ground state of the quadrupole model, we use similar ideas as in the SSH model MPS derived in Sec. II. In analogy to the particle-hole transformation of Eq. (7) on the B sublattice of the SSH chain, we define new physical modes a τ,x by performing a particle-hole transformation on the sublattices 3 and 4 while leaving the sublattices 1 and 2 unaffected, The vacuum |Ω of the new modes, satisfying a τ,x |Ω = 0 for τ = 1, . . . , 4, contains exactly 2N x N y physical particles as required for the quadrupole model at half-filling. Thus, the particle-hole transformation allows us to express the quadrupole PEPS in terms of a separate and parity-even local tensor for each lattice site. Parityevenness is required in order to ensure that the state is independent of the order of contractions in the network [56].
The quadrupole PEPS has bond dimension D Quad = 2 and is constructed from four local tensors A Due to the close relation between the SSH and the quadrupole models, we are guided in our ansatz for the local tensors of the quadrupole PEPS by the MPS tensors of Eq. (11). This can be most easily seen from Fig. 6(b). We obtain the quadrupole model by coupling neighboring sites both horizontally and vertically according to the pattern of an SSH chain. In the PEPS tensors, the couplings of the horizontal SSH chains are transmitted by the left and right virtual fermions associated to parameters α x , β x analogous to Eq. (11). Similarly, the vertical SSH chains are implemented using the top and bottom virtual fermions associated to the parameters α y , β y . Using the analogy to Eq. (11), we therefore obtain local PEPS tensors whose non-vanishing elements are given in terms of five real parameters γ, α x , α y , β x and β y as All the other tensor elements are equal to zero. The phases in Eq. (29) were chosen such as to ensure that there is a flux π through every plaquette. As sketched in Fig. 6(b), the parameters β x and β y represent the coupling of the physical leg to the virtual legs corresponding to the horizontal and vertical bonds pointing into the unit cell, respectively. Similarly, α x and α y control the coupling of the physical leg to the virtual legs pointing out of the unit cell.
The PEPS of Eq. (29) has a global U(1) symmetry analogous to Eq. (12) for the SSH MPS. Indeed, the local tensors are invariant under a combination of U(1) rotations of the physical and virtual legs given by for sites on the sublattices τ = 1, 2, and for sites on the sublattices τ = 3, 4. Here, U (ϕ) = 1 0 0 e iϕ is the U(1) rotation acting on a single fermion. Eq. (30) implies that the virtual legs for the local tensors on the sublattices τ = 1, 2 and τ = 3, 4 transform as holes and particles. Hence, each pair of virtual legs associated with the same nearest-neighbour bond transforms oppositely under the U(1) rotation, such that the PEPS is invariant under the physical part of Eq. (30). The charge associated with this symmetry is such that the U(1) symmetry ensures that the state lies exactly at half-filling of the lattice just as in the onedimensional case.

Parent Hamiltonian
For all values of the parameters γ, α x , α y , β x and β y , the PEPS from Eq. (29) can be expressed as a Gaussian TNS for free fermions. In Appendix H, we show that the PEPS with parameters α ≡ α x = α y and β ≡ β x = β y is the unique ground state of an extended version of the quadrupole model H PEPS with C 4 symmetric hoppings, with a staggered chemical potential that breaks C 4 symmetry and with an additional next-to-nearest neighbour hopping between sites on the same sublattice τ with amplitude t τ . The pattern of next-to-nearest neighbour hoppings is shown in Fig. 6(c). The couplings of H PEPS depend on the parameters of the PEPS as Similarly to the parent Hamiltonian H MPS of the SSH model MPS of Eq. (11), H PEPS describes different phases depending on the values of the parameters α, β and γ. When two out of the three parameters vanish, the system is in an atomic insulator state: If α = β = 0, the particles are localized on the orbitals τ = 3, 4, whereas they are localized on the orbitals τ = 1, 2 if α = γ = 0 or β = γ = 0. In contrast, when α = 0 and β, γ = 0, all hoppings of the parent Hamiltonian vanish except for the nearest-neighbour hopping within unit cells. Similarly, when β = 0 and α, γ = 0, the only non-zero hopping is the nearest-neighbour hopping between unit cells. In both cases, the system is in a dimerized phase with a staggered chemical potential. Setting β = γ and α = γ, respectively, we recover the trivial and OAI dimerized phases of the quadrupole model with vanishing chemical potential. Finally, when all three parameters are nonzero, the system has a non-vanishing nearest and nextnearest neighbour hopping as well as a finite staggered chemical potential.

Corner states and entanglement spectrum
The characteristic (d − 1)-dimensional gapless edge states of a conventional d-dimensional TI are reflected in the state's bulk ES [22]. Similarly, the ES of the quadrupole model in its dimerized OAI phase hosts gapless corner states as long as the virtual cut is compatible with the protecting symmetries [24]. In order to further characterize the different phases described by the PEPS from Eq. (29), we therefore study the ES of the state defined on a torus w.r.t. a rectangular subsystem A Lx×Ly of L x and L y unit cells in the horizontal and vertical direction.
For β = 0, the PEPS from Eq. (29) is in a dimerized phase where the system splits into four-site plaquettes shifted from the unit cell by one site in the horizontal and vertical direction. For α = γ, we obtain the OAI dimerized phase of the quadrupole model, whereas for α = γ there is a non-zero staggered chemical potential that breaks the symmetries protecting the OAI phase. Due to the dimerization, with open boundaries the system has SSH chains with a staggered chemical potential at the edges, and corner sites which decouple from each other and the bulk. Correspondingly, the singleparticle ES of the PEPS w.r.t. the subsystem A Lx×Ly has three distinct contributions from the bulk, the edges and the corners as sketched in Fig. 7. The bulk consists of (L x − 1)(L y − 1) plaquettes decoupled from the rest of the system, which contribute 4(L x − 1)(L y − 1) nonentangled modes with levels |λ bulk | = 1 (in the dimerized limit) to the single-particle ES.
On the other hand, the corner and edge sites belong to plaquettes crossed by the boundary of the subsystem A Lx×Ly . In Appendix H 3 we show that the four boundaries provide 4(L x − 1 + L y − 1) entangled modes with levels ±λ edge , where Since the number of these entangled modes grows linearly with the size of the edge, the boundary leads to an area law term in the entanglement entropy. Furthermore, each corner site hosts exactly one mode corresponding to a level ∓λ corner with Here, the negative sign holds for the top right and bottom left corner sites on the sublattices τ = 1 and τ = 2, respectively, whereas the positive signs applies to the top left and bottom right corner sites on the sublattices τ = 3 and τ = 4, respectively. For α = γ, the PEPS of Eq. (29) describes the OAI phase of the quadrupole model. Indeed, in this case the edges have entanglement levels λ edge = 1/ √ 2 and the four corners have degenerate levels λ corner = 0 corresponding to maximally entangled corner modes. On the other hand, for α = 0, both λ edge = 1 and λ corner = 1 such that the TNS describes an atomic state.
Finally, if the OAI dimerized phase is perturbed by a small non-zero value β = 0, one can check numerically that the corner modes acquire a finite splitting λ corner = 0 and are no longer perfectly localized at the corner sites. This confirms that the TNS with β = 0 is not in the OAI phase of the quadrupole model as expected from the breaking of the mirror symmetries and C 4 symmetry by the next-nearest neighbour hopping and the chemical potential.

C. 3D chiral hinge PEPS from dipole pumping
The PEPS of Eq. (29) with parameters α = α x = α y and β = β x = β y describes a dimerized dipole pumping cycle of the quadrupole model if α, β and γ follow the parametrisation φ pump which we found for the dimerized charge pumping cycle of the SSH model MPS. The coupling constants of the parent Hamiltonian H PEPS derived from the parametrisation φ pump are shown in Fig. 8(a). They differ from those of the dipole pumping Hamiltonian of Ref. [3] only by a factor of 1/ √ 2 for the hopping amplitudes, which does not affect the topology of the interpolation. The single particle ES of the PEPS along φ pump is shown in Fig. 8(b). In the first half of the cycle, the system is in a dimerized phase with each unit cell decoupled from the rest of the system. Hence, the single-particle ES contains only the bulk bands with λ = ±1. However, in the second half of the cycle for 0 < t < π, the contributions from the edges and corners can be clearly distinguished. The four corner modes connect the bands with λ = ±1 and are degenerate in the OAI dimerized phase of the quadrupole model obtained for t = π/2. In the three-dimensional second-order TI of Ref. [2] obtained by the identification of the time t along the dipole pumping cycle with the momentum k z in the third direction, these corner modes generate the chiral modes localized at the hinges. Following the steps described in Sec. III, we may use the PEPS of Eq. (29) along the interpolation φ pump to define a three-dimensional PEPS for the second order hinge TI. Moreover, we can construct a topologically trivial three-dimensional state from the interpolation φ triv from Eq. (17). These PEPSs have a finite bond dimension D Quad = 2 in the x and y directions. In the hybrid realmomentum space where the third dimension corresponds to the momentum k z or time t, the states have a finite bond dimension equal to one also in the third direction.
On the other hand, due to the non-locality of the inverse FT, in real space their bond dimension D z,τ in the third direction for sites on the sublattice τ for τ = 1, . . . , 4 grows with the system size N z . As in Sec. III E, we can estimate D z,τ from the ES of a column of sites on the sublattice τ w.r.t. the subsystem A L of the first L sites.
The number of entangled levels in the single-particle ES of a column of sites on the sublattice τ = 1 of the three-dimensional PEPSs is shown in Fig. 9(a). As we can see from the smallest system size N z = 28, for both φ pump and φ triv the number of entangled modes is identical to the two-dimensional case from Sec. III when replacing N y with N z . Moreover, the spectra displayed in  As we show in Appendix G, for φ pump this number of entangled modes is the maximal number compatible with the symmetries of the quadrupole PEPS, here the U(1) symmetry from Eq. (30), and the mirror symmetry M xy of the local tensor on the sublattice τ = 1. Indeed, the latter causes the decoupling of one superposition of the left and down virtual legs, and similarly for the up and right virtual legs. For φ triv , there is an additional reduction of the number of entangled modes since the left and down virtual legs decouple due to the vanishing parameter β = 0 along the trivial interpolation (see Appendix G).
Analogously to Sec. III E, we therefore conclude that the bond dimension D z,τ in the third direction grows exponentially with N z for τ = 1, 2, 3, 4. Due to the mirror symmetries of the quadrupole model, D z,τ has the same value as the vertical bond dimension of the twodimensional PEPSs obtained from cyclic interpolations of the SSH model. The increase in spatial dimensionality therefore does not cause an increase of the bond dimension in the (d + 1) st direction.

V. CONCLUSION
In this article, we showed how to use charge pumping to define TNSs for (d + 1)-dimensional conventional or higher-order TIs starting from TNSs of TIs in d space dimensions. To that end, we constructed a Gaussian fermionic MPS for the SSH model with bond dimension D SSH = 2 in d = 1 dimension, and a Gaussian fermionic PEPS for the topological quadrupole model with bond dimension D Quad = 2 in d = 2 dimensions. We proved that these TNSs have local gapped parent Hamiltonians with up to next-nearest neighbour hopping, and thereby showed that they describe the SSH model along a charge pumping cycle and the quadrupole model along a dipole pumping cycle, respectively. We employed these TNSs to construct a two-dimensional PEPS for a Chern insulator and a three-dimensional PEPS for a chiral hinge HOTI, respectively. The (d + 1)-dimensional TNSs inherit the finite bond dimension D SSH and D Quad in the first d dimensions, respectively. In a hybrid coordinate system where the (d + 1) st dimension corresponds to momentum, the (d+1)-dimensional TNSs have a finite bond dimension in this direction. In contrast, we showed that the bond dimension in the (d + 1) st direction grows exponentially in a real-space coordinate system.
Our results suggest several directions for future work. On one hand, it would be interesting to study if a realspace PEPS for the Chern insulator with a polynomi-ally growing bond dimension can be found by truncating the Schmidt values of the real-space column in our construction. Such a result could potentially provide insight into the physical origin for the obstructions preventing the existence of chiral PEPSs with a finite bond dimension. On the other hand, we expect the TNSs constructed here to be useful for finite-size simulations despite their growing bond dimension, since the bond dimension is finite in all but one direction. By their local nature, they could be employed as the building block of interacting (d + 1)-dimensional TIs obtained by Gutzwiller projection or parton constructions [57][58][59][60]. Eigenstates and thermal states of free-fermion systems are given by Gaussian states which satisfy Wick's theorem. They are hence fully characterised by their covariance matrix (CM) of two-point correlation functions [51]. In this appendix, we review the definitions of the complex and real CM for pure and mixed states in Appendix A 1. We summarise the relation with the entanglement spectrum in Appendix A 2 and provide the concrete expression for the CM of a Gaussian state parametrised as the exponential of a fermion bilinear in Appendix A 4.

Definitions
We consider a generic system of N fermionic DOFs with annihilation operators a j and creation operators a † j for j = 1, . . . , N . In a pure or mixed state of this system, its covariance matrix is defined as where µ, ν = 1, . . . , 2N and χ = (a 1 , . . . , a N , a † 1 , . . . , a † N ) is the mode vector. The blocks R and Q of dimension N × N in Eq. (A1) are anti-Hermitian and anti-symmetric, respectively, i.e. R † = −R and Q T = −Q, such that G is anti-Hermitian. For a generic mixed state, G satisfies the inequality GG † ≤ 1 4 1 and its eigenvalues come in complex conjugate pairs ± i 2 |λ j | with 0 ≤ |λ j | ≤ 1 for 1 ≤ j ≤ N . For pure states, GG † = 1 4 1 such that its eigenvalues are given by |λ j | = 1.

Relation to entanglement spectrum
We frequently consider the restriction of a pure quantum state |ψ to a subsystem A of the entire system, which is described by the reduced density matrix ρ A = TrĀ[|ψ ψ|] obtained by tracing over the DOFs in the complementĀ of A. The many body entanglement Hamiltonian H Ent w.r.t. this partition is given by the logarithm of the reduced density matrix as [23] where Z = Tr e −HEnt . If |ψ is a Gaussian state with CM G, then ρ A is Gaussian with its CM G A = (G µν ) µ,ν∈A given by the restriction of G to the modes of A. In this case, the manybody entanglement Hamiltonian is a bilinear function of the fermionic mode operators defined by a square matrix referred to as the single-particle entanglement Hamiltonian. The eigenvalues β j of the single-particle entanglement Hamiltonian are related to the eigenvalues ± i 2 |λ j | of G A as |λ j | = tanh βj 2 , where 1 ≤ j ≤ L and L is the number of modes in the subsystem A [52]. We therefore refer to the collection {|λ j |} 1≤j≤L (A5) as the single particle ES of ρ A . The single particle ES can also be computed from the Majorana CM Γ of Eq. (A3) for the Gaussian state |ψ . One obtains the levels |λ j | for 1 ≤ j ≤ L directly as the singular values of the block (Γ µν ) µ,ν∈A associated with the DOFs in A. On the other hand, the levels can be computed indirectly from the off-diagonal block (Γ µν ) µ∈A,ν∈Ā describing the correlations between A and its complementĀ. Indeed, the two blocks are coupled by the constraint ΓΓ T = 1 originating in the purity of the state |ψ . Concretely, the singular values of the offdiagonal block are given by [63]. This allows to infer the single-particle entanglement energies for weakly entangled modes |λ j | ≈ 1 with an improved numerical accuracy.

Gaussian projections and Schur complements
In the construction of GfTNSs to be discussed in Appendix B, we will frequently encounter the following situation: Let us take a total system of N + M fermionic modes with mode operators a j , a † j for j = 1, . . . , N + M , and consider the subsystem defined by the M last modes j = N + 1, . . . , N + M . Let |Q be a Gaussian pure state of the total system and |ω a Gaussian pure state of the subsystem of the last M modes. Then, the projection which is a state of only the first N modes, is again a Gaussian state. The Majorana CMs Γ ψ , Γ Q and Γ ω of the three Gaussian states are given by the expression in Eq. (A3), where the expectation values are taken in the states |ψ , |Q and |ω , respectively. Γ ψ can be computed from Γ Q and Γ ω using a Schur complement coming from the Gaussian integration over the last M modes [51,[64][65][66]. Indeed, we can write the anti-symmetric Majorana CM Γ Q of the un-projected state of the total system as where the block A of size 2N ×2N refers to the Majorana modes corresponding to a j , a † j for j = 1, . . . , N , and the block D of size 2M × 2M refers to the Majorana modes corresponding to a j , a † j for j = N + 1, . . . , N + M . The block B of size 2N ×2M encodes the correlations between the first N and last M modes. In terms of these blocks and the CM Γ ω of size 2M ×2M , the CM of the projected state is Mathematically, this is the Schur complement of the block corresponding to the last M modes.

Parametrisation of Gaussian states
We consider a normalised Gaussian state with even fermionic parity which, anticipating Appendix B, we will denote by |ψ = |Q . |Q is parametrised in terms of the anti-symmetric complex matrix M of dimension N × N as [51,67] where |Ω is the fermionic vacuum and N is a normalisation factor. An relevant case is given by states whose parametrisation matrices have non-zero entries only in the first row and column (or more generally, in the n th row and column). In particular, we will see in Appendix E and Appendix H that the local tensors of the MPS of Eq. (11) and the PEPS of Eq. (29) define Gaussian states of this form. In this case, M ij M kl a † i a † j a † k a † l = 0 for all 1 ≤ i, j, k, l ≤ N since all non-vanishing terms contain the factor (a † 1 ) 2 = 0 ((a † n ) 2 = 0 in the general case). Therefore, the series expression of the exponential in Eq. (A10) terminates after first order such that leading to the CM

Appendix B: Construction of GfTNS
Gaussian fermionic tensor network states (GfTNS) describe free fermion systems. They have the advantage that the contraction of the network can be performed in terms of covariance matrices, allowing efficient computations even for large systems. In this appendix, we review the construction of Gaussian TNSs in one and two spatial dimensions. We begin in Appendix B 1 by illustrating the construction of a TNS via fiducial states using the simple example of a bosonic MPS, where we do not need to take care of fermionic signs in tensor products. In Appendix B 2, we move to the case of fermionic physical and virtual particles. In Appendix E, the concepts introduced here are illustrated using the pedagogical example of the SSH model MPS from Sec. II C.

Fiducial state approach
We recall that the construction of a TNS can be performed in several different but equivalent ways. Within the most well-known approach, the local information on the state is encoded in the local tensor with physical and virtual legs. For example, for a one-dimensional bosonic MPS |ψ bMPS the local tensor at position x takes the form A[x] ix lxrx , where i x is the physical index, l x is the left virtual index and r x is the right virtual index. The global state is obtained by contracting, for each nearest-neighbour bond, the two virtual legs associated with this bond which belong to neighboring tensors. This is done by identifying and summing over the corresponding virtual indices. For example, the contraction of two neighboring MPS tensors gives lx+1rx+1 . More formally, this can be expressed as follows: First, we create a total maximally entangled state in the virtual layer whose role it is to implement the contraction of bonds. The total maximally entangled state is the tensor product over all nearest-neighbour bonds of a maximally entangled state of the two virtual particles for this bond. In our MPS example, the maximally entangled state for the bond between sites x and x + 1 is Secondly, at each site we translate the local tensor into a local projection map, which maps the virtual particles onto the physical particle at this site. The representation matrix of the local projection map is given by the local tensor. For instance, the MPS local projection map at site x isÂ The bond between sites x and x + 1 is contracted by applying the projection mapsÂ[x] andÂ[x + 1] to |ω x,x+1 , givinĝ The global state is then obtained by contracting all bonds, i.e. by applying the product of all local projection maps to the total virtual maximally entangled state. For instance, the MPS on a chain with N x sites and periodic boundaries takes the form where the lattice site indices are x = 0, . . . , N x − 1. This is a state of only the physical particles. For GfTNSs, we will deal with fermionic particles and therefore follow a slightly different, but equivalent, approach to construct the TNSs using so-called fiducial states (see Refs. [47,67] for pedagogical introductions). In this approach, we consider fiducial states instead of local projection maps, which contain only creation operators. The fiducial state on each lattice site lies in the joint Hilbert space of the physical and virtual particles on this site. Its basis coefficients are given by the entries of the local tensor. For example, the local fiducial states of the MPS are Hence, the fiducial states are equivalent to the local projection maps or local tensors; in particular, they contain all the local information about the TNSs. The total state is obtained by projecting the tensor product of all local fiducial states on the total virtual maximally entangled state. For the MPS, we have which is equivalent to the expression in Eq. (B4). We are now ready to generalize the fiducial state formalism to fermionic physical and virtual particles.

Fermionic particles
We consider a one-or two-dimensional lattice system of free fermions with f fermionic modes per lattice site x, which are associated to the physical creation operators a † τ,x with τ = 1, . . . , f . A Gaussian fermionic TNS for this system with physical dimension 2 f and bond dimension 2 ξ is obtained by associating ξ complex virtual fermionic modes with each physical lattice site x and nearest-neighbour direction α, where α = L, R for onedimensional MPSs with left and right nearest neighbour bonds and α = L, U, R, D for two-dimensional PEPSs on the square lattice with left, up, right and down nearestneighbour bonds [33,67]. We denote the creation operators for these virtual modes by b † α,j,x where j = 1, . . . , ξ labels the different modes per bond and lattice site. For each lattice site, there are hence n modes = f + 2ξ modes for a fermionic MPS and n modes = f + 4ξ modes for a two-dimensional fermionic PEPS. The mode operators needed for the construction of the SSH MPS are discussed in the first paragraph of Appendix E.
The local information about the TNS is contained in the local fiducial states |Q x , introduced in Eq. (B5), which are equivalent to the local tensors A[x] used in the main text. We can easily translate between the two approaches, since the basis elements of the local fiducial state are by definition equal to the local tensors (cf. Eq. (B5)). For fermions, we write |Q x = Q x |Ω where Q x is a polynomial of creation operators which acts on the vacuum to create the fiducial state. For fermionic MPSs and PEPSs in one and two dimensions we have [49,56] A where |i with i = 0, . . . , 2 f − 1 is a basis for the Fock space associated with the physical mode operators {a † τ,x } 1≤τ ≤f , and |l with l = 0, . . . , 2 ξ − 1 is a basis of the Fock space associated with the left virtual mode operators {b † L,j,x } 1≤j≤ξ on site x (similarly |u for the up mode operators {b † U,j,x } 1≤j≤ξ , |r for the right mode operators {b † R,j,x } 1≤j≤ξ and |d for the down mode operators {b † D,j,x } 1≤j≤ξ ). For the SSH model MPS, the fiducial states obtained thus are given in Eq. (E1).
A Gaussian fermionic TNS has the property that all local fiducial states |Q x satisfy Wick's theorem. In this case, the global physical state |ψ is also Gaussian [33,51]. We denote by Γ Q the Majorana CM of the product of all fiducial states, x Q x |Ω , also referred to as total fiducial state. From now onwards, we consider GfTNSs with parity-even local tensors as discussed in Sec. II B, whose local fiducial states therefore have an even number of physical and virtual fermions. The maps Q x can thus be expressed as in Eq. (A10) as the exponential of a quadratic form of the physical and virtual creation operators on the site x, which are contained in the last n modes entries (χ x ) n modes +m 1≤m≤n modes of the mode vector of Eq. (B7) [51,67]. Concretely, the local fiducial state is parametrised as with an anti-symmetric square matrix M x of dimension n modes . For the SSH model MPS, this matrix is given in Eq. (E2).
In order to illustrate the concepts introduced above, let us construct a Gaussian maximally entangled state of the virtual fermions for each nearest-neighbour bond x x . Let us denote by α and α the type of virtual fermion involved in the bond x x on the site x and x , respectively. For instance, if x = x +x, then α = R and α = L. A fermionic maximally entangled state for this bond is then given by [33] Due to the fermionic anti-commutation relations, this expression is not symmetric under exchange of (α , x ) and (α , x ); we say that the bond points from the initial site x to the final site x . The state |ω x x is a Gaussian state satisfying Wick's theorem. To see this, we first consider the simple case of a one-dimensional MPS with ξ = 1 virtual fermion per nearest-neighbour bond and lattice site. In this case, the virtual maximally entangled state from Eq. (B13) for the bond x, x +x becomes where the bond points from site x to site x +x. This is a Gaussian state of the form of Eq. (A11) parametrised by the anti-symmetric matrix M = (i/2)σ 2 . Hence, its complex CM G can be computed from Eq. (A12) and is given by where σ 2 = 0 −i i 0 denotes the second Pauli matrix. According to Eq. (A3), the corresponding Majorana CM Γ is obtained by conjugation with the matrix S, and is found to be where σ 1 = 0 1 1 0 denotes the first Pauli matrix. For a two-dimensional TNS with ξ = 1, the real CM of the virtual state for the vertical nearest-neighbour bond x, x −ŷ is also given by Eq. (B16). In this case, we choose the bond to be oriented downwards from x to x −ŷ to comply with our ordering of the virtual legs as L, U, R, D. For TNS with more virtual fermions, ξ > 1, the virtual maximally entangled state from Eq. (B13) is a tensor product of multiple states of the form of Eq. (B14). Hence, the corresponding real CM is a direct sum of multiple copies of the Γ from Eq. (B16). We denote by Γ ω the Majorana CM of the total virtual maximally entangled state ⊗ x x |ω x x . Since the total virtual maximally entangled state is a tensor product over all bonds, Γ ω is a direct sum of multiple copies of the Γ from Eq. (B16).
As explained in Appendix B 1, the global physical state |ψ is obtained from the constituents introduced above by projecting the total fiducial state onto the virtual maximally entangled state, We recall from Appendix A 3 that for Gaussian states, this projection can be formulated in terms of CMs and gives rise to a Schur complement. This is the approach we take in the following.
Let us see how we can apply the general Schur complement expression of Eq. (A9) in order to compute the Majorana CM Γ |ψ for the physical state from the CMs Γ Q and Γ ω defined above for the fiducial and virtual maximally entangled states. We introduce the symbols p and v to collectively refer to all Majorana mode operators for the physical and virtual fermions. Therefore, (Γ Q ) pp and (Γ Q ) vv denote the blocks of the CM of the total fiducial state that describe the reduced state of only the physical and virtual degrees of freedom (DOFs), respectively. On the other hand, the blocks (Γ Q ) pv = − (Γ Q ) T vp encode the correlations between physical and virtual fermions. In Eq. (B17), we are projecting the fiducial state of physical and virtual modes onto a maximally entangled state of the virtual modes in order to obtain a state of only the physical modes. According to Eq. (A9), this is represented by the Schur complement of the virtual block (vv) given by (B18) Note that the CM on the LHS of this equation tracks every physical mode independently. Hence, its dimension is proportional to the system size and becomes large for our cases of interest. We will now specialize Eq. (B18) to translation invariant states, where we can achieve a massive reduction of the size of the physical CM down to the number of Bloch bands independent of the system size.

Appendix C: Translation invariant GfTNSs
In this appendix, we specialize the formalism reviewed in Appendix B to translation invariant GfTNSs, where the contraction of the network can often be performed analytically. In Appendix C 1, we compute the CM of a translation invariant GfTNS, which we then use in Appendix C 2 to construct a parent Hamiltonian for the state.

Covariance matrix
For translation-invariant Gaussian TNSs, we introduce the Fourier transform (FT) of the physical and virtual mode operators as for all τ = 1, . . . , f , j = 1, . . . , ξ, α = L, R for MPSs and α = L, U, R, D for two-dimensional PEPSs. Here, the FT in two spatial dimensions with N x and N y sites in the horizontal and vertical direction, respectively, position vector x = (x, y) and momentum vector k = (k x , k y ), is given by The momenta in the horizontal (k x ) and vertical (k y ) direction take values k x = 2πj Nx with 0 ≤ j ≤ N x − 1 and k y = 2πj Ny with 0 ≤ j ≤ N y − 1. The FT for a single spatial direction is analogous.
We define the FT of the mode vector χ x from Eq. (B7) to be Therefore, k , a 2,k , . . . , a f,k , b L,1,k , b L,2,k , . . . , b R,1,k , . . . , mixes mode operators at momenta k and −k similarly to a Nambu spinor. Due to the translation invariance of the GfTNS, the FT brings the CMs of the physical state, the total fiducial state and the total virtual maximally entangled state into a block diagonal form. We denote the Majorana CM of the total fiducial state w.r.t. the Fourier transform of the mode operators by where the last equality defines the Majorana CMΓ Q (k) restricted to the block of momentum k. Note that the size ofΓ Q (k) is given by twice the number of Bloch bands, 2×n modes , and therefore no longer grows with the system size. Analogous statements hold for the physical CM Γ |ψ with the Fourier block matrixΓ |ψ (k) of size 2f , and for the CM of the total virtual maximally entangled state Γ ω with Fourier blockΓ ω (k) of size 2(n modes − f ).
Since the TNS is translation invariant, the local fiducial state |Q x is the same on every unit cell. Hence, the Fourier CM of the total fiducial state is localised at momentum k = 0,Γ Here, Γ Qx is the CM of the local fiducial state |Q x on a single site of dimension 2n modes . It has the block structure where the real anti-symmetric blocks A and D of dimension 2f and 2(n modes − f ) describe the physical and virtual subspaces, respectively, whereas the block B encodes the coupling between physical and virtual modes.
On the other hand, the CM of the total virtual maximally entangled stateΓ ω (k) has a non-trivial momentum dependence since the maximally entangled states from Eq. (B13) connect different unit cells. It is a direct sum of the contributions from the different spatial directions.
For a single virtual fermion with ξ = 1, the Majorana Fourier CM of the horizontal bonds oriented from left to right is the Fourier transform of the matrix Γ from Eq. (B16). It reads [40] Γ This matrix is written in the basis of the FT of the Majorana operators constructed from the complex modes where we have omitted the index j since ξ = 1. The contribution from the vertical bonds is given by Eq. (C7) with k x → −k y , since the vertical bonds are oriented downwards and hence in the direction of negative k y . Therefore, the Fourier CM of the virtual bonds is anti-hermitian and satisfies the identi- We are now in a position to compute the physical Majorana CM in momentum space by taking the FT of Eq. (B18). In momentum space, the RHS simplifies significantly due to the expression for the fiducial state CM given in Eqs. (C5) and (C6). We thus find [33] The physical state is well-defined if the matrix inversion can be carried out, i.e. if the determinant is not equal to zero. We emphasize that the matrices in Eq. (C8) have a constant size given by the number of Bloch bands, such that Eq. (C8) can typically be evaluated analytically. In Appendix E, we use Eq. (C8) to evaluate the Bloch CM of the SSH model MPS from Sec. II A, and in Appendix H, we use Eq. (C8) to evaluate the Bloch CM of the quadrupole model PEPS of Sec. IV B.

Parent Hamiltonian
A translation invariant Gaussian fermionic TNS |ψ has an infinite number of parent Hamiltonians for which it is an exact ground state. For any non-negative scalar function (k) ≥ 0 on the Brillouin zone, is a parent Hamiltonian for the Gaussian fermionic TNS |ψ [33,40]. Here, The properties of the parent Hamiltonian H depend on the dispersion function (k). If (k) > 0 is strictly positive throughout the Brillouin zone, H is gapped. Moreover, if all matrix entries of the product (k)Γ |ψ (k) are polynomials in e ±ikx and e ±iky , the parent Hamiltonian is strictly local.
A natural, but not unique, choice for the dispersion function is given by (k) = q(k) from Eq. (C9). Indeed, since D andΓ ω (k) are anti-hermitian and of even dimension (see Eqs. (C6) and (C7)), it follows that q(k) is real. If moreover q(k) is strictly positive throughout the Brillouin zone, implying that the PEPS has exponentially decaying real-space correlations, the parent Hamiltonian H q is gapped and strictly local with all terms acting on at most 2ξ successive unit cells [40,47].
Appendix D: GfTNSs with conserved particle number In the main text, we consider Gaussian fermionic TNSs with a conserved particle number: the ground states of both the SSH model and the quadrupole model lie at half filling. The TNSs are written in a basis related to the physical basis by a staggered particle-hole conjugation (cf. Eq. (7) for the SSH MPS and Eq. (28) for the quadrupole PEPS). Hence, the U(1) symmetry of the local tensors, which imposes the conservation of the physical particle number, also takes a staggered form (cf. Eq. (12) for the SSH MPS and Eq. (30) for the quadrupole PEPS). In this appendix, we rephrase this U(1) symmetry in the language of fiducial states, and show that it enforces many vanishing elements for the CM of the local fiducial state. These are equivalent to the vanishing of the off-diagonal block Q = 0 of the complex CM of a state with conserved particle number (see Appendix A 1), but expressed in the basis after the staggered particle-hole transformation. We focus on the one-dimensional case since the computation for two-dimensional GfTNSs is analogous.
For one-dimensional Gaussian fermionic MPS, we consider a U(1) symmetry of the local tensor of the general form where η τ,x , η L,j,x , η R,j,x ∈ {±1} for τ = 1, . . . , f and j = 1, . . . , ξ. Here, whose many-body basis representation is given in Eq. (D1). Here, each individual operatorÛ in the product acts on exactly one fermion. For instance, the operator acting on the physical fermion τ is given bŷ and similarly for the virtual fermions. This agrees with the matrix representation of the U(1) rotation acting on a single mode given in Eq. (D2). We observe that the U(1) operator of the physical fermion τ from Eq. (D4) satisfies Extending this to all of the modes in the mode operator of Eq. (B7), we find that where 1 ≤ µ, ν ≤ 2n modes . Here, we collected all the parameters η into a vector of length 2n modes using the same ordering as in the mode vector of Eq. (B7), η x = η 1,x , η 2,x , . . . , η f,x , η L,1,x , η L,2,x , . . . , η R,1,x , . . .
Therefore, they vanish unless (η x ) µ = (η x ) ν . This shows that the symmetry of Eq. (D1) forces the vanishing of half the elements of the complex CM of the local fiducial state, where 1 ≤ µ, ν ≤ 2n modes .

Appendix E: Covariance matrix for SSH MPS
In this appendix, we illustrate the formalism of GfTNS introduced in Appendix B and C using the example of the MPS from Eq. (11) describing charge pumping in the SSH model. After expressing the state as a GfTNS in Appendix E 1, we demonstrate the computation of its Bloch CM and parent Hamiltonian in Appendix E 2.

Expression as GfTNS
Our goal is to express the SSH charge pumping MPS, which is already fully defined by Eq. (11) in the language of local tensors, as a GfTNS using the formalism from Appendix B. Since the MPS has physical dimension 2 and bond dimension 2, it is described by f = 1 physical fermion and ξ = 1 virtual fermion per nearestneighbour bond and lattice site (cf. the first paragraph of Appendix B 2). On the A sublattice, the annihilation operators for the physical, left virtual and right virtual modes are a A,x , b L,A,x and b R,A,x , respectively, and similarly for the B sublattice. Following the recipe given above, we need to find the local fiducial states |Q A,x = Q A,x |Ω and |Q B,x = Q B,x |Ω that match the local tensors A and B from Eq. (11) at each unit cell x. The link between the local tensors and fiducial states is then given by Eq. (B10) stating that the local tensor entries are the basis coefficients of the fiducial state. For example, the tensor element A 1 01 = β tells us that the local fiducial map Q A,x contains a term βa † A,x b † R,A,x . Performing this matching for every non-zero tensor entry, we find that the fiducial maps are given by In a second step, we want to write the local fiducial states as Gaussian states satisfying Wick's theorem in order to express the MPS as a GfTNS. Since the local tensors are parity even, the local fiducial states |Q A,x and |Q B,x can be parametrised as in Eq. (B12) using the exponential of anti-symmetric coefficient matrices M A,x and M B,x . In the case of the SSH pumping MPS, this is very simple, since the fiducial states we derived in Eq. (E1) have the form of Eq. (A11), allowing us to directly read off M A,x and M B,x (imposing antisymmetry).
We find that the coefficient matrices are given by Here, we defined the quotients a = α/γ and b = β/γ of the parameters used in Eq. (11), and we absorbed the remaining factor γ into the normalisation constant N in Eq. (A11). The complex CMs of these local fiducial states are computed from M A,x and M B,x using Eq. (A12), and are then transformed to the Majorana representation using Eq. (A3). One finds that the Majorana CM Γ Q A,x is (E3) where we introduced the shorthand notation c ≡ 1 + 4a 2 + 4b 2 , and the basis is the Majorana basis derived from (a A,

Bloch CM and parent Hamiltonian
Now that we have expressed the MPS from Eq. (11) as a GfTNS, we want to use the power of the formalism introduced in Appendix C to derive its CM and parent Hamiltonian on a chain with periodic boundary conditions. We will compute the Bloch CM by evaluating Eq. (C8), and from there obtain the parent Hamiltonian via Eq. (C10).

a. CM for a unit cell
The expression for the Bloch CM from Eq. (C8) is valid for a translation-invariant GfTNS. In particular, the CM Γ Qx with its blocks A, B and D from Eq. (C5) is the Majorana CM of the fiducial state of a unit cell, not a single site. In order to proceed, we therefore need to derive the fiducial state of an entire unit cell by contracting the virtual bond within a unit cell. In the language of GfTNSs, this is described by projecting the fiducial states Q A,x Q B,x |Ω of the two sites in one unit cell onto the virtual maximally entangled state connecting them, This projection is a special case of Eq. (A7). Hence, the CM Γ Qx of the resulting state is given, as in Eq. (A9), by the Schur complement of the of the block of size 4×4 corresponding to the Majorana mode operators constructed Here, A and D are the blocks on the diagonal of the direct sum Γ Q A,x ⊕ Γ Q B,x corresponding to the Majorana modes derived from Γ Qx has the block form where A and D are real anti-symmetric blocks of size 4 × 4, and B is a real block of size 4 × 4. We find that the diagonal blocks A and D are of the form where r and s are parameters. Specifically, the physical and virtual blocks are Here, the denominator is a consequence of the matrix inverse in Eq. (E5). In addition, the block containing the coupling between physical and virtual fermions is We can now directly compute the Fourier Majorana CMΓ |ψ (k x ) of the physical state defined by the MPS on a chain with N x sites and periodic boundary conditions. Γ |ψ (k x ) is given by the Schur complement in Eq. (C8), where the Fourier Majorana CMΓ ω (k x ) for the virtual bonds is given by Eq. (C7), and A, B and D are given in the previous section.
For the SSH pumping MPS, Eq. (C8) is an matrix equation of size 4 × 4 since there are two physical Bloch bands. The matrix inverse can be evaluated analytically using the special parametrisation Z (1) (r, s) from Eq. (E7). Indeed, the Fourier Majorana CM for the virtual bonds from Eq. (C7) can be written asΓ ω (k x ) = Z (1) (0, e ikx ). One easily checks that Using these identities, the matrix inverse in Eq. (C8) can be performed by hand, and we evaluate the determinant q(k x ) from Eq. (C9) as Unless γ = 0 and |α| = |β|, q(k x ) is strictly positive such that the MPS is well-defined everywhere except for these parameter values.
After the matrix inversion in Eq. (C8), the remaining matrix multiplications in Eq. (C8) can be performed using a computer algebra system. We thus find that the Fourier Majorana CMΓ |ψ (k x ) of the physical state is again of the form of Eq. (E7), with parameters c. Parent Hamiltonian Using Eq. (C10), we can now find a Bloch parent Hamiltonian H for the MPS from its Majorana Fourier CM which we computed in Eq. (E13). In order to gain a physical understanding of the parent Hamiltonian, we express it in terms of the original complex physical modes before the particle-hole transformation of Eq. (7) given byâ after the FT of Eq. (C1). We find where the functions s(k x ) and r(k x ) are defined in Eq. (E13).
In order for H to be gapped and strictly local, we need to find a dispersion function (k x ) which is strictly positive such that (k x )s(k x ) and (k x )r(k x ) are polynomials in e ±ikx . As explained in Appendix C 2, a natural choice is (k x ) = q(k x ). Indeed, q(k x ) computed in Eq. (E12) cancels the denominator of Eq. (E13). Since q(k x )s(k x ) and q(k x )r(k x ) contain the factor e ±ikx up to second order, the parent Hamiltonian H q contains hopping terms between up to second-nearest neighbour unit cells.
Due to the special structure of Eq. (E8), we can in fact obtain a more short-ranged parent Hamiltonian from the choice proportional to q(k x ), which is in turn proportional to the denominator of r and s from Eq. (E13b). (k x ) is strictly positive for all parameter values that lead to a well-defined state. The factor 1/(a 4 + b 4 + 1) is a normalisation ensuring that the parent Hamiltonian matches Eq. (3) if the MPS is given by the parametrisation φ pump of Eq. (15). We see that the parent Hamiltonian H with Bloch representation w.r.t. the original complex physical modes has hopping only up to nearest-neighbour unit cells.
Appendix F: Column covariance matrix of real-space (d + 1)-dimensional TNS from charge pumping of d-dimensional TNS with conserved particle number In Sec. III C, we introduced the tensors A col and B col describing the real-space Chern PEPS restricted to a column of sites on the A and B sublattice at positions {(x, y)} 0≤y≤Ny−1 . A col and B col are defined by the application of the inverse FT of Eq. (22) to the physical and horizontal virtual legs of a column of the hybrid Chern PEPS at positions {(x, k (j) y )} 0≤j≤Ny−1 . In this appendix, using the representation of the SSH pumping MPS as a GfTNS from Appendix E, we compute A col and B col explicitly in terms of their CMs. In order to demonstrate the generality of the result, we will consider a general (d + 1)-dimensional TNS constructed from charge pumping of a d-dimensional TNS, which is assumed to have a conserved number of physical particles, and hence possess a U(1) symmetry of the form discussed in Appendix D.
Thus, let |ψ d (t) be the Gaussian fermionic TNS in d spatial dimensions along a cyclic interpolation parametrised by the time t ∈ (−π, π]. With the same basis as in Eq. (B7), we collect the physical and virtual mode operators for one unit cell x ∈ Z d of |ψ d (t) into the mode vector of length 2n modes . The physical and virtual mode operators now depend on the time t along the interpolation. As explained in Appendix B, the d-dimensional TNS is defined by its Gaussian local fiducial state Q x (t)|Ω , which is characterised by its complex CM G x (t) of dimension 2n modes . We assume that the TNS has a conserved number of particles, such that the local fiducial state Q x (t)|Ω of |ψ d (t) has a U(1) symmetry of the form discussed in Appendix D. This symmetry determines which physical and virtual modes correspond to holes and particles: for each 1 ≤ µ ≤ n modes , an entry (η x ) µ = 1 or (η x ) µ = −1 in the vector η x from Eq. (D6) indicates that the mode µ has a particle-or hole-like character, respectively. Note that η x does not depend on the time t, such that the holeor particle-like character of the modes remains unchanged along the interpolation.
We can now move to the hybrid (d + 1)-dimensional TNS, which is defined by Eq. (20) of the main text. From Sec. III B we recall that the local fiducial state (equivalent to the local tensor) of the hybrid state at the position (x, k (j) d+1 ) is given by Q x (t (j) )|Ω containing the modes χ x (t (j) ) from Eq. (F1). In particular, due to the tensor product in the (d + 1) st direction, virtual fermions in this direction are not needed.
We can now easily write the CM of one column of the hybrid (d + 1)-dimensional TNS, given by the sites at positions {(x, k (j) d+1 )} 0≤j≤N d+1 −1 . Indeed, due to the absence of virtual legs in the direction d + 1, the contraction of the bonds in this direction of the hybrid column amounts to a trivial tensor product in the language of local tensors. In terms of fiducial states, this corresponds to a direct sum of CMs. Hence, the complex CM G hybrid Here, G hybrid x is written in the basis {(χ x ) µ (t (j) )} 1≤µ≤2n modes ,0≤j≤N d+1 −1 of all operators for the physical and virtual modes in the first d directions in the column.
We now consider the (d + 1)-dimensional real-space state restricted to a column {(x, x d+1 )} 0≤x d+1 ≤N d+1 −1 , which is obtained by applying the inverse FTF in direction d + 1 to the physical and virtual legs of the hybrid column. The complex CM G col x of this at position x is therefore given by This expression can be further simplified due to the constraint of Eq. (D8) imposed on G x (t (j) ) by its U(1) symmetry, implying G x (t (j) ) µ,µ = 0 unless (η x ) µ = (η x ) µ . Indeed, we may thus define the matrix that mixes elements of the complex CM of the ddimensional state at times t (j) and −t (j) according to whether the modes µ and µ transform as particles or holes, respectively. The 2n modes -dimensional blocks of the column CM G col x are then given by the FT of this matrix, This expression is explicitly invariant under real-space translations x d+1 → x d+1 + 1 acting on both the physical modes and the virtual modes in the first d directions.
Hence, the inverse FT of Eq. (F4) guarantees the translation invariance of the column CM G col x in the direction d + 1.
Finally, we want to investigate which form the global U(1) symmetry related to particle number conservation takes for G col x . This result will be used in Appendix G. The U(1) is inherited from the invariance of each fiducial state Q x (t)|Ω under the U(1) symmetry of Eq. (D3), which is independent of t. From the expression of the latter in second quantization (cf. Eq. (D4)), we see that the generator of the global U(1) symmetry of the column state is whose expectation value vanishes.
Appendix G: Disentangled modes in single-particle ES at fixed particle number In this appendix, we show that as stated in Eq. (25), the maximal number of entangled modes compatible with the U(1) symmetry of the SSH model MPS in the ES of the Chern PEPS column state is given by max{3L, N y }, where L is the number of sites in the subsystem. We proceed by deriving a lower bound on the number of disentangled modes in the ES of a state with a conserved particle number (and therefore an upper bound on the number of entangled modes). We first consider a generic state in Sec. G 1, before specialising to the fiducial state of a real-space column of the (d + 1)-dimensional TNS in Sec. G 2.

Insulator with filling fraction q
We consider a non-interacting system of N fermionic DOFs with creation and annihilation operators a † j , a j for j = 1, . . . , N . We define a bi-partition of the system into the subsystem A and its complementĀ, where the DOFs of A are described by the first N A modes j = 1, . . . , N A . Let |ψ be a pure state of this system with a conserved particle number and filling fraction q, such that the total number of occupied modes in |ψ is qN . We denote by n λ=1 A the number of entanglement levels with the value λ = 1 in the single-particle ES of |ψ restricted to A. We now want to show that this number is bounded below by the filling fraction as Proof: Let H = N i,j=1 h ij a † i a j be a non-interacting flat-band Hamiltonian whose ground state is |ψ , where h is a Hermitian matrix of dimension N ×N . The occupied modes in |ψ are given by the qN orthogonal eigenstates u (k) of h with energy -1, i.e.
where k = 1, . . . , qN . W.r.t. the bi-partition into A,Ā, each basis state u (k) falls into exactly one of the following three categories (assuming that the eigenstates are ordered accordingly): (1) For k = 1, . . . , m 1 , the states satisfy u (k) i = 0 for i = N A +1, . . . , N such that the corresponding occupied modes are composed of DOFs of the subsystem A. According to Theorem 1 of Ref [68], each such state leads to an entanglement level λ = 1 in the singleparticle ES of |ψ restricted to A. Hence, m 1 ≤ n λ=1 A .
(2) For k = m 1 + 1, . . . , m 2 , u The numbers m 1 , m 2 are assumed to be maximal in the sense that no linear combination of eigenstates of the category (3) lies either purely in A or purely inĀ.
In the next paragraph, we will show that the number qN − m 1 of states from the categories (2) and (3) is no larger than the number N − N A of DOFs inĀ. This proves the claim, since the number of states from the different categories can then be estimated as leading to Eq. (G1). As a final step, we need to show that qN − m 1 ≤ N −N A . This follows from the linear independence of the with coefficients µ k . The scalar product withũ (l) shows that 0 = µ l for l = m 1 + 1, . . . , m 2 : Indeed, the orthogonality of the u (k) implies i>N Aũ where the sum runs only over category (3). Due to the maximality of m 1 , m 2 as defined above just below point (3), we then have 0 = qN k>m2 µ k u (k) . Otherwise, this would be a linear combination purely in A since the part inĀ vanishes by assumption. We thus also get 0 = µ l for l = m 2 + 1, . . . , qN , proving that qN − m 1 ≤ N − N A .
2. Real-space column of (d + 1)-dimensional TNS We now apply Eq. (G1) to the case where |ψ is the fiducial state of one real-space column of the (d + 1)dimensional pumping TNS, whose CM is computed in Eq. (F5). The column state has N = n modes N d+1 degrees of freedom, where N d+1 is the number of sites in the direction d+1 and n modes the number of physical and virtual particles per lattice site of the d-dimensional TNS. Below we will show that in a suitable basis, the column state has a conserved particle number qN = N d+1 n η− , where n η− is the number of values −1 in the first n modes entries of the vector η x from Eq. (D6). We consider the single-particle ES of the column state w.r.t. the subsystem A L of the first L sites 0 ≤ x d+1 ≤ L − 1 of the column, which has N A L = Ln modes DOFs. By Eq. (G1), the number n λ=1 A L of entanglement levels with the value λ = 1 in this spectrum is lower bounded as We now show that with a suitable particle-hole transformation we can find a single-particle basis w.r.t. which the column state has a conserved particle number N d+1 n η− . In Eq. (F5) the state is expressed in a basis where the CM G col x has a non-vanishing off-diagonal block Q corresponding to superconducting correlations, such that only the parity of the particle number is conserved (cf. Eq. (A1)). Below Eq. (F5) we note that half the entries of G col x vanish, namely those with (η x ) µ = (η x ) µ . A new basis, in which the complex CM has only a diagonal block R and hence a fixed particle number, is created as follows: For those modes 1 ≤ m ≤ n modes with (η x ) m = −1, we exchange annihilation and creation operators, since the latter have (η x ) m+n modes = 1. This is a particle-hole transformation corresponding to the mapping for all physical modes with η τ,x = −1 while leaving modes with η τ,x = 1 unchanged, and similarly for the virtual modes.
The number of particles in the column state is given by the expectation value from Eq. (F8) of the generator of its global U(1) symmetry. Note that this is not the physical particle number, but rather the number of particles in the system composed of the physical legs and the virtual legs in the first d spatial dimensions. Under the transformation of Eq. (G6), the particle number expectation value from Eq. (F8) transforms as explained in Eq. (13). We thus find that in the new basis of modes the particle number is N d+1 n η− as claimed above.

a. Chern PEPS
Let us apply these results to the Chern PEPS derived from the SSH pumping MPS. For simplicity, we restrict ourselves to a column A col of A sites with n modes = 3 DOFs per site and η A = 1 and η L,A = η R,A = −1 (cf. the discussion below Eq. (13)). Therefore, n η− = 2 and Eq. (G5) becomes This shows that Eq. (25) gives the maximal number of entangled modes compatible with the U(1) symmetry of the SSH model MPS. For the PEPS defined by the trivial cycle φ triv from Eq. (17), the discussion above can be refined: since β = 0 throughout the interpolation, all right virtual modes decouple from the column tensor such that there trivially are N y entanglement levels |λ| = 1. Let us investigate if the U(1) symmetry causes additional decoupled levels in the system of the coupled physical and left virtual particles. This system has only two DOFs per site, namely the physical leg with η A = 1 and the left virtual leg with η L,A = −1. In this case, Eq. (G5) gives a trivial lower bound for the number of decoupled modes with levels λ = 1, n λ=1 A L ≥ max{2L − N y , 0} = 0 for L ≤ N y /2. (G8) Hence, the U(1) symmetry does not cause any additional decoupled modes, and the number of entangled modes is given by 2L as discussed in the main text.

b. Chiral hinge PEPS
The discussion for the three-dimensional chiral hinge PEPS from Sec. IV C is analogous, where a column of sites on the sublattice 1 has n modes = 5 DOFs per site with η 1 = 1 and η L,1 = η U,1 = η R,1 = η D,1 = −1. Therefore, n η− = 4 and Eq. (G5) implies n λ=1 For the mirror-symmetric case with couplings α = α x = α y and β = β x = β y , this bound can be refined. Indeed, by defining the linear combinations b LD ±,1 = (b L,1 ± b D,1 )/ √ 2 and b U R ±,1 = (b R,1 ± b R,1 )/ √ 2, the local fiducial state |Q [1] from Eq. (H1) below can be written as Therefore, when considering only one column of sites on the sublattice 1, two virtual fermionic modes decouple from the local tensor on each site. Effectively, the remaining coupled system therefore has only n modes = 3 DOFs per site with η 1 = 1 and η LD −,1 = η U R −,1 = −1, such that the number of additional disentangled modes due to the U(1) symmetry can be estimated as above for the Chern PEPS. Hence, with the identification N y → N z , the bound on the number of disentangled modes is given by Eqs. (G7) and (G8) for the PEPSs derived from φ pump and φ triv , respectively.

Appendix H: Quadrupole PEPS as GfTNS
In this appendix we apply the formalism of GfTNSs to the quadrupole model pumping PEPS of Eq. (29). In Appendix H 1, we show how the state can be expressed as a GfTNS. This allows us to compute its Bloch CM and a Bloch parent Hamiltonian on the torus in Appendix H 2. These sections are completely analogous to Appendix E for the SSH pumping MPS. We therefore refer the reader to this appendix for a detailed explanation of each step. Finally, in Appendix H 3 we discuss the ES of the PEPS when the parameters are chosen such that the state represents the OAI dimerized phase of the quadrupole model.

Expression as GfTNS
The quadrupole model pumping PEPS is defined in Eq. (29) of the main text in the language of local tensors. Here, we want to re-express this state in the formalism of GfTNSs. We recall that each unit cell consists of 2 × 2 lattice sites, and that the PEPS has physical dimension 2 and bond dimension 2. This corresponds to f = 1 physical fermion per lattice site and ξ = 1 virtual fermion per nearest-neighbour bond and lattice site, represented in Fig. 6(b) by blue and red circles, respectively.
As explained in Appendix B, the local PEPS tensors A [τ ] on the four sublattices τ = 1, 2, 3, 4 from the main text correspond to local fiducial states |Q [τ ] , whose basis coefficients are given by the local tensors (see Eq. (B10)). We write a τ for the annihilation operator of the physical fermion and b L,τ , b U,τ , b R,τ , b D,τ for the annihilation operators of the left, up, right, down virtual fermions on the sublattice τ (we dropped the unit cell index due to the translation invariance). Let |Ω denote the vacuum annihilated by all these operators. Applying Eq. (B10), we see that the local fiducial states on the four sublattices derived from the local tensors of Eq. (29) are These fiducial states are of the form of Eq. (A11) with only zero-and second-order terms in the creation operators. Hence, the fiducial states are Gaussian and can be parametrised as in Eqs. (A10) and (B12) with antisymmetric coefficient matrices (H2d) Here, we defined the quotients a x = α x /γ, a y = α y /γ, b x = β x /γ and b y = β y /γ of the parameters from Eq. (29). We absorbed the remaining factor γ into the normalisation constant N in Eq. (A11).
We have thus successfully written the PEPS from Eq. (29) as a GfTNS. Before proceeding, we choose the orientation of the virtual bonds as follows: All horizontal bonds are oriented from left to right and all vertical bonds are oriented from top to bottom. The order of the individual fermions within these blocks is indicated by the labels next to each circle.

Bloch CM and parent Hamiltonian
Having expressed the quadrupole model pumping PEPS as a GfTNS, we now want to use the formalism from Appendix C to compute the Bloch CM and a Bloch parent Hamiltonian for the state on a torus. We will compute the Bloch CM by evaluating Eq. (C8), and from there obtain the parent Hamiltonian via Eq. (C10).

a. CM of unit cell
The Majorana CM Γ Qx in Eq. (C8) refers to the fiducial state of a unit cell, not that of a single site. We therefore need to compute Γ Qx from the fiducial states for each individual sublattice given in Eq. (H1). This is done by contracting the four virtual bonds within one unit cell. In the language of fiducial states, we project the tensor product |Q [1] ⊗ |Q [2] ⊗ |Q [3] ⊗ |Q [4] of the fiducial states for a unit cell on the product of the maximally entangled states of the four virtual bonds connecting the lattice sites within this unit cell. This is analogous to the computation leading to Eq. (E13), and we refer the reader to Appendix E 2 a for details. Here we will only state the result for the Majorana CM Γ Qx .
The fiducial state of one unit cell, after contraction of the virtual bonds within the plaquette, describes four physical fermions and two virtual fermions per direction left, up right, down (see Fig. 10). We collect the four physical, horizontal virtual and vertical virtual fermions Eq. (C9) as q(k) = 2 8 1 + a 4 + b 4 + a 2 b 2 cos k x + a 2 b 2 cos k y 4 1 + 2 √ 2a 2 + 2a 4 + b 4 4 .

(H12)
This is strictly positive unless γ = 0 and |α| = |β|, such that the PEPS is well-defined everywhere except for these parameter values.

c. Parent Hamiltonian
We are now in a position to find a parent Hamiltonian for the physical state on a torus, using Eq. (C10) and the result for the Bloch Majorana CMΓ |ψ from Eq. (H13). To get an intuitive understanding of the Hamiltonian, we find it useful to express it in terms of the original fermionic modes before the particle-hole transformation of Eq. (28). Their FT is related to the FT of the new modes asâ τ,k = a τ,k for τ = 1, 2 andâ τ,k = a † τ,−k for τ = 3, 4. With respect to the original fermionic modes, Eq. (C10) is given by where a summation over τ, τ is implied in the second line. As discussed in Appendix H 2 c, the properties of the parent Hamiltonian are determined by our choice of dispersion relation (k). By setting (k) = q(k), we would obtain a gapped parent Hamiltonian with coupling between up to fourth-nearest neighbour unit cells. However, we can find a more short-ranged parent Hamiltonian due to the special form of the fiducial state CM using Eq. (H5). Indeed, by choosing the dispersion function (k) = 1 + a 4 + b 4 + a 2 b 2 cos k x + a 2 b 2 cos k y a 4 + b 4 + 1 (H15) which is proportional to the denominator of r, s and u from Eq. (H13), we obtain a second-nearest neighbour parent Hamiltonian with Bloch representation where σ 0 is the identity matrix of dimension two.

ES in topological quadrupole phase
When b = 0, the system described by the PEPS of Eq. (H2) splits into decoupled four-site plaquettes shifted from the unit cell by one site in both directions, corresponding to the OAI dimerized phase of the quadrupole model if a = 1. The Majorana CM Γ Plaquette describing the physical state of one such decoupled plaquette takes the form of Eq. (H5) with As an application, we will derive the ES contributions from the edges and corners given in Eqs. (34a) and (34b), respectively. We begin with the four corners. The Majorana CM Γ corner,τ for a single corner site on the sublattice τ is given by the corresponding block of dimension 2 on the diagonal of Γ Plaquette . Specifically, Γ corner,τ = 0 −r r 0 (H18) for τ = 1, 2, 3, 4. We now transform Γ corner,τ to the basis of the original complex fermionic modes before the particle-hole transformation of Eq. (28). Then, the CM of the corner site has a vanishing off-diagonal block Q * corner,τ = 0 and a diagonal blockR * corner,τ = ∓ i 2 λ corner with λ corner = r, where the negative and positive signs hold for the sublattices τ = 1, 2 and τ = 3, 4, respectively. Using the expression for r from Eq. (H17a) and a = α/γ, we obtain the formula for the corner ES level given in Eq. (34b) with one level per corner.
This corresponds to the formula given in Eq. (34a).