Out-of-time ordered correlation functions for the localized $f$ electrons in the Falicov-Kimball model

TL;DR

This paper calculates exact out-of-time correlation functions for localized $f$ electrons in the Falicov-Kimball model, revealing insights into quantum chaos and scrambling within a mean-field framework.

Contribution

It provides an exact evaluation of OTOCs for localized electrons in the Falicov-Kimball model using dynamical mean-field theory, highlighting regimes of quantum chaos.

Findings

Different regimes of quantum chaos and scrambling identified.

OTOCs show similarities with disordered quantum spin models.

Winding numbers of block Toeplitz matrices characterize chaos regimes.

Abstract

We provide an exact evaluation of the out-of-time correlation (OTOC) functions for the localized $f$-particle states in the Falicov-Kimball model within dynamical mean-field theory. Different regimes of quantum chaos and quantum scrambling are distinguished by the winding numbers of the block Toeplitz matrices used in the calculation. The similarities of these fermionic OTOCs and their logarithmic derivatives for time evolution with the OTOCs for quantum spin models with disorder are also discussed.

Equations

$$F(r,t)=⟨[W_0(t),V_r(0){]}^{†}[W_0(t),V_r(0)]⟩,$$

$$L(r,t)=\frac{\mathrm{d}\ln F(r,t)}{\mathrm{d}t}=\frac{1}{F(r,t)}\frac{\mathrm{d}F(r,t)}{\mathrm{d}t},$$

$$H=-\sum _{ij}{t}_{ij}{d}_i^{†}{d}_j+\sum _i{H}_i$$

$$H_i=U{n}_{id}{w}_{i1}-\mu n_{id}+(E_f-\mu )w_{i1}$$

$$w_{i1}=n_{if}=f_i^{†}{f}_i$$

$$G_d(\omega )$$

$$={\int }_{-\infty }^{+\infty }\mathrm{d}\epsilon \frac{\rho (\epsilon )}{\omega +\mu -\Sigma (\omega )-\epsilon },$$

$$G_d(\omega )$$

$$=w_0{\mathcal{G}}_0(\omega )+w_1{\mathcal{G}}_1(\omega ).$$

$$\rho (\omega )=\frac{1}{t^{\ast }\sqrt{\pi }}{\mathrm{e}}^{-(\omega /t^{\ast }{)}^2}.$$

$$G_{d,\mathbf{k}}(\omega )=\frac{1}{\omega +\mu -\Sigma (\omega )-{\epsilon }_{\mathbf{k}}},$$

$$f_i(t)$$

$$f_i^{†}(t)$$

$${\mathcal{U}}_i(t,0)$$

$${\tilde{\mathcal{U}}}_i(t,0)$$

$$G_{ijf}^c(t,t^{\prime })=-\mathrm{i}⟨{\mathcal{T}}_c{f}_i(t)f_j^{†}(t^{\prime })⟩={\delta }_{ij}{G}_f^c(t,t^{\prime }),$$

$$G_f^{>}(t)$$

$$=-\mathrm{i}{w}_0{\mathrm{e}}^{-\mathrm{i}(E_f-\mu )t}$$

$$\times {\det }_{[0,t]}\Vert \delta (t_1-t_2)-{\mathcal{G}}_0(t_1-t_2)U\Vert$$

$$G_f^{<}(-t)$$

$$\times {\det }_{[0,t]}\Vert \delta (t_1-t_2)+{\mathcal{G}}_1(t_1-t_2)U\Vert .$$

$${\mathcal{G}}_0(t^{\prime }-t^{\prime \prime })$$

$$\times \mathrm{Im}\frac{1}{\omega +\mu -\lambda (\omega )}$$

$${\mathcal{G}}_1(t^{\prime }-t^{\prime \prime })$$

$$\times \mathrm{Im}\frac{1}{\omega +\mu -U-\lambda (\omega )}$$

$$G_f^r(t)$$

$$=-\mathrm{i}\Theta (t){\mathrm{e}}^{-\mathrm{i}(E_f-\mu )t}\{w_0{\det }_{[0,t]}\Vert \mathbf{I}-{\mathcal{G}}_{0}U\Vert$$

$$+w_1{\left({\det }_{[0,t]}\Vert \mathbf{I}+{\mathcal{G}}_{1}U\Vert \right)}^{\ast }\}$$

$$=G_0^f(t)+G_1^f(t)$$

$$G_f^r(\omega )$$

