A new class of SYK-like models with maximal chaos

We investigate a model closely related to both the original Sachdev-Ye-Kitaev (SYK) model and the N\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \mathcal{N} $$\end{document} = 1 supersymmetric SYK model. It consists of N real Majorana fermions and M auxiliary bosons with Yukawa interactions. We consider the large N and M limit and keep the ratio M/N fixed. The model has two branches characterized by the conformal dimensions of fields, which we compute as a function of the ratio M/N. One of the branches contains the supersymmetric saddle for M = N. As we take the limit M/N → ∞ both branches coincide and we obtain the same conformal dimensions as SYK. Furthermore, we determine the Lyapunov exponent of the model and find maximal chaos independent of M/N.

In particular the SYK model is a (nearly) conformal field theory (CFT) in the infrared, and is assumed to have a nearly anti de sitter (AdS) dual in this regime [10][11][12]. For these low energies the model can be described by a Schwarzian [13], which also appears on the bulk side in the AdS 2 dilaton gravity.

JHEP01(2019)166
In this paper we consider a particular model closely related to the N = 1 supersymmetric extension of SYK. Instead of having an equal number, N , of fermions and bosons we consider the case where we have M bosons and N fermions and study its behaviour as a function of the ratio M/N .
In section 2 we will introduce the model and discuss in more detail the relation to the (supersymmetric) SYK model. Afterwards, in section 3, we consider the effective action. We derive the equations of motion and consider the solutions at strong coupling. We find two families of solutions that we label by their conformal dimensions at M = N (rational or irrational). Afterwards we elaborate on the behaviour of the conformal dimensions at different regimes, e.g. as M/N → ∞ we obtain the same behaviour as in the original SYK model. In section 4 we compare the entropy of both solutions and determine that the rational solution is the dominant saddle for all M/N .
In section 5 we compute the Lyapunov exponent and find that is independent of M/N due to a subtle cancellation.

The model
The model consists out of N Majorana fermions obeying {ψ i , ψ j } = δ ij and M (auxiliary) bosons. We will use indices a, b to denote the bosons and i, j, k for the fermions (no ambiguity will arise). The Lagrangian is given as follows: where ψ denote the Majorana fermions and φ the bosons. The coupling C aij is defined to be antisymmetric in the last two indices, which are contracted with the Majorana fermions. In [23] a similar term was studied as a perturbation upon the ordinary SYK model. The fermions are dimensionless, whereas the bosons φ and couplings C have dimension of E 1/2 . Notice that we have two parameters M and N . We are interested in taking the limits of both M and N going to infinity but keeping M/N fixed. In other words we have that M = αN for some fixed α. From now on we will always assume that two a indices are summed up to M whilst the other i, j, k, . . . are summed up to N . Lastly, we let the coupling be disordered averaged by the following distribution: Here J has the dimension of energy and is larger than zero. We can now compute some basic one-loop diagrams for both the fermions and the bosons. We show the one-loop corrections to the two point functions in figure 1, which are proportional to some power of M/N (that can easily be checked). In fact one can check that any boson loop adds a factor of M N and each fermion loop

Relation to SYK
Let us first examine the relation to the original SYK model [1,24] with Hamiltonian To check the similarity we start by plugging in the algebraic equation of motion for φ a back into the Lagrangian. The equation of motion is found to be φ a = i 2 C a ij ψ i ψ j . After plugging it into eq. (2.1) we obtain the Hamiltonian: This is also the presentation that one can see in [23]. We can then use the antisymmetry in the last two indices of C aij and the commutation relations of the Majorana fermions to rewrite this to: where we defined the constant E 0 = − 1 16 C 2 aij (recall that a is summed to M and i, j up to N ). Comparing now to the standard SYK Hamiltonian, eq. (2.4), we find: (2.7) The notation indicates that the asymmetry on the right hand side is only in i, j, k and l, which in turn of course means that J ijkl is completely asymmetric. The above expression for the J coupling shows us that the model is essentially obtained by performing a Hubbard-Stratonovich (HS) transformation on SYK. Of course apart from this HS transformation we have also chosen a different distribution (see eq. (2.2)) compared to SYK. This means that J ijkl are no longer the independent Gaussian variables and this is the cause of the differences between the models. Figure 2. This is the leading diagram that contributes to interactions between replicas. In the figure the top three lines would be associated with a different replica index than the other three below. One can check that this figure is proportional to 1/N .

