Dynamics of quasiperiodically driven spin systems

TL;DR

This paper investigates the complex dynamics of a spin system driven by Fibonacci-sequence kicks, revealing classical chaos with zero Lyapunov exponent and quantum fractality, along with ergodic thermalization when interacting with another spin.

Contribution

It introduces a detailed analysis of Fibonacci-driven spin dynamics, deriving invariants, and contrasting classical chaos with quantum fractality and thermalization behaviors.

Findings

Classical phase space shows chaos with zero Lyapunov exponent.

Quantum states exhibit fractal fluctuations in spin operator means.

Interaction induces ergodic dynamics and infinite temperature thermalization.

Abstract

We study the stroboscopic dynamics of a spin-$S$ object subjected to $\delta$-function kicking in the transverse magnetic field which is generated following the Fibonacci sequence. The corresponding classical Hamiltonian map is constructed in the large spin limit, $S \rightarrow \infty$. Upon evolving such a map for large kicking strength and time period, the phase space appears to be chaotic; interestingly, however, the geodesic distance increases linearly with the stroboscopic time implying that the Lyapunov exponent is zero. We derive the Sutherland invariant for the underlying $SO(3)$ matrix governing the dynamics of classical spin variables and study the orbits for weak kicking strength. For the quantum dynamics, we observe that although the phase coherence of a state is retained throughout the time evolution, the fluctuations in the mean values of the spin operators exhibit fractality which is also present in the Floquet eigenstates. Interestingly, the presence of an interaction with another spin results in an ergodic dynamics leading to infinite temperature thermalization.

Figures

Equations

$$\hat{H}(t)={\omega }_0{\hat{S}}_z+\lambda {\hat{S}}_x\sum _{n=-\infty }^{\infty }\delta \left(t-\sum _n{T}_n\right),$$

$$\text{T+, T−, T+, T+, T−, T+, T−, T+, ⋯}.$$

$${\hat{\mathcal{F}}}_n=e^{-i{T}_n{\hat{S}}_z}{e}^{-i\lambda {\hat{S}}_x}.$$

$${\hat{S}}_x^{n+1}$$

$${\hat{S}}_y^{n+1}$$

$${\hat{S}}_z^{n+1}$$

$$J_{n}~{}=~{}\left(\begin{array}[]{ccc}\cos T_{n}&-\sin T_{n}~{}\cos\lambda&\sin T_{n}~{}\sin\lambda\\ \sin T_{n}&\cos T_{n}~{}\cos\lambda&-\cos T_{n}~{}\sin\lambda\\ 0&\sin\lambda&\cos\lambda,\end{array}\right),$$

$${\hat{\mathcal{F}}}_1=e^{-i\lambda {\hat{S}}_x},{\hat{\mathcal{F}}}_2=e^{-i2T_0{\hat{S}}_z}{e}^{-i\lambda {\hat{S}}_x},$$

$${\hat{\mathcal{F}}}_{m+2}={\hat{\mathcal{F}}}_{m+1}{\hat{\mathcal{F}}}_m,$$

$$\left(\begin{array}[]{c}s_{x}^{n+1}\\ s_{y}^{n+1}\\ s_{z}^{n+1}\end{array}\right)~{}=~{}J_{n}~{}\left(\begin{array}[]{c}s_{x}^{n}\\ s_{y}^{n}\\ s_{z}^{n}\end{array}\right),$$

$$U_m=e^{i{\alpha }_m{\hat{n}}_m\cdot \vec{\sigma }}=\cos {\alpha }_m{\mathbb{I}}_2+i\sin {\alpha }_m{\hat{n}}_m\cdot \vec{\sigma },$$

$$I_s=x_m^2+x_{m+1}^2+x_{m+2}^2-2x_m{x}_{m+1}{x}_{m+2}-1$$

$$R_{m+2}=R_{m+1}{R}_m.$$

$$R_m=e^{i{\varphi }_m{\hat{e}}_m\cdot \vec{T}},$$

$$R_{ij}^m$$

$$J_1=e^{-i\lambda T_x}\text{and}{J}_2=e^{-i2T_0{T}_z}{e}^{-i\lambda T_x},$$

$$I_s=-[\sin T_0\sin (\lambda /2){]}^2.$$

$$J_m=e^{-i{\beta }_G^m[\lambda T_x+(2T_0/{\beta }_G)T_z]/\sqrt{5}},$$