$$+w_1{\left({\det }_{[0,t]}\Vert \mathbf{I}+{\mathcal{G}}_{1}U\Vert \right)}^{\ast }\},$$

$$A_f(\omega )=-\frac{1}{\pi }\mathrm{Im}{G}_f^r(\omega ).$$

$$C_{0,1}^f(\omega )=1\mp U{\mathcal{G}}_{0,1}^t(\omega )$$

$$=1\mp U\mathrm{Re}{\mathcal{G}}_{0,1}^r(\omega )\mp \mathrm{i}U\tanh \frac{\beta \omega }{2}\mathrm{Im}{\mathcal{G}}_{0,1}^r(\omega ),$$

$${\check{C}}_{\alpha }^f(\omega )=-{\mathrm{e}}^{\mathrm{i}\omega W}{C}_{\alpha }^f(\omega ),$$

$$H=\frac{1}{2}\sum _{ij}{V}_{ij}{n}_i{n}_j-\mu \sum _i{n}_i,$$

$$\frac{\mathrm{d}{c}_i(t)}{\mathrm{d}t}$$

Full text

I Introduction

Long-time coherent quantum many-body dynamics is important for condensed matter physics, quantum information science, optical physics, and high-energy physics [1, 2]. First, it is connected with the spreading of quantum information, its detection, and decoding. These dynamics can be examined in a system of qubits whose Hamiltonian governs unitary dynamics where the initial state of one qubit becomes nonlocal and is distributed across a system by the unitary dynamics—the so-called quantum information scrambling phenomenon. The rate with which information becomes nonlocal is usually quantified by the out-of-time-ordered correlators [3, 4]

[TABLE]

where $W_{0}$ and $V_{r}$ are local operators at positions 0 and $r$.

Out-of-time-ordered correlators (OTOC) were proposed by Larkin and Ovchinnikov [5, *larkin_1969:1200] as a measure of importance of quantum motion with respect to the quasiclassical approximation in superconductivity. Later, OTOCs were used to describe the evolution of the overlap of perturbed and unperturbed states as indicators of the quantum motion stability in chaotic and regular systems [7].

For chaotic systems, the initially local operator $W_{0}(t)$ expands, under Heisenberg dynamics, with a speed called the butterfly velocity and overlaps later with the operator $V_{r}$ resulting from the growth of the OTOC function inside the emergent light-cone in space-time. The growth rate of the OTOC is influenced by the interplay between quantum fluctuations and local scrambling, and depends on the size of the system and on its temperature. For thermal quantum systems with a large number of degrees of freedom, the initial classical exponential growth rate is characterized by the Lyapunov exponent bounded by the temperature $\lambda_{L}\leq 2\pi k_{B}T/\hbar$ —this bound is a bound on the growth rate of chaos [8, 9]. However, for quantum dissipative systems, e.g. systems in contact with a thermal environment, we find that instead of the short time exponential growth, they exhibit a long time decay with a decay rate that mirrors the classical Lyapunov exponent [10, 11].

In order to estimate the growth rate of an OTOC, in particular, for exponential growth, it is convenient to consider the logarithmic derivative of the OTOC [12, 13, 14, 15]

[TABLE]

which is equal to the Lyapunov exponent for “fast scramblers” obeying the chaos bound, and can provide detailed behavior of the OTOC at its wavefront, even when it is not exponentially growing. Relaxation dynamics in many-body systems can be characterized by different time scales $\tau$, each with their own temperature and wave-vector dependencies, which affect the scrambling, and result in an anomalous behavior for the Lyapunov exponents [16].

For quantum information scrambling, one is interested in the OTOCs of the local two-state operators—that is, for the qubits. As a rule, spin models with a direct interaction between the spins are used to describe their unitary dynamics, see e.g. [1] and references therein. On the other hand, in real systems, the qubits are strongly localized and their dynamics, in many cases, are governed by the indirect interaction through other degrees of freedom, which often become retarded and long range.

One of the simplest models with two-level states interacting through (nonlocal) band electrons is the Falicov-Kimball model (FKM) [17], which considers localized $f$-states which locally interact with itinerant $d$-electrons. It was shown by Brandt and Schmidt [18] and Kennedy and Lieb [19] that the ground state phase diagram of the FKM can be obtained from the Ising-like (lattice gas) model with an effective interaction between classical spins (occupation of $f$ states) after integrating out the fast $d$-electrons, and the cases of 1D and 2D lattices have been examined in detail [20, 21, 22, 23].