Relation to supersymmetric SYK
The N = 1 supersymmetric SYK model was introduced in [22], the Lagrangian density is given by: There are two important differences compared to the model described in eq. (2.1). The first important aspect is that there are N bosons, which is the same as the number of fermions (which has to be true for supersymmetry). Secondly the coupling C ijk in the supersymmetric case has to be completely antisymmetric. Note that the equal number of bosons and fermions is also necessary for the antisymmetry in the coupling. In other words, starting from eq. (2.1) we can obtain the supersymmetric model by setting M = N and making the coupling completely antisymmetric. It is precisely when the coupling is completely antisymmetric (and hence M = N ) that the Lagrangian is invariant under supersymmetry transformations.

Effective action and conformal dimensions
To find the effective action we will follow the standard procedure of averaging over the disorder in C aij by using the replica trick (see appendices in [19,25]). As in the usual SYK case we will assume replica diagonal matrices. To justify this we have to compare log Z and log Z since assuming replica diagonal matrices corresponds to evaluating the latter instead of the former. The usual argument (see [25]) is to consider diagrams that are in log Z but not in log Z. The leading diagram belonging to the former but not the latter is shown in figure 2 and as can be verified it is suppressed by 1 N M . Thus in the large N limit these contributions will be subdominant. Using replica symmetry, the result of the disorder average becomes:

JHEP01(2019)166
For the last term we introduced brackets below the sum to indicate that j, k sum up to N whilst a sums up to M . We now introduce bilocal fields for both the fermions and bosons as follows: We insert them into the partition function by Lagrange multipliers. Afterwards we are only left with Gaussian integrals for both the fermions and bosons. Completing these leads to: On the left hand side we divided out a factor of N , but could just as well have taken out M (recall that M/N is fixed).
Let us now vary with respect to G φ and G ψ to obtain the self energies: These equations can also be obtained using the melonic structure of the Feynman diagrams at large N and M , just as in ordinary SYK. The Schwinger-Dyson equations are obtained by varying with respect to the Σ (we assume time translation symmetry and go to Fourier space):

Two saddle points
In order to solve the above equations we have to assume the strong coupling limit βJ 1. This implies that in eq. (3.6) we can ignore the first terms on the right hand side. Hence

JHEP01(2019)166
we can write the equations as follows (we have Fourier transformed back to time): We then use the following (conformal) form for the two point functions: To obtain conditions on the conformal dimensions we plug these into the saddle point equations above, eq. (3.7) and eq. (3.8). Afterwards we Fourier transform using [22,24]: Some other useful relations for Γ functions are . (3.14) After plugging this all in we obtain the following relations: These relations (for M = N ) have also been derived in [22]. By comparing the frequency dependent parts we obtain the first condition on the conformal dimensions: As a side note, under this condition the saddle point equations have the conformal symmetry, very analogous to the original SYK model:

JHEP01(2019)166
where f (τ ) a smooth function (in one dimension Conf(R) ∼ = Diff(R)). To obtain results for finite temperature we use this symmetry with f being the exponential map for example. Coming back to eq. (3.15), we can obtain another constraint by taking the quotient, which yields the (transcendental) equation: This result, for M = N , is also obtained in [22], although it contains some typos. In [23] it is also shown for M = N , albeit in a different form. The second condition, eq. (3.18), can also be recast to an equation for ∆ φ using eq. (3.16): First we solve eq. (3.18) for M being equal to N . This case overlaps with supersymmetric SYK (as commented upon in the introduction) and we find the same solutions as in [22]. The first solution is given by: We label this solution as the "rational" solution. In the supersymmetric model this solution is the one that preserves supersymmetry. In that case, the supersymmetric Ward identity G φ = ∂ τ G ψ , together with eq. (3.9) implies ∆ φ = ∆ ψ + 1 2 [22], obviously obeyed by eq. (3.20).
There is another solution with positive conformal dimensions, it is however irrational:

21)
As one can easily check this does not satisfy ∆ φ = ∆ ψ + 1 2 and hence would break supersymmetry. A similar situation arises in [26] where there are also two solutions, one preserving and one breaking the supersymmetry.

