Tunable Quantum Chaos in the Sachdev-Ye-Kitaev Model Coupled to a Thermal Bath

TL;DR

This paper investigates how coupling SYK models to a thermal bath affects quantum chaos indicators like Lyapunov exponent and butterfly velocity, revealing tunable chaos and thermalization behaviors.

Contribution

It introduces a model of SYK chains coupled to a thermal bath, demonstrating controllable quantum chaos properties and spatially dependent thermalization effects.

Findings

Lyapunov exponent decreases with stronger bath coupling

Butterfly velocity transitions from 2T to linear T dependence

Spatial dependence of chaos metrics when only chain ends are coupled

Abstract

The Sachdev-Ye-Kitaev (SYK) model describes Majorana fermions with random interaction, which displays many interesting properties such as non-Fermi liquid behavior, quantum chaos, emergent conformal symmetry and holographic duality. Here we consider a SYK model or a chain of SYK models with $N$ Majorana fermion modes coupled to another SYK model with $N^2$ Majorana fermion modes, in which the latter has many more degrees of freedom and plays the role as a thermal bath. For a single SYK model coupled to the thermal bath, we show that although the Lyapunov exponent is still proportional to temperature, it monotonically decreases from $2\pi/\beta$ ($\beta=1/(k_BT)$, $T$ is temperature) to zero as the coupling strength to the thermal bath increases. For a chain of SYK models, when they are uniformly coupled to the thermal bath, we show that the butterfly velocity displays a crossover from a $\sqrt{T}$-dependence at relatively high temperature to a linear $T$-dependence at low temperature, with the crossover temperature also controlled by the coupling strength to the thermal bath. If only the end of the SYK chain is coupled to the thermal bath, the model can introduce a spatial dependence of both the Lyapunov exponent and the butterfly velocity. Our models provide canonical examples for the study of thermalization within chaotic models.

Equations

$$H_{\chi }=\frac{1}{4!}\sum _{ijkl}^N{J}_{ijkl}{\chi }_i{\chi }_j{\chi }_k{\chi }_l,$$

$$\overline{J_{ijkl}}=0,\overline{J_{ijkl}^2}=\frac{3!J^2}{N^3}.$$

$$G(i\omega {)}^{-1}=-i\omega -\Sigma (i\omega ),\Sigma (\tau )=J^{2}G(\tau {)}^3.$$

$$G(\tau )=b\frac{\text{sgn}(\tau )}{|\tau {|}^{\frac{1}{2}}},$$

$$J^2{b}^4=\frac{1}{4\pi }.$$

$$G_{\beta }(\tau )=b\frac{\text{sgn}(\tau )}{{\left|\frac{\beta }{\pi }\sin \left(\frac{\pi \tau }{\beta }\right)\right|}^{\frac{1}{2}}}.$$

$$\frac{1}{N}\mathcal{F}({\tau }_1,{\tau }_2,{\tau }_3,{\tau }_4)=\frac{1}{N^2}\sum _{j,k=1}^N⟨{\mathcal{T}}_{\tau }{\chi }_j({\tau }_1){\chi }_j({\tau }_2){\chi }_k({\tau }_3){\chi }_k({\tau }_4)⟩-G({\tau }_{12})G({\tau }_{34}).$$

$$F_{AB}(t)=\text{Tr}[y{A}^{†}(t)y{B}^{†}(0)yA(t)yB(0)],\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{where}\text{leavevmode}\text{leavevmode}y=e^{-\frac{1}{4}\beta \hat{H}}.$$

$$\frac{S}{N}=-\frac{1}{2}\log \det ({\partial }_{\tau }-\Sigma )+\frac{1}{2}\int d{\tau }_{1}d{\tau }_2\left[\Sigma ({\tau }_1,{\tau }_2)G({\tau }_1,{\tau }_2)-\frac{J^2}{4}G({\tau }_1,{\tau }_2{)}^4\right].$$

$$\frac{S}{N}\propto {\int }_0^{\beta }d\tau \frac{1}{2}\left[({ϵ}^{\prime \prime }{)}^2-{\left(\frac{2\pi }{\beta }\right)}^2({ϵ}^{\prime }{)}^2\right].$$

$$\frac{S}{N}=-\frac{\#}{J}\int d\tau \left\{\tan \left(\frac{\pi f}{\beta }\right),\tau \right\},\text{where}\{f,\tau \}\equiv \frac{f^{\prime \prime \prime }}{f^{\prime }}-\frac{3}{2}{\left(\frac{f^{\prime \prime }}{f^{\prime }}\right)}^2.$$

$$H=H_{\chi }+H_{\psi }+H_c,$$

$$H_c=\frac{1}{4}\sum _{ijkl}{u}_{ijkl}{\chi }_i{\chi }_j{\psi }_k{\psi }_l,$$

$$\overline{J_{ijkl}}=0,\overline{J_{ijkl}^{{}^{\prime }}}=0,\overline{u_{ijkl}}=0,$$

$$\overline{J_{ijkl}^2}=\frac{3!J^2}{N^3},\overline{J_{ijkl}^{{}^{\prime }2}}=\frac{3!J^2}{N^6},\overline{u_{ijkl}^2}=\frac{2!u^2}{N^5}.$$

$$G_{\chi }(\omega {)}^{-1}=-i\omega -{\Sigma }_{\chi }(\omega ),\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{\Sigma }_{\chi }(\tau )=J^2{G}_{\chi }^3(\tau )+u^2{G}_{\psi }^2(\tau )G_{\chi }(\tau ).$$

$$G_{\psi }(\omega {)}^{-1}=-i\omega -{\Sigma }_{\psi }(\omega ),\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{\Sigma }_{\psi }(\tau )=J^2{G}_{\psi }^3(\tau ).$$

$$G_{\chi }(\tau )=a\frac{\text{sgn}(\tau )}{|\tau {|}^{\frac{1}{2}}},\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{G}_{\psi }(\tau )=b\frac{\text{sgn}(\tau )}{|\tau {|}^{\frac{1}{2}}},$$

$$a^4{J}^2+u^2{a}^2{b}^2=\frac{1}{4\pi },\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{b}^4{J}^2=\frac{1}{4\pi }.$$

$$G_{\chi }(\tau )=a{\left[\frac{\pi }{\beta \sin \frac{\pi \tau }{\beta }}\right]}^{\frac{1}{2}}\text{sgn}(\tau ),\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{G}_{\psi }(\tau )=b{\left[\frac{\pi }{\beta \sin \frac{\pi \tau }{\beta }}\right]}^{\frac{1}{2}}\text{sgn}(\tau ).$$

$$F_{\chi \chi }(t_1,t_2)=\text{Tr}[y{\chi }_i(t_1)y{\chi }_j(0)y{\chi }_i(t_2)y{\chi }_j(0)],\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}y=\rho (\beta {)}^{\frac{1}{4}}.$$

