Imaginary time crystal of thermal quantum matter

TL;DR

This paper introduces the concept of an imaginary time crystal in thermal quantum matter, revealing a new phase characterized by crystalline order in imaginary time and unique temperature-dependent behaviors, supported by quantum Monte Carlo simulations.

Contribution

It demonstrates the existence of an imaginary time crystal phase in thermal quantum systems, a novel state arising from bath-induced interactions and exhibiting distinctive properties.

Findings

Presence of crystalline order in imaginary time space.

Temperature-oscillating behavior of physical observables.

Distinct ground-state and thermal properties from conventional systems.

Abstract

Spontaneous symmetry breaking is responsible for rich quantum phenomena from crystalline structures to superconductivity. This concept was boldly extended to the breaking of time translation, opening an avenue to finding exotic phases of quantum matter with collective time modulation and correlation. Here we report that a thermally open quantum ensemble manifests in the dual space of imaginary time with crystalline ordering due to a bath-induced retarded interaction. Exact quantum Monte Carlo simulations are performed to show that this imaginary time crystal phase exhibits characteristic ground-state and thermal properties absent in conventional quantum manybody systems, especially the striking temperature-oscillating behavior of its physical observables.

Equations

$${\hat{\rho }}_s=\frac{1}{Z_s}{e}^{-\beta H_s-S_R},$$

$$S_R=-{\int }_0^{\beta }d\tau {\int }_0^{\beta }d{\tau }^{\prime }\sum _i(n_i(\tau )-n_0)D(\tau -{\tau }^{\prime })(n_i({\tau }^{\prime })-n_0).$$

$$D(\tau -{\tau }^{\prime })$$

$$F(\tau -{\tau }^{\prime })$$

$$\hat{H}_{i}=\left(\begin{array}[]{cc}\Sigma_{1}(k,\omega)&\omega_{d}\\ \omega_{d}&\Sigma_{2}(k,\omega)\\ \end{array}\right)\,,$$

$$S$$

$$=$$

$$D(\tau -{\tau }^{\prime })$$

$$\propto$$

$$D(\tau-\tau^{\prime})=-\frac{\alpha}{\sigma^{2}}\big{[}1-\frac{(|\tau-\tau^{\prime}|-\tau_{0})^{2}}{\sigma^{2}}\big{]}e^{-\frac{(|\tau-\tau^{\prime}|-\tau_{0})^{2}}{2\sigma^{2}}}\,.$$

$$\int d\tau d{\tau }^{\prime }\delta n(\tau )D(\tau -{\tau }^{\prime })\delta n({\tau }^{\prime })=\sum _{ij}m({\tau }_i)V({\tau }_i-{\tau }_j)m({\tau }_j),$$

$$V({\tau }_i-{\tau }_j)=\alpha e^{-\frac{(|{\tau }_i-{\tau }_j|-{\tau }_0{)}^2}{2{\sigma }^2}}.$$

Full text

I Introduction

The phases of quantum matter are, in general, categorized according to symmetries and their spontaneous breaking patterns. Recently, a striking possibility was put forward Wilczek (2012), dubbed “time crystal”, where the interacting particles spontaneously organize themselves into a “time” coherent state that breaks temporal translational symmetry. Temporal periodicities spontaneously emerge analogous to spatially ordered crystals. This phase has raised considerable interest Shapere and Wilczek (2012); Li et al. (2012); Wilczek (2013); Sacha (2015); Else et al. (2016); Khemani et al. (2016); Yao et al. (2017); Syrwid et al. (2017); Khemani et al. (2017); Russomanno et al. (2017); Gong et al. (2018); Huang et al. (2018); Sacha and Zakrzewski (2018); Iemini et al. (2018); Chew et al. (2019), and was subsequently proven forbidden in thermodynamic equilibrium states Bruno (2013); Watanabe and Oshikawa (2015). A variant of the phase was proposed  Else et al. (2016); Khemani et al. (2016); Yao et al. (2017) by considering certain non-equilibrium condition, and has quickly been observed in experiments Choi et al. (2017); Zhang et al. (2017).