Arbitrary M and N
Let us now vary the ratio M/N and find the conformal dimensions as a function of this ratio. We solve eq.  branches (or families) of solutions, labelled by their behaviour at M = N . The rational solution was also found in [23]. When M/N becomes large the rational and irrational flow to the same point. This can be understood by considering the defining equations eq. (3.18) and eq. (3.19). When one takes the limit of M/N going to infinity there is only one solution left: This is an interesting point, since we find that the conformal dimension of the fermion is exactly that which one finds in the original SYK model [24]. In particular they find that the conformal dimension for an arbitrary number of interactions q SYK is given by 1 q SYK . Let us consider our model with q interactions as described in appendix A. We derive there the relation (q − 1) = 1 2 q SYK and furthermore that the solution for ∆ ψ in the general q case (eq. (A.5)): We find thus in the M/N → ∞ limit that we obtain the same behaviour as SYK, and conclude that the two models have the same infrared fixed point.
The behaviour for small M/N can be understood by taking the small M/N limits in the defining equations. This is equivalent to considering the limit N/M going to infinity in eq. (3.19). Consider the following two limits:

JHEP01(2019)166
Applying these in eq. (3.19) shows us that for small M/N we either need to consider the case where ∆ φ is very small or the case where it goes to one. The latter corresponds to the rational family. For the irrational ∆ φ 1 case we find from eq. (3.19) that for small M/N it behaves as: The corresponding dimensions for the fermion can of course be found by eq. (3.16). Lastly let us investigate the rational ∆ ψ for small M/N . Observing figure 3, we assume that ∆ ψ is small and consider eq. (3.15) and eq. (3.18). We can then solve as follows: The solution corresponding with ∆ ψ going to zero might be understood from figure 1, since in this limit the right diagram corresponding with the self energy corrections to the fermion vanishes. We then expect the fermion to reduce to its free propagator, i.e. with ∆ ψ = 0. Such an argument can however not be made for the case where instead ∆ φ goes to zero.

Dominant saddle
In this section we will determine what is the dominant saddle by comparing the entropies of both solutions for arbitrary M/N . For the case M/N = 1 we can do the computation analytically for the rational case as we will see below. For the computation we will follow [22] and use the model for a q-interaction (meaning a vertex with one boson and q − 1 fermions, with q odd), see appendix A for an overview of the changes. The free energy becomes: Now we derive with respect to q (we continue the values of q to the reals) such that we don't have to evaluate the first terms. We take the fields to be on-shell such that we only need to explicitly take the partial derivative of the last term where the Gs are now the finite temperature versions, obtained by conformal symmetry mentioned before (eq. (3.17)):

