A generalization of Sachdev-Ye-Kitaev

The SYK model: fermions with a q-body, Gaussian-random, all-to-all interaction, is the first of a fascinating new class of solvable large N models. We generalize SYK to include f flavors of fermions, each occupying Na sites and appearing with a qa order in the interaction. Like SYK, this entire class of models generically has an infrared fixed point. We compute the infrared dimensions of the fermions, and the spectrum of singlet bilinear operators. We show that there is always a dimension-two operator in the spectrum, which implies that, like in SYK, there is breaking of conformal invariance and maximal chaos in the infrared four-point function of the generalized model. After a disorder average, the generalized model has a global O(N1) × O(N2) × … × O(Nf) symmetry: a subgroup of the O(N) symmetry of SYK; thereby giving a richer spectrum. We also elucidate aspects of the large q limit and the OPE, and solve q = 2 SYK at finite N.


Introduction
The Sachdev-Ye model [1], as recently revived and simplified by Kitaev [2], possesses, for large N , three remarkable properties: conformal invariance in the infrared, solvability, and maximal chaos. While there are models that contain some of these properties, SYK is the first to have all three, as was recognized by Kitaev in a series of incredibly insightful seminars [2].

JHEP02(2017)093
Broadly speaking, until recently two classes of large N theories have been studied: matrix models and vector models, in which the dynamical variables transform in the adjoint or fundamental representation of a local or global SU(N ) or O(N ) symmetry group, respectively. Matrix models are closely related to string theories [3][4][5][6][7][8], with the most concrete realization being the duality between supersymmetric gauge theories and string theory in Anti-de Sitter space [9]. N = 4 super Yang-Mills is conformally invariant, and at large 't Hooft coupling the bulk gravity has black holes, so it should be maximally chaotic. However, it is not easily solvable. Vector models also have a long history and recently have been shown to be dual to interesting gravity theories. The critical O(N ) vector model is conformally invariant and solvable, and the bulk dual is higher spin Vasiliev theory [10,11]. However, it is integrable for large N , so it is not likely to be chaotic. Roughly speaking, matrix models are too difficult to be explicitly solvable, while vector models are too simple to have the same rich properties. One would like a model that lies in between: one that is sufficiently complicated to be chaotic, while still simple enough to allow for direct analytic calculations for strong coupling. SYK is such a model.
At large N , the dominant Feynman diagrams for matrix models are planar diagrams, whereas the dominant diagrams for vector models are bubble diagrams. The SYK model is dominated by a new class of Feynman diagrams, which have been referred to as sunset, or watermelon, diagrams. The SYK model may be just one example out of a much broader and new class of models. Past studies of large N models have been extremely fruitful for understanding both quantum field theories and string theories. One may hope that the study of SYK-like models will also prove productive.
SYK is a quantum mechanics model, living in 0+1 dimensions. While two dimensional CFTs have been extensively studied and categorized, one dimensional CFTs have not. In fact, it has been argued that 0 + 1 dimensional CFTs with nontrivial dynamics do not actually exist [12]. SYK confirms this: the four-point function breaks the SL(2, R) conformal invariance [2,13,14], consistent with holographic studies of AdS 2 [15]. It appears that in one dimension a theory can at best only be "nearly" conformally invariant. In SYK the breaking of conformal invariance, to leading order in 1/N , is confined to a single dimension-two operator appearing in the OPE, so the power of conformal invariance is still largely applicable.
The SYK model consists of N 1 Majorana fermions χ i , with a q-body Hamiltonian with quenched disorder, H = i 1 ,...,iq J i 1 ,...,iq χ i 1 χ i 2 · · · χ iq . (1.1) The model has qualitatively similar properties for any choice of even q ≥ 4. The couplings J i 1 ,...,iq are independently chosen from a Gaussian, O(N ) invariant, distribution with zero mean and a variance proportional to J 2 N 1−q . When evaluating observables, say correlation functions, a disorder average is performed at the end of the calculation. Of course at large N the model is self-averaging: randomly chosen, but fixed, J i 1 ,...,iq give the same results as disorder averaged J i 1 ,...,iq . One can alternatively think of the J i 1 ,...,iq as nearly static free bosonic fields; at leading order in 1/N , this gives the same connected correlation functions [16], and furthermore, allows one to gauge the O(N ) symmetry [17]. To leading order JHEP02(2017)093 in 1/N the fermions are non-interacting, and the two-point function of the fermions satisfies a simple integral equation which can be explicitly solved near the infrared fixed point. The fermions start with dimension 0 in the UV, and flow to dimension ∆ = 1/q in the IR. After the disorder average, the dynamics is invariant under an O(N ) global symmetry, χ i → O ij χ j , with OO T = 1, much like a vector model. The bilinear, primary, fermion operators, singlets under O(N ), are schematically N i=1 χ i ∂ 2n+1 τ χ i . In the UV, these operators have dimension 2n + 1. In the IR, the dimensions receive an order-one shift for small n, and approach 2∆ + 2n + 1 asymptotically for large n. The standard AdS/CFT dictionary relates the dimensions of CFT single-trace operators for matrix theories, or bilinear singlet operators for vector models, to the masses of particles in the bulk dual. This would imply that the SYK dual has a tower of particles in the bulk, with masses, in units of the AdS radius, roughly spaced by two. This spectrum differs from N = 4/AdS 5 × S 5 duality where for large 't Hooft coupling only a small number of massless modes survive, or vector model/Vasiliev duality, where a tower of massless modes appears in the bulk. In [14] it was argued that the bulk dual of SYK might be a string theory with the string scale comparable to the AdS radius, and thus non-local or stringy. But what the dual of SYK is, and the extent to which it is nonlocal, remains an open problem.
The goal of this paper is to generalize the SYK model. We would like to understand how large the class of such models is, and which features are generic and which are special to SYK. This paper will not add anything new to the bulk interpretation of SYK, but the dual bulk theory, whatever it is, should be able to incorporate this more general class of models.
Two seemingly important ingredients in SYK are: (a) 0 + 1 dimensions, where the fermions are dimensionless, thereby ensuring that any product of fermions is a relevant perturbation, and (b) quenched disorder, which plays an important role in the solvability at large N . The generalization we explore is one in which there are f flavors of fermions, χ a i , where i = 1 . . . N a and a = 1 . . . f , with a Hamiltonian, where I is a collective site index, and the subscript on the fermion is the site while the superscript is the flavor. The number of sites, N a , for each fermion, as well as the order of the interaction, q a , can depend on the flavor a, as long as N a /N remains finite as N = a N a → ∞.
In section 2 we derive the Schwinger-Dyson equation for the two-point functions of the fermions, and find that the model (1.2) generically has an IR fixed point. While the IR dimension for SYK (1.1) was ∆ = 1/q, for (1.2) a set of f transcendental equations determine the dimensions ∆ a . In the limit of large q a , these have simple analytic solutions. Furthermore, for large q a one only needs to sum a particular subset of Feynman diagrams, thus yielding an explicit expression for the two-point function and the spectral function.
In section 3 we study the spectrum of composite operators. After the disorder average, the generalized model (1.2) has an O(N 1 ) × O(N 2 ) × · · · × O(N f ) symmetry. The singlet bilinear operators are Na i=1 χ a i ∂ 1+2n τ χ a i for any a ∈ {1, . . . , f }. So we expect there to be f towers of operators. We derive equations determining the IR dimensions of these operators.

