Local Criticality, Diffusion and Chaos in Generalized Sachdev-Ye-Kitaev Models

Download as pdf or txt
Download as pdf or txt
You are on page 1of 38

Prepared for submission to JHEP

Local criticality, diffusion and chaos in generalized


Sachdev-Ye-Kitaev models
arXiv:1609.07832v2 [hep-th] 29 May 2017

Yingfei Gu,1 Xiao-Liang Qi1 and Douglas Stanford2


1
Department of Physics, Stanford University,
382 Via Pueblo Mall, Stanford, CA, USA
2
School of Natural Sciences, Institute for Advanced Study,
1 Einstein Dr, Princeton, NJ, USA

E-mail: [email protected], [email protected], [email protected]

Abstract: The Sachdev-Ye-Kitaev model is a (0 + 1)-dimensional model describing Ma-


jorana fermions or complex fermions with random interactions. This model has various
interesting properties such as approximate local criticality (power law correlation in time),
zero temperature entropy, and quantum chaos. In this article, we propose a higher di-
mensional generalization of the Sachdev-Ye-Kitaev model, which is a lattice model with
N Majorana fermions at each site and random interactions between them. Our model can
be defined on arbitrary lattices in arbitrary spatial dimensions. In the large N limit, the
higher dimensional model preserves many properties of the Sachdev-Ye-Kitaev model such
as local criticality in two-point functions, zero temperature entropy and chaos measured by
the out-of-time-ordered correlation functions. In addition, we obtain new properties unique
to higher dimensions such as diffusive energy transport and a “butterfly velocity” describ-
ing the propagation of chaos in space. We mainly present results for a (1 + 1)-dimensional
example, and discuss the general case near the end.

ArXiv ePrint: 1609.07832


Contents

1 Introduction 1

2 Review of the SYK model 3

3 The generalized SYK model 5


3.1 Definition of the chain model 5
3.2 The effective action and the saddle point 7
3.3 Fluctuations of the collective fields and four-point functions 10
3.4 Symmetry breaking and pseudo-Goldstone mode 12

4 The OPE region and the Energy Transport 13


4.1 The OPE at finite p 14
4.2 The h = 2 contribution for small p 15
4.3 Energy transport 16

5 Chaos and the butterfly velocity 19

6 General construction and higher dimensions 22


6.1 Generalized SYK model on higher dimensional lattices 23
6.2 Models with global symmetry 25
6.3 General q 27

7 Conclusion and discussion 27

A Diagrammatic derivation for four-point functions 28

B Summation trick and the prefactor 32

C Diffusion and the butterfly velocity in general dimensions 33

1 Introduction

Holographic duality, also known as the anti-de-Sitter space/conformal field theory (AdS/CFT)
correspondence, refers to the duality between quantum many-body systems in d spatial di-
mensions and gravitational theories in (d + 1)-dimensional asymptoticaly anti-de Sitter
geometries [1–3]. Various pieces of evidence suggest that holographic duality is a generic
phenomenon that applies beyond the super Yang-Mills theories where the original conjec-
ture was proposed. However, it is difficult to find concrete models for which the duality can
be verified directly by comparing bulk and boundary calculations. A very interesting devel-
opment in this direction is the Sachdev-Ye-Kitaev model[4], which is a (0 + 1)-dimensional

–1–
model describing random interactions between N Majorana fermions. When N is large
and the temperature is low, the model has an emergent approximate time reparameteriza-
tion symmetry which is weakly broken. In this limit, two-point functions and four-point
functions can be computed. The results suggest that the model has a weakly coupled
holographic dual, which includes dilaton gravity in an approximate AdS2 geometry[5–10],
weakly coupled to an infinite number of matter fields[11].
The SYK model is a modification of a quantum spin model proposed by Sachdev and
Ye more than 20 years ago[12] (which was also related to holographic duality in Ref. [13]).
In the large N and low temperature limit N  βJ  1 (with β the inverse temperature
and J the average coupling strength), the behavior of the model is controlled by a large
N saddle point, with the fermion two-point function G(τ, τ 0 ) = N1 j hχj (τ )χj (τ 0 )i playing
P

the role of a semi-classical “order parameter” with small fluctuations suppressed by N1 . At


low temperature, G(τ, τ 0 ) has a power law dependence on τ − τ 0 in the infared, suggesting
that the low energy dynamics of this model might be conformally invariant [14]. In fact,
the (0 + 1)-d conformal symmetry is only approximate, and in fact the low-temperature
dynamics are dominated by the specific way in which the symmetry is broken [4, 11].
A particularly interesting type of four-point function is the out-of-time-order corre-
lation function (OTOC) F (t) = hχj (t)χk (0)χj (t)χk (0)i, which has been proposed as a
measure of chaos in quantum systems[15–20] (see Ref. [21–24] for experimental propos-
als for measuring OTOC). Physically, the decrease of this four-point function measures
the increase of the size1 of the operator anticommutator {χj (t), χk (0)}, which indicates
how sensitive the system is to an initial perturbation created by acting with the fermion
operator χj (0). The exponential time dependence of the connected part of the OTOC
F (t)conn. ∝ eλL t defines an inverse time scale λL which can be considered a quantum
analog of Lyapunov exponent. Ref. [25] proved a general upper bound λL 6 2π β (for a
regularized form of OTOC with imaginary time evolution), which is saturated for theories
with an Einstein gravity dual. Interestingly, the SYK model also saturates this upper
bound [4, 11, 26]. Many other aspects of the SYK model (and a similar model for complex
fermions) have been investigated recently[8, 10, 26–32].
Given the interesting properties of the (0 + 1)-dimensional SYK model, it is natural
to look for higher dimensional cousins. In this paper, we propose a family of higher di-
mensional variants of the SYK model, which remain solvable in the large N limit. Our
model can be defined on an arbitrary discrete lattice in arbitrary spatial dimensions. There
are N fermions on each site with an SYK Hamiltonian. Different sites have independent
SYK couplings, and neighboring sites are coupled by random four-fermion terms. Using
the same techniques as those applied to the SYK model, we can study two-point functions
and four-point functions in space-time. The spatial locality of the model allows us to study
transport properties and propagation of quantum chaos in space-time. We find that the
disorder-averaged two-point functions vanish between different lattice sites, and have the
same local critical behavior as in the SYK model within each site. Our model also has the
1
At infinite temperature, the “size” of the anticommutator simply means its 2-norm. At finite tem-
perature, it is generalized to the thermal expectation value h{χj (t), χk (0)}2 iβ . Note that we study the
anticommutator instead of the commutator because the operators are fermionic.

–2–
same zero temperature entropy at each site as the SYK model. Correlation between differ-
ent sites begins at the level of four-point functions. Similar to the SYK model case, one can
consider the fermion four-point function as being mediated by a series of collective fields,
with the leading contribution coming from energy fluctuations. In our generalized models,
the four-point function allows us to study the dynamics of collective fields in space-time.
In particular, we find a diffusive dynamics of energy density, which means this model de-
scribes a diffusive strongly correlated metal phase. The OTOC can also be studied, which
now has both spatial and temporal dependence. We show that at low temperature the
Lyapunov exponent still saturates the chaos bound 2π β . The propagation of quantum chaos
in space-time can be characterized by a “butterfly effect velocity” x = vB t, which means
that the anticommutator {χjx (t), χjy (0)} becomes significant at t = |x − y|/vB + (const).
Interestingly, the diffusion constant D and butterfly velocity vB in our generalized SYK
2
vB
model satisfy a simple relation2 D = 2πT , which realises a bound conjectured on diffusion
in incoherent metal and agrees with the holographic calculation on incoherent black hole
[35–37].
The remainder of the paper is organized as follows. In section 2 we will briefly review
the properties of the original SYK model. In section 3, we define the generalized SYK
model in higher dimensions. For concreteness, we work on a (1 + 1)-dimensional example
and study its correlation functions and thermodynamic properties in detail. In section 4
we study the dynamics of collective fields in this system based on an operator-product
expansion of fermion four-point functions. In particular, we show that one of the collective
fields describes the diffusion of energy in this system. In section 5, we study the OTOC
of the (1 + 1)-dimensional model and obtain the Lyapunov exponent and the butterfly
velocity. In section 6 we discuss the general form of the model in generic dimensions and
graphs. In section 7 we end the paper with a summary and discussion of further topics.

2 Review of the SYK model

In this section, we briefly review some basic facts about the SYK model. The SYK model
consists of N Majorana fermions χj , j = 1, 2, . . . , N with a random four-fermion interaction
[4] X
H= Jjklm χj χk χl χm , {χj , χk } = δjk (2.1)
16j<k<l<m6N

where {Jjklm } are independent random couplings with zero mean Jjklm = 0, and variance
of individual elements is given by 3!1 N 3 Jjklm
2 = J 2 . Here, the symbol J without indices
is a dimensionful constant that sets the scale of the Hamiltonian, and the factors 3!1 N 3
are for later convenience. The interaction is all-to-all, so that there is no spatial locality,
and this model should be considered as a (0 + 1)-d quantum mechanical system. The
model is solvable at large N , and exhibits holographic behavior at strong coupling N 
βJ  1. In this limit, the (finite temperature) two-point function[12, 14] G (τ1 , τ2 ) :=
2
See Ref. [33, 34] for other discussions on diffusion and butterfly effect in solid state systems.

–3–
1 PN
N j=1 hTτ χj (τ1 )χj (τ2 )iβ has emergent PSL2 (R) symmetry3

πτ12 −2∆
 
∆ βJ
G (τ1 , τ2 ) = b sin , 0 6 τ12 < β (2.2)
π β
 
1 1
b= − ∆ tan(π∆), ∆ = 1/4
π 2
where we denote τ1 − τ2 by τ12 for simplicity here and below. One can also compute the
four-point function in the holographic limit[4, 11]. The interesting piece of the four-point
function is the connected part, which begins at order N1
N
1 1 X
F(τ1 , τ2 , τ3 , τ4 ) := 2 hTτ χj (τ1 )χj (τ2 )χk (τ3 )χk (τ4 )iβ − G(τ1 , τ2 )G(τ3 , τ4 ). (2.3)
N N
j,k=1

At leading order of 1/βJ and 1/N , when we consider a configuration of imaginary (Eu-
clidean) times with the ordering τ1 > τ2 > τ3 > τ4 , then the correlator F factorizes:
" ! ! #
F(τ1 , τ2 , τ3 , τ4 ) 8 πτ12 πτ34
= βJ −1 − 1 + O(1/βJ) (2.4)
G(τ12 )G(τ34 ) παK β tan πτβ12 β tan πτβ34
where αK ≈ 2.852. The factorization indicates that the OPE of two fermion operators
is dominated by the conserved quantity in this model — the Hamiltonian itself, and the
fluctuation of energy determine the four point function in this configuration.
Quantum chaos can be diagnosed by the out-of-time-ordered correlator (OTOC). The
OTOC is calculated by starting with the Euclidean correlator with a “j-k-j-k” configuration
of times, e.g. τ1 > τ3 > τ2 > τ4 , and then giving the times τ1 and τ2 large real-time parts.
Making the particular choice that all four points are evenly spaced on the imaginary time
circle, we have [11]
 
3β β β 8 π 2π
F( + it, + it, , 0) = βJ (1 − cosh t) + O(1/βJ) (2.5)
4 4 2 παK 2 β
which exhibits an exponential
 growth
 for real time t greater than the dissipation time

td ∼ β, i.e., F ∼ −βJ exp β t for t  β. This growth was recognized as a signature
of chaos[15, 16, 18], moreover, the corresponding growth exponent λL = 2π β saturates the
chaos bound proposed in Ref. [25].
The two-point function and four-point function shown above can be calculated via
standard Feynman diagrams with large N simplification. Another approach is to use
the disorder-averaged effective action. After introducing a pair of auxiliary bi-local fields
“Green’s function” G(τ1 , τ2 ) and “self-energy” Σ(τ1 , τ2 ) and carrying out the disorder av-
erage of the partition function, the fermions can be integrated out, leaving [4, 38]:
Z
Z = DGDΣ exp (−N Seff [G, Σ])