JHEP01(2019)166
The integral can then be computed straightforwardly (using the periodicity in the τ variables): where C is a constant independent of β. The constant C is a diverging quantity independent of q contributing to the ground state energy, but will not contribute to the zero temperature entropy, similar to the scenario in [22]. It is important to note that apart from the overall factor, the M/N dependence is also in A q−1 B (eq. (A.3)) and the conformal dimension ∆ ψ (figure 3). The zero temperature entropies S R and S I are labelled by their rational or irrational origin, see eq. (3.20) and eq. (3.21) respectively, and given by: where then i ∈ {I, R}, depending on the branch we consider. We will call the integrands (including the constants) also i(q) and r(q). The expressions for A q−1 B are shown in eq. (A.3) and are also dependent on the conformal dimensions. Since we don't know the analytical expressions for ∆ ψ at arbitrary M/N we will solve this problem numerically and afterwards fit these as a function of q. In particular we proceed as follows: • Fix a value for M/N • Solve ∆ R ψ and ∆ I ψ numerically • Fit these conformal dimensions as a power series in 1 q • Find the irrational integrand i(q) and the rational integrand r(q) After we have the integrands we can compare both the entropies by integration. There are a couple remarks regarding the power series for the ∆ ψ . First of all, we already obtain good fits for the numerical solution when going up to order 1/q 4 . Secondly, for large q the leading power of all the numerically computed ∆ ψ behave as 1 2q . In order to illustrate the described process we will first work it out for a particular value M/N ≈ 10 −1 (to be explicit M/N = 0.03 . . .). As it turns out the other values of M/N follow qualitatively the same behaviour. Afterwards we have a short discussion for M/N = 1, in which case there is some more analytical control for the rational case. We have done the above procedure for a range of values of M/N , in particular at M/N being of order 10 −2 , 10 −1 , 1, 10 1 and 10 2 . As it turns out the behaviour is very similar at all of these values, so let us just describe one of them, 10 −1 .
Having fixed the value of M/N and proceed with the next step, to solve the equation for ∆ ψ as a function of q (see eq. (A.5)): . We only show a small amount of these points such that the fits can be clearly distinguished. As can be seen, the power series up to 1/q 4 describes the results with high accuracy.
This equation is first solved numerically for both branches. Afterwards we fitted the data with the following power series: The results can be seen in figure 4. Note that we plotted only a small region of q for clarity, we solved all the cases from q = 3 at least up to q = 300. Apart from this we only show some of the data points, such that the fits can be clearly seen.
As mentioned before, we can see that both the ∆ ψ have a leading order term close to 1 2q . To check whether the leading order term is in fact equal to 1 2 we can check if the fitted value becomes closer to 1 2 for larger q. More precisely, we repeat the fitting procedure for different ranges of q, every time with a higher starting point, up to q = 300. The endpoint we fix at q = 415, which is large enough to obtain good fits for all the different starting points. In this way the behaviour for large q, captured by the leading order 1 q , will become more dominant. Let us call the fitted 1 q coefficient of the irrational branch f 1 such that: We follow the outlined procedure for this irrational case and show the results in figure 5.
As can be seen the leading order coefficient asymptotes to 1 2 for larger starting points of q. This means that the leading order behaviour in q of the ∆ ψ is captured by . This plot shows the difference between the fitted 1/q coefficient for ∆ I ψ , called f 1 (eq. (4.9)), and 1 2 as a function of the starting point of fitting. The endpoint of all the fits is fixed at q = 415, such that there are enough data points to fit for all the different starting points of q. It can be seen that the fitted coefficient asymptotes to 1 2 for increasing q ranges. The dots are the fitted coefficients, whilst the line approximates them well as 10 −2 q 3 .
means that in the strict q → ∞ limit the two branches coincide, and we can already expect their entropies to be the same in that limit.
Let us now compute the integrands from eq. (4.5), which we first do numerically and then compare those with the obtained fits. We show the results in figure 6. For clarity we show again only a subset of the computed q values. It becomes clear that i(q) > r(q) (which also holds outside the plotted range of q), i.e. the slope of the irrational entropy is always larger than the rational one. We can then expand the integrands for large q to find also the leading behaviour of the actual entropies. This yields (we plugged in the expression for A q−1 B): The fact that the first terms are the same is no coincidence, it arises from the equality of the first two terms in eq. (4.7), (the leading order is 1 2q ). Now we can find the entropies behaviour at large q: The first thing to notice is that we find in the strict q → ∞ that the entropies coincide, which we expected since then both the branches coincide. Apart from this we see that for q 1 (where the above approximation is valid) the entropy of the rational branch dominates the one of the irrational branch. But we can extend this now to arbitrary q, using the behaviour of the derivatives of the entropies (i(q) and r(q)). The numerical results shown in figure 6 show that i(q) > r(q). This means that the slope of the rational entropy is always smaller than the irrational one and hence we have S R > S I for all (finite) q, meaning that the rational branch is the dominant saddle.