In addition, the FKM possesses exact solutions for infinite dimensional lattices [24, 25] within the dynamical mean-field theory (DMFT) [26, 27, 28]. Despite the $f$-electron occupation number being a conserved quantity (an integral of motion), the $f$-electron spectral function is a complicated object which displays nontrivial features [29, 30, 31]. The local $d$-electron density-density OTOC for the FKM was calculated and its enhancement near the metal-insulator transition point was obtained [32]. But, since the $f$-electron density-density OTOC for the FKM is exactly equal to zero, one must look at OTOCs constructed from $f$-electron creation and annihilation operators, which describe transitions between local $f$-states, similar to the spin-flip operators for spin systems.

In this article, we investigate OTOC functions constructed from $f$-electron creation and annihilation operators (state-flip operators for the two-level system) for the FKM on an infinite dimension hypercubic lattice. First, we introduce the model and recall the derivation of the $f$-electron DOS. Next, we consider the lattice gas model and obtain analytic results for the single-particle DOS and OTOCs for any lattice, and then examine the infinite-dimensional limit. Finally, we calculate OTOCs for $f$-electrons for the FKM and provide numerical results.

II Model Hamiltonian and single-particle propagators

The Falicov-Kimball model (FKM) [17] is the simplest model for strongly correlated electron systems and describes a system of the itinerant $d$-electrons which interact with localized particles. Its Hamiltonian

[TABLE]

contains a hopping contribution from the itinerant electrons (kinetic energy) and a local part

[TABLE]

which includes the interaction of strength $U$ between the itinerant and localized particles, Here, $n_{id}=d_{i}^{\dagger}d_{i}$ is the local density of $d$ electrons and $w_{i1}$ is the occupation of site $i$ by the localized particles, which can take two values, 0 or 1, and describes a two-level state subsystem without specifying its internal dynamics. As a rule, it is associated with the localized $f$-electron state, i.e., it is the expectation value of the localized electron density operator

[TABLE]

and its internal dynamics is determined by the $f$-electron creation $f_{i}^{\dagger}$ and annihilation $f_{i}$ operators. These local $f$-electron densities are conserved in the model, as they commute with the Hamiltonian.

The FKM possesses an exact solution in the limit of infinite dimensions $D\to\infty$, after scaling the nearest-neighbor hopping by $t_{ij}=t^{*}/2\sqrt{D}$, and can be solved exactly with DMFT [28, 25]. Below, we shall put $t^{*}=1$ as the energy scale.

In DMFT, the self-energy is local [26, 27], and the lattice model is mapped onto an effective impurity problem with the same self-energy and some auxiliary self-consistent $\lambda$-field [24]. Both the self-energy $\Sigma(\omega)$ and $\lambda$-field $\lambda(\omega)$ are found from solutions of the self-consistent set of equations

[TABLE]

where $G_{d}(\omega)$ is the $d$-electron Green’s function for the impurity problem and for the FKM it becomes

[TABLE]

Here, $w_{1}=n_{f}=\langle f^{\dagger}f\rangle$ and $w_{0}=1-n_{f}=\langle ff^{\dagger}\rangle$ are average occupations of the localized $f$-states on the lattice and $\rho(\omega)$ is the unperturbed density of states, which for the infinite-dimensional hypercubic lattice is a Gaussian

[TABLE]

Finally, for the $d$-electron propagator we find

[TABLE]

where $\varepsilon_{\mathbf{k}}=-\frac{t^{*}}{\sqrt{D}}\sum_{\alpha}\cos k_{\alpha}a$ is the unperturbed band energy.

Finding the solution for the $f$-electron propagator is more challenging [29, 33, 25]. We follow the derivation in [31]. The time evolution of the $f$-operators is given by

[TABLE]

where the evolution operators are

[TABLE]

Based on this, the contour-ordered single-particle propagator is defined to be

[TABLE]

where the time ordering is performed on the Keldysh contour and the Coulomb interaction $U$ acts only on the upper (forward) branch of contour. It is a local quantity because there is no hopping for the $f$ electrons. In particular, the greater $G_{f}^{>}(t)$ and lesser $G_{f}^{<}(t)$ Green’s functions can be expressed in terms of functional determinants of Toeplitz matrices calculated on the time interval $[0,t]$ with $t>0$ (for details, see [31])

[TABLE]

and [$G_{f}^{<}(t)=-\{G_{f}^{<}(-t)\}^{*}$]