JHEP02(2017)093
We prove that for any choice of parameters: f , N a 's, q a 's, there is always a dimension-two operator in the spectrum. In SYK, the dimension-two operator is responsible for both the breaking of conformal symmetry in the four-point function and for maximal chaos. The same properties hold for the generalized model.
An instructive case to study is the generalized model with all N a equal to N/f and all q a equal to q. It has the symmetry O(N/f ) × · · · × O(N/f ). The spectrum contains a tower identical to that of SYK with a qf body interaction, along with a new tower that appears with a degeneracy of f − 1. Indeed, this model is similar to SYK with N fermions and a qf body interaction, but the full O(N ) symmetry is broken and consequently more singlet operators exist, allowing for a richer model.
In appendix A we consider the path integral for the generalized model. This provides an alternate way of computing the correlation functions, with the saddle point giving the Schwinger-Dyson equations for the two-point functions, and the leading 1/N fluctuations about the saddle giving the four-point function. In appendix B we consider (1.2) with an additional scalar; a special case of this includes supersymmetric SYK [19]. Finally, appendix C solves SYK for q = 2 at finite N . SYK for q = 2 is like N fermions with a random mass matrix; the randomness makes it nontrivial, though it is less interesting than q ≥ 4. This appendix can be read independently of the rest of the paper.
2 Two-point function

SYK
Let us recall the SYK model [2]. 1 It contains N Majorana fermions with the anticommutation relation {χ i , χ j } = δ ij . The action is, where the coupling J i 1 ,...,iq is totally antisymmetric and, for each i 1 , . . . , i q , is chosen from a Gaussian ensemble. The two-point function of the J i 1 ,...,iq is taken to be, At leading order in 1/N , (2.2) is equivalent to the simpler normalization, The particular scaling with N in the choice (2.3) is in order to obtain a nontrivial large N limit, while the other factors are for convenience. One can consider SYK for any even q ≥ 2, JHEP02(2017)093 with q = 4 being the prototypical case [2]. At q = ∞ there are some simplifications [14]. The case q = 2 is simplest, and is equivalent to an O(N ) vector fermion with a random mass matrix, although in many ways it is qualitatively different from the q > 2 models. We solve the q = 2 SYK at finite N in appendix C. At zero coupling, the Euclidean two-point function T χ i (τ )χ j (0) ≡ G(τ )δ ij is given by, where the factor sgn(τ ) (sgn(τ ) = 1 for τ > 0 and sgn(τ ) = −1 for τ < 0) accounts for the fermion anticommutation. To leading order in 1/N , the Schwinger-Dyson equations for the two-point function drastically simplify and are given by, The first of these, (2.5), is the standard equation expressing the two-point function in terms of the one-particle irreducible self-energy Σ(ω). The second equation, which is written in position space, is a special feature of SYK (see figure 1). At leading order in 1/N , the only diagrams being that survive are nested sunset diagrams; all others are suppressed by some power of 1/N . These equations can be combined into a single integral equation; however, an analytic solution to this equation is not known. At strong coupling, |Jτ | 1 (equivalently, the infrared limit), one can drop the iω in (2.5), to get, One can verify that, is a solution to (2.7) provided one takes,