$${\varphi }_m=-\frac{{\beta }_G^m}{\sqrt{5}}\bar{\lambda }\text{and}{\hat{e}}_m=\frac{1}{\bar{\lambda }}(\lambda ,0,2T_0/{\beta }_G),$$

$$\left(\begin{array}[]{c}s_{x}^{m}\\ s_{y}^{m}\\ s_{z}^{m}\end{array}\right)=\left(\begin{array}[]{c}e_{1}e_{3}(1-\cos\phi_{m})-e_{2}\sin\phi_{m}\\ e_{2}e_{3}(1-\cos\phi_{m})+e_{1}\sin\phi_{m}\\ \cos\phi_{m}+e_{3}^{2}(1-\cos\phi_{m})\end{array}\right).$$

$$(s_x^m,s_y^m)=(e_1{e}_3(1-\cos {\varphi }_m),e_1\sin {\varphi }_m),$$

$$|\Theta ,\Phi ⟩=(1+|z{|}^2{)}^{-S}{e}^{z{\hat{S}}_{+}}|S,-S⟩,$$

$$|\psi (m)⟩={\hat{\mathcal{F}}}_m|\psi (0)⟩.$$

$$|\phi ⟩=\frac{1}{\sqrt{2S+1}}\sum _{l=1}^{2S+1}{e}^{il\phi }|l⟩,$$

$$p(\phi )=|⟨\phi |\psi (m)⟩{|}^2.$$

$$X_{\Omega }=\sum _{m=1}^N{x}_m{e}^{i2\pi \Omega m},$$

$$I_q^{\nu }=\sum _{m=1}^{2S+1}|{\chi }_{\nu }(m){|}^{2q},I_q=\frac{1}{2S+1}\sum _{\nu }{I}_q^{\nu }\sim (2S+1{)}^{-{\tau }_q},$$

$$D_q=\underset{S\to \infty }{\lim }\frac{S_q^{\nu }}{\log (2S+1)},S_q^{\nu }=-\sum _{m=1}^{2S+1}\frac{|{\chi }_{\nu }(m){|}^{2q}}{(q-1)}.$$

$$\hat{H}(t)$$

$${\hat{H}}_0$$

$${\hat{\mathcal{F}}}_1^I=e^{-iT{\hat{H}}_0}{e}^{-i\lambda {\hat{S}}_x^A},{\hat{\mathcal{F}}}_2^I=e^{-iT{\hat{H}}_0}{e}^{-i\lambda {\hat{S}}_x^B}.$$

$$|{\psi }_{AB}(0)⟩=|\Theta ,\Phi {⟩}_A\otimes |\Theta ,\Phi {⟩}_B.$$

Full text

I Introduction

In recent years quasiperiodic systems have attracted a lot of interest in various contexts ranging from quasicrystals quasi_lattice ; q_crys_book ; q_crys ; q_crys1 ; q_crys2 to localization-delocalization transitions Inguscio ; hanggi ; Huse ; Bloch ; Bloch1 , multifractality mult_frac ; kohmoto ; mult_frac1 ; frac , topological phases q_topo , and so on. The creation of a quasiperiodic potential in one dimension using bichromatic optical lattices Inguscio has led to the realization of the well known Aubry-André model AA which has been studied extensively both theoretically hanggi ; Huse and experimentally Bloch ; Bloch1 ; AA_light , in the context of observing localization phenomena, particularly many-body localization in interacting systems MBL_rev ; abanin_rmp . In the presence of a periodic drive such many-body localized states exhibit drive-induced delocalization and thermalization of isolated interacting quantum systems Bloch_drive ; abanin_drive , and its connection with the underlying chaotic dynamics and random matrix theory has been explored Sray_drive . As an extension, several interesting questions can be addressed related to the dynamical behavior of quantum systems under a quasiperiodic drive. One such issue is the emergence of steady states in quasiperiodically driven interacting quantum systems a_sen ; a_dutta and most interestingly, its connection with spectral properties and random matrix theory which is related to ergodicity.