[TABLE]

Here

[TABLE]

and

[TABLE]

are the time ordered (causal) impurity Green’s functions for the empty and filled $f$ states, respectively, using the boundary condition for the step function, given by $\Theta(0)=0$ [34]. The retarded Green’s functions is defined in the usual way via

[TABLE]

with a Fourier transform

[TABLE]

which defines the single-particle DOS as

[TABLE]

The long time asymptotic of the Toeplitz determinants in (19) can be obtained by using Szegő’s theorem and the Wiener-Hopf approach, see [31] and references therein. It is determined by the winding number (index) of the characteristic functions $C_{0}^{f}(\omega)$ and $C_{1}^{f}(\omega)$, which are given by

[TABLE]

respectively. Here, $\mathcal{G}_{0,1}^{t}(\omega)$ is the Fourier transform of the time-ordered impurity Green’s function in (17) and (18). We choose both time indices $t^{\prime}$ and $t^{\prime\prime}$ to be on the top branch of the Keldysh contour, and $\mathcal{G}_{\alpha}^{r}(\omega)$ is the retarded Green’s function. According to Szegő’s theorem, the long time asymptotic is determined by the logarithm of this function, which is not an analytic function. The imaginary part of $C_{\alpha}^{f}(\omega)$ crosses the zero line only at $\omega=0$ and the winding number is determined by the sign of the real part at this point. Note that we can not calculate the winding number via its strict mathematical definition; instead we use Nyquist (Cole-Cole) type plots that illustrate their topology along with the changes of the sign of the imaginary part for negative values of the real one.

At half filling and for $U<0.866$, when both $1-U\operatorname{Re}\mathcal{G}_{0}^{t}(\omega=0)>0$ and $1+U\operatorname{Re}\mathcal{G}_{1}^{t}(\omega=0)>0$, the winding number is equal to zero (see Appendix A). Its logarithm is well defined and the corresponding determinants decay exponentially at large times for high temperatures, while they behave as a power-law for low temperatures (with a singular power law at $T=0$ only). This behavior can be seen in the change of the asymptotics of the logarithmic derivative from a constant to a $c/t$ dependence (Fig. 1). The corresponding DOS displays a Lorentzian peak for high temperatures, which is replaced by a power law singularity in the limit as $T\to 0$ in agreement with the seminal results of Brandt and Urbanek [29].

On the other hand, for $U>0.866$, when both $1\mp U\operatorname{Re}\mathcal{G}_{0,1}^{t}(\omega=0)<0$, the logarithm of the characteristic function is multivalued with a branch cut, and the winding number is equal to $-1$ (see Appendix A). In this case, the characteristic function (22) is replaced by

[TABLE]

where $W=\Delta t\to 0$ is the time discretization, and the long time asymptotic becomes a more complicated nonexponential one with oscillations of the logarithmic derivatives (Fig. 2). Such oscillations are strongest at high temperatures and decrease their amplitude and increase their period as the temperature is lowered. Now, the corresponding DOS possesses a gap at low temperatures due to freezing of thermal fluctuations [29] and the presence of strong-correlations creating a Mott gap; the gap region is filled in by thermal excitations for increasing temperatures.

Szegő’s theorem and the Wiener-Hopf approach work better for smooth functions. In our case, the functions (17) and (18) are smooth at high temperatures, whereas at low temperatures, the Fermi distribution $f(\omega)$ becomes step-like and one must employ results for singular integral equations [35, 36]. Thus, at $T=0$, the long time asymptotic becomes a power law in the metallic phase for $U<U_{c}=\sqrt{2}$ and rapidly goes to zero for the Mott insulator ($U>U_{c}$).

III Lattice gas

As mentioned above, the FKM can be projected onto a classical lattice gas (Ising) model with effective interactions. So, we first consider the simple lattice gas model

[TABLE]

where $n_{i}=c_{i}^{\dagger}c_{i}$ is the occupation of site $i$ by a spinless Fermi-type particle, with values $n_{i}=1$ or [math]. This model corresponds to a classical Ising model with $S_{i}^{z}=n_{i}-\frac{1}{2}=\pm\frac{1}{2}$.

Despite the occupation operator $n_{i}$ being a conserved quantity, $[n_{i},H]=0$, one can still examine nontrivial behavior for the dynamics and for the density of states for such fermions. From the equations of motion, we have that

[TABLE]

This then yields the time evolution of the annihilation and creation operators, respectively,

