Contents
 easter egg
 Energy diffusion and the butterfly effect in inhomogeneous Sachdev-Ye-Kitaev chains
 Yingfei Gu, Andrew Lucas and Xiao-Liang Qi Department of Physics, Stanford University, Stanford, CA 94305, USA
 yfgu@stanford.edu, ajlucas@stanford.edu, xlqi@stanford.edu July 27, 2019
 Abstract: We compute the energy diffusion constant D, Lyapunov time τ\textscl and butterfly velocity v\textscb in an inhomogeneous chain of coupled Majorana Sachdev-Ye-Kitaev (SYK) models in the large N and strong coupling limit. We find D⩽v2\textscbτ\textscl from a combination of analytical and numerical approaches. Our example necessitates the sharpening of postulated transport bounds based on quantum chaos.

## 1 Introduction

A few years ago it was noted that many experimentally realized “strange metals” seem to be characterized by Drude “transport time” of order [1]. As this time scale was proposed to be the “fastest possible” time scale governing quantum dynamics [2, 3], it was conjectured by Hartnoll [4] that these strongly interacting strange metals remained metallic due to a fundamental bound on diffusion: . In strongly interacting systems without quasiparticles, many-body quantum chaos provides a natural velocity and time scale for such a diffusion bound [5, 6]. In some large- models, it is natural to define a Lyapunov time and butterfly velocity as [7, 8]

 ⟨V(x,t)W(0,0)V(x,t)W(0,0)⟩T∼1−1Nexp[tτ\textscl−xτ\textsclv\textscb], (1)

for rather general Hermitian operators and . is an analogue of the Lyapunov time from classical chaos: it describes the rate at which quantum coherence is lost. Similarly, is called the butterfly velocity: it governs the speed of chaos propagation, and is somewhat analogous to a state-dependent Lieb-Robinson velocity [9]. Since we now have a velocity and time scale which can be defined in any strange metal, [5, 6] noted that was naturally related to in some simple holographic settings. A natural question to ask is whether Hartnoll’s conjecture becomes

 D≳v2\textscbτ\textscl. (2)

Originally (2) was observed to hold for charge diffusion constant [5], with theory-dependent O(1) prefactors, but there are now multiple known counterexamples [10, 11, 12]. It is more compelling that (2) hold for the energy diffusion constant: as argued in [11, 13, 14], characterizes the loss of quantum coherence, a process related to quantum “phase relaxation” which should also characterize energy fluctuations and diffusion. Furthermore, additional evidence for an energy diffusion bound of the form (2) has arisen in holographic models in the low temperature limit [6, 15], and models of Fermi surfaces coupled to gauge fields [13]: at weak coupling, [16, 17] have also proposed a relation between diffusion and chaos.

Much of the recent literature focuses on “homogeneous” models of disorder: these can crudely be thought of as models where momentum is not a conserved quantity (hence, by Noether’s Theorem, microscopic translation symmetry has been broken), yet the effective equations governing transport remain spatially homogeneous. However, transport coefficients can be sensitive to how translation symmetry has been broken [14], so it is important to test the robustness of any transport bound in inhomogeneous models. Such a test was performed for charge diffusion in a family of holographic models [10], and the inequality in (2) was found to be reversed. In this paper, we test (2) for the energy diffusion constant in an inhomogeneous system: an inhomogeneous analogue of a Sachdev-Ye-Kitaev (SYK) chain of Majorana fermions. We will describe this solvable model of a “strange metal” without quasiparticles in more detail in Section 2. Our main result is that in this model the energy diffusion constant is upper bounded by chaos:

 D⩽v2\textscbτ\textscl. (3)

We will prove this in the limit where the inhomogeneity is parametrically slowly varying, and provide examples where the ratio is arbitrarily small, in Section 3.