One way to generate a quasiperiodic drive is by perturbing the system under consideration following the Fibonacci sequence; such a sequence has a rich mathematical structure giving rise to an invariant of the corresponding dynamical systems pandit_inv ; kohmoto_inv ; Sutherland . This way of generating a quasiperiodic drive can provide an alternate way to study the quasiperiodic structures observed in Fibonacci lattices kohmoto_inv ; Sutherland ; Ksaha . The Fibonacci drive can also be generated for a series of incommensurate frequencies known as metallic means, a common example of which is the golden ratio, $\beta=(\sqrt{5}+1)/2$. These driving protocols can give rise to the realization of a strange non-chaotic attractor Ruelle leading to a fractal-like dynamics which has been theoretically studied for dynamical maps with quasiperiodicity fractal ; pikovsky ; lai ; ramaswamy as well observed experimentally SNA_exp . The evolution of a spin-$1/2$ system under quasiperiodic perturbation reveals various interesting dynamical behaviors and temporal correlations orland ; lebowitz ; geisel_quasi . A quasiperiodic drive can also lead to slow relaxation to non-equilibrium steady states which has been investigated for interacting spin systems in the presence of a disordered magnetic field potter_mag . It is interesting to investigate if there is any underlying fractality in the eigenmodes of the time evolution operator even for non-interacting systems under a quasiperiodic drive and the fate of such critical states in the presence of interactions.

In this work we first consider a simple model of a non-interacting spin of magnitude $S$ (which has a well-defined classical limit) subjected to a quasiperiodic drive. Our objective is two-fold. First, we want to study the dynamics of the corresponding classical system and analyze the fluctuations around the mean values of the spin components. Second, we want to study the quantum dynamics and spectral properties of the time evolution operator, particularly its change in the presence of interaction. The paper is organized as follows. In Sec. II we describe the time-dependent Hamiltonian and the details of the quasiperiodic driving protocol. We then study the dynamics of the corresponding classical system by suitably taking the limit $S\rightarrow\infty$ in Sec. II.1. This is followed by the derivation of an invariant for this dynamical system in Sec. II.2. The quantum dynamics of a finite spin $S$ and the spectral properties of the time evolution operator are discussed in Sec. II.3. The effects of interaction are studied in Sec. III. Finally, we summarize our results, discuss possible experiments which can test our results, and conclude in Sec. IV.

II Driven spin model in transverse magnetic field

The dynamics of a spin-$S$ particle in the presence of a time-dependent magnetic field can be described by the Hamiltonian

[TABLE]

where the first term in the Hamiltonian represents the effect of a magnetic field with strength $\omega_{0}$ applied along the $\hat{z}$ direction (we have absorbed the gyromagnetic ratio in the definition of $\omega_{0}$), and the second term represents a $\delta$-function kicking due to a transverse magnetic field in the $\hat{x}$ direction with strength $\lambda$. Here $\hat{S}_{x,y,z}$’s denote the spin angular momentum operators and $T_{n}$ is the time lapse between the $(n-1)$-th and $n$-th kicks. Here we consider the case $T_{n}=T_{0}(1\pm\epsilon)\equiv T_{\pm}$, i.e., the time lapse can take two values $T_{+}$ or $T_{-}$ which follows the Fibonacci sequence:

[TABLE]

We will scale the energy (time) by $\omega_{0}$ $(1/\omega_{0})$ and set $\hbar=1$ throughout this paper.

The time evolution operator describing the dynamics of such system between the $(n-1)$-th and $n$-th kicks is given by

[TABLE]

Note that for a periodic drive in which all the $T_{n}$’s are equal, Eq. (3) reduces to the usual Floquet operator. Following the time evolution of any operator between the $(n-1)$-th and $n$-th kicks under the Floquet matrix $\hat{\mathcal{F}}_{n}$ (namely, $\hat{A}_{n+1}=\hat{\mathcal{F}}_{n}^{\dagger}\hat{A}_{n}\hat{\mathcal{F}}_{n}$) Floquet , we obtain the Heisenberg equations of motion for the spin operators as follows,

[TABLE]

This is a linear map for the spin operators of the form: $(\hat{S}_{x}^{n+1},\hat{S}_{y}^{n+1},\hat{S}_{x}^{n+1})=J_{n}(\hat{S}_{x}^{n},\hat{S}_{y}^{n},\hat{S}_{x}^{n})$, where the transfer matrix $J_{n}$ can be written as

[TABLE]

where the $T_{n}$’s are given in Eq. (2). We will consider the case $\epsilon=1$, so that Eq. (3) becomes

[TABLE]

where $\hat{\mathcal{F}}_{1}$ represents kicking the spin-S object by a magnetic field in the $\hat{x}$ direction, and $\hat{\mathcal{F}}_{2}$ corresponds to time evolution of the system under $\hat{S}_{z}$ for time interval $2T_{0}$ followed by another kick in $\hat{S}_{x}$.

