Bulk Localised Transport States in Infinite and Finite Quasicrystals via Magnetic Aperiodicity

Robust edge transport can occur when particles in crystalline lattices interact with an external magnetic field. This system is well described by Bloch's theorem, with the spectrum being composed of bands of bulk states and in-gap edge states. When the confining lattice geometry is altered to be quasicrystaline, then Bloch's theorem breaks down. However, we still expect to observe the basic characteristics of bulk states and current carrying edge states. Here, we show that for quasicrystals in magnetic fields, there is also a third option; the bulk localised transport states. These states share the in-gap nature of the well-known edge states and can support transport along them, but they are fully contained within the bulk of the system, with no support along the edge. We consider both finite and infinite systems, using rigorous error controlled computational techniques that are not prone to finite-size effects. The bulk localised transport states are preserved for infinite systems, in stark contrast to the normal edge states. This allows for transport to be observed in infinite systems, without any perturbations, defects, or boundaries being introduced. We confirm the in-gap topological nature of the bulk localised transport states for finite and infinite systems by computing common topological measures; namely the Bott index and local Chern marker. The bulk localised transport states form due to a magnetic aperiodicity arising from the interplay of length scales between the magnetic field and quasiperiodic lattice. Bulk localised transport could have interesting applications similar to those of the edge states on the boundary, but that could now take advantage of the larger bulk of the lattice. The infinite size techniques introduced here, especially the calculation of topological measures, could also be widely applied to other crystalline, quasicrystalline, and disordered models.