[TABLE]

where the evolution operators are given by

[TABLE]

Now, for the single-particle Green’s function we find [$G_{ij}(t,t^{\prime})=\delta_{ij}G_{i}(t,t^{\prime})$]

[TABLE]

where we have used the fact that $n_{i}(t)=n_{i}$ is a conserved quantity, which allows us to replace the time-ordered exponents by ordinary ones.

In the case of a nearest-neighbor interaction, one has

[TABLE]

where $Z$ is the coordination number (number of nearest neighbors for site $i$). Calculating the local density of states, we obtain

[TABLE]

Here $w_{i}(s)$ are statistical weights of the configurations with $s$ occupied neighbors of site $i$, and

[TABLE]

Then, for the homogeneous phase, we have

[TABLE]

Now, we examine the OTOCs for the lattice gas model. For the square modulus of the anticommutator, we have

[TABLE]

Using $c_{i}^{\dagger}(t)c_{i}(t)=n_{i}$, $c_{j}(0)c_{j}^{\dagger}(0)=1-n_{j}$, and

[TABLE]

we find both local and nonlocal terms

[TABLE]

Similarly, another OTOC gives us

[TABLE]

which are also nonlocal. One can see, that these OTOCs expand only with the distance of the interaction range $V_{ij}$ and represent a direct instantaneous interaction. In the limit $Z\to\infty$ for the nearest-neighbor interactions (34), the nonlocal contributions vanish as $1/Z$.

On the other hand, for the square of the anticommutator, we have only local contributions given by

[TABLE]

while for nearest-neighbor interactions (34), we find

[TABLE]

One can see, that the expression (43) for the OTOC constructed by the square of the anticommutators is similar to the one for the retarded Green’s function (33), both having the square of the exponent (evolution operator). It can potentially be considered as a higher moment of a (quasi)probability.

IV OTOCs for the Falicov-Kimball model

Similar to the lattice gas model, we can examine different OTOCs for the $f$ electrons of the FKM.

IV.1 Anticommutator squared magnitude

The first one is the squared magnitude of the anticommutator of the creation and annihilation operators at different times $[f_{i}(t),f_{j}^{\dagger}(0)]_{+}$, which satisfy

[TABLE]

For the FKM, the product of the creation and annihilation operators at the same times in the first two terms give either $P_{i}^{+}=f_{i}^{\dagger}(t)f_{i}(t)$ or $P_{i}^{-}=f_{i}(t)f_{i}^{\dagger}(t)$, which are conserved quantities, and do not depend on time. These terms yield static correlation functions. But, the last two terms of the time evolution of the operators are determined by (10) and (11) and, as a result, we find

[TABLE]

where in the last term we have an evolution on the doubly folded Keldysh contour (Fig. 3).

For the local function, $i=j$, we have

[TABLE]

The nonlocal function $(i\neq j)$ satisfies

[TABLE]

In the infinite-dimensional case, it contains two evolution operators (12) which act on the site $i$ and do not touch the site $j$. As a result, we can carry one of the evolution operators through the $f_{j}^{\dagger}$ and they mutually cancel. Finally, we obtain

[TABLE]

On the other hand, in the finite $D\to\infty$ case, one has to replace the evolution operators (12) and (13) for the local impurity problem by the evolution operators for the lattice, which makes it much more difficult to find a solution and we do not consider it further here. This does give nonlocal $1/Z$-corrections, but, in contrast to (40), they are governed by the effective long-range retarded interactions.

IV.2 Other trivial anticommutator functions

We can also consider an anticommutator of two creation $[f_{i}^{\dagger}(t),f_{j}^{\dagger}(0)]_{+}$ or two annihilation $[f_{i}(t),f_{j}(0)]_{+}$ operators and construct another OTOC. The square of such anticommutators yield OTOC functions which are equal to zero. We have

[TABLE]

For the squared magnitude, similar to (48), we find that

[TABLE]

V Square of the anticommutator

For the square of anticommutator, we have

[TABLE]

The last two terms are equal to zero and we obtain

[TABLE]

The evolution operators do not change the occupation of the $f$ states and the nonlocal contributions with $i\neq j$ contain products of two creation and two annihilation operators at the same sites, which vanish. Hence, in the $D=\infty$ limit, we have only local OTOCs that are nonzero. They are given by

[TABLE]

This expression can be rewritten in a form similar to the one for the $f$-electron propagator [31]