Starting from two such $SU(N)$ matrices, $\hat{\mathcal{F}}_{1}$ and $\hat{\mathcal{F}}_{2}$, the successive Floquet operators in a Fibonacci sequence can be generated using the recursion relation,

[TABLE]

where the initial matrices $\hat{\mathcal{F}}_{1}$ and $\hat{\mathcal{F}}_{2}$ are given in Eq. (6). We would like to point out that at a Fibonacci time step $m$, the stroboscopic time is given by $n=F_{m}$, where $F_{m}$ is the $m$-th Fibonacci number. For large $m$, the stroboscopic time increases exponentially as $n\sim e^{\beta_{G}m}$ where, $\beta_{G}=(\sqrt{5}+1)/2$ is the golden ratio. The advantage of using such a recursion relation is that one can numerically obtain the steady state of the system after a very long time scale. Henceforth, we will adopt Eq. (7) to study the dynamics for a very large duration in stroboscopic time.

II.1 Classical Dynamics

We will first discuss the dynamics of the corresponding classical system. The classical limit of such a spin system can be obtained by considering the large spin limit, $S\rightarrow\infty$. Then the spin variables $\hat{S}_{x,y,z}$ can be classically described by the components of a spin vector $\vec{S}\equiv(S_{x},S_{y},S_{z})$. We scale the spin operators by the magnitude S to obtain the classical spin variables, $s_{i}=S_{i}/S$ which follow the commutation relations $[s_{i},s_{j}]=i\epsilon_{ijk}s_{k}/S$. In the limit $S\rightarrow\infty$, the commutators vanish and the variables become classical. Thus using Eqs. (4), the stroboscopic time evolution of the corresponding classical spin variables in between consecutive kicks can be described by the following linear Hamiltonian map,

[TABLE]

where the transfer matrix $J_{n}$ is given in Eq. (5). By evolving Eq. (8) stroboscopically in time we obtain the trajectories on a unit sphere as shown in Fig. 1 for different driving parameters. We observe that for a small driving strength $\lambda$ and time period $T_{0}$ the trajectories are regular and precess over time as depicted in Fig. 1 (a). Over a small time scale such a regular trajectory is plotted in the projected plane of $s_{x}-s_{y}$ in Fig. 2 (a). However, for large values of $\lambda$ and $T_{0}$, the dynamics is no longer regular and eventually covers the whole surface of the sphere shown in Fig. 1 (b). It is interesting to note that the transfer matrix in Eq. (5) is unimodular and therefore its eigenvalues have the form: $1,e^{\pm i\varepsilon}$, $\varepsilon$ being the eigenphase. As a result, the Lyapunov exponent always turns out to be zero blumel ; schack ; kurths . To further illustrate this, we compute the growth of the geodesic distance on the Bloch sphere. We start from the initial point $s_{i}=(s_{x},s_{y},s_{z})=(0,0,1)$ and evolve it under successive kicking following Fibonacci sequence. The resulting trajectory in $s_{x}-s_{y}$ plane is depicted in Fig. 2a. The geodesic distance between the initial point $s_{i}$ and the time-evolved point $s_{f}$ is given by $d=\cos^{-1}({\vec{s}}_{i}\cdot{\vec{s}}_{f})$. In Fig. 2b we have plotted $d$ as a function of the stroboscopic time $n$; we see a linear growth which indicates that the Lyapunov exponent is zero in this case.

II.2 Fibonacci sequence of $SO(3)$ matrices and Sutherland invariant

We first discuss the case of $SU(2)$ matrices multiplied according to a Fibonacci sequence for which Sutherland found an invariant; then we will discuss how this invariant generalizes to the case of $SO(3)$ matrices. Starting with two $SU(2)$ matrices $U_{1}$ and $U_{2}$, we generate a Fibonacci sequence of matrices defined by the recursion relation $U_{m+2}=U_{m+1}U_{m}$. Let us parametrize $U_{m}$ as,

[TABLE]

where $0\leq\alpha_{m}\leq\pi$, $\hat{n}_{m}$ is a unit vector, $\mathbb{I}_{2}$ denotes the $2\times 2$ identity matrix, and $\vec{\sigma}=(\sigma_{x},\sigma_{y},\sigma_{z})$ denotes the Pauli matrices. Defining $x_{m}=\frac{1}{2}\text{tr}(U_{m})=\cos\alpha_{m}$, Sutherland showed that the quantity