Statistical field theory shows that a $d$-dimensional quantum thermal ensemble is equivalently described by the path integral representing the partition function in $d+1$ Euclidean space where time is imaginary. In the field theory dual description, the thermal distribution weighted by the well-known temperature dependent factor is fully captured by integration over auxiliary imaginary time (“iTime”) ranging from 0 to $\beta=1/k_{B}T$ ($k_{B}$ the Boltzmann constant). In this manner, finite-temperature physics is mapped to the dynamics of Euclidean action. Motivated by such a profound relation between temperature $T$ and iTime, Wilczek had speculated the novel possibility of spontaneous translational symmetry breaking in iTime in the end of his original paper, dubbed as “imaginary time crystal” (iTC) Wilczek (2012). One may wonder what class of thermal ensembles can show crystalline order in iTime space, whether the “time” analogue of phonons exist and what novel effects should manifest in ground state and thermal response functions.

Here, we find a class of open quantum ensembles that shows the iTC phase. A key feature of our mechanism is a non-monotonic bath-induced retarded interaction, whose dependence is on the relative time in two body scattering. That satisfies the two conditions what previously might have seem conflicting, i.e., being time dependent but fully translational-invariant. As we will show in the following, integrating out certain engineered local thermal baths yields such a retarded interaction, which favors crystalline patterns in iTime with periods determined by the position of interacting potential minimum. Thermodynamic quantities of the iTC exhibit a striking behavior of oscillation in temperature, which provides an experimental observable effect for this novel phase.

II Model and method.

We study hard-core bosonic models in one-dimensional (1D) lattice of $L$ sites as well as a $L\times L$ square lattice with system Hamiltonian $H_{s}=\sum_{\langle ij\rangle}-J[a_{i}^{\dagger}a_{j}+a_{j}^{\dagger}a_{i}]$ where $J$ denotes the hopping amplitude and $\langle ij\rangle$ indicates a pair of adjacent sites. We focus on the half-filling case throughout this paper. On each site $i$ the density operator of the hard-core boson $n_{i}$ additionally couples to a local bath, as shown in Fig. 1(a). We assume the total system (system+bath) is in thermodynamic equilibrium with temperature $T$. The partition function takes the form $Z_{tot}={\rm Tr}_{s}{\rm Tr}_{b}e^{-\beta H_{tot}}$ (${\rm Tr}_{s/b}$ denotes tracing over the system/bath degree of freedom). We first integrate out the bath degrees of freedom, which yields a retarded density-density interaction term in iTime. Thus the open system can be described by the reduced density matrix

[TABLE]

where $H_{s}$ is the system Hamiltonian and $Z_{s}={\rm Tr}_{s}\hat{\rho}_{s}$. $S_{R}$ describes effective action of the translational-invariant onsite bath-induced retarded interaction:

[TABLE]

where $n_{0}=\frac{1}{2}$ is the average particle density, and the site-independent kernel function $D(\tau-\tau^{\prime})$ depends on the bath properties. For instance, a harmonic oscillator bath with frequency $\omega_{0}$ gives rise to the kernel function $D(\tau-\tau^{\prime})\sim e^{-\omega_{0}|\tau-\tau^{\prime}|}+e^{-\omega_{0}(\beta-|\tau-\tau^{\prime}|)}$ for large $\beta$. We will show later that the iTC phase can exist for a broad class of kernel functions; here, we chose:

[TABLE]

where $|\tau|\leq\beta$, and $\alpha$ is the retarded interaction strength.

The model with the retarded interactions defined in Eq. (3) can be more than a toy model of merely academic interest. In general, the interaction induced by baths composed of free bosons in thermal equilibrium is usually attractive, while a fermionic bath-induced interaction is more complex as sometimes it can be oscillating between attraction and repulsion (e.g.the Ruderman-Kittel-Kasuya-Yosida (RKKY) interaction Ruderman and Kittel (1954); Kasuya (1956); Yosida (1957)). Now we need a RKKY-like interaction, but with oscillating decay in iTime instead of real space. Consider $D(\tau)=\sum_{m}\frac{e^{i\omega_{m}\tau}}{i\omega_{m}-\Sigma_{b}(i\omega_{m})}$ with $\omega_{m}=\frac{2\pi m}{\beta}$ the Matsubara frequency and $\Sigma_{b}(i\omega_{m})$ the self-energy of the bath. In general, an imaginary part in self-energy leads to an exponential decay of the quasi-particle in real time (finite lifetime), which corresponds to an oscillation in the imaginary time after analytic continuation. Motivated by Kozii and Fu’s recent work about an effective non-Hermitian description for fermionic quasiparticles Kozii and Fu (2017), In Appendix A, we propose a similar microscopic bath model that can induce the oscillating retarded interaction as shown in Eq. (3). Moreover, by analogy with interaction in conventional crystals (e.g. Van der Waals force), we should emphasize that the time-oscillating feature in the interaction is not a necessary condition for the iTC phase. In Appendix B, we have chosen a different form of retarded interaction and numerically shown that the phase can actually exist as long as the retarded interaction is nonlocal and non-monotonic in itime with at least one minimum, whose position defines the “time” lattice constant of iTC.