Thus, we do not expect (2) to be generically true in disordered strange metals: holding as a strict inequality up to a finite O(1) prefactor which may be theory-dependent111For example, [13] found in their model, in a clean theory. In the SYK chain that we study, the numerical prefactor is 1, and in SYK-like holographic models the constant ranges from 0.5 to 1 [15]. Examples where the prefactor can be arbitrarily small can be found in [6]. but robust to disorder. Furthermore, in a generic, nearly translation invariant dimensional field theory without a global symmetry, the natural diffusion constant is the diffusion of energy. As noted in [6], this diffusion constant will be parametrically large due to the fact that translations are only weakly broken [14]. Hence, we now have examples of strange metals where is either much larger or much smaller than . It is still an interesting open question whether transport properties of strange metals are related to quantum chaos. Any such relation may be restricted to particular diffusion constants and/or models. For example, it may be the case that diffusion constants of translation invariant field theories are related to chaos [13, 15]. We hope that variational methods, related to those developed in [18, 19, 20] for hydrodynamic and holographic models, may be useful in providing rigorous bounds on transport and chaos in disordered strange metals.

For the remainder of the paper, we set .

## 2 The Inhomogeneous SYK Chain

The SYK model is a strongly interacting large- model in spacetime dimensions. It was introduced a long time ago as a model of disordered quantum magnets [21, 22]; it was revived more recently [23] due to its possible connection to holography [24, 25], which is a toy model for quantum gravity [26, 27, 28, 29, 30, 31]. Although it has since been shown that the SYK model does not admit a simple holographic dual [32], it does share many fascinating properties of holographic theories, including being “maximally chaotic” [33].

The model that we introduce is a generalization of the SYK chain model developed in [34].222See [35, 36] for another way of adding a spatial dimension to the SYK model. Consider a one-dimensional lattice of sites, where on each lattice site there exist Majorana fermions (), obeying the standard commutation rules . The Hamiltonian of these fermions is

 H=L∑x=1⎛⎝∑i

where the couplings and are all assumed to be independent Gaussian random variables drawn from a distribution with zero mean and following variance:

 ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯Jijkl,xJi′j′k′l′,x=3!N3J20,x\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}ii′\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}jj′\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}kk′\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}ll′% \mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}xy,¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯J′ijkl,xJ′i′j′k′l′,x=1N3J21,x\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}ii′\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}jj′\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}kk′\mbox% {   {1 0 -.18 1}\hbox to 0.0pt{δ}}ll′\mbox{   {1 0 -.% 18 1}\hbox to 0.0pt{δ}}xy, (5)

The important (and only) difference comparing to [34] is that we do not assume the variances and take the same value for each . We are interested in the thermodynamic limit .

One can show that the replica-diagonal333Off-diagonal sectors in replica space do not contribute at the orders in that we study. partition function, at inverse temperature , can be written as a path integral over two bilocal fields:

 ¯¯¯¯Z=∫DGDΣexp[−NSeff[G,Σ]] (6)

with the Euclidean time action

 Seff.=L∑x=1 ⎧⎪⎨⎪⎩−logPf(∂τ−Σx)+12β∫0d2τ[ΣxGx−J20,x4G4x−J21,x4G2xG2x+1]⎫⎪⎬⎪⎭ (7)

The “Green’s function” and “self energy” are functions of two time variables. The product in above formula is abbreviation for ; and are similar products. and show up to exactly quadratic order in (7) because the random couplings and were Gaussian random variables. As the manipulations in this section are essentially identical to [32, 34], we only present the few steps where they differ in an important way (due to the absence of translational invariance on average).

It is convenient to rewrite the interaction term (-term) into the following form:

 ∑x(J20,xG4x+J21,xG2xG2x+1) = ∑x{(J20,x+J21,x+J21,x−12)G4x+12G2x[J21,x(G2x+1−G2x)+J21,x−1(G2x−1−G2x)]} (8)

If one chooses and such that for each ,

 J2≡J20,x+J21,x+J21,x−12 (9)

is a constant independent of , then the effective on-site coupling is easily seen to be -independent. The saddle point equations of become

 G−1x(τ1,τ2) =Σx(τ1,τ2)−\mbox{   {1 0 -.18 1}\hbox t% o 0.0pt{δ}}′(τ1−τ2), (10a) Σx(τ1,τ2) =J20,xGx(τ1,τ2)3+J21,xGx+1(τ1,τ2)2+J21,x−1Gx−1(τ1,τ2)22Gx(τ1,τ2), (10b)