$$F_{\psi \psi }(t_1,t_2)=\text{Tr}[y{\psi }_i(t_1)y{\psi }_j(0)y{\psi }_i(t_2)y{\psi }_j(0)],\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}y=\rho (\beta {)}^{\frac{1}{4}}.$$

$$F_{\chi \chi }(t_1,t_2)=\int d{t}_{3}d{t}_4{K}_{R,\chi \chi }(t_1...t_4)F_{\chi \chi }(t_3,t_4),$$

$$K_{R,\chi \chi }(t_1...t_4)=3J^2{G}_{R,\chi }(t_{13})G_{R,\chi }(t_{24})G_{lr,\chi }(t_{34}{)}^2+u^2{G}_{R,\chi }(t_{13})G_{R,\chi }(t_{24})G_{lr,\psi }(t_{34}{)}^2,$$

$$G_{R,\chi }(t)=\sqrt{2}a\theta (t){\left[\frac{\pi }{\beta \sinh \frac{\pi t}{\beta }}\right]}^{\frac{1}{2}},\text{leavevmode}\text{leavevmode}\text{leavevmode}{G}_{lr,\chi }(t)=a{\left[\frac{\pi }{\beta \cosh \frac{\pi t}{\beta }}\right]}^{\frac{1}{2}},$$

$$G_{R,\psi }(t)=\sqrt{2}b\theta (t){\left[\frac{\pi }{\beta \sinh \frac{\pi t}{\beta }}\right]}^{\frac{1}{2}},\text{leavevmode}\text{leavevmode}\text{leavevmode}{G}_{lr,\psi }(t)=b{\left[\frac{\pi }{\beta \cosh \frac{\pi t}{\beta }}\right]}^{\frac{1}{2}}.$$

$$F_{\chi \chi }(t_1,t_2)=\frac{e^{-h\frac{\pi }{\beta }(t_1+t_2)}}{{\left[\cosh \frac{\pi }{\beta }{t}_{12}\right]}^{\frac{1}{2}-h}},$$

$$1=\frac{\Gamma (\frac{5}{2})\Gamma (\frac{1}{2}-h)}{\Gamma (\frac{3}{2})\Gamma (\frac{3}{2}-h)}\times 4\pi \left(J^2{a}^4+\frac{u^2{a}^2{b}^2}{3}\right).$$

$$-h=1+\frac{1}{2}{k}^2-\frac{\sqrt{k^4+4k^2}}{2},\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{with}\text{leavevmode}\text{leavevmode}k=\frac{u^2}{J^2}.$$

$${\lambda }_L=\frac{2\pi }{\beta }\left(1-\frac{\sqrt{k^4+4k^2}-k^2}{2}\right).$$

Full text

1 Introduction

As a concrete solvable model displaying maximal quantum chaos, non-Fermi liquid behavior and holographic duality, the Sachdev-Ye-Kitaev (SYK) model proposed by Kitaev Kitaev2 based on the early work of Sachdev and Ye SY has been widely studied recently both on the field-theory side Kitaev2 ; Comments ; spectrum1 ; spectrum2 ; spectrum3 ; Liouville ; Liouville2 ; SYK new ; SYK new2 ; SYK new3 and the gravitational side Comments ; bulk Yang ; bulk spectrum Polchinski ; bulk2 ; bulk3 ; bulk4 ; bulk5 ; syk-bh ; SYK g new1 ; SYK g new2 ; SYK g new3 . On the field theory side, the model contains $N$ Majorana fermions $\chi_{i}\ (i=1,...,N)$ with random four-fermion interaction terms. The Hamiltonian is

[TABLE]

where the normalization of $\chi_{i}$ is given by anticommutation relation $\{\chi_{i},\chi_{j}\}=\delta_{ij}$. The coupling constants $\{J_{ijkl}\}$ are antisymmetric with respective exchanging of indices. Their mean vanishes and the variance is finite as

[TABLE]

The interaction is all-to-all, thus the model is usually viewed as a (0+1)-d quantum mechanics model. The specific choice of $N$ scaling in Eq. (2) leads to an elegant large-$N$ structure, with which the two-point and four-point correlation functions can be calculated analytically Comments . The exact solution in large-$N$ limit displays a number of very intriguing properties, as we briefly summarized below.

First, it displays an emergent conformal symmetry and the non-Fermi liquid behavior. To the leading order of $\frac{1}{N}$ expansion, the imaginary-time two-point function $G(\tau)\delta_{ij}\equiv\langle\mathcal{T}_{\tau}\chi_{i}(\tau)\chi_{j}(0)\rangle$ satisfies a simple form of Schwinger-Dyson equations as

[TABLE]

In the infrared limit, one drops the $-i\omega$ term, and the Schwinger-Dyson equations can be solved by using following ansatz

[TABLE]

which leads to

[TABLE]

A Fourier transformation of the Green’s function shows the divergence of spectral function as $1/\omega^{1/2}$ at low energy, which signals a non-Fermi Liquid behavior. The form of the Green’s function also indicates that the system acquires an emergent conformal symmetry in the IR limit with the scaling dimension of $\chi$ equals to $1/4$, which can also be seen from the fact that the fix-point action only has the four-fermion interaction terms. With the help of this conformal symmetry and by utilizing conformal mapping $\tau\rightarrow\tan\frac{\pi\tau}{\beta}$, the low temperature Green’s function (with $\beta J\gg 1$) is then given by

[TABLE]

Second, it shows maximally chaotic behavior. With the knowledge of the two-point Green’s function at finite temperature, one can further calculate the four-point Green’s function of the SYK model within the large-$N$ expansion. At the leading order $\mathcal{O}(N^{0})$, the four-point function is given by its disconnected part, as the connected part starts at order $\mathcal{O}(\frac{1}{N})$

[TABLE]

With the analytical continuation of $\mathcal{F}(\tau_{1},\tau_{2},\tau_{3},\tau_{4})$ to $\mathcal{F}(\frac{3\beta}{4}+it,\frac{\beta}{4}+it,\frac{\beta}{2},0)$, one can calculate the out-of-time-ordered correlation (OTOC) function for $\chi_{i}$ and $\chi_{j}$. The definition for a regularized OTOC for operator $A$ and $B$ is:

[TABLE]