Path integral Quantum Monte Carlo (QMC) simulation is a stochastic numerical method to study the equilibrium properties of quantum many-body systems based on importance sampling of the world line configurations in Euclidean spacetime. Following the algorithm proposed in Ref. Cai et al. (2014), we generalize the QMC simulation with worm-type updates Prokof’ev et al. (1998); Pollet et al. (2007) to study the open systems whose distribution functions [Eq. (1)] deviate from the Boltzmann distribution. The key point is to evaluate the integrals resulting from the retardation and include them into the QMC acceptance ratio during the updates of sampling. Notice that since the retarded interactions [Eq. (2)] only appear in the diagonal parts (under the Fock basis) of the effective actions, they do not cause extra minus sign problem. Thus our QMC simulations can be considered as numerically exact since our “system” Hamiltonian (of hard-core bosons) is also positive definite. Retarded interactions have also been implemented in other different QMC algorithms Assaad and Lang (2007); Weber et al. (2017).

III Groundstate properties.

The bosonic models with retarded interactions can be solved numerical exactly by Quantum Monte Carlo simulations (see Methods). We first focus on the ground state of the total system (system+bath). In finite size scaling, the inverse temperature is chosen as $\beta=L(2L)$ for 1D (square) lattice, and the ground state in the thermodynamic limit is approached in the limit $L\rightarrow\infty$. In the 1D case without retardation, the ground state is a quasi-superfluid ground state characterized by algebraic decaying correlation functions in both space and iTime. As the retarded interaction is switched on, spatial-temporal configurations with periodically oscillating world lines are favored since the interacting energies are lowered in these patterns. The competition between the hopping and the retarded interactions will lead to novel quantum phases.

The crystalline order can be identified from the density correlation functions in space and iTime: $S(r,\tau)=\frac{1}{L}\sum_{i}\langle(n_{i}(0)-n_{0})(n_{i+r}(\tau)-n_{0})\rangle$. We first focus on the onsite correlations with $r=0$. As shown in Fig. 5(a), for small $\alpha$, $S(0,\tau)$ still decays algebraically with $\tau$, but is accompanied by oscillations. While for sufficiently large $\alpha$ a long-range order is built up and a characteristic time period emerges, indicating a spontaneous breaking of translational symmetry in iTime, characterized by the order parameter $m^{2}_{\tau}=\frac{1}{\beta}\int_{0}^{\beta}d\tau S(0,\tau)\cos 2\pi\omega_{c}\tau.$ A finite-size scaling of $m_{\tau}$ is shown in Fig.5(c), from which we can find that for large $\alpha$, $m_{\tau}$ extrapolates to a finite value in the thermodynamic limit, indicating a crystalline order in iTime. The existence of a quantum phase transition can be further verified by the correlation length $\xi_{\tau}$ along iTime ($S(0,\tau)\sim e^{-|\tau|/\xi_{\tau}}$ for sufficiently large $\tau$), which can be derived from the structure factor (Fourier transformation of $S(0,\tau)$) Sandvik (2010). The normalized correlation length $\xi_{\tau}/\beta$ as a function of $\alpha$ for different system size and $\beta$ is plotted in Fig. 5(e), where we find a crossing point at $\alpha_{c}^{\tau}=0.23(1)$, indicating a scaling invariant critical point. In the iTC phase, $\xi_{\tau}/\beta$ linearly scales with the system size, and extrapolates to a finite value in the thermodynamic limit (inset of Fig. 5(e)).