[TABLE]

is independent of $m$ Sutherland .

We now consider $SO(3)$ matrices denoted as $R_{m}$. Starting with two such matrices $R_{1}$ and $R_{2}$, we generate a Fibonacci sequence using the recursion relation

[TABLE]

Let us parametrize $R_{m}$ as follows,

[TABLE]

where $0\leq\phi_{m}\leq 2\pi$, $\hat{e}_{m}$ is a unit vector and $\vec{T}=(T_{x},T_{y},T_{z})$ are the generators of $SO(3)$ matrices. One can show that the matrix elements of $R_{m}$ given by $R^{m}_{ij}$ and the components of $\hat{e}_{m}$ given by $e_{i}^{m}$ are related as

[TABLE]

where $\epsilon_{ijk}$ is the totally antisymmetric matrix with $\epsilon_{123}=1$. Now using the standard mapping between the spin-1/2 and spin-1 representations of the angular momentum group, $\alpha_{m}=\phi_{m}/2$ and $\hat{n}_{m}=\hat{e}_{m}$ one can obtain the Sutherland invariant $I_{s}$ for a given $SO(3)$ matrix $R_{m}$. [There is a subtlety here: since the same $SO(3)$ matrix $R_{m}$ corresponds to two different $SU(2)$ matrices, $U_{m}$ and $-U_{m}$, whose traces (divided by 2) are given by $x_{m}$ and $-x_{m}$, one has to check at each step of the Fibonacci sequence which of the two possible values of the trace gives the correct value of the Sutherland invariant in Eq. (10).]

As an example we consider the transfer matrix $J_{n}$, given in Eq. (5) which is a $SO(3)$ matrix. Starting with two such matrices $J_{1}$ and $J_{2}$ given by

[TABLE]

we find that the corresponding Sutherland invariant is

[TABLE]

We have checked numerically up to a large Fibonacci step $m\sim 1000$ that $I_{s}$ remains constant during the time evolution.

It is known that the $m$-th Fibonacci number, given by $F_{m}=\left(\beta_{G}^{m}-\frac{(-1)^{m}}{\beta_{G}^{m}}\right)/\sqrt{5}$, quickly approaches the value $\beta_{G}^{m}/\sqrt{5}$ as $m$ increases. We then find that

[TABLE]

which has been derived under the approximation that $\lambda$, $T_{0}\ll 1$ so that the commutators arising from $[T_{i},T_{j}]$ do not grow much within a small time scale. This leads us to define

[TABLE]

where $\bar{\lambda}=\sqrt{\lambda^{2}+4T_{0}^{2}/\beta_{G}^{2}}$. We will now study what happens when $R_{m}$ acts on the column $(s_{x},s_{y},s_{z})=(0,0,1)$ as numerically shown in Fig. 2 (a). Using Eq. (13), we see that

[TABLE]

We then see that

[TABLE]

where $e_{1}=\lambda/\bar{\lambda}$ and $e_{3}=2T_{0}/\beta_{G}\bar{\lambda}$. We see that the point $(s_{x}^{m},s_{y}^{m})$ describes an ellipse as shown by the solid line in Fig. 2 (a) whose center lies at $(e_{1}e_{3},0)$ and the lengths of the axes are $2e_{1}e_{3}$ and $2e_{1}$ in the $\hat{x}$ and $\hat{y}$ directions respectively.

We would like to mention here that quasiperiodic driving of $SU(2)$ matrices has been studied in detail in Refs. a_sen and a_dutta . It has been found that the long-time behavior of the system depends to some extent on the value of the Sutherland invariant $I_{s}$ a_sen . The behavior is particularly simple near $I_{s}=0$ and $-1$ (which are respectively the maximum and minimum possible values of $I_{s}$). Near $I_{s}=0$, the trajectory is given by a circle on the Bloch sphere a_dutta . This is similar to the behavior of $SO(3)$ matrices discussed above, namely, the point moves on an ellipse when $T_{0}$ and $\lambda$ are small which corresponds to a very small value of $I_{s}$ according to Eq. (15).

II.3 Quantum Dynamics

We now analyze the quantum dynamics of a spin-$S$ object governed by the time-dependent Hamiltonian in Eq. (1). To this end we construct the initial wave function $|\psi(0)\rangle$ from a spin coherent state given by radcliff

[TABLE]