J2
Z  
1 4
Seff [G, Σ] = − log Pf (∂τ − Σ) + dτ1 dτ2 Σ(τ1 , τ2 )G(τ1 , τ2 ) − G(τ1 , τ2 ) . (2.6)
2 4
!
3 af +b πτ a b
This symmetry acts as f → cf +d
where f = tan β
, and ∈ SL2 (R).
c d

–4–
The large N prefactor implies that the problem is essentially classical. The saddle point
reproduces the Schwinger-Dyson equations that can also be derived via Feynman diagrams:
1
G(iω) = , Σ(τ ) = J 2 G(τ )3 (2.7)
−iω − Σ(iω)

where we assumed time translation symmetry to simplify the equations.


One can also obtain the connected part of the four point function by considering
quantum fluctuations around the saddle point. Among all the quantum fluctuations, there
is a special class induced by the reparametrization of the time circle f ∈ Diff(S 1 ), which
contributes the leading piece at strong coupling βJ  1. The effective description of
this part turns out to be well-captured by a local action proportional to the Schwarzian
derivative[4, 11]. Remarkably, the same form of effective action also appears in the AdS2
Einstein-dilaton theory[5–10].
One can also derive the thermodynamic properties from the large-N saddle point free
energy:

J2
 Z  
F 1 1 4
= − log Pf (∂τ − Σ) + dτ1 dτ2 Σ(τ1 , τ2 )G(τ1 , τ2 ) − G(τ1 , τ2 ) (2.8)
N β 2 4
γ
= U − S0 T − T 2 + . . . (2.9)
2
In the second line we write the free energy in a low temperature expansion,4 where U ≈
−0.0406J is the ground state energy, S0 ≈ 0.232 is the zero temperature entropy[4, 38],
and γT = cv = 16πα √K
2βJ
≈ 0.396
βJ is the specific heat[11]. The entropy term can be derived
by inserting the conformal saddle point solution (2.2) in the effective action. The specific
heat can be derived from knowledge of the leading (in 1/βJ) correction to the conformal
saddle, but the energy requires the exact (numerical) finite βJ solution.

3 The generalized SYK model

In this section, we will present a simple way to generalize the SYK model to higher di-
mensions while keeping the solvable properties of the model in the large-N limit. For
concreteness of the presentation, in this section we focus on a (1 + 1)-dimensional example,
which describes a one-dimensional array of SYK models with coupling between neighbor-
ing sites. It should be clear how to generalize, and we will discuss more details of the
generalization to arbitrary dimensions and generic graphs in section 6.

3.1 Definition of the chain model


In (1 + 1)-d, our model describes a coupled array of SYK model sites, as shown in figure 1.
Each site containins N  1 Majorana fermions with SYK interactions drawn independently
for each site. Each pair of neighboring sites are then further coupled via random four
4
Starting at T 3.77 , this expansion is expected to also involve non-integer powers given by the dimensions
of irrelevant operators in the model.

–5–
0
Jjklm
Jjklm

j m j l
k l
k m

Figure 1: A chain of coupled SYK sites: each site contains N  1 fermion with SYK
interaction. The coupling between nearest neighbor sites are four fermion interaction with
two from each site.

fermion interactions with two of the fermions from each site. The Hamiltonian has the
following form:
 
M
X X X
0
 
H= 
 Jjklm,x χj,x χk,x χl,x χm,x + Jjklm,x χj,x χk,x χl,x+1 χm,x+1 

x=1 16j<k<l<m6N 16j<k6N
16l<m6N
(3.1)

where x labels the lattice sites, and {χj,x }j=1,2,...,N ; x=1,2,...,M are the Majorana fermion op-
erators satisfying anti-commutation relations and periodic boundary condition: {χj,x , χk,y } =
δxy δjk , χj,0 ≡ χj,M . N is the number of Majorana fermions on each site and M is the num-
ber of the sites, or equivalently, the length of the chain. In this expression, we restrict the
range of indices in the sum such that each term only appears once. The first term describes
the on-site SYK interaction, while the second term is the nearest neighbor random four
fermion coupling. The random couplings {Jjklm,x } and {Jjklm,x 0 } are drawn independently
for each value of x, from a distribution with zero mean and variance defined in the following
way:
1 3 2
0
Jjklm,x = Jjklm,x = 0, N Jjklm,x = J02 , N 3 Jjklm,x
02 = J12 . (3.2)
3!
The normalization, especially the factors of N , are chosen to make the large N limit uniform
and ensures the dimension 1 coupling constants J0 and J1 represent the average strength
of the thermal bath seen by each fermion field.
Comparing to the original SYK model, our model clearly has spatial locality. One
can view our model as either M coupled SYK sites, or equivalently as a big SYK site
with N M Majorana fermions but with inhomogeneous coupling strength—the non-local
couplings (between sites |x − y| > 1) are suppressed. As will be discussed in section 6, it is
not essential to have the 2-2 coupling between neighboring sites. Introducing more generic
couplings such as χj,x χk,x+1 χl,x+2 χm,x+3 is straightforward as long as the coefficients of
these terms are all independent variables. However, in this section we will focus on the 2-2
coupling case for simplicity.
A first observation of the model is that the Hamiltonian doesn’t contain any quadratic
term, so that the free propagator is diagonal not only in flavor indices but also in spatial

–6–
j, x j, x

Figure 2: The leading order diagrams only connect fermions with same flavor and spatial
coordinate under random average of disorder fields (dashed line).

coordinate hTτ χj,x (τ1 )χk,y (τ2 )ifree = 12 δjk δxy sgn(τ12 ). Furthermore, at leading order the
interaction vertices under random average only contribute to the diagonal part, see figure 2.
Therefore, the dressed Green’s function hTτ χj,x (τ1 )χk,y (τ2 )i is also diagonal in both flavor
and spatial coordinates 5 . In addition, the diagonal Green’s functions hTτ χj,x (τ1 )χj,x (τ2 )i
are independent of j due to an SO(N ) symmetry of the model after averaging over disorder.
Therefore we have hTτ χj,x (τ1 )χk,y (τ2 )i = δjk δxy N1 N
P
l=1 hTτ χlx (τ1 )χlx (τ2 )i.

3.2 The effective action and the saddle point

α α

α α

β β
(a) Replicon diagonal ∼ N

(b) Off-diagonal ∼ 1/N 2

Figure 3: Replicon diagonal v.s. off-diagonal contributions to the partition function Z n :


solid lines connects same replica index. Different replica indices can be connected only by
dashed lines, which are disorder fields. The replicon off-diagonal diagram on the right is
suppressed by 1/N 3 compared to the diagonal diagram on the left.

In this section, we employ the large N effective action approach to analyze the model.
In principle, to study the quenched problem, one should introduce n replicas and analyze
the disorder averaged partition function Z n . However, as in the original SYK model,
our chain model self-averages at large N . One can verify this by checking the replicon
off-diagonal contribution to the averaged partition function Z n . Figure (3) shows the
leading replicon off-diagonal diagram, which is suppressed by 1/N 3 relative to the connected
diagonal diagrams. We will be interested in correlators at order 1/N , so we can assume
n
Z n = Z and directly work with the disorder averaged partition function and correlators.
5
For the specific 2-2 interaction we choose in the chain model, the two-point functions connecting different
sites vanish even without averaging over disorder, because of a Z2 fermion parity conservation on each
site. However, even in the more general models that we will discuss in section 6, the cross-site two-point
functions still vanish after averaging over disorder, as a consequence of local SO(N ) symmetry on each site
after random average.

–7–
In other words, quenched equals annealed at the order we work, so we can study annealed
quantities.
To derive the effective action, we begin by integrating the partition function Z[{Jjklm , Jjklm 0 }]
0
over disorder {Jjklm } and {Jjklm } with Gaussian distributions. This will produce eight-
fermion interaction terms that are non-local in time:
Z
Z = Dχ exp(−S[χ])
   4
M Z β N Z β N
X X 1 1 X
S[χ] = dτ χj,x (τ )∂τ χj,x (τ ) − dτ1 dτ2 J02  χj,x (τ1 )χj,x (τ2 )
 0 2 8N 3 0
x=1 j=1 j=1
 2  2 
N
X X N 
+J12  χj,x (τ1 )χj,x (τ2 )  χj,x+1 (τ1 )χj,x+1 (τ2 )  . (3.3)

j=1 j=1

The expression of the averaged partition function can be further simplified by introducing
bilocal auxiliary fields Gx (τ1 , τ2 ) and Σx (τ1 , τ2 ). Σx (τ1 , τ2 ) are Lagrange multipliers which
impose the constraints
N
1 X
Gx (τ1 , τ2 ) = χj,x (τ1 )χj,x (τ2 ). (3.4)
N
j=1

Integrating out the fermion fields χjx (τ ) one obtains the effective action of the bilocal
fields:
Z
Z = DGDΣ exp(−N Seff [G, Σ])
M  β
J02
Z 
X 1
Seff [G, Σ] = − log Pf (∂τ − Σx ) + dτ1 dτ2 Σx (τ1 , τ2 )Gx (τ1 , τ2 ) − Gx (τ1 , τ2 )4
2 0 4
x=1
J12

2 2
− Gx (τ1 , τ2 ) Gx+1 (τ1 , τ2 ) .
4
(3.5)

All terms except the last term ∝ J12 describe decoupled SYK models at each site, and the
J12 term couples the bilocal fields Gx (τ1 , τ2 ) on neighboring sites. In the large N limit, the
saddle point of this action determines the leading order behavior of the two-point function
Gx (τ1 , τ2 ). The saddle point equation δS δSeff
δG = 0 and δΣ = 0 produce the Schwinger-Dyson
eff

equations (assuming time translation symmetry):


1
Gx (iω) = (3.6)
−iω − Σx (iω)
J12
 
2 2 2 2

Σx (τ ) = Gx (τ ) J0 Gx (τ ) + Gx−1 (τ ) + Gx+1 (τ ) (3.7)
2

–8–
with Gx (τ ) = Gx (τ + τ1 , τ1 ) and G(iω) its fourier transformation, and similarly for Σx (τ )
and Σx (iω). This set of equations can be equivalently derived via Feynman diagrams:

x = x + x x (3.8)
 
x x−1 x+1
x J12 x−1 x+1
x x = J02 +

+ ,

(3.9)
2

x x x

where the thick lines represent dressed Green’s functions and the gray disk represents the
self-energy.
Notice that the chain model under disorder average is translation invariant and also
1
has a translation invariant free propagator Gfree
x (τ ) = 2 sgn(τ ). Therefore we consider
translation invariant solutions of the Schwinger-Dyson equations:

Gx (τ ) = Gs (τ ), Σx (τ ) = Σs (τ ); (3.10)
1
Gs (iω) = , Σs (τ ) = (J02 + J12 )Gs (τ )3 . (3.11)
−iω − Σs (iω)

Comparing to equation (2.7), one sees the Schwinger-Dyson equations reduce to exactly p the
same form as those of a (0 + 1)-d SYK model with the coupling constant J = J0 + J12 . 2

Therefore, we can directly apply the SYK model results in Ref. [4, 11]. In particular, we
immediately know that the solution in the conformal limit N  βJ  1 has the following
familiar form:

πτ12 −2∆
 
s ∆ βJ
G (τ1 , τ2 ) = b sin , 0 6 τ12 < β (3.12)
π β
 
1 1
b= − ∆ tan(π∆), ∆ = 1/4.
π 2
p
Here and below, we will denote the effective coupling J02 + J12 in our model as J.
In summary, the two-point function in our chain model is local in space and power-law
decaying (at low temperature) in time, a behavior known as local quantum criticality[39,
40]. This saddle point solution (and the finite βJ corrections to it) is the starting point
for studying other properties of the model.
For example, inserting the solution into the action, we derive the saddle point approx-
imation to the partition function. This gives the order N term in the free energy:

J2 s
 Z  
F 1 s 1 s s 4
= − log Pf (∂τ − Σ ) + dτ1 dτ2 Σ (τ1 , τ2 )G (τ1 , τ2 ) − G (τ1 , τ2 ) .
NM β 2 4
(3.13)

Here we have used space translation symmetry to simplify the notation. (M is the number
of lattice sites, i.e. the spatial volume of the system.) Using the saddle point solution
one can see that thepfree energy density NFM agrees exactly with that in the SYK model
Eq. (2.8) with J = J02 + J12 . Therefore, in the large N limit our model has exactly the