and they admit an -independent approximate solution: , with

 Gs(τ) =bΔ(βJ\mbox{   {1 0 -.18 1}\hbox t% o 0.0pt{π}}sin\mbox{   {1 0 -.18 1}\hbox to 0.0pt{π}}τβ)−2Δ,0⩽τ<β (11a) b =1\mbox{   {1 0 -.18 1}\hbox to 0.0pt{π}}(12−Δ)tan(\mbox{   {1 0 -.18 1}\hbox to 0.0pt{% π}}Δ),Δ=14,

which becomes exact at (conformal) limit. The system also has a uniform specific heat per site . Thus, as in [34], this saddle point is identical to the -dimensional SYK model of [23, 32] at coupling constant . If the choice (9) is not made, then the saddle point equations do not admit a homogeneous solution, and it is unclear what the effective theory is.

In the strong coupling limit, , and long wavelength limit, the physics of interest to us is governed by the fluctuations induced by reparametrization modes , which act as . To quadratic order of the infinitesimal fluctuations, and leading order in expansion, the effective action for the fluctuations has a simple form in Fourier space: defining , , we find

 Seff.=1256\mbox{   {1 0 -.18 1}% \hbox to 0.0pt{π}}∑xy∑nϵn,x|n|(n2−1)(α|n|βJ\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}xy+Cxy)ϵ−n,y, (12)

where all -dependence is contained in the tridiagonal matrix

 Cxy=13J2⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝⋱ −J21,x−1 0 0−J21,x−1 J21,x−1+J21,x −J21,x 00 −J21,x J21,x+J21,x+1 −J21,x+10 0 −J21,x+1 ⋱⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠, (13)

and is a constant determined by numerics [32]. A few more steps of this derivation are contained in Appendix A. The long wavelength limit alluded to earlier is the regime when the eigenvalues of are not larger than , which (as we will see) do exist even for the disordered matrix. The derivation of this effective action is identical to [34]: in this previous work, was translation invariant and so (12) was written in momentum space, where the matrix becomes diagonal.

By writing in the form

 Cxy =DTxzΛzwDwy =13J2⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝⋱ −1 0 00 1 −1 00 0 1 −10 0 0 ⋱⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠T⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝⋱ 0 0 00 J21,x 0 00 0 J21,x+1 00 0 0 ⋱⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝⋱ −1 0 00 1 −1 00 0 1 −10 0 0 ⋱⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠ (14)

we immediately recognize that it is positive definite and can be interpreted as the first-order finite-difference discretized version of the differential operator

 Cxy∼13J2(−ddxJ1(x)2ddx)discretized. (15)

The interpretation of as an approximate differential operator becomes exact when varies slowly. Letting denote spatial averages over , suppose that

 E[J21,xJ21,y]∼f(|x−y|M) (16)

with a (large) integer, and a non-zero function for O(1) argument. To leading order in , the low-lying spectrum of the discrete operator will be identical to the continuum differential operator (15).

In our model, we will take to be an arbitrary function of , simply constrained to (otherwise as defined in (9) would be negative). The properties of the matrix will then depend on the inhomogeneity that we encode through -dependent .

## 3 Diffusion and Chaos

From the effective action (12), we are able to extract the thermal response functions. The procedure is identical to [34] and a diffusion pole is found in the energy density () two-point function:

 ⟨Tttx,nTtty,−n⟩T∼(|ωn|\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}xy+2πJαCxy)−1≡(|ωn|\mbox{   {1 0 -.18 1}% \hbox to 0.0pt{δ}}xy+˜Cxy)−1. (17)

Upon proper analytic continuation to real time, we interpret (17) as having diffusive poles (on the negative imaginary axis) whenever is an eigenvalue of . is analogous to a tight-binding-model hopping matrix. If commutes with a discrete translation operator, then we expect plane wave eigenstates, the lowest-lying of which will have an eigenvalue in a chain of length . If is random, then strictly speaking all eigenstates of at fixed are localized in the continuum. However, because is analogous to a discretized differential operator (15) which has an exact delocalized zero mode, the low-lying spectrum of will look diffusive on length scales . In other words, the localization length grows faster than [37], and the smallest nontrivial eigenvalue scales as in this case as well. Hence, the diffusion constant is finite.

In fact, the lowest-lying non-trivial eigenvectors of are well approximated by plane waves:

 ux∼eiqx,q=±2\mbox{   {1 0 -.1% 8 1}\hbox to 0.0pt{π}}L, (18)

