The Fibonacci quasicrystal: case study of hidden dimensions and multifractality
TThe Fibonacci quasicrystal: case study of hidden dimensions andmultifractality
A. Jagannathan
Laboratoire de Physique des Solides,Universit´e Paris-Saclay, 91405 Orsay,France (Dated: January 1, 2021)
The distinctive electronic properties of quasicrystals stem from their long range struc-tural order, with invariance under rotations and under discrete scale change, but withouttranslational invariance. This review introduces the multiple manifestations of multifrac-tality in the Fibonacci chain, a 1D paradigm for quasicrystals. Quasiperiodic systemshave been recently of interest furthermore for their topological properties, describedin this review for the simplest 1D case. Some important generalizations of the basictight-binding models, and closely related models for other quasiperiodic systems arealso discussed.
CONTENTS
I. Background 1II. Geometric properties of the Fibonacci chain 2A. Substitution method 21. Inflation matrix 3B. Cut-and-project method 41. Structure factor 52. Conumbering scheme 6C. The characteristic function method 6III. Models and methods. Exact results 6A. Diagonal and off-diagonal Fibonacci models 7B. Multifractal energy spectra 8C. Gap labeling and topological indices 9D. Trace map method 10E. Log-periodic oscillations 12F. The wavefunction for E = 0 12G. Chern numbers. Bulk-edge correspondence 14IV. Approximate methods 15A. Perturbative methods 15B. Approximate renormalization group 161. RG for hopping model (1) 162. RG for hopping model (2) 173. RG for diagonal model 17V. Perturbative RG theory of multifractal spectra andstates 17A. Multifractality of the energy spectrum 18B. Gaps, stable gaps and topological numbers 19C. Multifractality of wavefunctions 201. Hamiltonian in the conumber basis 212. Energy-position symmetry 213. Multifractal exponents for wavefunctions 214. Spectrally averaged dimensions 22VI. Dynamical properties 22A. The diffusion exponent β E = 0 26B. Landauer approach 26 C. Kubo-Greenwood approach 27D. Many body metallic and insulating states 27VIII. Perturbations. Disorder and boundary effects 28A. Finite disorder and approach to Andersonlocalization 28B. The proximity effect 29IX. Generalized Fibonacci models 29A. Phonon models 29B. Mixed Fibonacci models 30C. Interference and flux dependent phenomena 31X. Related quasiperiodic models 31A. Other substitutional chains 32B. Product lattices 32XI. Interactions and quasiperiodicity 33A. Heisenberg and XY chains 33B. Many body localization 34C. Anomalous diffusion properties 34XII. Experimental systems 34XIII. Summary and Outlook 35Acknowledgments 36References 36 I. BACKGROUND
The Fibonacci chain is a one-dimensional quasiperiodicstructure that is closely related to the three dimensionalicosahedral quasicrystals discovered by Schechtman et al(Shechtman et al. , 1984). The study of electronic proper-ties of quasicrystals thus logically begins with the studyof electrons in a 1D Fibonacci chain. Tight-bindingHamiltonians for the Fibonacci chain have been exten-sively investigated. Properties of the energy spectra andcritical states of the chain have been studied using a vari-ety of methods, exact solutions, perturbation theory andnumerical analysis. Aperiodic Schrodinger operators in a r X i v : . [ c ond - m a t . s t a t - m ec h ] D ec general and quasiperiodic systems in particular are anactive subject of mathematical physics, and are reviewedin (Damanik et al. , 2015).Quasiperiodic Hamiltonians are of growing interest fortheir nontrivial topological properties (Huang and Liu,2018, 2019; Verbin et al. , 2013). Like the well-knownAubry-Andr´e-Harper (AAH) model with which it is of-ten compared, the one-dimensional Fibonacci quasicrys-tal possesses a higher dimensional “parent” system, fromwhich it inherits topological characteristics. Thus, inthe Fibonacci chain, there appear topologically-protectedboundary states equivalent to the edge states of the two-dimensional Integer Quantum Hall Effect (Verbin et al. ,2015). Unlike the AAH model, which is critical at a singlepoint in its parameter space, the Fibonacci models whichwe discuss here are critical for all values of the strengthof the quasiperiodic modulation. Multifractal states areomnipresent in the phase diagram of quasicrystals. Evenif a complete theoretical description is still lacking forthese systems, there have been a number of recent workson the spatial characteristics of multifractal states in qua-sicrystals (Kalugin and Katz, 2014; Mace et al. , 2017) andsome of their physical manifestations such as the proxim-ity effect (Rai et al. , 2019b), or in dynamical phenomenaincluding growth of entanglement (Mac´e et al. , 2019).We also show how disorder affects electronic propertiesof the quasicrystal. Although the focus will be primarilyon single particle properties, we will briefly discuss inter-acting quasiperiodic systems. In particular, many bodylocalization (MBL) due to quasiperiodic potentials (Iyer et al. , 2013) has been an active research topic. One ofthe questions concerns the differences in critical behav-ior, if any, from MBL due to random potentials (Khemani et al. , 2017). It is becoming possible to study a numberof models experimentally with cold atoms in optical po-tentials, possibly extending to realizations of generalizedFibonacci problems (Singh et al. , 2015).The review is organized as follows: Sec.II on geometri-cal aspects presents details of the structure of Fibonaccichains, approximant structures and the structure factor.Sec.III presents the basic models and lists some impor-tant exact results. Sec.IV introduces the perturbativerenormalization group technique and its qualitative pre-dictions. Sec.V discusses in detail some properties ofthe pure hopping model. Sec.VI discusses wave packetdynamics. Sec.VII presents some results for chain con-ductances. Sec.VIII discusses aspects of disorder andboundary effects. Sec.IX discusses Hamiltonians thatcombine diagonal and off-diagonal quasiperiodic modula-tions. Sec.X presents related 1D models and extensionsto higher dimensions. Sec.XI lists results for interactingquasiperiodic systems. Sec.XII presents experimental Fi-bonacci systems. Sec.XIII concludes with a summary andoutlook. II. GEOMETRIC PROPERTIES OF THE FIBONACCICHAIN
The Fibonacci chain is a 1D quasicrystal according tothe revised definition of the IUCr (IUCr, 1992). The def-inition states that a quasicrystal, like a periodic crystal,is a material having a sharp diffraction pattern composedof Bragg peaks. This section describes how to generatethe Fibonacci chain and reviews some of its geometri-cal and structural properties. We will introduce severalmethods, each of which is helpful in its way for a betterunderstanding of electronic properties in this system.
A. Substitution method
The substitution method explicitly introduces the no-tion of scale invariance of the quasicrystal, later usedin the renormalization group transformation. The Fi-bonacci substitution rule, σ , acts on the two letters A and B and transforms these as follows: σ : (cid:40) A → ABB → A. (1)Letting the substitution act repeatedly on the letter B generates a sequence of words C n = σ n ( B ) of increasinglength, as shown for the first few members in Table 1.These chains are finite approximants of the Fibonaccichain, which is obtained in the limit n → ∞ . It is easyto see that the lengths of the words are equal to theFibonacci numbers F n , defined by the recursion relation F n = F n − + F n − with F = F = 1. The ratio of twoconsecutive Fibonacci numbers tends to the golden meanas n → ∞ , F n − F n − = τ n lim n →∞ τ n = τ (2)where τ n are the rational approximants of the goldenmean τ = (1 + √ /
2. The lengths of the chains aregiven by F n ∼ τ n in the large n limit. n C n F n However, where a d dimensional periodic crystal needs at most d reciprocal lattice vectors to index all of its Bragg peaks, a d dimensional quasicrystal requires a greater number D > d ofreciprocal lattice vectors.
FIG. 1 Illustration of inflation transformations, going fromthe C (top) to C (middle) to C (bottom) chain. A (resp.B) tiles are shown by grey (resp. black) rectanglesTable 1. The first six approximants built using thesubstitution rule Eq.1 Inflation/deflation of tiles
The substitution methodshows the hierarchical relations between the chains andsuggests that problems on the chain could be tackledusing renormalization group methods. Consider a 1Dtiling of A and B tiles such that the ratio of their lengths l A /l B = τ . The approximant chains C n correspond toa series of finite tilings, which can be transformed intoone another by so-called inflation and deflation opera-tions. Using the substitution 1 in reverse, one goes froma chain of F n tiles to a chain of F n − tiles. Note that thiscorresponds to a “site decimation” process which elimi-nates a certain subset of sites. Rescaling all the tiles by afactor of τ n restores the length of the chain to its originalvalue, as illustrated in Fig.1. The infinite chain is invari-ant under inflation/deflation – i.e. the FC has a discretescale invariance. Concatenation
From Fig.1 (or Table 1) one sees thatthe n th chain C n can be obtained by the concatenationof two shorter chains C n − and C n − . This propertywill be useful later for the transfer matrix method III.D.Repeating this operation, one obtains a relation betweenthe n th chain and n − n − C n = C n − ⊕ C n − ⊕ C n − (3)where the symbol ⊕ here denotes concatenation (joiningthe chains in the specified order from left to right). Thisrecursive construction holds also for the energy spectrumin perturbative RG, as we will later see in Sec.V.
1. Inflation matrix
Let N ( n ) A and N ( n ) B be the number of occurrences of Aand B in an approximant chain C n . From the substitu-tion rule 1, it is easy to see that N ( n +1) A = N ( n ) A + N ( n ) B and N ( n +1) B = N ( n ) A , with the initial condition N (0) A = 0and N (0) B = 1. This relation can be put in a matrix form, (cid:34) N ( n +1) A N ( n +1) B (cid:35) = (cid:20) (cid:21) (cid:34) N ( n ) A N ( n ) B (cid:35) (4) where the 2 × M . The eigenvalues of M are λ = τ , and λ = − τ − .The corresponding eigenvectors are { τ, } and { , − τ } .The first eigenvector gives the relative frequencies of theA- and B-tiles, by virtue of the Perron-Frobenius theo-rem. The ratio N A /N B tends to τ when n → ∞ . Formore details on symbolic substitutions, see (Baake andGrimm, 2013a).At this point, it is useful to make a brief digressionrelating to other types of self-similar binary chains gen-erated by the substitution method.1. The so-called silver mean chain can be obtainedfrom repeated application of the rule σ Ag : (cid:40) A → AABB → A. (5)It can be easily checked, by writing the inflationmatrix for this case, that the eigenvalues satisfy theequation λ − λ − λ = ( √ B → A, A → A n B ( n ≥ α of an n th degree monic polynomialequation with integer coefficients, such that α isgreater than 1, while all other roots are of modulusless than 1).2. Consider next the substitution rule B → A, A → ABBB , giving rise to a self-similar aperiodic struc-ture which is however not a quasicrystal. Thisstructure is non-Pisot – it can be checked that thetwo eigenvalues of the inflation λ = (1 ± √ / Fluctuations
Consider a subsystem of N letters fromthe infinite Fibonacci chain. For values of N which arenot Fibonacci numbers the number of A’s and B’s in thesample δ N = N A − N B fluctuate around the mean value δ . The behavior of δ N for large N is described in termsof η , the so-called “wandering exponent” (Luck, 1993) as δ N − δ ∼ N dη (6)where d is the space dimension. In periodic systems, fluc-tuations are sub-extensive (are due to boundary effects)and η = 0. Fluctuations are similarly negligible in thethermodynamic limit for the which are Pisot structures:including the Fibonacci chain, and all the well-knownquasiperiodic systems in higher dimensions, such as thePenrose 2D and 3D tilings, Ammann-Beenker tiling. Forthese structures η = 0, and as a consequence their diffrac-tion spectrum is pure point, i.e. consists solely of Braggpeaks (Godr`eche and Luck, 1992).It is instructive to compare this with geometry fluctu-ations in random systems. Consider a random sequenceof the letters A and B for some fixed probability p A and p B = 1 − p A of the letters. The law of large numbers leadsto the exponent η = 1 / η = ln | λ | / ln λ (Godr`eche and Luck, 1992; Luck, 1993).This property of bounded geometrical fluctuations ofthe quasicrystal is the principal reason for electronicstates in quasicrystals being relatively more extended,as compared with, for example, the critical states at themetal-insulator transition in disordered structures. B. Cut-and-project method
The cut and project method is a means of generatingthe quasicrystal starting from a higher dimensional pe-riodic lattice. Here this is a 2D square lattice, and thequasicrystal is had by projection onto the physical axiscalled x (cid:48) in Fig.2. To be selected, a point must lie withinthe red strip S of slope given bytan θ = ω (7)where the notation ω = 1 /τ = ( √ − / x (cid:48) axis (the physical axis), horizontal bonds and verticalbonds project onto the red and blue intervals respectively.The y (cid:48) axis is called the perpendicular (or internal) space.Expressing all lengths in units of the square latticeparameter a , vertices of the square lattice are located at (cid:126)R mn = m(cid:126)x + n(cid:126)y where m and n are integers. The x (cid:48) -coordinate and the coordinate along the perpendiculardirection y (cid:48) are given by x (cid:48) = m cos θ + n sin θy (cid:48) = − m sin θ + n cos θ (8)upto an overall shift. To be selected, the point mustsatisfy the condition 0 ≤ y (cid:48) < W , where W = sin θ +cos θ is the cross-section of the strip S .After projection on the x (cid:48) axis, the spacing be-tween nearest neighbors can have two values: cos θ =1 / √ ω or sin θ = ω/ √ ω , corresponding to thelengths of A and B tiles of the previous section. Translational symmetries, phonons and phasons
The ori-gin of the x (cid:48) y (cid:48) axes is arbitrary: i.e. shifting the strip per-pendicularly to itself does not result in a new quasiperi-odic structure. This statement requires a clarificationof what is meant by equivalence between quasiperiodictilings. Tilings are said to be equivalent (or locally iso-morphic) if the same sequences of tiles – which may beof arbitrarily large size – can be found in both struc-tures, and with the same frequency of occurrence. Tosay it differently, one can make the two tilings overlapout to arbitrarily large distances, by suitably translatingone with respect to the other.The higher dimensional representation of the FC showsthat the structure has a translational symmetry with re-spect to displacements parallel to the physical space, butalso in the perpendicular direction. This leads theoreti-cally to the possibility of having two kinds of Goldstonemodes : phonons and phasons. Phonons are 1D vibra-tional modes in the usual sense. To visualize phasonmodes, one can look at Fig.3 which shows how the pro-jected structure would change if the lattice were to beslightly displaced in the direction perpendicular to thestrip. One of the vertices ( m, n ) in the figure moves outof the strip while, simultaneously, the vertex ( m − , n +1)enters the strip. The net result is a small discontinu-ous jump of one site, while the other points in its neigh-borhood remain unaffected. This so-called phason jumpcorresponds to exchanging the A and B tiles around agiven vertex, AB → BA. The shift produces a new chainstructure which is equivalent to the old one. In acceptedterminology, a phason mode is a coherent excitation ofthe perfect Fibonacci chains, corresponding to long wave-length fluctuations. In contrast, when phason flips areintroduced in a random uncorrelated fashion all along achain, this would give rise to a geometrically disorderedchain. In practice, however, spontaneous phason flipspresumably have a significant energy cost, and are un-likely to be excited at low temperatures.
Approximants
Periodic approximants of the FC aregenerated by taking a rational slope for the strip S . If theslope is chosen to be the ratio of two successive Fibonaccinumbers tan θ n = F n − /F n − then the projected chain has a repeating structure con-sisting of F n − and F n − tiles of type A and B respec-tively, and the total number of tiles is F n . It is eas-ily checked that one obtains the same approximant se-quences already seen in Table 1. For example, for n = 3the strip has a slope equal to , and one obtains a peri-odic repetition of the motif ABA . FIG. 2 Schema of the cut-and-project method. Selected points (joined by a broken line) of a 2D square lattice are projectedonto the x (cid:48) axis, giving the binary quasiperiodic sequence of red and blue tiles. The infinite selection strip S is colored red.FIG. 3 Illustration of a phason flip or local rearrangement of the tiling when the lattice is slightly displaced with respect tothe selection strip. The old and new structures are obtained by projecting the blue(red) broken lines.
1. Structure factor
As we have said at the beginning of this section, thedefining characteristic of a quasicrystal is that it has apure point diffraction pattern, i.e. Dirac delta-functionpeaks of the structure factor. The indexing of peaks pro-ceeds as for crystals however it requires a set of D recip-rocal lattice vectors where D is larger than the spatialdimension d . This distinguishes the quasicrystal from aperiodic lattice where the number of reciprocal latticevectors is, of course, equal to d . The peaks of the struc-ture factor of the FC occur for q-vectors given by linearcombinations of two reciprocal lattice vectors.It is rather easy to see that the structure factor of theperfect FC has only Bragg peaks, in the higher dimen-sional representation. The position of the peaks of thestructure factor of the FC can be easily deduced fromthose of the square lattice. The reciprocal space of thesquare lattice is given by (cid:126)G hk = 2 π ( h(cid:126)x + k(cid:126)y ). When pro-jected onto the physical direction, q (inclined at angle θ with respect to q x ), these correspond to the peaks of thequasicrystalline structure factor. Thus, peaks occur forwave-vectors q which are indexed by two integers h and k q = 2 π ( h cos θ + k sin θ )= 2 π √ ω ( h + ωk ) (9)where the last line shows that the peaks lie at positionsgiven by a linear combination of two incommensuratewave vectors.Fig.4a) shows some of the positions obtained by pro-jecting the vertices of the reciprocal space lattice (shownin black) onto the q -axis. In Eq.9 note that h and k cantake all possible integer values, resulting in a dense distri-bution of peaks of the structure factor along the q axis. Inpractice the observable peaks of the structure factor arehowever far fewer, for most of the peaks have intensitieswhich are too small to observe. This occurs because theintensities of the peaks depend on a form factor, namely,the Fourier transform of the selection strip. The strip isdescribed by the function S ( y (cid:48) ) = 1 for 0 ≤ y (cid:48) ≤ W and0 elsewhere. As a result, peak intensities are modulatedby the function |S ( q (cid:48) ) | = 4 sin ( W q (cid:48) / /q (cid:48) (10)where the perpendicular reciprocal space coordinate, q (cid:48) = πa ( − h sin θ + k cos θ ). The S function is akin to the Airyfunction for diffraction through a slit in optics, and itis shown in Fig.4b. The resulting variation of the peakintensities of the FC are shown in Fig.4c. A discussion ofdiffraction patterns of aperiodic systems from a mathe-matical perspective can be found in (Baake and Grimm,2013b). On patterns and their probabilities
The cut-and-projectmethod provides a convenient way to compute the prob-ability of a given pattern to occur in the FC. These prob-abilities are proportional to the length of a corresponding“acceptance zone” in perpendicular space. To illustratethe method, consider the probability of occurrence of a 1-letter pattern: the B tile. To be selected, the correspond-ing vertical bond, of projected length cos θ must lie withinthe window. The acceptance zone is W − cos θ = sin θ ,and the probability of B is sin θ/W . The probabilityof finding an A tile is proportional, by a similar ar-gument, to cos θ/W . The ratio of the probabilities is p ( A ) /p ( B ) = cot θ = τ , as expected.Consider a 2-letter pattern, such as AA. This patterncorresponds, in the 2d square lattice, to 2 consecutivehorizontal bonds sandwiched between two vertical bonds.The probability of the pattern AA is given by (cos θ − sin θ ) /W ≈ . AB and BA haveequal probabilities of ≈ .
2. Conumbering scheme
Sire and Mosseri (Sire and Mosseri, 1990) showed thatthere is an alternative ordering of sites for systems ofco-dimension 1, i.e. having D = d + 1 namely only oneperpendicular dimension. This is of course true in thepresent case of the FC. The so-called conumber of thesite i is given by its position along the y (cid:48) axis. Thisis illustrated in Fig.5 for a periodic approximant of 13sites. Notice that in this numbering scheme, the twonearest neighbors (in real space) of the site c are c =Mod[ c ± F n − , F n ].This numbering orders the sites according to their localenvironments as follows: • sites with conumbers 1 < c < F n − have an B-tileto the left and an A-tile to their right. • sites with F n − + 1 < c < F n − have A-tiles bothon the left and on the right. • the remaining F n − sites have an A-tile to the leftand an B-tile to their right.As we saw in Sec.II.A the FC is a self-similar structure.Its inflation symmetry is coded by a hierarchical struc-ture of the conumbering indices. Each of the three groupsof sites listed above is, in turn, composed of three sub-groups having the same properties but on a bigger length scale. In particular, consider the central site, accordingto the conumber scheme, for an approximant for which n is a multiple of 3. This site always has the same lo-cal environment (A-tiles on both sides) under inflationof n → n −
3, all the way down to n = 1. The conum-bering scheme will be useful later in Sec.V for a compactrepresentation of the spectrum and states of the hoppingmodel. C. The characteristic function method
The characteristic function method, based on the circlemap for generating aperiodic sequences, is a convenientway to generate aperiodic chains. In particular, it pro-vides a useful parametrization of chains in terms of aphason angle φ . The j th letter of the Fibonacci chain isobtained via the characteristic function χ j χ j = sign [cos(2 πj/τ + ϕ ) − cos( π/τ )] (11)where ϕ is a constant. The j th letter is A if χ j = 1 and B if χ j = 0. The Fibonacci approximants C n definedearlier are found by suitable choice of 0 ≤ ϕ < π . When ϕ is varied, phason flips occur, a single flip at a time. Inthis way, one can generate a family of F n + 1 chains oflength F n , corresponding to different values of ϕ (see FigII.C).Alternative forms of the characteristic function can befound in the literature, as for example χ j = [( j + 1) τ ] − [ jτ ] − m >
0, where [ X ] stands for the integer partof X . In this case χ j = 1 stands for A and χ j = 0 for B. III. MODELS AND METHODS. EXACT RESULTS
Much of the extensive literature on electronic proper-ties of the FC is devoted to the study of tight-bindingmodels of the form H = (cid:88) n (cid:15) n c † n c n − ( t n c † n +1 c n + h.c. ) (12)where the (cid:15) n are the site energy at the nth site, while t n is the hopping amplitude between sites n and n + 1.We will assume that these parameters are defined by lo-cal rules, i.e., the values depend on the environment ofeach site. The first papers on this family of Hamiltoni-ans appeared very shortly before the discovery of qua-sicrystals (Kohmoto et al. , 1983b; Ostlund et al. , 1983),and were followed by many other groundbreaking pa-pers in the next few years (Evangelou, 1987; Kalugin et al. , 1986; Kohmoto and Oono, 1984; Kohmoto et al. ,1987; Luck, 1989; Ostlund and Pandit, 1984; Sutherlandand Kohmoto, 1987a; Tang and Kohmoto, 1986). Whenwritten for approximant chains with periodic boundary FIG. 4 a) Plot of the reciprocal lattice of the square lattice, showing the ( q x , q y ) and ( q, q (cid:48) ) axes and the projections of somerepresentative points. b) Plot of the form factor corresponding to the selection strip c) Plot of relative intensities of the firstfew peaks for q > q values are shown in units of (2 π/a ).FIG. 5 The n = 5 ( N = 8) approximant showing the realspace labels of sites (along the chain) as well as their alterna-tive labeling using their conumbers (perpendicular direction). conditions, the Hamiltonian Eq.12 is invariant with re-spect to translations of the selection strip in the 2D space(Sec.II.B). A. Diagonal and off-diagonal Fibonacci models
The Hamiltonian of Eq.12 is termed mixed, as therecan be spatially modulations in both diagonal and off-diagonal terms. However, two important simple casescontain all the essential new physics, namely:1.
Off-diagonal case
Site energies are assumed to beconstant i.e. (cid:15) n = (cid:15) , while the hopping amplitudes t n can can take the value t A or t B , according toa Fibonacci sequence. This is also referred to asthe pure-hopping Fibonacci Hamiltonian since theconstant energy term can be dropped by a suitable FIG. 6 Successive phason flips in the n = 4 approximant( N = 5, the first and last sites are equivalent under a trans-lation) as φ is varied in Eq.11. redefinition of the energy. H = − (cid:88) n t n c † n +1 c n + h.c. (13)Absorbing t B in the definition of the units of energy,leaves as sole parameter the ratio ρ = t A /t B , whichcontrols all the properties of the chain. We willwithout loss of generality, henceforth assume bothamplitudes to be positive, since these signs can bechanged by a gauge transformation.2. Diagonal case
Here the quasiperiodicity is assumedto be present in the diagonal term, while thehopping amplitudes are assumed to be uniform, t n,n +1 = t , for all values of n : H = (cid:88) n (cid:15) n c † n c n − t (cid:88) n c † n +1 c n + h.c. (14)where the onsite potentials (cid:15) n take on two dis-crete values (cid:15) A and (cid:15) B , according to a Fibonaccisequence. As in model 1 above, there is only onenontrivial parameter in this model, and it dependson the energy difference, ε = ( (cid:15) A − (cid:15) B ) /t .It will be instructive to compare the properties of theabove models with those of the quasiperiodic Aubry-Andre-Harper (AAH) model (Aubry and Andr´e, 1980;Gordon et al. , 1997; Harper, 1968). This model is equiv-alent to a tight-binding problem of an electron hoppingin a 2D square lattice subjected to a magnetic field, fora flux per plaquette equal to φ = ωφ , where φ = h/ e is the flux quantum. The resulting quasiperiodic AAHHamiltonian is of the form H = (cid:88) n tc † n +1 c n + h.c. + 2 V cos( 2 πnτ + φ ) c † n c n (15)where the strength of the onsite potential energy 2 V depends on the hopping amplitude along the directiontransverse to the chain and φ is a phase. A well-knownduality transformation takes the AAH model Eq.15 canbe transformed into a Hamiltonian of the same form butwith the exchange t ↔ V . When V = t the model isself-dual. This is the critical AAH model, having manyproperties in common with the Fibonacci model, as de-scribed below. B. Multifractal energy spectra
Spectra can be classified into three types : continu-ous spectra associated with extended states (as in peri-odic solids), pure point spectra associated with localizedstates (as in disordered solids) and singular continuousspectra associated with multifractal states. If one definesthe scaling of the integrated density of states, N ( E ) (de-fined as the fraction of states of energy equal to or lessthan E ) in the vicinity of the energy E by N ( E + ∆ E ) − N ( E ) ∼ ∆ E α (16)then the three cases correspond to α = 1, α = 0 and0 < α <
1. The wavefunctions typically associated withthis last type of fractal spectrum are “critical” states –neither extended nor localized.For the AAH model, the nature of the spectrum de-pends on the parameter
V /t : the spectrum is continuousfor
V /t <
1, pure point for
V /t >
V /t = 1.For the Fibonacci models in Eqs.12 and 13 the energyspectrum is singular continuous as soon as there is ape-riodicity, however small (Bellissard et al. , 1989b; Delyon et al. , 1985; S¨uto, 1989). The situation is analogous tothat of the Anderson model for 1D disordered metals,where the critical value for localization in 1D for disor-der strength is zero.It is noteworthy that the spectra in all three cases arepure spectra – meaning purely singular continuous, or pure point, or absolutely continuous. In general, how-ever, models may have spectra with several different com-ponents, and there may be one or more mobility edgesseparating different regions. This is the case for theAnderson model in 3D for disorder strengths which aresmaller than the critical value, for example. In one di-mension also there can be mobility edges. This occurs ingeneralized Harper models, where the potential energydepends, for example, on two incommensurate wave vec-tors, as discussed by Hiramoto and Kohmoto (Hiramotoand Kohmoto, 1989, 1992), and reconsidered recently in,for example (Das Sarma and Xie, 1988; Ganeshan et al. ,2015; Liu et al. , 2015).Fig.7 shows three typical forms of the densities of statesfor the three models. The left hand figure shows theDOS for the off-diagonal model for the case t A = 0 . t B = 1 and periodic boundary conditions. This modelhas a chiral symmetry – for each solution | ψ (cid:105) of en-ergy E , one has a solution | ψ (cid:48) (cid:105) of energy − E such that (cid:104) n | ψ (cid:48) (cid:105) = ( − ) n (cid:104) n | ψ (cid:105) . The spectrum is therefore symmet-ric around 0, as can be seen in Fig.7. The horizontalaxis represents the dimensionless energy E/t B . The bandstructure is seen to be composed of three main clusters,each of which comprises three sub-clusters, and so on.Each level broadens into a band when periodic bound-ary conditions are assumed. The process of subdivisioninto smaller bands continues as one considers larger andlarger approximant chains.The spectrum for the diagonal model for (cid:15) A (cid:54) = (cid:15) B isvery similar. The spectrum, which is asymmetric, hastwo main clusters. These clusters are in turn composedof three sub-clusters, which trifurcate into three clustersand so on. Part b) of Fig.9 shows the IDOS computednumerically for an approximant chain of N = 144 sites.The parameters were taken to be (cid:15) A = − (cid:15) B = 2 t , withperiodic boundary conditions. For clarity, the origin ofenergy of this curve has been shifted with respect to thecurve shown in a). Finally, for comparison, in part c) ofthe figure, we show the N ( E ) for the AAH model com-puted at criticality V = 1, with τ = τ n in the cosineterm, for a system of N = 144 sites with periodic bound-ary conditions.The figures in Fig.7 illustrate the characteristic fea-tures of multifractal structures, namely local power lawsingularities of the DOS, N ( E +∆ E ) − N ( E ) ∼ ∆ E − α ( E ) .As the system size is increased, band widths shrink, andscale with different exponents α . To fully describe thismultifractal spectrum, one must specify the full set ofexponents α and their densities f ( α ). This can be donenumerically by standard methods of multifractal analy-sis (Halsey et al. , 1986). The method consists of definingthe “partition function”Γ n ( q, τ ) = (cid:88) E (1 /F n ) q (∆ n ( E )) τ (17)where ∆ n ( E ) is taken to be the width of the energy band FIG. 7 Densities of states dN/dE for three models a) the off-diagonal Fibonacci model Eq.13 with t A /t B = 0 . ε = 4 and c) AAH model Eq.15 at criticality. associated to the energy level labelled E . One determines τ as a function of q by requiring that Γ ∼ n →∞ . The function f ( α ) is the Legendre transform of τ ( q ),given by α q = d ( q − d q dq (18)and f ( α q ) = qα q − τ q (19)The function f ( α ) gives the fraction of sites around whichthe DOS scales with the power α . f ( α ) is typically a con-vex curve extending between the extremal values α min and α max . For the periodic crystal when the spectrum iscontinuous, this curve reduces to just two points: α = 1describing the interior of the band, and α = 1 /
2, due tovan Hove singularities at the band edges, as at the bot-tom of the band where dN ( E ) ∼ ( E − E min ) / . In thequasiperiodic case, singularity strengths vary, dependingon the energy. Fig.8, taken from (Rudinger and Piechon,1998) shows f ( α ) values (indicated by crosses) computednumerically for t A /t B = 0 .
2. Exact expressions can beobtained for scaling exponents α ( E ) for two special en-ergies, as we will explain later. C. Gap labeling and topological indices
In one dimensional problems, the IDOS N ( E ), definedas the fraction of states of energy less than E , is alsoequal to the number of changes of sign (nodes) of thewavefunction per unit length. One can thus introducea “wave number” k ( E ) = N ( E ) / E . For the n th periodic approximant,the spectrum consists of F n bands, each corresponding towave vectors k j , such that the IDOS between two bands j and j + 1 has the value j/F n .The gap labeling theorem (Bellissard et al. , 1989a)states that, of the class of models given by Eq.12, theIDOS N ( E ) within the gaps must take values given by N ( E ) = p + q τ (20) FIG. 8 f ( α ) computed numerically (crosses) for the off-diagonal Fibonacci model for ρ = 0 .
2. The solid line cor-responds to an analytic expression obtained using the tracemap (see text) (figure reprinted from R¨udinger and Pi´echon(Rudinger and Piechon, 1998)) where p and q (not to be confused with the multifractalparameter q defined in the last section) are integers. Thepair of indices (p , q) label all of the possible gaps, butthe gap labeling theorem does not specify whether ornot a given gap is actually opened. It suffices to specifyonly one integer, q, as p is then fixed by the condition0 ≤ N ( E ) ≤
1. q is a topological index, it is robust underperturbations. It is a winding number that describes, forexample, the variations of the edge modes in a systemwhich has interfaces (see Sec.III.F).Fig.9 shows the IDOS for each of the three modelsplotted as a function of the dimensionless energy
E/t .The horizontal lines indicate values of IDOS given byEq.20 for gaps corresponding to values of − ≤ q ≤ et al. , 2012;Verbin et al. , 2013). In an elegant experiment using aphotonic waveguide array, the topological equivalence ofthe Fibonacci and Harper models was explicitly shownby (Verbin et al. , 2015).The gap widths for a given label q differ from one0 FIG. 9 Integrated DOS N ( E ) plotted versus the energy E (in dimensionless units) for three topologically equivalent modelsa) off-diagonal model Eq.13 ( ρ = 0 .
7) b) diagonal model Eq.14 ( ε = 4) and c) critical AAH model Eq.15. In all the plots,energies were shifted and normalized such that the total band widths are equal to 1 (system size N = 144, periodic boundaryconditions).FIG. 10 Gap widths plotted on a log scale versus the topo-logical index q for the n = 16 approximant in the off-diagonalmodel. Red data points stand for transient gaps (see Sec.V)(figure reprinted from (Mac´e et al. , 2017)) model to another, and they tend to decrease with q, al-though not in monotonic fashion. This is shown in Fig.10which shows the result of gap widths plotted against q,for the n = 16 approximant. These were computed forthe pure-hopping model using an approximate renormal-ization group, as outlined in Sec.IV. The black (red) datapoints correspond to stable (transient) gaps respectively(these are gaps which tend to finite (zero) width as n tends to infinity). D. Trace map method
The trace map analysis is a powerful technique whichhas led to many results for quasiperiodic Hamiltonians.The starting point is the definition of transfer matricesrelating the wavefunction amplitudes on the n + 1th siteto the amplitudes on sites n and n −
1. From Eq.12 onehas the tight-binding equations( E − (cid:15) n ) ψ n + t n ψ n +1 + + t n − ψ n − = 0 (21) This relation can be re-expressed in terms of a 2 × (cid:20) ψ n +1 ψ n (cid:21) = (cid:20) ( E − (cid:15) n ) t n − t n − t n (cid:21) (cid:20) ψ n ψ n − (cid:21) = (cid:32) n (cid:89) k =2 T k,k +1 (cid:33) (cid:20) ψ ψ (cid:21) = M n (cid:20) ψ ψ (cid:21) (22)where we have introduced the local transfer matrix T k,k +1 , and the global transfer matrix M n which is aproduct of n − M n = .....T , T , .1. Diagonal model.
The local transfer matrix dependson the onsite energy (cid:15) n and the amplitudes for hop-ping onto the sites to the left and to the right ofsite n . For the diagonal model, there are only twodifferent possible transfer matrices, namely T A = (cid:20) ( E − (cid:15) A ) t −
11 0 (cid:21) ,T B = (cid:20) ( E − (cid:15) B ) t −
11 0 (cid:21) (23)Let us now consider the approximant chain oflength N = F n , and let x n = Tr M n be the half-trace of the transfer matrix. For an energy E to cor-respond to an allowed (normalizable) wavefunction,the half-trace of M N must satisfy the condition | x n | ≤
1. Thanks to the concatenation propertyof chains already mentioned in Sec.II.A, namely C n +1 = C n ⊕ C n − , the global transfer matrices forsuccessive approximants satisfy (Kohmoto et al. ,11983b; Ostlund et al. , 1983) M n +1 = M n − M n (24)Given Eq.24, it can be shown that the half-tracessatisfy a three term recursion relation (Kohmoto et al. , 1983b): x n +1 = 2 x n x n − − x n − (25)with the initial conditions x − = 1 , (diagonal model) x = ( E − (cid:15) B )2 , x = ( E − (cid:15) A )2 (26)To determine which energies belong in the spec-trum one computes the iterates x n Using Eq.25and checks to see whether they remain boundedand within the interval ( − , Off-diagonal model.
A similar set of relations ob-tains in the case of the off-diagonal model. Onestarts by introducing three different transfer ma-trices, T AA , T AB and T BA , for the three bond con-figurations AA,AB and BA which are possible inthe FC. These matrices are defined by T AA = (cid:20) E/t −
11 0 (cid:21) ,T AB = (cid:20) E/t − t B /t A (cid:21) ,T BA = (cid:20) E/t − t A /t B (cid:21) (27)The problem of writing the global transfer matrixcan be simplified (Kohmoto et al. , 1987), as thehopping on Fibonacci chains can be described withjust two matrices T AA and T AB T BA . If we renamethese transfer matrices as T B and T A respectively,the global transfer matrices for approximant chainscan be written exactly as in the diagonal case.The recursion relations Eq.24 and Eq.25 thereforehold also for the off-diagonal model. The initialconditions for the half traces in the off-diagonalcase are x − = 12 (cid:18) t B t A + t A t B (cid:19) , (off-diagonal model) x = E t B , x = E t A (28)The recursion relation for the traces Eq.25 constitute adynamical system in a three-dimensional space. Defining the variables x = x n − , y = x n and z = x n +1 , Eq.25maps a given point as follows: x → x (cid:48) = yy → y (cid:48) = zz → z (cid:48) = 2 yz − x (29)One of the invariants of the dynamical system Eq.29 isthe quantity (Kohmoto et al. , 1983b; Kohmoto and Oono,1984; Kohmoto et al. , 1987): I = x + y + z − xyz −
1= 14 ( (cid:15) A − (cid:15) B ) (diagonal model)= 14 (cid:18) t A t B − t B t A (cid:19) (off-diagonal model) (30)where the last two equations were written using the ini-tial conditions for the diagonal and off-diagonal modelrespectively. For a proof of the invariance of I , a so-calledFricke character, see (Baake et al. , 1993). Under the dy-namical map, points move on the surface I = const , inthe three dimensional space. Details of the form of thesesurfaces for different values of the parameters and of dif-ferent kinds of orbits are given in (Kalugin et al. , 1986).Orbits which escape to infinity, such that lim n →∞ x n isinfinite, correspond to energies which are not in the spec-trum. Numerically, this is found to be the case of almostall energies, consistent with the fact that the spectrumhas a Lebesgue measure of zero. Periodic orbits withlim n →∞ x n ≤ Itcan be shown, by tracing orbits for successive periodicapproximants, that the band structure has a self-similarstructure – is a Cantor set.Kohmoto, Kadanoff and Tang obtained exact resultsfor two cases where the trace map leads to a periodic orbit(Kohmoto and Banavar, 1986; Kohmoto and Oono, 1984;Kohmoto et al. , 1987; Ostlund et al. , 1983). By consid-ering the linearized map around these special points ofthe spectrum, they found the scaling exponents for thecorresponding bands in terms of the “escape rates” of thedynamical map. The band widths are given by ∆ ∼ ω (cid:15) ,where λ is the minimal eigenvalue of the linearized map(see below).1. Solution for the band center. The trace map has aperiodic orbit consisting of the six-cycle (0 , , a ) → ( − a, , → (0 , − a, → (0 , , − a ) → ( a, , → (0 , a, → (0 , , a ), where a = √ I + 1. The scalingexponent for this band can be expressed in terms ofthe eigenvalue (cid:15) of the linearized equation around There are in principle two other more “exotic” possibilities iii)aperiodic and bounded orbits and iv) recurrent orbits where thepoint returns to the allowed region and for which lim n →∞ x n ≤ α ctr (Kohmoto and Oono, 1984) is α ctr = ln τ / ln (cid:15) (cid:15) = [ (cid:112) I ) + 2(1 + I )] (31)2. Solution for band edges. These correspond to two-cycles of the trace map, ( a, b, b ) → ( b, a, a ) → ( a, b, b ) where a = J + √ J − J and b = J −√ J − J where J = [3 + √
25 + 16 I ]. The scalingexponent for these bands is expressed in terms of (cid:15) , the eigenvalue of the linearized map, as follows α edge = ln τ / ln (cid:15) (cid:15) = [8 J − (cid:112) (8 J − − / et al. , 1987) that the α ctr and α edge values correspond to the extremal values,namely, α min and α max . In fact, however, R¨udinger andPi´echon (Rudinger and Piechon, 1998) showed that thisis not always the case, by analyzing the trace map inthe vicinity of a 4-cycle that governs scaling for IDOSvalue N ( E ) = and . Their analysis showed thatthe maximal value of α occurs at these points when ρ issmall. Another conjecture concerned the possibility thatthe spectrum becomes monofractal when ρ = ρ c (Zhong et al. , 1995). This conjecture was based on the fact that α ctr is equal to α edge for the hopping ratio ρ c ≈ . ρ = ρ c , α / isdifferent from – larger than – the other two exponents.Thus the spectrum is not a monofractal. An approximateanalytical expression for f ( α ) derived in (Rudinger andPiechon, 1998) is shown in Fig.8, along with the numer-ical data for ρ = 0 . E = 0 wave-function in the hopping model is found by the trace mapcalculation to be | ln ρ | / ln τ , in good agreement with nu-merical calculations of the wavefunction (Kohmoto andBanavar, 1986).Discussions of trace maps and their dynamical proper-ties can be found in the reviews by (Baake et al. , 1993;Damanik et al. , 2015). Generalizations of the trace mapfor mixed Hamiltonians with two (and more) parameterswill be briefly discussed in Sec.IX. As one might expect,these generalized models have a larger parameter space,with more possibilities for the spectra. They can admit,for example, extended states for special energies. E. Log-periodic oscillations
Discrete scale invariances can lead to log-periodic oscil-lations in thermodynamic properties according to a gen-eral argument (Derrida et al. , 1984). Within the for-
FIG. 11 (top) Log-log plot showing the numerically com-puted IDOS versus energy for the hopping model in the vicin-ity of E = E min (for ρ = 0 . N ), showing the main period ofthe oscillation and some smaller periods. malism of renormalization group, Gluzmann and Sor-nette (Gluzman and Sornette, 2002) considered an ob-servable f ( x ) in a system close to criticality, where f ( x ) = µ − f ( γx ) under a renormalization transforma-tion. They showed that f ( x ) has a power law scaling“decorated” by a log-periodic function f ( x ) = x m P ( x ).The power is given by m = ln µ/ ln γ and the period of theoscillations is log γ . This is indeed what one observes forthe IDOS N ( E ) in the Fibonacci model. Fig.11a) showsin a log-log plot the IDOS versus energy in the bottomof the band. The points are obtained by numerical di-agonalization, and the straight dashed line indicates theaverage IDOS. Fig.11b) shows the fluctuations of ln(N)around the average value. The period of the oscillations,close to 1, corresponds to the inflation factor which inthis case is γ = τ for the side bands (for a description ofthe renormalization transformation of the Hamiltoniansee Sec.IV). The main period and some smaller ones, in-dicating a fractal structure can be seen. F. The wavefunction for E = 0 This subsection discusses an exact solution for one ofthe wavefunctions of the hopping model, whose multi-fractal properties can be computed as a function of ρ .The wavefunction at the band center E = 0 (shown in3 FIG. 12 Center state wavefunction | ψ i ( E = 0) | plotted ver-sus site index i (numerically computed for n = 12 chain with ρ = 0 . Fig.12 for a finite approximant) can be determined ex-actly, by a recursive construction. It turns out to be aparticularly simple form of the critical states proposedby Kalugin and Katz (KK) for the ground state of tilings(Kalugin and Katz, 2014). These authors argued that fora family of quasiperiodic Hamiltonians that includes allthe standard ones, the ground state | ψ (cid:105) can be writtenas a product of two factors. The amplitude on site i isgiven by (cid:104) i | ψ (cid:105) = C ( i ) e κh ( i ) (33)where κ is a (real) constant. The prefactors C ( i ) dependon the local configuration of the atoms around site i .Long range correlations between sites are given by theexponent h , called the height field , which is the integralof a quasiperiodic function. The E = 0 wave functionof the Fibonacci chain has the form of Eq.33 with theparticular choice C ( i ) = ± i belongs. This solution for the 1D chainis a relatively tractable case study to illustrate some ofthe properties of the KK eigenstates which are of coursemore complex in higher dimensional quasicrystals.From the Hamiltonian Eq.13, the following relation ob-tains for the E = 0 wave function amplitudes on i and i + 2 t i +1 ψ i +2 + t i ψ i = 0 (34)There are two independent E = 0 solutions, one for eachsublattice. The two being equivalent in the limit of theinfinite FC let us henceforth consider the even sub-latticesolution, for sites i = 2 m . There are three possibilities Note that the expression 33 can be considered as a generaliza-tion of the usual Bloch form for the wavefunctions in a periodiclattice, which can be written as a product of a periodic function u n ( x ) (in 1D) and an exponential e ih ( x ) where h ( x ) = kx is theintegral of the Bloch wave vector k . FIG. 13 Definition of the local arrow function and its integral,the height function (figure reprinted from (Mace et al. , 2017)) for the bond configurations between sites i and i + 1,namely, AA, AB or BA . The three cases are given bythe relation ψ m +1 = − ρ A ( m ) ψ m (35)where ρ = t A /t B and the A (for arrow ) function is de-fined locally according to the configuration of the bondsbetween the two sites A ( m ) = +1 ( AB ) A ( m ) = − BA ) A ( m ) = 0 ( AA ) (36)Fig.13 shows the arrow function for the even sites of asmall chain segment. The figure shows the arrows corre-sponding to the three bond configurations which are pos-sible using the conventions: → (rightarrow) for the bondsequence AB , ← (left arrow) for the bond sequence BA and no arrow for the bond sequence AA .Repeating the recursion relation Eq.35, one obtains ψ m = ( − m ρ h ( m ) ψ = ( − m e κh ( m ) (37)where h ( m ) = (cid:80) m A ( j ), κ = ln ρ and ψ was set equalto 1. This expression is of the KK form, Eq.33, with aconstant prefactor on all sites. The function h in the ex-ponent is an integral of a quasiperiodic function, A ( m ).Fig.13 shows the height function for the first few evensites. As the length of the chain gets larger, the heightfunction fluctuates more and more. Fig.14 shows theheight function calculated for a long segment of the Fi-bonacci chain. The properties of the wave function canbe determined when the distribution of heights is known.One can show by explicit calculation, that the wave func-tion is multifractal, and express all of its generalized di-mensions D φq in terms of ρ . This can be done by in-troducing inflation matrices to relate the heights in the n th and the n − n limit, theheight distribution P ( h ) satisfies a diffusion equation asa function of t (the number of inflations), as follows P ( t ) ( h ) ∼ √ πDt exp (cid:18) − h Dt (cid:19) (38)4 FIG. 14 The height function plotted for a long chain, showingits fluctuations at large distances (figure reprinted from (Mace et al. , 2017)) where the “diffusion coefficient” is given by D = √ .We refer the reader to (Mace et al. , 2017) for details.Since P ( h ) is a symmetric distribution around h = 0,the resulting form of the wave function has a symme-try between peaks and valleys, as can be seen in Fig.12.The typical value of h in a chain of N = τ n sites isgiven by the standard deviation of the gaussian aftera “time” t = n , h typ ∼ √ Dn . This leads to a typ-ical value of ψ which falls off with the chain length, N = τ n , as ψ (2 N ) ∼ e − cst √ D ln N . The spatial de-cay of the wavefunction is thus faster than power lawbut slower than exponential. In contrast, note that,for the randomly disordered off-diagonal model, a sim-ilar argument gives the wavefunction at E = 0 to be ψ (2 N ) ∼ e − cst √ DN (Economou and Soukoulis, 1981;Inui et al. , 1994; Theodorou and Cohen, 1976). This isa stretched exponential function, decaying much fasterthan the E = 0 wavefunction of the Fibonacci chain.The fractal dimensions of ψ can be exactly computed.These quantities are deduced from the scaling of the mo-ments of the wavefunctions, which are defined as follows.Let the q -weight of the wavefunction ψ be defined by: χ q ( ψ, R ) = (cid:80) i ∈R | ψ i | q (cid:0)(cid:80) i ∈R | ψ i | (cid:1) q (39)where the sums run over all sites in a given region R .The q -weight is a measure of the fraction of the presenceprobability contained inside region R .Consider a sequence of regions R n whose radius growsto infinity as n → ∞ . The q th fractal dimension, D ψq ( ψ ),is the scaling of the q -weight with the volume of the re-gion: D ψq ( ψ ) = lim n →∞ − q − χ q ( ψ, R t )log Ω( R n ) (40)where Ω is the number of sites inside region R . As wehave already seen for the DOS in Sec.V.A, one can thencompute the Legendre transform of the fractal dimen- FIG. 15 f E =0 ( α ) spectrum as given by Eq.40 for the E = 0wave function for different values of the hopping ratio ρ . Notethe symmetry of f E =0 around its maximum value. (Repro-duced with permission from (Mace et al. , 2017)). sions, f E =0 ( α ). This function gives the fraction of sitesfor which the wavefunction scales with the power α . f E =0 ( α ) curves obtained by the exact calculation fordifferent values of ρ are plotted in Fig.15. As the figureshows, in each case the α values lie within a finite inter-val, indicating that ψ is multifractal. As ρ approaches1, the support of the function shrinks to a single point, α = 1, corresponding to the extended state. Anotherpoint to note: thanks to the symmetry of the heights dis-tribution, the function f E =0 ( α ) is symmetric around itsmaximum, as was already observed in numerical studiesof this wave function (Evangelou, 1987; Fujiwara et al. ,1989). In other words, as the system size increases, theminima/maxima of ψ scale to zero/infinity in the samemanner.The transmission coefficient for this E = 0 wave func-tion can be calculated exactly. It is interesting to notethat this quantity (which will be discussed in sectionSec.VII) has the exact value of 1 (i.e. transmission is perf ect ), for certain sites along the chain, for arbitrarilylarge distances. G. Chern numbers. Bulk-edge correspondence
The close connection between topological phases andquasicrystals which are described in a higher dimensionalspace has been pointed out by many authors (Huangand Liu, 2018, 2019; Kraus et al. , 2012; Verbin et al. ,2013). These show that, for 1D quasicrystals, there are2D Chern numbers and that there are topologically pro-tected boundary states similar to those in a 2D quantumHall system.It is expected that edge modes should be present in thequasicrystal, analogous to those in the AAH model. Justas changing the arbitrary phase φ in the AAH modelEq.15 leads to tuning the edge mode energy, one cantune edge modes in FC approximants by varying the ar-bitrary phase ϕ in Eq.11. This is seen in Fig.16, which5shows for an 89-site chain the energy levels as a functionof the parameter 0 ≤ ϕ ≤ π (the figure was renderedsymmetric with respect to π by shifting the angle by ϕ = − ωπ ( N + 1)). It can be seen that the levels remainflat as φ is varied until a sudden phason flip occurs some-where along the chain. There are, in all, N such flips inthe interval. The label of the gap q gives the number ofgap-crossings of the states, which are seen most clearlyin the main gaps of the spectrum in Fig.16.Experimental studies of the hopping Hamiltonian us-ing a polaritonic cavity modes to detect eigenmodes andtheir energies can be performed, as shown in (Baboux et al. , 2017; Tanese et al. , 2014). In the experiments,nearest neighbor cavities were spaced so as to be linkedby a strong or a weak coupling, following the Fibonaccisequence. Chains of given length N , one for each φ -value,were fabricated in a Fabry-Perot geometry, resulting inedge modes located at the central mirror symmetric po-sition. For a given gap of label q, the corresponding edgemode was observed to cross the gap q times, confirmingthat this is indeed a winding number. The sign of q de-termines the sense of the gap crossing (from the upper tolower edge, or vice versa). FIG. 16 Plot of the energy spectrum versus φ for energy levelsof the off-diagonal model in an open 89-site chain ( ρ = 0 .
7) .The number of gap-crossings of states is most easily countedfor the largest gaps ( | q | ≤ One can use the winding property of states to “pump”charge adiabatically across the chain, as was done in ex-periments in photonic quasicrystals by Kraus et al (Kraus et al. , 2012; Verbin et al. , 2013). The “jumping” of statesfrom one edge to the other as the φ parameter is changed,is shown in Fig.17 for the q = 2 gap states. This topolog-ical pumping can have consequences for physical proper-ties. We will see later that the winding number q appearsin quantities such as the induced superconducting orderparameter as discussed in Sec.VIII.R¨ontgen et al studied real space bond configurationsfor the appearance of edge modes in (R¨ontgen et al. ,2019). They showed that the edge modes are linked tolocal “resonators” – this term used by them to denoteclusters of bonds which are symmetric under reflection. FIG. 17 Plot of the states on either side of the q = 2 gapfor six values of φ between 0.2 π and 0.8 π showing the upper(orange) and lower (blue) band edge states. One sees thesestates evolve from extended in-band states, becoming local-ized edge states, the edge states exchanging positions, andfinally becoming in-band states again The localization length and energy of the edge state de-pend on the resonator “size” (i.e. the distance out towhich they possess reflection symmetry).Chern numbers have, interestingly, also been observedby a light diffraction experiment (Dareau et al. , 2017).In the experiment, a digital micromirror device was usedto realize a set of approximant chains of fixed length anddifferent values of φ in Eq.11. The behavior of the diffrac-tion peak at different wave vectors k was shown to dependon the associated topological number q. IV. APPROXIMATE METHODSA. Perturbative methods
Many different kinds of perturbative calculations havebeen done to study electronic properties of the FC, bothin real and reciprocal space.Luck (Luck, 1989) carried out a perturbation expan-sion for the diagonal model in terms of the Fourier com-ponents of the potential. Taking the onsite energies tobe (cid:15) A = V and (cid:15) B = − V by a suitable shift of the origin,one can compute the Fourier transform of the potential, V G N ( k ), and the structure factor, S N ( k ) = V | G N ( k ) | .In the thermodynamic limit the structure factor can beshown to have power law singularities at a dense set ofreciprocal space vectors k (as seen in Sec.II.B). In thevicinity of each of the peaks, one has S ( k ) dk ∼ | k − k | α (41)6For a general potential the singularity strength α canvary depending on the peak, whereas for the quasicrystal α = 1 for all peaks according to the arguments we pre-sented in Sec.II. Luck showed for the general case thatthe width of the gap opened at the unperturbed energy (cid:15) ( K = k /
2) is related to α , the singularity at k = k .Specifically, for the gap where the IDOS N ( E ) = k / π ,the gap width is given by∆ ∼ V β (42)where β = 2 / (2 − α ). In the case of the quasicrystal, β = 1 for all the gaps. This analysis shows that, for weakquasiperiodic potentials, the plateaux of the IDOS arerelated in a natural way to the module of wave vectors.It should be noted that this perturbation theory does notconverge, even for arbitrarily weak potentials, as pointedout by Kalugin et al (Kalugin et al. , 1986), because ofthe nature of the Fourier module of the quasicrystal –consisting of a dense distribution of peaks. Nevertheless,the indexing of gaps using this method is robust. Thegap labeling theorem (Bellissard et al. , 1989a) providesthe rigorous justification that the indexing continues tohold for arbitrarily strong potentials.Other perturbative approaches have been proposed.for the off-diagonal model, a real space perturbationtheory considers the effects of weak quasiperiodicity(Rudinger and Sire, 1996; Sire and Mosseri, 1989). An-other approach, for the diagonal model, which starts froma strong atomic limit, and treats hopping as a perturba-tion was introduced by Barache and Luck (Barache andLuck, 1994). B. Approximate renormalization group
This subsection describes the main ideas behind a per-turbative real space RG method due to Niu and Nori (Niuand Nori, 1986, 1990) and Kalugin, Kitaev and Levitov(Kalugin et al. , 1986; Levitov, 1989). This approach hasbeen extremely fruitful for describing a great number ofstatic and dynamic properties of the FC. In this sectionwe will describe the basic notions of this RG for the off-diagonal and the diagonal models.
1. RG for hopping model (1)
We begin with details of the RG for the off-diagonalmodel Eq.13, where t B > t A . The hopping ratio ρ = t A /t B lies in the range 0 ≤ ρ ≤
1. Recall that onecan assume that both t A and t B are positive for, if not,the solutions can be found from our model by a suitablemapping - or “local gauge transformation” of the wave-functions. The goal of RG is to obtain a description ofthe spectrum and states perturbatively in ρ .For ρ = 0, the chain breaks up into disconnectedgroups of sites which can be classified as follows: atom sites – the sites sandwiched between A bondsmolecule sites – pairs of sites linked by a B bondThe spectrum, in this limit, consists of only three dis-crete degenerate levels: the E = 0 level of the atoms,and E = ± t B for the molecular bonding/antibondinglevels. For a chain of N = F n sites, the degeneracy ofthe E = 0 level is given by the number of atoms, F n − .The degeneracy of the E = ± t B levels is given by thenumber of molecules, F n − , in the chain. For small non-zero ρ (cid:28)
1, these three levels split into three clustersof levels, the molecular bonding (+m) and antibonding(-m) bands and the atomic (0) cluster. The separationsbetween these three clusters for small ρ is roughly t B . Inperturbation theory, the three clusters do not mix, andcan be treated as three independent systems for the cal-culation of the effective Hamiltonians. Molecular RG
Consider the Hamiltonian H for a Fi-bonacci chain of length F n and consider the lowest molec-ular bonding level located at the energy − t B and havinga wavefunction which is non-zero on the molecule sites. Itis easy to check that the chain formed by molecules is pre-cisely the n − et al. , 1986;Niu and Nori, 1986, 1990) that, upto an overall con-stant shift, the new effective Hamiltonian H (cid:48) is againa Fibonacci hopping Hamiltonian, with the renormal-ized hopping amplitudes t (cid:48) A and t (cid:48) B . The old chain andthe new chain after decimation of atoms are indicatedin Fig.18, along with the two new hopping amplitudes.The renormalized hopping amplitudes are given to lowestorder in ρ by t (cid:48) A = zt A t (cid:48) B = zt B (43)where z = ρ/
2. Note that the hopping ratio is unchangedto lowest order under RG, as the new weak and stronghopping amplitudes satisfy t (cid:48) A /t (cid:48) B = ρ (cid:48) = ρ . To sum-marize, the effective Hamiltonian for the bonding set oflevels is, up to a global shift, that of the n − E = − t B is split into 3 levels, which canbe labeled − +, − −− , separated by gaps of width t (cid:48) B .A similar analysis shows that the effective Hamiltonianfor the antibonding levels “+” is an identical FC of F n − sites with hopping amplitudes given by t (cid:48) A = zt A (weak)and t (cid:48) B = − zt B (strong). The original level located at E = − t B is split into 3 levels, labeled ++ +0 and + − . Atomic RG
Consider the central band located around E = 0. The levels around E = 0 correspond to wavefunc-tions which are largest on the atom sites, of which thereare F n − . It is easy to check that the new chain formedby these atom sites is nothing but the n − H (cid:48) depends on two new renor-malized hopping amplitudes. The new strong and weak7 FIG. 18 Illustration of the molecular deflation rule: the fifthapproximant is transformed to the third.FIG. 19 Illustration of the atomic deflation rule: here thefifth approximant is transformed to the second. bonds, t (cid:48)(cid:48) A and t (cid:48)(cid:48) B , are given to lowest order in ρ by t (cid:48)(cid:48) A = zt A t (cid:48)(cid:48) B = zt B (44)where z = ρ . As for the molecular RG, the hoppingratio is preserved, since the new weak and strong hop-ping amplitudes satisfy ρ (cid:48)(cid:48) = t (cid:48)(cid:48) A /t (cid:48)(cid:48) B = ρ . To summarize,the effective Hamiltonian for the atom set of levels is theHamiltonian of a chain of F n − sites and with renormal-ized hoppings. The original level is therefore split into 3levels, labeled 0+, 00 and 0 − . These are separated bygaps of width t (cid:48)(cid:48) B .This process can be repeated until one reaches thethree first chains. The result for the clustering structureis a succession of trifurcations as illustrated in Fig.20a).One can reverse the process, alternatively, and track eachband as it splits into three sub-bands when n increasesby 2 or by 3. Doing this, one sees that each of the F n levels of the spectrum follows a unique path under suc-cessive renormalizations (termed renormalization path)of the form { c c c .... } where c i can take the three val-ues 0 , ±
1. This trifurcation scheme also holds for thecritical AAH model (Hiramoto and Kohmoto, 1989).
2. RG for hopping model (2)
Consider the FC in which the strong bond correspondsto A-bonds. Again, one can assume, as for the precedingmodel that the amplitudes are positive without loss ofgenerality, with t A (cid:28) t B . The perturbation theory isnow carried out in powers of t B /t A . Inspection of theFC shows that in the limit t B = 0 the chain breaks upinto diatomic and triatomic molecules. These moleculesgive rise to five energy levels E/t A = ±√ , ± ,
0. Itcan be shown that the new effective Hamiltonians withineach of the five bands is a Fibonacci hopping Hamiltonian
FIG. 20 Cluster and subclustering structures for the threedifferent models described in Sec.IV (figure reprinted from(Niu and Nori, 1990)) (Niu and Nori, 1990) of the type discussed above (modelA). Thus each of the five levels trifurcates, and continuethereafter to trifurcate under successive RG steps. Thisis illustrated in Fig.20b.
3. RG for diagonal model
For the diagonal model Eq.14, a perturbation theoryin t/ ( (cid:15) A − (cid:15) B ) shows once again the recursive structure ofthe energy spectrum. For t = 0, one has two isolated lev-els, E = (cid:15) A (degeneracy F n − ) and E = (cid:15) B (degeneracy F n − ). For small non-zero t , one can compute the neweffective Hamiltonians in perturbation theory. It is foundthat these are given, once again, by two hopping param-eters, one strong and one weak. Thus after one RG step,we are led back to the Fibonacci hopping model, leadingto a splitting into three levels, and thereafter, with eachsuccessive RG step, trifurcations. This is illustrated inFig.20c. V. PERTURBATIVE RG THEORY OF MULTIFRACTALSPECTRA AND STATES
In this section we review the predictions of perturba-tive RG described above for multifractal properties ofthe hopping model. Sec.V.A shows how RG methods in-troduced in the last section can be applied to the caseof the hopping model to compute spectral properties asdone in (Pi´echon et al. , 1995; Zheng, 1987). Sec.V.Cand wave functions in (Mac´e et al. , 2016; Pi´echon, 1996;Thiem, 2015a). In Sec.VI we describe the use of RG tostudy the quantum dynamics of wave packets, following8the approach used in (Pi´echon, 1996; Thiem, 2015a).
A. Multifractality of the energy spectrum
The RG method described in the preceding sectionshowed that the spectrum of a chain of number of sitesequal to F n can be mapped, after one RG step, to thespectra of two shorter chains F n − and F n − . We nowapply this to obtain quantitative information on the spec-tral properties, following Zheng (Zheng, 1987) and Be-nakli et al (Pi´echon et al. , 1995). This can be formallyexpressed by the relation H n = ( zH n − − t s ) (cid:124) (cid:123)(cid:122) (cid:125) bonding levels ⊕ ( zH n − ) (cid:124) (cid:123)(cid:122) (cid:125) atomic levels ⊕ ( zH n − + t s ) (cid:124) (cid:123)(cid:122) (cid:125) antibonding levels + O ( ρ )(45)to lowest order in ρ . Thus, given the first three spectra W ( n ) corresponding to n = 0 , ,
2, one can construct allthe n > n = 0 and n = 1 forexample, the approximant chains have only one hoppingamplitude, t B or t A . Applying periodic boundary condi-tions the spectrum for n = 0 is the band − t B < E < t B ,with a band width of ∆ (0) = 2 t B . The spectrum of the n = 1 chain is a narrower band − t A < E < t A . The n = 2 chain is an alternating sequence of t A and t B , thusthe spectrum has two bands separated by a gap as shownin Fig.21. One can now proceed to construct the spec-trum for n = 3. W (3) is composed of two (bonding andantibonding) lateral bands and one central (atom) band.The side bands are simply W (1) multiplied by the factor z and translated in energy by ± t B . The central band is W (0) multiplied by the factor z . This procedure can beused to construct all the spectra shown in the figure.In each RG step, the band widths are reduced by thefactors z or z . The resulting band width ∆ of a givenlevel in the n th approximant depends on the sequenceof RG steps that were taken. A simple example is thetwo bands at the top and bottom edges of the spectrum,which remain always molecular, throughout successiveRG transformations. They have the RG paths 111 ... and111 ... having ∼ n/ α edge , defined through∆ edge ∼ F − /αn = ω n/α are given by∆ edge = z n t A α edge = log( ω ) / log( z ) (46)to lowest order in ρ . The second simple case concernsthe atomic level at E = 0, which has a RG path of 000 ... having n/ n a multiple of 3). Thus itsband width and its scaling exponent α ctr are given by∆ ctr = z n t A α ctr = log( ω ) / log( z ) (47) FIG. 21 Schematic view of the recursive construction of spec-tra using RG. The first three spectra for n = 1 , n = 3 (two bands) are shown at the top of the figure.The n th spectrum is obtained by the union of the n − n − z and z ,and shifted as described in Eq.45. The labels g and g ∗ referto transient and stable gaps, see text (figure reprinted from(Rudinger and Piechon, 1998)) Comparison with the exact results of Kohmoto et al forthese two exponents obtained using the trace map, onesees that these values of α represent the first terms of anexpansion in ρ . Other levels have a mixture of atom andmolecular RG, so that the scaling is given in general by z n m z n a where n m ( a ) is the number of molecular(atomic)RG steps in its RG path. These numbers are not indepen-dent as they must satisfy the condition n ∼ n m + 3 n a .It is useful now to introduce the variable x = n m /n asa measure of the degree to which a given RG path hasmolecular character. For very long chains, x is a contin-uous variable in the interval [0 , / x = 0 corresponds to the energy in the middle of thespectrum. The maximum value x = 1 / x versus the level index, for the n = 12 approximant. Notethat, in general, many different levels can share a givenvalue of x (cid:54) = 0, whereas the value x = 0 corresponds to asingle state which occurs only in every third chain ( n amultiple of 3).The exponent α can be computed for a given value of9 FIG. 22 x ( E a ) plotted versus the level index a for each ofthe energy levels E a of the n = 12 chain. Lines are drawn toguide the eye. given x , and is given by α = − log F n log ∆= n log ω ( n m log z + n a log z )= log ωx log z/z / + log z / (48)Let the number of levels scaling with a power α be N ∼ F f ( α ) n , thus defining the function f ( α ). This is justthe number of levels corresponding to a given value of n a and n m , which is 2 n m n m ! / ( n m + n a )!. With the use ofStirling’s approximation this leads to N ( x ) ∼ F g ( x ) n with g ( x ) = 1log ω ( x x + 1 + x x )+ 1 − x − x )) (49)Eqs.48 and 49 determine the function f ( α ), and describethe multifractal scaling of the spectrum. The exponentfor any band can be computed if its RG path is specified.The expressions for α ctr and α edge are the leading termsin an expansion in small ρ of the exact formulae Eqs.31and 32 obtained by Kohmoto et al. We note finally thataccording to the results in Eqs.48 and 49 there is a spe-cial point when ρ = 1 /
8. At this point z = z / , and thespectrum becomes mono-fractal, because all bands scalein the same way according to our perturbation theory.However, in fact, as we have mentioned in the discus-sion at the end of Sec.III.D, exact results show that thespectrum is not a monofractal for any ρ .The scaling of the total band width, B ( n ) = (cid:80) F n j ∆ ( n ) j with the system size can now be determined. From Eq.45,one can deduce a recursion equation relating the totalband width of the n th chain to those of the n − n − B ( n ) = 2 zB ( n − + zB ( n − (50) Defining the exponent b by B ( n ) ∼ F − bn for large n , fromEq.50 one sees that b must satisfy the equation1 = 2 zω − b + zω − b (51)Recursive relations can be written likewise for all of themoments of the density of states (i.e. the inverse bandwidth ∆ − j ). Recall that the generalized dimensions forthe spectrum are the exponents corresponding to the q thmoment of the DOS. Following the thermodynamical for-malism introduced before, one defines the partition func-tion Γ ( n ) ( q, τ ) = F n − q (cid:80) F n j ∆ − τj , where τ q = D q ( q − ( n ) ( q, τ ) = 2 ω qn z τ Γ ( n − ( q, τ ) + ω qn z τ Γ ( n − ( q, τ ) (52)For each q , the corresponding τ value is obtained by re-quiring that Γ be stationary. This results in the condition1 = 2 ω q z (1 − q ) D q + ω q z (1 − q ) D q (53) Relations for generalized dimensions
The Hausdorff di-mension D F of the spectrum is given by D , which sat-isfies the equation 2 z D F + z D F = 1 (54)Although derived in the limit of small values ρ (cid:28)
1, thisrelation nevertheless gives a rather good value even forrelatively large values of ρ – one obtains D F = 0 .
76 for ρ = 0 . The information dimension, D , enters in aninequality for the diffusion exponent, discussed in Sec.VI.The exponent D is also of special interest, in particular,for dynamics, as will be shown later in Sec.VI. One sees,from Eqs.53 and 51 that the band width exponent b isrelated to the D q via D δ = 1 / (1 + δ ). B. Gaps, stable gaps and topological numbers
In the construction of spectra with the RG recursionscheme, it becomes apparent that two kinds of gaps ap-pear in the spectra of approximant chains: there are tran-sient gaps and stable gaps. To understand these notionsconsider the spectra of the first few approximants shownin Fig.21: one sees that the spectrum for n = 2 has agap labeled g which disappears for n = 3 and n = 4,reappearing as a smaller gap for n = 5, and going to zeroas n tends to infinity. This is an example of a transientgap. Stable gaps are descendants of the gaps marked g ∗ whose widths remain finite. Writing recursion relations Strictly speaking, this approach is valid only for strong quasiperi-odic modulations. However, these results remain pertinent evenfor moderate to weak quasiperiodicity, and calculations on finitechains show that, when the gaps are opened, they persist for all ρ , closing only in the periodic case. P ( g ), one can show thatit has a power law form P ( g ) ∼ g − (1+ D F ) (55)The limiting value as n tends to infinity of the two maingaps, g ∗ , and the width of the spectrum ∆ ∗ have alsobeen computed in terms of ρ (see (Pi´echon et al. , 1995)for details). Gap labeling
Given the recursive structure of the spec-trum, it is easy to determine the labels for each of thegaps of the system of N = F n levels. For the j th plateauof the IDOS N ( E ) = N j /N , one must solve the rela-tion N j = Mod[q F n − , F n ] to obtain q. Stable gaps andtransient gaps are indicated by different colors for theirq labels in Fig.23 (black for stable and red for transient).As we mentioned before in Sec.III.C, the widths of gapstend to decrease with q. However as Fig.10 where thelogarithm of gap width is plotted versus q shows, thedependence is non-monotonic. The gap widths have log-periodic oscillations in q, and one sees a self-similar struc-ture. It is also interesting to note that the smallest valuesof q up to a certain maximum correspond to stable gaps.For q above a certain value, there are only transient gaps. FIG. 23 Recursive construction of the spectra of approximantchains showing the gap structure. Labels give the values of qdetermined according to Eq.20 in red for transient, and blackfor stable gaps (overbars indicate negative sign).
C. Multifractality of wavefunctions
The RG construction of the energy spectrum inSec.V.A has its parallel for the construction of wavefunc-tions. Consider the wavefunction ψ ( E ) for an allowedenergy E . If the energy is located in a side-band, thesupport of ψ ( E ) is concentrated on the molecular sites.If the energy is located in a side-band, the support of ψ ( E ) is concentrated on atom sites. Under RG, the ini-tial chain is transformed to a shorter chain, and the site FIG. 24 (top) Ground state probabilities ψ i plotted versussite index i . (middle) E = 0 state probabilities ψ i plottedversus site index i . (bottom) Probabilities ψ i plotted versussite index i for state a = 34. Data obtained by diagonalizationfor N=144 chain with PBC for ρ = 0 . i maps to site i (cid:48) of the new chain. One can introduce, aswe did for the energy recursion relations, wavefunctionrenormalization factors λ and λ , and write | ψ ( n ) i ( E ) | = λ | ψ ( n − i (cid:48) ( E (cid:48) ) | if E is atomic | ψ ( n ) l,r ( E ) | = λ | ψ ( n − i (cid:48) ( E (cid:48) ) | if E is molecular (56)where E (cid:48) denotes the energy after renormalization. In thecase of molecular RG (second line), there are two sites l (left) and r (right) forming the molecule correspondingto the site i (cid:48) . The wavefunction renormalization factorsare given by λ ≈ / λ ≈
1, to lowest order in ρ .The higher order corrections are important to keep for acorrect description of multifractality, as shown in (Mac´e et al. , 2016). Given the RG path of a state, with the helpof Eqs.56, one can reconstruct the corresponding ψ ( E ).To illustrate the different structures that are obtained,let us consider some examples. The RG path of the low-est level of the spectrum is { ... } , and the wavefunc-tion constructed recursively has the largest amplitudeson pairs of sites which derive from molecules on a largerlength scale, which derive from still bigger molecules andso on. The top part of Fig.24, shows the wavefunction1for a chain of 144 sites ( n = 11) as computed numeri-cally for a value of ρ = 0 .
2. One sees here the charac-teristic double-peak composed of double-peak structureon several scales. The figure shows deviations from theRG construction that we outlined: there are for examplesmall peaks which we ignored in our lowest order approx-imation. There are asymmetries in the amplitudes of thedouble peaks. These are due to higher order terms in ρ , which in our example, is not particularly small. Themiddle figure shows the wavefunction in the same chainfor the energy E = 0, whose RG path is { , , , ... } . Thepeaks are this time primarily localized on two atom sites.Finally the bottom figure shows an example of a ran-domly chosen mixed wavefunction that has both atomicand molecular components in its construction.
1. Hamiltonian in the conumber basis
In preparation for the discussion of wavefunctions, wediscuss the representation of the Hamiltonian in the con-umbering basis introduced in Sec.II.B.2. In this basis,the hopping Hamiltonian takes the form of a T¨oplitz ma-trix, where the non-zero elements lie at the distance F n − and F n − from the principal diagonal (Sire and Mosseri,1990). For example, the Hamiltonian for the n = 5 chain( N n = 8) can be written as follows H = . . . t A . t B . .. . . . t A . t B .. . . . . .t A . t B t A . . . . . t A .. t A . . . . . t A t B . . t A . . . . .. t B . . t A . . . .. . t B . . t A . . . As can be read off directly from the matrix, sites num-bered 1 through 3 and 6 through 9 have a strong bond t B and thus form molecules in the first RG step, while sites4 and 5 have weak bonds to either side and are atomsites. Within the two groups of F n − molecular sites,under a second RG step, one further has a subgroupinginto F n − molecules and F n − atoms of “second genera-tion”. This sub-grouping occurs in the middle block aswell. The conumbering scheme thus automatically classi-fies sites according to the same rules as the band struc-ture.
This remark has its importance for the discussionof wavefunctions which follows next.
2. Energy-position symmetry
As we have said, the construction of states and of thespectrum follow the same schema. Classifying sites ac-cording to their conumber corresponds, in fact, exactly to way energies are ordered. This leads to a remark-able approximate symmetry between states and energies.Fig.25a shows an intensity plot of the numerically calcu-lated (for an 89 sites chain with periodic boundary condi-tions) values of ψ ( E a ) j plotted against the conumber j for each allowed energy E a . Fig.25b shows the result forthe wavefunctions after four RG steps. Two observationscan be made: i) the similarity of the RG-constructed andthe numerical data is manifest and ii) the figures show areflection symmetry with respect to the diagonal, i.e. if i and a represent the position and the energy, then | ψ i,a | = | ψ a,i | (57)to lowest order in ρ . This striking symmetry betweenspatial and spectral variables holds for sufficiently small ρ .
3. Multifractal exponents for wavefunctions
The recursion relations for wavefunctions are analo-gous to those presented for the spectrum. These rela-tions involve the two different rescaling factors in therecursion formulae, λ and λ . Thus all the wavefunctionsare multifractal, in the same way that the spectrum ismultifractal due to the presence of two distinct shrinkingfactors z and z . Defining the fractal dimensions D ψq forthe wavefunction ψ ( E ) by χ nq ( E ) = (cid:88) i | ψ ni ( E ) | q ∼ (cid:18) F n (cid:19) ( q − D ψq ( E ) (58)For q = 2 the quantity χ ( E ) is nothing but the inverseparticipation ratio (IPR). For a given system size, theinverse of the IPR provides an indication of the spatialspread of the wavefunction. The scaling of the IPR as afunction of the system size determines whether or not astate is localized. The exponent D ψ ( E ) is an often usedindicator helping to locate the metal-insulator transitionin, for example, the 3D Anderson model. To recall, thevalue D ψ ( E ) = 1 indicates that the state E is extended,while D ψ ( E ) = 0 characterizes a localized state. In-termediate values 0 < D ψ ( E ) < n = 12 approximant. It is interesting to comparethis plot with Fig.22. The IPR of levels appears, in fact,to track their x values. This correlation is explained bythe fact that molecular states are more delocalized thanatom states, therefore, the higher the number of molecu-lar RG steps, given by the variable x = n m /n the largerthe IPR of the state. To be sure, this argument is onlyvalid for sufficiently small ρ .The full set of exponents D ψq ( E ) can be computed inthe perturbative RG scheme (Mac´e et al. , 2016). To low-est order, these are given by the value of x (which mea-sures the extent to which a given state is of molecular2 FIG. 25 Upper figure: Intensity plot of wavefunctions as afunction of the conumber (x-axis) versus the energy level in-dex (y-axis). Lower figure: the first few steps of the geomet-rical construction of the wavefunctions according to pertur-bation theory Eq.56 (reproduced with permission from (Mac´e et al. , 2016)) type), as follows D ψq, ( E ) = − x ( E ) log 2log ω + O ( ρ ) (59)This simple result for very small ρ gives a monofractal,since the fractal dimensions do not depend on q at leadingorder. For larger ρ , however, the higher order correctionsshow that the wavefunctions are indeed multifractal (seeMace et al (Mac´e et al. , 2016)). For small ρ , Eq.59 showsthat the smaller x , the smaller the fractal dimension.Themost extended states, according to Eq.59 are those atthe spectrum edges, where x = . The equation alsopredicts that D ψq, ( E ) = 0 upto higher order corrections for the level in the center of the band where x = 0. Thisindicates that the state is localized or close to being local-ized. However, as the exact calculations in the precedingsection Sec.III.F showed, it is in fact a critical state. Thediscrepancy is corrected by including the higher ordercorrections which are lacking in Eq.59 (as was done in(Mac´e et al. , 2016)). FIG. 26 a) Inverse participation ratios χ for all states ofa n = 12 chain ( ρ = 0 .
4. Spectrally averaged dimensions
In certain contexts, one may need to know, not thebehavior of a single eigenstate, but the average behaviorof wavefunctions close to a certain energy (such as theFermi energy). One can for example wish to determinethe average value of the fractal dimension D ψ ( E ), withinsome energy interval ∆ E . This averaged exponent occursfor example in some rigorous inequalities for dynamicalquantities as described in the next section. The averagedwavefunction exponents could also be relevant for otherphysical properties, such as the Kondo screening of impu-rities. These exponents can be calculated by consideringa generalization of the χ -function in Eq.58 to include asum over all energies. The results for the averaged di-mensions D ψ ( E ), for two different values of ρ , are shownin Fig.27 (taken from (Mac´e et al. , 2016)). VI. DYNAMICAL PROPERTIES
The quantum diffusion of wavepackets is determinedby the spectral properties of the underlying Hamiltonian.Extended wavefunctions and continuous spectra typicallylead to ballistic motion, where the particle moves with awell-defined group velocity. Singular continuous spectrasuch as that of the FC result in more complex behaviors.Some of the characteristics of wave packet diffusion arediscussed in this section.3
FIG. 27 The averaged fractal dimensions of the wavefunctions D ψq as a function of q , for ρ = 0 . , .
5. Dots: numericalresults, solid line: theoretical predictions.
A. The diffusion exponent β The mean square displacement in a time t of an elec-tron starting from the site i of the chain is given by d ( i , t ) = (cid:88) j ( i − i ) p ( i, i , t ) (60)where p ( i, i , t ) is the probability to be on site i at time t for a given starting position i , normalized such that (cid:80) i p ( i, i , t ) = 1. The time dependence of d ( i , t ) is inprinciple extremely complex, due to the multifractalityboth of the density of states and of the wavefunctions.The exponent β ( i ) describes its leading long time be-havior after smoothing fluctuations, i.e. d ( i , t ) ∼ t β ( i ) (61)for sufficiently long times t . In the simplest cases, val-ues of the exponent are well-known: for electrons in aperiodic crystal there is ballistic propagation with a con-stant group velocity. The system has translational invari-ance and β ( i ) = β independently of the initial position.For localized electrons, d ( t ) tends to a constant at largetimes, and β = 0. For standard diffusion, as in the caseof electrons in a moderately disordered crystal, β = 1 / β values depend on the choice of origin i ofthe wavepacket. One can define an effective averaged dif-fusion exponent β by considering the averaged quantity d ( t ) = (cid:104) d ( i , t ) (cid:105) where the brackets denote the averageover initial positions i . Then d ( t ) ∼ t β defines the glob-ally averaged value of the diffusion exponent β .Numerical results for the values of β in the FC for dif-ferent values of w = t A /t B are shown by circular symbolsin Fig.28 (taken from Thiem and Schreiber (Thiem andSchreiber, 2013)). One sees that the diffusion exponentincreases in value monotonically with t A /t B and reachesthe expected value of 1 in the limit of the periodic chain.Data for d dimensional product lattices (see Sec.X) alsoincluded in the figure, show that the exponent β does notdepend on the dimensionality in this class of models. FIG. 28 Plots of the diffusion exponent β as a function of thehopping ratio for the Fibonacci chain and its d -dimensionalgeneralizations as described in Sec.X. (reproduced with per-mission from (Thiem and Schreiber, 2013)) It can be noted that this exponent enters in manyphenomenological theories of transport in quasicrystals.From a dimensional analysis the diffusivity D ∼ d ( τ ) /τ and the conductivity (using the Einstein relation) shouldscale as τ β − , where τ is some characteristic cut-off timefor diffusion. For sub-diffusive motion, i.e. β < , theelectrical conductivity σ would then decrease with in-creasing τ . This behavior is the opposite of what onewould observe in metals, where conductivity increaseslinearly with τ according to the Drude formula. In realquasicrystals, in fact, experiments show that conductiv-ity decreases when structural quality of the sample isimproved for example by annealing (Mayou et al. , 1993). B. RG method for dynamical quantities
Studies of dynamics using the approximate RG schemehave been carried out by Abe and Hiramoto (Abe and Hi-ramoto, 1987), Pi´echon (Pi´echon, 1996) and Thiem andSchreiber (Thiem and Schreiber, 2012). These studiesuse a recursive approach to computing d ( i , t ) as well asall of the generalized moments of the displacement, d q ( i , t ) = 1 /N (cid:88) j ( i − i ) q p ( i, i , t ) (62)To lowest order, one can write two different recursive re-lations for the probability p ( n ) , depending on whether theinitial site is an atomic or a molecular site. Following theapproach in (Pi´echon, 1996), let the initial site be labeled i , and i refer to a site of the n th chain. In the renormal-ized chain, similarly, the initial site has the label i (cid:48) and i (cid:48) refers to a site of the new chain. The length scale renor-malization factors is either ω or ω depending on thetype of RG transformation (atom or molecule), while thecorresponding energy-time renormalization factors are z and z . The relations between probabilities defined on the4old and new chains can be stated as follows : p ( n ) ( i, i , t ) | ato ≈ ω n p ( n − ( i (cid:48) , i (cid:48) , zt ) p ( n ) ( i, i ± , t ) | mol ≈ ω n p ( n − ( i (cid:48) , i (cid:48) , zt ) (63)In the second relation, the initial site in the case of themolecule are further labeled with a ± standing for theleft and right atoms of that molecule. The first of theserelations says that the probability to go from site i tosite i in a time t in the n th chain is reduced by a factor ω with respect to the probability to go from site i (cid:48) to site i (cid:48) (a distance shorter by ω ) in a time zt in the n − d ( t ) ∼ t β ato with β ato = ln ω / ln z
2. In the opposite situation of a pure molecular diffu-sion process, one obtains another exponent β mol = ln ω / ln z Thiem and Schreiber defined an average value of β using the relative fraction of atom and molecular sites,as follows: β = τ − τ + 1 β ato + 2 τ + 1 β mol (64)They noted that this value fits the dynamics well forsmall ρ , up to about ρ = 0 .
02 (when the approximationsmade for the wavefunctions become inadequate).Using the recursion relations in Eq.63, Pi´echon derivedrecursion formulae for the moments d q ( i , t ), and for theiraverages d q ( t ) = (cid:104) d q ( i , t ) (cid:105) , as follows d nq ( t ) = ω − q ) n d n − q ( zt ) + 2 ω − q ) n d n − q ( zt ) (65)Let us suppose there is a fixed point solution of theprobability, p ∗ , in the limit n → ∞ . Eq.65 implies aself-consistency condition : p ∗ ( r, t ) = ω p ∗ ( ω r, zt ) + 2 ω p ∗ ( ω r, zt ) (66)Note that the probability depends on two scaling fac-tors leading to a multiscaling property: the function p ∗ ( r, t ) cannot be written as a function of a single vari-able. Consider now the probability p ( r = 0 , t ). Letting p (0 , t ) ∗ ∼ t − γ and comparing Eq.66 with Eq.53 one seesthat γ = D , one of the spectral dimensions introducedearlier.More generally, if one assumes that higher momentsof the diffusion distance scale as d q ( t ) ∼ t qσ q , then one finds that σ q = D − q . This result shows that the dif-fusion exponents σ q ≥ D , thus satisfying the Guarneriinequalities (Guarneri, 1993) for all values of q . This for-malism, and similar conclusions are applicable to otherquasiperiodic chains such as the AAH model studied by(Evangelou and Katsanos, 1993; Ketzmerick et al. , 1997). FIG. 29 Comparative plot of diffusion exponent β , general-ized spectral dimensions and wavefunction dimensions com-puted numerically for different values of the hopping ratio( w = t A /t B ) (reproduced with permission from (Thiem andSchreiber, 2013)) We have seen already in Sec.V.A that in perturbationtheory the spectrum was approximately monofractal for z = z / . This simplification occurs for the wavepacketdynamics as well – the exponent σ q = D F for all q , where D F = ln ω / ln z . In this case, the dynamics can be ex-pected to be simple diffusion with a single exponent How-ever this is only true to the extent of the approximationsmade (see the caveat based on the trace map analysisSec.III). Relations between exponents and inequalities
Guarneri(Guarneri, 1993) derived an inequality stating that β ≥ D , where D is the earlier mentioned information di-mension of the spectrum. Another inequality derivedby Ketzmerick et al (Ketzmerick et al. , 1997) involvesthe average wavefunction exponent D ψ and it statesthat β ≥ D /D ψ . Both these inequalities are satisfiedin the FC, as can be seen from Fig.29. For the diag-5onal Fibonacci model, (Damanik, 2006; Damanik andTcheremchantsev, 2007) showed that β must satisfy cer-tain bounds that depend on the onsite energies (cid:15) A and (cid:15) B . C. Autocorrelation function exponent
As for the mean square distance, the time dependentcorrelation function in quantum systems with Cantorspectra is expected to have a power law decay, fallingoff as t − δ at long times. Ketzmerick, Petschel and Geiselhave argued (Ketzmerick et al. , 1992) that the exponent δ should be equal to D . Several authors (Ketzmer-ick et al. , 1992; Thiem and Schreiber, 2012; Zhong andMosseri, 1995) have studied the behavior of the smoothedautocorrelation function given by, C ( i , t ) = 1 t (cid:90) t P ( i , i , t (cid:48) ) dt (cid:48) (67)which gives the integrated probability up to time t forthe particle to be found at the initial position. Av-eraging over all initial positions, one obtains C ( t ) =(1 /N ) (cid:80) i C ( i, t ) ∼ t − δ (cid:48) . This quantity is easier to fitthan the autocorrelation function. Note that δ (cid:48) may bedifferent from δ due to logarithmic corrections comingfrom the integral in Eq.67. Fig.30 taken from (Thiemand Schreiber, 2013) shows the smoothed autocorrela-tion function for different values of t A /t B (the variable w in the figure). The exponent should tend to the expectedvalue 1 in the periodic limit when the spectrum becomescontinuous ( t A = t B or (cid:15) A = (cid:15) B ). D. Log-periodic oscillations
Hiramoto and Abe observed that there are oscillationssuperposed on top of the average power law behavior ofthe rms distance d ( i , t ). These oscillations were studiedin more detail by Lifshitz and Even-Dar Mandel (Lifshitzand Even-Dar Mandel, 2011). A similar oscillatory be-havior is seen for the return probability p (0 , t ). Fig.31taken from Thiem (Thiem, 2015b), shows ln d ( t ) plottedversus ln t/t B . Several periods of oscillations of d ( t ) canbe clearly seen in the plot.The empirical form that was used to fit these datais d ( i , t ) ∼ t β ( i ) (1 + αe if ln t ). The frequencies f aredifferent, depending on the time scale that is considered,and they follow a hierarchical rule. The shortest timescorrespond to the fastest oscillations, whose frequency In practice, however, as pointed out in Yuan et al (Yuan et al. ,2000), it is hard to get convergence in numerical studies, thusmany early numerical works obtained values significantly smallerthan 1 even for the periodic case.
FIG. 30 Log-log plot of the smoothed autocorrelation func-tion for different values of the hopping ratio ( w = t A /t B ),showing the fit to power law. (reproduced with permissionfrom (Thiem and Schreiber, 2013))FIG. 31 depends on the nature of the initial site i – whetherit is atom or molecule. At the shortest time scale theoscillations have a period (in dimensionless units) of 2 π for an initial site of type atom, and π when the initial siteis molecule type. Going to longer times the oscillationshave frequencies which are smaller by factors z and z respectively. Thiem has given a quantitative account ofthese oscillations in terms of the perturbative RG theory(Thiem, 2015b). They stem, at each length scale of theRG process, from resonances due to the molecular energy6level splitting. The author notes that, on the other hand,such oscillations are not observed for the quasiperiodiccritical AAH model (Thiem, 2015b). This may be dueto the difference of the potentials – since the potentialenergy is a continuous-valued function in the AAH model,there are no molecular clusters in its RG scheme, andthus no characteristic resonance frequencies. VII. TRANSPORT
Some frequently asked questions concern the resistivityof a quasicrystal: how does it depend on the sample size,disorder, temperature, on external magnetic fields, etc? The following studies of transport properties of 1DFibonacci chains help to shed light on these questions.
A. Exact result for E = 0 The E = 0 state transmission coefficient in the hoppingmodel can be calculated from the exact solution given inSec.III.F. This transmission coefficient is proportional tothe zero temperature conductivity at half-filling, whenthe Fermi energy of the system is E F = 0. Let us con-sider a FC of length L = 2 n attached to periodic chains(input and output chains) at either end. The transmis-sion coefficient is given by (Beenakker, 1997; Economouand Soukoulis, 1981) T n = 4( x n + x − n ) (68)where x n = | ψ (2 n ) /ψ (0) | . Using the expression eq.37one obtains (Mace et al. , 2017) T n = 1cosh ( κ ( h ( n ) − h (0)) (69)where h the height function, it will be recalled, dependssolely on the geometry. This expression shows that thetransmission T n = 1 – there is perfect transmission –when the heights h ( n ) and h (0) are equal. This type of“intermittent” transparency occurs for sites separated bydistances that can tend to infinity.One can also compute the harmonic mean of the trans-mission over a chain of L + 1 sites. The harmonic meanwhich is more adapted than the arithmetic mean for sys-tems with very large fluctuations is defined by (cid:104)T (cid:105) L = (cid:32) L − L (cid:88) i T − i (cid:33) − (70)The scaling of (cid:104)T (cid:105) L with system size L can be expressedanalytically using the height distribution Eq.38, whichyields (cid:104)T (cid:105) L ∼ (cid:32) LL φ (cid:33) − (71) FIG. 32 Mean E = 0 transmission coefficient as a function oflength L of the Fibonacci chain. The dashed and continuouslines show the analytical prediction of Eq.71 and numericalresults for n = 27 ( N= 196418 atoms) respectively where the exponent φ can be expressed in terms of κ ,the golden mean τ and the maximal eigenvalue of a cer-tain product of generalized inflation matrices (see (Mace et al. , 2017) for details). This expression predicts a powerlaw decay of the mean transmission as a function of dis-tance L , in accord with the numerical data, as seen inFig.VII.A. Note that there is no contradiction betweenthis power law and the fact that the chain can be trans-parent for very long distances. The apparent paradox isexplained by noting that the power law represents an average behavior, while for a given system there arelarge fluctuations around the mean (as already observedby Sutherland and Kohmoto (Sutherland and Kohmoto,1987b)), whereas in a given chain, the transmission co-efficient can be (and is) equal to 1 for certain sites. B. Landauer approach
The Landauer formula relates the resistivity ρ ( n ) tothe transfer matrix T n of a system of n sites as follows ρ n = 14 (cid:0) T Tn T n − (cid:1) (72)where T T denotes the transpose of the transfer matrix.Using the trace map techniques described in Sec.III,Sutherland and Kohmoto (Sutherland and Kohmoto,1987b) studied the behaviors near the band edge and theband center. They showed that the resistance grows nofaster than a power law of the system size. They pointedout as well that there is a wide distribution of powersgoverning its growth with system size, and this leads tovery large fluctuations of the resistance. They specu-lated, finally, that this behavior is qualitatively also to beexpected for other energies in the band. Their conjectureas to power law behavior of the resistivity was proved byIochum and Testard (Iochum and Testard, 1991).7For comparison with measurements, it is pertinent toconsider the average resistance where the average is takenover states lying within some energy interval. The energyinterval chosen should depend on factors such as the tem-perature, or energy-scale corresponding to inelastic scat-tering, disorder conditions, etc. For a system of length L and some appropriately chosen ∆ E , ρ ( L ) = 1 n ( E )∆ E (cid:88) E ρ ( L, E i ) (73)defines an average resistance which was studied for theFC by Das Sarma and Xie using the Landauer formalism(Das Sarma and Xie, 1988). The result is a power lawbehavior for the resistivity, ρ ∼ L − a . The authors notedthat the power law holds for other fillings, provided thatthe Fermi level is not close to a large gap. C. Kubo-Greenwood approach
Sanchez et al (S´anchez et al. , 2001; S´anchez and Wang,2004) developed an RG approach for the conductivitystarting from the Kubo-Greenwood formula: σ ( µ, ω, T ) = 2 e ¯ hπm V (cid:90) dE f ( E ) − f ( E + ¯ hω )¯ hω Tr (cid:0) Im G + ( E + ¯ hω )Im G + ( E ) (cid:1) (74)where V is the volume of the system, G + ( E ) is the re-tarded one-particle Green’s function, f is the Fermi-Diracdistribution, with Fermi energy µ for temperature T . Ap-plying the RG method to the Kubo-Greenwood formulathey could study the AC and DC conductivities for verylarge systems (S´anchez and Wang, 2004). In the hoppingproblem, for half-filling, Sanchez and Wang found thatthe scaling exponents of the conductivity and of the den-sity of states have a similar dependence. Fig.33 showsthe parameter b (defined by σ ∼ b − n/ ) and the param-eter d (defined by DOS ∼ d − n/ ) plotted against t A /t B .The authors concluded that these results show the exis-tence of an Einstein relation, σ ∼ ( dN/dE ) D , whereinthe conductivity and density of states are related by thediffusivity D . D. Many body metallic and insulating states
Even in the absence of electron-electron interactions,many body effects can play an important role in qua-sicrystals. Varma et al considered many body effects –this amounts to taking into account the Pauli principle– on the conductance in Fibonacci chains (Varma et al. ,2016). They computed the effective localization length Λfor various band fillings, using a formalism developed byKohn in his theory of the insulating state (Kohn, 1964).This localization length, which is related to the real partof the conductivity tensor, is more sensitive to spectralgaps, and to transport properties than the more famil-iar single particle localization length defined in terms of
FIG. 33 Plot of the conductivity exponent b (circles, for thedefinition see text) and the density of states exponent d (fol-lowing the analytic formula (Kohmoto et al. , 1987), dashedline) as a function of the hopping ratio γ = t A /t B (figurereproduced from (S´anchez and Wang, 2004) with permission the decay of the envelope of a given wavefunction. Ina conductor, this many body localization length scalesto infinity, while in an insulator, this quantity saturateswith the system size. Varma et al showed that in thehopping model, there are special values of the band fill-ing corresponding to an insulating state. These fillingscorrespond to values of IDOS of ω and ω correspond tothe positions of the large gaps. This is shown in Fig.34where Λ − is plotted versus inverse system size 1 /L . Thefillings and correspond to a metallic state. FIG. 34 Scaling of the many body localization length Λ withsystem size L , for different band fillings. One sees metallicbehavior (Λ decreasing with L ) for fillings of 1/2 and 1/4,and insulating behavior (Λ constant with L ) for fillings g − and g − where g is the golden mean VIII. PERTURBATIONS. DISORDER AND BOUNDARYEFFECTS
Disorder can be expected to play an important role forthe electronic transport in real quasicrystals for the fol-lowing reason. We have seen in the previous section that,for a perfect FC, transport at T = 0 can be describedby different power laws, and can either scale towards aninsulating or to a metallic state in the thermodynamiclimit, depending on the Fermi energy. Adding disorderof any form – be it structural defects or chemical sub-stitutions or phonons – introduces a cut-off timescale, τ , whose value depends on disorder. One can thereforehave a conductance which is finite in a quasicrystal thatis weakly disordered. For increasing disorder, one expectsthat the quantum interference phenomena which lead, ina perfect quasicrystal to the multifractality of the den-sity of states and the multifractality of states, would beprogressively suppressed. Combined, these effects mightcontribute to improving transport as disorder increases,however, these are still open questions.The problem of a single impurity in a FC was studiedby several authors. The trace map method was used tocompute localization lengths of an impurity state repre-sented by a delta-function potential in (Naumis, 1999).The effect of structural defects, namely a single phasondefect (in which a pair of bonds is locally exchanged, asfor example AB → BA has been studied (Naumis andArag´on, 1996). This study concluded that the presenceof a single impurity affects all of the states and leads toan increase of the fractal dimension of the spectrum. Aweak form of structural disorder was considered in (Vel-hinho and Pimentel, 2000), by allowing randomness inthe substitution rules for building chains. The conclu-sion reached in this and a later study (Huang and Huang,2004) is that this type of disorder is irrelevant, in thatthe Lyapunov exponents of states were not changed. Tomodify the critical states of the pure system, the disordermust break some symmetries of the Fibonacci Hamilto-nian as in the models we discuss in the following section. A. Finite disorder and approach to Anderson localization
We now focus on the effects of adding a finite bulk dis-order to the two models under discussion, Eqs.13 and 14.Quite general rigorous arguments show that, in one di-mensional models an infinitesimal disorder leads to local-ized states (Delyon et al. , 1985). This has been confirmedby a variety of numerical studies (Liu and Riklund, 1987;Naumis, 1999). Das Sarma and Xie (Das Sarma and Xie,1988) studied the effect of randomness in a Fibonacci QCusing a scattering model for a system in which the scat-terers of constant height were placed in a Fibonacci se-quence of two distances a and b . The system was coupledto leads and conductance computed using the Landauer formula G = e ¯ h T / (1 − T ), where T is the transmis-sion coefficient of F n site system. They reported that,while disorder (in the positions of the scatterers) even-tually localizes all states, small disorder did not changethe physics qualitatively. Introducing a shuffling of thesequence of Kronig-Penney type scatterers leads as wellto localization of states and consequently an exponentialdecay of the conductance.There is no doubt that sufficiently large disorderstrength leads to strong localization. Interestingly, how-ever, the approach to localization can be quite compli-cated and state-dependent. It has been shown in (Jagan-nathan et al. , 2019) that there are interesting crossoverphenomena going from the pure system to the localizedsystem as the strength of disorder is increased. The au-thors studied the change of critical states of the purehopping model when hopping amplitudes are randomlyperturbed from their initial values t B and t A = ρt B .While most states tend to become more and more lo-calized as the disorder strength is increased, some statesgo the other way initially, delocalizing before starting tolocalize like the other states. This can be seen in Fig.35which shows the IPRs for the pure chain (blue) and theaveraged IPRs for weakly disordered chains (red). Manyof the IPR values increase under adding disorder, i.e. thewavefunction gets more localized. However, for certainstates, it can be seen that the IPR value decreases forweak disorder, indicating delocalization. This result isrobust under change of system sizes and occurs for stateshaving a particular kind of RG path. FIG. 35 Log(IPR) plotted for each of the states for a N = 144Fibonacci chain without disorder (blue dots) and with weakdisorder (orange dots). Note the reduction of the IPR forcertain states (for example α = 11 and α = 24) which showre-entrant localization. The complete set of IPRs can be described in terms ofscaling functions – i.e. for a given state one can collapsethe data for different system sizes L and different disorderstrengths W onto a single curve. The finite size analysisin (Jagannathan et al. , 2019) showed that despite theirdifferent approaches to localization, ALL states are de-scribed by a single critical exponent, ν . The value of φ depends on the ratio t A /t B and was found numericallyfrom finite size scaling plots as in Fig.36 – ν = 0 .
53 for ρ = 0 .
33. The re-entrant behavior of the IPR can be ex-plained in terms of the perturbative RG theory (see (Ja-gannathan et al. , 2019; Jagannathan and Tarzia, 2020)for details). However the phenomenon seems to be moregeneric, and the re-entrance behavior is observed for the9diagonal model, as well as for other generalizations of theHamiltonian.
FIG. 36 Data for the normalized IPR of several low lyingstates plotted versus the scaling variable
W L /ν showing thedata collapse for different disorder and system sizes (from(Jagannathan and Tarzia, 2020) reproduced with permission) B. The proximity effect
It is well-known that it is possible to induce super-conducting correlations in a non-interacting conductingsystem by coupling it to a superconductor – the so-calledproximity effect. The proximity effect provides a way toobserve experimentally the unique properties of criticalstates in the FC.Connecting the FC to a bulk superconductor and usinga mean field theory, Rai et al (Rai et al. , 2019b) computedthe distribution of ∆ i , the induced local superconductingorder parameter on site i . There are very large spatialfluctuations of ∆ i , due precisely to critical states. Fig.37(top) shows the profile of the order parameter (OP) as afunction of the site number. One sees here the character-istic multifractal properties reflected in the variations ofthe order parameter. Fitting the average curve obtainedby changing the phason angle parameter φ , one sees thatthe OP decays as a power law in the distance from theN-S interface (Rai et al. , 2020). The power, which varieswith t A /t B , is expected to depend both on the exponentof the density of states at the Fermi level and the aver-aged fractal dimensions of the wavefunctions near E = 0.Along with the other states, edge modes contribute tothe induced order parameter on each site. In fact, if onecycles through chains as a function of the phason angleparameter φ in Eq.11, the induced order parameter ∆ i at a given site oscillates and the periods are just thetopological numbers of gaps close to the Fermi energy(Rai et al. , 2019a). This can be seen in the Fig.38 whichshows the variations of the OP at the midpoint of thechain as a function of the phason angle φ . Two differentchains are shown to emphasize that the basic periods donot change when going from smaller to larger systems,only additional periods appear. The lower figure showsthe power spectrum of the oscillations, and the periodswhich are present in the curve of ∆ mid . One sees the FIG. 37 The superconducting order parameter ∆( i ) in a Fi-bonacci chain placed in contact with a BCS superconductorat both ends, plotted versus position showing power law decay( t A /t B = 0 .
8, decay exponent=0.6) periods 4 , , , ,
6, – these correspond precisely to the q values of the largest gaps near the Fermi energy (here, E F = 0). IX. GENERALIZED FIBONACCI MODELSA. Phonon models
Phonon modes in a quasiperiodic chain can be studiedby considering the set of equations Hψ n = K n − ψ n − + K n ψ n +1 − ( K n − + K n ) ψ n = Eψ n (75)where ψ n denotes the displacement of atom n of mass m with respect to its equilibrium position, and the cou-plings K n can take two values, K A or K B . One canconsider, alternatively, another version of the model inwhich the masses m can vary, and the couplings are con-stant. The operator H generalizes the discretized Lapla-cian operator, and its eigenvalues E = mω yield thefrequencies of phonon modes. The phonon problem istackled using the methods we have seen already for theclosely related electron problem (Kohmoto et al. , 1983a;Lu et al. , 1986; Luck and Petritis, 1986; Nori and Ro-driguez, 1986; Ostlund et al. , 1983; ? ). The trace mapequation is the same as in the electron problem, namelyEq.25, and the same kind of analysis applies. The spec-trum of energies E = mω has a Cantor-set structureas can be seen in Fig.39 (taken from (Luck and Petri-tis, 1986)) which shows the IDOS versus energy (thesequantities are denoted in the figure by H and z ). In con-trast to the electronic case, for phonons, the scaling isnon-uniform as a function of energy, and gaps becomevery small as E tends to zero. At the lowest frequen-cies, the integrated density of states, which is plotted inFig.39 looks almost indistinguishable from that of a pe-riodic chain. Like the periodic chain, the IDOS has a0 FIG. 38 (top) The superconducting order parameter ∆( mid )at the mid-point of a Fibonacci chain versus the phason angle φ showing oscillations. (bottom) Fourier spectrum showingsome main periods present in the oscillations. (reproducedwith permission from (Rai et al. , 2019a)) van Hove singularity, IDOS ∼ √ E , for small E . Thisbehavior seems to indicate at first glance that the longwavelength Goldstone modes are robust with respect tothe quasiperiodic modulation. However Luck and Petritispresented a rigorous argument to show that the spectrumdoes not have any absolutely continuous component, evenfor frequency tending to zero.The gap labeling theorem is seen to hold, as expected,and the letters A, B and C indicate some important gapshaving gap labels given by the three smallest Fibonaccinumbers (Luck and Petritis, 1986). Luck and Petritisshowed that the upper edge of the spectrum is describedby a new 6-cycle a → − b → − a → b → − a → − b .This observation was then used to show that, close tothe upper edge, the quantity 1 − N ( E ) where N is theIDOS has a power law modulated by log-periodic oscil-lations. Ashraff and Stinchcombe (Ashraff and Stinch-combe, 1989) computed the dynamic structure factor forthe Fibonacci chain. Quantum diffusion properties in FIG. 39 (from (Luck and Petritis, 1986) reproduced with per-mission) this model have been studied (Kohmoto and Banavar,1986; Lifshitz and Dar Mandel, 2011) and found to shareelectronic properties of multifractal structure and logperiodic oscillations. For more information on phononmodes in quasicrystals we refer the reader to (Janssen et al. , ????).
B. Mixed Fibonacci models
The term mixed models is used to denote a generalmember of the family of models Eq.12, where diago-nal and off-diagonal quasiperiodic modulations are bothpresent. These are relevant for experiments, for exam-ple, as real systems can be expected to have both formsof quasiperiodicity.Many of the techniques, including the powerful transfermatrix analysis can be extended to mixed models. Maciaand Dominguez-Adame (Macia and Dominguez-Adame,1996) considered a mixed model in which the A and Bsites have onsite energies of α or β , following a Fibonaccisequence. The hopping amplitudes are assumed to havetwo possible values, t AB = t BA and t AA = γt AB . Theinitial step consists of defining the basic transfer matri-ces. Choosing, without loss of generality, units such that β = − α and t AB = 1 one obtains four different transfermatrices in this model as follows : X = (cid:20) ( E + α ) −
11 0 (cid:21) , Y = (cid:20) ( E − α ) /γ − /γ (cid:21) ,W = (cid:20) ( E − α ) − γ (cid:21) , Z = (cid:20) ( E − α ) −
11 0 (cid:21) (76)Macia and Dominguez-Adame showed that after renor-malization the global transfer matrix in this model has astructure identical to the one for the Fibonacci sequencefor the diagonal model, Eq.22. This is achieved by defin-ing blocks of sites via T A = ZY X and T B = W X . They1showed that in finite chains the energies of certain trans-parent states (with transmission coefficients of unity) areof the form E ( k ) = ± (cid:112) α + 4 cos ( kπ/N ) (77)with | α | ≤ k an integer such that N τ = kπ with k = 0 , , ... . The states with such energies E ( k ) were con-firmed as being extended states via a multifractal analy-sis. In the limit of infinite system lengths their transmis-sion coefficient is less than unity, but remains finite.Two such states are shown in Figs.40, having transmis-sion coefficients of about 0.59 and 0.74. These states arenon-periodic (appearances to the contrary). Other crit-ical states of extended and self-similar character whichare fully transparent were seen to exist as well (bottomfigure). It is interesting to note that the distribution ofcharge does not provide an indication of the transmissiv-ity of states: the charge distribution in the self-similarstates are less homogeneous than those shown in Figs.40,yet their transmission coefficients are higher. This andother mixed models, and many of their electronic prop-erties are described in a review article (Maci´a, 2005).In a different approach using perturbation theory, Sireand Mosseri have investigated a model for approximantchains in which t i takes two values, t A or t B , accordingto the Fibonacci sequence. The onsite potentials V weretaken to be V AA = − λ/ V AB = V AB = λ/
2. Themodel thus has two parameters λ and ρ . By consider-ing n th generation approximant chains, Sire and Mosserishowed that there are gap closings and quasi-extendedstates for two different families of solutions. In partic-ular, they showed that as n → ∞ there are Bloch-typeextended wavefunctions that can be described in terms ofa wave-vector k ∈ (0 , π ). The energies of such states aredistributed throughout the band, with the band edgescorresponding to k = 0 and k = π (Sire and Mosseri,1990). Similar conclusions as to the existence of suchextended Bloch-type states were reached by Kumar andAnanthakrishna (Kumar, 2017; Kumar and Ananthakr-ishna, 1987). Extended states with periodic envelopesmay exist in mixed systems can exist even for certaindisordered cases, as reported in (Huang and Gong, 1998). C. Interference and flux dependent phenomena
Transmission properties of chains of Aharonov-Bohmrings of two different sizes, and connected in a Fibonaccisequence, have been studied (Chakrabarti et al. , 2003).The transport properties now become flux dependent.It was observed that transmission decreases as a powerlaw in the number of rings and that there are resonantstates for specific flux values. An RG analysis using theLandauer formalism and the trace map method showsthat the transmission coefficient possesses a self-similarstructure (Nomata and Horie, 2007).
FIG. 40 Extended state in a F approximant chain for t A /t B = 2, α = 0 .
75 and E = − .
25 b)Extended state for t A /t B = 2, α = 0 . E = 5 / F approximant with t A /t B = 2, α = 0 . E = E (1160)(from (Macia and Dominguez-Adame, 1996)reproduced with permission) X. RELATED QUASIPERIODIC MODELS
Related models of particular interest include a class ofaperiodic 1D chains which can be obtained by generaliz-ing the substitution rules we have introduced in Sec.II.A.Higher dimensional lattices obtained by direct products2of chains have also been studied.
A. Other substitutional chains
Although we have focused on a single quasiperiodicsystem, described by the golden mean τ , many of themethods used are generalizable to other irrationals. Thenature of the irrational number (algebraic or not) is ofprimary importance for the geometric properties and inconsequence for the electronic properties as well. As wehave said, quasicrystals are a special class of structuresbased on Pisot numbers. The so-called metallic meanswhich are solutions of the equation x − nx − et al. , 2000). A study of the multifractalexponents for the central E = 0 state for metallic meanchains is done in (Mace et al. , 2017). Energy spectraof generalized Fibonacci-type quasilattices having self-similar as well as quasiperiodic structure were studied in(Fu et al. , 1997). A gap labeling theorem is shown toexist in these cases. FIG. 41 Schema of 2D product lattices hosting two differ-ent types of hopping models. a) direct product Hamiltonianwith t A (long bonds)and t B (short bonds) b) Labyrinth model(hopping along one of the diagonals of each plaquette) (from(Thiem and Schreiber, 2013) reproduced with permission) Other well-known aperiodic (though not quasiperi-odic) systems obtained by substitution rules are theThue-Morse, the period-doubling and Rudin-Shapiro se-quences. We refer to the review in (Maci´a, 2005) for adiscussion of their electronic properties.
B. Product lattices
The 1D chain can be used as the basis for extensionsto arbitrary dimensions d . Taking d = 2, for example,the direct product C n × C n of two Fibonacci approxi-mants aligned along x and y axes forms a 2D lattice of squares and rectangles (Lifshitz, 2002). As can be seenin Fig.41a) its connectivity is that of the square lattice.The energy spectrum and wavefunctions for tight-bindingmodels on these direct product lattices have been stud-ied. In the pure hopping vertex model, electrons can hopalong the two directions with amplitudes t A or t B . TheHamiltonian is separable into two independent Fibonaccichain problems. The energies are the sum of two 1D ener-gies E ij = E i + E j , and the corresponding wavefunctionsgiven by the product ψ ij ( x, y ) = ψ i ( x ) ψ j ( y ) where E i and ψ i are solutions to the 1D problem. The propertiesof the spectrum depend on the value of ρ . The spectrumof the d = 2 product lattice is purely singular continuousfor ρ < ρ c where ρ c ≈ . ρ > ρ c , the spectrum has a continuous part (Mandel andLifshitz, 2008; Sire, 1989; Thiem and Schreiber, 2013).The labyrinth model (Sire et al. , 1989), a 2D variantbased on the direct product of chains (see Fig.41b), alsohas properties derivable from the 1D solutions. The gen-eralized dimensions describing multifractal properties ofwavefunctions in d dimensional product lattices was in-vestigated in (Thiem and Schreiber, 2011; Yuan et al. ,2000). The exponents for d -dimensions are simply pro-portional to the 1D exponents, D ψ,dq = dD ψ, q . Dy-namical exponents have been computed for these lattices(Thiem and Schreiber, 2012; Zhong and Mosseri, 1995).The return probability exponent shows a d dependence, γ d = dγ (1) . The diffusion exponent β is expected ac-cording to theory to be constant as the dimensionality d increases and this is indeed found numerically, as can beseen in Fig.28. The autocorrelation function exponent δ depends on the dimensionality, and for higher d increasesfaster as a function of the modulation strength parameter w = t A /t B – see Fig.42. A complete account of the elec-tronic properties of such d -dimensional generalizations isgiven in (Thiem and Schreiber, 2013). FIG. 42 Plot of the autocorrelation function exponent δ asa function of the strength of quasiperiodic modulation (thevariable w = t A /t B corresponds to ρ of our text) for productlattices of dimensions 1,2 and 3 (from (Thiem and Schreiber,2013) reproduced with permission) d -dimensional tilings The spatial connectivity of the3product lattice systems is simple in that they have anunderlying average periodic structure: d -dimensional hy-percubic lattices. The electronic properties of such prod-uct lattices are seen to be “inherited” from the parentchains. The situation is different for other quasiperi-odic tilings. For 2D and 3D tilings, such as the Penrosetiling, spectra and wavefunctions remain difficult to com-pute analytically, with the exception of the ground state(Mace et al. , 2017). As for the Fibonacci quasicrystal,nontrivial topological properties are to be expected inthese higher dimensional cases. The possibility of higherorder topological insulators based on the Penrose and oc-tagonal tilings has been discussed, for example, in (Chen et al. , 2020; Fulga et al. , 2016). XI. INTERACTIONS AND QUASIPERIODICITY
The topic of interacting quasiperiodic systems requiresa separate review. This section will be restricted to givingbrief outlines of some of the main contributions, alongwith a non-exhaustive list of references.The effects of quasiperiodic perturbations in interact-ing fermionic chains was investigated by Vidal et al forcontinuum models using renormalization group in (Vidal et al. , 1999, 2001). Considering in particular the case ofmetallic mean chains, they found that there was a metal-insulator transition (Vidal et al. , 2001) for repulsive inter-actions. Hiramoto (Hiramoto, 1990) did a Hartree-Fockanalysis to study the effect of weak interaction U. Thestudy showed that the singular continuous single parti-cle spectrum persists in the presence of interactions, incontrast to the critical Harper model where the singularcontinuous behavior is destroyed by U.In a study of the Hubbard model on a Fibonacci chainby weak-coupling renormalization group and density ma-trix renormalization group methods (Hida, 2001), Hidashowed that, for the diagonal Fibonacci model, weakCoulomb repulsion is irrelevant in the sense of RG andthe system will behave as a free Fibonacci chain. Forstrong Coulomb repulsion the system becomes a Mott in-sulator and, in the spin sector, can be modeled in termsof a uniform Heisenberg antiferromagnetic chain. For theoff-diagonal case, he obtained a Mott insulator with a lowenergy sector that could be described in terms of a Fi-bonacci antiferromagnetic Heisenberg chain. Gupta et al(Gupta et al. , 2005) studied the DC electrical conductiv-ity for half filling, using Hartree-Fock mean field theory,to see the interplay of interactions and quasiperiodicity.They concluded that, while each of these factors taken in-dividually tend to decrease the conductivity, there maybe enhancement of the conductivity due to competitionbetween them.The evolution of multifractality in an interactingfermion chain was studied in (Mac´e et al. , 2019). Con-trary to naive expectations, these authors found that adding repulsive interactions did not lead to enhanceddelocalization. Fig.43 shows the half chain von Neu-mann entropy plotted against time for different strengthsof the quasiperiodic potential (controlled by a parame-ter h , with h = 0 being the periodic case, and h = 1the strongly quasiperiodic chain). For a periodic chain,(black curve) the entanglement entropy grows as a power-law in the time, S ( t ) / ∼ t /z . The exponent z increasesas the strength of the quasiperiodicity h is increased.This further confirms that the free Fibonacci chain is in-termediate between the free delocalized state and an An-derson localized state, from the point of view of its trans-port properties. One also sees log-periodic oscillationssuperposed on the power-law behavior, especially visiblefor the yellow curve (strong quasiperiodicity limit). FIG. 43 Entanglement entropy as a function of time for dif-ferent values of the strength of quasiperiodic modulation h (black curve: periodic chain, yellow chain: strongly quasiperi-odic chain). The inset shows the dependence of the power z on h (reproduced with permission from (Mac´e et al. , 2019)) A. Heisenberg and XY chains
The properties of spin chains with Fibonacci couplingshave been investigated by a number of methods. Hidaanalyzed the spin 1/2 Heisenberg model using densitymatrix renormalization group finding that quasiperiodicmodulations are relevant in this case and that the groundstate is in a different universality class from that of theXY chain (or free particle) problem (Hermisson, 2000;Hida, 1999). The entanglement entropy S of aperiodiccritical chains was studied by Igloi et al (Igl´oi et al. ,2007). For these chains, the half chain entanglement in-creases as S ∼ c ln L + cst , where L is the chain lengthand c is the central charge (=0.5 for the periodic case).They found that for Fibonacci XY chains the quasiperi-odic modulation is marginal, in that the central charge c depends on ρ = J A /J B (the ratio of spin-spin couplings).For Fibonacci Heisenberg chains, the quasiperiodicity is4strongly relevant and the prefactor is given by c (0) ≈ . B. Many body localization
The question of many body localization due toquasiperiodic potentials was raised in (Iyer et al. , 2013;Khemani et al. , 2017). These studies asked whether theMBL transitions are different in the presence of quasiperi-odic potentials compared to random potentials, and if soin what ways. The MBL transition in the quasiperiodicAAH model, which can be experimentally realized in in-teracting boson and fermions cold atoms systems, hasbeen studied, and shown to lead to a new type of “non-random” universality class (Khemani et al. , 2017). It isinteresting to ask whether there are any significant dif-ferences between many body localization in AAH and inFibonacci chains. Details of the transition have been dis-cussed in (Mac´e et al. , 2019; Varma and ˇZnidariˇc, 2019).
C. Anomalous diffusion properties
Settino et al (Settino et al. , 2020) have studied dy-namics in interacting aperiodic many body systems in-cluding the FC . They showed that, for the onsite Fi-bonacci model, the singular continuous spectrum for thenon-interacting problem remains and induces an anoma-lous dynamics. Lo Gullo et al (Lo Gullo et al. , 2017)studied aperiodic discrete time quantum walk problems,relevant in quantum computing, and which could be re-alized using optical fibers (Nguyen et al. , 2020). Theycomputed the energy spectra and the spreading of aninitially localized wave packet for different cases, find-ing that in the case of Fibonacci and Thue-Morse chainsthe system is superdiffusive, whereas for Rudin-Shapiro,another substitutional chain, it is strongly subdiffusive.They propose that the different dynamics are linked tothe nature of the spectra in the two cases – singular con-tinuous for the former, discrete for the latter.
XII. EXPERIMENTAL SYSTEMS
Fibonacci sequences occur naturally in 3D icosahe-dral quasicrystals and also in dodecagonal quasicrystals.These are structures which are based on the golden mean.Fig.44a which shows an STM image of copper adatomsdeposited on an icosahedral AlPdMn quasicrystal (Mc-Grath et al. , 2012) provides a striking illustration ofthis connection. As the height profile in Fig.44b shows,the distances between columns is a Fibonacci sequence.These rows of aperiodically spaced layers are coupled tothe bulk, and the resulting Hamiltonians are likely tobe fairly complicated. For experimental investigations ofthe 1D model, it is therefore useful to fabricate artificial
FIG. 44 a) STM image of a 100 × et al. , 2012) systems in order to study the Fibonacci chain , and wewill now describe a few of these.The off-diagonal tight-binding Fibonacci model hasbeen experimentally realized in a polaritonic gas in aquasi-1D cavity (Baboux et al. , 2017; Tanese et al. , 2014).Some of the theoretical predictions for the energies andthe eigenmodes of this system were observed. The dis-crete scale invariance of the spectrum and the gap la-beling theorem were thus experimentally verified. Theexistence of gap states was checked, and their spatialstructure mapped out. Most interestingly, their topolog-ical winding numbers could be experimentally measuredby varying the phason angle φ .Optical waveguides fabricated by using a femtosecondlaser beam to inscribe quasiperiodic spatial modulationin a bulk glass have been studied in (Kraus et al. , 2012).Injecting light into these allowed to study propagatingand localized modes, and to demonstrate the topologi-cally protected edge modes predicted by theory. Topo-logical pumping of photons was demonstrated in (Verbin5 et al. , 2015), where the topological equivalence betweenthe Harper and the Fibonacci models was experimen-tally verified. Related subjects of recent investigations,outside the scope of this review, are higher order topolog-ical insulators using quasicrystals (Chen et al. , 2020) andtopological quantum computation using Fibonacci anyonchains (Chandran et al. , 2020; Feiguin et al. , 2007).Merlin et al(Bajema and Merlin, 1987; Merlin et al. ,1985) fabricated semiconductor superlattices following aFibonacci sequence and studied their properties usingRaman spectroscopy. The samples were composed of twotypes of films : 27 nm-thick layers of GaAs and 43 nm-thick layers of GaAlAs (such that the ratio of thicknessesis close to the golden mean). These layers were piled ontop of each other in a quasiperiodic sequence along the z direction. The Raman frequency shifts were comparedfor periodic and for quasiperiodic structures. In addi-tion to the acoustic phonons in this system, plasmon-polariton modes are argued to play an important rolein (Albuquerque and Cottam, 2003), where their Ramancross-sections for the samples are discussed. Hawrylaket al (Hawrylak et al. , 1987) have computed the plas-mon spectrum for such superlattices, described its scalingproperties and explicitly computed the f ( α ) spectrum.Metamaterials in the nanoscale also provide a widearray of possibilities for aperiodic structures. Opti-cal transmission spectra of photonic band-gap Fibonacciquasiperiodic nanostructures composed of both positive(SiO2) and negative refractive index materials are dis-cussed by de Medeiros et al (de Medeiros et al. , 2007).In another direction, there are possibilities to makequasiperiodic sequences in 1D biomaterial (DNA-based)systems and study the consequences for transport (Albu-querque et al. , 2005).Magnetic multilayers composed of Fe/Cr layers andstudied by MOKE (magneto-optic Kerr effect) and FMR(ferromagnetic resonance) should display self-similarmagnetization versus curves and interesting thermody-namic signatures (Bezerra et al. , 2001). Another promis-ing system is composed of epitaxially grown layers of Feand Au using ultra-high-vacuum vapor deposition (Suwa et al. , 2017) which were theoretically predicted to haveanomalous magnetoresistance (Machado et al. , 2012).Fibonacci nanowire arrays have been studied recentlyby Lisiecki et al who discuss their magnonic proper-ties and interest for possible applications (Lisiecki et al. ,2019).Quantum dots can be used to make artificial crystalsand quasicrystals (Kouwenhoven et al. , 1990). They havealready been used to study magnetotransport in a pe-riodic crystal. It may therefore be possible to studytransport in artificial Fibonacci chains made with quan-tum dots. On a quite different length scale, microwavepropagation in dielectric resonators have been used withsuccess to simulate the tight-binding model for graphene(Bellec et al. , 2013). Preliminary work shows that this is also an extremely versatile system in which to studyelectronic properties of Fibonacci chains Cold atoms in optical potentials constitute a particu-larly fertile ground to realize quasiperiodic models andstudy their properties under controlled conditions. Anumber of recent theoretical studies have thus lookedat generalizations of the tight-binding models that arerelevant to cold atom experiments. The discrete-valuedFibonacci potential is more difficult to realize experi-mentally, compared to the Harper model, which can berealized by applying an incommensurate laser potential(Fallani et al. , 2007; Lye et al. , 2007). Recently, how-ever, Singh et al (Singh et al. , 2015) proposed a means ofrealizing generalized Fibonacci models on chains basedon the cut-and-project method, in a 2D optical lattice.If realized, this would provide opportunities to study ex-perimentally multifractal states and probe the multiscaledynamics in the Fibonacci quasicrystal.
XIII. SUMMARY AND OUTLOOK
This review has relevant for 1D quasiperiodic tight-binding models. We have introduced some of the mainconcepts and techniques for the 1D Fibonacci Hamilto-nian, an important model from a fundamental viewpoint.This family of models is also interesting from the pointof view of applications : in electronic devices but alsofor their mechanical properties, use in phononics, optics,or magnetic systems to name only a few. Although wefocused here on electronic properties, there are many in-teresting and closely related problems for electromagneticwave modes in aperiodic media. The unique optical re-flectivity properties of aperiodic multilayers suggest ap-plications as mirrors with omnidirectional reflectivity forall polarizations of incident light over a wide range ofwavelengths (Axel and Peyri`ere, 2010), and more gener-ally in nanodevices (Macia, 2012).One of the principal characteristics of the quasicrystalis the existence of multifractalities as a function of theenergy, of the space coordinates and of temporal corre-lations. We have discussed these properties with explicitcalculations and for specific examples. We have men-tioned a few consequences of these multifractal states forphysical properties : transport, disorder induced local-ization/delocalization and the proximity effect.The models discussed here can hopefully serve asa starting point for investigations of more complexquasiperiodic systems. Much more could be done, exper-imentally, to investigate the distinctive dynamical prop-erties of quasiperiodic chains. The effects of interactionsis another subject for future investigations. Last but not F. Pi´echon and F. Mortessagne, private communication.
ACKNOWLEDGMENTS
I gratefully acknowledge numerous discussions withMichel Duneau, Pavel Kalugin, Jean-Marc Luck, Nico-las Mac´e and Fr´ed´eric Pi`echon.
REFERENCES
Abe, S., and H. Hiramoto (1987), Phys. Rev. A , 5349.Albuquerque, E., and M. Cottam (2003), Phys. Rep. ,225.Albuquerque, E. L., M. S. Vasconcelos, M. L. Lyra, and F. A.B. F. de Moura (2005), Phys. Rev. E , 021910.Ashraff, J. A., and R. B. Stinchcombe (1989), Phys. Rev. B , 2670.Aubry, S., and G. Andr´e (1980), Ann. Israel Phys. Soc (133), 18.Axel, F., and J. Peyri`ere (2010), Journal of Physics A: Math-ematical and Theoretical (3), 035005.Baake, M., and U. Grimm (2013a), Aperiodic Order , Chap. 4.Baake, M., and U. Grimm (2013b),
Aperiodic Order , Ency-clopedia of Mathematics and its Applications, Vol. 1 (Cam-bridge University Press).Baake, M., U. Grimm, and D. Joseph (1993), I. J. Mod. Phys.B , 1527.Baboux, F., E. Levy, A. Lemaˆıtre, C. G´omez, E. Galopin,L. Le Gratiet, I. Sagnes, A. Amo, J. Bloch, and E. Akker-mans (2017), Physical Review B (16), 161114.Bajema, K., and R. Merlin (1987), Phys. Rev. B , 4555.Barache, D., and J. M. Luck (1994), Phys. Rev. B , 15004.Beenakker, C. W. J. (1997), Rev. Mod. Phys. , 731.Bellec, M., U. Kuhl, G. Montambaux, and F. Mortessagne(2013), Phys. Rev. B , 115437.Bellissard, J., B. Iochum, E. Scoppola, , and D. Testard(1989a), Comm. Math. Phys. , 527.Bellissard, J., B. Iochum, E. Scoppola, and D. Tes-tard (1989b), Communications in Mathematical Physics (3), 527.Bezerra, C., E. Albuquerque, and M. Cottam (2001), PhysicaA: Statistical Mechanics and its Applications (1), 341.Chakrabarti, A., R. A. R¨omer, and M. Schreiber (2003),Phys. Rev. B , 195417.Chandran, A., F. J. Burnell, and S. L. Sondhi (2020), Phys.Rev. B , 075104.Chen, R., C.-Z. Chen, J.-H. Gao, B. Zhou, and D.-H. Xu(2020), Phys. Rev. Lett. , 036803.Damanik, D. (2006), Philosophical Magazine (6-8), 883,https://doi.org/10.1080/14786430500199286.Damanik, D., M. Embree, and A. Gorodetski (2015), Prog.Math. Phy. , 307.Damanik, D., and S. Tcheremchantsev (2007), J. Amer.Math. Soc. , 799.Dareau, A., E. Levy, M. B. Aguilera, R. Bouganne, E. Akker- mans, F. Gerbier, and J. Beugnon (2017), Phys. Rev. Lett. , 215304.Das Sarma, S., and X. C. Xie (1988), Phys. Rev. B , 1097.Delyon, F. m. c., Y.-E. L´evy, and B. Souillard (1985), Phys.Rev. Lett. , 618.Derrida, B., C. Itzykson, and J. Luck (1984), Commun.Math.Phys. , 115.Economou, E. N., and C. M. Soukoulis (1981), Physical Re-view Letters (9), 618.Evangelou, S. (1987), J. Phys. C , L295.Evangelou, S., and D. Katsanos (1993), J. Phys. A (Math) , L1243.Fallani, L., J. E. Lye, V. Guarrera, C. Fort, and M. Inguscio(2007), Phys. Rev. Lett. , 130404.Feiguin, A., S. Trebst, A. W. W. Ludwig, M. Troyer, A. Ki-taev, Z. Wang, and M. H. Freedman (2007), Phys. Rev.Lett. , 160409.Fu, X., Y. Liu, P. Zhou, and W. Sritrakool (1997), Phys.Rev. B , 2882.Fujiwara, T., M. Kohmoto, and T. Tokihiro (1989), Phys.Rev. B , 7413.Fulga, I. C., D. I. Pikulin, and T. A. Loring (2016), Phys.Rev. Lett. , 257002.Ganeshan, S., J. H. Pixley, and S. Das Sarma (2015), Phys.Rev. Lett. , 146601.Gluzman, S., and D. Sornette (2002), Phys. Rev. E ,036142.Godr`eche, C., and J. M. Luck (1992), Phys. Rev. B , 176.Gordon, A. Y., S. Jitomirskaya, Y. Last, and B. Simon(1997), Acta Math. (2), 169.Guarneri, I. (1993), Europhysics Letters , 729.Gupta, S., S. Sil, and B. Bhattacharyya (2005), Physica B , 299.Halsey, T., M. Jensen, and L. Kadanoff (1986), Physical Re-view A (2).Harper, P. (1968), Proc. Phys. Soc. A , 874.Hawrylak, P., G. Eliasson, and J. J. Quinn (1987), Phys.Rev. B , 6501.Hermisson, J. (2000), J. Phys. A , 57.Hida, K. (1999), J. Phys. Soc. Jpn. , 3177.Hida, K. (2001), Phys. Rev. Lett. , 1331.Hiramoto, H. (1990), J. Phys. Soc. Jpn. , 811.Hiramoto, H., and M. Kohmoto (1989), Phys. Rev. B ,8225.Hiramoto, H., and M. Kohmoto (1992), Intl. J. Mod. Phys.B , 281.Huang, D., and D. Huang (2004), Phys. Rev. B , 205124.Huang, H., and F. Liu (2018), Phys. Rev. Lett. , 126401.Huang, H., and F. Liu (2019), Phys. Rev. B , 085119.Huang, X., and C. Gong (1998), Phys. Rev. B , 739.Igl´oi, F., R. Juh´asz, and Z. Zimbor´as (2007), EurophysicsLetters (EPL) (3), 37001.Inui, M., S. A. Trugman, and E. Abrahams (1994), Phys.Rev. B , 3190.Iochum, B., and D. Testard (1991), J. Stat. Phys. , 715.IUCr, (1992), Acta Crystallographica A48 , 928.Iyer, S., V. Oganesyan, G. Refael, and D. A. Huse (2013),Phys. Rev. B , 134202.Jagannathan, A., P. Jeena, and M. Tarzia (2019), Phys. Rev.B , 054203.Jagannathan, A., and M. Tarzia (2020), Eur. Phys. J. B ,46.Janssen, T., G. Chapuis, and M. de Boissieu (????), in Ape-riodic Crystals , Chap. 4 (Oxford University Press). Kalugin, P., and A. Katz (2014), J. Phys. A , 315206.Kalugin, P. A., K. A. Yu., and L. S. Levitov (1986), JETP (2), 410.Ketzmerick, R., K. Kruse, S. Kraut, and T. Geisel (1997),Physical Review Letters (September), 10.1103/Phys-RevLett.77.1413, arXiv:9604005 [quant-ph].Ketzmerick, R., G. Petschel, and T. Geisel (1992), Phys. Rev.Lett. , 695.Khemani, V., D. N. Sheng, and D. A. Huse (2017), Phys.Rev. Lett. , 075702.Kohmoto, M., and J. Banavar (1986), Physical Review B ,563.Kohmoto, M., L. P. Kadanoff, and C. Tang (1983a), Phys.Rev. Lett. , 1870.Kohmoto, M., K. L.P., and C. Tang (1983b), Physical ReviewLetters (3), 1870.Kohmoto, M., and Y. Oono (1984), Physics Letters ,145.Kohmoto, M., B. Sutherland, and C. Tang (1987), PhysicalReview B (3), 1020.Kohn, W. (1964), Phys. Rev. , A171.Kouwenhoven, L. P., F. W. J. Hekking, B. J. van Wees, C. J.P. M. Harmans, C. E. Timmering, and C. T. Foxon (1990),Phys. Rev. Lett. , 361.Kraus, Y. E., Y. Lahini, Z. Ringel, M. Verbin, and O. Zil-berberg (2012), Phys. Rev. Lett. , 106402.Kumar, V. (2017), in Journal of Physics: Conference Series ,Vol. 809 (IOP Publishing) p. 012023.Kumar, V., and G. Ananthakrishna (1987), Phys. Rev. Lett. , 1476.Levitov, L. (1989), J. Phys. (France) , 707.Lifshitz, R. (2002), J. Alloys and Compounds , 186.Lifshitz, R., and S. Dar Mandel (2011), Phil. Mag. , 2792.Lifshitz, R., and S. Even-Dar Mandel (2011), Phil. Mag. ,2792.Lisiecki, F., J. Rych(cid:32)ly, P. Ku´swik, H. G(cid:32)lowi´nski, J. W. K(cid:32)los,F. Groß, N. Tr¨ager, I. Bykova, M. Weigand, M. Zelent,E. J. Goering, G. Sch¨utz, M. Krawczyk, F. Stobiecki,J. Dubowik, and J. Gr¨afe (2019), Phys. Rev. Applied ,054061.Liu, F., S. Ghosh, and Y. D. Chong (2015), Phys. Rev. B ,014108.Liu, Y., and R. Riklund (1987), Phys. Rev. B , 6034.Lo Gullo, N., C. V. Ambarish, T. Busch, L. Dell’Anna, andC. M. Chandrashekar (2017), Phys. Rev. E , 012111.Lu, J. P., T. Odagaki, and J. L. Birman (1986), Phys. Rev.B , 4809.Luck, J., and D. Petritis (1986), Journal of Statistical Physics , 289.Luck, J. M. (1989), Phys. Rev. B , 5834.Luck, J.-M. (1993), Europhysics Letters , 359.Lye, J. E., L. Fallani, C. Fort, V. Guarrera, M. Modugno,D. S. Wiersma, and M. Inguscio (2007), Phys. Rev. A ,061603.Mace, N., A. Jagannathan, P. Kalugin, R. Mosseri, andF. Piechon (2017), Physical Review B (4), 045138.Mac´e, N., A. Jagannathan, and F. Pi´echon (2016), PhysicalReview B (20), 205153.Mac´e, N., A. Jagannathan, and F. Pi´echon (2017), in Journalof Physics: Conference Series , Vol. 809 (IOP Publishing)p. 012023.Mac´e, N., N. Laflorencie, and F. Alet (2019), SciPost Phys , 050.Machado, L. D., C. G. Bezerra, M. A. Correa, C. Chesman, J. E. Pearson, and A. Hoffmann (2012), Phys. Rev. B ,224416.Maci´a, E. (2005), Reports on Progress in Physics (2), 397.Macia, E. (2012), Rep. Prog. Phys. , 036502.Macia, E., and F. Dominguez-Adame (1996), Physical Re-view Letters , 2957.Mandel, S. E.-D., and R. Lifshitz (2008),Philosophical Magazine (13-15), 2261,https://doi.org/10.1080/14786430802070805.Mayou, D., C. Berger, F. Cyrot-Lackmann, T. Klein, andP. Lanco (1993), Phys. Rev. Lett. , 3915.McGrath, R., H. Sharma, J. Smerdon, and J. Ledieu (2012),Phil. Trans. Roy. Soc. A , 2930.de Medeiros, F., E. Albuquerque, and M. Vasconselos (2007),Surface Science , 4492.Merlin, R., K. Bajema, R. Clarke, F. Y. Juang, and P. K.Bhattacharya (1985), Phys. Rev. Lett. , 1768.Naumis, G. G. (1999), Phys. Rev. B , 11315.Naumis, G. G., and J. L. Arag´on (1996), Phys. Rev. B ,15079.Nguyen, D., T. Nguyen, and R. e. a. Khrapko (2020), Sci.Rep. , 7156.Niu, Q., and F. Nori (1986), Phys. Rev. Lett. , 2057.Niu, Q., and F. Nori (1990), Physical Review B (16),10329.Nomata, A., and S. Horie (2007), Phys. Rev. B , 115130.Nori, F., and J. P. Rodriguez (1986), Phys. Rev. B , 2207.Ostlund, S., and R. Pandit (1984), Phys. Rev. B , 1394.Ostlund, S., R. Pandit, D. Rand, H. J. Schellnhuber, andE. D. Siggia (1983), Phys. Rev. Lett. , 1873.Pi´echon, F. (1996), Physical Review Letters , 4372.Pi´echon, F., M. Benakli, and A. Jagannathan (1995), Physi-cal Review Letters (26), 5248.Rai, G., S. Haas, and A. Jagannathan (2019a), Phys. Rev. B , 165121.Rai, G., S. Haas, and A. Jagannathan (2020), Phys. Rev. B xxx , xxx.Rai, G., S. W. Haas, and A. Jagannathan (2019b), PhysicalReview B , 165121.R¨ontgen, M., C. V. Morfonios, R. Wang, L. Dal Negro, andP. Schmelcher (2019), Phys. Rev. B , 214201.Rudinger, A., and F. Piechon (1998), J. Phys. A (Math.) ,155.Rudinger, A., and C. Sire (1996), J. Phys. A (Math) ,3537.S´anchez, V., L. A. P´erez, R. Oviedo-Roa, and C. Wang(2001), Phys. Rev. B , 174205.S´anchez, V., and C. Wang (2004), Phys. Rev. B , 144207.Settino, J., N. W. Talarico, F. Cosco, F. Plastina, S. Manis-calco, and N. Lo Gullo (2020), Phys. Rev. B , 144303.Shechtman, D., I. Blech, D. Gratias, and J. W. Cahn (1984),Phys. Rev. Lett. , 1951.Singh, K., K. Saha, S. A. Parameswaran, and D. M. Weld(2015), Phys. Rev. A , 063426.Sire, C. (1989), Europhysics Letters , 483.Sire, C., and R. Mosseri (1989), J. Phys. (France) , 3447.Sire, C., and R. Mosseri (1990), J. Phys. (France) , 1569.Sire, C., M. Remi, and J. Sadoc (1989), J. Phys. (France) , 3463.Sutherland, B., and M. Kohmoto (1987a), Physical ReviewB , 5877.Sutherland, B., and M. Kohmoto (1987b), Phys. Rev. B ,5877.S¨uto, A. (1989), J. Stat. Phys. , 525. Suwa, T., S. Tomita, N. Hosoito, and H. Yanagi (2017), Ma-terials , 1209.Tanese, D., E. Gurevich, F. Baboux, T. Jacquemin,A. Lemaˆıtre, E. Galopin, I. Sagnes, A. Amo, J. Bloch, andE. Akkermans (2014), Physical Review Letters (14),146404.Tang, C., and M. Kohmoto (1986), Physical Review B (3),2041.Theodorou, G., and M. H. Cohen (1976), Physical Review B (10), 4597.Thiem, S. (2015a), Phil. Mag. (11), 1233.Thiem, S. (2015b), Philosophical Magazine (11), 1233.Thiem, S., and M. Schreiber (2011), Eur. Phys. J. B , 415.Thiem, S., and M. Schreiber (2012), Phys. Rev. B , 224205.Thiem, S., and M. Schreiber (2013), J. Phys. Cond. Mat. ,075503.Varma, V. K., S. Pilati, and V. E. Kravtsov (2016), Phys.Rev. B , 214204.Varma, V. K., and M. ˇZnidariˇc (2019), Phys. Rev. B , 085105.Velhinho, M. T., and I. R. Pimentel (2000), Phys. Rev. B , 1043.Verbin, M., O. Zilberberg, Y. E. Kraus, Y. Lahini, and Y. Sil-berberg (2013), Phys. Rev. Lett. , 076403.Verbin, M., O. Zilberberg, Y. Lahini, Y. E. Kraus, and Y. Sil-berberg (2015), Phys. Rev. B , 064201.Vidal, J., D. Mouhanna, and T. Giamarchi (1999), Phys.Rev. Lett. , 3908.Vidal, J., D. Mouhanna, and T. Giamarchi (2001), Phys.Rev. B , 014201.Yuan, H. Q., U. Grimm, P. Repetowicz, and M. Schreiber(2000), Phys. Rev. B , 15569.Zheng, W. M. (1987), Phys. Rev. A , 1467.Zhong, J., J. Bellissard, and R. Mosseri (1995), J. Phys.Cond. mat. , 3507.Zhong, J., and R. Mosseri (1995), J. Phys. Cond. Mat.7