For sufficiently large $\alpha$, it is interesting to notice that the crystalline pattern not only emerges in iTime, but also in real space, which can be identified by the equal-time correlation function $S(r,0)$ as well as its structure factor $S(Q)=\frac{1}{L}\sum_{r}e^{iQr}S(r,0)$. In the case of half-filling, the spatial crystalline pattern is a charge-density-wave (CDW) state with a periodicity twice that of the lattice, as shown in Fig. 5(b). A finite size scaling of the order parameter $m_{s}=\sqrt{S(Q=\pi)}$ is shown in Fig. 5(d). Since there is no direct nearest-neighboring repulsive interaction in our system Hamiltonian, the CDW order is formed by an effective repulsive interaction induced by the retardation: the retarded interaction (3) favors oscillating world lines in iTime, and thus requires the hard-core bosons being separated as far as possible to avoid blocking the oscillations of each other. Therefore at half-filling the tendency towards CDW order is expected in real space. Similar to the analysis above, we also plot the normalized correlation length $\xi_{s}/L$ in real space as function of $\alpha$ in Fig.5(f), where a continuous phase transition is found at $\alpha=0.245(5)$. Above this value, long-range CDW correlations emerge in real space.

To further investigate the properties of the iTC, we calculate the single particle Green’s function: $G(r,\tau)=\langle a_{i}^{\dagger}(\tau)a_{i+r}(0)\rangle$. The equal-time part $G(r,0)$ is shown in Fig. 3(a), which decays algebraically with $r$ (accompanied by an oscillation) even for large $\alpha$, indicating a quasi-superfluid order in the iTC. This order can be further verified from the superfluid density, which can be calculated as $\rho_{s}=\frac{L}{\beta}\langle W^{2}\rangle$, where $W$ is the winding number denoting the net times by which the world-lines wrap around the 1D lattice. As shown in the inset of Fig. 3(a), $\rho_{s}$ grows monotonically with $\alpha$, indicating that the retarded interactions enhance the mobility of bosons. Now we turn to the temporal part of $G(r,\tau)$, and focus on its $k=0$ component: $G(k=0,\tau)=\frac{1}{L}\sum_{r}G(r,\tau)$. Fig. 3(b) shows that $G(k=0,\tau)$ decays exponentially with $\tau$, indicating a finite single-particle charge gap. This agrees with the results of the compressibility $\kappa=\frac{\beta}{L}(\langle N^{2}\rangle-\langle N\rangle^{2})$ with $N$ the total particle number. Fig. 3(c) show that $\kappa$ decays to zero for large $\alpha$, thus the iTC is incompressible.

In summary, the iTC phase induced by the retarded interaction (3) spontaneously breaks the translational symmetries in both real space and iTime. Besides these intriguing properties, this phase is unique in the sense that it is incompressible but possess quasi-long range superfluid order, so it is fundamentally different from the conventional bosonic ground states in any closed quantum many-body systems: e.g. the (quasi-)superfluidity, Mott insulator, supersolid, Bose glass or Bose metal. The nature of this phase can be qualitatively understood in the strong coupling limit $\alpha\gg J$, where the retarded interaction dominates in the partition function. For each site, the retarded interaction (3) favors configurations with alternating $n(\tau)$ that oscillates between 0 and 1 with a frequency $2\pi\omega_{c}$, which is reminiscent of the soliton solutions in the BCS pairing model Mukhin (2009); Galitski (2010), similar states have recently been proved to be metastable in a wide class of closed quantum many-body systems Mukhin and Galimzyanov (2019). Since a world line needs to be continuous, the temporal-spatial configurations minimizing the retarded interacting action (2) are highly degenerate. As shown in the Fig. 3(d), these configurations can be understood as a synchronous movement of the bosons in iTime, in the sense that different bosons hop simultaneously towards the same direction; These synchronous movements can reduce the possibility of the collisions of world lines of different hard-core bosons, and at the same time, minimize the retarded interacting energy. In this hugely degenerate manifold, those configurations with nonzero winding number contribute to the nonzero superfluid density (e.g. the 2nd figure in Fig.3(d)). The incompressible nature can also be understood in such a picture: an extra particle added into the half-filled system will inevitably block or impede the synchronous movement of other bosons due to its hard-core nature, and thus will increase the retarded interacting energy and give rise to a finite charge gap.