which can be verified numerically. See Appendix B for more comments on this equation. The eigenvalue of such a will be , with an effective diffusion constant. In the large limit, we may compute by solving the following differential equation as :

 (19)

The constant is computed in Appendix B:

 1D=E[1D(x)]. (20)

This equation has a straightforward physical interpretation. Because the specific heat in our model is -independent to leading order [34], is proportional to the thermal conductivity. One can approximate our inhomogeneous chain by joining together homogeneous SYK chains of length . Within each of these chains, the thermal conductivity is proportional to the diffusion constant . Joining together these segments of length , we find a resistor network: hence, the thermal resistivity spatially averages. This leads to (20). In Appendix B, we argue that (20) is applicable even beyond the large limit, under some assumptions which work relatively well in practice numerically, so long as finite size effects are small.

Now we study the butterfly velocity, defined by out-of-time-ordered correlation functions of spatially separated operators. In order to extract , we study the (properly regularized) connected out-of-time-ordered correlation function. One finds [34], in the region that

 1N2∑i,j⟨χi,x(t),χj,y(0)χi,x(t),χj,0(0)⟩T,connected∼e2% \mbox{   {1 0 -.18 1}\hbox to 0.0pt{π}}βt(2\mbox% {   {1 0 -.18 1}\hbox to 0.0pt{π}}β+˜Cxy)−1. (21)

Comparing to (1), we observe that in this model, as in the usual SYK model,

 τ\textscl=β2\mbox{   {1 0 -.18 1}\hbox to 0.0pt{π}}. (22)

This matrix inverse is the discrete analogue of the Green’s function

 2\mbox{   {1 0 -.18 1}\hbox to 0.0pt{π}}βG(x;y)−ddx(D(x)dG(x;y)dx)=\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}(x−y). (23)

In the long range disorder limit, when at each point

 M≫√2\mbox{   {1 0 -.18 1}\hbox to 0.0pt{π}}βD(x), (24)

the solution of this equation is exponentially decaying [10]:

 G(x;y)∼e−|x−y|/v\textscbτ\textscl (25)

with

 1v\textscb=√2\mbox{   {1 0 -.18 1}\hbox to 0.0% pt{π}}βE[1√D(x)]. (26)

This equation can be derived by noting that, away from the points , we can change coordinates from to , defined by . One then finds

 2\mbox{   {1 0 -.18 1}\hbox to 0.0pt{π}}βD(s)G=d2Gds2, (27)

and when varies slowly, one can straightforwardly write

 G=exp[−∫ds√2\mbox{   {1 0 -.18 1}\hbox to 0% .0pt{π}}D(s)β]=exp[−∫dx√2% {1 0 -.18 1}\hbox to 0.0pt{π}D(x)β]. (28)

Hence, we find that and are not equal. Using the Cauchy-Schwarz inequality it is straightforward to conclude

 v2\textscbτ\textscl=E[1√D(x)]−2⩾E[1D(x)]−1=D. (29)

The physics at play here is essentially the same as in holography, where charge diffusion was shown to obey a similar inequality, for the same reasons [10].

We have numerically computed and in SYK chains of finite length with periodic boundary conditions. is found by averaging the two smallest non-vanishing eigenvales of . is found by computing the typical value of for .444We must also normalize the value of found numerically by a factor very close to unity, to account for finite size effects. This factor depends only on . So long as , we find that agrees with the “resistor chain” prediction (20) for any , so long as , to within about 0.1% residual error (which is possibly a numerical finite size effect). Indeed, the derivation of (20) in Appendix B does not rely on the assumption that is slowly varying, so this is not surprising.