–9–
same zero temperature entropy per fermion S0 ≈ 0.232, and specific heat cv ≈ 0.396
βJ per
fermion. It should be noted that the zero temperature entropy is extensive.
We should remark here that the exact agreement on thermodynamic properties with
the SYK model only holds at leading order in the 1/N expansion. This is similar and
related to the discussion of the two-point function, where the agreement with the SYK
model only holds at leading order. We will see below that 1/N effects, which are the
quantum fluctuations around the large N saddle point, are different in our chain model
than in the SYK model.

3.3 Fluctuations of the collective fields and four-point functions


The fermion four-point function can be determined from the fluctuations about the saddle
point just discussed. At large N the saddle is sharp, and the connected four point function
is small, of order N1 . This is determined by the Gaussian fluctuations of bilocal fields
Gx (τ1 , τ2 ) and Σx (τ1 , τ2 ) about the saddle.
It is convenient to expand about the saddle using variables gx , σx defined by

Gx (τ1 , τ2 ) = Gs (τ1 , τ2 ) + |Gs (τ1 , τ2 )|−1 gx (τ1 , τ2 ),


Σx (τ1 , τ2 ) = Σs (τ1 , τ2 ) + |Gs (τ1 , τ2 )|σx (τ1 , τ2 ), (3.14)

where we have rescaled the fluctuation fields gx , σx by prefactors |Gs |−1 and |Gs | for con-
venience. 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.
Using Eq. (3.4), the connected, averaged four-point function of the fermions (defined anal-
ogously to Eq. (2.3)) can be written as the connected two-point function of gx (τ1 , τ2 ):
1
Fxy (τ1 , τ2 ; τ3 , τ4 ) = hGx (τ1 , τ2 )Gy (τ3 , τ4 )i − hGx (τ1 , τ2 )ihGy (τ3 , τ4 )i
N
= |Gs (τ12 )|−1 |Gs (τ34 )|−1 hgx (τ1 , τ2 )gy (τ3 , τ4 )i. (3.15)
F
More precisely, N is the connected part of the fermion four point function. We will see
that F is of order one at large N , so the connected correlator is of order N1 .
To compute the hgx gy i correlator, we expand the effective action to second order in
the fluctuation fields g, σ, which leads to
Z
1 X
δSeff [g, σ] = − d4 τ σx (τ1 , τ2 )Gs (τ13 ) · |Gs (τ34 )| · Gs (τ42 ) · |Gs (τ21 )|σx (τ3 , τ4 )
4 x
!
Z X1 3J 2 X
+ d2 τ σx (τ1 , τ2 )gx (τ1 , τ2 ) − gx (τ1 , τ2 )Sxy gy (τ1 , τ2 ) . (3.16)
x
2 4 x,y

The spatial kernel Sxy is a tight-binding hopping matrix

J12
(δx,y±1 − 2δx,y ) ,
Sxy = δx,y + (3.17)
3J 2
p
where we continue to use the notation that J = J02 + J12 .

– 10 –
It is straightforward to integrate out σx and obtain a quadratic action for gx alone.
We define Ke as the (symmetrized) four-point function kernel of the SYK model[4, 11]:

e (τ1 , τ2 ; τ3 , τ4 ) = 3J 2 Gs (τ13 ) · |Gs (τ34 )| · Gs (τ42 ) · |Gs (τ21 )|.


K (3.18)

The effective action of gx is


3J 2
Z X h i
δSeff [g] = d4 τ gx (τ1 , τ2 ) Ke −1 (τ1 , τ2 ; τ3 , τ4 ) δxy − Sxy δ(τ13 )δ(τ24 ) gy (τ3 , τ4 ),
4 x,y
(3.19)

which determines the fermion four-point function:


1 1
Fxy (τ1 , τ2 ; τ3 , τ4 ) = hgx (τ1 , τ2 )gy (τ3 , τ4 )i
N |G (τ12 )Gs (τ34 )|
s

1 1 2  e −1 −1
= K − S . (3.20)
N |Gs (τ12 )Gs (τ34 )| 3J 2
Here we have written the correlator in a compact matrix form.
Comparing Eq. (3.20) with the four-point function of the SYK model, the only differ-
ence is with the spatial kernel S. Replacing Sxy by δxy (i.e. taking J1 = 0) reproduces the
SYK result, as expected.
The behavior of the four-point function (3.20) can be analyzed by diagonalizing the
kernel (Ke −1 − S). First of all, due to translation symmetry it is straightforward to do a
spatial Fourier transformation to (x − y) and define the Fourier component
1 1 1 2 h e −1 i−1
Fp (τ1 , τ2 ; τ3 , τ4 ) = K − s(p)δ (τ13 ) δ (τ24 )
N N |Gs (τ12 )Gs (τ34 )| 3J 2
2J 2
with s(p) = 1 + 12 (cos p − 1) . (3.21)
3J
Then one can diagonalize the temporal kernel K e in the same way as in the SYK model. We
write the antisymmetric eigenfunctions Ψh,n (τ1 , τ2 ) where n labels the fourier mode for the
sum of the two times, and h specifies the dependence on the difference of the times. Writing
the corresponding eigenvalues of Ke as k(h, n), the four-point function can be expressed as

1 1 1 2 X k(h, n)
Fp (τ1 , τ2 ; τ3 , τ4 ) = Ψh,n (τ1 , τ2 ) Ψ∗ (τ3 , τ4 )
N N |G (τ12 )G (τ34 )| 3J 2
s s 1 − s(p)k(h, n) h,n
h,n
(3.22)

(where we have used the fact that the symmetrized kernel is Hermitian). The only difference
from the original SYK model is the factor of s(p) 6 1 in the denominator.6 The details of
the eigenvectors and eigenvalues of the temporal kernel, k(h, n) and Ψh,n have been worked
out in Ref. [11], and we will use them in the following sections. It should be noted that
Eq. (3.22) returns to the SYK result at zero momentum p = 0, as one expects from the
form of the effective action.
6
It is interesting to note that the same type of modification was found in a different generalization of
the SYK model in appendix G of [11].

– 11 –
τ1 τ2 τ
Gx (τ1 , τ2 )
0
β

x x
(a) Green’s function (b) Space-time picture for chain model

Bz

(c) Analogy to ferromagnetic spin chain with pinning field

Figure 4: (a) Green’s function Gx (τ1 , τ2 ) is a function of two imaginary time variables, each
defined on the imaginary time circle. It transforms covariantly under the reparametrization
field fx ∈ Diff(S 1 ). (b) The space-time picture for the chain model. The Schwinger-Dyson
equation at conformal limit has global reparametrization symmetry, but the conformal
solution spontaneously breaks Diff(S 1 ) to PSL2 (R). Moreover, the UV term −iω in (3.11)
breaks the emergent reparametrization and lifts the Goldstone modes to pseudo-Goldstone
modes. (c) The situation in the SYK chain model is similar to a ferromagnetic spin chain
with a small pinning field, where the SU(2) symmetry is “almost spontaneously” broken
to U(1), leading to a pseudo-Goldstone mode.

3.4 Symmetry breaking and pseudo-Goldstone mode


Although the general discussion above is sufficient for calculating four-point functions,
it is helpful to gain more physical understanding by analyzing symmetry properties of
the saddle point solution. As in the SYK model, our effective action (3.5) admits an
approximate reparametrization symmetry of time in the IR limit ω → 0, where one can
ignore the −iω in the first term. The reparameterization transformation is defined by
∆
f ∈ Diff(S 1 ) : Gx (τ1 , τ2 ) → f 0 (τ1 )f 0 (τ2 ) Gx (f (τ1 ), f (τ2 )). (3.23)
The symmetry Diff(S 1 ) is broken spontaneously to PSL2 (R) by the solution in equa-
tion (3.12). If this symmetry were exact, the spontaneous symmetry breaking would have
produced infinite number of Goldstone modes corresponding to the spatially dependent
reparametrizations fx ∈ Diff(S 1 )/ PSL2 (R). The −iω term in the effective action (3.11)
plays the role of a small explicit symmetry breaking field, which selects the solution in Eq.
3.12 and turns the Goldstone bosons to pseudo-Goldstone bosons. This is similar to the
situation in a ferromagnetic spin chain with small external magnetic field (see figure 4).
The effect of this symmetry breaking term is small at large βJ, which means that the
pseudo-Goldstone bosons have large fluctuations and make the most important contribu-
tions to the long-wavelength dynamics. In particular, they are responsible for the diffusion
of energy that we will analyze in section 4 and the chaos characteristics we will study in
section 5.

– 12 –
From the view point of the effective action, the pseudo-Goldstone modes are those
fluctuations along the “nearly-flat” direction around the saddle point in the potential of
Seff [G]. As we have discussed, these fluctuations correspond to the spontaneously and ex-
plicitly broken reparametrization symmetry, and therefore can be parametrized by residue
target space Diff(S 1 )/ PSL2 (R) at each point in space. More explicitly, fx acts on the
saddle point in the following way
∆
fx (τ ) ∈ Diff(S 1 ), Gs (τ1 , τ2 ) → Gfx (τ1 , τ2 ) := fx0 (τ1 )fx0 (τ2 ) Gs (fx (τ1 ), fx (τ2 )). (3.24)

For small deformations of time fx (τ ) = τ + x (τ ), the quadratic action for x (τ ) can be


determined by diagonalizing the kernel Ke and using (3.19). Building on the results of [11],
we will find below that to quadratic order in the spatial momentum p, this leads to the
action
√ !
1 X 2αK 2 2 J12 2
S= n,p n (n − 1) + 2 p |n|(n2 − 1) −n,−p , (3.25)
256π n,p βJ 3J
M Z
1 X β i( 2πn τ −xp)
fx (τ ) = τ + x (τ ), n,p = √ dτ e β x (τ ). (3.26)
M x=1 0

The first term is familiar from the SYK model. It is local in time, and can be interpreted as
a quadratic approximation to the Schwarzian derivative action at each site. The coefficient
1
βJ tells us that at large βJ limit the reparameterization fields are soft modes, due to the
approximate reparameterization symmetry. The second term describes a simple coupling
of the reparameterization modes at different sites, but with a nonlocal form as a function
of time.7 We will see that together, these two terms determine both the energy diffusion
dynamics and chaos behavior.8
In the following two sections, we will analyze two different regions of the four-point
function with different ordering of the time variables τ1 , τ2 , τ3 , τ4 as shown in figure (5).
The four-point functions with ordering τ1 > τ2 > τ3 > τ4 , which we will sometimes refer
to as j-j-k-k order, determines the dynamics of collective fields in our model, while the
order τ1 > τ3 > τ2 > τ4 which we will refer to as j-k-j-k determines OTOC after analytic
continuation, and characterizes chaos in our system.

4 The OPE region and the Energy Transport

As we have seen from the two-point functions and four-point functions, single fermion fields
χjx (τ ) do not propagate spatially in our model. The only fields that have nontrivial dy-
namics are the collective bilocal fields gx (τ1 , τ2 ), which are singlet in the flavor SO(N ) sym-
metry. In the four-point function Fxy (τ1 , τ2 ; τ3 , τ4 ), if we denote τ = 21 (τ1 + τ2 − τ3 − τ4 ),
and take the limit |τ2 − τ1 |, |τ4 − τ3 |  |τ |, we are effectively taking an operator product
7
An interesting question is what is the full non-linear form of the effective action which generalizes the
Schwarzian action in (0 + 1)-d case, which we will leave for future work.
8
Notice that this action has three zero modes, n = −1, 0, 1 for each point in space. These are set to zero
by the quotient in Diff(S 1 )/ PSL2 (R).

– 13 –
imaginary time imaginary time

τ1 τ1
τ2 τ3

τ3 τ2
τ4 τ4
space space
x x

(a) j-j-k-k order (b) j-k-j-k order

Figure 5: Two regions of the four-point function are illustrated. Factorization in the
configuration at left gives the propagating bilocal fields. The configuration at right can be
continued to the OTOC which diagnoses chaos.