[TABLE]

where the integration and time ordering are performed over the doubly folded Keldysh contour (Fig. 3). Now the time-dependent field $U_{c}^{[0,t],[0,t]}(t^{\prime})$ is nonzero and equal to $U$ on the two forward branches for $t^{\prime}\in[0,t]$ and is zero otherwise. In the same way as in [31], we can derive the final expression for the OTOC:

[TABLE]

In a discretized form, $t_{m}=mW$ with $t_{0}=0$ and $t_{N}=t$, the first determinant looks like

[TABLE]

where

[TABLE]

is the lesser Green’s function,

[TABLE]

is the greater one, and for the equal times $t^{\prime}=t^{\prime\prime}$ we choose $\Theta(0)=0$ [34], which gives

[TABLE]

Similar expressions can be written for the second determinant in (57).

One might expect these results to be similar to those of the $f$-particle propagator. Indeed, they are. We find that the long time asymptotic of the determinants in (68) or $F_{\alpha}(t)$ in (58) are determined by the Szegő-Widom theorem for block Toeplitz matrices [37], which employ characteristic functions

[TABLE]

Now, the imaginary part of each characteristic function can become zero either at one frequency, $\omega=0$, or at two frequencies, $\omega=0$ and when $1\mp U\operatorname{Re}\mathcal{G}_{0,1}^{r}(\omega)=0$. We examine this OTOC in detail in the next section.

VI Results and Discussion

At half filling $n_{d}=n_{f}=0.5$, the contributions $F_{0}(t)$ and $F_{1}(t)$ in (58) are equal and the total OTOC $F(t)$ is real. But, the partial contributions from the sites with empty $f$ particle states, $F_{0}(t)$, and with occupied $f$ particle states, $F_{1}(t)$, display intriguing time dependencies. For $U<0.44$, the real part of both characteristic functions (72) is positive at frequencies where the imaginary part vanishes and the winding number is zero (see Appendix B). The OTOC function (57) displays exponential decay at large times for high temperatures, which is replaced by a power law as $T\to 0$ (Fig. 4). The corresponding logarithmic derivative

[TABLE]

approaches a constant value at large times for high temperatures, with a decay rate within the “chaos bound” $|L_{0}(t)|/(2\pi T)<1$ 111Despite the chaos bound of [8, 9] is applicable for the short timescale, we compare the decay rate with this temperature bound to elucidate thermal effects in relaxation., whereas it displays $1/t$ dependence for low temperatures, which confirms the power law asymptotics of the OTOC. Consequently, the corresponding “density-of-states” 222Now, DOS is defined by the real part of OTOC because (52) omits factor “$-\mathrm{i}$” in contrast to the one-particle Green’s function (14)

[TABLE]

displays a Lorentzian peak at high temperatures which is replaced by a power law singularity as $T\to 0$.

Next, for $0.44<U<1.235$, the real part of the characteristic function becomes negative at $\omega=0$, the logarithm of the characteristic function is now a multivalued function with a branch cut and the winding number is equal to $-1$ (see Appendix B). Now, the characteristic function (72) must be replaced by

[TABLE]

where $W=\Delta t\to 0$ is the time discretization, and the long-time asymptotic becomes more complicated (Fig. 5) with oscillations of the logarithmic derivatives $L_{0}(t)$ at high temperatures similar to the one for the single-particle propagators for $U>0.866$ (Fig. 2). With decreasing temperature, these oscillations are replaced by $1/t$ asymptotics, because we are still in the metallic phase $U<\sqrt{2}$. Now the OTOC’s “density-of-states” displays a two peak structure separated by a pseudogap at high temperatures with the development of the power law singularity in the center of the pseudogap as $T\to 0$. The sharp peak in the “density-of-states” is observed at low temperatures up to $U$ values where the “density-of-states” becomes negative, which reflects an increase in the quantum effects for these cases (such as Kirkwood-Dirac quasiprobabilities (KDQ) [40, 41, 42], which can also become negative).