Symmetries and excitations. The iTC phase simultaneously breaks the continuous and discrete translational symmetry in temporal and spatial directions. The spacial and temporal spontaneous symmetry breakings lock with each other, leaving an unbroken subgroup of intertwined space-iTime translational symmetry. The unit cell of the iTC is spanned by the primitive basis vectors $\mathbf{b}_{1/2}$ with a rhombic symmetry as shown in Fig. 3(e), as opposed to a “decoupled” direct product of vectors in $\tau$ and $x$ directions with a rectangular symmetry. Due to the fact that the time continuum is locked with the discrete spatial lattice we started with, small translational fluctuations of the atomic positions along the primitive basis vectors do not generate gapless “phonons”, but excitations with finite energy cost due to the discreteness of the spatial translational symmetry in our model. The space-iTime locking avoids the dangerous gapless modes in renormalization, so remarkably stabilizes the iTC phase.

IV Temperature effect of iTime crystal.

The most striking effect of iTC phase can be found at finite temperature. In conventional systems, the temperature dependence of a physical quantity could be either monotonic or non-monotonic with one or several peaks. In contrast, an iTC could exhibit exotic thermal behavior that its physical observable depends on temperature in an oscillating way. This intriguing property can be understood as follows: the crystal structure in iTime will perfectly fit at certain temperatures satisfying the condition that $\beta$ is a multiple of the period of iTC; otherwise, for general $\beta$ the mismatch between them will introduce defects or distortions, which impede the synchronous movement of the bosons, thus are detrimental to crystalline order in both iTime and space. As an example, in Fig. 3(f), we plot the CDW ordering parameter as a function of $\beta$, which shows an oscillating behavior with peaks at the positions of integer multiples of the iTC period. This abnormal thermal behavior can be considered as an experimental diagnostic of the iTC phase.

Extension to square lattices. Up to now, we have focused on 1D in the spatial part, while a generalization of our method to 2D square lattices is straightforward. Similar to the 1D case, we find that the temporal crystalline order also emerges for sufficiently large retarded interactions (left inset of Fig. 4). For a fixed $L$ and $\beta$ there seems to be a jump of the iTC order parameters $m_{\tau}^{2}$ as a function of $\alpha$ (Fig. 4), indicating that it is likely a 1st order quantum phase transition. However, due to the finite size effect and the limited system size in our simulation, we cannot preclude other possibilities. Another important difference between 1D and 2D is that in the latter case, there exists true superfluidity long-range order characterized by the condensate fraction $\rho_{c}=\frac{1}{L^{2}}\sum_{r}G(r,0)$. In Fig. 4 (b), we plot $\rho_{c}$ as a function of $\alpha$, where we also find a discontinuity occurs at $\alpha_{c}=0.405(5)J$. We find that even in the iTC phase, the superfluid order still persists (right inset of Fig. 4).

V Conclusion and outlook

We introduce a class of quantum many-body systems with oscillating retarded interactions in iTime and performed numerically exact QMC simulations. Our study finds an iTC phase with exotic zero-temperature properties and abnormal thermal behavior. To conclude this paper, we wish to outline some differences between our results and other related studies recently. First, our model is conveniently presented by a Euclidean action, so it is beyond the scope of Hamiltonian problems and does not contradict the “no time crystal ground state” theorem Bruno (2013); Watanabe and Oshikawa (2015). Second, even though it corresponds to a “time”-dependent problem, the “time” dependence enters our model through the relative time difference in terms of two-body interactions, in contrast to the early Floquet models whose time-dependence is explicitly through overall time  Else et al. (2016); Yao et al. (2017). This difference fundamentally changes the symmetry and a continuous symmetry model opens the wonder of what the analogue of Goldstone modes in spatial crystals, i.e, the ubiquitous and all important phonons, will be after the spontaneous breakdown of time translational symmetry. The potential of discovering “time” phonons is now open for further work. Finally, a thermodynamic physical system is found to support our “theoretical” time-crystal model. Standard statistical physics is applied to show the correspondence of our imaginary-time quantum mechanical model with retarded interaction (dependent on relative time!) to the ensemble of hard-core bosons embedded in engineered baths. This report therefore extends the study of time crystal to the Euclidean spacetime corresponding to thermal ensemble systems.