where $\Theta$ and $\Phi$ are the polar and azimuthal angles, respectively, representing the orientation of the classical spin vector of magnitude $S$, and $z=e^{-i\Phi}\tan(\Theta/2)$. The time-evolved state after the $m$-th Fibonacci kick is given by

[TABLE]

To compute the distribution of the relative phases of the time-evolved state we first construct the phase state given by oberthaler ,

[TABLE]

where $\varphi=\varphi_{0}+2\pi l^{\prime}/(2S+1)$ and $l^{\prime}\in[1,2S+1]$. We choose $\varphi_{0}=-\pi$ so that the relative phase lies in the range $-\pi$ to $\pi$. The phase distribution can be obtained by projecting $|\psi(m)\rangle$ onto the phase state as given by

[TABLE]

In Fig. 3 (a) we have shown the time evolution of $p(\varphi)$, and in Fig. 3 (b) the snapshots of $p(\varphi)$ at different times are plotted. We observe the phase distribution $p(\varphi)$ remains a highly peaked function, however, its peak position changes under stroboscopic evolution.

It indicates that the phase diffusion does not take place even for large kicking strengths, and the phase coherence is not lost during the time evolution. We also compute the expectation values of the spin operators, i.e., $\langle\hat{S}_{z}\rangle$ where the average $\langle.\rangle$ is taken with respect to the time-evolved state $|\psi(m)\rangle$. We have verified that under quantum dynamics $\langle\hat{S}_{i}\rangle$’s are in agreement with the classical variables $s_{i}$ obtained from the dynamical map. The time evolution of the mean values of these operators exhibit fluctuations similar to the peak of $p(\varphi)$ in Fig. 3.

To understand the nature of the fluctuations in the dynamics of spin variables, we compute the cumulative sum of the Fourier transform of $\langle\hat{S}_{z}\rangle$ given by pikovsky ; lai ; ramaswamy ; orland ; lebowitz ,

[TABLE]

where $x_{m}=\langle\hat{S}_{z}\rangle$ computed at the $m$-th Fibonacci step, and $\Omega$ is the frequency. From the power spectrum, $X_{\Omega}$ vs $\Omega$, we observe that several frequency modes are present in the fluctuations which confirms that the dynamics is not at all regular even for small values of $\lambda$ and $T_{0}$ for which the spin variables are seen to undergo a precessional motion in Fig. 1 (a).

The corresponding dynamics in the $\text{Re}[X_{\Omega}]-\text{Im}[X_{\Omega}]$ plane exhibits a fractal-like structure as depicted in Fig. 4 (a). The degree of such fractal motion can be quantified using the relation $|X_{\Omega}|^{2}\sim N^{\beta}$, where the exponent $\beta=1\,(2)$ signifies random (regular) paths respectively, whereas $\beta\neq 1,2$ corresponds to a fractal path in the $\text{Re}[X_{\Omega}]-\text{Im}[X_{\Omega}]$ plane which can be observed from the logarithmic plot in Fig. 4b. We have computed $\beta$ for different values of $\Omega$, i.e., $\beta(\Omega=0.2)=1.436$ and $\beta(\Omega=0.4)=1.1$, for $\lambda=\pi/50$ and $T_{0}=\pi/100$. We have checked that for different choices of $\lambda$ and $T_{0}$ the behavior remains qualitatively similar.

Next, we investigate the spectral properties of the Floquet operator, $\hat{\mathcal{F}}_{m}$, after a sufficiently large Fibonacci step $m$. Due to unitarity, we have $\hat{\mathcal{F}}_{m}|\chi_{\nu}\rangle=e^{i\varepsilon_{\nu}}|\chi_{\nu}\rangle$, where $\varepsilon_{\nu}\in[-\pi,\pi]$ and $|\chi_{\nu}\rangle$ are the eigenphase and eigenvector corresponding to the $\nu$-th Floquet eigenmode. We compute the moments of the eigenstates $|\chi_{\nu}\rangle$ of the Floquet operator $\hat{\mathcal{F}}_{m}$ at the $m$-th Fibonacci step using the relation given by

[TABLE]

where $\chi_{\nu}(m)=\langle\chi_{\nu}|\alpha_{m}\rangle$, $|\alpha_{m}\rangle$ being the computational basis. The exponent $\tau_{q}$ is related to the fractal dimension $D_{q}$ as, $\tau_{q}=D_{q}(q-1)$ geisel1 ; moessner ; Meo ; lfsantos . The fractal dimension can equivalently be computed from lfsantos ; bogomolny ; fradkin