which diagnoses the quantum butterfly effect Kitaev1 ; bh1 ; bh2 ; bh3 and has been applied to various models recently otoc cft ; otoc anyon ; OTOC-Boson Hubbard ; OTOC-Keldysh ; OTOC-Luttinger ; OTOC-MBL2 ; OTOC-MBL3 ; OTOC-MBL4 ; OTOC-MBL5 ; OTOC-MBL6 ; OTOC-Ruihua MBL1 ; OTOC-Yao ; OTOC-classical ; OTOC-quantum ; OTOC-quantum channel ; OTOC-protocal ; OTOC new . There are also experimental measurements of OTOC in different systems OTOC-exp1 ; OTOC-exp2 . The OTOC in many chaotic systems has a universal behavior as $F_{AB}(t)\sim c_{0}-c_{1}\exp(\lambda_{L}t)$ at a time scale called the scrambling time $t_{s}$. $\lambda_{L}$ defines the Lyapunov exponent for a quantum system. Assuming $t_{s}\gg t_{d}$ ($t_{d}$ is defined as the decay time for two point function), it is proved that $\lambda_{L}$ is bounded by $2\pi/\beta$, and the models with holographic duality is believed to saturate the bound Kitaev2 ; bh1 ; bh2 ; bh3 ; prove . For the SYK model, one can explicitly show that $F(t)\equiv F_{\chi_{i}\chi_{j}}(t)\sim-\beta J\exp(2\pi t/\beta)$ when $t\gg t_{d}=\beta$. The Lyapunov exponent extracted from it is exactly $2\pi/\beta$.

Thirdly, a gravitational model shares the same effective action as the low-energy theory of the SYK model. In the SYK model, by using the replica trick to treat the disorder average, a non-local action can be deduced after introducing bi-local fields $G(\tau_{1},\tau_{2})$ and $\Sigma(\tau_{1},\tau_{2})$:

[TABLE]

The saddle point equations for $G(\tau_{1},\tau_{2})$ and $\Sigma(\tau_{1},\tau_{2})$ (by assuming the translational invariance in time) are the same as shown in Eq. (3). Using saddle-point solutions, one can show there is a finite zero-temperature entropy which is another signature of a non-Fermi Liquid resulted from the large degeneracy of ground states numerics wenbo . Expanding $G$ and $\Sigma$ around their saddle point solutions gives the low-energy effective action for the SYK model. Note that the low energy physics of the SYK model is dominated by the reparametrization modes due to the emergent conformal symmetry, and a reparametrization to the time variable $\tau\rightarrow f(\tau)$ acts on the $G(\tau_{1},\tau_{2})$ field as $G(\tau_{1},\tau_{2})\rightarrow(f^{\prime}(\tau_{1})f^{\prime}(\tau_{2}))^{\frac{1}{4}}G(f(\tau_{1}),f(\tau_{2}))$. For small deformation $\tau\rightarrow\tau+\epsilon(\tau)$, the effective action can be approximated by an action of the fluctuation $\delta_{\epsilon}G$, which has an elegant form as the Schwarzian action:

[TABLE]

For finite transformation $\tau\rightarrow f(\tau)$, the action can be written as

[TABLE]

$\#$ is some constant number given in Ref.Comments . This Schwarzian effective action roots in the conformal symmetry of the SYK model and can give rise to the maximal Lyapunov exponent. Interestingly, on the gravity side, the same action appears for the dilaton gravity theory in two-dimensional near anti-de Sitter spacetime (NAdS2) bulk Yang .

Motivated by the aforementioned intriguing properties of the SYK model, recently there appear many interesting generalizations of the SYK model numerics wenbo ; Yingfei1 ; Yingfei2 ; generalization 1 ; wenbo susy ; susy2 ; no disorder1 ; no disorder2 ; no disorder3 ; no disorder4 ; no disorder5 ; no disorder6 ; no disorder7 ; no disorder8 ; no disorder9 ; no disorder10 ; Altman ; thermal transport ; high-D1 ; high-D2-con ; yyz condensation ; sk jian ; generalization 2 ; transition1 ; our ; Balent ; thick ; new gsyk1 ; no disorder new ; susy3 . The models with $U(1)$ symmetry are studied numerically for both bosons and fermions in Ref.numerics wenbo by exact diagonalization, which shows the nearly degeneracy of the many-body spectrum at low energy for the fermionic model and the spin-glass state for the bosonic model. There are also efforts toward a more precise holographic description for which a supersymmetric version of the SYK model is studied wenbo susy ; susy2 ; susy3 . By coupling the SYK models Yingfei1 ; Yingfei2 , an lattice model with spatial degree of freedom can be realized, where the butterfly velocity and the diffusion constant can be calculated and are compared with the holographic result. The issue of chaotic to non-chaotic transiton is considered with several different generalization of the SYK models Altman ; sk jian ; yyz condensation ; our ; Balent . Among all these generalizations, the non-Fermi Liquid phases given by the SYK-type models all possess a maximal chaotic behavior and Lyapunov exponent saturating to $2\pi/\beta$, with only one exception that contains time dependent interaction Comments ; thick .

In this paper, we consider a SYK model or a chain of SYK models with $N$ Majorana fermion modes coupled to another SYK model with $N^{2}$ Majorana fermion modes. The latter contains many more degrees of freedom and can be viewed as as a kind of thermal bath. Our model is time independent and this model is also solvable in the large-$N$ limit. In view of the three properties mentioned above, we will show that, on one hand, this model still displays an emergent conformal symmetry and a non-Fermi liquid behavior, but on the other hand, it is not maximally chaotic. In fact, we will show that the Lyapunov exponent of this model can be tuned to any value between zero and $2\pi/\beta$. Whether this model also has a gravitational dual remains unclear.

In section 2, we discuss a single SYK model coupled to the thermal bath. We show analytically illustrate that $\lambda_{L}$ for the small system monotonically decreases from $2\pi/\beta$ to zero as the coupling strength to the thermal bath increases. In section 3, we consider a chain of SYK models. When the chain is uniformly coupled to the thermal bath, the butterfly velocity displays a crossover from $\sqrt{T}$-dependence at relatively high temperature to a new linear $T$-dependence at low temperature. If only the end of the SYK chain is coupled to the thermal bath, we found that there is a spatial dependence of both the Lyapunov exponent and the butterfly velocity.

2 $(0+1)$-d SYK Model Coupled to a Thermal Bath

2.1 The Model and the Two-Point Functions

In this section, we introduce a generalized version of the SYK model that contains a small SYK model with $N$ Majorana fermions (denoted by $\chi_{i}$, $i=1,2,...,N$) to a large SYK model with $N^{2}$ Majorana fermions (denoted by $\psi_{i}$, $i=1,2,...,N^{2}$). We will call the two subsystem as SYKχ and SYKψ, respectively. In fact, the number of Majorana fermions in the large SYK model does not have to be $N^{2}$, and it can be generally $N^{\alpha}$ as long as $\alpha>1$ such that it dominates over $N$ in the large $N$ limit. Here we choose the number of Majorana fermions in the larger cluster to be $N^{2}$ just for concreteness. The Hamiltonian of the coupled SYK system is written as:

[TABLE]

where both $H_{\chi}$ and $H_{\psi}$ take the same form as Eq. (1), with $\{J_{ijkl}\},\{J^{{}^{\prime}}_{ijkl}\}$ being the random couplings in $H_{\chi}$ and $H_{\psi}$, respectively. $H_{c}$ is the coupling Hamiltonian defined as