Our model can be derived from open quantum systems that have been studied extensively, such as quantum optomechanical resonators, QED cavities of interacting light and atoms, trapped ions Barreiro et al. (2011) and Rydberg atoms Bernien et al. (2017)) coupled to synthetic thermal baths, etc. The crucial ingredient to realize the proposed idea is to look for a system that has retarded interaction and allows the system-bath coupling to be easily tuned. The above main conclusions are firmly supported by numerical exact QMC simulations of bosonic models with oscillating retarded interactions in relative iTime. Our method can be straightforwardly generalized to other bosonic and spin systems with retarded interactions, in one or higher dimensions. The competition between time-oscillating interaction and quantum-thermal fluctuations is expected to give rise to other forms of nontrivial temporal-thermal orderings that are completely absent in the conventional Hamiltonian systems. A generalization to fermionic models may be even more interesting since it opens up new possibilities of finding novel “non-fermi liquid”. It is also interesting to study the real time counterpart of our model, which may be derived by integrating out a non-equilibrium bath. This may provide a new perspective for current time crystal researches since such an effective action has the full continuous (rather than discrete) symmetry in temporal translation.

Acknowledgment. We appreciate insightful discussions with F. Wilczek and V. Galitski. This work is supported in part by the National Key Research and Development Program of China (Grant No. 2016YFA0302001), NSFC of China (Grant No. 11674221 and No.11745006), Shanghai Rising-Star Program, Eastern Scholar Professor of Distinguished Appointment Program and Project of Thousand Youth Talents (Z. C.) and by AFOSR Grant No. FA9550-16-1-0006, MURI-ARO Grant No. W911NF-17-1-0323 through UC Santa Barbara, and NSF China Overseas Scholar Collaborative Program Grant No. 11429402 sponsored by Peking University (W.V. L.).

Note: During the completion of this manuscript, we notice the preprint by K.B. Efetov Efetov (2019) on thermodynamic time crystal. The subject is related but the result and models are completely different and independent.

Appendix A A microscopic model of fermionic bath

In this section, we provide a microscopic model of the bath that can induce the retarded interaction as shown in the main text. The bath is composed of fermionic quasiparticles in interacting or disordered quantum many-body systems. Recently, Kozii and Fu  Kozii and Fu (2017) proposed an effective non-Hermitian Hamiltonian to described such quasiparticles, where the effect of interaction or disorder can be encoded into the quasiparticle’s self-energy, whose imaginary part represents a finite lifetime of the quasiparticles. For a quasiparticle with a given momentum, its lifetime can be determined by the electron-electron (or electron-phonon) scattering, or the scattering by impurities. Even though a complete microscopic Hamiltonian of an interacting fermionic system should be Hermitian, the finite lifetime of the quasiparticle give rise to the non-Hermitianity of the effective Hamiltonian description of the quasiparticles.

Motivated by Ref. Kozii and Fu (2017), we propose a microscopic model for the bath of quasi-particles, which can be effectively described by a two-band fermionic non-Hermitian Hamiltonian as: $H=\sum_{i}\mathbf{\Psi}_{i}^{\dagger}\hat{H}_{i}\mathbf{\Psi}_{i}$, where $\mathbf{\Psi}_{i}=[\Psi_{i,1},\Psi_{i,2}]^{T}$, and $\Psi_{i,1(2)}$ is an annihilating operator of the bath quasiparticle at site $i$ with band index $1(2)$, and

[TABLE]

where $\omega_{d}$ is interband hopping amplitude, and we assume that the bath at each lattice site is independent of others. $\Sigma_{1(2)}(k,\omega)$ is the bath quasi-particle self-energy, whose real and imaginary part correspond to its dispersion relation and decay rate, respectively. Since the bath under consideration here is spatially local, thus we can neglect its dispersion and only keep its imaginary part. For simplicity, we only keep the constant terms in the imaginary part, thus the self energy can be approximated as $\Sigma_{1}(k,\omega)=\Sigma_{2}(k,\omega)\approx-i\Gamma$. In Ref.Kozii and Fu (2017), a microscopic derivation of $\Gamma$ from an electron-phonon interacting model has been discussed.