JHEP02(2017)093
The Fourier transform of G(τ ), given in (2.8), is useful in verifying this, where we defined, . (2.11) What is special to SYK is that the IR Schwinger-Dyson equations (2.7) are invariant under reparameterization of time, ). Therefore, although (2.8) is at zero temperature, we can easily construct the finite-temperature two-point function by mapping the real line to a circle [2,[35][36][37].

A generalization of SYK
The model we introduce is a generalization of SYK (2.1). It contains f flavors of fermions, with N a fermions of flavor a, each appearing q a times in the interaction, so that the Hamiltonian couples q = f a=1 q a fermions together. We continue to let the subscript on the fermion χ a i denote the site i ∈ {1, . . . , N a }, while the superscript a will now denote the flavor a ∈ {1, . . . , f }. Explicitly, the action is, where I is a collective index, I = i 1 , . . . , i q 1 , . . . , j 1 , . . . , j q f . The coupling J I is antisymmetric under permutation of indices within any one of the f families, and is drawn from a Gaussian distribution, where the disorder average J I J I is given by 14) It will be convenient to make the following definitions, The class of models (2.12) for large N is characterized by f − 1 independent continuous parameters 0 < κ k < 1, as well as the q k , which can be any positive integers provided that their sum is even. After the disorder average, SYK (2.1) has O(N ) symmetry, while in the generalized model (2.12) the symmetry is broken to the subgroup O( In the SYK model (2.1), one can generalize the action to have multiple interaction terms, with different q, each coming with its own independent disorder J i 1 ...iq . This sum JHEP02(2017)093 of SYK Hamiltonians is just as solvable as SYK, however it is not especially interesting since in the IR the term with smallest q will be dominant. In the model (2.12), one can also consider generalizing the action to include sums of interaction terms. However, now the IR can be more interesting, since there can be multiple terms (with the same total q) that are equally important in the IR.

Two-point function
The free two-point function for each χ a i , is again given by (2.4). Away from the UV it will continue to be the case that the two-point function is diagonal in flavor and site space. Denoting the two-point function for the flavor k fermion by G k (τ ), the self-energy for the flavor k fermion is (see figure 2 Making use of (2.14), (2.15), this simplifies to, An alternative way to obtain (2.17) is by performing the replica trick to do the disorder average, introducing mean fields, integrating out the fermions, and taking the large N saddle point; see appendix A. For one flavor, (2.17) reduces to the SYK expression for the self-energy (2.6). 2 The factors are as follows. The factor ( qa!) 2 comes from the square of the prefactor of the interaction term in (2.12). There is a factor of q k from the number of contractions with the ingoing fermion, and another q k with the outgoing fermion, and a (q k −1)! from the contraction of the remaining flavor k fermions amongst themselves. There is also a factor of qa! from contractions of flavor a fermions, for all the other flavors.

JHEP02(2017)093
We first determine whether there is an IR fixed point and if so, what is the IR dimension ∆ k of the fermions of flavor k? In the IR, the two-point function should take the form, To find the normalization b k and dimension ∆ k , we first insert the above ansatz into (2.17) and take the Fourier transform, 3 Inserting (2.19) and (2.18) into the IR limit of (2.5), Σ k (ω)G k (ω) = −1, gives, . (2.21) The first equation is just the statement that the IR dimension of the coupling J I is zero.
Simplifying the second gives, Equating all the (2.22), for k ranging from 1 to f , gives f − 1 equations. Combined with (2.20), for any given choices of κ k and q k , we have a set of f equations for the f unknown dimensions ∆ k . 4 These equations have simple solutions in the limit of q a 1, as we show in the next section.

Large q k
If the number of fermions of flavor k appearing in the interaction (2.12) is large, q k 1, then from (2.20) we know that ∆ k 1. Let us assume q k 1 for all k. In this limit, (2.22) simplifies to b qa a = 1 2 κ k Q k ∆ k , with the solution, (2.23) Eq. (2.23) shows that a hierarchy in the q k 's for different flavors, or in the κ k 's, will lead to a hierarchy in the ∆ k 's. The smallness of the dimensions ∆ a suggests one should be able to solve for the twopoint function at all energies. This was done for SYK at large q in [14]. Here we perform JHEP02(2017)093 an analogous computation for the generalized model in the large q k limit. The two-point functions are taken to be, Taking the Fourier transform, for which we use the shorthand F , Inverting to get G k (ω) −1 , (2.5) allows us to identify, where in the second equation we have done an inverse Fourier transform of the first. Combining with (2.17) gives, We have such an equation for every k. Thus, we can express g a in terms of g k for any a, k, Summing (2.27) for all k and using (2.28) gives, where the rescaled J is kept finite in the large q k limit. The solution to (2.29) is easily derived for finite temperature, namely with g a (τ ) = g a (τ + β) [14], where v is defined implicitly in terms of J . At zero temperature (2.30) becomes, which combined with (2.28) gives, where ∆ k is given by (2.23). Having the exact solution, we can take the IR limit J|τ | 1 to find the individual normalizations of the two-point function (2.

Graphical solution
The two-point function in the large q limit can alternatively be found by summing an appropriate set of Feynman diagrams. We will show how this works in SYK. Due to large q combinatorics, the diagrams that appear most often are like the ones shown in figure 3(a), rather than those in figure 3(b). The Feynman diagrams that are summed at large q can be characterized as those diagrams that, under a single vertical cut, break up into two tree diagrams. The self-energy can therefore be found recursively, as shown in figure 3(c). The equation corresponding to figure 3(c) is, where G 0 (ω) is the free two-point function (2.4). Rearranging (2.34) gives where the second equation is the inverse Fourier transform of the first. Letting we get, which is (2.27) for one flavor.

Spectral function
The frequency space two-point function follows from (2.24), (2.32), Introducing a Schwinger parameter, and performing the τ integral in eq. (2.38) gives, The spectral function (as defined by eq. (C.4)) for the flavor k fermion is therefore, where ∆ k is given by (2.23). Since ∆ k 1, this is sharply peaked around small λ. If there is only one flavor, then ∆ = 1/q. This spectral function is for q 1. For q = 2, the SYK spectral function is instead a Wigner semicircle (C.5).

Effective action
We have so far discussed the model (2.12) directly in terms of the fermions, finding the two-point function at large N through study of Feynman diagrams. It is useful to also consider the path integral approach. Employing the replica trick, one can carry out the disorder average, and then integrate out the fermions after the introduction of new fields For one flavor, this reduces to the effective action for SYK [2] (see [1] for an analogous expression for the SY model, and [18] for the Dirac fermion version of SYK). The large N saddle point of the action gives the Schwinger-Dyson equations for the two-point function found previously from Feynman diagrams. In particular, varying S eff with respect toG k (τ 1 , τ 2 ), and assuming time-invariance, gives (2.17), while varying with respect tõ Σ k (τ 1 , τ 2 ) yields (2.5) for each flavor. The saddle point solutions are denoted by G k (τ 1 , τ 2 ) and Σ k (τ 1 , τ 2 ). To leading order in 1/N , the free energy is given by the saddle of (2.42),

JHEP02(2017)093
Following [14], one can differentiate with respect to J to get, For large q a , G a (τ ) was found explicitly in section 2.3. Also, since the partition function only depends on βJ, it follows that J∂ J = β∂ β . Thus for large q a , where v is defined in terms of J in (2.30). Up to the choice of normalization of the variance of the disorder, J I J I , this is the same as for SYK with N fermions and a f a=1 q a body interaction. So the entropies are also the same. In order to see a distinction, one must study the 1/N corrections.

Dimensions of composite operators
We begin by reviewing and adding some detail to the computation in [2] for the IR dimensions of the SYK composite operators. The primary O(N ) invariant bilinear operators are, where the coefficients d nk are chosen so that the operators are primary. For instance, The general form of d nk will not be important for us. Note that operators that are schematically of the form i χ i ∂ 2n τ χ i are not primary. We would like to compute the overlap between the state created by the composite operator O n acting at time τ 0 , and two fermions at times τ 1 and τ 2 , respectively. In other JHEP02(2017)093 Figure 4.
(a) The diagrams being summed to compute the three-point function This can be done iteratively, as shown in (b) (see eq. (3.4)), with the kernel shown in (c) adding rungs to the ladder. In the IR we can simplify (b) to get (d).
words, the three-point function, χ i (τ 1 )χ i (τ 2 )O(τ 0 ) , which we will denote by v(τ 0 ; τ 1 , τ 2 ). If the two fermions just propagated without interacting with each other, this would be found by Wick contractions, is the first diagram that appears in figure 4(a). We must also include a sum over all the ladder diagrams in figure 4(a). One can perform the sum by solving the equation (see figure 4 where the kernel is the operator that adds a single rung (see figure 4(c)), where τ ab ≡ τ a −τ b . Letting the composite have dimension h, in the IR the solution to (3.4) will take the form of conformal three-point function, For h > 2∆, the term G 0 χχO is much smaller than (3.6) in the IR, τ 12 1, so we can drop it in (3.4). Thus, (3.4) .7) gives [2], where ψ(∆) was defined in (2.11). A plot of g(h) is given in figure 5. One can see that there is a tower of h's for which g(h) = 1. For large h the solutions to g(h) = 1 are approximately h ≈ 2∆ + 2n + 1.
There are solutions to g(h) = 1 for h < 2∆ as well. These solutions immediately follow from the h > 2∆ solutions due to the symmetry g(h) = g(1 − h). However, they do not cor- Figure 6. The diagonal (a) and off-diagonal (b) components of the kernel (3.13), (3.14) for two flavors with q 1 = q 2 = 3. The coloring scheme is the same as in figure 2.
respond to dimensions of composite operators. Recall that dropping the first term in (3.4) was justified in the IR only for h > 2∆. (For h < 2∆, it is instead justified in the UV). Knowing the dimensions of the "single-trace" operators, one can say something about the bulk dual of SYK. The AdS/CFT dictionary relates the dimensions of single-trace operators to the masses of bulk fields, for AdS 2 . So the dual of SYK has a tower of particles in the bulk, one for each solution to g(h) = 1 for h > 2∆. For large integer n, these have approximate masses m ≈ 2∆ + 2n + 1.

Generalized model
We now generalize the calculation to the model (2.12). The operators (3.1) now have a superscript O a n to account for the different flavors, and the kernel is now a matrix in flavor space, K mn (τ 1 , τ 2 , τ 3 , τ 4 ), where m denotes the flavor of the incoming fermions on the left at times τ 1 , τ 2 and n denotes the flavor of the outgoing fermions on the right at times τ 3 , τ 4 , see figure 6. The off-diagonal component K kl has flavor k propagators along the rails, while the rung consists of q k − 1 flavor k propagators, q l − 1 flavor l propagators, and q a flavor a propagators for all a = k, l, where the combinatorial factor in front is, 5 The diagonal components of the kernel are similar, but with slightly different propagator powers and combinatorial factors. Using (2.14), (2.15) and simplifying we get the diagonal 5 The factor of N qa a , for a = l, k, comes from the site index summation within the rung. For flavors k, l, there are only q k − 1, q l − 1 propagators in the rung, so those give factors of N There is then an additional factor of N l because the Feynman diagrams are built by adding the kernel to the left (see figure 4); so the l index will get summed over.
As in SYK, we must find the eigenvectors and eigenvalues of the kernel. Letting g be an eigenvalue, and v a (τ 12 ) the components of an eigenvector, where the first term on the left is from the diagonal term in the kernel (3.13) and the second term is from the off-diagonal terms (3.14). Inserting the propagator (2.18) and evaluating gives, In order to eliminate the dependance on ω, we must choose, In SYK we know that (3.8) is a special case of (3.6) with 2α = 2∆−h. Similarly here, we let

JHEP02(2017)093
where the diagonal and off diagonal components ofK are, . (3.23) In getting from (3.18) to (3.22) we made use of the product of normalizations of the propagators b qa a given in (2.21). If there is one flavor,K 11 reduces to (3.9). The next step is to find all h for which there is an eigenvalue g ofK (3.22) that equals 1. This is in principle straightforward: for any q a , κ a in (2.12) ones solves (2.20), (2.22) to find the IR dimensions ∆ a of the fermions, then for fixed h one finds the f eigenvalues of K, and then for each of those eigenvalues solves for all h such that the eigenvalue is equal to one. Aside from some special cases, we can not write a general and explicit answer for the h's. However, it is easy to see that there will always be a dimension 2 operator in the spectrum. For h = 2 (3.23) simplifies to, Inserting this into (3.22), one can easily verify that the following vector is an eigenvector ofK with eigenvalue one. Verifying this requires using q a ∆ a = 1, and nothing else. Perhaps surprisingly, it is not even required that the ∆ a are actual dimensions: one does not need to impose (2.22). The dimension-two operator is important: it leads to the breaking of conformal invariance and to maximal chaos; we will comment more on it in the next section. Another universal feature (for any number of flavors greater than one) is the seeming presence of a dimension-one operator. Inserting ρ(h = 1) = −1 into (3.22), one can verify that there are f − 1 eigenvectors ofK that have eigenvalue one. For any k ∈ {2, . . . , f }, such an eigenvector has two nonzero components, In fact, verifying this requires no assumptions on ∆ a . The presence of these dimensionone operators suggests a symmetry. In fact, this symmetry is simple to see from the effective action (2.42). 6 One can rescaleG qaG a (τ 1 , τ 2 ), for any a = 1, while leaving the IR limit of (2.42) invariant.

Equal q a , κ a
A simple and instructive case is when all the q a are equal to some q, and all the κ a are equal, for all flavors a. The dimensions ∆ a are then, by symmetry, all equal to ∆ = 1 f q . The matrixK in (3.22) factorizes, where ρ(h) is given by (3.23) and is independent of the flavor. The eigenvalues ofK are thus, where I = i k = f n k + p − 1, where n k ranges from 1 to N 1 , and gets rid of the flavor index on the fermions. There are now f N 1 sites; however, this is not the same as SYK with a qf body interaction, since the interactions are not all-to-all, being restricted to occur between particular qf sets of sites.

Four-point function
Having found the dimensions of the bilinear singlet operators O a n , the next step is to compute their OPE coefficient. The OPE between two fermions will include all the O a n and their descendants, and will take the form, where c a,b n are the OPE coefficients and C n (τ 12 , ∂ τ 1 ) = 1+. . . is fixed by conformal invariance. The OPE coefficient can be extracted by computing the three-point function between the two fermions and O, which in section 3.1 was labelled as v(τ 0 ; τ 1 , τ 2 ) and satisfied eq. (3.4). Thinking of K as a matrix with indices (τ 1 , τ 2 ), (τ 3 , τ 4 ), the formal solution of (3.4) is, where we have generalized (3.4) to account for multiple flavors. Notice that when we computed the dimensions in section 3.1, we were allowed to drop the G 0 χχO term in (3.4), arguing it was unimportant in the IR. However, for finding the OPE coefficients one is interested in the UV, τ 12 1, so this term is essential. We will also be interested in the four-point function. Defining the bilocal, and proceeding formally, one can perform a double OPE expansion on the four-point function, 7 g a (τ 1 , τ 2 )g b (τ 3 , τ 4 ) = 1 N a N b n,e c a,e n c b,e n C n (τ 12 , ∂ τ 1 )C n (τ 34 , ∂ τ 3 ) 1 |τ 13 | 2hn,e . (3.33) The right-hand side is a sum of conformal blocks, given by hypergeometric functions of the conformally invariant cross ratio, n,e c a,e n c b,e n x hn,e 2 F 1 (h n,e , h n,e , 2h n,e , x) , x = τ 12 τ 34 τ 13 τ 24 .
(3.34) This is similar to two-dimensional CFTs, except here we have one cross-ratio instead of two. At large N , the leading and first subleading in 1/N pieces of the four-point function are, one uses the kernel (3.13), (3.14) to add rungs to the ladder. Summing all the ladder diagrams,

SYK
The technical challenge in evaluating (3.37) explicitly comes from inverting 1 − K. Recall the procedure used in SYK. One first finds a complete basis of eigenvectors of the kernel. This turns out to be given by (3.6) with h ranging over even positive integers h = 2, 4, 6, . . ., as well as h = 1/2 + is where s > 0 [13,14,38]. One then projects (3.37) onto this basis and performs the sum/integral over the discrete and continuous tower of h's to find (3.34) with OPE coefficients c n [14], , where α 0 (q) = 2πq (q−1)(q−2) tan π q , (3.38) where g(h) is given by (3.9) and h n are the solutions of g(h n ) = 1. Eq. (3.38) is for h n > 2. There is an additional complication that occurs for the h = 2 block. One can notice that g(h = 2) = 1, and since h = 2 is part of the basis of eigenvectors used to invert 1 − K, this causes the four-point function to diverge in the conformal limit. The h = 2 block must therefore be treated outside the conformal limit. Moving slightly away from the IR, the eigenvalue g(h = 2) gets slightly shifted away from 1, and so the h = 2 block gives a finite but large, and non-conformal, contribution to the four-point function. Since its prefactor is dominant, its growth controls the behavior of the finite temperature out-of-time-order four-point function used to probe chaos [39,40]. The growth of the h = 2 block occurs with a Lyapunov exponent 2πT that saturates the chaos bound [41]. 8

JHEP02(2017)093
In section 3.1.1 we found that the generalized model always contains a dimension two operator; assuming its OPE coefficient doesn't vanish, this implies that in the IR the generalized model, like SYK, both breaks conformal invariance and is maximally chaotic.
To compute the four-point function (3.37) for the generalized model with generic q a and κ a , one would need to repeat the procedure used for SYK, accounting for the additional complexity of having flavor. However, in the case that all the q a are equal and all the κ a are equal, it is simple to find the four-point function, and this is the case we focus on.

Equal q a , κ a
If all the q a are equal to q, and all the κ a = 1/f , then the kernel matrix (3.13), (3.14) factorizes into a flavor-space matrix and a function of times, we diagonalize (3.37) in flavor space, forming O T F O, to find, (3.36). Both of the terms appearing in F ab are similar to what occurs in SYK, so we can write the answer, where, Here α 0 (q) is given by (3.38) and ρ(h) is given by (3.23) with fermion dimension ∆ = 1/2q. To be clear, the h n appearing in (3.45) and (3.46) are the solutions of (2q − 1)ρ(h n ) = 1 and −ρ(h n ) = 1, respectively. Recall that we found in section 3.1.2 that with two flavors

JHEP02(2017)093
(with q 1 = q 2 = q, κ 1 = κ 2 ), the spectrum of bilinear composite operators contains two towers: a tower that matches the 2q body SYK tower, and a new tower, see figure 7. The OPE coefficients c S n are for the 2q body SYK tower. Note that (3.45) is for h n > 2; as discussed before, the contribution of the h = 2 block diverges in the conformal limit.
The OPE coefficients c A n are for the new tower. Notice that this vanishes for the h = 1 operator. The OPE coefficients c 1 n and c 2 n , in terms of (3.45), (3.46), are given by, where, for simplicity of presentation, rather than writing c a,b n , we have explicitly separated the two towers. 9 A more intuitive way to think about this four-point function is to define the symmetric and antisymmetric combinations of the bilocals (3.32), The symmetric correlator probes only the SYK tower, and matches the SYK 2q body four-point function. The antisymmetric correlator, probes only the new tower. One can also reproduce (3.50), (3.51) from the path integral picture, see appendix A.1.

Discussion
The SY model [1] involves all-to-all interactions between spins in some representation of SU(M ), H = N i,j=1 M µ,ν=1 J ij S µ i ν S ν j µ , with Gaussian-random couplings J ij . Writing the spins as products of two fermions, this becomes a four-fermion interaction. One of the key realizations of [1] was that the model is solvable in the double scaling limit, N → ∞, M → ∞, M/N → 0. It was recognized in [2] that a simpler model is one that avoids spins altogether and goes directly to the fermions, H = J ijkl χ i χ j χ k χ l . There is then a clear generalization to a model with a q-index coupling and a q-body interaction [2].
In this paper, we have made another straightforward generalization, involving f flavors of fermions with N a sites for each flavor and a f a=1 q a body interaction. Perhaps It was recognized in [14] that the SYK model simplifies in the limit q 1. Here we pointed out that in the q 1 limit, only a particular subset of Feynman diagrams need to be summed. For any even q ≥ 4, the SYK model has qualitatively similar properties. In the generalized model introduced in this paper, there are more parameters to vary, and one may wonder if there are corners of parameter space which either lead to simplifications or qualitative differences.
We have only begun exploring the parameter space, focusing on the symmetric case of an equal number of sites for each flavor, as well as interaction orders q a that are independent of the flavor. The main qualitative difference, as compared to SYK with a qf body interaction, is more singlet operators resulting from a symmetry that is a subgroup of the O(N ) symmetry of SYK. One feature we found, that holds for any choice of parameters, is the presence of a dimension-two bilinear singlet operator in the infrared. Another was the presence of a dimension-one operator; however, in the symmetric case considered, its OPE coefficient vanished.
Nontrivial and solvable models are both rare and valuable. It is now clear that the class of SYK-like models is much larger than just the SY model. Just how large this class is, if there are further generalizations, and the precise characterization of the Feynman diagrams, at each order in 1/N , are all still open problems. We may hope that exploring this structure will provide guidance towards understanding the dual string theory, if there is one.

JHEP02(2017)093
where α is the replica index, α ∈ {1, . . . , M }, a is the flavor, a ∈ {1, . . . , f }, i a is the site index, i a ∈ {1, . . . , N a }, and I is a collective site index, I = i 1 , . . . , i q 1 , . . . , j 1 , . . . , j q f , and P [J I ] is the probability distribution for the J I (2.13). Doing the Gaussian integral over the disorder, (A.1) becomes, Having done the disorder average, we see that there is a O(N 1 ) × O(N 2 ) × · · · × O(N f ) symmetry. We thus introduce the collective fields, by inserting delta functions, whereΣ αβ a (τ 1 , τ 2 ) acts as a Lagrange multiplier. We insert into (A.2) such a delta function for each replica index pair α, β and each flavor a. This gives, Integrating out the fermions gives

B Model with a scalar
In this appendix we consider a model with a boson field. It is a slight variant of (2.12) and has the action, . . , i q 2 , . . . , j 1 , . . . , j q f , and q 1 = 1, and q = f a=2 q a . This has a similar interaction as (2.12), but the first flavor is with a boson instead of a fermion. The boson is taken to be auxiliary, having UV dimension [φ] = 1/2, so the coupling J I has dimension 1/2. We restrict to only one boson, q 1 = 1, in order to ensure that the interaction is relevant. We also require q be even. 10 The only technical distinction between finding the IR dimensions for (B.1) compared with (2.12) is that the boson propagator is symmetric in time. The ansatz for the IR boson propagator is, where ψ(∆) is defined in (2.11), while for the IR fermion propagator it is, for k ≥ 2. In the IR, one drops the free propagator appearing in the Schwinger-Dyson equation (2.5), so for both the boson and the fermions one has Σ k (ω)G k (ω) = −1. The self-energy is again given by (2.16), with the disorder average normalization given in (2.14).
Repeating the several steps in section 2.2 gives the following equations, 10 A supersymmetric variant of SYK was introduced in [19] (see also [42]). A more minimal supersymmetric SYK, with only the interaction J ijk φiχjχ k , is being studied in [43,44]. This interaction would be a special case of (B.1) with one fermion flavor with q2 = 2.

C Random mass matrix fermions
The simplest SYK model is for q = 2: fermions with a random mass matrix. For this case, all computations can be performed exactly, without the restriction of being near the fixed points or working at large N . In this appendix, we solve the q = 2 model. For infinite N this is trivial, while for finite N it is slightly more involved but follows from standard matrix model techniques. One should keep in mind that the q = 2 case has multiple features that are not representative of SYK at larger q; in particular, it is not chaotic.

C.1 Infinite N
The Schwinger-Dyson equations for the two-point function (2.5), (2.6) are integral equations for general q, but become a simple quadratic equation for q = 2. The solution is One can also find this directly by summing non-crossing rainbow diagrams (see figure 9), where G 0 (ω) is the bare propagator (2.4), while C n are the Catalan numbers, 11 We thank Yu Nakayama for sharing his results and explaining the supersymmetric model.

JHEP02(2017)093
The Catalan numbers are the number of different ways n + 1 factors can be completely parenthesized; here the parentheses are the rainbows. Summing (C.2) gives (C.1). One can also write (C.1) as, where the spectral function ρ(λ) is the Wigner semi-circle, with support for |λ| < 2J, At finite temperature, the frequencies in (C.1) should be viewed as the Matsubara frequencies, ω n = (2n + 1)π/β. Taking the discrete Fourier transform of (C.4) gives, where 0 < τ < β. In the limit of zero temperature βJ 1, we can evaluate the integral to obtain, where L 1 is the modified Struve function, and I 1 is the modified Bessel function. While both L 1 and I 1 grow exponentially, the two-point function (C.7) decays monotonically with |Jτ |. The combination is sometimes denoted by M 1 ≡ L 1 − I 1 . We can do a strong coupling expansion of (C.7), where we see that the first term matches what was expected from the IR limit of the Schwinger-Dyson equations (2.8), (2.9).

Comments
One comment is that in summing the Feynman diagrams giving (C.2) it is important to work at finite temperature. Each of the diagrams individually has IR divergences: the Fourier transform of any of the individual terms in (C.2) will diverge in the limit of β → ∞. Of course, one could have chosen to regulate the IR divergence in some way other than working at finite temperature. However, finite temperature is natural. The point is just that the dimensionless expansion parameter is βJ. Another comment is that aside from the implicit appearance of β in the Matsubara frequencies, (C.2) has no explicit β dependance. This is a property that is special to q = 2. For SYK with q ≥ 4, one could solve the Schwinger-Dyson equations perturbatively around weak coupling, giving an expansion of the form, (ω n β) 2l , (C.9) with some coefficients g kl . One can derive recursion relations for g kl , but we have not found a way of solving them.

JHEP02(2017)093
We have been considering Majorana fermions. One can instead study q = 2 with Dirac fermions, At leading order in 1/N , this gives the same two-point function (C.1). Working with Dirac fermions gives slightly more flexibility, as one can introduce a chemical potential. With no chemical potential, as in (C.10), one is at half-filling. Explicitly, consider a single free Dirac fermion H = ω 0 c † c. (Adding a chemical potential just corresponds to adding to (C.10) such a term for each fermion, with chemical potential µ = −ω 0 .) The finite-temperature two-point function is trivially, Since we are at finite temperature, fields have the time range 0 < τ < β. The two-point function is a function of the difference between two times, and so naturally has the range −β < τ < β. However, from (C.11) we see that for 0 < τ < β, G(τ − β) = −G(τ ). We can thus restrict to 0 < τ < β. The filling fraction Q is defined as the expectation value of the occupation number, (C.12) We can choose the filling fraction by choosing ω 0 . It is clear that for any finite temperature, if ω 0 = 0, then there is no energy cost to being in the state |1 versus |0 , and so the filling fraction is 1/2. Note that the limits of T → 0 and ω 0 → 0 do not commute. If we set T = 0 at finite ω 0 (including ω 0 = 0), then we get zero filling: from (C.11), G(τ ) = θ(τ ).
Finally, the q = 2 model sums rainbow diagrams. There are many other models that sum rainbow diagrams. For instance, two-dimensional QCD [45] has the same rainbow diagrams, where the fermions are the quarks, and the disorder lines are the gauge field propagators. Also, the recently studied three-dimensional U(N ) k Chern-Simons theory coupled to scalars or fermions also sums rainbow-like diagrams [46,47]. A simple large N quantum mechanics model that sums rainbow diagrams is the Iizuka-Polchinski model [48] (see also, [16,49]). The IP model has a harmonic oscillator in the adjoint representation of U(N ) plus a harmonic oscillator in the fundamental representation of U(N ), coupled through a trilinear interaction. In the limit that the mass of the adjoint goes to zero, this is essentially the same as the model eq. (C.10), at leading order in 1/N . The reason we say essentially the same is because in the IP model the fundamental is effectively at zero filling. In other words, its free two-point function is θ(τ ) as opposed to 1 2 sgn(τ ), and correspondingly, the infinite N two-point function after summing the rainbow diagrams is only the first term, I 1 , in (C.7). In addition, at subleading orders in 1/N , differences will arise between the model eq. (C.10) and the IP model. This is because the adjoint propagator will receive quantum corrections, whereas the J ij J ij "propagator" in eq. (C.10) is always a constant.

Four-point function
It is simplest to write the four-point function in frequency space. This is defined as, Written as a series in 1/N , ijkl + . . . . At leading order in 1/N there is a disconnected piece, At first subleading order in 1/N , the four-point function is a sum of ladder diagrams, like SYK for general q. However for q = 2 there is an extreme simplification, since the rungs only contain the disorder lines. Since there is no momentum exchange, summing the ladders in frequency space simply involves summing a geometric series, which gives, It is only necessary to establish the first term in the s-channel piece. The other term, as well as the t and u channels, follow from antisymmetry. Through a Fourier transform and analytic continuation of (C.14), one finds there is no exponential growth in the out-of-timeorder four-point function [16], and so the q = 2 model is not chaotic.

C.2 Finite N C.2.1 Dirac fermion
We now compute the two-point function at finite N for the random mass matrix fermion, (C.10). For fixed coupling J ij , this is just N free fermions with mass matrix J ij , so the nontrivial part is to perform the disorder average. Specifically, Consider first the trivial case of N = 1. This is just a fermion with a random mass. Then (C.15) reduces to (C.4) with a spectral function,

JHEP02(2017)093
So at N = 1 the spectral function is a Gaussian, while at N = ∞ it is the Wigner semicircle (C.5). The two-point function for N = 1 can also be found by summing Feynman diagrams, (2n − 1)!! (JG 0 (ω)) 2n . (C. 18) Writing the double factorial as a Gaussian integral, and interchanging the integral and the sum, we recover (C.17). Explicitly performing the integral gives the two-point function in terms of the complimentary error function, , ω > 0 .
(C. 19) In the zero temperature limit we also find, We now move on to the case of general N , using the method of orthogonal polynomials to evaluate (C.15). We can write (C. 16 We now take linear combinations of the columns of (C.22), transforming it into a matrix with i, j element, φ j (λ i ), where φ j (λ i ) is a polynomial with lowest element 1 and highest element λ j i . The determinant (C.22) remains invariant under these operations. We can write the determinant as a sum of permutations of the integers from 0 to N − 1, (C.23) We choose the φ n such that, dλ φ n (λ)φ m (λ) e −λ 2 /2J 2 = f n δ nm . The φ n will be proportional to the Hermite polynomials, defined as, Now to evaluate (C.15), note that, and so we find that the spectral function is, Evaluating the sum gives, (C.30) where we have used that f k =J 2k+1 √ 2π k!, which follows from (C.26), (C.27). A plot of (C.30) is shown in figure 10.
An alternative way to write the two-point function is to perform the integral over λ in (C.4) before evaluating the sum over k appearing in the spectral function (C.29). After the introduction of a Schwinger parameter, the integration over λ yields a Laguerre polynomial. Using that the sum of the Laguerre polynomials is an associated Laguerre polynomial N −1 k=0 L k (x) = L 1 N −1 (x), we find The procedure is now similar to the Dirac case. We can write the determinant as a sum of permutations of the even integers from 0 to N − 2, The φ n are the same as in the Dirac case. The partition function now involves just the even normalization constants, For evaluating the two-point function, note that since the eigenvalues come in pairs, The two-point function is thus, This is similar to (C.31), except it involves a sum only over the even Laguerre's. The spectral function is, (C.46)  Figure 11. Plot of the spectral function (C.46) at N = 50 for the q = 2 Majorana SYK. This differs from the random mass matrix fermion spectral function in the region of small λ, see figure 10; the distinction goes away at infinite N .
This is similar to the spectral function for the Dirac fermion (C.30), except for the addition of the last term in (C.46) that is 1/N suppressed relative to the first two (and the trivial distinction that occurs at order 1/N betweenJ andJ), see figure 11.
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.