[TABLE]

where $u_{ijkl}$ are also random coupling numbers. $\{J_{ijkl}\},\{J^{{}^{\prime}}_{ijkl}\},\{u_{ijkl}\}$ are antisymmetric random variables with zero mean

[TABLE]

and their variances are

[TABLE]

The reason for making this specific choice of $N$ dependence of the variances will be clear soon because it gives a nice structure for the theory in the large-$N$ limit.

As shown in Fig. 1 (a), we use straight line and wavy line to denote their Green’s functions $G_{\chi}(\tau)=\langle\mathcal{T}_{\tau}\chi_{i}(\tau)\chi_{i}(0)\rangle$ and $G_{\psi}(\tau)=\langle\mathcal{T}_{\tau}\psi_{i}(\tau)\psi_{i}(0)\rangle$, respectively, where $\mathcal{T}_{\tau}$ denotes time-order operator in the imaginary time $\tau$. Note that the Green’s functions are diagonal in term of the fermion indices because of the disorder average. In the large-$N$ limit, the Green’s functions can be be determined by diagrammatic method and by solving the Schwinger-Dyson equations. For the SYKχ, one would expect two contributions to its self-energy $\Sigma_{\chi}(\tau)$, as given by two diagrams in Fig. 1 (b). The choices of the orders of $N$ in Eq. (15) ensures that the same two diagrams are both of order $\mathcal{O}(N^{0})$, which is the lowest order in $\frac{1}{N}$ for the Green’s function. From these diagrams we can obtain the Schwinger-Dyson equations for the SYKχ system as

[TABLE]

For the SYKψ system, in the large-$N$ limit, one may also expect two contributions for its self-energy $\Sigma_{\psi}(\tau)$, as shown in Fig. 1. However, the choices in Eq. (15) results in the suppression of the second diagram by $\frac{1}{N}$, thus it can be neglected at the leading order. This large-$N$ structure physically makes sense, since physically the properties of the larger system should not be affected by the small system at the leading order of $\frac{1}{N}$. Then we have the Schwinger-Dyson equations for the SYKψ system:

[TABLE]

In the strong coupling limit, we first drop the $-i\omega$ term, and the Green’s functions obey the form as the original single SYK model

[TABLE]

where the coefficients $a$ and $b$ are determined by

[TABLE]

Note that the coefficients are different from the original SYK model, which is the key to the discussion below. Consequently, the finite temperature version of the Green’s functions are

[TABLE]

2.2 Four Point Functions and the Tunable Lyapunov Exponent

The Lyapunov exponent is defined from four point function OTOC. To calculate the four-point function, one can solve the self-consistency equations for the four point functions in the real time to first obtain the asymptotic behavior in the chaos limit as shown in Ref.Kitaev1 ; Comments , from which one can extract the Lyapunov exponent perfectly. We use $F_{\chi\chi},F_{\psi\psi}$ and $F_{\chi\psi}$ to denote the four point functions of four $\chi$’s, of four $\psi$’s and of two $\chi$’s and two $\psi$’s. Let us consider the following OTOC in the real time,

[TABLE]

[TABLE]

where fermions are separated by a quarter of the thermal circle, at the lowest order $\mathcal{O}(N^{0})$, the four point function is given by the disconnected part, and what we are interested in is the sub-leading part that describes chaotic behavior. Below, we use $F(t_{1},t_{2})$ to refer to the sub-leading part of the four point function without further indication. As in the original SYK model, they are determined by the ladder diagrams in the large-$N$ limit.

We use self-consistency equations to determine the asymptotic behavior of $F_{\chi\chi}$, i.e., the function $F_{\chi\chi}$ is an eigenfunction of $K_{R,\chi\chi}$ with eigenvalue one:

[TABLE]

where $K_{R,\chi\chi}(t_{1}...t_{4})$ is the retarded kernel for evaluating the four point function. We shall first analyze the structure of the retarded kernel. One finds that the kernel consists of two parts at lowest order given by the first and second diagram in Fig. 2 (a). Other possibilities, like the third diagram in Fig. 2 (a), is in fact of the order $\frac{1}{N^{2}}$. At leading order for the connected part of four-point function, it can be neglected. From the diagrams, one can obtain the retarded kernel $K_{R,\chi\chi}(t_{1}...t_{4})$ as

[TABLE]

where $G_{R}(t)$ is the real time retarded correlator and $G_{lr}(t)\equiv iG(it+\frac{\beta}{2})$ is the Wightmann correlator. They are given by

[TABLE]

[TABLE]

By taking an ansatz form of

[TABLE]

and substituting Eq. (24) - (26) into Eq. (23) while requiring that the eigenvalue of the kernel $k_{R}(h)=1$, one finds

[TABLE]

The value of the parenthesis in the above equation equals to $1/4\pi$ in the original SYK model, while here the value depends on $u$ explicitly and implicitly via coefficients $a$ and $b$. One can solve $h$ from Eq. (19) and Eq. (28), the result is

[TABLE]

One can see that Since $\lambda_{L}=-h\frac{2\pi}{\beta}$, we find that

[TABLE]

When $k=u^{2}/J^{2}=0$, the two systems decouple, and the Lyapunov exponent of $F_{\chi\chi}$ recovers $2\pi/\beta$, which is simply the maximumly chaotic value in the SYK model. However, as one increases the interaction between the two systems, i.e., $u^{2}/J^{2}$ increases, the Lyapunov exponent decreases. When $u^{2}/J^{2}\rightarrow\infty$, $\lambda_{L}\rightarrow 0$. Fig. 3 shows one of the central results of this work where we plot the dependence of $\lambda_{L}$ on $u^{2}/J^{2}$. For small $u^{2}/J^{2}$, the Lyapunov exponent decreases linearly with $u^{2}/J^{2}$:

[TABLE]

The above calculation deals with the four point function in the real time and in the chaos limit directly. One can also first calculate the exact four point function in the imaginary time and then continue it to the real time, which will lead to the same result Wenbo .

If we take the Lyapunov exponent $\lambda_{L}$ as a measurement of the chaos in one system, then this result tells us that we can tune chaotic behavior in one SYK system by changing the magnitude of its interaction with a much larger system. The underlying mechanism is based on the hierarchy of the scrambling times between these two systems. In the SYK model, the scrambling time $t_{s}$ is proportional to $\beta\log N$, where $N$ is the number of Majorana fermions in the system. The large difference in the size of the two subsystems results in a large difference in the scrambling time. At the time scale when the SYKχ system enters the chaos region, the SYKψ system has not scrambled at all. Thus, through the coupling between the two systems, the chaotic behavior of the SYKχ system is weaken by the larger system SYKψ.