Other M/N values
The procedure outlined in the previous section we have followed also for the other values of M/N around 10 −2 , 10 1 and 10 2 . As it turns out the behaviour is always (qualitatively) similar. We indeed find that the rational branch is the dominant saddle for all the M/N values. The larger the value of M/N the smaller the difference between the entropies becomes. This can be deduced from figure 3, where the branches flow to the same point as we increase M/N . Let us in this section mention the M/N = 1 case, where we have more analytical control for the rational branch. Let us first consider the rational branch. In this case we can always solve the exact dependence of the conformal dimension on q (see appendix A): Note that this coincides with the leading order found for other M/N values. The above expression means that the integrand in eq. (4.5), r(q), can be exactly computed: and hence we can compute the zero temperature entropy for the rational case: To fix the integration constant C we will consider the limit q → ∞. We can follow exactly [22], section II.C. There the results in a large q expansion are obtained: Note that v is an integration constant related to βJ [22]. From the expressions above it becomes clear that the q → ∞ limit reduces to free fermions. It also allows us to fix the constant C since: where in the last line we used the above observation that it should reduce to a free fermion entropy. We will see that the irrational entropy also behaves similarly as C 2 + f 1 q and hence for q → ∞ we find that it has the same integration constant, i.e. C 2 = C.
Just as in the M/N = 1 cases, we can't solve analytically the q dependence of the irrational solution. We did manage to find a good fit by ∆ I ψ = − 1 2q + 1 q−1 , which matches the numerical results well. Note it also reduces to 1 2q for large q. In figure 7 we plot the numerical results for both the rational and irrational cases. To conclude which of the entropies is bigger (i.e. which is the dominant saddle) we will investigate again the integrands as a function of q, see figure 8. From this plot we can see that the irrational integrand is (as in all the M/N cases) bigger than the rational one and as q increases their difference decreases. We can see that the irrational integrand (i(q)) is larger than the rational one (r(q)).
We now investigate the behaviour for large q. By using the approximate solution (see figure 7) we can obtain an expression for i(q) at large q. This is done by using the approximate solution in the integrand and expanding: From this we can find the leading order behaviour in the difference between the entropies: where the integration constants cancelled each other, such that only q dependent terms remain (recall that as q → ∞ we should obtain the free fermion entropy). We also used the rational entropy behaviour for large q, which can be obtained from eq. (4.15): We can then conclude that for large q: S I < S R . Further more since the slope of S R is always smaller than that of S I (since i(q) > r(q)) we find that this conclusion holds for any q. So, as for the other values of M/N , we find that the rational branch is the dominant saddle.