[TABLE]

$S_{q}^{\nu}$ is the Rényi entropy corresponding to the $\nu$-th eigenvector $|\chi_{\nu}\rangle$. In Figs. 5 (a) and (b) we have plotted $\log I_{q}$ and the average Rényi entropy $S_{q}=\sum_{\nu}S_{q}^{\nu}/(2S+1)$ respectively as a function of $\log(2S+1)$ for different values of $q$. From the slope of the linear fitting as shown in Fig. 5 (a) we obtain $\tau_{q}$ which is plotted versus $q$ in Fig. 5 (c) exhibiting a nontrivial behavior which confirms the existence of fractality in the eigenvectors. The fractal dimension is obtained from the slope of the dashed line at large $q$, $D_{q}\sim 0.6$ as shown in Fig. 5 (c). These plots have been generated after a sufficiently large number of Fibonacci steps, say, $m\sim 30$, after which the behavior does does not change with $m$ as shown in Fig. 5. We also observed that the qualitative behavior remains the same for different choices of the driving parameters $\lambda$ and $T_{0}$.

Further, we explored the fractality in the eigenspectrum by computing the local number of eigenstates $\Delta N_{\varepsilon}$ within the interval $\Delta\varepsilon$ around the eigenphase $\varepsilon$. Typically $\Delta N_{\varepsilon}$ follows the relation $\Delta N_{\varepsilon}\sim(\Delta\varepsilon)^{\alpha}$ as depicted in the $\log-\log$ plot in Fig. 5 (d) for $\varepsilon\sim 0$, i.e., at the center of the band. From the slope of the linear fitting in Fig. 5 (d) we obtain $\alpha\sim 0.69$ which signifies fractality in the Floquet spectrum kohmoto ; geisel1 ; geisel ; pandit .

III Effects of interactions

To study the effects of interactions, we consider a system of two interacting spins and the $\delta$-function kicking is applied to both the spins following a Fibonacci sequence. Our goal is to study the effect of interactions on the quantum dynamics governed by the quasiperiodic drive. The Hamiltonian describing this system is given by

[TABLE]

This is a simple generalization of the non-interacting model discussed above, where $J$ is the strength of the interaction between the two spins, and the last term represents kicking applied to the two spins following a Fibonacci sequence. The Floquet operators for the first two Fibonacci step are therefore given by

[TABLE]

The subsequent matrices in the Fibonacci sequence are then generated using Eq. (7). By diagonalizing the Floquet matrix $\hat{\mathcal{F}}_{m}^{I}$ we obtain the eigenphases $\varepsilon_{\nu}$ and the corresponding eigenvectors $|\chi_{\nu}\rangle$ (see Sec. II.3). We first compute the spacing between the successive eigenphases, i.e., $\delta_{\nu}=\varepsilon_{\nu+1}-\varepsilon_{\nu}$, and calculate the distribution of the $\delta_{\nu}$’s following the procedure described in Ref. Haake in order to keep the normalization $\int P(\delta)d\delta=1$ and mean $\int\delta P(\delta)d\delta=1$. In Figs. 6 (a) and (b) we have plotted the normalized quasienergy spacing distribution $P(\delta)$ for small and large values of $J$ respectively. We note that in the presence of interactions, the fractality of the system vanishes and level repulsion sets in. For large $J$, $P(\epsilon)$ resembles the Wigner surmise and corresponds to the Gaussian orthogonal ensemble of random matrix theory as shown in Fig. 6 (b); this indicates a possible thermalization of the system for large values of the interaction strength $J$ sray .

We compute the moments of the Floquet eigenstates using Eq. (25) followed by a calculation of the exponent $\tau_{q}$ which is related to the fractal dimension $D_{q}$. In Fig. 6 (c) we have plotted $\log I_{q}$ vs $2\log(2S+1)$ for different values of $q$. From the linear fitting we obtain the slope $\tau_{q}$ which we have plotted in Fig. 6 (b) as a function of $q$ and compared with the case $J=0$. We note that in contrast to the non-interacting case, we find for finite values of the interaction strength $J$ that $\tau_{q}\sim q$; this implies $D_{q}\sim 1$ indicating an ergodic nature of the Floquet eigenstates.

We further elucidate this fact from the wave packet dynamics. We construct an initial wave function from the product of two spin coherent states given by

[TABLE]