Now let us turn to the OTOC $F_{\psi\psi}(t_{1},t_{2})$ of four $\psi$ fermions. As illustrated in the last section, at the lowest order of $1/N$, the SYKψ system should not be affected by the small system. We will find this is indeed the case here. One would find that the leading order of the connected part $F_{\psi\psi}(t_{1},t_{2})$ is at order $\mathcal{O}(1/N^{2})$, and the only contribution to the kernel at this order is shown in Fig. 4 (a). From the self-consistency equation

[TABLE]

and by taking the ansatz form

[TABLE]

with the eigenvalue $k_{R}(h^{\prime})=1$, one finds

[TABLE]

Solving this one obtains $h^{\prime}=-1$, which means that the Lyapunov exponent $\lambda_{L}$ is $2\pi/\beta$. The Lyapunov exponent of $F_{\psi\psi}$ is the same as the one in the original SYK model, which matches our previous expectation.

If we want to understand the behavior of $F_{\chi\chi}$ at next order $\mathcal{O}(\frac{1}{N^{2}})$ precisely, we need to take into account the $\mathcal{O}(\frac{1}{N})$ correction to the correlators, and analyze the kernel structure as before. Some useful hint can be obtained by simply looking at the third diagram in Fig. 1 (a), which is the $1/N^{2}$ order contribution to the kernel $K_{R,\chi\chi}$. It contains the four-point function $F_{\psi\psi}$ as a propagator inside (see Fig. 4 (b)), and since the Lyapunov exponent of $F_{\psi\psi}$ reaches the maximal value $2\pi/\beta$, we anticipate that at order $1/N^{2}$, $F_{\chi\chi}(t_{1},t_{2})$ will also display a Lyapunov exponent $2\pi/\beta$. This gives an exponential growth contribution as $\frac{1}{N^{2}}\exp(\frac{2\pi}{\beta}t)$ which reveals the maximally chaotic behavior at longer time when both systems are scrambled. Similarly, since the lowest order of the connected part of OTOC $F_{\psi\chi}$ is at $\mathcal{O}(\frac{1}{N^{2}})$ and it also contains $F_{\psi\psi}$ as an inner propagator, we anticipate that $F_{\psi\psi}$ will also grow as $\exp{\frac{2\pi}{\beta}t}$ at order $\frac{1}{N^{2}}$.

From this point of view, one can see different chaotic behavior at different time scale. As far as the small system is concerned, there are two scrambling time scales, with the first one at the shorter time and the second one at the longer time. The small system is not maximally chaotic at the first scrambling time, because of the non-chaotic environment at that time scale. And it is maximally chaotic at the second scrambling time at longer time scale, because of the maximally chaotic environment at that time scale.

2.3 The Effective Action

In this section we will discuss an alternative effective action derivation for the change of the Lyapunov exponent in a perturbative manner. This will be useful for the later discussion of the SYK chain coupled to the environment.

The derivation of the effective action utilizes the original fermion path integral formalism and carries out the disorder average. By introducing bi-local fields $G_{\chi}(\tau_{1},\tau_{2}),G_{\psi}(\tau_{1},\tau_{2})$ and Lagrange multiplier fields $\Sigma_{\chi}(\tau_{1},\tau_{2}),\Sigma_{\psi}(\tau_{1},\tau_{2})$ that sets $G_{\chi}(\tau_{1},\tau_{2})=\frac{1}{N}\Sigma_{i}\chi_{i}(\tau_{1})\chi_{i}(\tau_{2})$ and $G_{\psi}(\tau_{1},\tau_{2})=\frac{1}{N^{2}}\Sigma_{i}\psi_{i}(\tau_{1})\psi_{i}(\tau_{2})$, one obtains the nonlocal action as

[TABLE]

The saddle point equations of the fields are the same as in Eq. (16) and Eq. (17). We expand the fields around their saddle point solution as

[TABLE]

where $g(\tau)$, $g^{\prime}(\tau)$, $\sigma(\tau)$, $\sigma^{\prime}(\tau)$ denote the fluctuation of the fields. By keeping up to the quadratic order, one obtains the effective action

[TABLE]

Integrating out $\sigma$ and $\sigma^{\prime}$, and to the leading order of $\frac{u^{2}}{J^{2}}$, one obtains

[TABLE]

in which

[TABLE]

and

[TABLE]

$K^{\prime}$ does not have $J$ and $u$ dependences since the coefficient $b$ in $G_{\psi,s}(\tau)$ will cancel out the $J^{2}$ prefactor, while $K$ has $J$ and $u$ dependences.

The exact four point function can be derived by first diagonalizing the kernels $K$ and $K^{\prime}$, with their eigenvalues denoted by $k(h,n)$ and $k^{\prime}(h,n)$ and the eigenfunction corresponding to $(h,n)$ denoted by $\Psi_{h,n}$. For the SYKψ system, at the leading order, one only needs to consider the first term in the effective action. The connected part of the imaginary time four-point function for SYKψ system can be written as

[TABLE]

However, for $h=2$, one has $k^{\prime}(h,n)=1$, which leads to the vanishing of effective action for $h=2$ modes as well as a divergence of the four-point function in above equation. The divergence comes from the reparametrization modes which are soft modes in the conformal limit. For them, one must consider the correction away from the conformal limit to obtain a meaningful result. On the other hand, these modes are the dominant part of the four-point function in the long time limit, and result in the exponential growth with the Lyapunov exponent $2\pi/\beta$ for the OTOC. This is the same story as the original SYK model, as well as the SYKψ system. However, for the SYKχ system, the interaction with the bath makes the difference. If we integrate out $g^{\prime}$ in the effective action, at the leading order, the effective action for $g$ is

[TABLE]

With this, one obtains the expression for the connected part of the imaginary time four-point function

[TABLE]

For $h=2$, $k^{\prime}(h,n)=1$, but there is no divergence due to $\frac{u^{2}}{J^{2}}\neq 0$. Since the divergence is removed, the long time behavior of the four-point function should be a result of both $h=2$ and $h\neq 2$ modes. In principle, the asymptotic behavior of OTOC in the chaos limit contains contributions from the $h\neq 2$ modes. However, for $u^{2}/J^{2}\ll 1$, the reparametrization modes still compose the dominant part of the four-point function, thus we can first focus on the $h=2$ part, and then consider the correction from $h\neq 2$ part perturbatively. We will show below that this leads to the same result as in previous section for the perturbative regime $u^{2}/J^{2}\ll 1$.

Concentrating on the $h=2$ part, we need to write an explicit action for the reparametrization modes. For $u=0$, the SYKχ and SYKψ systems are decoupled, thus we should have two sets of reparametrization mode $\epsilon$ and $\epsilon^{\prime}$, corresponding to the reparametrization of $G_{\chi}$ and $G_{\psi}$, respectively. For finite $u$, since $G_{\chi}$ and $G_{\psi}$ are coupled through Eq. (16), they should be reparametrized together. However, since we are working in the $u^{2}/J^{2}\ll 1$ regime, we shall still keep two sets of reparametrization modes, and treat the $u^{2}$ term as the first order perturbation that couples the two sets of modes together. Explicitly, the two sets of reparametrization modes are introduced as:

[TABLE]

with functions $f_{n}$ defined in Ref.Comments and $\int f_{n}f_{m}\propto|n|(n^{2}-1)\delta_{m,-n}$. The eigenvalues of $K^{\prime}$ for the reparametrization modes are approximated by $1-|n|\frac{\sqrt{2}\alpha_{K}}{2\pi\beta J}$, while the eigenvalues of $K$ for the reparametrization modes are approximated by $\frac{J^{2}}{J^{2}+u^{2}}-|n|\frac{\sqrt{2}\alpha_{K}}{2\pi\beta J}$, in which $\alpha_{K}$ is a numerical factor Comments . The effective action is then given by:

[TABLE]

where we neglected a term proportional to $N$ for quadratic term of $\epsilon^{\prime}$. To the lowest order for $\epsilon_{n}^{\prime}$, the action is just Schwarzian derivative while for $\epsilon_{n}$ there is an additional term $\propto|n|(n^{2}-1)$.

The $h=2$ contribution to the leading order of the connected part four-point function of the SYKχ system can be read from this effective action as :

[TABLE]

where $\tau=(\tau_{1}+\tau_{2}-\tau_{3}-\tau_{4})/2$. By analytical continuation with $\tau_{12}=\tau_{34}=\beta/2$ and $\tau=it$, we obtain the growing term in the long time as

[TABLE]

We see that the Lyapunov exponent is still $2\pi/\beta$, because we have not taken the $h\neq 2$ corrections into account yet. Here we have two small parameters $1/\beta J$ and $u^{2}/J^{2}$, which are both of order $\epsilon$. Suppose $F_{\chi\chi}$ can be written in the form as

[TABLE]

It is analyzed in Ref.Comments that the $h\neq 2$ part gives a growing term $\frac{6\pi}{\beta}te^{\frac{2\pi}{\beta}t}$. By matching the above expansion, we can obtain

[TABLE]

where the first term is the first order correction to the Lyapunov exponent away from the conformal limit, while the second term is due to the interaction with the thermal bath. In the conformal limit $\beta J\rightarrow\infty$, the Lyapunov exponent is given by

[TABLE]

which is consistent with the result Eq. (30) and Eq. (32) in the sec. 2.2.

3 (1+1)-d SYK Models Coupled to a Thermal Bath

In this section we consider the (1+1)-dimensional generalization of the previous model. By generalizing to a one dimensional chain model, we can study richer physics including the spatial dependence of OTOC, the energy transport property and so on. There are several proposals for generalizing the SYK model to higher dimensions Yingfei1 ; Yingfei2 ; high-D1 ; high-D2-con , here we mainly follow Ref. Yingfei1 ; Yingfei2 and focus on the generalization to one dimensional chain model. Two different configurations of chain model will be discussed in this section. In both cases, we first prepare an one dimensional SYK chain, with $N$ Majorana fermions $\chi_{i,x}$ on each site, and prepare a SYKψ system with $N^{2}$ Majorana fermions $\psi_{i}$. In the first configuration, we use the SYKψ system as a global bath, and couple every site of the SYKχ chain to the SYKψ system uniformly. In the second configuration, we use the SYKψ system as a local bath and attach the SYKψ system to one end of the SYKχ chain. The first model preserves the translational symmetry along the chain, while the second model does not.

3.1 (1+1)-d SYK Chain Coupled to a Global Thermal Bath

In this subsection, we study the first configuration as illustrated in Fig. 5. To maintain the translational symmetry, here we choose the interaction strength with the $N^{2}$ fermions SYKψ bath to be uniform. If the interaction strength is chosen different, this will be a new type of inhomogeneous SYK chain model, which can be analyzed using similar method as in Ref. Yingfei2 .

The Hamiltonian of the system is

[TABLE]

where $H^{c}$ is defined as in Eq. (13). The random couplings in each term are $\{J_{jklm,x}\}$, $\{J^{{}^{\prime}}_{jklm,x}\}$, $\{u_{ijkl}\}$, $\{\tilde{J}_{ijkl}\}$ in order, with

[TABLE]

As we do not want the bath SYKψ system to be affected by the existence of the chain at leading order of $\frac{1}{N}$, we demand the number of sites $M$ satisfy $\frac{M}{N}\rightarrow 0$ in the large-$N$ limit.

Using the replica method, after doing the disorder average and introducing the bi-local fields, one arrives at the effective action as (in below we omit the terms for the SYKψ system):

[TABLE]

This effective action shows that the translational invariant saddle point solutions should satisfy:

[TABLE]

[TABLE]

where $a$ and $b$ are the factors in front of the Green’s functions $G^{s}(\tau)$ and $G_{\psi}(\tau)$ as in Eq. (18). The value of $b^{2}/a^{2}$ can be tuned by changing the variance $\tilde{J}^{2}$ of the random couplings $\{\tilde{J}_{jklm}\}$. Here we set $b^{2}/a^{2}=1$ for the simplification of the equations by changing the value of $\tilde{J}^{2}$. We define

[TABLE]

then the Eq. (57) becomes:

[TABLE]

In the conformal limit $N\gg\beta J\gg 1$, one has the saddle point solution

[TABLE]

Expanding the fields around the saddle point solutions as

[TABLE]

[TABLE]

one can expand the effective action to quadratic order as

[TABLE]

The spatial kernel $S_{xy}$ is defined as

[TABLE]

with its Fourier transformation into momentum space being

[TABLE]

Defining $\tilde{K}$ as the symmetrized four-point function kernel of the SYK model

[TABLE]

and by integrating out the $\sigma_{x}$ fields, the effective action can be written as

[TABLE]

One can read the four point function with two fermions on the $x$ site and two fermions on the $y$ site from this action, which is

[TABLE]

By Fourier transforming $(x-y)$ to the momentum space, one obtains the four point function in momentum space as

[TABLE]

One can further get the behavior of the four point function from this expression by diagonalizing the kernel $\tilde{K}$ in the same way as in the original SYK model. The kernel has eigenvalues $k(h,n)$ and eigenfunctions $\Psi_{h,n}(\tau_{1},\tau_{2})$ which have been worked out in the Ref. Comments . With these eigenvalues and eigenfunctions, an exact expression for the four-point function in the momentum space is given by

[TABLE]