By performing the standard procedure of path integral for fermions N.Nagaosa (1999), the partition function can be expressed as $Z=\int\mathfrak{D}[\bar{\psi}]\mathfrak{D}[\psi]e^{-S[\bar{\psi},\psi]}$, where $\psi_{i}$ is the Grassmann number satisfying $\Psi|\psi\rangle=\psi|\psi\rangle$ with $|\psi\rangle$ the fermonic coherent state (we neglect the site and band indexes here for simplicity), and $\mathfrak{D}[\psi]$ denotes the functional integral over $\psi(\tau)$. The effective action of the bath can be written as:

[TABLE]

where we performed the Fourier transformation: $\psi(\tau)=\frac{1}{\beta}\sum_{m}e^{-i\omega_{m}\tau}\psi(i\omega_{m})$; $\bar{\psi}(\tau)=\frac{1}{\beta}\sum_{m}e^{i\omega_{m}\tau}\psi(i\omega_{m})$, and $\omega_{m}=\pi T(2m+1)$ is the Matsubara frequency.

We further assume that the system particles (i.e., hard-core bosons) couple to the bath quasiparticles locally via its density operator $n_{i}$: $H_{sb}=\sum_{i,\sigma}n_{i}(\alpha_{i,\sigma}+\alpha_{i,\sigma}^{\dagger})$. By integrating out the bath degrees of freedom, one can obtain the bath-induced retarded interaction between the system particles: $S_{ret}=\sum_{m}D(\omega_{m})n(\omega_{m})n(-\omega_{m})$, where the kernel function takes the form of $D(\omega_{m})\propto\frac{1}{\det[-i\omega_{m}\hat{1}-\hat{H}]}=\frac{1}{(\omega_{m}-\Gamma)^{2}+\omega_{d}^{2}}$ with $\Gamma$ and $\omega_{d}$ defined the same as in Eq.(4).

By performing the Fourier transformation again and take the zero temperature limit, we obtain $S_{ret}=\int d\tau d\tau^{\prime}D(\tau-\tau^{\prime})n(\tau)n(\tau^{\prime})$, with

[TABLE]

which is the kernel function shown in the main text.

Appendix B An alternative model for iTC

In the main text, we choose an oscillating decaying retarded interaction in imaginary time. However, this oscillating decaying feature of the retarded interaction is not a necessary condition for the iTC phase, which is found to exist in a much broader circumstances. That is, it can exist as long as the retarded interaction is non-monotonic in $|\tau-\tau^{\prime}|$, and has at least one minimum whose position defines the lattice constant of iTC along $\tau-$direction, similar to the conventional crystals in space.

In this section, we will numerically verify this point by showing the existence of the iTC phase with a different retarded interaction with the kernel function:

[TABLE]

The physics with such a complex form of the kernel function can be easily understood by making analogy with the conventional spatial crystals, where the interacting objects are point-like. However, in the iTC, the interacting objects are continuous lines (world-line). Thus to make such an analogue, one needs first to do a dual transformation, given world line configuration can be mapped to a particle picture, i.e., the gas of instantons and anti-instantons (also referred in literature as to kinks and anti-kinks) with the same interacting energy N.Nagaosa (1999):

[TABLE]

where $m(\tau_{i})=\pm 1$ indicates a kink or antikink with the position $\tau_{i}$. The interacting energy $V(\tau-\tau^{\prime})$ between a pair of kinks can be directly evaluated by performing integral over $\tau$ and $\tau^{\prime}$, which gives rise to

[TABLE]

Obviously, for those point-like objects (kinks and anti-kinks), this kind of interaction favor a crystalline pattern in $\tau$-direction with lattice constant $\tau_{0}$, which corresponds to the iTC pattern in the world-line configurations.

To incorporate the quantum fluctuations in this picture, we perform QMC simulations of the model that is similar to what is discussed in the main text (hopping+retarded interaction), the only difference being that the kernel function replaced by Eq. (9). In Fig. 5(c), we plot the onsite temporal density-density correlations, which shows a long-range crystalline order in $\tau$-direction for a strong retarded interaction ($\alpha=6J$), indicating a spontaneous translational symmetry breaking in imaginary time. The lattice constant of the iTC is determined by the characteristic itime length (the minimum position) of the retarded interaction ($\tau_{0}$). This example indicates that the oscillating feature of the retarded interaction is not necessary for the iTC phase, and hence proves that the phase can be observed in much boarder situations.