expansion (OPE) of the fermion fields χj,x (τ1 )χj,x (τ2 ) and similarly for χk,y (τ3 )χk,y (τ4 ).
The four-point function becomes a sum over the propagators of different collective fields
that appear in the OPE. Roughly speaking, different collective fields correspond to different
families of eigenvalues in Eq. (3.22).
The eigenvalues and eigenfunctions of the temporal kernel have been studied in Ref. [4,
11]. In the strong coupling limit βJ  1, the two-point functions take the conformal form
(3.12). This correlation function is covariant under Möbius transformations (PSL2 (R)) of
the time coordinate, which is used in diagonalizing the temporal kernel. It turns out[4, 11,
26] that the eigenvalues of the kernel are given by the simple formula
3 π 1 1
k(h, n) → k(h) = − tan (h − ), h ∈ {2, 4, 6, . . . , + iR+ } = I. (4.1)
2h − 1 2 2 2
The label set I consists of a discrete series and a continuum[11, 26]. For small momentum
p2 . βJJ 2 , we will see that βJ
1
corrections to this formula are important. First, however,
1
we study the case of finite p where the conformal large βJ limit is straightforward.

4.1 The OPE at finite p


The fact that the eigenvalues (4.1) don’t depend on n is a consequence of the conformal
symmetry of the eigenvalue problem at large βJ. It makes it possible to do the sum over
n in Eq. (3.22), which gives particular hypergeometric functions of the cross ratio

sin πτβ12 sin πτβ34


η= . (4.2)
sin πτβ13 sin πτβ24

The sum and integral over the remaining eigenvector index h ⊂ I can then be done following
Ref. [11] by deforming the contour of integration over the continuum portion of I. As the
contour is pushed to infinity, one encounters two sets of poles: one set conveniently cancels
the contribution from the discrete h ⊂ I and the other gives the answer (in the region

– 14 –
τ1 >τ2 >τ3 >τ4 where η < 1)

Γ(h)2 h
 
4π X h − 1/2 k(h)
Fp = −Gs (τ12 )Gs (τ34 ) Res η 2 F1 (h, h, 2h, η)
3 π tan πh/2 1 − s(p)k(h) Γ(2h) h=hm (p)
m=0
(4.3)
The only difference in (4.3) from the regular SYK case [11] is that we have the factor of
s(p). Here 2 F1 (h, h, 2h, η) is the hypergeometric function, which is regular and approaches
constant 1 at small η. The values {hm (p)} that we sum over are the infinite set of positive
solutions of the pole equation 1 − k(h)s(p) = 0. They are approximately double-integer
spaced, approaching hm (p) → 32 + 2m for large m.
Taking the residues, we have
X
Fp = −Gs (τ12 )Gs (τ34 ) c2m (p)η h 2 F1 (h, h, 2h, η) (4.4)

h=hm (p)
m
4π h − 1/2 Γ(h)2 1
c2m (p) = . (4.5)
3k 0 (h) tan πh/2 Γ(2h) s(p)2 h=hm (p)

τ12 τ34
If one considers the OPE limit J −1  τ12 , τ34  τ = τ13 ' τ24 , we have η ' β ,
(π sin πτ
β
)2
and
 −2hm (p)
s s
X X h (p)−2∆ hm (p)−2∆ 2 β πτ
Fp ∝ G (τ12 )G (τ34 ) c2m (p)η hm (p) ' τ12m τ34 cm (p) sin ,
m m
π β
(4.6)

with ∆ = 14 the dimension of the fermion field. For fixed p, the dimension hm (p) labels the
time PSL2 (R) representation of the operator φm appearing in the OPE.
It is also useful to consider the fourier transform of Fp back to position space.9 One can
show that c2m (p) is analytic in a strip surrounding the real p axis, so the fourier transform
will decay exponentially in x. More precisely, c2m (p) has both a pole and a branch cut
for complex p. The pole is at the location where hm (p) becomes an even integer, and the
tan(πhm /2) in the denominator diverges. The branch cut is due to the multivaluedness of
the solution to k(hm ) = s(p)−1 and occurs at the location where k 0 (hm (p)) = 0. For large
m, one can show that both singularities are at a distance of order log m from the real p
axis, so the contribution of the m mode will decay in space as e− log(m)|x| for large x and
large m. Here x is measured in lattice units, so the modes only propagate a few lattice
sites. However, note that the decay factor in the exponential, log(m)|x| ∼ log(h)|x|, grows
less rapidly with h than in a conformally invariant (1 + 1) dimensional theory, where we
would have h|x|.

4.2 The h = 2 contribution for small p


For small momentum p, the first solution to k(hm )s(p) = 1 approaches h0 (0) = 2. Because
of the factor of tan πh2m in the denominator of Eq. (4.5), this leads to a large OPE coefficient
9
Here we consider only the modes m > 1 that have smooth behavior for p → 0. We will consider the
m = 0 (i.e. h ≈ 2) contribution in detail below.

– 15 –
c20 (p) that diverges like p−2 as p → 0. It will be important to understand how this is cut
off at finite βJ.
The would-be divergence can be traced to the h = 2 ⊂ I eigenfunctions of K. e In
this section we study the h = 2 piece in more detail, which requires taking into account
the non-conformal corrections to k(h, n). The correction can be derived by including the
1/βJ correction to the conformal form of the two-point function, as has been discussed in
Ref. [11]. In the following we will apply the SYK results to our generalized model.
To start, we note that due to time translation symmetry, the eigenfunctions Ψ2 (τ1 , τ2 )
τ1 +τ2
−i 2π n
for h = 2 have the form Ψ2 (τ1 , τ2 ) = Ψ2,n (τ12 )e β 2 . In the conformal limit βJ  1,
Ψ2,n (τ12 ) has the following explicit form:

γn sin n πτβ12 πτ12


Ψ2,n (τ12 ) = fn (τ12 ), fn (τ12 ) := − n cos n . (4.7)
2 sin πτβ12 tan πτβ12 β

3
In the above expression, γn2 = π2 |n|(n 2 −1) is the prefactor to normalize Ψ2,n . In the first order

perturbation theory, one can use these zeroth order (conformally covariant) eigenfunctions
to obtain the first order corrections to the eigenvalues, which are given in Ref. [11] by

2 αK
k(h = 2, n) ' 1 − |n| + ... (4.8)
βJ
where αK ≈ 2.852. Summing over the eigenfunctions with this improved formula for the
eigenvalue, we find the h = 2 contribution

Fp,h=2 (τ, τ12 , τ34 ) 8X e−iωn τ 1


s s
' √ fn (τ 12 )fn (τ 34 )
G (τ12 )G (τ34 ) π n 1 − s(p)(1 − α 2 |n| |n|(n2 − 1)
K βJ )

16J X e−iωn τ 1
'√ 2 2π|n|
fn (τ 12 )fn (τ 34 ). (4.9)
2αK n |n|(n − 1) β + Dp2


Here ωn = β n is the Bosonic Matsubara frequencies. In the second step we have expanded
J2
s(p) in the long wavelength limit s(p) ' 1 − 3J12 p2 , and defined a constant D which, as will
be shown below, is the energy diffusion constant of the system:

2πJ 2
D= √ 1 . (4.10)
3 2JαK
Notice that Eq. (4.9) has a smooth p → 0 limit, which reduces to the corresponding formula
of the SYK model in Ref. [11].

4.3 Energy transport


In subsections 4.1 and 4.2 we discussed contributions to the fermion four-point functions
that are naturally organized in terms of the OPE. The operator product of two fermion
fields χj,x (τ + τ 212 ), χj,x (τ − τ 212 ) generates an infinite family of collective fields, including
the reparameterization field (h = 2 contribution) and the local conformal fields φm . In the
strong coupling limit βJ  1, the dominant contribution for small p is the h = 2 piece,

– 16 –
χj χk χj χk
1
∼ c2m (p) ∼ −iω+Dp2

χj χk χj χk
(a) OPE from the conformal part (b) OPE from the correction of h =
2 parts

Figure 6: (a) The conformal limit of the four-point function corresponds to collective fields
that are locally critical and short-range correlated in space (see equation (4.6)). (b) The
non-conformal corrections of the four-point function corresponds to the reparameterization
field, which has a diffusive dynamics in space-time.

which has diffusive behavior as can be seen from Eq. (4.9) and (3.25), also see figure 6 for
an illustration. It is natural to relate the h = 2 contribution to energy transport properties,
since the time reparameterization field is related to diffeomorphisms.
To understand the relation to energy transport more explicitly, we start from p = 0,
where the h = 2 four-point function (4.9) reduces to the SYK result. As shown in Ref.
[11], Fp=0,h=2 (τ, τ 12 , τ 34 ) has a factorized form in the j-j-k-k order region (τ1 >τ2 >τ3 >τ4 ),
which satisfies the following equality[9]:

β2
 s  s 
1 ∂G (τ12 ) ∂G (τ34 )
Fp=0,h=2 (τ, τ12 , τ34 ) = . (4.11)
N N cv ∂β ∂β
Here cv is the specific heat per fermion, such that
N M cv ∂hHi
2
=− = hH 2 i − hHi2 ≡ hH 2 iconn. (4.12)
β ∂β
is the thermal fluctuation of the total energy. Here H is the Hamiltonian of the chain
model. Thus Eq. (4.11) can be rewritten as
1 (∂Gs (τ12 )/∂β) (∂Gs (τ34 )/∂β) 2
Fp=0,h=2 (τ, τ12 , τ34 ) = M hH iconn.
N (∂hHi/∂β)2
∂Gs (τ12 ) ∂Gs (τ34 ) 2
=M hH iconn. (4.13)
∂hHi ∂hHi
Eq. (4.13) implies that the h = 2 piece of the fermion four-point function at p = 0 de-
scribes the thermal fluctuation of the two-point function Gs (τ12 ) induced by the fluctuation
of total energy. If we define the energy density of the chain model as T 00 (x), and define

its Fourier component as T 00 (p) = M −1/2 x T 00 (x)e−ipx , we have T 00 (p = 0) = H/ M .
P

It is natural to generalize Eq. (4.13) to finite (small) momentum and express the energy
density correlation function as
Fp,h=2 (τ, τ12 , τ34 ) N c2v Fp,h=2 (τ, τ12 , τ34 )
hTτ T 00 (−p, τ )T 00 (p, 0)iconn. ' ∂Gs (τ ∂Gs (τ
= . (4.14)
NM2 12 ) 34 ) β 4 ∂Gs (τ12 ) ∂Gs (τ34 )
∂hHi ∂hHi ∂β ∂β

– 17 –
Eq. (4.14) is expected to hold in the limit τ12 , τ34 → 0 (which guarantees the j-j-k-k order
for any finite τ ). In this limit, by a Taylor expansion of Fp,h=2 (τ, τ12 , τ34 ) in τ12 and τ34
one obtains

Fp,h=2 (τ, τ12 , τ34 ) J X |n|(n2 − 1)e−iωn τ  2πτ12 2  2πτ34 2


' √ .
Gs (τ12 )Gs (τ34 ) 9 2αK n 2π|n|
+ Dp 2 β β
β

πα
Evaluating ∂β G(τ ) for small τ and using cv = √K
16 2βJ
and Eq. (4.14), we find
 2 2 
β ωn
N cv X |ωn | 4π 2 − 1
hTτ T 00 (−p, τ )T 00 (p, 0)iconn. = 2 2
e−iωn τ (4.15)
β n
|ω n | + Dp

with ωn = 2πn τ
β . This equation directly gives the Matsubara correlator C00 (p, iωn ). By
10 τ
analytically continuing iωn → ω + iδ in C00 (p, iωn ) and subtracting a contact term[41],
one obtains the retarded energy-density correlation function in the (p, ω) space11 :

Dp2
 
R N cv −iω 3 N cv
C00 (p, ω) = − 1 + O(ω ) ' − . (4.16)
β −iω + Dp2 β −iω + Dp2

Here we have taken the small frequency limit and omitted the higher order terms in ω.
Eq. (4.16) tells us that energy density perturbations in our system satisfy a diffusion
equation with diffusion constant D given by Eq. (4.10). From energy density correlations
one can also derive the thermal conductance:
N cv Dω 2
 
iωβ R
κ0 (ω, p) = Re C 00 (ω, p) = . (4.17)
p2 ω 2 + (Dp2 )2