A. Motivation
In crystalline materials, e.g. condensed matter or cold atoms in optical lattices, the standard picture according to band theory is that a system is either an insulator or metal [1,2]. During the 1980s, this picture began to change with the discovery of topological states of matter [3][4][5]. For example, topological Edge States (ESs) can occur when a charged particle in a crystal interacts with an external magnetic field [6,7]. The physics of charged particles in a two-dimensional crystalline lattice with an applied strong magnetic field is a well-studied problem for both the single-particle [8][9][10][11][12][13][14] and many-body [15][16][17][18] regimes. There have also been numerous experimental realisations and proposals [7,[19][20][21][22].
Recently, there has been renewed interest in adding a magnetic field to scenarios involving quasicrystalline lattices [37][38][39][40][41][42][43]. In quasicrystals, the concepts of bands and band-gaps are difficult to consistently define, since Bloch's theorem is not enforceable without approximations to the overall structure. While recent results have confirmed the presence of ESs in a magnetic field [38,43] and studied the appearance of higher-order topological states [44][45][46] in quasicrystals, there have been few tangible differences from their study in periodic systems.
In this paper, we will show that the now standard picture of insulators, metals and topological insulators with surface states is not the full story for quasiperiodic systems. When the confining potential is quasicrystalline Figure 1. Illustration of the different types of in-gap states in (a) finite and (b) infinite quasicrystals when a uniform magnetic field (blue arrow) is present. In the case of a finite quasicrystal with N sites (a), the green arrow depicts potential transport across a conventional ES, which forms at the boundaries of the lattice. The quasicrystal may also permit the formation of BLT, whose potential transport is depicted by the red arrow. In the case of a truly infinite quasicrystal (b), we no longer have boundaries and hence no ES, but we retain the BLT state and it's supported transport.
there are two competing not necessarily commensurate length scales from the magnetic field and quasiperiodic lattice. This results in a magnetic aperiodicity which directly leads to the observation of Bulk Localised Transport (BLT) states. As BLT states arise from the magnetic aperiodicity, they are significantly different from previous states found in the internal sections of fractal lattices [47][48][49]. These fractal lattices are almost entirely composed of effective hard edges with no discernible bulk, meaning these 'internal' ESs are truly just conventional ESs on an unconventional lattice. The BLT states are not an artefact of effective edges introduced through an impurity or set of dislocations [13,[50][51][52][53][54] from the quasicrystalline lattice. We will also note that the BLT states appear to be of a different character than the higher-order zerodimensional corner modes found in quasicrystals [44], as these are entirely bound to corners on the hard boundary of a finite system. The BLT states do share many of the properties of ESs, but are entirely localised within the bulk, as illustrated in Fig. 1. Importantly, the BLT states support transport (or currents) in the same way as the well-studied ESs, which could make them useful for future applications utilising BLT.

B. Overview
Before detailing the methods and calculations, we overview the main results of our work. This section is intended to review our work, with references to the rest of this paper where we discuss the approach and results in detail.
Quasicrystals and magnetic fields. We envisage a system of charged particles existing on a quasicrystalline lattice under the influence of a uniform perpendicular magnetic field. This system is illustrated in Fig. 1. The Hamiltonian of this system is well described by the Hofstadter vertex model, introduced in Sec. II, which modifies the standard Hofstadter model to the vertex model of a quasicrystalline tiling, illustrated in Fig. 3. Solving for the finite size spectrum and states of the Hofstadter vertex model is relatively straight forward. However, defining bands and in-gap states is difficult due to the breakdown of Bloch's theorem. We define the in-gap states, and hence the bands, via topological measures that are non-zero for in-gap states and outlined in Sec. III.
Infinite size quasicrystals. As the states we will consider are supported by the bulk of the lattice, it is interesting to consider if they are retained in the spectrum of the infinite Hamiltonian. To consider the infinite-size quasicrystalline lattice, we will utilise an infinite-size error controlled algorithm; which abandons square truncation to allow for the interaction of sites outside a finite patch of the lattice to be accounted for. The infinite-size algorithm is described in Sec. II B and can be utilised to calculate the spectrum of other infinite-dimensional operators [55]. We also detail a new extension of the infinitesize algorithm in order to calculate topological measures, as described in Sec. III A. This allows for states with ingap characteristics in the infinite-size quasicrystal (or any crystal or aperiodic lattice) to be identified.
BLT states. We first illustrate BLT states in Sec. IV by considering in detail the Hofstadter vertex model of the quasicrystalline Ammann-Beenker (AB) tiling. As already stated, the BLT states are peculiar as they are ingap, but entirely localised within the bulk of the lattice, with no component on the edge. Examples of BLT states are given for the AB tiling in Fig. 5(c) and 8, clearly illustrating the bulk nature of these in-gap states. The BLT states are also shown to be in-gap from their nonzero topological measures. Interestingly, the proportion of states that are of BLT type in a finite system converges to a non-zero value with increasing system size, as shown in Fig. 7. This is contrary to the regular ESs; which become a vanishingly small proportion of the states with increasing system size. As expected from the arguments illustrated in Fig. 1, we find that BLT states exist in the infinite-size spectrum, and are of the same character as those in the finite system, see Fig. 10 for examples. We confirm that BLT states exist in the spectra of other quasicrystals for both finite and infinite cases in Sec. V. For example, we consider 5-fold, 7-fold, 10-fold, and 12fold examples in Fig. 16. We also show examples of BLT states in Fig. 17 where the rotational symmetry is broken and severe defects are included. In summary, we find that BLT states are ubiquitous in quasicrystals in magnetic fields, even in the presence of extreme deformations to the lattice.
Supported bulk transport. Key to the future consideration and application of BLT states is the fact that they support transport like any other in-gap state. This results in BLT which is independent of the existence or form of any edges. In this way, BLT can be considered to be robust against perturbations along the edge, much like regular edge states are robust against perturbations in the bulk. We show explicitly that transport is supported by the BLT states of the AB Hofstadter vertex model in both the finite and infinite size cases in Figs. 18 and 19 respectively. We also show that the location of this transport can be varied due to the BLT states being supported on different parts of the bulk. BLT will also exist in cases of extreme deformations to the system, as BLT states are retained. We also show in Fig. 20 that the BLT states of other quasicrystals also support transport along them as would be expected. The BLT discussed and exhibited in this paper is not in direct correspondence with the regular transport of an electron in a magnetic field. On a lattice, the cyclic motion of the electron is truly in correspondence to the cyclic motion around a single tile of the lattice. BLT is a larger collective effect of the quasicrystalline lattice and magnetic field which goes beyond the circling of a single tile but it is still, of course, cyclic in nature.
Origin of BLT states. We find through example toy models in Sec. VII that the BLT states are a direct result of the interplay of the length scales between the magnetic field and quasicrystal. We coin the term magnetic aperiodicity in relation to these competing length scales. The origin of this interplay can be motivated by considering a particle looping around the individual distinct tiles of the quasicrystalline lattice, as illustrated in Fig. 2. Each distinct tile has a different area, which is incommensurate in relation to the areas of the other distinct tiles. When a particle circulates one of these tiles, it gains a phase that is dependent on the flux, which in turn is dependent on the area of the tile, as illustrated in Fig. 2. This means when a particle circulates around each distinct tile it gains a phase that is incommensurate with the phase associated with circulating the other distinct tiles, as the areas are incommensurate to each other and the magnetic field is uniform. It is the convolution (or interference) of these incommensurate phases that is the interplay of the length scales of the magnetic field and quasicrystalline lattice. With the end result being the Figure 2. Origin of magnetic aperiodicity, showing two rhombic prototiles of a quasicrystalline tiling. In (a,c), a particle (denoted by the purple circle) is initially localised to one corner of the prototile. Each tile has a unique area of A1 and A2, and their fraction A1/A2 is irrational. After some time T in (b,d), the particle circulates around the prototile and acquires a unique phase of either θ1 or θ2. The fraction θ1/θ2 will again be irrational due to the presence of incommensurate areas, leading to a magnetic aperiodicity on the full quasicrystalline tiling.
generation of BLT states and support of BLT.

C. Terminology
As this work straddles the fields of condensed matter, quantum simulators, and spectral computations from mathematical physics, there are certain terminologies we need to ensure are defined consistently to avoid confusion between readers of different fields.
First, we will refer to a lattice as being periodic or crystalline as long as a unit cell can be defined with an associated Brillouin zone. This allows for the use of Bloch's theorem to calculate the band structure of the lattice. Note, this does not exclude the presence of edge states, as they are dependent on the boundary. Even in the case of open boundaries, bands can still be calculated. Therefore, a periodic lattice is defined independent of the boundary conditions, which would usually be considered as periodic, infinite or open. The strict definition of a lattice itself can be considered to only apply to periodic systems. However, the definition of a lattice in condensed matter physics is more general and can be interpreted as defining a group of discrete connected points. Throughout this work, we will follow this convention and refer to lattices as being any group of discrete connected points.
We will often refer to bulk states and in-gap states. Bulk states are all states that are allowed in the system with real quasimomentum according to Bloch's theorem. In-gap states are all other solutions to the Schrödinger equation for the lattice Hamiltonian, which are the complex quasimomentum solutions [13]. In general, we do not have access to the quasimomentum from numerical Figure 3. Construction of the Ammann-Beenker (AB) tiling, using the (a) incommensurate square and rhombus as prototiles. The aperiodic tiling is generated from these prototiles, leaving no gaps. Here, we take a circular cutoff in tiling space to show a (b) finite sample of AB tiling and preserve rotational symmetry with respect to the origin (centre of the tiling). Note, the enforcement of rotational symmetry here is arbitrary and plays no role in the formation of BLT in this study, as we show in Sec. V. The corresponding vertex model (c) is then defined by setting bonds as the edges of tiles, and the lattice sites as the intersection of tile edges. In this example, the total number of lattice sites is N = 185.
approaches. As Bloch's theorem breaks down in quasicrystals, we must turn to alternative methods to define if a state of the spectrum is in-gap. For this we turn to the topological measures of the Bott index and local Chern marker, as discussed in Sec. III A.
The meaning of an edge state must also be clearly defined. In spectral problems of infinite-size operators studied by the applied mathematics community, one of the main problems tends to be the removal of spectral pollution. Spectral pollution refers to a set of states in the spectrum which are not actually part of the infinite-size spectrum. These typically manifest in the form of edge states due to finite size effects. While the edge states do not exist in the spectrum of infinite-size operators, they are physical states of finite-size system, with distinct observable properties. Note, the BLT states outlined in this paper are not spectral pollution, and are in fact part of the spectrum of the infinite-size operator, which we will show in Sec. IV B. This means that not all in-gap states are spectral pollution, as is usually thought.

A. Hofstadter Vertex Model
We will consider lattices generated from the vertex model of aperiodic tilings. A 2D tiling is a countable family of closed sets (prototiles) which covers the entire 2D plane without any gaps or overlaps [26,56,57]. Aperiodic tilings are a subclass of tilings that exhibit long-range order, but no short-range translational invariance. Finding tiles that enforce quasiperiodicity is not a simple task, and the initial aperiodic tiling patterns contained thousands of distinct tiles [26]. Penrose discovered an aperiodic tiling requiring only a few rhombic tiles [58]. Since then, there has been a multitude of aperiodic tilings discovered with a variety of non-crystalline rotational symmetries [26,59]. We will focus on an AB vertex model as an example, which has 8-fold rotational symmetry and may be generated from an incommensurate rotation and projection of the 4D hypercubic lattice [60][61][62] (see Appendix B). We illustrate the quasicrystalline AB vertex model and its construction from an aperiodic tiling in Fig. 3. However, our results are not specific to this 8fold tiling, and we will show that BLT can occur in other quasicrystalline lattices in Sec. V.
The vertex model takes the aperiodic tiling and considers a lattice site to exist at each vertex and a bond to be present along the edges of tiles [63][64][65][66], as shown in Fig. 3(b) to (c). We will consider the vertex model of the AB tiling with a perpendicular constant magnetic field, as depicted in Fig. 1. The single-particle Hamiltonian is then where θ jk is the Peierls phase [67] due to the magnetic field between sites j and k, j, k is the sum over all N vertices/sites connected by an edge, and |j the state of a particle occupying site j. We consider the Landau gauge A(r) = Bxŷ = (φ/A)xŷ, with the magnetic field strength B, flux φ (measured in terms of the flux quantum φ 0 = 2π) and penetrating area A. We will take the area A to be that of the square tile of the AB tiling (the qualitative results, including the observation of BLT states, are independent of the choice of A). Units of = e = 1 are considered throughout this work, and we will work in units of energy J. The Hamiltonian (1) is well-understood when applied to periodic systems [8] and can even result in similar physics when applied to some quasicrystals [38][39][40]43].

B. Infinite-Size Algorithm
To rigorously probe the spectral properties of the infinite tiling directly, we will use a set of new computational techniques for infinite-dimensional spectral problems [68]. We begin with a description of computing the spectrum, with further details presented in Appendix A. As an example, the results of the infinite algorithm are shown for an infinite square lattice with a Hamiltonian of Eq. (1) in Fig. 4. The algorithm perfectly replicates the fractal Hofstadter butterfly usually generated through the consideration of periodic boundary conditions [8]. The removal of ESs, or in this case (of the infinite tile) spectral pollution, can be seen by the difference between Figs. 4(a) and (b). As discussed previously, spectral pollution is the terminology used to describe spurious eigenvalues present between parts of the essential spectrum.
We will utilise a new algorithm [55], developed by one of the co-authors, which allows for the calculation of the spectrum of a full infinite-dimensional operator with error control. This algorithm is general in its applications and will be of use to physical scenarios other than that considered here. For example, extensions to unbounded operators and partial differential operators can be found in Ref. [69]. For the present paper, the algorithm is of particular use since (i) the aperiodic nature of quasicrystals makes it a considerable challenge to approximate the spectrum of the full infinite-dimensional operator without finite-size effects and (ii) the approximation error can be computed and reaches effectively zero from a physical standpoint, as detailed below. There has also been an approach developed to obtain the exact solutions of quasicrystals of infinite size through the use of a superspace [70]. This method is specifically constructed to handle quasicrystalline problems by converting them to higherdimensional periodic problems and could be applied in the scenarios discussed in this work. However, we find the algorithm utilised here to be highly efficient at handling the problem of a quasicrystal of infinite size in a uniform magnetic field.
In infinite dimensions, the Hamiltonian H can be represented by an infinite Hermitian matrix,Ĥ = {Ĥ ij } i,j∈N , which acts on l 2 (N), the space of square summable sequences. A suitable ordering of the sites (e.g. ordering by positional radius from an origin) leads to a matrixĤ with finitely many non-zero entries in each column. In other words, there exists a function f : N → N such thatĤ ij = 0 if i > f (j), thus describing the sparsity ofĤ. Sparse Hamiltonians are a subclass of operators that are dealt with in Ref. [55] by considering the function where P m denotes the orthogonal projection onto the linear span of the first m basis vectors and σ 1 denotes the smallest singular value of the corresponding rectangular matrix. The rectangular truncation P f (n) (Ĥ − z)P n corresponds to including all of the interactions of the first n sites (the first n columns ofĤ) without needing to apply boundary conditions (see, for example, Fig.  1 of Ref. [71]). This is in sharp contrast to standard methods that typically take a square truncation of the matrixĤ (corresponding to a truncation of the tile) with a boundary condition. This difference allows us to prove convergence, provide error control, and also lends itself to adaptive computations of the full infinitedimensional operator. Physically, F n (z) is the squareroot of the ground state energy of the folded Hamiltonian P n (Ĥ − z) * (Ĥ − z)P n . In Ref. [55], it is shown that F n (z) converges down to the distance of z to the spectrum of H (uniformly on compact subsets of C) as n → ∞ which, together with a local optimisation routine, leads to the computation of the spectrum and approximate states with error control as n → ∞. Algorithmic steps are provided in Appendix A.

A. Topological Measures
Probing topological invariants in quasicrystalline systems is difficult, due to the ill-defined Brillouin zone and the breakdown of Bloch's theorem. For two-dimensional lattices subject to an external perpendicular magnetic field, the topological invariant of each respective band is its Chern number [6,7]. Differences between Chern numbers of the bands above and below a band gap are then equivalent to the number of ESs that appear within the band gap via the bulk-boundary correspondence. To obtain such topological invariants, integrals are normally taken across the first Brillouin zone of the system [6,72]. This explicitly means that the underlying lattice needs to be crystalline for these integrals to be generic properties of the bulk. However, there are measures that are independent of the boundary and have been shown to be equivalent to properties of the Chern number. These are the Bott index [73] and the local Chern marker [74], which we will utilise here.

Finite Systems
The Bott index is a spectral quantity defined for each individual state, related to the commutativity of two ma-trices [73,75]. It is defined for the nth eigenstate as withV n x/y being the projected position operators, for the nth state cumulative projectionŝ where |m is the mth eigenstate (we follow the usual convention of listing states in order of increasing eigenvalues) and the unitary diagonal position operators arê Note, thex S /ŷ S are position operators which must be scaled between 0 and 1. By calculating the Bott index, we can measure obstructions to the formation of a maximally localised Wannier basis that span the occupied states [76,77]. In order to find localised Wannier states, it is usually necessary to find continuous and periodic logarithms of families of unitaries [78]. This is discussed in detail in Ref. [79], where they also show that such logarithms can not be defined when the Chern numbers associated to the occupied states are non-zero. In essence the problem comes down to the definition of the fractional power of a matrix, and issues around the choice of a branch cut in the complex plane to define the required logarithm [79,80]. Note, this logarithm is not the same as the one utilised in the definition of the Bott index. If at a given energy the Bott index is non-zero, then any state with that energy must be in-gap, as the occupied band of states below must have a non-zero Chern number associated to it. The Bott index has proven useful in disordered and quasicrystalline systems, where bands cannot be defined through Bloch's theorem [41][42][43]77]. However, it is computationally expensive, as it requires the logarithm of a matrix whose size scales with the lattice. The Bott index does not lend itself to the large system sizes of this work and is, in general, ill-defined in the infinite-size case. Therefore, we will for the most part turn to an alternative measure in the local Chern marker.
The topological invariant can also be projected into the real physical space of the system, and for the Chern number this is known as the local Chern marker [74,[81][82][83]. Unlike the Bott index, the local Chern markers are defined on every single site j of the lattice for the nth eigenstate as with A c a reference area of the lattice and x n =Q nxP n ,ŷ n =P nŷQn , wherex /ŷ are the position operators. We will take A c to be the area of the square tile of the AB lattice. The local Chern marker has already been used to distinguish topological states in quasicrystals and disordered systems [38,84]. For large lattice sizes, we find that the local Chern marker is a more efficient way of distinguishing in-gap states since we do not need to compute a matrix logarithm for each state. Moreover, it can be extended to infinite systems. We will consider an effective Chern marker C n for a given state |n , which takes the integer of maximal counts in the distribution of C n j , which we find to be in agreement with the Bott index for crystals and quasicrystals.

Infinite Systems
Associated with the Hamiltonian H is a projectionvalued measure, E, whose existence is guaranteed by the spectral theorem [85] and whose support is the spectrum Sp(H). This diagonalises H, even when there does not exist a basis of normalisable eigenfunctions (recall that we are working in an infinite-dimensional Hilbert space): In finite dimensions and for compact Hamiltonians, E consists of a sum of Dirac measures, located at the eigenvalues, whose values are the corresponding projections onto eigenspaces. More generally, however, there may be a continuous component of the spectrum and spectral measure. Generalisations of the spectral projectors in Eq. (6) can be given in terms of E aŝ where we now label over energy values E, which also covers the possibility of continuous spectra.
The key ingredient that allows approximations of E to be computed is the formula for the resolvent In Ref. [86], it is shown how to compute the action of the resolvent with error control via the rectangular truncations P f (n) (Ĥ−z)P n . Using this, we compute a smoothed approximation of E via convolution with a rational kernel K for smoothing parameter > 0. Taking z = x + i , the classical example of this is Stone's formula which corresponds to convolution with the Poisson kernel As ↓ 0, this approximation converges weakly (in the sense of measures) to E. However, for a given truncation size, if is too small the approximation becomes unstable due to approximating the sum of Dirac masses that correspond to the spectral measure of the truncation of H. There is an increased computational cost for smaller , which typically requires larger truncation parameters.
Since we want to approximate spectral properties without finite-size effects, it is advantageous to replace the Poisson kernel with higher-order rational kernels developed in Ref. [87]. This allows a larger for a given accuracy, thus leading to a lower computational burden. The reason for using rational kernels is that, through a weighted distribution of resolvents, we can recover a generalised Stone's formula which converges with mth order of convergence in [87].
Here, c.c. denotes taking the adjoint, the constants α j and a j can be found in Appendix A, and * represents convolution. With this in hand, and for a given energy value E, we can write down smoothed generalisations of the spectral projectors in Eq. (6) aŝ where I denotes the identity operator. Note that the convolution [K * E] is a bona fide operator-valued function, and so the above definition makes sense. Finally, we definex and the smoothed infinite-dimensional local Chern marker on a basis site j up to energy value E as In the examples that follow, we will take = 0.05 and the 6th order kernel described in Appendix A, where we also give further algorithmic details. We will consider the same definition of an effective Chern marker, C n , for the infinite size as we do for the finite size.
As an example of this new infinite size algorithm for a topological measure, the results for the effective Chern marker for infinite size are shown for a square crystalline lattice with a Hamiltonian of Eq. (1) in Fig. 4. The Hofstadter butterfly has no states with a non-zero effective Chern marker present (as expected), apart from a handful of states which, for a given resolution, have a vanishing Chern marker as the algorithm converges. The algorithm developed here to compute topological properties for infinite size systems is of wide applicability and could be used to study other quasiperiodic, aperiodic, or even periodic lattices with complex structure.

B. Radial Measure
To distinguish between BLT states and ESs in the gap, we will define a radial measure of with ρ n j = |ψ n j | 2 /max |ψ n | 2 being the rescaled probability density if this quantity is at least 0.75max |ψ n | 2 , and ρ n j = 0 otherwise. Here, ψ n is the wavefunction of the nth eigenstate, N f is the number of elements for which ρ n j > 0, and r j is the jth site normalised radial coordinate with 0 ≤ r j ≤ 1. Therefore, L n is defined such that for every state 0 ≤ L n ≤ 1, with ESs having L n ∼ 1 and BLT states L n < 1. The radial measure gives the degree of which the density profile is localised towards the lattice centre. Note that, for regular bulk states, we also have that L n < 1. However, we distinguish between bulk states and in-gap BLT states via the topological measures already discussed.

C. Transport Properties
The most significant physics that emerges due to the presence of in-gap states is their resulting transport properties. For ESs, this means that transport is supported along the edge of the system. The BLT states characterised in this paper support transport along localised regions within the bulk of the lattice. For the finite lattice, this requires the evolution of the current state ψ(t 0 ) under the time evolution such that the final state is In our calculations, we will use a Trotter decomposition of the evolution unitary into discrete time steps. Note, that our Hamiltonian is always time-independent and we do not drive the system in any way. For the infinite size lattice, we cannot just apply the time evolution to a finite truncation, since we want to avoid finite-size effects. For a holomorphic function g, Cauchy's integral formula yields where γ is a closed contour looping once around the spectrum. Transport properties are computed via the choice g(z) = exp(−izt). The contour integral is computed using quadrature and approximations of the resolvent (H − z) −1 via rectangular truncations as above (see Ref. [86] for details and an example of fractional diffusion on a quasicrystal). In particular, the rectangular truncation of the Hamiltonian is chosen adaptively through a posteriori error bounds. This allows us to perform rigorous computations with error control that are guaranteed to be free from finite-size or truncation/discretisation effects, directly probing the transport properties of the infinite lattice. This is difficult to achieve via other methods such as truncating the tile since it can be difficult to predict how large the truncation needs to be a priori.

IV. BULK LOCALISED TRANSPORT STATES
As already discussed, for Hamiltonian (1)  in general two states that can be found: ordinary bulk states and in-gap ESs. In quasicrystals composed of single tiles, such as the Rauzy tiling, it is also observed that there are again two states, bulk states and in-gap ESs [38]. For Hamiltonian (1), this is all that would usually be expected unless there were perturbations or defects present. These perturbations could include, for example, the introduction of impurities, which can have in-gap states bound to them [13,[50][51][52][53][54], or the presence of internal hard edges, as in fractal lattices [47][48][49].
In quasicrystalline lattices constructed from multiple incommensurate tiles, a single peculiar in-gap bulk state was recently observed in Ref. [43]. However, the origins of this state were not known and it was introduced as a potential one-off peculiarity. In this section, we show that the state observed in Ref. [43] turns out to be one of many BLT states.
Examples of all three possible states, i.e. a bulk state, ES and BLT state are shown in Fig. 5 for a finite-size AB tiling vertex model with 1273 sites, along with their corresponding local Chern marker. Note, in this figure, and other figures shown later, we saturate the colour maps of the local Chern marker distributions to a range that shows the variation in values within the bulk. This is necessary because of large divergences in the local Chern marker near edges of the system, which occur in order for the sum of all local Chern markers in a given state to be zero [38,74]. In the bulk of the lattice, these fluctuations will generally be small, allowing for the effective Chern marker of a state to be visualised under a suitable range for the colourmap.

A. Scaling of the In-Gap States
First, we split our spectra into the ES, BLT state and bulk state components. We split the spectra by using the radial measure of Eq. (16), combined with the effective Chern marker of the state. We can see that if we vary the flux, as shown in Fig. 6, then the number of ESs or BLT states at any given flux can vary. There is also a characteristic dip in the in-gap states at a single flux within 0 ≤ φ/φ 0 ≤ 1 as is the case for ESs in a crystalline lattice.
There is an interesting trend in the BLT states compared to the ESs, in that the number of BLT states as a proportion of the total number of states is increasing with system size. This is contrary to what is expected for normal ESs in the gap, which will decrease as a proportion of the spectrum with system size, due to the area of the bulk increasing in size faster than the edge increases in linear size. By considering larger system sizes for the finite system in Fig. 7, we show that, as expected, the proportion of ESs tends towards zero. However, the proportion of BLT states increases with system size and appears to converge towards a non-zero number. In fact, the values near convergence are ∼ 20 − 40% of the total states, which is considerable and seen across a broad range of flux values.
The fact that BLT states do not make up a small percentage of the total states is a sign of their truly bulk nature. These are states that scale with the size of the bulk, rather than the edge. It is then reasonable to expect that their location in the lattice and density configuration will be as varied as that of the ordinary bulk states themselves. In other words, they can support transport over large and varied sections of the lattice bulk. This is in stark contrast to the transport supported by ESs, which is necessarily located along boundaries. It is important to mention here that the observed convergence with system size is no guarantee that the BLT states will be preserved for the infinite size lattice, as the usual thermodynamic limit approach to the infinite size is ill-advised in quasicrystals due to their aperiodic nature. We can then ask an intriguing question -are the BLT states preserved for the infinite-size quasicrystal?

B. A Zoo of Bulk Localised Transport States
We can now report that there is in fact a rich and varied zoo of BLT states for the AB vertex model in a magnetic field, for both finite and infinite systems. The BLT states are also prevalent in a variety of other quasicrystals, which we will discuss in Sec. V.
We first fix the flux to φ/φ 0 = 0. 69 and give examples of the varied structure of BLT states in Fig. 8. To find the large variety of BLT states realisable at this flux, we simply need to extend the system size from the previous consideration in Fig. 5. The examples shown in Fig. 8 are  for 7753 sites in the 8-fold lattice. We show that not only is the original BLT state retained, as shown in Fig. 8(a), but we also realise BLT states in different regions of the lattice, both far into the bulk and nearer the boundary of the system. The form of the states shown in Fig. 8 further explains the trend of the proportion of BLT states as we increase the system size shown in Fig. 7. As the quasicrystal becomes larger, there are more regions of the lattice to which the BLT states can localise to. We also show the local Chern markers for each state in Fig. 8, where it can be seen that the effective Chern marker of the bulk is non-zero, as the majority of sites in the bulk have a non-zero local Chern marker.
BLT states are not a peculiarity of a single flux, as was shown in Fig. 6, and we show a range of example BLT states at different flux in Fig. 9. All the states shown in Fig. 9 have a non-zero Bott index and effective Chern marker. It is clear from the examples of Fig. 9 that the BLT states are not an artefact of a single region or a subset of regions of the lattice. Instead, they appear throughout the system, with their location dependent on the flux. This hints that their origins are due to an interplay of the constant magnetic field with the quasiperiodicity of the lattice, which we will explore further in Sec. VII. The appearance of BLT states is as diverse as the usual bulk states observed on the quasicrystal, reflecting the quasiperiodic nature of the system. In Fig. 9(d,e), examples are shown which have localisation to a single site, or a rotationally symmetric set of single sites, with other examples showing intricate localised structures. A key property of the BLT states appears to be their ability to appear anywhere in the quasicrystalline lattice. Their ability to be localised at differing regions of the lattice could prove useful in future applications, as their position is not restricted and could be tuned without altering the lattice geometry.
The truly intriguing question for BLT states is whether they survive in the infinite-size system? We confirm that the BLT states can indeed exist in the infinite-size quasicrystal by applying the method of Sec. II B. A set of example states are shown in Fig. 10, all with a nonzero effective Chern marker. The similarity of these BLT states to the finite-size results already discussed is striking. Again, it is also shown that the region where BLT states localise is flux-dependent. The persistence of BLT states into the infinite-size is an important difference for the usual edge states of Hamiltonian (1). From periodic systems, we would expect all states in the infinite size to be in the bulk bands, with no in-gap states. However, with the BLT states preserved, a quasicrystal in a magnetic field can have states with properties that are usually considered to be related to being in-gap, even for the infinite-size lattice without defects or boundaries.

C. Prevalence of the Bulk Localised Transport States
The perfect, periodic fractal nature of the Hofstadter butterfly is not preserved for quasicrystalline systems. For quasicrystals, the energy-flux plane still contains a self-similar structure, which is, in general, aperiodic. The structure of the energy-flux plane for quasicrystals has been previously studied [37,43], and will now be utilised to show the prevalence of BLT states throughout the single-particle 'phase diagram' of Hamiltonian (1).
In Fig. 11 we show the finite-size Hofstadter butterfly for two lattice sizes, with the colour map being the effective Chern marker C n for each state and the corresponding radial measure of Eq. (16) for every in-gap state (i.e. non-zero C n ). From this figure, it can be seen that there is a significant number of states that are in-gap and within the bulk, i.e. with L n < 1 and C n = 0. We also show the energy-flux plane for the infinite-size quasicrys- Figure 11. Spectra of the quasicrystalline lattice, where each state is labelled according to the effective Chern marker C n and radial measure L n for (a,c) N = 1273 and (b,d) N = 3041 sites. The spectra (a) has similar properties to what is observed in the square lattice, but with a more intricate and quasiperiodic structure. For the larger lattice (b), one also observes that the total number of topological states throughout the spectra increases. For visualisation purposes, states with zero C n are plotted in light, translucent gray, in order to display only the locality of different in-gap states in (c,d). tal in Fig. 12. It is clear that BLT states are present for the infinite system through large ranges of flux, with states possessing non-zero effective Chern marker being present even after ESs are removed from the system. The regions where BLT states are present in the infinite-size quasicrystal map well to those present in the finite size.
The aperiodic Hofstadter butterflies in both the finite and infinite size show that the BLT states are not a single set of peculiar states limited to the examples shown in the previous section. Instead, the BLT states are present in the majority of parameter space for Hamiltonian (1), with BLT states even dominating the spectrum for particular ranges of the flux. We can go to larger values of Figure 13. Spectra of the quasicrystalline lattice, where each state is labelled according to (a) the effective Chern marker C n and (b) the radial measure L n across a larger range of φ/φ0 for N = 1273. As expected, the Hofstadter butterfly is aperiodic and has a rich, self-similar structure, with many different regions hosting BLT states.
the flux, as shown in Fig. 13, and observe even more BLT states. The prevalence of the BLT states and their variety could make the utilisation of their supported BLT particularly interesting.

V. BULK LOCALISED TRANSPORT STATES IN OTHER QUASICRYSTALS
We have seen how prominent BLT states are within the spectra of the AB vertex model. We now show that BLT states can also populate the spectra of other kinds of quasicrystals, both with and without rotational symmetries, as long as magnetic aperiodicity, i.e. the interplay of the magnetic field with the incommensurate areas of the quasicrystal, is retained.
We plot the infinite size Hofstadter butterfly for the 5fold lattice in Fig. 15, labelled according to the effective Chern marker. This again demonstrates the removal of conventional ESs from the spectra of the 5-fold vertex model and the retention of BLT states across a broad range of flux. In other words, the appearance of BLT states is not limited to particular values of the magnetic field or the exact geometry of the lattice, and they can dominate the spectra of quasicrystals with any rotational symmetries or structure, as long as magnetic aperiodicity is retained.
We now show this by considering the BLT states present in vertex models of 5-fold, 7-fold, 10-fold and 12fold rhombic tilings, which are all deduced from projections of higher-dimensional cubic lattices (see Appendix  In these examples, we also now have rhombic tilings with more than two prototiles, namely the 7-fold and 12-fold lattice. As the global rotational symmetry of the rhombic quasicrystal increases, the number of prototiles will also increase to ensure that no gaps are left in the tiling. The Hamiltonian is still that described by Eq. (1).
In Fig. 16, we plot several example BLT states for each Hofstadter vertex model in the infinite size. This clearly illustrates the appearance of similar-looking in-gap states to the ones observed in the 8-fold AB tiling. As expected, the structure of states on other vertex models have the  rotational symmetry consistent with the tiling itself, but preservation of the global symmetry is not important to the formation of BLT states (we will show explicit examples of this later in this section). The appearance and structure of BLT states in different quasicrystals is as rich as the examples shown for the AB tiling. All the states shown in Fig. 16 have non-zero effective Chern marker and, hence, will support transport along them in the bulk. These results confirm that BLT states as defined in this work are a property of magnetic aperiod- icity, independent of the underlying details of the lattice structure.
To demonstrate that BLT states are not a consequence of rotational (or reflection) symmetry we consider the AB and the 5-fold cases with two large arc segments removed, as shown in Fig. 17. This is a rather extreme deformation to the bulk of the system, similar to the consideration of substantially deformed boundaries for regular in-gap ESs. BLT states for different flux are also shown in Fig. 17, and all these states again have non-zero effective Chern marker. While the global rotational symmetry has been broken, the BLT states can still form and are prevalent throughout the spectrum, showing that they are truly a property of the bulk. The BLT states can even localise around regions with local symmetry, mimicking the BLT states observed for the non-deformed quasicrystal. Furthermore, we also observe that these BLT states may become much more highly localised to certain regions of the lattice than what is seen with global rotational symmetry. Importantly, they can also appear in this extremely perturbed lattice as states extended and connected in both domains, as shown in Fig. 17(b,f,g), allowing for transport around the deformations, akin to the transport of regular edge states around perturbations.

VI. BULK LOCALISED TRANSPORT
For any in-gap topological state, one of the most interesting properties is their support of transport. In the case of crystalline lattices in constant magnetic fields, this is usually considered by launching a state (or particle) along the edge of the system, and observing the robust transport of a component of the state around the boundary of the system. If a boundary is instead formed within the lattice, e.g. between a topological and nontopological region, then transport can also be supported along such features due to the presence of in-gap states bound to the interface.
The BLT states found in this work are another type of in-gap state, and they should be fully expected to support transport along their locality. Indeed, we find that a quasicrystal in a magnetic field can also support longlived BLT within the bulk, as depicted in Fig. 1 18 and 19 respectively for the AB tiling. It is clear, that a component of the initial state populates the BLT state, allowing for BLT to be supported. Note, here, we do not attempt to load into the BLT states, meaning we loose population into the other states that overlap with the initial state. For experimental scenarios, it would be prudent to instead use state preparation schemes or shortcuts to adiabaticity to load efficiently into the BLT states [88][89][90].
The transport shown for the infinite-size lattice in Fig. 19 is remarkable. For the infinite size, the usual transport along the edge is not present due to the lack of said edge. For the periodic infinite size, we can still have transport carrying states due to the presence of an interface or impurities. However, in the quasicrystal, an interplay of the magnetic field and quasiperiodicity can localise the particle to the bulk in the infinite size without any alterations to the lattice. These infinite size states are then as capable of carrying transport along them as any edge state in the finite system.
We also show that this is again not a fluke of the AB Hofstadter vertex model by showing examples of BLT for Figure 19. Infinite-size time evolution of BLT states at time t = 0J −1 (first column), t = 200J −1 (second column), t = 500J −1 (third column) and t = 1000J −1 (fourth column). The first two rows correspond to the BLT states in Fig. 18 (φ = 0.69φ0), but now computed on the full infinite tiling. The third row corresponds to an excitation of the BLT state at φ = 0.2φ0 on Fig.  10(c). Note, that our Hamiltonian is always time independent and we do not drive the system in any way. Figure 20. Infinite-size time evolution of BLT states for other quasicrystals at time t = 0J −1 (first column), t = 200J −1 (second column), t = 500J −1 (third column) and t = 1000J −1 (fourth column). The first row corresponds to the 5-fold BLT state in Fig. 20(b), with φ = 0.69φ0. The second row corresponds to an excitation of the 12-fold BLT state on Fig. 20(g), with φ = 0.4φ0. some other quasicrystals in Fig. 20. These results are for the infinite size lattice and are again remarkable in their form. We show examples of BLT for the infinite tilings with 5-fold and 12-fold symmetries. Again, we observe transport supported by the BLT states which is confined within the bulk of the lattice and possible in systems regardless of the details of the boundary.
A key property of BLT could be the possibility of using its varied location throughout the lattice itself. In a crystalline finite system, the transport is only along the edge, or any internal edges imposed by the lattice structure and/or defects. However, with BLT we have shown that it is possible to have transport supported in multiple locations within the lattice for any given flux. By varying the flux, we can also tune where the BLT is supported. BLT can then give a degree of control for transport to occur within the bulk. By showing that this can be done in the infinite size, we have shown that BLT is independent of the size of the lattice, or the exact form of the boundary. This is in stark contrast to in-gap ESs, whose presence can in some cases be entirely reliant on the geometry of the edge [91].

VII. FORMATION OF BULK LOCALISED TRANSPORT STATES IN A TOY MODEL
So far, we have demonstrated that BLT states and their supported BLT are prevalent throughout the spectra and lattice of quasicrystals in magnetic fields. In this section, we will characterise how these states form through a toy model. First, we will consider a simple but misguided toy model that shows the BLT states are not formed by effective edges from the varying local coordination number. These are equivalent to dislocations being present in the system. We will then describe a simple toy model which shows that the BLT states arise due to the interplay of the quasiperiodic lattice, and mainly the irrational areas present, with the constant magnetic field.

A. Dislocation Toy Model
At first glance, it is easy to think that the BLT states must form only due to the quasicrystalline nature of the lattice alone. This would be through effective edges being formed in the system via the local aperiodic variation in the coordination number for each site. In many ways, this would be similar to edge states being bound to a dislocation or defect in the lattice structure. We, therefore, consider a toy model on a square lattice where we have a central region with lattice constant l 2 = 13τ 10 , where τ is the golden ratio, and an outer region of lattice constant l 1 = 1. The lattice constants l 1 and l 2 are incommensurate and lead to a dislocation along the boundary of these two lattices as shown in Fig. 21(a). We couple all sites along this dislocation that are within 1.5l 1 of each Figure 21. Zoomed in portions of the square lattice toy models used to look at the formation of BLT states. In (a), we show a square lattice with period l2 = 13τ 10 embedded in a larger square lattice of period l1 = 1. Connections, coloured red, are generated between the boundaries of both lattices that are less than a threshold of 1.5l1 apart, resulting in a coordination number dislocation. In (b), we show a square lattice with lattice spacing l that has a positional disorder applied to a selection of sites. The sites coloured yellow are randomised from their original position in a small radius of 0.25l units. The bottom row of figures (c,d) show the bottom right corners of the same lattices in order to emphasise the incommensurate areas between different unit cells more clearly.  Some examples of the states along the dislocation are shown in Fig. 22. The density profiles of these states look similar to those of the BLT states but with a crystalline 4-fold rotational symmetry, as would be expected for this toy model. They also have a corresponding non-zero Bott index and are, therefore, in-gap and will support transport. However, if we look at the local Chern marker, we observe a striking difference to the BLT states in the quasicrystal. The states along the dislocation are a product of the change of the local Chern marker across the dislocation. This in itself is not a surprise and is usually the reason why in-gap states appear on internal edges in these systems [47][48][49]. It is clear though that this is not how the BLT states form, as we do not see any change in the local Chern marker across the interface for quasicrystalline BLT states. Therefore, even though at first glance, the variation in the coordination number appears to map well to how the BLT states are formed, it cannot be their true origin.

B. Magnetic Aperiodicity Toy Model
The formation of BLT states is due to an interplay of the magnetic field with the aperiodicity of the quasicrystalline structure. This is to be expected from the results in Sec. IV B, where the location of the BLT states was shown to be largely dependent upon the applied magnetic field strength. Through another toy model on an originally square lattice, we can show that the magnetic aperiodicity introduced into the system is responsible for forming BLT states. In this model, we take a subset of lattice sites on a square perimeter and vary their location to a small degree, while retaining the constant coordination number and connectivity of the square lattice. An example of this toy model is shown in Fig. 21(b), where we also enforce a 2-fold rotational symmetry for comparison. This toy model then has a small region where there is a disorder in the area of tiles. This would then map to a disorder in the flux penetrating tiles in this region, as we have cells with irrational areas to each other.
Some examples of the states along the magnetically disordered region are shown in Fig. 23. The states and their corresponding local Chern markers look very similar to the BLT states observed throughout this work. As expected, the states also have a 2-fold rotational symmetry and show clearly that the formation of BLT states is due to the interplay of the constant magnetic field with the self-similar aperiodic structure of the lattice. This would also explain why BLT states were not observed before in quasicrystals composed of a single tile such as the Rauzy tiling [38], as there is then no magnetic aperiodicity introduced into the system. Note, the Rauzy tiling does have a local variation in the coordination number, further supporting the conclusions from the dislocation toy model. Finally, it is worth noting that if periodic boundary conditions are applied to the hard edges of these toy models, then the BLT states remain persistent in the spectra. In other words, this again demonstrates that the formation of BLT states is not linked to the presence or geometry of the system's boundary.
To observe BLT states, it is therefore key to have the aperiodic nature of the lattice compete with the magnetic field. The simplest way to do this is to have a lattice constructed from building blocks that have incommensurate areas to one another, like the two prototiles of the AB tiling. We would then expect that BLT states are present in all quasicrystalline lattices which do not have a single, or multiple commensurate building blocks. We would also expect that the presence of BLT states should be independent of the vector potential, as long as the choice does not destroy the interplay with the incommensurate building blocks of the lattice.

VIII. CONCLUSIONS
We have shown that the conventional picture of insulators, metals, and topological insulators with surface states is not the full story for quasiperiodic systems. When the quasiperiodic nature of the quasicrystalline lattice interacts with the magnetic field, it is possible for states to form that are unique in their character. These states are localised to the bulk, but are in-gap and support transport along them in the bulk. The BLT states are far from a peculiarity and can exist throughout the spectra of 2D quasicrystals in magnetic fields due to a magnetic aperiodicity. Magnetic aperiodicity is not an artificial construct and is natural in the majority of 2D quasicrystals, due to the incommensurate nature of the building blocks of the lattice. For the case of the vertex models studied in this work, this incommensurate nature arises from the incommensurate areas of the prototiles.
BLT states are even sometimes the dominant in-gap state for large finite system sizes. We have confirmed that the BLT states can exist in a variety of regions within the bulk of the lattice, with their position being dependent on the flux, or equivalently the magnetic field strength. Through the use of a new numerical technique, we have also shown that BLT states are present in the infinite-size lattice. This is quite remarkable as they are in-gap states, with corresponding non-zero topological measures. By exciting regions within the bulk of the lattice, we have also shown that BLT is supported in both the finite and infinite size.
A future direction of interest could be to realise BLT states of quasicrystals in an experimental setting. While there are current developments of realising quasicrystalline problems in cold atoms [29,32], a promising setting for current realisation of the BLT states could be photonic lattices [92]. Photonic lattices allow a high degree of controllability in the lattice geometry and are favourable for realising the BLT states due to a large number of sites being possible. The lattice is usually etched into the two-dimensional plane of a fused silica crystal, with time then being transposed onto the third dimension of the crystal. Photonic lattices are good simulators for single-particle physics, and synthetic gauge fields can be realised via helical waveguides [93], a Floquet step-like approach [22,94], or a strain across the lattice [95]. Furthermore, there has also been a recent proposal to realise a Penrose quasicrystal in a synthetic vector potential using helical waveguides [35]. As the key to the realisation of BLT states is the magnetic aperiodicity, we would expect that any of these current techniques could potentially be utilised to probe the physics of BLT states.
The discovery of a zoo of BLT states in quasicrystals opens a number of intriguing open problems and future research directions. One possible line of work is to consider the nature of the spectra for the infinite-size tiling. Usually, we would consider these states to all be ordinary bulk states, but as we have shown this is not necessarily the case, even without any impurities or true disorder in the system. Another interesting question is how the BLT states would appear or alter the physics in the presence of interactions. The study of the physics of two-dimensional quasicrystals including interactions is an emergent topic [96][97][98]. With the recent advances in many-body numerical techniques [99][100][101], the physics of BLT states and their ramifications in quasicrystalline lattices in the many-body regime could soon be probed.
The reliance on magnetic aperiodicity also initiates an interesting set of questions. This is not a sole property of quasicrystals and could be incorporated into crys-tal structures like that considered briefly to generate an aperiodic Hofstadter butterfly in Ref. [43], or extensions of the toy model of magnetic aperiodicity considered in this work. This approach could be used to design effective regions in the system where BLT states are desired to localise, through the appearance of incommensurate phases. These regions could then be tuned to support or not support transport through changes in the applied magnetic field strength. This would potentially allow the design and control of specific regions supporting localised transport within the lattice bulk, with interesting applications.
The understanding of the physics of electronic-like states in quasicrystals is at an early stage of its development, especially concerning potential applications to quantum technologies. There is much work to be done to realise the potential of these fascinating and complex structures. This work has shown one possible exotic behaviour of these systems in magnetic fields; the existence of BLT states and their supported transport. Future work to understand the applications of BLT states to quantum problems, like the applications considered for typical edge states [3,6,102,103], could be particularly fruitful.
Extensions to unbounded operators and partial differential operators can be found in Ref. [69].
Recall that in our setting, the Hamiltonian H can be represented by an infinite Hermitian matrix,Ĥ = {Ĥ ij } i,j∈N and we are given a function f : N → N such thatĤ ij = 0 if i > f (j), thus describing the sparsity of H. For z ∈ R, the key quantity to compute is where P m denotes the orthogonal projection onto the linear span of the first m basis vectors and σ 1 denotes the smallest singular value of the corresponding rectangular matrix. The function F is an upper bound for the distance of z to the spectrum Sp(H), and converges down to this distance uniformly on compact sets as n → ∞. There are numerous ways to compute F n , such as standard iterative algorithms or incomplete Cholesky decomposition of the shifts P n (Ĥ − z)P f (n) (Ĥ − z)P n (see the supplementary material of [55] for a discussion). The other ingredient is a grid of points G n = {z (n) 1 , ..., z (n) j(n) } ⊂ R providing the wanted resolution r n over the spectral region of interest.
The algorithm is sketched in Algorithm 1, whereF n denotes the described suitable approximation of F n (which can be computed in parallel). The simple idea of the method is a local search routine. IfF n (z) ≤ 1/2, we search within a radiusF n (z) around z to minimise the approximated distance to the spectrum. This gives the set M z which is our best estimate of points in the spectrum near z. The output is then the collection of these local minimisers. The algorithm's output, Γ n (H), converges to the spectrum Sp(H) of the full infinite-dimensional operator as n → ∞ (for suitable r n → ∞). Note that this convergence is free from edge states. Moreover, the error bound of the algorithm satisfies sup z∈Γn(H) dist(z, Sp(H)) ≤ E n (A2) and the output E n converges to zero as n → ∞ (proven in Ref. [55]). For a desired accuracy δ > 0, we simply increase n until E n ≤ δ. In our numerical experiments we chose δ to be smaller than the required spectral resolution. Finally, the output V n consists of the approximate states corresponding to the output Γ n .
Algorithm 1 Computation of spectrum and the associated approximate states with error control. The computation ofF n can be performed in parallel.
1: For z ∈ Gn, approximate Fn(z) to accuracy (2rn) −1 from above. Call the approximationFn(z) and assume it takes values in (2rn) −1 Z. 2: For z ∈ Gn, let vn(z) denote the approximation of the right-singular vector of P f (n) (Ĥ − z)Pn corresponding to the smallest singular value. Otherwise, set Mz = ∅.

Computing Spectral Measures and Local Chern Markers
In this section, we will assume access to a routine that approximates the action of the resolvent (H − z) −1 on a vector with error bounds. In the scenario of the current paper, this can be done through the rectangular truncations P f (n) (Ĥ − z)P n and solving the resulting overdetermined linear system in the least squares sense. The residual converges to zero as n → ∞ and can be used to provide the needed error bounds through an adaptive selection of n (see [86,Th. 2

.1]).
We use the high-order kernel machinery developed in Ref. [87], where the following definition is made.
Definition 1 (mth order kernel) Let m ∈ N and K ∈ L 1 (R). We say K is an mth order kernel if: (ii) Zero moments: We set K (·) = −1 K(·/ ). In Ref. [87], theorems were proven on the rates of approximating a probability measure µ by the convolution K * µ. For mth order kernels, under suitable conditions, mth order convergence in holds (i.e. the error scales as m up to logarithmic factors). These rates can be carried over to approximating the projection-valued measure E described in Sec. III A. High-order kernels can be constructed using rational functions as follows. Let {a j } m j=1 be distinct points in the upper half plane and suppose that the constants {α j } m j=1 satisfy the following (transposed) Vandermonde system: Then the kernel is an mth order kernel [87], and we have the following generalisation of Stone's formula (A5) As a natural extension of the Poisson kernel, whose two poles are at ±i, we consider the choice We then determine the residues by solving the Vandermonde system in Eq. (A3). The first six kernels are explicitly written down in Table I. With this in hand, and for a given energy value E, we compute the smoothed spectral projections in Eq. (13) using the trapezoidal rule. The quantitiesx E andŷ E can be computed via successive applications of the relevant projectors. This is outlined in Algorithm 2, which computes the local Chern markers over a grid of energy values of spacing ∆E. In practice, the algorithm has two levels of parallelism. We can compute resolvents in parallel across different energy values E j and we can perform the algorithm in parallel for different sites indexed by i. Algorithm 2 Computation of local Chern markers. There are two levels of parallelism: the computation of resolvents and the entire algorithm for different sites indexed by i.
Input:Ĥ, f , m ∈ N (order of kernel), > 0 (smoothing parameter), ∆E (energy or spectral spacing), i (site index), Ac (reference area of the lattice), L (lower bound for Sp(H)) and M ∈ N (number of energy values). where K is the mth order kernel from

Computing Transport Properties
Finally, we will discuss the computation of transport properties. Given an initial wavefunction ψ 0 , we wish to compute where γ is a closed contour looping once around the spectrum. Suppose that the spectrum is located in an interval [a, b] ⊂ R. We take γ to be a rectangular contour split into four line segments: two parallel to the imaginary axis with real parts a − 1 and b + 1 and two parallel to the real axis with imaginary parts ±η (η > 0). Along these line segments we apply Gaussian quadrature with enough quadrature nodes for the desired accuracy (the number of nodes can be found by bounding the analytic integrand). Suppose that the weights and nodes for the quadrature rule applied to the whole of γ are {w j } N j=1 and {z j } N j=1 . Then the approximation of (A7) is given by  Table I. The numerators and residues of the first six rational kernels with equispaced poles (see Eq. (A6)). We give the first m/2 residues because the others follow by the symmetry αm+1−j = αj.
The (H − z j ) −1 ψ 0 are computed using the adaptive method outlined in Sec. A 2, which can be performed in parallel across the different quadrature nodes. We also reuse these computed vectors for different times t (numerically, this requires η to not be too large and suitable N can be selected for a finite interval of desired times t).

Appendix B: Construction of Aperiodic Tilings
The quasicrystalline tilings considered in this work are deduced from 2D cut-and-project sets of a D dimensional lattice Z. With this method, the key idea is to map lattice points in R D → R 2 . We consider rhombic tilings in this work, which requires that the Z are hypercubic structures. The set of basis vectors of Z are then simply the permutations of 1 in R D zero vectors. In order to produce the tiling, we first define a rotation on the points V ∈ Z relative to the origin 0 where V is defined by the basis vectors of Z, R is an incommensurate rotation operator and W is a transformed position. For hypercubic lattices, the columns of R naturally define the transformed basis vectors, leading to the following constraints on R where C i is the ith column of R in vector form and |R| = 1 for a rotation operator. This is effectively a statement that the normalisation and orthogonality of basis vectors should be invariant under rotation. We can now define two unique subspaces of the R D superspace of Z. First, we have the tiling space T ∈ R 2 , which is the 2 dimensional projection of the transformed points W . In other words, the points in T are defined by the first two elements of W . Second, we have the internal space I ∈ R D−2 , which is the D − 2 dimensional projection of W . Similar to before, a point in I is defined by the next D − 2 elements of W , after the points T .
To now form the tiling, we also define a cutoff for points in I. This is done by taking the projected unit cell of Z as a bounding volume for points in I. The reason this is necessary is due to the fact that the projection of all W to T will lead to a dense set of points, with no quasicrystalline structure and many overlapping tiles. Therefore, the final tiling in T now only accepts points whose dual in I is bounded by the projected unit cell, ensuring no tiles overlap.
For each tiling generated from a D dimensional hypercube, there may also exist a family of isomorphism classes. These are tilings with equivalent or up to a 2×D rotational symmetry, but with different local properties and frequency of tiles. These can be found by simply offsetting projected points in I by a small vector s.