The would-be divergence from $h=2$ is removed by $s(p)<1$, similar as in the section 2.3. Similarly, when the divergence has been removed, the long time behavior of the four point function comes from both the $h=2$ and $h\neq 2$ parts. Following the same procedure as in sec. 2.3, because for $p^{2}\ll 1$ and $u^{2}/J^{2}\ll 1$, $s(p)\approx 1$, the four-point function in the long time and energy transport properties in the long wavelength limit are still mainly governed by the $h=2$ reparametrization modes. It is reasonable for us to first focus on the $h=2$ part, and then take into account of the corrections coming from $h\neq 2$ parts perturbatively when necessary. For the $h=2$ modes, we make reparametrizations to the saddle point solution

[TABLE]

with small deformations $f_{x}(\tau)=\tau+\epsilon_{x}(\tau)$. Using the same method as in sec. 2.3, we derive the effective action of the chain model in terms of the reparametrization modes:

[TABLE]

where $\epsilon_{n,p}$’s are the reparametrization modes by Fourier transforming $\epsilon_{x}(\tau)$:

[TABLE]

Thus the $h=2$ contribution to the four-point function can be derived as ($\tau=(\tau_{1}+\tau_{2}-\tau_{3}-\tau_{4})/2$):

[TABLE]

in which

[TABLE]

$D$ is the diffusion constant of the system. To see this, one needs to consider the energy transport property of the system. The diffusive dynamics in space-time is governed by the reparameterization fields, which are the same as the ones that govern the long time behavior of fermion four-point function. The retarded energy-density correlation function in the $(p,\omega)$ space can be derived by utilizing the result of fermion four-point function in the OPE limit $|\tau_{2}-\tau_{1}|,|\tau_{4}-\tau_{3}|\ll\tau$. By the same method as in Ref. Comments Yingfei1 , we get the expression for the retarded energy-density correlation function $C^{R}_{00}(p,\omega)$:

[TABLE]

The denominator $(-i\omega+\frac{1}{\tau_{u}}+Dp^{2})$ comes from a diffusion-type equation in real space

[TABLE]

where $\phi(x,t)$ is the energy density function. From this equation, one finds the physical meaning of $D$ and $\tau_{u}$, where $D$ is the diffusion constant, while the $\frac{1}{\tau_{u}}\phi$ term will lead to an exponential decay of the energy density function, with $\tau_{u}\sim\frac{J}{u^{2}}$ being the characteristic time scale for the decay. Physically, this characterizes the process of energy leaking into the bath SYKψ system through interaction, which is similar to the process that the energy flow crossing the black hole horizon being absorbed into the black hole, and can no longer be extracted by an exterior observer.

After the discussion of the energy transport property, we now turn to discuss the behavior of the OTOC in the chaos limit. The OTOC as a function of both space and time contains two majorana fermions from position $x$ and two from position [math]:

[TABLE]

from which we can extract both the Lyapunov exponent and the butterfly velocity. At the leading order $\mathcal{O}(N^{0})$, the four-point function is given by a disconnected part $-G(\frac{\beta}{2})^{2}$, and the $\mathcal{O}(\frac{1}{N})$ part of the function can be derived by first doing analytical continuation of the Eq. (74) with $\tau_{12}=\tau_{34}=\beta/2$ and $\tau=it$, and then Fourier transforming back to real space. By analytical continuation, we obtain the growing term in the long time as

[TABLE]

We see that it is exponential growing in the long time. In the same spirit as in sec. 2.3, we need to consider the $h\neq 2$ contributions to get the modified Lyapunov exponent. By same method, after taking into account the first order correction of $p^{2},u^{2}/J^{2}$, as well as $1/(\beta J)$, one finds

[TABLE]

where

[TABLE]

$F(x,t)$ is given by Fourier transforming Eq. (80) back into real space

[TABLE]

The residue at the pole $p^{*}=i\sqrt{\frac{\frac{2\pi}{\beta}+\frac{1}{\tau_{u}}}{D}}$ dominates the integral for large $x$. One obtains

[TABLE]

We find that for large $x$, the first order correction to the Lyapunov exponent exactly cancels each other, leaving a Lyapunov exponent $2\pi/\beta$ unchanged as in Ref. Yingfei1 . However, this will not hold for small $x$ and the Lyapunov exponent will be tuned. The butterfly velocity is given by

[TABLE]

Defining $\tau_{L}=\frac{\beta}{2\pi}$, we have

[TABLE]

This relationship is very interesting. On one hand, it still satisfies

[TABLE]

which is proposed to be a bound for the butterfly velocity vb bound ; vb energy1 ; vb energy2 ; thermal transport ; holographic review ; Sachdev new ; Blake . It is found to held in a SYK chain Yingfei1 and holographic theories vbbh .

On the other hand, the scaling law of $v_{B}^{2}/D$ shows a crossover behavior as temperature increases, as shown in Fig. 6. As we mentioned, $\tau_{u}$ is the time scale which governs the decay of the energy density function, while it can also define a temperature scale that is proportional to $u^{2}/J$. When the temperature of the system is much higher than this temperature scale, i.e., $k_{B}T\gg u^{2}/J$, Eq. (85) tells

[TABLE]

which is the same temperature scaling law as in Ref. Yingfei1 ; Yingfei2 . However, when the temperature of the system is far below that temperature scale, i.e., $k_{B}T\ll u^{2}/J$, Eq. (85) tells

[TABLE]

This scaling law is similar as that in the Landau’s Fermi liquid theory, in which $v_{B}$ is replaced by Fermi velocity $v_{F}$, which is a constant in low temperature, while $D$ scales as $\frac{1}{T^{2}}$. However, the difference is that here $D$ is a constant that does not depend on the temperature, while $v_{B}$ scales as $T$.

To conclude this subsection, we find that, by coupling the one dimensional SYK chain to a larger SYK bath, the butterfly velocity can not only acquire different values, but also acquire different temperature dependence. From Eq. (84), it can be seen that when we increase the interaction strength $u$ with the bath, the butterfly velocity decreases, which physically means that the information leaks into the bath. Meanwhile, the dependence of the butterfly velocity on temperature transits from $v_{B}\propto\sqrt{T}$ to $v_{B}\propto T$.

3.2 (1+1)-d SYK Chain Coupled to a Local Thermal Bath

In this subsection, we study another configuration, in which we couple the SYKψ system to one end of the SYK chain as shown in Fig. 7. Different from the model in sec. 3.1, this model is no longer translational invariant. In this configuration, our focus will be the spatial dependence of both the Lyapunov exponent and the butterfly velocity. We call the site which is nearest to SYKψ as the site $1$, the second nearest as the site $2$, and so on.

The Hamiltonian of the system can be written as

[TABLE]

The random couplings in each term are $\{J_{ijkl,x}\},\{J^{{}^{\prime}}_{ijkl,x}\},\{J^{{}^{\prime}}_{ijkl,0}\},\{\tilde{J}_{ijkl}\}$, with

[TABLE]