For p = 0 this formula reduces to κ0 = N cv D = Cv D, with Cv = N cv the specific heat per


site.

In summary, in this section we analyzed the fermion four-point function in the limit
τ12 , τ34 → 0, which determines the behavior of SO(N ) singlet collective fields in the chain
model. The only collective field with nontrivial spatial dynamics is the time reparameteri-
zation field, the dynamics of which describes diffusion of energy in this disordered system,
with a temperature independent diffusion constant. Although there are infinite number of
other collective fields in this system, all of them are locally critical and decay exponentially
in space with order 1 correlation length.
10
One may concern that the fourier coefficient in Eq. 4.15 looks divergent at large |ω| which makes the
analytic continuation less rigorous. However, the expression in the denominator only keeps the leading terms
in (|ωn |/J) expansion. In general, it includes all higher power terms, i.e. Dp2 + |ωn | + i>2 ci (|ωn |/J)−hi
P

where ci are coefficients that can be determined by numerical calculations, and hi are the negative roots of
eigenvalue equation k(h) = 1 (see Eq. 4.1 and Ref. [4, 11]). We start from i = 2 in the summation since the
first negative root h1 = −1 corresponds to the |ωn | term. For i > 2, the roots are generally non-integers
with spacing roughly 2: h = −2.77354, −4.67946, −6.63197, . . .. These higher power terms regulate the
large |ω| behavior for the coefficients, but won’t contribute the low energy physics discussed later.
11
We would like to thank Subir Sachdev for pointing out an error in our earlier version.

– 18 –
5 Chaos and the butterfly velocity

Another interesting aspect of the model is the OTOC, which have been studied as a diag-
nostic of chaos. In addition to the exponent λL that determines the exponential growth of
(anti)commutators in the large N dynamics within a single site, the spatial locality in our
generalized model allows us to study the dynamics of chaos in space. An important param-
eter is the speed vB [16, 42] defined as the speed of growth of the “filled light cone” that
marks the region where operators have large (anti)commutators with an initial operator.
It has the interpretation of the speed at which the butterfly effect spreads in space, and
has also been related to the Lieb-Robinson velocity [43] recently[44, 45]. In this section we
will evaluate vB .

t1 = t + i 3β
4

t3 = i β2

t2 = t + i β4

t4 = 0

Figure 7: Double Keldysh-Schwinger contour with operators equally spaced in imaginary


time. (Note that here we use the convention that the real/imaginary part of t is the
real/imaginary time, which is different from our convention in previous section.)

The OTOC can be computed by analytic continuation of the imaginary time ordered
correlation function. A convenient special case to study is the “regularized OTOC,” where
the operators are equally spaced in imaginary time as shown in figure 7:
N  
1 X 3β β β
F (x, t) = 2 χj,x (t + i )χk,0 (i )χj,x (t + i )χk,0 (0)
N 4 2 4 β
j,k=1
N
1 X
= 2 Tr (rχj,x (t)rχk,0 (0)rχj,x (t)rχk,0 (0)) , (5.1)
N
j,k=1

β
with r = ρ(β)1/4 = e− 4 H /Z 1/4 . The regularization does not change the qualitative be-
havior of the OTOC. For fixed x we expect the following behavior for F (x, t): at early
times it should be approximately equal to the disconnected product −G( β2 )2 (the minus
sign is because we have fermions), and at late times it should be small, indicating a large
anticommutator between χj,x (t) and χk,0 (0). The butterfly velocity vB is defined as the
rate at which the region where F is small expands outwards as we increase t. Of course, to
actually see that F becomes small, we would have to sum all orders in the N1 expansion.
We can only compute the first N1 term, but we assume that the exact F becomes small at
around the time where this correction becomes comparable to the order one disconnected
contribution −G( β2 )2 .

– 19 –
In the imaginary time configuration, the N1 part of the four point function is given by
F
−N , where F is the function studied in the previous section. More precisely, we studied
the spatial fourier transform Fp . To compute the N1 term in F (x, t), we will continue Fp
to the configuration in Eq. (5.1), and then finally fourier transform back to position space.
We will see that there are some subtleties involved in getting an expression for Fp that is
accurate for the small momenta that dominate this fourier transform.
To begin, we will warm up by studying the case where p2  βJJ 2 so that we can use
1
the conformal limit of the kernel. From Eq. (4.2), we have that the cross ratio of the times
in the configuration (5.1) is
2
η= . (5.2)
1 − i sinh 2πt
β
For large t, this is small and purely imaginary, but we begin the continuation from t = 0,
where η is greater than one. The four point function and this continuation are discussed
in detail for the SYK model in Ref. [11]. The only difference in our case is that we have
to insert a factor of s(p). After the contour manipulation and expansion of the relevant
hypergeometric functions [11], one finds that the growing term at small η is

Fp (η) 4 (h−1/2) kR (1 − h) Γ( 12 − h2 )Γ(h − 12 )
1−h
∼− 0 (1 − h) (−iη) (5.3)
G(τ12 )G(τ34 ) 3 tan(πh/2) s(p)kR 21−h Γ( h2 )


h=h∗ (p)
cos π( 14 − h
2) 3 1 + 3s(p)
kR (1 − h) = k0 (h) = , h∗ (p) = . (5.4)
cos π( 14 + h
2)
2h − 1 2

(h∗ (p)−1)t
Using Eq. (5.2), we see that this implies an exponential growth Fp ∼ e β , so we
have a momentum-dependent chaos exponent λL (p) = 2π β [h∗ (p) − 1].
J2
This exponent is largest at small p, where we have h∗ (p) ≈ 2 − 2J12 p2 . This leads to
a p−2 divergence in the prefactor coming from the tan(πh/2) in the denominator of (5.3).
As before, this divergence can be traced to the h = 2 eigenfunctions, which we can treat
more accurately by directly continuing Eq. (4.9) to the configuration (5.1). In appendix
B, we show that the growing term after this continuation is

t √
Fp,h=2 (t) J eβ 2αK
β 2
∼ − · 2π 2
, α≡ . (5.5)
G( ) α β + Dp 4π
2

This expression has a smooth p → 0 limit. However, notice that the exponential growth

t
is e β independent of p, whereas we have argued that at finite p the exact answer should

(h (p)−1)t
be modified to e β ∗ . Also, at finite βJ, we expect a modification of the growth
1
exponent proportional to βJ .
Both of these modifications must come from the sum over h 6= 2 eigenfunctions. We
can perturbatively compute both of them as follows. We have two small parameters, p2
1
and βJ , and we consider both to be first order quantities, of order . If we expand in , we
expect to find
Fp (t) 1 λL ()t 2π
∼− e , B = b1  + b2 2 + ..., λL = (1 − λ1  + ...). (5.6)
G( β2 )2 B() β

– 20 –

t
At order −1 , we have 1
b1 e
β . Comparison with (5.5) then determines b1 . At order 0 we

t
expect a term 2πλβb1 te
1
coming from the small shift in λL . Since this is at order 0 , we can
β

compute it in the theory with p = 0 and βJ = ∞. Indeed, Ref. [11] did find such a term

t
in this limit, which for our case evaluates to 6π β te
β . Matching to the expected term, we

find λ1 = 3b1 , giving the formula


 
Fp (t) 1 2π [1−3b(p)]t α 2π 2
∼− eβ , b(p) = + Dp . (5.7)
G( β2 )2 b(p) J β

1
The exponent in Eq. (5.7) is correct to order p2 and order βJ and the inverse of the prefactor
is correct to the same accuracy. One can check that when βJ = ∞ it agrees with the exact
result Eq. (5.3) to the claimed order in p2 .
We now have an expression that is sufficiently accurate at small p, so we can do the
final step and transform Eq. (5.7) to position space to compute F (x, t). Approximating
the discrete fourier transform as an integral, we have

1 ∞ dp eipx 2π
Z
F (x, t) [1−3b(p)]t
β 2
=1− eβ + ... (5.8)
−G( 2 ) N −∞ 2π b(p)
q

The integrand has a pole at momentum p = i βD , which dominates the behavior for large
x, leading to12
F (x, t) a1 2π (t−x/vB ) 2 2πD
'1− eβ , vB = . (5.9)
−G( β2 )2 N β

This is the main result of this section, giving the butterfly velocity vB . The order N1 term
competes with the order one term, indicating that the anticommutator has become large,
β
when t = 2π log aN1 + vxB . The formula for vB given in Eq. (5.9) agrees with the relation
identified in holographic theories, see Ref. [37].
It is remarkable to find a simple relation between the diffusion constant and butterfly
velocity of this type. In this model, it is a consequence of the fact that the same repa-
rameterization degrees of freedom are responsible both for energy diffusion dynamics and
the OTOC chaos behavior. This is a property that the model shares with conventional
holographic theories, where the gravitational field in the bulk determines both of these
phenomena. It would be interesting to work out whether vB 2 = 2πD continues to hold at
β
1
higher orders in βJ . At least naively, because other modes besides reparameterizations will
become important, one would not expect this equality to persist beyond infinite βJ.
There is an interesting subtlety in the fourier transform that we glossed over above.
In fact, the pole dominates only for x  vβJ Bt
. This means that for x . vJB log N , the pole
analysis will not be correct at the time when the anticommutator becomes large. For such
12
Some constants that appear in the following equations are
1 3
β 2 J2 βJ 2

3βJ 2 J 3αK
a1 = √ , a2 = √ , a3 = √ , a4 = .
2αK J1 2παK J1 2βJ 4πJ12

– 21 –
Figure 8: We sketch two level sets of the function F (x, t). The red/dotted line corresponds
to F = 1 − N1 , and the blue/dashed line corresponds to F = 12 . The blue/dashed curve
is the butterfly cone: operators above this location have large anticommutators with the
operator at the origin. The gray solid line marks the transition between the behavior (5.9)
to the right, and (5.10) to the left.

x, we can approximate (5.8) another way by replacing b(p) in the denominator by b(0) and
doing the Gaussian integral. This leads to

F (x, t) a2 2π 2
'1− √ e β (1−a3 )t−a4 x /t . (5.10)
−G( β2 )2 N t

which is accurate for x  vβJBt


. In this region, we find that the “butterfly cone” is rounded
out, see figure 8. It is rather striking that these two different regions of behavior charac-
terized by (5.10) and (5.9) were also identified in the analysis of stringy corrections to the
holographic computation of F (x, t) in [46].
We will close this section with one further comment. A surprising feature of the large

t
x behavior (5.9) is that the growth as a function of time is given by e β , despite the finite
1
p and βJ corrections to λL in Eq. (5.7). These corrections, which are both negative for
real momenta, cancel against each other when we evaluate at the pole at imaginary p. It
1
would be interesting to know if this persists at higher orders in βJ . Note that in the small
1
x region (5.10), the growth as a function of time is decreased by a βJ correction, with the
same coefficient as in the original SYK model.

6 General construction and higher dimensions

In previous sections, we have focused on a (1 + 1)-dimensional chain example of the gener-


alized SYK model, but our construction actually applies to generic dimensions and more
general forms of interactions. In this section, we will discuss the general form of our model.
We will start from the simple case of higher dimensional regular lattices and then discuss
the even more general cases beyond that.

– 22 –
u
Juuwz

w
Jxxyy

z x y

Figure 9: An example of two-dimensional square lattice model. The Hamiltonian could


contain all kinds of random four-fermion terms, for example, Jjklm,uuwz χj,u χk,u χl,w χm,z
and Jjklm,xxyy χj,x χk,x χl,y χm,y as shown in the figure.

6.1 Generalized SYK model on higher dimensional lattices


In the chain model example, the different sites are only coupled by a 2-2 coupling with two
fermions from each site. This restriction is chosen for simplicity, which is not necessary.
In general, our model can be defined on arbitrary graphs, including higher dimensional
translation invariant lattices and non-translation invariant graphs (see figure (9) for an
illustration on square lattice). We denote the set of sites in the graph as Γ, and label the
sites by x, y, z.... On each site there are N Majorana fermions {χj,x }j=1,...,N . The fermions
are coupled via random four-fermion interactions
X X
H= Jjklm,xyzw χj,x χk,y χl,z χm,x , (6.1)
x,y,z,w∈Γ j,l,k,m