For $U>0.866$, the imaginary part of the characteristic function (72) crosses the zero line at frequencies where $1\mp U\operatorname{Re}\mathcal{G}_{0,1}^{r}(\omega)=0$, but for $U<1.235$ the winding does not change because the real part is still positive at these points, whereas for $U>1.235$ the real part becomes negative and the winding number is equal to $-2$ (see Appendix B) resulting in a more complicated long time behavior of the OTOCs and their logarithmic derivatives with damped oscillations (Fig. 6). Now, the OTOC’s DOS is more robust, displays weak temperature dependence, and possesses negative values in a wide frequency interval. For large values of interaction $U$, the quantum correlations become dominant and the “flip” of an $f$-state produces strong disturbance for the itinerant $d$-electrons, which can subsequently lead to local heat fluctuations. Such quantum heat fluctuations can be associated with negativity of quasiprobability distribution, e.g. Kirkwood-Dirac quasiprobability, which, on the other hand, can be a signature of suppressing of the decoherence and of quantum chaos and scrambling [43, 44].

However, the sum rule (integral over frequency) for OTOC’s DOS is conserved and is equal to one:

[TABLE]

which means that the appearance of negative DOS at small frequency values leads to a DOS enhancement elsewhere. For the first moment we obtain

[TABLE]

which equals zero at half filling $n_{d}=n_{f}=0.5$, and one can proceed similarly and systematically to determine higher-order moments.

It should be noted, that the Szegő-Widom theorem for block Toeplitz matrices [37] is applicable for smooth functions. For this case in the fermionic FKM, the functions are continuous at high temperatures. At low temperatures, both the matrix elements in (68) and characteristic functions (72) develop step-like features from the Fermi factors and the OTOC’s long time behavior is now governed by localization singularities given by the Anderson orthogonality catastrophe [45, 35, 34]. These cases of three different winding numbers for the characteristic function are also distinguished by the behavior of the logarithmic derivative at low temperatures (see Appendix C, Fig. A.6).

At high temperatures, the oscillatory behavior of the logarithmic derivative is similar to the one seen in [12, 13] for quantum spin models behind the light cone after the wave front passes. This analogy is plausible, because we consider a local OTOC function, which is always behind the wave front and at high temperature the effects of the Fermi statistics become negligible and the electrons behave like Pauli particles. Nevertheless, the localized $f$ particles act as impurity scattering centers for the itinerant $d$ electrons and the destruction of the oscillations of logarithmic derivative for large Coulomb interaction $U$ is similar to the one observed in [13] with the an increase of disorder, when one observes the post wavefront ripples damped by disorder. This static disorder can protect the fast scrambling, which can be thought of as a breaking of the “bound on chaos”, when at long times we have $|L_{0}(t)|>2\pi k_{B}T/\hbar$, which, is also manifested by the negative OTOC “density-of-states” $A(\omega)$, see Fig. 6. According to Yunger Halpern et al. [41], OTOCs can be considered as moments of the Kirkwood-Dirac quasiprobabilities which implies that negative values of the “density-of-states” $A(\omega)$ are signatures of quantum information scrambling.

In cases where the calculated OTOC function decays with time, the question arises “where has the corresponding quantum information gone?” For lattice models with direct interactions between spins, it spreads to the nearby sites due to quantum information scrambling. For the FKM, there is no direct interaction between $f$ particles, only an effective retarded one mediated through the itinerant $d$ electrons, and in addition to quantum scrambling, quantum information can be transferred to the $d$ electrons. Exactly how this occurs is a project for future investigations.

VII Conclusions

We have examined how an indirect retarded interaction between qubits can affect the scrambling of quantum information. We considered the dynamics of the localized two-level states ($f$-electrons) interacting with the itinerant $d$-electrons in the Falicov-Kimball model and found exact results for the OTOCs constructed via state-flip operators ($f$-electron creation and annihilation operators) within dynamical mean field theory. Despite the classical (Ising) nature of our two-level states (the occupation $n_{f}$ is a conserved quantity), their dynamics is nontrivial and is governed by the winding of the characteristic functions for the corresponding Toeplitz matrices. With an increase of the local interaction potential $U$ between $f$ and $d$ electrons, the logarithmic derivative of the OTOCs start to display oscillations at high temperatures similar to those seen in quantum spin models. On the other hand, at low temperatures Fermi statistics and screening effects related to the orthogonality catastrophe become important, which strongly affects the OTOC’s behavior. With a further increase of the interaction $U$, we find a fast damping of these oscillations due to the effects of disorder and, on the other hand, makes the OTOC’s “density-of-states” (moments of the Kirkwood-Dirac quasiprobabilities) become negative, which is another manifestation of the quantum information scrambling.