As shown in Figure 1, we see that while agrees very well with the “hydrodynamic” theory, agrees with the hydrodynamic theory at large , while partly approaching from above as . This behavior is not surprising: as becomes shorter, the four-point function (21) begins to “self-average” over the inhomogeneity in a manner analogous to diffusion. Nevertheless, holds as a strict inequality in the inhomogeneous systems that we have studied numerically, even when (no correlations among . Our violation of the relation is not limited to the regime of long range disorder.

So far, the discrepancies between and are only on the order of a few percent in our numerical data. Yet (29) implies that there is no possible upper bound on diffusion due to quantum chaos. Defining a natural probability measure on as

 p(X)dX≡E[Θ(X+dX−D(x))Θ(D(x)−X)], (30)

we estimate that if

 p(X→0)∼Xa,with −12

then but .555Since our analytic calculation of requires (24), and disorder is correlated over sites, we estimate that in a chain of length the minimal value of scales as . Requiring that requires that , up to some dimensionless constant. Hence, we conclude that the inhomogeneous SYK chain with but finite only strictly exists in a somewhat subtle thermodynamic limit with and taken to simultaneously, making sure to obey . We have looked for this parametric breakdown of the relationship between and in smaller chains with .

As shown in Figure 2, we see qualitative agreement (but quantitative disagreement) with our hydrodynamic predictions for and (accounting for finite size effects). As , we observe that the ratio becomes dependent on the length of the chain, and decreases for the longer chain. This provides evidence that even when , this inhomogeneous SYK chain may have but in the thermodynamic limit.

Finally, let us comment on the low temperature limit (while, of course, taking as well such that ). At small enough temperature, keeping the inhomogeneity fixed, (24) will break down. Figure 2 suggests that the large limit is not required to obtain substantial deviations from . A more interesting subtlety that arises at very low temperature is the difference between periodic inhomogeneity and random inhomogeneity. For periodic inhomogeneity, one can diagonalize the matrix as a periodic tight-binding hopping matrix in a larger unit cell, and at low enough temperatures, the solution to (23) can be well-approximated by considering only the physics of the lowest band. In this regime, one will recover . As the period of the periodic inhomogeneity grows longer, the temperature above which decreases. We expect that for random inhomogeneity (where the eigenstates do not form bands, but are in fact localized) one finds at all finite temperatures, and provide some numerical evidence for this in Figure 3.

## 4 Outlook

We have presented a modification of the SYK chain model of [34], in which there is an upper bound on the diffusion constant: . As we pointed out in the introduction, this suggests that there is no (simple) generic bound relating transport and quantum chaos in all strange metals.

One might ask whether our violation of (2) could be found in a “homogeneous” model that does not rely on explicit translation symmetry breaking in the low energy effective description. In the SYK chains that we have studied, this is easy to accomplish at leading order in , because there is no difference betwen averages over annealed disorder vs. quenched disorder. Hence, consider the partition function

 Z=∫DJ21,xP[J21,x]ZSYK[J21,x], (32)

with the partition function defined in (6), and a translation invariant function where have support on some finite domain. We may think of as either a partition function for the “slow” dynamical variables , or as accounting for certain correlated non-Gaussian fluctuations in the random variables and of the microscopic Hamiltonian (4).666In this latter case, and are no longer completely independent random variables. The reason for this is that after integrating over on a single site, we generically find . Since the joint probability distribution of the does not factorize, we conclude that these random variables are no longer independent. Hence, integrating over fluctuations in restores translation invariance, but necessarily introduces non-trivial correlations between the and . So long as for each choice of , by linearity, this inequality will remain true even in the homogeneous model (32) after averaging over the ensemble of random couplings.

## Acknowledgements

We thank Mike Blake for helpful comments on a draft of this paper. AL was supported by the Gordon and Betty Moore Foundation’s EPiQS Initiative through Grant GBMF4302. YG and XLQ are supported by the National Science Foundation through the grant No. DMR-1151786, and by the David & Lucile Packard Foundation.

## Appendix A Derivation of Effective Action

It is convenient to expand about the saddle using renormalized variables defined by

 Gx(τ1,τ2) =Gs(τ1,τ2)+|Gs(τ1,τ2)|−1gx(τ1,τ2), Σx(τ1,τ2) =Σs(τ1,τ2)+|Gs(τ1,τ2)|σx(τ1,τ2), (33)

where we have rescaled the fluctuation fields by prefactors and . It should be noticed that although the saddle point is uniform in space and translation invariant in time, the fluctuation fields have generic space-time dependence.

Now we expand the effective action to second order in the fluctuation fields , which leads to

 Seff.[g,σ]≈Sseff −14∫d4τ∑xσx(τ1,τ2)Gs(τ13)⋅|Gs(τ34)|⋅Gs(τ42)⋅|Gs(τ21)|σx(τ3,τ4) +∫d2τ(∑x12σx(τ1,τ2)gx(τ1,τ2)−3J24∑x,ygx(τ1,τ2)Sxygy(τ1,τ2)). (34)

The spatial kernel is a tight-binding hopping matrix

 Sxy=\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}xy+13J2[(−J21,x−J21,x−1)\mbox{   {1 0 -.18 1% }\hbox to 0.0pt{δ}}xy+\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}x,y−1J21,x+\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}x,y+1J21,x−1]=\mbox{   {1 0 -.18 1}\hbox to 0.0% pt{δ}}xy−Cxy,, (35)