here in the sum over x, y, z, w, one can always define an order of the sites, and avoid
duplication due to different order of x, y, z, w. For fixed sites x, y, z, w, we can further
restrict the range of indices j, k, l, m in the sum to avoid duplication due to the different
order of j, k, l, m.13 This definition makes sure that each 4-fermion term only appears once
in the Hamiltonian. {Jjklm,xyzw } are random couplings which are completely independent,
with mean and variance
2 1 2
Jjklm,xyzw = 0, Jjklm,xyzw = J , (6.2)
N 3 xyzw
where Jxyzw is fixed in the large N limit. This model includes SYK model and the chain
model discussed above as two special cases. The SYK model corresponds to taking Jxyzw =
J uniform for all x, y, z, w, which gives a completely non-local Hamiltonian (or equivalently
it can be treated as the case Γ only contains a single site). The chain model is obtained
by setting Jxyzw to zero except {Jxxxx }x∈Γ and {Jxxyy }x∈Γ for y = x+1.
Remarkably, this general model can still be solved in the same way as the chain model.
From the Feynman diagrams shown in figure 10 one can directly see that the disorder
averaged two-point function is still diagonal between different sites, even if there is no
13
For example: if x = y = z = w, we restrict 1 6 j < k < l < m 6 N ; and if x = y 6= z = w, we restrict
1 6 j < k 6 N and 1 6 l < m 6 N , as we did in the chain model.

– 23 –
Jxyzw x y
y
z w Jxyzw
x z x
x y
w
(b) Ladder in generalized
(a) A watermelon diagram
model

Figure 10: (a) A typical “watermelon” diagram for the general model. Just like SYK
model and the (1 + 1) − d chain model, in the general model the coupling Jxyzw is also
diagonal, such that only fermions with the same flavor and spatial coordinate are con-
nected under random average of disorder fields (dashed line). (b) A typical ladder in the
generalized model. One needs to sum over all possible z and w in the middle to get the
ladder kernel Kxy .

symmetry reason for it to vanish. In path integral approach, we can write down the
effective action of the same collective fields Gx (τ1 , τ2 ) and Σx (τ1 , τ2 ):
Z Y
Z= DGx DΣx exp(−N Seff [G, Σ])
x∈Γ
X 1
Z 
Seff [G, Σ] = − log Pf (∂τ − Σx ) + dτ1 dτ2 Σx (τ1 , τ2 )Gx (τ1 , τ2 ) ,
2
x∈Γ
X  Z 
2
− Cxyzw Jxyzw dτ1 dτ2 Gx (τ1 , τ2 )Gy (τ1 , τ2 )Gz (τ1 , τ2 )Gw (τ1 , τ2 ) ,
x,y,z,w∈Γ
(6.3)

where Cxyzw are combinatorial factors which depend on how many sites in xyzw coin-
1 1
cide. For example Cxxxx = 2·4! , Cxxyy = 2·2!·2! for x 6= y.14 In the large N limit,
the corresponding saddle point equation always admits a translation invariant solution
Gx (τ1 , τ2 ) = Gs (τ1 , τ2 ), which is identical to the two-point function of SYK model with
effective coupling
8 X
J2 = 2
Cxyzw Jxyzw . (6.4)
|Γ|
xyzw∈Γ

Here we denote |Γ| as the total number of sites.


Similar to the chain model, we can expand Gx (τ1 , τ2 ) around the saddle point by
defining Gx (τ1 , τ2 ) = Gs (τ12 ) + |Gs (τ12 )|−1 gx (τ1 , τ2 ). Expanding the effective action to
14 Pk
In general, a group of sites xyzw defines a partition (n1 , n2 , ..., nk ) of 4: i=1 ni = 4 with ni the
multiplicity of site i. For example, xxxx corresponds to partition (4) and xxyy corresponds to the partition
(2, 2). The combinatorial factor is given by Cxyzw = 21 Q 1ni !
i

– 24 –
quadratic order of gx (τ1 , τ2 ) leads to

3J 2 X
Z  
δSeff [g] = d4 τ gx (τ1 , τ2 ) Ke −1 (τ1 , τ2 ; τ3 , τ4 )δxy − Sxy δ(τ13 )δ(τ24 ) gy (τ3 , τ4 ).
4 x,y
(6.5)

The effective action has the same form as equation 3.19 except for a different spatial kernel
Sxy :

4 X 2
Sxy = (Cxyzw +Cxzyw + Cxzwy + Czxyw + Czxwy + Czwxy )Jxyzw . (6.6)
3J 2
z,w∈Γ