The expectation values of the corresponding spin observables $\langle\hat{S}_{x,y,z}^{A,B}\rangle$ can be computed from the time-evolved wave function $|\psi_{AB}(m)\rangle=\hat{\mathcal{F}}_{m}^{I}|\psi_{AB}(0)\rangle$. We observe that for small values of the interaction strength $J$, $\langle\hat{S}_{x,y,z}^{A,B}\rangle$ saturates to different nonzero values, whereas for large $J$, all the spin observables decays to zero as depicted in Figs. 7 (a) and (b) respectively. We compute the reduced density matrix corresponding to either of the spins A and B from the relation

[TABLE]

where Tr(.) represents partial tracing with respect to spin B or A. In Fig. 8 we have shown the structure of $\hat{\rho}_{A}^{m}$ after Fibonacci step $m=50$ which, we have checked, is sufficient to obtain the steady states. For small values of $J$, $\hat{\rho}_{A}^{m}$ contains both the diagonal and off-diagonal entries in the eigenbasis of $\hat{S}_{z}$ indexed by $m_{1}$, $m_{2}$ which can be observed in Fig. 8 (b). On the other hand, for large values of $J$, $\hat{\rho}_{A}^{m}$ becomes completely diagonal as is evident from Fig. 8 (c), with equally weighted entries, i.e., $\rho_{A}^{m}(m_{1},m_{1})\sim 1/(2S+1)$. Such an observation indicates that in the presence of interactions the coherent state picture is lost and leads us to conclude that the system approaches a diagonal ensemble and thermalizes to infinite temperature sray ; fazio ; nandkishore .

IV Conclusion

To summarize, we have studied the dynamics of a spin-$S$ object which has a well-defined classical limit and is subjected to quasiperiodic kicking following the Fibonacci sequence. By evolving the corresponding classical Hamiltonian map of the spin variables we obtained the phase portraits which, for small kicking strength $\lambda$ and time period $T_{0}$, exhibits regular orbits which precess over an unit sphere. Interestingly, for increasing $\lambda$ and $T_{0}$ the dynamics appears to be chaotic; however, the Lyapunov exponent vanishes. We have calculated the Sutherland invariant which constrains to some extent the dynamics governed by the transfer matrix with $SO(3)$ symmetry. Fluctuations of the classical counterpart of the spin dynamics exhibit a fractal structure which is verified by its Fourier spectrum analysis. It turns out that for an initially chosen spin coherent state, the phase coherence is retained during the time evolution under Fibonacci drive even for large $\lambda$ and $T_{0}$ indicating q classical-quantum correspondence. Fractality in classical dynamics is observed from the spectral analysis. More interestingly, the fractality is also present in the internal structure of the Floquet matrix governing the full quantum dynamics which has been investigated from the scaling of the Rényi entropy, as well as from the moments of the Floquet eigenvectors and the quasienergy spectrum. Finally, we have considered two such spin-$S$ objects interacting with each other and driven quasiperiodically. We have shown that in the presence of the interaction, the fractal behavior vanishes and level repulsion sets in the Floquet quasienergy spectrum. In the dynamics, we observed that for increasing interaction strength $J$, the average values of the components of spin operators for both the spins saturate to zero and the reduced density matrix for either of the spins becomes diagonal. The emergence of the diagonal ensemble with equally weighted diagonal elements indicates thermalization of the system to infinite temperature and corresponds to the microcanonical ensemble of statistical mechanics.

In conclusion, the existence of critical Floquet eigenstates exhibiting fractality in a quasiperiodically driven non-interacting spin model and its crossover to ergodic dynamics in the presence of interactions constitutes the central result of our work. Kicked spin models have already been realized in experiments considering the angular momentum of an atom in a suitable hyperfine state q_state_exp ; q_state_exp1 ; the kicking can be generated by using a short magnetic pulse k_top . The dynamics of such kicked systems can also be investigated in circuit QED experiments Cqed . The signature of fractality and its change to ergodic behavior in the dynamics can be found in the experiments by measuring the discrete time Fourier amplitude spectrum of time varying physical observable such as the spin variables in our study as well from the dimension measurement SNA_exp . The models discussed in our paper can thus be realized experimentally and our results can be tested in similar experiments.

Acknowledgments

S.R. thanks Anandamohan Ghosh for useful discussions. D.S. thanks DST, India for Project No. SR/S2/JCB-44/2010 for financial support. S.R. acknowledges the financial support and hospitality provided by IISER Kolkata for this work.