Chaos
In this section we will investigate the chaos or Lyapunov exponent of the model as a function of the ratio M/N . We will first review shortly the basics of such a computation and then move on to our model. The main tool for quantifying quantum chaos are so called Out of Time Order Correlators (OTOC) [5,[27][28][29][30]. For a more elaborate review of chaos and calculating these correlators see chapter 8 in [17], the first section of [8] and a discussion in [25].

JHEP01(2019)166
From a quantum mechanical point of view we can take two arbitrary Hermitian operators V and W and consider the commutator [W (it), V (0)] (with real time t ∈ R). The argument of the operator is imaginary since we consider it to be Euclidean time, as will be the case for our operators later on. The commutator describes the influence of small changes of V on later measurements of W (or the other way around). One particular indicator of these effects of chaos, which we will also use, puts the operators on the thermal circle [8]: where the brackets denote the thermal trace, the precise factors of β will not be important for us. For late enough times t (to be precise, between the dissipation and scrambling time [8]), quantum chaos dictates that this correlator will grow exponentially. By considering all the terms that arise in the above correlator one can show [8] that the exponential growth of the correlator arises due to the exponential behaviour of the related correlator: These out of time order correlators F (t) are usually studied in the context of quantum chaos, and we will use these as well. Schematically the OTOC eq. (5.2) behaves as [8,17,29]: The exponent λ L is called the Lyapunov exponent and it quantifies the chaos of the system. In the coming section our goal is to extract this Lyapunov exponent from the OTOCs. In general one can follow two approaches. The most obvious one is to compute the full four point function and continue these Euclidean correlators to real time. An easier option, however, is to consider the so called retarded kernel and its eigenfunctions [1,24,29]. In the context of ladder diagrams, kernels are the operators that add one more ladder to the diagram. For the OTOCs it has to be the retarded kernel due to the complex time contours specified by OTOCs similar to eq. (5.1). For a review of this procedure including the complex time contours, ladder diagrams and the application to ordinary SYK see [17].
The key idea of this procedure is to consider an exponentially growing OTOC on which the kernel(s) are acting. Under the assumption of this exponential growth one can find that it is precisely the eigenfunctions of the kernel with eigenvalue one that govern the chaotic behaviour. More intuitively, the growth rate of OTOC is determined by the demand that adding another ladder should not change the total sum. In the rest of this section we explain this procedure in more detail.

Retarded kernels
Let us now turn to our model and consider the four point functions (or OTOCs) that we want to compute. In [31,32] the chaos is also calculated for similar circumstances and we will comment upon this method at the end. t2 t4 (c) K 12 kernel Figure 9. Here we show the relevant kernels for the chaos computation. The subscripts of the kernels denote that they are elements of a matrix. The total matrix acts on a vector consisting of diagrams which starts with either two ψ or two φ lines.
We will consider the four point functions ψψψψ , ψψφφ , φφψψ and φφφφ . This is because acting with kernels on these diagrams will result in mixing between them and hence we can not consider them separately. The explicit OTOCs we will consider are of the form: where y is defined as y 4 = ρ(β). Diagrammatically these OTOCs are four point functions (ladder diagrams) with an arbitrary large amount of rungs. The other combinations of ψ and φ listed above have similar expressions and are denoted by F ψφ , F φψ and F φφ . The two subscripts of F ij denote the two incoming and two outgoing species, respectively. Let us now consider all the (retarded) kernels necessary for our model, which we draw in figure 9. Note that there is no kernel K 22 since there is no such interaction in the Lagrangian. It then becomes clear that acting with the K 12 and K 21 kernels causes mixing between the four point functions.
To get expressions for them we need first the necessary propagators in the diagrams. For the horizontal propagators we need the retarded ones (due to the complex time contours, see [17]): Recall that the arguments are imaginary since we consider complex Euclidean time. We can then use the finite temperature two point functions from eq. (4.3) (at τ > 0) to find: And similarly for φ:

JHEP01(2019)166
Lastly we need the ladder rung (lr) propagator, 1 which is obtained by simply continuing the Euclidean propagator τ → it + β 2 : Here x denotes ψ or φ and b x denotes A or B respectively. The form of this propagator is the same for fermions and scalars since we only need to consider τ > 0 here. We can then write down the expressions for the kernels. Note that each vertex gets a factor i from inserting it on a Lorentzian time fold in the contour and apart from this we also give K 11 and K 21 an additional minus sign due to the ordering of the contour (see also [31,32]). The resulting form of the kernels is: 2 The times t ij = t i − t j are shown in figure 9.

Integral matrix equation
Now that we have obtained the retarded kernels we go back to our four out of time order correlators. All together they obey an integral matrix equation as shown in figure 10, this is a generalization of the one particle version seen for example in [17]. In the figure we have put all the OTOCs in a four component vector seen on the very left (and right) side. These are exactly the OTOCs we named F ψψ , F ψφ , F φψ and F φφ before. Our (drawing) conventions are such that for the very left vector the times t 1 and t 2 are on the top left and bottom left of each four point function, respectively. The first vector on the right hand side denotes the free contributions to the four point functions. Clearly F ψφ and F φψ don't have these since there is no such free propagator. The matrix consists out of the retarded kernels discussed above and depicted in figure 9. Note that the matrix product in the last term also has an implicit convolution (which we will explicitly compute later on).
Quantum chaos implies that for late enough times these OTOCs will show exponentially growing behaviour, as discussed shortly in the introduction of this section. So let us make the following exponential growth ansatz: where i, j can denote ψ or φ and f ij denote functions of the time difference. Under the assumption of exponential growth the matrix equation figure 10 will simplify due to suppression of the zeroth order contributions. In fact, as one can easily check, the free diagrams will exponentially vanish compared to the exponential growth of the other terms. We are then left with the following equation: The F s now obey the ansatz eq. (5.10) and the matrix multiplication still involves the convolutions. However, we see that the problem can in fact be reduced to two identical problems because of the block diagonal structure. Hence we don't have to refer to the outgoing lines of the OTOCs (the second subscript of the F ij ) and consider simply: This leads to the following equations: where we have now explicitly written out the convolutions. The two equations are mixed and can be combined to give: F ψ (t 1 , t 2 ) = (K 11 * F ψ )(t 1 , t 2 ) + (K 12 * (K 21 * F ψ ))(t 1 , t 2 ) .

Lyapunov exponents
To actually solve the integrals we need to find the functions f i (t 12 ) in eq. (5.10) such that eq. (5.12) is satisfied. We take the following form of the functions, similar to [17,24,31,32]: Here the C i denote non-zero real constants, and we have h as the free exponential growth parameter. The Lyapunov exponent can be found by λ L = − 2πh β . The crucial integral for the computations is as follows 3 (5.17) Using the above identity we can then calculate the following integrals, reminiscent of eigenvalue equations:

JHEP01(2019)166
Note that here A and B are the coefficients of the (retarded) propagators, for which we only have an expression for A 2 B, see eq. (3.15). We can now use the above k ij along with eq. (5.18) in the integral equation eq. (5.15), to get: Of course, one could also have used the eigenfunctions (eq. (5.18)) first in eq. (5.13) and afterwards solved the mixing. Either way the equation resulting from the chaos regime is: We pick then some fixed M/N (which fixes ∆ ψ and A 2 B) and numerically solve this equation for the Lyapunov exponent λ L = − 2 π h β , which yields the solution h = −1. As it turns out, for h = −1 all the M/N dependence drops out and we find in fact maximal chaos for all values of M/N : Motivated by these numerical results we analytically checked whether k 11 + k 21 k 12 = 1 for h = −1. To do so we use the identities from eq. (3.13) and also the following: . (5.27) Using these identities and the expressions for A 2 B, we obtain the following simplified expressions, valid at h = −1: Hence even though k 11 and k 12 k 21 individually depend on M/N (since ∆ ψ does), the combined result exactly cancels. Lastly let us shortly mention another method of obtaining the chaos, outlined in [31,32]. In these articles the approach is to take the matrix of kernel eigenvalues, the k ij , and diagonalize it. Afterwards one of the eigenvalues is set to one. Let us consider this matrix: for which the resulting eigenvalues are k ± = 1 2 k 11 ± k 2 11 + 4 k 21 k 12 . The growing behaviour is found when k + = 1, which amounts to k 11 + k 21 k 12 = 1, consistent with our method.

JHEP01(2019)166 6 Discussion
In this article we have investigated new SYK-like models with M bosons and N fermions. The parameter M/N determines the behaviour of the model and for M = N our model is related to the supersymmetric SYK model. We have found that there are two families of solutions in the model, distinguished by their conformal dimensions which we plotted as a function of M/N in figure 3. For M = N the rational solution coincides with the supersymmetric solution found in [22]. Another interesting regime is the M/N → ∞ limit. We have shown that the solution of the fermionic conformal dimension in this limit is given by: This shows us that in this limit our model and SYK have the same behaviour, i.e. they have the same infrared fixed point. We have shown that the rational branch is the dominant saddle for arbitrary values of M/N (and arbitrary q). For the generic values of M/N we have shown this by first solving the problem numerically and then fitting it (with high accuracy). When M = N it can be done in a more analytical manner.
Apart from this we investigated the Lyapunov exponent and found it to be λ L = 2πT , independent of M/N . This is due to some non trivial cancellations in the M/N dependences. Due to the maximal chaos the model has a holographic interpretation as a black hole. It would be interesting to understand the role of M/N in this holographic description. Concretely it would be interesting to find the Schwarzian for this model, in particular the coefficient in front of the Schwarzian action, related to the heat capacity, and its dependence on M/N . We leave this for future research.

A The model for a q-point interaction
In this appendix we will shortly show how the model and some results change when we consider an interaction vertex of degree q, so an interaction with one boson and q − 1 fermions. The integer q is supposed to be odd, but later we will continue it to arbitrary real values. The model we consider in the main text has q = 3. When we apply this, the coupling i 2 C aij φ a ψ i ψ j goes to i (q−1)! C ai 1 i 2 ...i q−1 φ a ψ i 1 . . . ψ i q−1 . Integrating out the bosons would lead to a Hamiltonian with a vertex containing 2(q − 1) fermions. One can thus relate this to the interaction parameter q SYK in [24] by (q − 1) = 1 2 q SYK . We further assume that q N .

JHEP01(2019)166
The disorder average is now chosen as follows: Once again we can take the conformal form (see eq. (3.9)) for the two point functions. By following the computations done for the q = 3 case we find that the conformal symmetry is present under the condition that (compare to eq. (3.16)): The equations for the constants in the two point functions, A q−1 B, yield (compared to eq. (3.15)): . (A.4) The resulting transcendental equation for the conformal dimensions reads: One may check that the rational value ∆ ψ = 1 2q is a solution at M = N . Furthermore as we take the limit M/N → ∞ we see that we find the solution: . (A.6) Using the above found relation (q − 1) = 1 2 q SYK we see that this is equal to 1 q SYK .
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.