The discussion so far does not rely on translation symmetry, and applies to general
graphs. If the graph is a d-dimensional lattice and the coupling has translation sym-
metry, the spatial kernel will also have the same symmetry, so that Sxy = S(x − y)
where the labels x, y should be considered as d-dimensional vectors. In this case the
spatial kernel can be diagonalized by Fourier transformation. If Sxy is short-ranged, the
Fourier transformation s(~ p) is a smooth function of p~, which can be expanded at small p
p) ' 1 − i=1,2,...,d ai p2i . In the same way as in the chain model, we obtain an energy
P
as s(~
diffusion, and ai determines the diffusion constant of energy along the i-th direction.
Following the same approach as in the chain model case, we can also study the OTOC
in the general model. Similar to the (1 + 1)-d case, one finds a Lyapunov exponent 2π β
saturating the chaos bound, and a butterfly velocity vB (for translation invariant systems).
v2
When the spatial kernel Sxy is short-ranged, the universal relation D = 2πTB
still holds.
More details of the higher-dimensional calculation is given in Appendix C.

6.2 Models with global symmetry


As we learned from Ref. [12, 27, 29], a complex fermion version of the SYK model can be
defined, which has similar properties such as local critical two-point functions. The Hamil-
tonian of the complex fermion at chemical potential µ = 0 is H = j,k,l,m Jjklm c†j c†k cl cm
P

with 1 6 j < k 6 N , 1 6 l < m 6 N . Jjklm are also independent random variables. Since
one can always write complex fermion operators in Majorana operators, the complex SYK
model can be viewed as a Majorana model with 2N Majorana fermions cj = 12 (χj1 + iχj2 ).
The coefficients in this Majorana model are different from that of an SYK model because
of the U(1) symmetry requirement.
It is natural to generalize the single-site complex fermion SYK model to models de-
fined on a general graph with a generic global symmetry. The most general form of the
Hamiltonian is given by

X N
X X
P
H= Jjklm,xyzw ηPabcd χj,a,x χk,b,y χl,c,z χm,d,w (6.7)
x,y,z,w∈Γ j,l,k,m=1 a,b,c,d;P

with j = 1, 2, ..., N , a = 1, 2, ..., L labels N L Majorana fermions {χj,a,x } at each site


x ∈ Γ. The index a carries a representation of a global symmetry group. (For Majorana

– 25 –
operators the representation has to be real. In other words, χi,a,x is transformed under a
subgroup of SO(L).) ηPabcd for each P is an invariant rank-4 tensor in the symmetry group.
For example, for symmetry group U(1) ' SO(2), there are two possible invariant tensors
η0abcd = δ ab δ cd , η1abcd = ab cd . 15 We assume the couplings Jijkl
P are independent variables

with

P

P
2 1 P
2
Jijkl,xyzw = 0, Jijkl,xyzw = 3 Jxyzw . (6.8)
N
Once the Hamiltonian is defined one can try to study the effective action in the same
way as before. Now we have to introduce the matrix bilocal field
1 X
Gab
x (τ1 , τ2 ) = hχj,a,x (τ1 )χj,b,x (τ2 )i (6.9)
N
j

and the corresponding Lagrange multiplier Σab x (τ1 , τ2 ). The effective action after integrating
out fermions is
X   1Z 
ab ab ab
Seff [G, Σ] = − log Pf ∂τ δab − Σx + dτ1 dτ2 Σx (τ1 , τ2 )Gx (τ1 , τ2 )
2
x∈Γ
X Xh 2
P P
− Cxyzw Jxyzw
x,y,z,w∈Γ P
Z 
P 0 bb0 cc0 dd0
· dτ1 dτ2 ηabcd ηaP0 b0 c0 d0 Gaa
x (τ ,
1 2τ )G y (τ ,
1 2τ )G z (τ ,
1 2τ )G w (τ1 2 (6.10)
, τ )

Although the effective action is more complicated, the large N saddle point approxi-
mation still apply, at least if we keep L finite in the large N limit. 16 The saddle point
condition of this effective action gives the Schwinger-Dyson equations of G and Σ. In
general there may be saddle points where Gab and Σab have off-diagonal matrix elements,
which we don’t know how to solve analytically. However, there always exists a diagonal
saddle point solution Gab ab s
x (τ1 , τ2 ) = δ G (τ12 ), for which the effective action reduces to the
SYK model with a coupling
8 X X P 2 P
J2 = P
Cxyzw Jxyzw P
ηabcd ηabcd . (6.11)
|Γ|L x,y,z,w
P

Therefore one can always take an expansion around this saddle point and study its stability
using our knowledge about SYK model solution. If this saddle point is stable, these models
have the same locally critical two-point functions as the SYK model, but have different
1 ab
N fluctuations since the fluctuation gx (τ1 , τ2 ) is a matrix field. This will lead to different
four-point functions. We will leave more systematic investigation of these models to future
works.
15
A natural way to define ηPabcd is to write ηPabcd = fab
Pm Pm Pm
fcd , with fab the Clebsch-Gordan coefficients
that maps the representation carried by ab to the representation P . However, it should be noted that the
choice may be redundant. In other words, due to the anti-commutation of Majorana fermion operators,
one may not need all representations P to expand an invariant tensor.
16
The suppression of inter-replica coupling by large N still applies to in this limit, so that we expect the
disorder averaged effective action with a single copy is still meaningful.

– 26 –
6.3 General q
Another type of generalization one can consider is to include interactions that involve q
fermions at a time, rather than just four. For q > 4, the detailed analysis of the chain
model will be very similar to q = 4, with a few minor changes in numerical coefficients.
The relation vB2 = 2πD/β will remain correct at large coupling. One can also consider

including terms with different values of q in the same model. In general, the terms with
higher values of q are more RG-irrelevant, so the terms with the lowest values of q will
dominate in the infrared. Nevertheless, the subleading terms can have interesting effects.
Suppose we have terms with q1 , q2 , with q1 < q2 , and we use the same dimensionful coupling
J for both interactions. Then a perturbative analysis of the Schwinger Dyson equations
indicates that for τ  β, τ J  1 we will have (omitting coefficients)
 
sgn(τ ) sgn(τ ) sgn(τ ) 2 q2
G(τ ) ∝ 2 + 2 + + ... p= +2 −1 . (6.12)
|Jτ | q1 |Jτ | q1
+1 |Jτ |p q1 q1

The first term is the naive conformal limit in the pure q1 theory. The second term is
the correction to the conformal limit, again in the pure q1 theory; this term leads to the
αK correction to the kernel that eventually gives the action for the reparameterization
modes. The third term is the new feature of the theory with both q1 and q2 . It reflects the
contribution of the irrelevant q2 -fermion operator, which affects the correlator at quadratic
order. The dimension of this operator near the IR is ∆ = qq21 , and when ∆ < 32 we have
p < q21 + 1, so this third term will dominate over the second in (6.12). In such a case, the
analysis of the reparameterization action would have to be redone, following e.g. appendix
D of [9]. It would be interesting to compute the energy and chaos dynamics in the resulting
model.

7 Conclusion and discussion

In this paper, we generalize the (0 + 1)-dimensional SYK model of Majorana fermions to


higher spatial dimensions. The generalized model retains many interesting properties of
the SYK model, such as local criticality, extensive zero temperature entropy and maximal
chaos. On top of that, the spatial locality of our generalized model leads to many new
physical properties. We find that single Majorana fermions in our model do not propagate
between different sites, but collective modes made by pairs of fermions have nontrivial
spatial dynamics. In particular, the most important collective mode in the low-energy-
long-wavelength limit is the time reparameterization field, the dynamics of which describes
the diffusion of energy in this system, with a temperature independent diffusion constant
D. This result tells us that our model describes a strongly correlated diffusive metal. The
dynamics of the same reparameterization field also determines OTOC of fermion operators,
from which we can also study the butterfly effect in this model. Our result shows that
chaos spreads in space with a “butterfly velocity” vB . Remarkably, at strong coupling, the
diffusion constant and the butterfly velocity satisfies a simple relation D = vB 2 /2πT , in

consistency with the proposal in the literature about incoherent metals[35–37].

– 27 –
Our model pointed out a new class of solvable interacting lattice models in condensed
matter physics. Usually, solvable models are mapped to weakly interacting theories such
as mean field theories, so that they are not “chaotic”, while the interesting phenomena in
chaotic systems cannot be studied in solvable models. The generalizations of SYK model
is a rare example of solvable but still chaotic systems. Therefore this model provides an
interesting platform for studying various properties of strongly correlated systems, such as
thermalization, entanglement propagation, dissipative transport, etc. It is also natural to
ask whether further generalizations of these models allow us to investigate the possibility
of many-body localization and phase transition between localized and delocalized phases.
From the perspective of holographic duality, the generalized SYK models might be
considered as models that are dual to some kind of incoherent black hole (see[37] and
references therein), but the details of this duality require further work. At strong coupling,
the models do share a key property with conventional holographic systems, which is that a
single set of degrees of freedom dominate and describe both energy diffusion and the chaos
behavior. In a holographic theory, the relevant degrees of freedom are the bulk gravitational
field. Here, it is a reparameterization of time that can vary from place to place. It would
be interesting to derive the action (3.25) for these degrees of freedom from a subset of the
metric degrees of freedom on some black hole background, in a similar way as the derivation
of (0+1)-dimensional Schwarzian action from the Einstein-dilaton theory in approximately
AdS2 background [7, 9]. It would also be interesting to understand the relationship of
this action to recent work on hydrodynamic actions [47, 48]. Another natural question is
whether there are higher-dimensional translation invariant generalizations of SYK models
which are dual to weakly coupled gravity theories in the bulk.

Acknowledgment

We would like to thank Mike Blake, Richard Davison, Luca V. Delacrétaz, Wenbo Fu,
Tarun Grover, Sean Hartnoll, Alexei Kitaev, Juan Maldacena and Steve Shenker for helpful
discussions. We especially acknowledge Subir Sachdev for helpful discussions and comments
on the draft. This work is supported by the National Science Foundation through the grant
No. DMR-1151786 (YG and XLQ). D.S. is supported by the Simons Foundation grant
385600.

A Diagrammatic derivation for four-point functions

In this section, we present a diagrammatic derivation of the fermion four-point functions.


The four-point functions are the leading correlation functions which couple different sites.
In the language of collective field Gx (τ1 , τ2 ), the four-point function comes from its quantum
fluctuations. The four-point functions are essential for understanding transport properties
and also the OTOC measure of chaos, as discussed in section 4 and 5.

– 28 –
x, τ1 y, τ3

...

x, τ2 y, τ4

Figure 11: Ladder diagrams connecting two sites far apart. To connect two sites with
distance n = |x − y|, one need at least n ladders. And the four-point function includes all
possible such diagrams and the partner terms with (τ3 ↔ τ4 ).

The connected four-point function is defined as


N
Fxy (τ1 , τ2 ; τ3 , τ4 ) 1 X
:= 2 hTτ χj,x (τ1 )χj,x (τ2 )χk,y (τ3 )χk,y (τ4 )i − Gsx (τ1 , τ2 )Gsy (τ3 , τ4 )
N N
j,k=1
(A.1)

We follow the main text to use Gs to denote the saddle point of bilocal field, which is the
Green’s function here. Note that the four-point function is non-vanishing after disorder
average only if the spatial and flavor indices appear in pairs, i.e. SO(N ) singlet.
Similar to the SYK model, in the large N limit, the leading contributions to the
connected piece of four-point functions are of order 1/N and consist of the ladder diagrams:

F∼ + + ... − − − ... (A.2)

where the thick lines are dressed Green’s functions solved in the previous sections. The
interaction vertices are paired by random disorder fields Jx and/or Jx0 . Comparing to those
in the SYK model[4, 11, 26], the ladder diagrams in the chain model have extra labels
indicating the spatial coordinates (figure 11). Each ladder can either couple fermions on
the same site, or bring two fermions at site x to the neighboring site x ± 1, as is shown in
figure 12. The on-site terms are contributed by both Jx and Jx0 terms in the Hamiltonian,
while the nearest neighbor terms are only from the Jx0 term.
In general, we can complete the summation of ladder diagrams (equation A.2) using
the Schwinger-Dyson equation:

= − + (A.3)

where the gray box represents the dressed interaction vertex. In terms of algebraic formula,
this equation is written as

F = F0 + KF =⇒ F = (1 − K)−1 F0 (A.4)

with F = , K= , F0 = −

– 29 –
x x x x x y

x x J y y J0 x y J0

x x x x x y
0 0
(a) J-J contraction (b) J -J contraction, (c) J 0 -J 0 contraction,
y =x±1 y =x±1

Figure 12: Three types of the “rungs”: type (a) is the same as the one appears in
SYK, induced by the interactions between fermions at same sites; type (b) comes from the
interactions between fermions at site x and nearest neighbor y = x ± 1, but the “rails”
carry the same site indices, therefore, the effect of interaction doesn’t propagate to next
site; type (c) also comes from the interactions between nearest neighbor sites, and the rails
get shifted by ±1.

In this equation, the ladder kernel K = Kxy (τ1 , τ2 ; τ3 , τ4 ) is treated as an operator acting
on functions of one spatial coordinate x and two time variables. There are two kinds of
nonzero matrix elements of K, given by the diagrams in figure 12 (a) (b) and those in (c)
respectively:

1, x 3, x 1, x 3, x
Kxx (τ1 , τ2 ; τ3 , τ4 ) = J + J0
2, x 4, x 2, x 4, x
= (3J02 + J12 )Gsx (τ13 )Gsx (τ34 )2 Gsx (τ42 ) =: K1 (A.5)

1, x 3, y = x ± 1
Kx,x+1 (τ1 , τ2 ; τ3 , τ4 ) = J0
2, x 4, y = x ± 1
= J12 Gsx (τ13 )Gsx (τ34 )Gsy (τ34 )Gsx (τ42 ) =: K2 (A.6)

For the translation invariant saddle point solution, Gsx (τ ) is independent from x, so that
K1 and K2 are only different by the coefficient in front. Therefore Kxy (τ1 , τ2 ; τ3 , τ4 ) has
the separable form

Kxy (τ1 , τ2 ; τ3 , τ4 ) = Sxy K(τ1 , τ2 ; τ3 , τ4 ) (A.7)


J2
with Sxy = δx,y + 12 (δx,y±1 − 2δxy ) (A.8)
3J
K(τ1 , τ2 ; τ3 , τ4 ) = 3J G (τ13 )Gs (τ34 )2 Gs (τ42 )
2 s
(A.9)

The spatial kernel Sxy is a simple tight-binding hopping matrix ( i.e. an identity matrix
plus a lattice Laplacian), and the temperal kernel K(τ1 , τ2 ; τ3 , τ4 ) is identical to that of the
(0 + 1)-d
p SYK model with coupling constant J. (One should be reminded that we denote
J = J02 + J12 .) The separable form of the kernel Kxy (τ1 , τ2 ; τ3 , τ4 ) allows us to directly

– 30 –
· ... · = ...

Figure 13: A trick[11] to simplify the expressions by using the symmetrized kernel.

apply the results in the SYK model[4, 11] to diagonalize the kernel and solve the four-point
functions in our model.
P
Assuming a formal diagonalization Kxy = h,n,p k(h, n, p)|Ψh,n,p ihΨh,n,p |, we can ex-
P 1
press the four-point function by the inner product: F = h,n,p |Ψh,n,p i 1−k(h,n,p) hΨh,n,p |F0 i
where |Ψh,n,p i = Ψh,n,p (τ1 , τ2 , x) is some antisymmetric eigenfunctions in time which n la-
bels the fourier mode for the sum of the two times, and h specifies the dependence on the
difference of the times, and p labels the fourier mode for space. Technically, one can further
simplify the calculation using the symmetrized kernel17 [11] :

1, x 3, y
K
e xy (τ1 , τ2 , τ3 , τ4 ) := = Sxy K(τ
e 1 , τ2 , τ3 , τ4 ) (A.10)
2, x 4, y

where the symmetrized temporal kernel is defined as K(τ e 1 , τ2 , τ3 , τ4 ) = 3J 2 Gs (τ13 )·|Gs (τ34 )|·
Gs (τ42 )·|Gs (τ21 )|. The simplification works by the following steps: (1) add two rungs (with
absolute value for convenience) to the ladders in F (one on the left18 , one on the right),
see figure 13; (2) we get a multiple of “curved boxes” , then express it as a power of the
symmetrized kernel; (3) in the end, sum over all ladders, which is now a geometric series:

∞ 1 3 ∞
...
X X X
2 s s e n )xy
3J Sxx0 |G (τ12 )|Fx0 y (τ1 , τ2 ; τ3 , τ4 )|G (τ34 | = 2 =2 (K
x0 n=1 2 4 n=1
(A.11)

where the factor of 2 comes from the counting for extra term with 3 ↔ 4. This is the
general formal expression for the connected Euclidean four-point function. We can proceed
by diagonalizing the spatial kernel via plane waves eipx , i.e., k(h, n, p) = s(p)k(h, n) and
write the general expression above in momentum space:

2 X k(h, n)
Fp (τ1 , τ2 ; τ3 , τ4 ) = 2 s s
Ψh,n (τ1 , τ2 )Ψ∗h,n (τ3 , τ4 )
3J |G (τ12 )| · |G (τ34 )| 1 − s(p)k(h, n)
h,n
(A.12)

where k(h, n) is the eigenvalue of temporal kernel. This agrees with equation 3.22 derived
from effective action.
17
Roughly speaking, this trick is used to avoid the computation of the inner product hΨ|F0 i
18
the extra factor 3J 2 Sxx0 is needed to construct the kernel

– 31 –
B Summation trick and the prefactor

In this appendix we determine the large t behavior of the contribution to the OTOC from
the h = 2 modes. We choose a special configuration:

τ12 = τ34 = β/2, | Re(τ )| < π/2 (B.1)

where τ = τ1 +τ2 −τ
2
3 −τ4
is the center of mass time separation, and the requirement | Re(τ )| <
β/4 ensures the out-of-time order. The sum over h = 2 modes in equation (4.9) is

Fp,h=2 (τ ) 32J X (−1)n/2 cos(n 2π


β τ) n
=√ 2πn . (B.2)
Gs ( β2 )2 2αK n>2 even
n2 −1 β + Dp2

We would like to continue τ = it and determine the large t behavior of this sum. To do
this we consider the following integral
i∞+0 2π
π cos(ω β τ )
Z
1 ω 1
I= dω · · 2 · 2πω (B.3)
2πi −i∞+0 2 sin(πω/2) ω − 1 β + Dp2

Convergence is guaranteed by | Re(τ )| < β/4. The integrand has poles at ω = 1 and
2, 4, 6. . . . on R+ . When we deform the contour to the right (see figure 14), we find

Im ω

Re ω
1 2 4 6 8

Figure 14: Contour deformation from dashed line to the blue and red circles surrounding
the poles on the positive real axis.

X (−1)n/2 cos(n 2π
β τ) n

π cos( β τ )
I=− 2πn − 2π . (B.4)
n2 − 1 β + Dp2 4 β + Dp2
n>2 even

The key point is that when we continue to large real time τ = it, the integral I remains
convergent and non-growing in time, so we must have

X (−1)n/2 cos(n 2π
β τ) n

π cos( β τ )
2πn = − + (non-growing). (B.5)
n2 − 1 β + Dp2 4 2π
β + Dp
2
n>2 even

This directly gives (5.5).

– 32 –
C Diffusion and the butterfly velocity in general dimensions

In this section, we sketch the computation relevant to the diffusion and the butterfly
velocity in general dimensional models with translation symmetry. We won’t derive the
exact formula for most general case, but will instead present key steps that determines the
diffusion constant and butterfly velocity.
In the mode with transnational symmetry, the four-point function has simple expres-
sion in terms of momentum eigenvalue:
2 X k(h, n)
Fp~( τ1 , τ2 ; τ3 , τ4 ) = Ψh,n (τ1 , τ2 )Ψ∗h,n (τ3 , τ4 ) (C.1)
3J 2 Gs (τ12 )Gs (τ34 ) 1 − k(h, n)s(~
p)
h,n

Notice here the momentum p~ represents a general high dimensional vector. Further restrict
p) ' 1 − j aj p2j .
P
the model to be local, we have a small p expansion for eigen-value s(~
As we have noticed in the SYK and SYK chain model, the contribution relevant to the
energy fluctuation arises from h = 2 modes, which corresponds to the reparametrization
fields (pseudo-Goldstone modes). At strongly coupling √ limit, the eigenvalue for h = 2
2αK
receives a correction depends on n. k(h = 2, n) ' 1 − βJ |n| + . . .. Therefore, the pole
that determines the diffusion constant has the simple form:
X 1
Fp~ ∝ √ + ...
2αK |n|
aj p2j )
P
n 1 − (1 − βJ )(1 − j
X 1
∝ √ + ... (C.2)
2|n|αK
aj p2j
P
n βJ + j

this is the formula in Matsubara frequency ωn = 2π


β n, after rotating to real time, we have
diffusion pole:
1 2παj J
P 2 , Dj = √ (C.3)
−iω + j Dj pj 2αK
Next is to compute the butterfly velocity, which can be extract from the OTO corre-
lation function. Again, we choose the time configuration which placed four time equally
spacing around imaginary time circle, i.e., we are computing:
N  
1 X 3β β β
F (~x, t) = 2 χj,~x (t + i )χk,0 (i )χj,~x (t + i )χk,0 (0)
N 4 2 4 β
j,k=1
N
1 X
Tr rχj,~x (t)rχk,0 (0)rχj,~x (t0 )rχk,0 (0)

= 2 (C.4)
N
j,k=1

Here the spatial coordinate ~x represents a high dimensional vector. Analogous to the
computation in the chain model, we goes to momentum space, and plug in the h = 2
eigen-functions, which is the leading contribution:

X (−1)n/2 cos(n 2π
β τ) n
p, t)h=2 ∝
F (~ 2πn (C.5)
n2 − 1 Dj p2j
P
n>2 even β + j

– 33 –
Using the trick in appendix B and analytic to real time, we have a formula in momentum
space:
1 2π
t
p, t)h=2 ∝ 2π P
F (~ 2
e β (C.6)
β + D p
j j j

Here we get a Lyapunov exponent λL from h = 2 contribution. In general, the exponent


will receive a 1/βJ correction as well as small p2 correction from h 6= 2 part. Following the
argument in section 5, we have
    
1  2π 2π X
2 
 3αK
p, t) ∝ 2π P
F (~ 2
exp 1 − c  + D p
j j t , c= √ (C.7)
β + j Dj pj
β β  2 2πJ
j

In general dimension, we need to compute the following fourier transformation to get the
butterfly velocity in xj direction:
    
dd p eipj xj
Z  2π 2π X
2 

F (xj , t) ∝ exp 1 − c  + D p
j j t (C.8)
(2π)d 2π 2
P
β + j Dj p j
β β 
j

Notice t ∼ β log N is a large parameter. Therefore, we can integral over all direction except
pj by saddle point approximation, where the saddle point is pk = 0 k 6= j:

eipj xj
Z     
−(d−1)/2 dpj 2π 2π 2
F (xj , t) ∝ t 2π exp 1−c + Dj pj t (C.9)
(2π)d β + Dj p j
2 β β

Then this fourier transformation essentially goes back to the 1-dimensional case, for which
we know that at large xj , the integral is dominated by the pole 2π 2
β + Dj pj = 0. The pole
determines the exponential decaying profile for F (xj , t) at real space,
 

F (xj , t) ∝ exp (t − |xj |/vB,j ) (C.10)
β

with butterfly velocity


2
vB,j = 2πT Dj (C.11)

References

[1] Juan Maldacena. The large-n limit of superconformal field theories and supergravity.
International journal of theoretical physics, 38(4):1113–1133, 1999.
[2] Edward Witten. Anti-de sitter space and holography. Advances in Theoretical and
Mathematical Physics, 2:253–291, 1998.
[3] Steven S Gubser, Igor R Klebanov, and Alexander M Polyakov. Gauge theory correlators
from non-critical string theory. Physics Letters B, 428(1):105–114, 1998.
[4] Alexei Kitaev. A simple model of quantum holography. talks given at the KITP Program:
Entanglement in Strongly-Correlated Quantum Matter, April 7 and May 27, 2015.
[5] R. Jackiw. Lower Dimensional Gravity. Nucl. Phys., B252:343–356, 1985.

– 34 –
[6] C. Teitelboim. Gravitation and Hamiltonian Structure in Two Space-Time Dimensions.
Phys. Lett., B126:41–45, 1983.
[7] Ahmed Almheiri and Joseph Polchinski. Models of ads 2 backreaction and holography. arXiv
preprint arXiv:1402.6334, 2014.
[8] Kristan Jensen. Chaos and hydrodynamics near ads 2. arXiv preprint arXiv:1605.06098,
2016.
[9] Juan Maldacena, Douglas Stanford, and Zhenbin Yang. Conformal symmetry and its
breaking in two dimensional nearly anti-de-sitter space. arXiv preprint arXiv:1606.01857,
2016.
[10] Julius Engelsöy, Thomas G Mertens, and Herman Verlinde. An investigation of ads2
backreaction and holography. arXiv preprint arXiv:1606.03438, 2016.
[11] Juan Maldacena and Douglas Stanford. Comments on the sachdev-ye-kitaev model. arXiv
preprint arXiv:1604.07818, 2016.
[12] Subir Sachdev and Jinwu Ye. Gapless spin-fluid ground state in a random quantum
heisenberg magnet. Physical review letters, 70(21):3339, 1993.
[13] Subir Sachdev. Holographic metals and the fractionalized fermi liquid. Physical review
letters, 105(15):151602, 2010.
[14] Olivier Parcollet and Antoine Georges. Non-fermi-liquid regime of a doped mott insulator.
Physical Review B, 59(8):5341, 1999.
[15] AI Larkin and Yu N Ovchinnikov. Quasiclassical method in the theory of superconductivity.
Soviet Journal of Experimental and Theoretical Physics, 28:1200, 1969.
[16] Stephen H Shenker and Douglas Stanford. Black holes and the butterfly effect. arXiv
preprint arXiv:1306.0622, 2013.
[17] Stephen H. Shenker and Douglas Stanford. Multiple Shocks. JHEP, 12:046, 2014.
[18] Alexei Kitaev. Hidden correlations in the hawking radiation and thermal noise. In talk given
at Fundamental Physics Prize Symposium, 2014.
[19] Daniel A Roberts and Douglas Stanford. Diagnosing chaos using four-point functions in
two-dimensional conformal field theory. Physical review letters, 115(13):131603, 2015.
[20] Stephen H Shenker and Douglas Stanford. Stringy effects in scrambling. Journal of High
Energy Physics, 2015(5):1–34, 2015.
[21] Brian Swingle, Gregory Bentsen, Monika Schleier-Smith, and Patrick Hayden. Measuring the
scrambling of quantum information. Physical Review A, 94(4):040302, 2016.
[22] Guanyu Zhu, Mohammad Hafezi, and Tarun Grover. Measurement of many-body chaos
using a quantum clock. Physical Review A, 94(6):062329, 2016.
[23] Norman Y Yao, Fabian Grusdt, Brian Swingle, Mikhail D Lukin, Dan M Stamper-Kurn,
Joel E Moore, and Eugene A Demler. Interferometric approach to probing fast scrambling.
arXiv preprint arXiv:1607.01801, 2016.
[24] Jun Li, Ruihua Fan, Hengyan Wang, Bingtian Ye, Bei Zeng, Hui Zhai, Xinhua Peng, and
Jiangfeng Du. Measuring out-of-time-order correlators on a nuclear magnetic resonance
quantum simulator. arXiv preprint arXiv:1609.01246, 2016.

– 35 –
[25] Juan Maldacena, Stephen H Shenker, and Douglas Stanford. A bound on chaos. arXiv
preprint arXiv:1503.01409, 2015.
[26] Joseph Polchinski and Vladimir Rosenhaus. The spectrum in the sachdev-ye-kitaev model.
arXiv preprint arXiv:1601.06768, 2016.
[27] Subir Sachdev. Bekenstein-hawking entropy and strange metals. Physical Review X,
5(4):041025, 2015.
[28] Yi-Zhuang You, Andreas WW Ludwig, and Cenke Xu. Sachdev-ye-kitaev model and
thermalization on the boundary of many-body localized fermionic symmetry protected
topological states. arXiv preprint arXiv:1602.06964, 2016.
[29] Wenbo Fu and Subir Sachdev. Numerical study of fermion and boson models with
infinite-range random interactions. arXiv preprint arXiv:1603.05246, 2016.
[30] Dionysios Anninos, Tarek Anous, and Frederik Denef. Disordered quivers and cold horizons.
arXiv preprint arXiv:1603.00453, 2016.
[31] Antal Jevicki, Kenta Suzuki, and Junggi Yoon. Bi-local holography in the syk model. arXiv
preprint arXiv:1603.06246, 2016.
[32] Dmitry Bagrets, Alexander Altland, and Alex Kamenev. SachdevYeKitaev model as
Liouville quantum mechanics. Nucl. Phys., B911:191–205, 2016.
[33] Igor L Aleiner, Lara Faoro, and Lev B Ioffe. Microscopic model of quantum butterfly effect:
out-of-time-order correlators and traveling combustion waves. Annals of Physics,
375:378–406, 2016.
[34] Brian Swingle and Debanjan Chowdhury. Slow scrambling in disordered quantum systems.
arXiv preprint arXiv:1608.03280, 2016.
[35] Sean A Hartnoll. Theory of universal incoherent metallic transport. Nature Physics,
11(1):54–61, 2015.
[36] Mike Blake. Universal charge diffusion and the butterfly effect in holographic theories.
Physical Review Letters, 117(9):091601, 2016.
[37] Mike Blake. Universal diffusion in incoherent black holes. arXiv preprint arXiv:1604.01754,
2016.
[38] A Georges, O Parcollet, and S Sachdev. Quantum fluctuations of a nearly critical heisenberg
spin glass. Physical Review B, 63(13):134406, 2001.
[39] Qimiao Si, Silvio Rabello, Kevin Ingersent, and J Lleweilun Smith. Locally critical quantum
phase transitions in strongly correlated metals. Nature, 413(6858):804–808, 2001.
[40] Thomas Faulkner, Hong Liu, John McGreevy, and David Vegh. Emergent quantum
criticality, fermi surfaces, and ads 2. Physical Review D, 83(12):125002, 2011.
[41] Giuseppe Policastro, Dam T Son, and Andrei O Starinets. From ads/cft correspondence to
hydrodynamics. ii. sound waves. Journal of High Energy Physics, 2002(12):054, 2003.
[42] Daniel A Roberts, Douglas Stanford, and Leonard Susskind. Localized shocks. arXiv
preprint arXiv:1409.8180, 2014.
[43] E. H. Lieb and D. W. Robinson. The finite group velocity of quantum spin systems.
Commun. Math. Phys., 28:251–257, 1972.

– 36 –
[44] Daniel A Roberts and Brian Swingle. Lieb-robinson and the butterfly effect. arXiv preprint
arXiv:1603.09298, 2016.
[45] Yichen Huang, Yong-Liang Zhang, and Xie Chen. Out-of-time-ordered correlator in
many-body localized systems. arXiv preprint arXiv:1608.01091, 2016.
[46] Stephen H Shenker and Douglas Stanford. Stringy effects in scrambling. arXiv preprint
arXiv:1412.6087, 2014.
[47] Felix M Haehl, Ramalingam Loganayagam, and Mukund Rangamani. Adiabatic
hydrodynamics: the eightfold way to dissipation. arXiv preprint arXiv:1502.00636, 2015.
[48] Michael Crossley, Paolo Glorioso, Hong Liu, and Yifan Wang. Off-shell hydrodynamics from
holography. arXiv preprint arXiv:1504.07611, 2015.

– 37 –

You might also like