Time crystals: a review

Formation of space crystals relies on periodic self-organization of atoms due to their mutual interactions. Under certain conditions atoms arrange themselves in a periodic lattice that manifests itself in a periodic behavior in space of a probability density for a measurement of a single particle (an electron or an ion). Strictly speaking such a state cannot be the ground state of a many-body atomic system because it breaks the translation symmetry (assuming a typical nondegenerate ground state). To see this, let us consider a solid state Hamiltonian,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\tvec}[1]{\bf{#1}} \newcommand{\la}{\langle} \displaystyle H=\sum_{i=1}^N \frac{\tvec p_i^2}{2m_i}+\frac12\sum_{i\ne j}^NU_{ij}(\tvec r_i-\tvec r_j) \label{solid_state_h}, \nonumber \end{align} \tag{ 1 }$$

that describes N interacting particles in a finite volume V with periodic boundary conditions. If we shift all positions $\newcommand{\tvec}[1]{\bf{#1}} \tvec r_i$ by the same vector $\newcommand{\tvec}[1]{\bf{#1}} \tvec R$, the Hamiltonian does not change because it depends on relative distances between particles only. It means that the system possesses continuous (i.e. vector $\newcommand{\tvec}[1]{\bf{#1}} \tvec R$ can be arbitrary) space translation symmetry; the corresponding (unitary) space translation operator $\newcommand{\tvec}[1]{\bf{#1}} {{\mathcal T}}_{\tvec R}$ commutes with H and eigenstates $\newcommand{\tvec}[1]{\bf{#1}} \psi_n(\tvec r_1, \dots, \tvec r_N)$ of H are also eigenstates of $\newcommand{\tvec}[1]{\bf{#1}} {{\mathcal T}}_{\tvec R}$,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\tvec}[1]{\bf{#1}} \newcommand{\la}{\langle} \displaystyle {{\mathcal T}}_{\tvec R}\psi_n(\tvec r_1,\dots,\tvec r_N) = {\rm e}^{{\rm i}\varphi}\psi_n(\tvec r_1,\dots,\tvec r_N). \label{ta_eigenval} \nonumber \end{align} \tag{ 2 }$$

Taking into account equation (2) and the fact that

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\tvec}[1]{\bf{#1}} \displaystyle {{\mathcal T}}_{\tvec R}\psi_n(\tvec r_1,\dots,\tvec r_N) = \psi_n(\tvec r_1+\tvec R,\dots,\tvec r_N+\tvec R), \nonumber \end{align} \tag{ 3 }$$

it is easy to show that the probability density for the measurement of a single particle must be uniform in space if a system is prepared in the ground state (or any other nondegenerate eigenstate),

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\tvec}[1]{\bf{#1}} \displaystyle \rho(\tvec r_1+\tvec R)&=\int {\rm d}(\tvec r_2+\tvec R)\dots {\rm d}(\tvec r_N+\tvec R) \nonumber \\ &\quad \times \left\vert \psi_n(\tvec r_1+\tvec R,\dots,\tvec r_N+\tvec R)\right\vert ^2 \nonumber \\ = \int &{\rm d}\tvec r_2 \dots {\rm d}\tvec r_N\left\vert \psi_n(\tvec r_1,\dots,\tvec r_N)\right\vert ^2 = \rho(\tvec r_1),\nonumber \end{align} \tag{ 4 }$$

and no discrete structure is visible. However, crystalline properties can be observed in the two-point correlation function,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\tvec}[1]{\bf{#1}} \displaystyle \rho_2(\tvec r_1,\tvec r_2)=\int {\rm d}\tvec r_3\dots {\rm d}\tvec r_N\left\vert \psi_n(\tvec r_1,\tvec r_2,\tvec r_3\dots,\tvec r_N)\right\vert ^2, \nonumber \end{align} \tag{ 5 }$$

because the symmetry does not forbid $\newcommand{\tvec}[1]{\bf{#1}} \rho_2(\tvec r_1, \tvec r_2)$ to be non-uniform in space. If, for a fixed $\newcommand{\tvec}[1]{\bf{#1}} \tvec r_1$, the two-point correlation function reveals periodic behavior as a function of $\newcommand{\tvec}[1]{\bf{#1}} \tvec r_2$, then spontaneous breaking of the continuous space translation symmetry to a discrete translation symmetry can be predicted. The correlation function $\newcommand{\tvec}[1]{\bf{#1}} \rho_2(\tvec r_1, \tvec r_2)$ tells us what is the probability density for a measurement of the next particle provided the first particle has been detected at $\newcommand{\tvec}[1]{\bf{#1}} \tvec r_1$. Thus, it is enough to measure a single particle in order to see if a crystalline structure emerges. The measurement can be intentional, i.e. performed by apparatus in a laboratory, or simply due to the coupling of a system to its environment and the resulting possible particle losses.

The breaking of the space translation symmetry is related to the localization of the centre of mass of a system (van Wezel and van den Brink 2007). In the centre of mass coordinate frame it is apparent that the centre of mass degree of freedom decouples from the relative positions' degrees of freedom, i.e. the Hamiltonian (1) can be written as

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\tvec}[1]{\bf{#1}} \displaystyle H=\frac{\tvec P^2}{2mN}+{\rm relative~degrees~of~freedom.} \nonumber \end{align} \tag{ 6 }$$

The ground state of a system corresponds to the total momentum $\newcommand{\tvec}[1]{\bf{#1}} \tvec P=0$ and is entirely delocalized in the configuration space. Measurement of particle positions leads to the localization of the system and to an emergence of crystalline structures. It should be stressed that if one performed the measurement not of individual positions of particles but rather the relative distances between them, such a measurement would also reveal crystalline properties of the system without breaking the continuous space translation symmetry.

When the continuous space translation symmetry is broken to a discrete symmetry one can ask about the lifetime of the symmetry broken state because it is no longer the system eigenstate. In the thermodynamic limit defined as $N\rightarrow\infty$, the volume $V\rightarrow\infty$ but the particle density $N/V={\rm const.}$, the energy of the symmetry broken state is infinitely close to the ground state energy and its lifetime exceeds easily thousands of years. It can be estimated assuming that the center of mass is described by a wave-packet localized on a length scale $\sigma\approx10^{-11}$m, then, the corresponding kinetic energy is $E_k=\frac{h^2}{2mN\sigma^2}$. If such a system is isolated, the quantum spreading of the wave-packet leads to a delocalization of a space crystal. However, in order see that the crystalline structure is blurred, the center of mass must be delocalized on a length scale of the order of the lattice constant of a crystal, i.e. $a\approx 10^{-10}$ m, that takes time of the order of $t=a/(h/mN\sigma)\approx 5\cdot 10^4$ years for $mN=1$ kg. Fortunately, jewelery of our grandmothers is not isolated from the outside world and 'Diamonds are forever' (due to the measurement process).

Kinetic energy E_k is extremely close to the ground state energy and about $10^{-26}$ smaller than energy of an optical photon in our example. It implies that breaking of the symmetry practically costs no energy and can be induced by an infinitesimally small perturbation and, therefore, it occurs spontaneously.

An alternative method to predict the spontaneous breaking of a space translation symmetry is to apply a symmetry breaking perturbation, calculate the ground state of a system for finite N, take the thermodynamic limit and finally turn off the perturbation which leads to the symmetry broken state (Anderson 1997, Kaplan et al 1989, Koma and Tasaki 1993).

If a time-independent many-body system is prepared in an eigenstate $\psi_n$ corresponding to energy eigenvalue E_n, the probability density for detection of particles at a fixed position in the configuration space obviously does not change in time. It is a direct consequence of continuous time translation symmetry of a system. Indeed, time-independent Hamiltonian H commutes with the time translation operator which is simply the evolution operator ${{\mathcal T}}_t={\rm e}^{-{\rm i}Ht}$ (we assume $\hbar=1$) and eigenstates of H are also eigenstates of ${{\mathcal T}}_t$, thus,

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \vert \psi_n(t)\vert ^2=\vert {{\mathcal T}}_t\psi_n(0)\vert ^2=\vert {\rm e}^{-{\rm i}E_nt}\psi_n(0)\vert ^2=\vert \psi_n(0)\vert ^2. \nonumber \end{align} \tag{ 7 }$$

In 2012 Frank Wilczek proposed (Wilczek 2012) that time translation symmetry can be spontaneously broken in an analogue way to space translation symmetry breaking in the formation of space crystals. He coined the term time crystal for that phenomenon. If it exists then a time-independent many-body system, prepared in the ground state, can switch to a periodic motion in time under an infinitesimally weak perturbation. Experimentally it could be observed as a periodic behaviour in time of the probability density for a measurement of a system at a fixed point of the configuration space. It means that switching from space to time crystals we have to exchange the role of space and time. In the space crystal case we expect periodic behaviour in space at a fixed instant of time (i.e. at the moment when we perform a measurement of a system) while in the time crystal case we fix the position in space and ask whether a detector clicks periodically in time. More formal definitions of time crystals may be formulated, see e.g. Else et al (2016) and Khemani et al (2016a)—we shall stick to this intuitive one.

The idea of time crystals was proposed in two variants. Shapere and Wilczek (2012b) showed that a classical system can reveal periodic motion in the lowest energy state, see also (Ghosh 2014), while Wilczek himself (Wilczek 2012) presented an idea of quantum time crystal. In the classical case if we ask whether a time-independent system can move and at the same time possess the lowest energy, the answer seems to be obviously No! Indeed, in order to find the lowest energy for a classical particle we have to find an extremal value of a Hamiltonian,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \frac{\partial H}{\partial p}=0, \quad \frac{\partial H}{\partial x}=0, \label{classical_c} \nonumber \end{align} \tag{ 8 }$$

but that means no motion of a particle is possible because the first condition in (8) implies that the Hamilton equation,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \dot x=\frac{{\rm d}x}{{\rm d}t}=\frac{\partial H}{\partial p}=0. \label{hamilton_e} \nonumber \end{align} \tag{ 9 }$$

However, if we assume the energy of a particle of the form

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle E=\frac{\dot x^4}{4}-\frac{\dot x^2}{2}, \label{energy_c} \nonumber \end{align} \tag{ 10 }$$

one sees that the lowest energy corresponds to particle motion with velocity $\dot x=\pm 1$. This apparent contradiction with the conclusion based on the Hamilton equation (9) can be resolved when we realize that the energy (10) cannot be converted to the Hamiltonian smoothly. That is, the Hamiltonian is a multi-valued function of the momentum with cusps corresponding precisely to energy minima at $\dot x=\pm 1$ where the Hamilton equations are not defined (Shapere and Wilczek 2012b).

Encouraged by the classical analysis we can now look for a quantum time-independent many-body system that in the ground state can spontaneously switch to periodic motion and reveal crystalline properties in the time domain. Shapere and Wilczek (2012a) proposed how to quantize the classical single-particle system described by energy (10). However, in the following we will concentrate on a more general problem of many-body systems with more conventional kinetic energy terms that could reveal periodic motion in measurements repeated many times on the same realization of a system. That is, we would like to consider a situation where a detector placed at a certain point in the configuration space reacts periodically because particles are returning to this point.

A potential example fulfilling our requirements seems to be a superconducting device where an external magnetic field induces current of Cooper pairs. However, the flow of Cooper pairs is uniform and if, at a certain position in space, we detect particles at a certain moment of time, the next detection event will not be correlated temporally with the previous one and no periodic crystalline structure in time will be observed (Wilczek 2012, Yamamoto 2015). The original idea of Wilczek was more involved. It is known that interacting particles can form spontaneously inhomogeneous structures in space. On the other hand a charged particle on a ring (one-dimensional (1D) problem with periodic boundary conditions) subjected to magnetic flux can reveal non-vanishing probability current along the ring in the ground state if the flux is properly chosen. Combining these two observations, it should be possible to observe a spontaneous process where a many-body system prepared in the ground state switches to periodic evolution where an inhomogeneous particles density moves around a ring. Wilczek concentrated on bosons interacting via attractive contact potential on an Aharonov-Bohm ring and we will elaborate on this system in a moment. Very soon another proposition (Li et al 2012b) was suggested that ions on a ring may spontaneously form a space crystal when kinetic energy of ions is much smaller than Coulomb potential energy between them. Such a Wigner crystal, in the presence of a magnetic flux, reveals periodic motion along the ring even if the system is initially prepared in the ground state. Both proposals were immediately criticized by Bruno (2013a and 2013b) who pointed out that such scenario is impossible, for replies see Li et al (2012a) and Wilczek (2013b). Soon after Patrick Bruno showed that under quite general conditions spontaneous breaking of continuous time translation symmetry is not possible in any time-independent system prepared in the ground state (Bruno 2013c). Before we present the arguments of Bruno and other researchers we first describe details of Wilczek idea (Wilczek 2012).

Let us consider N bosons on a ring of unit length with attractive contact Dirac δ-interactions and in the presence of a magnetic flux α (a ring problem in the presence of a magnetic flux is dubbed Aharonov–Bohm ring). The particle mass and $\hbar$ are assumed to be equal to unity and the parameter g₀ that determines the strength of attractive interactions is negative. The Hamiltonian of the system,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H=\sum_{i=1}^N\frac{(\,p_i-\alpha)^2}{2}+\frac{g_0}{2}\sum_{i\ne j}\delta(x_i-x_j), \label{h} \nonumber \end{align} \tag{ 11 }$$

possesses continuous time and space translation symmetries and consequently probability density corresponding to any eigenstate is invariant under any translation in time and any translation of all particles along the ring.

Let us assume for a moment that the magnetic flux $\alpha=0$ and let us apply the mean field approximation. In the mean field approach all bosons are supposed to occupy the same single particle state ϕ, i.e. form a Bose–Einstein condensate. Many-body eigenstates are then of the form of a product state $\phi(x_1)\phi(x_2)\dots\phi(x_N)$. In order to obtain the ground state within the mean field approximation one has to find the minimal value of the energy of the system within the Hilbert subspace spanned by product states that reduces to the solution of the Gross–Pitaevskii equation (Pethick and Smith 2002),

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \left(-\frac12 \partial_x^2+g_0(N-1)\vert \phi\vert ^2\right)\phi=\mu\phi. \label{gpe} \nonumber \end{align} \tag{ 12 }$$

If the attractive particle interactions are sufficiently strong, i.e. $g_0(N-1)<-\pi^2$, it is energetically favorable to group particles together and the mean field solution ϕ breaks the space translation symmetry and becomes inhomogeneous in space (Carr et al 2000). The ϕ function is given by the Jacobi elliptic function but when $\vert g_0\vert N\gg 1$, it is well approximated by a bright soliton solution

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \phi(x)\approx\cosh^{-1}\left[\frac{g_0(N-1)}{2}(x-x_{\rm CM})\right], \label{bright_s} \nonumber \end{align} \tag{ 13 }$$

where x_CM is a parameter that describes the centre of mass position (Carr et al 2000, Pethick and Smith 2002). Thus, the mean field approach predicts that bosons form a Bose–Einstein condensate where all particles occupy a localized wavefunction (13). When we return to the many-body description we obtain that the many-body ground state has the space translation symmetry but this state is strongly vulnerable to any perturbation. In fact to break the space translation invariance it is enough to measure a position of a single particle (Delande et al 2013).

Wilczek expected that in the presence of a properly chosen magnetic flux α, one would not only observe spontaneous localization of density of particles in the configuration space but such a density will also move periodically along the ring and the motion sustains forever in the limit when $N\rightarrow\infty$, $g_0\rightarrow 0$ but $g_0(N-1)={\rm const.}$, i.e. in the limit where the mean field prediction remains unchanged, see (13). This is actually false as immediately pointed out by Bruno (2013a). Probably, the simplest way to demonstrate it, is to switch to the centre of mass coordinate frame (Syrwid et al 2017). Then, the Hamiltonian (11) reads

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H=\frac{(P-N\alpha)^2}{2N}+{\rm relative~degrees~of~freedom}, \label{h_cm} \nonumber \end{align} \tag{ 14 }$$

where P is the centre of mass momentum, i.e. the total momentum of the system, which is a conserved quantity. The centre of mass and the relative positions' degrees of freedom decouple and eigenstates of the system are determined by an independent choice of the centre of mass momentum $P_j=2\pi j$ (where j is integer) and the relative degrees of freedom quantum numbers. The ground state corresponds to

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \frac{\partial H}{\partial P_j}=2\pi\frac{j}{N}-\alpha\approx 0. \label{zero_c} \nonumber \end{align} \tag{ 15 }$$

In the limit when $N\rightarrow\infty$, equation (15) can be fulfilled exactly. This is very bad news because it means there is no probability current related to the centre of mass degree of freedom if the system is prepared in the ground state. Wilczek idea relied on the assumption that if the flux α is chosen properly, quantized values of particle momenta do not allow equation (15) to vanish and once the bright soliton is formed in the process of spontaneous breaking of space translation symmetry it will move too. We see it is not the case in the crucial limit, i.e. where the total number of particles increases. The idea of Li and co-workers was similar (Li et al 2012b). Instead of particles interacting via an attractive contact potential, they envision ions on an Aharanov–Bohm ring. A spontaneous breaking of space translation symmetry results in a formation of a space crystal which has been supposed to move in the presence of a magnetic flux even if an experiment started with the ground state. Again, the same line of arguments leads to the conclusion that it is not possible.

The impossibility of realization of the time crystal idea was also considered by Nozières (2013) who investigated a superfluid ring in a magnetic field and presented arguments that a charge density wave cannot reveal rotation induced by diamagnetic currents. However, general analysis of the impossibility of spontaneous time translation symmetry breaking in the ground state was performed by Bruno (2013c) and later by Watanabe and Oshikawa (2015). Bruno considered a general many-body system on an Aharonov-Bohm ring that is subjected to a perturbation that rotates periodically along the ring. He showed that for the ground state that breaks rotational symmetry, the moment of inertia of the system is always positive and imposing rotation increases energy. Bruno also considered the thermal equilibrium state in the rotating frame, where the Hamiltonian becomes time-independent, and reached a similar conclusion. However, for a periodically driven system analysis of thermal equilibrium in a rotating frame should be performed with caution because it actually implies that a system is assumed to be in contact with a reservoir of rotating particles. Moreover, it was shown that spontaneous emission of photons may correspond to a jump of an electron in an atom upwards in energy if the process is described in the rotating frame (Delande and Zakrzewski 1998). This obstacle was overcome by Watanabe and Oshikawa (2015) who did not assume the presence of a symmetry broken perturbation but focused on analysis of a correlation function. They proved that the two-point correlation function does not reveal any time dependence when a many-body system is prepared in the ground state or thermal equilibrium state in the limit when volume V of a system goes to infinity.

In the case of space crystals, spontaneous space translation symmetry breaking is demonstrated in the limit when a number of particles N and V tend to infinity but the density of particles is constant. Then, any perturbation is sufficient to break the symmetry and a symmetry broken state lives forever. In the case of time crystals, it is not necessary to take $V\rightarrow\infty$ because we do not have to assume periodic or any other behaviour of a system in space. We may set $V={\rm const.}$, increase N but keep the product of the coupling constant g₀ and N fixed, see equation (11).

Evolution of the phase of a Bose–Einstein condensate described within ground canonical ensemble can be considered in the context of time crystals. In the ground canonical formalism, the order parameter of a Bose–Einstein condensate reveals periodic oscillations $\newcommand{\ra}{\rangle} \newcommand{\la}{\langle} \newcommand{\tvec}[1]{\bf{#1}} \la\hat\psi(\tvec r, t)\ra=\psi_0{\rm e}^{-{\rm i}\mu t}$, where $\hat \psi$ is a bosonic field operator and μ is a chemical potential (Castin and Dum 1998, Pethick and Smith 2002). Such a periodic time evolution can be measured provided the condensate is coupled to another condensate. If a system is strictly isolated, i.e. when a number of particles is conserved, there is no reference frame and no time dependence can be detected. Volovik (2013) analyzed a class of such systems (see also Nicolis and Piazza (2012), Wilczek (2013a), Castillo et al (2014) and Thies (2014)) and introduced two relaxation times: energy relaxation time $\tau_E$ and relaxation time $\tau_N$ related to a total number of particles in our example. If $\tau_N\gg\tau_E$, a system with an average number of particles N relatively quickly relaxes to a minimal energy state at fixed N and then slowly relaxes to the equilibrium state where no oscillations are present. Although in the intermediate time $\tau_N\gg t\gg\tau_E$ one can observe breaking of time translation symmetry, a system is not strictly in the equilibrium state. While $\tau_N$ can tend to infinity, the strict limit $\tau_N=\infty$ cannot be taken if one wants to observe the oscillation because then there is no reference frame to measure them.

Similar class of systems is considered by Else et al (2017), however, the authors assume isolated systems. U(1) symmetries that are spontaneously broken are not exact symmetries but they are related to effective Hamiltonians. Spontaneous breaking of U(1) symmetry results in an order parameter that, according to an effective Hamiltonian, oscillates forever. However, at infinitely long times a system approaches a thermal state determined by a full Hamiltonian where no oscillations are observed.

The results presented indicate that in a many-body system spontaneous breaking of time translation symmetry cannot be observed if a system is prepared strictly in the equilibrium state. Wilczek idea in its original version turned out to be impossible for realization but it became an inspiration to other researchers. Particular progress has been achieved for systems that break discrete time translation symmetry—the problem described in section 3. Before considering it let us consider the possibility to utilize excited states of many-body systems.

Let us analyze if the spontaneous breaking of time translation symmetry can be observed for a many-body system prepared in an excited eigenstate. In other words we would like to answer the question if a time independent many-body system prepared in an excited eigenstate is able to self-organize in time and switch to periodic motion under infinitesimally weak perturbation. We additionally require that such a situation should not be a theoretical issue only but should be realizable experimentally.

Analysis of Wilczek model, i.e. the system described by the Hamiltonian (11), Syrwid et al (2017) suggests that indeed the spontaneous breaking of continues time translation symmetry can occur for an excited eigenstate. For any chosen value of the magnetic flux α, the ground state of the system corresponds to such a total momentum P_j that the probability current associated with the centre of mass motion vanishes if $N\rightarrow\infty$, see equation (15). However, if the system is prepared in an eigenstate with total momentum $P_N=2\pi N$, the probability current related to the center of mass motion,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \frac{\partial H}{\partial P_N}=2\pi-\alpha, \label{current} \nonumber \end{align} \tag{ 16 }$$

does not vanish provided the flux $\alpha\ne 2\pi$. Thus, for $\alpha\ne 2\pi$, let us choose among all eigenstates with the total momentum equal P_N the eigenstate $\newcommand{\ra}{\rangle} \vert \psi_0\ra$ with the lowest energy. Although such an eigenstate is not the ground state, it realizes Wilczek's idea. That is, the probability density related to this eigenstate is invariant under time translation transformation but if attractive particle interactions are sufficiently strong, in the limit when $N\rightarrow\infty$, any perturbation can make the density inhomogeneous. The inhomogeneous density should rotate around the ring with the period $T=(2\pi-\alpha){\hspace{0pt}}^{-1}$ that is determined by equation (16) (remember that the ring has a unit length).

Numerical simulations (Syrwid et al 2017) have confirmed that the density-density correlation function,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \newcommand{\ra}{\rangle} \displaystyle \label{rho2} \rho_2(x,t)\propto\langle\psi_0\vert \hat\psi^\dagger(x,t)\hat\psi(x,t)\;\hat\psi^\dagger(x_1,0)\hat\psi(x_1,0)\vert \psi_0\rangle, \nonumber \end{align} \tag{ 17 }$$

is inhomogeneous as a function of x and reveals periodic rotation around the ring with the lifetime t_c that increases linearly with N if $g_0\rightarrow 0$ but $g_0(N-1)=\rm constant$, i.e. in the limit when the mean field prediction remains unchanged, see figure 1. The density-density correlation function (17) corresponds actually to the probability density for the measurement of a particle at space point x and at time t provided that at $t=0$ another particle was already detected at x₁. It illustrates the nature of a spontaneous process: in order to learn whether the symmetry is broken one has to perform measurement and a minimal possible information, i.e. information about the position of a single particle, is sufficient to break the symmetry.

Zoom In Zoom Out Reset image size

Measurement of the position of a single particle results in a localization of the centre of mass of the system on a length scale of the order of the bright soliton size $\sigma_{\rm CM}(0)\approx [g_0(N-1)]^{-1}=\rm const.$, see equation (13). The localized centre of mass probability density spreads, $\sigma_{\rm CM}(t)\propto t/N$, in the course of time evolution determined by a free particle like Hamiltonian (14). Thus, time moment t_c when $\sigma_{\rm CM}(t_c)$ becomes comparable to the length of the ring scales linearly with N. That agrees with the numerical results. If more particles are measured initially, the centre of mass probability density can be localized on a length scale $\sigma_{\rm CM}(0)\propto N^{-1/2}$. In such a case, the lifetime of periodic evolution of the symmetry broken state scales like $N^{1/2}$ (Syrwid et al 2017).

It is not simple to prepare in a laboratory a many-body system in a specific excited eigenstate. However, ultra-cold atomic gases constitutes a perfect playground for realization, control and detection not only many-body ground states but also collectively excited states (Pethick and Smith 2002). The described spontaneous breaking of continuous time translation symmetry in excited states (Syrwid et al 2017) can be observed in ultra-cold atoms trapped in a toroidal potential in the presence of an artificial gauge field (Goldman et al 2014). The latter can be realized, e.g. by imposing rotation of a thermal cloud during the evaporative cooling because then the system dissipates to the lowest energy state in the rotating frame and the Coriolis force mimics a magnetic field (Pethick and Smith 2002). Proper choice of an artificial gauge potential allows one to obtain the flux $\alpha=2\pi$ for which the ground state is the desire state $\newcommand{\ra}{\rangle} \vert \psi_0\ra$ with the total momentum $P_N=2\pi N$. Then, turning off the gauge potential leaves the system in the state $\newcommand{\ra}{\rangle} \vert \psi_0\ra$ which is no longer the ground state of the Hamiltonian with $\alpha=0$ and breaking of time translation symmetry under any weak perturbation is expected. This observation will be in contrast with the same experiment but performed for $g_0(N-1)>-\pi^2$ where no spontaneous breaking of space or time translation symmetry occurs.

To summarize this section, any time-independent many-body system has the continuous time translation symmetry. It implies that probability density corresponding to any eigenstate does not change in time. Moreover, no spontaneous breaking of time translation symmetry can be observed if a system is in the ground state or the thermal equilibrium state. However, it was illustrated with the help of Wilczek's model that the continuous time translation symmetry can be spontaneously broken if a system is prepared in an excited eigenstate. Such a phenomenon could be in principle realized in ultra-cold atoms laboratory. The idea of time crystals initiated a novel field of research and inspired many physicists (Chernodub 2013, Mendonca and Dodonov 2014, Robicheaux and Niffenegger 2015, Yamamoto 2015, Faizal et al 2016).

The first example of the spontaneous breaking of discrete time translation symmetry to another discrete symmetry was given in 2015 (Sacha 2015b). As an illustration of the effect a system of ultra-cold atoms bouncing on an oscillating mirror in the presence of the gravitational field is considered. In order to explain the phenomenon we have to first describe the corresponding single particle problem, a classical version of which, called a bouncer has been introduced by Pustyl'nikov as a model for Fermi acceleration (Zaslavsky 1970, Pustyl'nikov 1978). This famous model of classical chaos has been studied experimentally (see e.g. Pierański (1983)), its quantum version has been often studied (Dembiński et al 1993, Buchleitner et al 2002), moreover, it was realized experimentally for cold atoms bouncing on a mirror formed by an evanescent wave (Steane et al 1995).

A single particle bouncing on an oscillating mirror is described (in 1D approximation and in the frame oscillating with a mirror and in the so-called gravitational units) by the following quantum Hamiltonian (Buchleitner et al 2002),

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H=-\frac12 \partial_z^2+V(z)+\lambda z\cos(\omega t), \label{gpe_s15} \nonumber \end{align} \tag{ 21 }$$

where $V(z)=z$ for $z\geqslant 0$ and $V(z)=\infty$ for $z<0$. In the frame oscillating with a mirror, the position of a mirror is fixed (at $z=0$) but the gravitational constant depends periodically on time. In the absence of the mirror oscillations ($\newcommand{\la}{\langle} \lambda=0$) all classical trajectories of a particle are periodic with a period increasing with the energy of the particle. When mirror oscillations are turned on, classical motion becomes irregular but some of periodic orbits survive. They are stable resonant orbits living in regular parts of the classical phase space. There is 1:1 resonant orbit where a particle moves periodically with a period equal to the period of the mirror oscillations $T=2\pi/\omega$. There exist 2:1 resonant orbit where a particle bounces on a mirror with a period twice longer than that of the oscillations as well as higher order resonances $s:1$. Switching to the quantum description, a motion of a particle is described by Floquet eigenstates. It may be surprising but classical-like motion of a particle on resonant orbits can be also observed in the quantum world. For example, suitable choice of parameters results in a Floquet eigenstate that is represented by a localized wave-packet moving periodically along 1:1 resonant orbit. This is an example of the so-called non-spreading wave-packet motion that was discovered more than 20 years ago (Henkel and Holthaus 1992, Bialynicki-Birula et al 1994, Delande and Buchleitner 1994, Buchleitner and Delande 1995, Holthaus 1995, Zakrzewski et al 1995) for a review see Buchleitner et al (2002).

We will focus on the 2:1 resonance case (Sacha 2015b). A single localized wave-packet moving along classical 2:1 resonant orbit cannot form a Floquet eigenstate because it moves with a period twice longer than T. However, a superposition of two such wave-packets that move with the period 2T but after T exchange their positions can form a proper Floquet eigenstate and indeed there exists such a state, see figure 2. Two localized wave-packets can actually form two mutually orthogonal superpositions. Therefore, there exist two Floquet eigenstates, $u_1(z, t)$ and $u_2(z, t)$, consisting of two localized wave-packets moving along the 2:1 resonant trajectory. The eigenstates $u_1(z, t)$ and $u_2(z, t)$ correspond to nearly degenerate (modulo $\omega/2$) quasi-energies $E_2-E_1=\omega/2+J$. A tiny splitting J is related to a tunneling process. By superposing $u_1(z, t)$ and ${\rm e}^{-{\rm i}\omega t/2}u_2(z, t)$ we can eliminate one of the two localized wave-packets. The remaining wave-packet evolves along 2:1 resonant trajectory but at the same time slowly tunnels to a position of the other missing wave-packet—the full tunneling process is completed after a period $\pi/J$ which is very long as compared to T. It is worth noting that in order to observe the described resonant behaviour, the resonant condition does not need to be strictly fulfilled. That is, if the resonance condition corresponds to an unperturbed classical orbit with a sufficiently high energy, in the quantum description the nearest in energy unperturbed quantum states will form the Floquet eigenstates that we look for. In other words small changes of the driving frequency ω do not change the behaviour and no fine tuning is necessary because the motion of localized wave-packets is protected by local constants of motion related to regular parts of the classical phase space (Buchleitner et al 2002).

Zoom In Zoom Out Reset image size

When we consider N particles bouncing on an oscillating mirror we observe a similar behaviour if particles are bosons and they do not interact. That is, in the Hilbert subspace spanned by Fock states $\newcommand{\ra}{\rangle} \vert n_1, n_2\ra$, where n₁ and $n_2=N-n_1$ are numbers of particles occupying single particle Floquet states u₁ and ${\rm e}^{-{\rm i}\omega t/2}u_2$, respectively, the lowest quasi-energy eigenstate corresponds to $\newcommand{\ra}{\rangle} \vert N, 0\ra$. In order to describe particles interacting via δ-contact potential one has to consider the following many-body Floquet Hamiltonian

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H_F&=\int\nolimits_0^{2T}{\rm d}t\int\nolimits_0^\infty {\rm d}z\;\hat\psi\left[H(t)+\frac{g_0}{2}\hat\psi^\dagger\hat\psi-{\rm i}\partial_t\right]\hat\psi, \nonumber \\ &\approx -\frac{J}{2}\left(\hat c_1^\dagger \hat c_1-\hat c_2^\dagger \hat c_2\right)+\frac14 (U-2U_{12})\nonumber \\ &\quad\times \left[(\hat c_1^\dagger)^2\hat c_2^2+(\hat c_2^\dagger)^2\hat c_1^2+2\hat c_1^\dagger \hat c_1\hat c_2^\dagger\hat c_2\right]+{\rm const.}, \label{manyH_s15} \nonumber \end{align} \tag{ 22 }$$

where U and $U_{12}$ are equal to integrals over time and space of products of the probability densities of two evolving localized wave-packets multiplied by g₀. In equation (22) we restricted to the Hilbert subspace spanned by the previously described Fock states $\newcommand{\ra}{\rangle} \vert n_1, n_2\ra$ and the bosonic field operator $\hat\psi(z, t)\approx u_1\hat c_1+{\rm e}^{-{\rm i}\omega t/2}u_2\hat c_2$ where $\hat c_1$ and $\hat c_2$ are standard annihilation operators. Such a restriction is valid provided the interaction energy is very small and couplings to the complementary Hilbert subspace can be neglected. It will be the case because the interaction energy per particle will be of the order of the tunneling splitting J that is extremely small as compared to any other energy scale of the system. Eigenstates of the Hamiltonian (22) correspond to many-body Floquet states that are also eigenstates of the time translation operator ${{\mathcal T}}_T$ what is apparent when one realizes that there are two classes of the eigenstates of (22). Eigenstates from the first class are spanned by Fock states with only even occupations of the ${\rm e}^{-{\rm i}\omega t/2}u_2(z, t)$ mode while eigenstates from the other class by Fock states with the odd occupations only (Sacha 2015b).

Assuming attractive ($g_0<0$) interactions that are very weak, i.e. $N\vert U-2U_{12}\vert <J$, the ground state $\newcommand{\ra}{\rangle} \vert \psi_0\ra$ of (22) matches the non-interacting result, $\newcommand{\ra}{\rangle} \vert \psi_0\ra\approx \vert N, 0\ra$. However, when

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle N\vert U-2U_{12}\vert >J, \label{cond_s15} \nonumber \end{align} \tag{ 23 }$$

it is energetically favorable to collect all bosons in a single localized wave-packet and consequently the ground state of the Hamiltonian (22) is a Schrödinger cat-like state which is clear if one writes such a many-body Floquet eigenstate in another Fock basis $\newcommand{\ra}{\rangle} \vert \tilde n_1, \tilde n_2\ra$ where $\tilde n_1$ and $\tilde n_2=N-\tilde n_1$ are occupations of the first localized wave-packet, $\phi_1=(u_1+{\rm e}^{-{\rm i}\omega t/2}u_2)/\sqrt{2}$, and the other localized wave-packet, $\phi_2=(u_1-{\rm e}^{-{\rm i}\omega t/2}u_2)/\sqrt{2}$, respectively. Then, the many-body ground state reads

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \newcommand{\ra}{\rangle} \displaystyle \vert \psi_0\ra\approx\frac{\vert N,0\ra+\vert 0,N\ra}{\sqrt{2}}. \label{cat_s15} \nonumber \end{align} \tag{ 24 }$$

The corresponding single particle probability density,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \newcommand{\ra}{\rangle} \displaystyle \rho_1(z,t)&=\langle\psi_0\vert \hat \psi^\dagger(z,t)\hat \psi(z,t)\vert \psi_0\rangle \nonumber \\ &\approx\frac{N}{2}\left(\vert \phi_1(z,t)\vert ^2+\vert \phi_2(z,t)\vert ^2\right), \label{single_p_s15} \nonumber \end{align} \tag{ 25 }$$

is plotted in figure 3 at different time moments. The discrete time translation symmetry is preserved in time evolution of the many-body Floquet eigenstate $\newcommand{\ra}{\rangle} \vert \psi_0\ra$ but this state is extremely vulnerable to any perturbation. After a measurement of a position x₁ of a single particle, the symmetry is gone because the quantum state of the remaining particles immediately collapses to one of the terms in the sum (24), i.e. to $\newcommand{\ra}{\rangle} \vert N-1, 0\ra$ or $\newcommand{\ra}{\rangle} \vert 0, N-1\ra$ depending on a result x₁ of the measurement. Then, time evolution shows that the original discrete time translation symmetry has been broken, i.e. the system evolves with the period 2T. The resulting state $\newcommand{\ra}{\rangle} \vert N-1, 0\ra$ (or $\newcommand{\ra}{\rangle} \vert 0, N-1\ra$) is robust against any further perturbation—one can perform many measurements and still the period of the time evolution remains 2T, see figure 4. Tunneling time from the state $\newcommand{\ra}{\rangle} \vert N-1, 0\ra$ to $\newcommand{\ra}{\rangle} \vert 0, N-1\ra$ or vice versa, in the limit $N\rightarrow\infty$, $g_0\rightarrow 0$ but $g_0N=\rm const.$, increases like $e^{\alpha N}/N$ with a positive constant α and very quickly becomes so long that it is not-measurable (Sacha 2015b).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Not only the ground and first excited states of (22) possess a Schrödinger cat-like structure. With $N\rightarrow\infty$ more and more eigenstates come in pairs of cat states (Ziń et al 2008, Oleś et al 2010) and experimental preparation of any initial state where most of atoms occupy a single localized wave-packet will result in time evolution with the period 2T that practically never ends. This is in contrast to the same experiment performed for $N\vert U-2U_{12}\vert <J$ where one will observe tunneling of atoms to another localized wave-packet after time period of the order of $1/J$.

The presented example constitutes an illustration that spontaneous self-organization in time of a periodically driven many-body system is possible, i.e. spontaneous breaking of discrete time translation symmetry to another discrete symmetry may occur. In the following section we shall discuss a similar behaviour observed for a very different model—a system of driven interacting spins prepared in such a way that during a single period T all the spins are flipped (being thus brought to the same orientation after 2T). Then, the corresponding Floquet eigenstates are macroscopic Schrödinger cat-like states of two possible orientations of spins. A measurement of the direction of even a single spin leads to a collapse of any of these cat-like states to a short range correlated time crystal state.

Interestingly, the original presentation of discrete time translation symmetry breaking (Sacha 2015b) was not originally noticed. It took more than a year to rediscover that phenomenon in quite a different setting (Else et al 2016, 2017, Khemani et al 2016b, von Keyserlingk and Sondhi 2016b) involving, importantly, the disorder induced effects.

Recent years brought an intensive study of many-body localization (MBL)—a phenomenon which seems to break a common wisdom about disordered many body systems. The latter were supposed to generically thermalize if evolved from some nonstationary initial state. Starting from a seminal work of Basko et al (2006) it seems it does not have to be the case, a sufficiently strong disorder may bring a many body system to a localized, non-ergodic phase (Oganesyan and Huse 2007, Žnidarič et al 2008). The phenomenon received a lot of attention with literally hundreds of publications in last 10 years (for excellent early reviews see Huse et al (2014) and Nandkishore and Huse (2015)).

While mainly observed for one-dimensional spin chains, MBL is considered by now to be a generic phenomenon in disordered many-body systems. The latter are often studied using strong periodic driving—being a natural current extension of strongly periodically driven single particle systems intensively studied in the last millennium (consider e.g. laser-atom interactions or microwave resonance techniques). Until quite recently it has been a common understanding that periodic driving of a many body system must supply (on average) energy to the system and lead to heating (see e.g. D'Alessio and Rigol (2014), Lazarides et al (2014a) and Lazarides et al (2014b)) for arbitrary initial states—see, however, Abanin et al (2015), Chandran and Sondhi (2016), Kuwahara et al (2016) and Mori et al (2016). Since any periodically driven system may be described by Floquet time-periodic eigenvectors (Shirley 1965, Sambe 1973), the inevitable heating suggests a lack of possibility to prepare an initial state in a form of a single, or few Floquet eigenstates. Interestingly, this believe is wrong for simple driven single particle problems where specific localized examples—the so-called nonspreading wavepackets (Buchleitner et al 2002) were mentioned already in the previous section. They have been also demonstrated experimentally (Maeda and Gallagher 2004). Let us note that lifetime of states that are superpositions of few Floquet eigenstates is infinite. These wavepackets were studied, however, in relatively simple small systems while the heating argument was usually presented in the thermodynamic limit.

Numerical experiments revealed that in the presence of disorder MBL persists for driven systems provided the frequency of the driving is high enough (Abanin et al 2015, Lazarides et al 2015, Ponte et al 2015) so the long-time system behavior may be analyzed by an effective Hamiltonian resulting from proper time averaging over the period of the drive. Soon it has been realized that periodically driven systems can not only remain many-body localized but that one can distinguish different 'phases' by means of appropriate correlators as discussed for spin systems by Khemani et al (2016b). Among those, in the present context, important is the phase named as π-spin glass. While the original analysis (Khemani et al 2016b) uses (Jordan-Wigner transformation based) link to earlier works on driven Majorana edge modes for noninteracting systems (Jiang et al 2011, Bastidas et al 2012, Thakurathi et al 2013) let us stay within the language of Schrödinger cat-like delocalized Floquet states as discussed in the previous section.

The system considered (Khemani et al 2016b) is a driven spin model with a binary periodic drive corresponding to the unitary evolution over the period $T=t_1+t_2$ given by $\newcommand{\e}{{\rm e}} U=\exp(-{\rm i}t_2H_x){\rm exp}(-{\rm i}t_1H_z)$ with

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \label{hamkhe} H_z&=\sum_ih_i\sigma_i^z+\sum_iJ_z\sigma_i^z\sigma_{i+1}^z \nonumber \\ H_x&=\sum_i \left(J_i\sigma_i^x\sigma_{i+1}^x +J_z\sigma_i^z\sigma_{i+1}^z\right). \nonumber \end{align} \tag{ 26 }$$

Aiming at discussing physics in the many-body localization regime, the authors analyse statistics of quasi-energies (eigenvalues of U) looking in particular at the average $\bar r$ ratio between the smallest and the largest adjacent energy gaps (Oganesyan and Huse 2007): $\newcommand{\e}{{\rm e}} r_n = \min[\delta^\epsilon_n, \delta^\epsilon_{n-1}]/ \max[\delta^\epsilon_n, \delta^\epsilon_{n-1}],$ with $\newcommand{\e}{{\rm e}} \delta^\epsilon_n= \epsilon_n - \epsilon_{n-1},$ and $\newcommand{\e}{{\rm e}} \epsilon_n\in [0, 2\pi)$ are the quasi-energies. In the MBL phase, one expects $\bar r$ to be close to the Poisson limit $\bar r^{\rm Poisson}=2\ln 2 -1\approx0.386$ (Atas et al 2013). That allows one to choose the proper values for J_z and suitable distributions for random values of $h_it_1$ and $J_it_2$ in (26).

Observe that H_x combined with H_z enjoys similar symmetries to Ising model. The Floquet quasi-energy eigenstates are also eigenstates of parity $P=\Pi_i \sigma^z_i$. Observe also that for $h_it_1=\pi/2$ the role of the first term in H_z will be to flip the spins in the x direction. Neglecting for a moment the other terms in U, its action in two periods of the drive brings the spins to the same orientation. One then expects that many-body localized (when disorder and interactions are included) Floquet eigenstates will resemble Schrödinger cat states with frozen domain walls deep in the spin glass phase

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\dow}{\downarrow} \newcommand{\ra}{\rangle} \displaystyle \vert \pm\ra= \frac{1}{\sqrt{2}}\left(\vert \uparrow\downarrow\uparrow\dots \uparrow\ra_x\pm\vert \downarrow\uparrow\downarrow\dots\downarrow\ra_x\right). \nonumber \end{align} \tag{ 27 }$$

They will come in pairs with different parity being almost degenerate in energy (modulo $\pi/T$) revealing the existence of $\pi/T$ phases resulting from discrete symmetry breaking (that destroys vulnerable cats). Figure 5 presents the disorder averaged spectral function of the spin raising operator, $\sigma_i^+$ on an arbitrary site i in the Floquet basis

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \newcommand{\ra}{\rangle} \displaystyle \label{eq:khe} {{\mathcal A}}(\omega)=\frac{1}{2^L}\sum_{\alpha\beta}\langle\phi_\alpha(0)\vert \sigma^+_i\vert \phi_\beta(0)\ra\;\delta(\omega-(\epsilon_\alpha-\epsilon_\beta)), \nonumber \end{align} \tag{ 28 }$$

where $\newcommand{\e}{{\rm e}} \epsilon_\alpha$ are Floquet eigenvalues with corresponding vectors $\newcommand{\ra}{\rangle} \vert \phi_\alpha\rangle$ with L being the system size. A peak at $-\pi/T$ is robust against changes of interaction strengths, disappearing only above a critical interaction strength. The effective period doubling observed may be interpreted as a potential Floquet realization of a time crystal (Khemani et al 2016b).

Zoom In Zoom Out Reset image size

The same group presented a nice general group analysis of symmetry broken phases in Floquet systems (von Keyserlingk et al 2016, von Keyserlingk and Sondhi 2016b) in a similar setting. This time the authors chose to flip the role of x and z axes considering the model defined by $\newcommand{\e}{{\rm e}} U=U_2U_1=\exp(-{\rm i}t_2H_2)\exp(-{\rm i}t_1H_1)$ with the period $T=t_1+t_2$ and

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \label{ham2} H_1&=\sum_ih_i\sigma_i^x\nonumber \\ H_2&=\sum_iJ_i\sigma_i^z\sigma_{i+1}^z, \nonumber \end{align} \tag{ 29 }$$

where h_i may be random while J_i are uniformly drawn from $[J_z-\delta J, J_z+\delta J]$ (to obtain many-body localization). This is a 'minimal' spin system with J_z controlling the interactions. Again for $h_i t_1=\pi/2$ the action of $\newcommand{\e}{{\rm e}} \exp(-{\rm i}t_1H_1)$ corresponds to flipping the spins along, this time, z direction, i.e. $U_1=\Pi_r{\rm i}\sigma^x_r$. The second unitary is diagonal in z components. Clearly the quasi-energy eigenstates will be now cat-like states

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \newcommand{\ra}{\rangle} \displaystyle \vert \pm\rangle\propto \vert \{d_i\}, p=\pm\rangle= \frac{1}{\sqrt{2}}\left(\vert \{m_i\}\rangle_z \pm \vert \{-m_i\}\rangle_z\right), \label{eqcat} \nonumber \end{align} \tag{ 30 }$$

where d_i is the expectation value of $\sigma_i^z\sigma_{i+1}^z$ and p is the Ising parity of the state. While the argument is presented for $h_it_1=\pi/2$ the numerical evidence (Khemani et al 2016b) reveals that the phase obtained is stable and robust. It is shown that broken symmetry phases are stable against weak local deformations of Floquet drives and it is stressed that the order revealed by π-spin glass (i.e. Floquet time crystal) is always spatio-temporal and never purely temporal as related to general symmetry properties—for details see von Keyserlingk et al (2016). The authors stress that absolutely stable states are generally possible by a combination of Floquet periodicity, broken symmetries and MBL—on the other hand no proof is given that MBL is a necessary condition for observing absolutely stable phases. It is shown that signatures of stable time crystals may be obtained starting from generic short-range initial states.

In a number of papers (Else and Nayak 2016, Potter et al 2016, Po et al 2016, Roy and Harper 2016, von Keyserlingk and Sondhi 2016a, Harper and Roy 2017) a further classifications of possible phases in periodically driven system have been attempted with the stress on novel interesting topological opportunities offered by Floquet Hamiltonians. This goes beyond the modest scope of this review centered on time crystals so we only mention those papers for an interested reader.

An independent clarification of time crystals for spin systems has been made by Else et al (2016). They consider a general time-periodic Hamiltonian and formulate simple criteria for the occurrence of a discrete time translation symmetry breaking (DTTSB) in terms of Floquet eigenstate properties. A simple statement is that Floquet eigenstates cannot be short-range correlated. Note that, as nicely pointed out (Else et al 2016), this condition implies lack of MBL in Floquet systems exhibiting DTTSB if MBL is strictly understood as an integrable disordered system with a complete set of local integrals of motion (LIOMs) (Serbyn et al 2013, Huse et al 2014). Non local Floquet states can not be all the common eigenstates of the Hamiltonian and the quasi-local operators so the system of LIOMs (together with the Floquet Hamiltonian) cannot form a complete set of observables. As shown (Else et al 2016) the condition of non local character of Floquet eigenstates is equivalent to the statement that for any time instant, τ and any state $\newcommand{\ra}{\rangle} \vert \psi(\tau)\rangle$ with short-range correlations there exist an operator $\hat O$ with the average in state $\psi(\tau)$ being not time-periodic with period T despite T-periodicity of the Hamiltonian itself. Note that this property is operational, i.e. it suggests that DTTSB can be detected by measurements on such an operator.

The example studied, a simple spin system (very close to that considered by Khemani et al (2016b), von Keyserlingk et al (2016) and von Keyserlingk and Sondhi (2016b)) applies within a single period two operations, a global rotation around x axis by π in time $t_1=\pi/2$ with the unitary

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \label{eqelse0} U_1=\exp\left[{\rm i}t_1\sum_i \sigma_x\right], \nonumber \end{align} \tag{ 31 }$$

followed by an evolution for time t₂ in a disordered Hamiltonian H_MBL,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \label{eqelse1} U_2=\exp\left[-{\rm i}H_{\rm MBL}t_2\right], \nonumber \end{align} \tag{ 32 }$$

with $T=t_1+t_2$ and

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \label{eqelse2} H_{\rm MBL}= \sum_i \left(J_i\sigma_i^z\sigma_{i+1}^z+h_i^z\sigma_i^z+h_i^x\sigma_i^x\right). \nonumber \end{align} \tag{ 33 }$$

The parameters $J_i, \ h_i^z$, and $h_i^x$ are uniformly and randomly chosen, $J_i\in[J/2, 3J/2]$, $h_i^z\in[0, h^z]$, and $h_i^x\in[0, h]$ with, typically $J=h^z=1$ and $h\ll J$. At the prize of repeating the reasoning presented above, we mention that the system is soluble for $h=0$, the eigenstates of H_MBL being the product of eigenvectors of individual $\sigma_i^z$, $\newcommand{\ra}{\rangle} \vert \{m_i\}\rangle_z$, $\newcommand{\ra}{\rangle} \sigma_k^z\vert \{m_i\}\rangle_z=m_k\vert \{m_i\} \rangle_z$ with $m_k=\pm 1$. Explicitly

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\ra}{\rangle} \displaystyle H_{\rm MBL}\vert \{m_i\}\rangle_z = [E_2(\{m_i\}) +E_1(\{m_i\})]\; \vert \{m_i\}\rangle_z, \nonumber \end{align} \tag{ 34 }$$

with $E_2(\{m_i\})=\sum_i J_im_im_{i+1}$ and $E_1(\{m_i\})=\sum_i h_i^zm_i$. Suppose for simplicity that all $h_i^z=0$. For $t_1=\pi/2$ in (31) the action of U₁ is to flip all the spins $\{m_i\}\rightarrow\{-m_i\}$ while U₂ adds a global phase. Clearly $\newcommand{\ra}{\rangle} \vert \{m_i\}\rangle_z$ cannot be a Floquet eigenstate—but both Schrodinger-cat like combinations $\newcommand{\ra}{\rangle} \vert \pm\rangle = (\vert \{m_i\}\rangle_z\pm\vert \{-m_i\}\rangle_z$ are T-periodic Floquet eigenstates. Such cat states are fragile to perturbations and DTTSB selects one of the pair (at random). Both $\newcommand{\ra}{\rangle} \vert \{m_i\}\rangle_z$ and $\newcommand{\ra}{\rangle} \vert \{-m_i\}\rangle_z$ evolve periodically at exact resonance with a period 2T. This reasoning can be easily extended to a still soluble case of $t_1=\pi/2$ and $h_i^z\in[0, h^z]$ (Else et al 2016).

The most important feature is that those TTSB solutions are robust against perturbations such as non-zero h or t₁ leading to a rotation around x axis by an angle slightly different from π. In effect $\newcommand{\ra}{\rangle} \newcommand{\la}{\langle} \langle\sigma_z\rangle$ shows pronounced oscillations with the period of 2T (compare figure 6) revealing DTTSB. The numerical studies use time-evolving block decimation (TEBD) scheme (Vidal 2004) for a system of 200 spins and $h=0.3$. It is difficult to separate the effects due to interactions and those due to randomness since the latter is inherent in the model studied.

Zoom In Zoom Out Reset image size

Quite a similar system is studied by Yao and coworkers (Yao et al 2017) who stress the robustness of the observed effect against perturbations. They concentrate on small deviations, say , of the ordered phase rotation, i.e. $\newcommand{\e}{{\rm e}} t_1=\pi/2-\epsilon$ in U₁ (31) and actually draw a phase diagram in the interactions J_z- plane—see figure 7. The system they study is obtained from (31)–(33) by taking $J_i\in [0.8J_z, 1.2J_z]$ and $h_i^x=0$ with $t_2=1$ (we are grateful to Norman Yao for providing us with updated information in the parameters used). The data leading to panel figure 7(a) are obtained after averaging over 100 disorder realizations.

Zoom In Zoom Out Reset image size

The main claim of the paper is that there is a phase transition from DTTSB phase to trivial paramagnet which is of the Ising type. Importantly, Yao and coworkers (Yao et al 2017) consider not only the abstract spin system but discuss small systems with long range interactions providing quantitative predictions for controlled ions chain. In that case the unitary evolution over the period $T=t_1+t_2$ is similar to the one given above with interactions between nearest neighbours replaced by the power law decaying long range interactions: $U(T)=U_1(t_1)U_2(t_2)$ with $\newcommand{\e}{{\rm e}} U_1(t_1)=\exp(-{\rm i}t_1\sum_i \sigma_i^x)$ and

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \label{yaoeq} U_2(t_2)= \exp\left[-{\rm i}t_2\left(\sum_{i\ne j}\frac{J_{ij}}{r_{ij}^\alpha} \sigma_i^z \sigma_j^z +\sum_i h_i^z \sigma_i^z\right)\right]. \nonumber \end{align} \tag{ 35 }$$

Interestingly the DTTSB time crystal behaviour seems to be even more robust with long range interactions as also evident from additional material in supplementary information to that paper. It shows also a difference between a disordered and clean cases showing in the character of the Fourier transform of spin-spin correlation function as reproduced in figure 8.

Zoom In Zoom Out Reset image size

A detailed proposition of ion chain experiment has proven successful. The joint effort with experimental team of Chris Monroe resulted in the report of an experimental realization of discrete time-crystals (Zhang et al 2017). The ion chain consists of 10 $^{171}$Yb$^+$ ions, each providing two sublevels of ²${\rm S}_{1/2}$ state: $\newcommand{\ra}{\rangle} \vert F=1, m_F=0\rangle$ and $\newcommand{\ra}{\rangle} \vert F =0, m_F=0\rangle$ denoted as effective 1/2-spin states $\newcommand{\ra}{\rangle} \newcommand{\dow}{\downarrow} \vert \downarrow\ra_z$ and $\newcommand{\ra}{\rangle} \vert \uparrow\ra_z$, respectively. Operations on a single spin are performed using optical Raman transitions. Interactions between spins are due to spin-dependent optical dipole forces. The AC Stark shifts due to off-resonant tightly focused laser beam allow one to introduce programmable disorder addressing each ion. This is accompanied by $\pi/2$ pulses to transfer the disorder from $\sigma_i^z$ to $\sigma_i^x$. In effect the following Floquet unitary operator is realized:

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \label{ion} U(T)&=\exp(-{\rm i}H_3t_3)\exp(-{\rm i}H_2t_2)\exp(-{\rm i}H_1t_1) \quad {\rm with} \nonumber \\ H_1&=\Omega(1-\epsilon)\sum_i\sigma_i^y, \nonumber \\ H_2&=\sum_{i,j}J_{ij}\sigma_i^x\sigma_j^x, \nonumber \\ H_3&=\sum_i h_i \sigma_i^x. \nonumber \end{align} \tag{ 36 }$$

The Rabi rotation $\Omega=\pi/2t_1$ so for $\newcommand{\e}{{\rm e}} \epsilon=0$ spins are flipped (as in theory propositions above). Non-zero, controllable allows for checking the robustness of the Floquet time crystal formed. The long range interactions $J_{ij}=J_0/\vert i-j\vert ^\alpha$ with $\alpha=1.5$, site dependent disorder corresponds to h_i uniformly chosen in the interval $[0, W]$ (Zhang et al 2017). While each of the steps in (36) looks simple and straightforward the impressive experimental strategy has to be carried out to actually implement (36) with the desired precision, the interested reader is advised to consult the original work and the supplementary information (Zhang et al 2017). A single Floquet period lasts $T=t_1+t_2+t_3\approx(15+27+33)\mu s$ being limited by a time needed for sufficient interactions between spins as well as by the fact that disorder h_i is applied to the ions consecutively.

The system is initiated in the separable state of all spins pointing downwards in the x-direction $\newcommand{\ra}{\rangle} \newcommand{\dow}{\downarrow} \vert \Psi\rangle={\vert \downarrow\downarrow\dots\downarrow\rangle_x}$ with $\newcommand{\ra}{\rangle} \newcommand{\dow}{\downarrow} \vert \downarrow\rangle_x=\frac{1}{\sqrt{2}}(\vert \downarrow\rangle_z+ \vert \uparrow\rangle_z)$. After evolution with a variable number of Floquet periods the magnetization along x is measured yielding the time correlation function

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \newcommand{\ra}{\rangle} \displaystyle \langle\sigma_i^x(t)\rangle=\langle\Psi\vert \sigma_i^x(t)\sigma_i^x(0)\vert \Psi\rangle. \label{monroe_mag} \nonumber \end{align} \tag{ 37 }$$

A typical duration of the experiment is about 100 Floquet periods. For optimally chosen t₁ (i.e. for $\newcommand{\e}{{\rm e}} \epsilon=0$ in (36)) and in the absence of the disorder and interactions one expects the magnetization to restore after integer multiples of 2T because a single period evolution flips the spins. For $\newcommand{\e}{{\rm e}} \epsilon\ne 0$, but small, spins rotation reveals beating that is indicated by two peaks in the Fourier transform of the time evolution of the magnetization—the peaks are located symmetrically around $1/2T$, see figure 9(a). Imperfections like the presence of disorder introduce dephasing among the spins, figure 9(b). The crucial observation is that the presence of interactions restores the orderly behaviour, the Fourier transform of the correlation function reveals a single peak at half the Floquet frequency, compare figure 9(c). The self-organization of the system may be understood as a manifestation of time crystal behaviour. The initial product state in the experiment $\newcommand{\ra}{\rangle} \newcommand{\dow}{\downarrow} \vert \downarrow\downarrow\dots\downarrow\ra_x$ may be thought of as a superposition of two Floquet eigenstates that are Schrödinger cat-like states $\newcommand{\ra}{\rangle} \newcommand{\dow}{\downarrow} \vert \pm\ra=\frac{1}{\sqrt{2}}(\vert \downarrow\downarrow\dots\downarrow\ra_x\pm\vert \uparrow\uparrow\dots\uparrow\ra_x)$. The Floquet eigenstates themselves evolve with the period T but they are extremely fragile. Even if the system was prepared in one of the Floquet eigenstates, it would collapse to a product state under infinitesimally weak perturbation and starts evolving with the period 2T that indicates spontaneous breaking of the original discrete time translation symmetry.

Zoom In Zoom Out Reset image size

Back to back report in the journal Nature (Choi et al 2017) describes Floquet time-crystalline order observed in a disordered dipolar many body system—somehow as far from the controlled 10 ions chain as one could imagine. The system studied is an ensemble of nitrogen vacancy spin impurities in diamond. Altogether about 10⁶ impurities contribute to the signal. Each nitrogen vacancy has unit spin but a two-level system is isolated, by applying an external magnetic field, and it can be manipulated using microwave radiation. In effect the following Hamiltonian is realized:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle H(t)&=\sum_i \left[\Omega_x(t)\sigma^x_i+ \Omega_y(t)\sigma^y_i+\Delta_i\sigma^z_i\right] \nonumber \\ &\quad+\sum_{ij}\frac{J_{ij}}{r_{ij}^3}\left(\sigma^x_i\sigma^x_j + \sigma^y_i\sigma^y_j - \sigma^z_i\sigma^z_j\right), \nonumber \end{align} \tag{ 38 }$$

where $\sigma^\mu_i$ are Pauli spin-1/2 operators, $\Omega_x$ and $\Omega_y$ are Rabi frequencies of microwave drivings, $\Delta_i$ are disorder on-site fields which may be only qualitatively estimated, $J_{ij}$ are dipolar, direction dependent coefficients with couplings decaying as a third power of the distance between vacancies.

While the physical system is very different from the one discussed previously, the assumed experimental scheme is a bit similar as shown in the top of figure 10. Initially spins are polarized along x-axis in the initiation stage. Then during a single period T two strategies are repeated: driving along x with $\Omega_x$ for $\tau_1$ followed by driving along y with $\Omega_y$ for $\tau_2$ (the Floquet period is $T=\tau_1+\tau_2$). As compared with the former experiment here the disorder and interactions are not switched off/on periodically but are given by a particular system realization. Almost π pulses are realized with $\Omega_y$, $\theta=\Omega_y\tau_2\approx\pi$. A rotation along x with $\tau_1$ serves as an effective 'interaction' time. For relatively short $\tau_1=92$ ns if $\theta=\pi$ a clear 2T-signal is observed in magnetization (or its Fourier transform)—however deviations from π destroy this synchronization, compare figure 10(c). Only when $\tau_1$ time is significantly increased ($\tau_1=989$ns) the Fourier transform peak remains locked at $1/2T$ regardless of the detuning of θ from the optimal value equal π. Thus the system self-organizes into a time crystal like behaviour with 2T periods.

Zoom In Zoom Out Reset image size

It is worth stressing that the same work reports also the 3T- behavior that may be realized modifying the driving scheme is such a way that after 3T the system of spins rotates closely to the initial configuration.

Interestingly authors of Choi et al (2017) stress that their experiment is not performed in the MBL regime, see also Zeng and Sheng (2017), which seems to contradict the claim that MBL is a prerequisite for time crystal behavior (Else et al 2016, Khemani et al 2016b, von Keyserlingk et al 2016, Yao et al 2017). Moreover, a recent study of long-range interacting systems (Ho et al 2017) further elaborates on that point—linking the stability of time crystals in the system studied to quantum critical regime. The authors, using perturbative arguments for their 3D system, suggest, that there is a sharp crossover between a region where the discrete time crystal is stable due to quantum criticality and the regime where discrete time crystal eventually dies. Early, disorder free propositions of discrete time crystals (Sacha 2015b) obtain their stability via a semiclassical connection to the classical stability of motion in a stable resonant $s:1$ elliptical island. Recently another realizations of clean time crystals in cold atoms have been proposed (Huang et al 2017). A series of clean quasi one-dimensional models is shown to exhibit Floquet time crystalline behaviour, robust in the strongly interacting regime due to emergent integrals of motion in the dynamical systems studied. It will be interesting to consider the connection between the last two approaches. Last but not least it has been shown that clean Floquet time crystal behaviour may be obtained in the well known kicked Lipkin–Meshkov–Glick model (Russomanno et al 2017). The stability of the time-crystal phase can be in this model directly analyzed in the limit of infinite size, discussing the properties of the corresponding classical phase space structure. Thus the role played by disorder and MBL in time crystal phenomenon is by no means obvious. There exists the numerical evidence (compare figure 8) that for some spin systems disorder makes time crystal behaviour stronger, still it does not seem to be necessary.

The real area of speculation (or careful analysis) opens when one considers the thermodynamic limit. It is possible that exponentially slow heating may occur, under favorable circumstances in periodically driven many body systems (Abanin et al 2015, 2017). The resistance of driven systems to expected heating has been even verified experimentally (Bordia et al 2017). Thus one may consider that before the inevitable thermalization, even in the absence of MBL, for some quite long time one may observe, in a finite system, a time crystal behavior. Such a behavior is called a pre-thermalization and the corresponding time crystal properties have been considered in such a regime (Else et al 2017). The authors of diamond experiment (Choi et al 2017) argue that their system cannot be considered as 'pre-thermalizing' yet they observe the time crystal features in the experiment. It may well be, therefore that for quite large yet finite systems the time crystal behavior remains robust for sufficiently long times that other mechanisms such as decoherence (Lazarides and Moessner 2017) take over. Clearly more studies of robustness and limitations of time-crystal phenomenon in different settings is desirable (Moessner and Sondhi 2017).

In condensed matter physics it is often assumed that a space crystal is already formed and its properties are analyzed with the help of a space periodic Hamiltonian $\newcommand{\la}{\langle} H(x+\lambda)=H(x)$. For ultra-cold atoms in an optical lattice the space periodicity is imposed by external laser fields, so a space periodic Hamiltonian arises without the spontaneous symmetry breaking. Counterparts of space periodic systems but in the time domain are periodically driven systems where $H(t+T)=H(t)$. As mentioned already there is a vast literature on the properties of periodically driven (sometimes called Floquet) systems. We concentrate here on a simple question whether systems with time periodicity may reveal non-trivial crystalline structures in time.

Guo et al (2013) have considered a single particle in a perturbed harmonic oscillator potential that is periodically driven,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H=\frac{p^2}{2}+\frac{\Omega^2z^2}{2}+\frac{\nu}{2}z^4+2V_0 z^s\cos\omega t \label{guo_ham} \nonumber \end{align} \tag{ 39 }$$

and have shown that the resonant driving, i.e. when $\omega\approx s\Omega$ with integer $s\gg 1$, results in the band structure of the quasi-energy spectrum of the system. This result has been demonstrated within the rotating wave approximation (RWA). That is, applying the unitary transformation $U(t)={\rm e}^{{\rm i}(\omega t/s)\hat a^\dagger \hat a}$, where $\hat a, \ \hat a^\dagger$ are the usual harmonic oscillator annihilation and creation operators, and dropping fast oscillating terms, the effective Hamiltonian becomes

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H_{\rm RWA}&=\left(\Omega-\frac{\omega}{s}\right)\hat a^\dagger \hat a+\frac{3\nu}{4\Omega^2}\hat a^\dagger\hat a(\hat a^\dagger \hat a+1) \nonumber \\ & \quad+\frac{V_0}{(2\Omega)^{s/2}}(\hat a^{\dagger s}+\hat a^{s}). \label{HRWA_Guo} \nonumber \end{align} \tag{ 40 }$$

The key observation is that H_RWA displays a new discrete symmetry that is not visible in the original Hamiltonian (39), i.e. it commutes with a unitary operator ${\rm e}^{-{\rm i}2\pi\hat a^\dagger\hat a/s}$. It is convenient to introduce the radial $\hat r$ and the angular $\hat \theta$ operators via $\newcommand{\la}{\langle} \hat a={\rm e}^{-{\rm i}\hat \theta}\hat r/\sqrt{2\lambda}$ and $\newcommand{\la}{\langle} \hat a^\dagger=\hat r{\rm e}^{{\rm i}\hat \theta}/\sqrt{2\lambda}$ where $\newcommand{\la}{\langle} \lambda=3\nu/(4\Omega^2\vert \omega/s-\Omega\vert)$ is an effective Planck constant. Then, $\newcommand{\la}{\langle} H_{\rm RWA}=\lambda^{-1}\vert \frac{\omega}{s}-\Omega\vert \hat g$ where

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \hat g=\frac14(\hat r^2+\lambda-1)^2+\frac{\mu}{2}\left[\left(\hat r{\rm e}^{{\rm i}\hat\theta}\right)^s+\left(\hat r{\rm e}^{-{\rm i}\hat\theta}\right)^s\right], \label{hrwa_radial} \nonumber \end{align} \tag{ 41 }$$

with

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \mu=\frac{2\lambda V_0}{\vert \omega/s-\Omega\vert }\left(\frac{\Omega\vert \omega/s-\Omega\vert }{3\nu}\right)^{s/2}. \nonumber \end{align} \tag{ 42 }$$

Looking at (41) it is easy to realize that H_RWA is invariant under the discrete rotation $\hat \theta\rightarrow \hat\theta +2\pi/s$, eigenstates of the system in the angular representation are Bloch-like waves $u_m(\theta){\rm e}^{{\rm i}m\theta}$ where $u_m(\theta+2\pi/s)=u_m(\theta)$ and the spectrum splits into bands as illustrated in figure 11. The performed RWA is valid provided $\Omega(\frac{\omega}{s}-\Omega){\hspace{0pt}}^2<\frac32 \nu$ what indicates the necessity of the anharmonic term in the Hamiltonian (39) (Guo et al 2013).

Zoom In Zoom Out Reset image size

Guo et al called such a behaviour a phase space crystal. We will see that a phase space crystal is actually a time crystal. That is, any resonantly driven system can be reduced to a solid state-like Hamiltonian by means of the secular approximation approach (Lichtenberg and Lieberman 1992, Buchleitner et al 2002) and a periodic structure in the phase space is always reproduced in the time domain (Sacha and Delande 2016).

Let us first illustrate this conjecture on the problem of a single particle bouncing on an oscillating mirror in the presence of the gravitational field, see the Hamiltonian (21). It was said that a suitable choice of the system parameters allows one to realize single particle Floquet eigenstates where localized wave-packets evolve along classical resonant periodic orbits. We have already analyzed the case of the 2:1 resonance. Now let us consider a general $s:1$ resonance (Sacha 2015a). With a suitable choice of parameters it is possible to realize Floquet eigenstates that consist of s well localized wave-packets which evolve along $s:1$ resonant orbit each of them with a period s times longer than the period T of the mirror oscillations. They exchange their positions so that the Floquet eigenstates are periodic with the fundamental period T, see an example for $s=4$ in figure 12(a). There are s such Floquet eigenstates $u_1(z, t), \dots, u_s(z, t)$ and they are related to quasi-energies $E_1, \dots, E_s$ which are nearly degenerate modulo $\omega/s$. When $s\rightarrow\infty$ these quasi-energies start forming a quasi-energy band. Moreover, one can extract s individual localized wave-packets $\phi_1(z, t), \dots, \phi_s(z, t)$ by a proper superpositions of s Floquet eigenstates. We anticipate that when $s\rightarrow\infty$ these wave-packets become counterparts of Wannier states known in solid state physics that are localized in sites of a spatially periodic potential and are obtained from Bloch eigenstates by a proper superposition (Kohn 1959).

Zoom In Zoom Out Reset image size

Now we have all elements in order to show that our system can reveal crystalline structure in the time domain. In figure 12(a) an example of a Floquet eigenstate is presented in the configuration space for two different moments of time so one can observe at which direction each of the localized wave-packets is propagating. However, the configuration space is not the domain where we can see crystalline structures. Let us choose a point in the configuration space, e.g. close to the turning point of a classical particle, and plot how the probability density for the detection of a particle at this point changes in time, figure 12(b). We can see that each wave-packet arrives at this point periodically and in the same way. That is guaranteed by the Floquet theorem. Thus, the time domain is a proper 'space' where a crystalline structure emerges. This conjecture becomes clearer when one restricts the description of a particle to the Hilbert subspace spanned by s evolving localized wave- packets and calculates the corresponding quasi-energy of a particle,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle E&=\int\nolimits_0^{s2\pi/\omega}{\rm d}t\int\nolimits_0^\infty {\rm d}z\;\psi^*\left(H(t)-{\rm i}\partial_t\right)\psi \nonumber \\ &\approx -\frac12\sum\limits_{j=1}^{s}(J_{j}a_{j+1}^*a_j+{\rm c.c.}), \label{sr_s15} \nonumber \end{align} \tag{ 43 }$$

where $\psi(z, t)\approx \sum_{j=1}^sa_j\phi_j(z, t)$ was substituted (Sacha 2015a). Equation (43) is nothing but the kinetic energy term of the tight-binding model known from solid state physics. The localized wave-packets, $\phi_j$, play a role of Wannier states and

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle J_j=-2\int\nolimits_0^{s2\pi/\omega}{\rm d}t\int\nolimits_0^\infty {\rm d}z\;\phi_{j+1}^*\left(H(t)-{\rm i}\partial_t\right)\phi_j \nonumber \end{align} \tag{ 44 }$$

are tunneling amplitudes. They are related to tunneling of a particle between neighbouring wave-packets (but neighbouring in the time domain) see figure 12(b). All the tunneling amplitudes have the same absolute values $\vert J_j\vert =J$. There are also longer range tunnelings, e.g. between next-to-nearest neighbour wave-packets but they are at least two orders of magnitude weaker in a typical situation and, similarly as in the tight-binding model, can be neglected.

We have demonstrated that the description of a resonant motion of a particle bouncing on an oscillating mirror can be reduced to a tight-binding model known in condensed matter physics (Guo et al 2013, Sacha 2015a). We have already mentioned that such an approach can be applied to any resonantly driven dynamical system ranging from ultra-cold atomic gasses to a Rydberg electron perturbed by a microwave field (Sacha and Delande 2016). Indeed, assume that a single particle integrable classical system described by an unperturbed Hamiltonian H₀ is driven periodically, i.e. there is another term in a system Hamiltonian, $H=H_0+H_1(t)$, which $H_1(t+T)=H_1(t)$. It is convenient to describe classical motion in action-angle variables of the unperturbed system (Lichtenberg and Lieberman 1992). For example in the case of a 1D system, we can perform canonical transformation from Cartesian position and momentum to canonically conjugate action I and angle θ for which $H_0=H_0(I)$. Then, the unperturbed motion can be solved immediately, $I=\rm constant$ and $\theta=\frac{\partial H_0(I)}{\partial I}t+\theta_0$. The latter equation describes a position of a particle on its trajectory. Now, let us assume that we drive the system harmonically with frequency ω which is close to a multiple of frequency of an unperturbed motion,

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \omega\approx s\frac{\partial H_0(I_s)}{\partial I_s}, \nonumber \end{align} \tag{ 45 }$$

with an integer s. Then, in the moving frame,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \Theta = \theta-\frac{\omega t}{s}, \label{Th_th} \nonumber \end{align} \tag{ 46 }$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle P = I-I_s, \nonumber \end{align} \tag{ 47 }$$

are slowly varying. Averaging the original Hamiltonian over time we obtain an effective time-independent (secular) Hamiltonian (Lichtenberg and Lieberman 1992)

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H\approx\frac{P^2}{2m}+V_0\cos(s\Theta), \label{pendulum_sr15} \nonumber \end{align} \tag{ 48 }$$

where m and V₀ are the effective mass and the amplitude of the effective Hamiltonian respectively. Thus, in the moving frame any resonant motion with $s\gg 1$ reduces to a solid state problem of a particle in a space periodic potential (48) (Lichtenberg and Lieberman 1992, Guo et al 2013, Sacha 2015a, Guo and Marthaler 2016, Sacha and Delande 2016).

In order to obtain a quantum effective description one can either quantize a classical secular Hamiltonian (48) (Buchleitner et al 2002) or apply the quantum version of the secular approximation from the very beginning (Berman and Zaslavsky 1977). By restricting to the lowest energy band of a quantum version of (48) one ends up with a tight-binding model of the form of equation (43). Interestingly, it is also possible to consider higher energy bands of the secular Hamiltonian (48) and investigate multi-band physics (Guo et al 2013, Sacha 2015a, Guo and Marthaler 2016, Sacha and Delande 2016).

Note that if a wave-function of a particle in the moving frame, $\psi(\Theta)$, reveals some crystalline structure, in the laboratory frame this wave-function, $\psi(\theta-\omega t/s)$, will show crystalline properties also versus time t if one fixes position close to a classical orbit because the transformation (46) is linear (Sacha and Delande 2016). On the other hand there is no guarantee that such a wave-function will reveal periodic behaviour as a function of a Cartesian position in space for fixed time because canonical transformation between action-angle variables and Cartesian position and momentum is non-linear in general. This is the case for example for a particle bouncing on an oscillating mirror where in the time domain one observes a crystalline structure but does not see it in the configuration space, figure 12.

Energy eigenstates of a single particle in the presence of a space periodic potential are given by Bloch waves that are extended in space. However, when a potential is not strictly periodic because some disorder is present, eigenstates may become exponentially localized around different points of configuration space due to destructive interference between different multiple scattering paths (Anderson 1958). Regardless of the magnitude of the disorder, provided it is random and decorrelated, Anderson localization is inevitable in one-dimensional (1D) systems as well as for time-reversal invariant spinless two dimensional (2D) systems (Abrahams et al 1979, Müller and Delande 2011). In the three dimensional (3D) world the situation is more complicated, typically disordered system reveal the so-called mobility edge. Eigenstates with energies below the mobility edge are localized and those above are extended (Abrahams et al 1979, Müller and Delande 2011).

Anderson localization can be also observed in the momentum space being called then the dynamical localization (Fishman et al 1982, Haake 1991, Moore et al 1995, Stöckmann 1999, Lemarié et al 2009). Interestingly, periodically driven single particle systems can reveal Anderson localization in the momentum space despite the classical diffusive behaviour. The diffusion in the phase space is suppressed by quantum interference effects.

In the present section it will be shown that yet another kind of Anderson localization is possible—Anderson localization in the time domain due to the presence of disorder in time (Sacha 2015a), i.e. when the disorder is added on top of the periodically changing force.

4.2.1. Anderson localization in time crystals.

The simple idea to realize Anderson localization in the time domain is to add random vibrations to a mirror that oscillates with a period T in the bouncing particle problem that has been already described above. It is assumed that a bouncing particle is resonantly driven by a periodically oscillating mirror and $s:1$ resonance condition is fulfilled. Moreover, there is an additional term $H^{\prime}(t)$ in the system Hamiltonian that introduces a weak temporal disorder—$H^{\prime}(t)$ fluctuates in time but fulfills periodic boundary condition on a long time scale, $H^{\prime}(t+sT)=H^{\prime}(t)$. It leads to additional terms in the tight-binding energy (43) which now takes the form

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle E \approx -\frac12\sum\limits_{j=1}^{s}(J_{j}a_{j+1}^*a_j+{\rm c.c.})+\sum\limits_{j=1}^s\epsilon_j\vert a_j\vert ^2, \label{disord_s15} \nonumber \end{align} \tag{ 49 }$$

where

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle \epsilon_j=\int\nolimits_0^{s2\pi/\omega}{\rm d}t\int\nolimits_0^\infty {\rm d}z\;H^{\prime}(t)\;\vert \phi_j\vert ^2, \label{disord_s15aa} \nonumber \end{align} \tag{ 50 }$$

are random numbers that can belong to any distribution by a proper engineering of the fluctuating part, $H^{\prime}(t)$, of the Hamiltonian (Sacha 2015a). Energy (49) constitutes actually a 1D Anderson model which possesses exponentially localized eigenstates. In figure 13 an example of a localized eigenstate is shown in the case when $\newcommand{\e}{{\rm e}} \epsilon_j$'s belong to a Lorentzian distribution.

Zoom In Zoom Out Reset image size

The entire Hamiltonian is time periodic with a long period sT. It guarantees existence of Floquet eigenstates periodic also with a period sT. Anderson localization in the time domain requires that the localization length in time has to be much smaller than sT. If a detector is situated close to a classical turning point in the configuration space, it is expected to click with probability that is localized exponentially around a certain moment of time and such a behaviour is repeated periodically with a period sT. It is in a complete analogy to Anderson localization in a 1D space crystal with periodic boundary conditions in space (a ring topology). In figure 14 a comparison of Anderson localization in a time crystal and a 1D space crystal is presented. In the space crystal case when a disorder is present, a particle localizes around a certain space point and moving periodically around the ring, one observes periodically a localized density profile. The larger the ring, the larger the space crystal. In the time crystal case, probability for a measurement of a particle at fixed space point is exponentially localized in time and this behaviour is repeated periodically with a period sT and the higher $s:1$ resonance the larger the time crystal.

Zoom In Zoom Out Reset image size

4.2.2. Anderson localization in time without crystalline structures.

Anderson localization in configuration space takes place when a disorder is present in a space crystal. However, in order to observe Anderson localization it is not necessary to start with a space crystal. One can begin with a purely disordered time-independent potential characterized by a finite correlation length and Anderson localization can be still observed (Kuhn et al 2007). It turns out that an analogous situation takes place in the time domain that will be illustrated with a simple model (Sacha and Delande 2016).

Let us consider a single particle on a ring described by the following classical Hamiltonian,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H=\frac{p^2}{2}+V_0g(\theta)\,f(t), \label{h_sd16} \nonumber \end{align} \tag{ 51 }$$

where θ is an angle that denotes a position of a particle on a ring and V₀ is amplitude of the perturbation. The perturbation is assumed to be regular in space (for fixed time t) and disordered in time (for fixed θ). In the following it is assumed that $g(\theta)=\theta/\pi$ for $\theta\in[-\pi, \pi[$ but because $-\pi$ and π correspond to the same point on the ring, $g(\theta)$ has a discontinuity and its Fourier expansion reads,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle g(\theta)=\sum\limits_{n\ne 0}g_n{\rm e}^{{\rm i}n\theta}, \label{g_s15} \nonumber \end{align} \tag{ 52 }$$

where $g_n=i(-1){\hspace{0pt}}^{n}/(\pi n)$. The function $f(t+2\pi/\omega)=f(t)$ is periodic but between $t=0$ and $2\pi/\omega$ it performs random fluctuations, i.e.

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle f(t)=\sum\limits_{k\ne 0}f_k{\rm e}^{{\rm i}k\omega t}, \nonumber \end{align} \tag{ 53 }$$

where $f_k=f_{-k}^*$ are independent random variables. Switching to the moving frame, $\Theta=\theta-\omega t$ and $P=p-\omega$, one can see that new position Θ and momentum P are slowly varying variables if resonant condition is fulfilled, $P\approx 0$. Then, averaging the Hamiltonian over the fast time variable results in a time-independent effective Hamiltonian,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H \approx \frac{P^2}{2}+U(\Theta), \label{heff_sd16} \nonumber \end{align} \tag{ 54 }$$

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle U(\Theta) = V_0\sum\limits_{k\ne 0}g_k f_{-k} {\rm e}^{{\rm i}k\Theta}+{\rm constant}, \label{heff_sd16a} \nonumber \end{align} \tag{ 55 }$$

where a particle moves in the presence of a time-independent disordered potential $U(\Theta)$. This potential can be thought of as the coherent addition of resonant terms between spatial harmonics of the potential and the corresponding temporal harmonics of the disordered driving amplitude (Sacha and Delande 2016). Statistical properties of Fourier components of $U(\Theta)$ can be engineered by a proper choice of a distribution for f_k. The coefficients g_k drops like $1/k$ and in order to deal with the disordered potential characterized by a small correlation length, the drop of g_k has to be compensated by an increase of absolute values of f_k.

Assume that $\vert g_kf_{-k}\vert \propto e^{-k^2/(2k_0^2)}$ and ${\rm Arg}(\,f_{-k})$ are random numbers chosen uniformly in the interval $[0, 2\pi[$. With a suitable choice of the amplitude V₀ and the correlation length of the disordered potential, which is equal to $\sqrt{2}/k_0$, eigenstates of the Hamiltonian (54) are Anderson localized even for energies higher than the standard deviation of the disordered potential (Sacha and Delande 2016). These eigenstates when plotted versus time for a fixed space point in the laboratory frame show up Anderson localization—an example is presented in figure 15.

Zoom In Zoom Out Reset image size

The model (51) can be realized in ultra-cold atomic gases as suggested in Delande et al (2017). Instead of trapping atoms in a ring-shape potential one can prepare a 1D sawtooth-shape potential that is modulated in time in a random manner. It is assumed that ultra-cold atoms are initially prepared in a time-independent optical lattice potential in a Mott insulator phase (Pethick and Smith 2002) that results in a sequence of independent slices of an atomic cloud, see figure 16. Then, the optical lattice is turned off, a modulated sawtooth potential is turned on and atoms are kicked so that their velocity fulfills resonance condition with the modulated sawtooth potential. As a result of Anderson localization in time, the slices of the atomic cloud keep their shape—they do not spread even though they fly with energy much greater than the amplitude of the modulated sawtooth potential.

Zoom In Zoom Out Reset image size

Other systems that seem promising for realization of Anderson localization in time are Rydberg atoms perturbed by fluctuating microwave field (Giergiel and Sacha 2017). A highly excited electron experiences resonant linearly polarized microwave field that performs random fluctuations. An additional static electric field along the microwave polarization axis stabilizes a Kepler orbit an electron is traveling along. It results in an effective Hamiltonian of the form of equation (54) in the frame moving with an electron because other degrees of freedom of the 3D system are frozen due to the presence of the external fields. One can observe Anderson localization of an electron along a classical Kepler orbit that corresponds to Anderson localization in time if a detector is placed at, e.g. nucleus.

In the present section we have so far considered single particle systems. However, not only a single atom but a bunch of atoms can bounce on an oscillating mirror (Sacha 2015a). Resonant bouncing of a single particle is described by quasi-energy of the form of a tight-binding model, equation (43). In the many-body case a similar simplification can be used if the interaction energy per particle is not much greater than a width of the quasi-energy band of the system (43). Close to a $s:1$ resonance a system of ultra-cold bosons with contact interactions may be then described in a restricted Hilbert subspace spanned by Fock states $\newcommand{\ra}{\rangle} \vert n_1, \dots, n_s\ra$, where n_i denotes occupation of a localized wave-packet $\phi_i$ evolving along the $s:1$ resonant orbit. The corresponding many-body Floquet Hamiltonian reads

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H_F&=\int\nolimits_0^{s2\pi/\omega}{\rm d}t\int\nolimits_0^\infty {\rm d}z\;\hat\psi\left[H(t)+\frac{g_0}{2}\hat\psi^\dagger\hat\psi-{\rm i}\partial_t\right]\hat\psi, \nonumber \\ &\approx -\frac12\sum\limits_{j=1}^{s}(J_{j}\hat a_{j+1}^\dagger \hat a_j+{\rm h.c.}) +\frac12\sum\limits_{j=1}^s U_{ij}\hat a_j^\dagger\hat a_j^\dagger \hat a_j\hat a_j, \label{manyHsr_s16} \nonumber \end{align} \tag{ 56 }$$

where the bosonic field operator $\hat \psi\approx\sum_{j=1}^s\phi_j\hat a_j$ and $\vert J_{j}\vert =J$, see (43). The coefficients

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle U_{ij}=g_0\int\nolimits_0^{s2\pi/\omega}{\rm d}t\int\nolimits_0^\infty {\rm d}z\;\vert \phi_i\vert ^2\;\vert \phi_j\vert ^2, \nonumber \end{align} \tag{ 57 }$$

describe effective interactions between particles that occupy different localized wave-packets. Note that despite the fact the original interactions between particles are contact, the effective interactions can have long range character (Anisimovas et al 2015, Eckardt and Anisimovas 2015, Sacha 2015a, Guo et al 2016) because all localized wave-packets pass each other in evolution along a 1D resonant orbit. The on-site interactions are dominating i.e. $\vert U_{ii}\vert >\vert U_{ij}\vert$ for $i\ne j$.

For sufficiently strong repulsive ($g_0>0$) contact interactions, when $U_{ii}\gg NJ/s$, the ground state of the system (56) is a single Fock state $\newcommand{\ra}{\rangle} \vert N/s, \dots, N/s\ra$ where well defined numbers of atoms occupy each localized wave-packet, long-time phase coherence is lost, a gap opens between the ground state level and excited levels and consequently a Mott insulator phase is realized in the time domain (Sacha 2015a), see schematic picture in figure 17.

Zoom In Zoom Out Reset image size

It is quite natural now to join the results of the previous parts and consider an interacting many-particle system driven resonantly, with at the same time added perturbation that fluctuates in time but is repeated with some large period sT. That extends the analysis of section 4.2.1 on interacting particle systems. Combining with the many body language of the previous section we are led to a genuine example that may result in a many-body localization induced by temporal disorder. Work in this direction is in progress (Mierzejewski et al 2017).

It is difficult to predict at which direction, research on time crystals will develop. There are open questions and new ideas are still coming. One of very intriguing topics concerns a possibility of realization of multi-dimensional time crystals. While it is definitely not possible to augment dimensionality of time, one can envisage 1D time crystals but with properties of, e.g. 3D space crystals. One proposal was already suggested in the literature that we shortly describe now.

The problem of a single particle on a 1D ring subjected to time-fluctuating perturbation (51) was generalized to the 3D case (Delande et al 2017). Consider a particle on a 3D torus described by the following Hamiltonian:

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H=\frac{p_\theta^2+p_\psi^2+p_\phi^2}{2}+V_0g(\theta)g(\psi)g(\phi)\,f_1(t)\,f_2(t)\,f_3(t). \label{h3d_dms17} \nonumber \end{align} \tag{ 58 }$$

The function $g(x)$ is the same as in equation (52). Time-dependent functions $f_i(t)$ are periodic with periods $2\pi/\omega_i$ but between $t=0$ and $2\pi/\omega_i$ they perform random fluctuations, i.e. $f_i(t)=\sum_{k\ne 0}f_k^{(i)}{\rm e}^{{\rm i}k\omega_i t}$ where $f_k^{(i)}=f_{-k}^{(i)*}$ are random numbers. Each degree of freedom is assumed to be in a resonance with one of the time-periodic functions $f_i(t)$ so that in the moving frame where $\Theta=\theta-\omega_1 t$, $\Psi=\psi-\omega_2 t$ and $\Phi=\phi-\omega_3 t$ and for the conjugate momenta $P_\Theta=p_\theta-\omega_1\approx 0$, $P_\Psi=p_\psi-\omega_2\approx 0$ and $P_\Phi=p_\phi-\omega_3\approx 0$, the effective secular Hamiltonian reads,

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle H_{{\rm eff}}=\frac{P_{\Theta}^2+P_{\Psi}^2+P_\Phi^2}{2}+V_{{\rm eff}}(\Theta,\Psi,\Phi), \label{heff_dms17} \nonumber \end{align} \tag{ 59 }$$

with $V_{{\rm eff}}=V_0h_1(\Theta)h_2(\Psi)h_3(\Phi)$ where

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\la}{\langle} \displaystyle h_i(x)=\sum_{k\ne 0}g_k f_{-k}^{(i)}{\rm e}^{{\rm i}kx}, \label{Eq:defh_dms17} \nonumber \end{align} \tag{ 60 }$$

provided $\omega_i$ are irrational numbers. The effective potential in (59) is time-independent and it is a product of three independent disordered potentials along each degree of freedom. From the general theory of Anderson localization in 3D (Abrahams et al 1979) one expects the existence of the mobility edge. That is, for properly engineered time fluctuating functions $f_i(t)$ and for a fixed amplitude V₀, all eigenstates of (59) with energies below the mobility edge are Anderson localized and above it are not. In the laboratory frame, when one asks how the probability density for a measurement of a particle at some fixed value of, e.g. θ changes in time, it turns out it is localized around a certain moment of time provided the relevant Floquet eigenstate corresponds to energy of (59) below the mobility edge (Delande et al 2017).