Note that in Eq. (90), the meanings of the symbols $J$ and $u$ are different from their meanings in the last subsection. For the physical picture to be clear, we only kept two parameters $u$ and $J$, although we could have four different parameters in the Eq. (90). The Schwinger-Dyson equations are

[TABLE]

with

[TABLE]

Suppose the solutions in the conformal limit still obey the ansatz

[TABLE]

then one can get a set of coupled equations for the coefficients $b$ and $\{a_{x}\}$ as:

[TABLE]

With the asymptotic boundary condition $a_{x}\rightarrow\textrm{const}$ when $x\rightarrow\infty$, one can solve the coupled equations numerically to get the coefficients $b$ and $\{a_{x}\}$.

Our goal is to calculate the OTOC functions in the chaos limit, and extract the Lyapunov exponents and butterfly velocities from them. In this model, the OTOCs $F_{1}(t_{1},t_{2})$, $F_{2}(t_{1},t_{2}),...$ 222$F_{i}(t_{1},t_{2})$ is the $\frac{1}{N}$ piece of the OTOC with four fermions from site $i$, as defined before. are coupled to each other. For $F_{1}(t_{1},t_{2})$, the self-consistency equation involves $F_{2}(t_{1},t_{2})$:

[TABLE]

where

[TABLE]

and

[TABLE]

For $F_{x}(t_{1},t_{2})$ with $x\geq 2$, the self-consistency equation involve $F_{x-1}(t_{1},t_{2})$ and $F_{x+1}(t_{1},t_{2})$:

[TABLE]

with

[TABLE]

and

[TABLE]

These equations can be casted into a matrix form with the kernel matrix $\mathcal{K}=\{K_{ij}\}$, and by solving the matrix integral equation, one would get different modes with different Lyapunov exponents, similar as different momentum modes in sec. 3.1. The matrix integral equation can be solved by the ansatz solution as

[TABLE]

From this one can solve a series of eigen-modes, each with different $h$, which corresponds to different Lyapunov exponents. A technical detail is that the kernel matrix $\mathcal{K}=\{K_{ij}\}$ is non-Hermitian, thus the eigen-modes $\{\vec{\mathcal{F}_{h}}\}$ are not orthogonal with each other. To keep the eigen-modes orthogonal with each other, one needs to use a symmetrized kernel matrix $\tilde{\mathcal{K}}=\{\tilde{K}_{ij}\}$. After getting a series of different modes, the $\frac{1}{N}$ part of the OTOC $\langle\chi_{i,x}(t)\chi_{j,x}(0)\chi_{i,x}(t)\chi_{j,x}(0)\rangle_{\beta}$ with four fermions on the $x$ site can be written as a summation of different modes:

[TABLE]

where $\alpha=\frac{4\pi J}{3\sqrt{2}\alpha_{K}}$. The coefficients $\{\tilde{\mathcal{F}}_{hx}\}$ are the mode expansion coefficients derived in numerics by using the symmetrized kernel matrix. The weight of different modes in summation should be determined by the effective action, while one can also compare with the results in the last subsection, and see it should be approximated by $\frac{1}{\frac{2\pi}{\beta}+\alpha(h+1)}$. For the OTOC $\langle\chi_{i,x}(t)\chi_{j,y}(0)\chi_{i,x}(t)\chi_{j,y}(0)\rangle_{\beta}$ with two fermions on the $x$ site and two fermions on the $y$ site, it can be derived as

[TABLE]

We first consider the OTOC with four fermions from the same site. We give the plot of the behavior of $F_{x}(t)$ for the first five sites away from the SYKψ system in the Fig. 8 (a). In numerics, we set $\beta=2\pi$, and $J=u=100$. We can see that, the closer to the SYKψ system, the slower the growth of $F_{x}(t)$. To quantify the growth rate of $F_{x}(t)$, we fit $F_{x}(t)$ after the dissipation time $t_{d}\sim\frac{\beta}{2\pi}$ with function $e^{\lambda_{L}t}$. We get the local Lyapunov exponent $\lambda_{L}$ for each site, which characterizes the growth rate of the OTOC. The result is shown in Fig. 8 (b). We see that the Lyapunov exponents acquire a non-trivial spatial dependence. In our (0+1)-d model, the interaction with a much larger cluster of SYK system results in the decrease of Lyapunov exponent of the small system. In this (1+1)-d model, this influence depends on the distance between the large system and each site on the chain.

To look at the spatial dependence of the butterfly velocity, we will study the OTOC $F_{x_{c}-a,x_{c}+a}(t)$ with different $a=0,1,2,...$ . We fit it with the form:

[TABLE]

where we assume that $\lambda_{L}$ and $v_{B}$ only depend on the center position $x_{c}$. In numerics, we can check that this assumption works well for small $a$. The dependence of $v_{B}$ on $x_{c}$ is shown in Fig. 9. From the plot we can see that, the closer the center position $x_{c}$ is to the SYKψ system, the smaller the butterfly velocity. This again shows that the presence of the larger system can also slow down the propagation of quantum chaos through space. The spatial dependence of Lyapunov exponent and butterfly velocity also backward prove the usefulness of the concepts in characterizing the degree of quantum chaos even locally.

4 Summary and Outlook

In this paper, we proposed several generalizations of the SYK models in which the characteristics of quantum chaos can be tuned by varying the coupling strength with a thermal bath. Owing to the separation of the scrambling time scales between the small and the large SYK system, and utilizing the nice large-$N$ structure of the SYK model, our set-up provides an exactly solvable model to illustrate the physics of how to tune the Lyapunov exponent and the butterfly velocity. We find a new temperature dependence law for the butterfly velocity as $v_{B}\propto T$ below the temperature scale induced by coupling to the bath. We also use simple numerics to illustrate the spatial dependence of the Lyapunov exponent and the butterfly velocity, which emphasizes the usefulness of the quantities as tools to diagnose the extent of quantum chaos locally.

Since the Lyapunov exponent and the butterfly velocity are both tunable in our models, our work provides a broader platform to test various inequalities proposed for chaotic systems and to study the thermalization of chaotic systems. One possible generalization of our model is to impose an $U(1)$ charge symmetry, then one could study the thermoelectric transport properties of the system thermal transport and see how they should be modified by interaction with the bath.

Another interesting direction is to search for a holographic description of the system. Since the Schwarzian action of the original SYK model is shared by NAdS2 gravity system, it is natural to ask whether the modified action given in this paper can be also found in some gravitational systems. More generally, it is still an open question whether the Lyapunov exponent in gravitational systems can be tuned by similar physical process inspired by this work. Similar questions can also be addressed on the behavior of butterfly velocity $v_{B}$. We defer these questions for further study.

Acknowledgements.

We thank Wenbo Fu, Yingfei Gu, Chao-Ming Jian, Yi-Zhuang You and Xiao-Liang Qi for helpful discussion. This work is supported by NSFC Grant No. 11325418 and MOST under Grant No. 2016YFA0301600.