with defined in (13). Next we integrate out and obtain a quadratic action for alone. We define as the (symmetrized) four-point function kernel of the SYK model:

 ˜K(τ1,τ2;τ3,τ4)=3J2Gs(τ13)⋅|Gs(τ34)|⋅Gs(τ42)⋅|Gs(τ21)|. (36)

We have defined . The effective action of is therefore

 Seff.[g]=Sseff+3J24∫d4τ∑x,ygx(τ1,τ2)[˜K−1(τ1,τ2;τ3,τ4)\mbox{%    {1 0 -.18 1}\hbox to 0.0pt{δ}}xy−Sxy\mbox{   {1 0 -.18 % 1}\hbox to 0.0pt{δ}}(τ13)\mbox{   {1 0 -.18 1}\hbox to 0.0% pt{δ}}(τ24)]gy(τ3,τ4), (37)

The leading contribution comes from the soft modes which can be identified as induced by reparametrization . They can be interpreted as energy fluctuations. We can linearized these modes and write them in the fourier space: , . For such modes, we know from [32] that the kernel obeys

 (˜K−1−1)ϵn,x=α|n|βJϵn,x (38)

Following [32], we now find the following effective action for the linearized modes given in (12). The additional numerical prefactors in (12), including , come from properly normalizing the eigenvector .

## Appendix B Eigenvalues of Inhomogeneous Diffusion Matrix

The following argument is reminiscent of [10]. Let us define the matrix

 Qxy=\mbox{   {1 0 -.18 1}\hbox to 0.0pt{δ}}xyeiqx, (39)

and consider the limit . We postulate that the lowest lying eigenvalues of are proportional to :

 ˜Cxyuy=Deffq2ux. (40)

Multiplying on both sides by we obtain

 Qxz(DTSDQ−1)zw(Qu)wy=Deffq2(Qu)x. (41)

We now look for a series solution to this equation of the form

 Qu=u0+iqu1−q2u2+⋯. (42)

Such a series expansion is reasonable – we expect in any finite chain that the lowest few eigenvectors are delocalized, and have confirmed this numerically. At , the matrix is simply the identity, and so we must take

 (u0)x=1. (43)

At , we find the equation

 Qxz[(DTSDQ−1)zw(u0+iqu1)w]=0. (44)

is invertible, and hence we conclude that if is non-trivial, the left-most must act on a non-trivial vector. Because has a one-dimensional kernel we conclude

 c(u0)z=(SDQ−1)zw(u0+iqu1)w (45)

We may fix the constant as follows. First left-multiply by , and keep only first order terms in , to obtain

 cS−1zw(u0)w=(DQ−1)zw(u0)w+iqDzw(u1)w≈−iqQ−1zw(u0)w+iqDzw(u1)w. (46)

In the second equality, we have exploited the fact that is an eigenvector of of eigenvalue . Now, we perform a non-rigorous sleight-of-hand: at leading order in , we may treat . We may then left-multiply by , to remove only the second term of (46), and obtain

 c=−iqu0⋅u0u0⋅S−1u0. (47)

On physical grounds, this is the statement that the eigenvector looks a lot like a plane wave, and this seems to be true numerically. We now have the corrections to the eigenvector , and to leading order use the variational principle to ‘exactly’ compute the effective diffusion constant. We find

 (u0+iqu1)⋅DTSD⋅(u0+iqu1)=q2c2u0⋅S−1u0=Deffq2u0⋅u0 (48)

which gives us, for a periodic chain of sites:

 1Deff=1L∑1Dx. (49)

One way to rigorously obtain (46) is to extend a chain of length into an infinite chain by tiling the same repeatedly: . Using the discrete translation symmetry and choosing the appropriate definition of , one can demand , and then sum over only sites in (46), killing the right-most term. In practice, we have often found this tiling to be unnecessary numerically.