Our results were obtained at half filling with $n_{d}=n_{f}=0.5$, but our main conclusions remain valid in when $n_{f}+n_{d}=1$ holds, with the remark that now asymptotic behavior of each contribution $F_{0}(t)$ and $F_{1}(t)$ in (57) will be different with different windings of the corresponding Toeplitz matrices. On the other hand, in the doped case, with $n_{f}+n_{d}\neq 1$, our system will be in a metallic phase and the low temperature behavior will be dominated by the localization singularities given by the Anderson orthogonality catastrophe, whereas at high temperatures the results will be similar to the one for un-doped case.

Acknowledgements.

This work was supported by the Department of Energy, Office of Basic Energy Sciences, Division of Materials Sciences and Engineering, under Grant No. DE-FG02-08ER46542. J.K.F. was also supported by the McDevitt bequest.

Data Availability

The data that support the findings of this article are openly available [46].

Appendix A Winding of characteristic function for one-particle propagator

According to the Szegő’s theorem (see [31] and references therein), the long time asymptotic of the $N\times N$ Toeplitz determinant

[TABLE]

is equal to

[TABLE]

with

[TABLE]

Here, $C(\xi)=\sum_{n=-\infty}^{+\infty}c_{n}\xi^{n}$ with $\xi=\mathrm{e}^{\mathrm{i}\theta}$ lying on the unit circle. The requirement for $C(\xi)$ to be continuous on the unit circle, is that its winding number (or index) is equal zero, namely

[TABLE]

and we also require that the summation is well defined, or $\sum_{n=-\infty}^{\infty}|c_{n}|<\infty$ for the theorem to hold.

After discretization of the time interval $[0,t]$ by $t_{n}=n\Delta t$, with $\Delta t=t/N$, the matrix elements $c_{n}$ in Eq. (15) satisfy

[TABLE]

By choosing $n\Delta t=t^{\prime}$ and $\theta=\omega\Delta t$, in the limit $\Delta t\to 0$ we obtain for the determinant the result

[TABLE]

where

[TABLE]

and $C(\omega)$ is defined by Eq. (22). Now, the previous condition (81) for the winding number is replaced by a requirement of analyticity of the integrand $\ln C_{\alpha}^{f}(\omega)$. This means that the real part of $C_{\alpha}^{f}(\omega)$ must be positive when the imaginary part changes sign; in other case its logarithm will become a multi-valued function, requiring a branch cut. This criterion can be easily visualized by the Nyquist (Cole-Cole) type plots: at half filling and for $U<0.866$, the Nyquist plot does not wind around the origin—that is, the winding number is zero, see Fig. A.1, whereas, for $U>0.866$, the Nyquist plot does wind around the origin and the winding number is $-1$, see Fig. A.2.

Appendix B Winding of characteristic function for OTOC

Similar to the $f$-particle propagator, the long time asymptotic of the determinants in (68) or $F_{\alpha}(t)$ in (58) is determined according to the Szegő-Widom theorem for block Toeplitz matrices [37] by the same expressions (79)–(84), but with the characteristic functions $C_{0,1}(\omega)$, see Eq. (72),

[TABLE]

Now, the imaginary part of each characteristic function can become zero either at one frequency, $\omega=0$, or at two frequencies, $\omega=0$ and when $1\mp U\operatorname{Re}\mathcal{G}_{0,1}^{r}(\omega)=0$. At half filling and for $U<0.44$, the Nyquist plots for $C_{0,1}(\omega)$ do not wind around the origin, Fig. A.3, and the winding number is zero. For $0.44<U<1.235$, the characteristic functions $C_{0,1}(\omega)$ wind once around the origin as shown in Fig. A.4, and the winding number is $-1$, whereas for $U>1.235$ it winds twice, as shown in Fig. A.5,—here, the winding number is $-2$.

Appendix C Winding of the logarithmic derivative

The existence of three different winding numbers for the characteristic function is easily seen on the Nyquist plots for the logarithmic derivative (73) too, see Fig. A.6 for an example. For $U<0.44$, the logarithmic derivative smoothly evolves to some steady value. For $0.44<U<0.866$ and at high temperatures, $L_{0}(t)$ first displays evolution similar to what was seen in the previous case. It then is replaced by circulations, whereas at low temperatures we still observe a smooth evolution to a steady value. For $0.866<U<1.235$, we observe circulations of the logarithmic derivative at high temperatures, which are destroyed (become more chaotic) with decreasing temperature and are replaced by some complicated looped trajectories. For $U>1.235$, due to the presence of pseudogap and Mott gap, the trajectories become even more chaotic, Fig. A.6.