Subspace-thermal discrete time crystals from phase transitions between different n-tuple discrete time crystals
Hongye Yu Tzu-Chieh Wei

UTF8gbsn





C. N. Yang Institute for Theoretical Physics,
State University of New York at
Stony Brook, Stony Brook, NY 11794-3840, USA
Department of Physics and Astronomy, State University of New York at
Stony Brook, Stony Brook, NY 11794-3800, USA


C. N. Yang Institute for Theoretical Physics,
State University of New York at
Stony Brook, Stony Brook, NY 11794-3840, USA
Department of Physics and Astronomy, State University of New York at
Stony Brook, Stony Brook, NY 11794-3800, USA



             

§ ABSTRACT 
We propose a new Floquet time crystal model that responds in arbitrary multiples of the driving period. Such an n-tuple discrete time crystal is theoretically constructed by permuting spins in a disordered chain and is well suited for experiment implementations. Transitions between these time crystals with different periods give rise to a novel phase of matter that we call subspace-thermal discrete time crystals, where states within subspaces are fully thermalized at an early time. However, the whole system still robustly responds to the periodic driving subharmonically, with a period being the greatest common divisor of the original two periods. Existing theoretical analysis from many-body localization fails to understand the rigidity of such subspace-thermal time crystal phases. To resolve this, we develop a new theoretical framework from the perspective of the robust 2π/n quasi-energy gap. Its robustness is analytically proved, under a reasonable conjecture, by a new perturbation theory for unitary operators. The proof applies beyond the models considered here to other existing discrete time crystals realized by kicking disordered systems, thus offering a systematic way to construct new discrete time crystal models. We also introduce the notion of DTC-charges that allow us to probe the observables that spontaneously break the time-translation symmetry in both the regular discrete time crystals and subspace-thermal discrete time crystals. Moreover, our discrete time crystal models can be generalized to higher spin magnitudes or qudits, as well as higher spatial dimensions. 



Subspace-thermal discrete time crystals from phase transitions between different n-tuple discrete time crystals
    Tzu-Chieh Wei
    Accepted Sep 3 2024 to ApJ Letters

===============================================================================================================






§ INTRODUCTION


A time crystal was originally proposed <cit.> by Wilczek as a phase of matter that spontaneously breaks the continuous time translation symmetry (τ). Despite that breaking it in equilibrium was later proved <cit.> not possible, various models <cit.> that spontaneously break the discrete time translation symmetry have been proposed. These discrete time crystals (DTC) are realized in periodic-driven (Floquet) systems, with some observables responding with periods that are multiples of the driving period. Comprehensive reviews of  time crystals can be found in Refs. <cit.>.

A major feature and challenge of a DTC system is its stability of subharmonic response against perturbations. A typical Floquet system's state can be thermalized to infinite temperature, whereas a DTC system contains a mechanism that robustly breaks ergodicity, thus preventing thermalization. The robust ergodicity breaking can be induced in several different settings, such as many-body localization <cit.> (MBL), Floquet prethermalization <cit.>, and others <cit.>.


Recently, various DTC phases have been observed in experiments <cit.>. Despite several higher n-tuple DTC models <cit.> having been theoretically proposed, there are rarely experimentally realized DTC models that exhibit periods higher than doubling. In this work, we propose a new model for an arbitrary integer n-tuple DTC feasible for experimental realizations, based on permutations on disordered one-dimensional spin-1/2 chains. For n=2, our model is reduced to a recently proposed swapping Floquet DTC model <cit.>. As our n-DTC model is based on “kicking" a disordered spin chain, it shares similar properties to existing MBL-DTC models, such as kicked Ising models <cit.>.

A remarkable property of an MBL-DTC system is the robust
2π/n gap in its quasi-spectrum (note we absorb the gap dependence on driving frequency into the effective Hamiltonian for simplicity). For typical MBL-DTC systems <cit.> with size L, when local period-T perturbations are added, the deviations from the 2π/n are exponentially <cit.> small e^-O(L) in L. Most proofs for this property rely on corollaries from the MBL <cit.>, such as the existence of a quasi-local unitary relating the perturbed eigenstates and unperturbed ones. In this work, we offer a new proof from the perturbative perspective that does not rely on the system being MBL. In our proof, a specific form of eigenstates with minimum quasi-spectral gap not closing severely, together with a conjecture assuming effects from globally different states can be neglected, suffice for the robust 2π/n gap against k-local period-T perturbations (with k finite). Such properties appear not just in MBL systems, hinting that MBL may not be a necessary condition for this type of DTC. 



Indeed, during the phase transition among our n-DTCs, we identify a new form of DTC that is fully thermalized but still exhibits robust subharmonic oscillations for some observables. The thermalization saturates quickly but is restricted in subspaces of the whole Hilbert space and remains there for an expentially long time ∼ e^O(L). Such restriction is due to the symmetries of our model in the unperturbed case but is robust against perturbations even if they break the symmetries. The subharmonic oscillations can be viewed as a thermal state restricted in one subspace, cyclically jumping to another subspace without mixing with other subspaces, which can be observed by measuring the symmetry charges that divide those subspaces. We call this new phase a subspace-thermal DTC (ST-DTC). The ST-DTC is different from MBL-DTC: In MBL-DTC <cit.>, the system retains the memory of its initial state forever, whereas, in ST-DTC, the memory is lost at the early time due to the thermalization in subspaces. The ST-DTC is also different from prethermal-DTC: in the prethermal-DTC, the lifetime scale exponentially <cit.> with the driving frequency e^O(ω_0/J) (J is the local energy scale), whereas the lifetime of ST-DTC is not directly related with ω_0, but with the system size e^O(L) similar to MBL-DTC. In addition, the robust subharmonic response in the prethermal-DTC only holds for low-energy initial states <cit.>, whereas in the ST-DTC it holds for almost every initial state. 

Our perturbative analysis for the robustness of 2π/n gap can explain the existence of the ST-DTC well. In addition, the analysis offers a systematic way to construct potential n-DTC models, which is much broader than MBL systems, and is not limited by the dimension or interaction range.
Moreover, as shown in previous works <cit.>,  a non-DTC phase is inevitable during the phase transition between different DTCs. In this work, we show that the transition can be done through an ST-DTC instead of a fully thermalized phase. 


The rest of the paper is organized as follows: In Sec. <ref>, we briefly outline the proof for robust 2π/n gap in DTC systems and leave details to Appendix. <ref> and <ref>, where we introduce a new unitary perturbation theory and apply it to the Floquet operator. In Sec. <ref>, we present the n-tuple DTC model by permuting a disordered 1-D spin-1/2 chain. In Sec <ref>, we show that the subspace-thermal DTC emerges, oscillating with period n_G≡(n_1,n_2), during the phase transition between n_1-DTC and n_2-DTC. With perturbation, the thermalization is robustly confined in the perturbed subspaces due to the emergent symmetries, and one can observe the robust subharmonic oscillations by measuring DTC-charges. In Sec. <ref>, we briefly discuss the potential experimental realization of our models, such as in current noisy intermediate-scale quantum (NISQ) computers. We conclude  in
Sec. <ref>.




§ SKETCH OF PROVING THE ROBUST 2Π/N GAP OF DTCS

In this section, we briefly outline the proof for the robustness of 2π/n gap of disordered DTCs against period-T local perturbations, and refer the readers to Appendix. <ref> and <ref> for further details.

For a Floquet system driven by a periodic Hamiltonian H_0(t+T)=H_0(t), the evolution of the system is governed by U_0(t)≡𝒯 e^-i ∫_0^t H_0(t') d t' (where 𝒯 denotes time ordering), whose stroboscopic properties can be fully described by the Floquet operator U_F^0≡ U_0(T). We define the quasi-energy ε of the U_F^0 by writing the eigenvalues of U_F^0 as e^-i ε. We are interested in how the quasi-spectrum of U_F^0 changes in the presence of perturbations.

Supposing a periodic perturbation λV̂(t+T)=λV̂(t) is added, the new Floquet operator becomes

    U_F(λ)=𝒯exp(-i ∫_0^T(H_0(t)+λV̂(t)) d t).

The perturbed Floquet operator U_F(λ) can be alternatively <cit.> written as
U_F(λ)≡ U_F^0 U_λ,
where U_λ= 𝒯exp(-i ∫_0^T U_0^†(t) λV̂(t) U_0(t)d t)=: e^-iλ V and we have defined an effective Hermitian term λ V for U_λ. 



To analytically study the stability of an L-site DTC system, we must accurately obtain the perturbed quasi-energies up to O(L)-th order. Existing perturbation theories <cit.> either have convergence problems for unitary operators or are hard to analyze in arbitrarily high order (more detailed discussions are in Appendix <ref>). Thus, we require a new perturbation theory for unitary operators that can be explicitly calculated to arbitrary order and rigorously converged for small enough perturbations. 





 §.§ Unitary perturbation theory


The goal of the unitary perturbation theory is to find quasi-energies and eigenstates for the perturbed unitary operator U_F(λ)≡ U_F^0 U_λ. The unitaries can be generic, not necessarily the Floquet operators for periodic-driven systems that we focus on here. For convenience, we define effective Hamiltonians for these unitary operators e^-i H_F(λ) ≡ U_F(λ) , e^-i H_F(0) ≡ U_F^0 and e^-i λ V ≡ U_λ (the last was already defined above).  Then we obtain the relation

    e^-i H_F(λ) =e^-i H_F(0)  e^-i λ V .

Compared to the standard definition of effective Hamiltonians for Floquet systems, we absorb the T in effective Hamiltonians for simplicity. Suppose we are interested in an eigenstate |ψ_i⟩ with quasi-energy E_i of U_F(λ), which is close to the eigenstate |i⟩ with quasi-energy ε_i of the unperturbed unitary U_F^0. This is equivalent to solving

    e^i [E-H_F(0)]|ψ⟩= e^ i λ V |ψ⟩.

Directly applying conventional perturbation theories to this equation can lead to convergence problems. To remedy this, we introduce a novel trick by rewriting the equation as

    tan(E-H_F(0)/2) |ψ⟩= tan(λ V/2) |ψ⟩,

which is equivalent to Eq. (<ref>), as they share exactly the same solutions [If we also regard the infinity as one point for tanπ/2, we can construct a bijective map between e^i A and tanA/2.]. Since eigenvalues of tanλ V/2 are small when λ→ 0, we can apply conventional perturbation theories to Eq. (<ref>) without convergence problems. Below, we use a similar procedure to the Brillouin-Wigner perturbation theory <cit.>. To begin with, we look for a solution |ψ_i⟩ close to the original eigenstate |i⟩ and define the projector P_i≡|i⟩⟨i|, as well as  the associated matrix R_i ≡∑_m≠ i(E-ε_m/2)|m⟩⟨m|.
Applying R_i to the Eq. (<ref>), we obtain (1-P_i)|ψ_i⟩=R_i tan(λ V/2) |ψ_i⟩,
which leads to |ψ_i⟩=1/1-R_i tan(λ V/2)P_i|ψ_i⟩, which can be expanded as follows,

    |ψ_i⟩=(1+R_i tan(λ V/2)+⋯) P_i|ψ_i⟩,

with the convergence condition tanλ V/2 < tanE_i-ε_m/2, for m≠ i. Using the relation ⟨i|  tan(E-H_F(0)/2)| ψ_i|=⟩⟨i|tan(λ V/2) |ψ_i⟩, we obtain the expansion for the perturbed quasi-energy E_i:

    tan(E_i (λ)-ε_i/2)= ⟨i|tanλ V/2|i|+⟩
    ∑_m_1 ≠ i⟨i|tanλ V/2| m_1|⟨%s|%s⟩⟩ m_1|tanλ V/2| i/tanE_i (λ)-ε_m_1/2+…,
 
where we see that the j-th term on the RHS is of order O(λ^j). Note that the perturbed energy E_i appears on both sides of the equation, so it is hard to solve it directly. Instead, one can solve it recursively. Let E_i^(j)(λ) be the j-th order approximation of E_i(λ) (i.e., |E_i^(j)(λ)-E_i(λ)|∼ o(λ^j) ), then one can obtain E_i^(j+1)(λ) by setting E_i(λ)=E_i^(j)(λ) on the RHS of Eq. (<ref>), and repeat the procedure recursively. 

So far, we have introduced a new perturbation theory for unitary operators based on the “tangent trick" applied to the Brillouin-Wigner perturbation theory. This theory works for general quasi-spectral perturbation by unitary operators, including Floquet operators of periodic-driven systems, which paves the way for our perturbative analysis for DTC models. We present in Appendix <ref> more rigorous and detailed discussions, along with the unitary perturbation theory near degeneracy.






 §.§ Differences in quasi-energy corrections are O(L)-order

In this section, we discuss a class of DTC systems with solvable points, where the Floquet operator U_F exactly transforms one state |χ_m⟩ to another globally different orthogonal state |χ_m+1⟩, in a cyclic way. In addition, we also require the system to exhibit a “local spectral gap”: states that differ locally have non-vanishing quasi-energy gaps, which are usually induced by random interactions. We remark that these conditions are satisfied by most existing MBL-DTC models <cit.>. 

For an n-DTC, the state returns to itself after n periods of driving. One can verify that if U_F acts in this way for the whole Hilbert space, all {|χ_m⟩} form an orthogonal basis, and thus, we can divide the Hilbert space into dynamically disjoint subspaces. Within the α-th subspace sector that has period n, the action of U_F is

    U_F|χ_α,m⟩= e^-i φ_α,m+1|χ_α,m+1⟩,

where 1≤ m ≤ n for n-DTC, and |χ_α,n+1⟩≡|χ_α,1⟩. In such a case, the quasi-energies ε_α,j and eigenstates |Φ_α,j⟩ can be exactly solved as follows, 

    |Φ_α,j⟩   =1/√(n)∑_m=1^n e^i (θ_α,m+m j2π/n )|χ_α,m⟩,
    ε_α,j   =1/n∑_m=1^n φ_α,m + j 2π/n,

where θ_α,m=1/n∑_i=1^n (m-i n) φ_α,i and  U_F|Φ⟩=e^-iε|Φ⟩. Thus, the quasi-energies are separated with exact 2π/n spacing, and the eigenstates are all composed of an equal-weight superposition of |χ_α,m⟩'s, with different relative phases among them for different eigenstates. When |χ_α,m⟩'s are globally different from each other, the eigenstates {|Φ_α,j⟩} are generalized Schrödinger's cat states. We define l_α being the minimum distance among states in the α-th subspace, and define |χ_α,m⟩'s to be globally different if l_α is O(L). The distance can be defined in various way, e.g., Hamming distance for Z-basis states. For general states, one can define the distance of |χ_α,m⟩ and |χ_α,m'⟩ being l if ⟨χ_α,m|𝒱|χ_α,m'|$⟩ is non-zero only if𝒱contains interactions involving at leastlsites. In DTC systems,l_αis typicallyO(L).

 A conjecture. Before proceeding to include perturbations, we introduce a conjecture we need in the proof. Suppose when a finite-body period-Tperturbation of strengthO(λ)is added to a (sufficiently strong) disordered DTC system defined above, the original eigenstate|Φ_α,i(0)⟩is perturbed to|Φ_α,i(λ)⟩. For a sufficiently small but non-vanishing perturbation strengthλ>0, we conjecture that there exists a finite and sufficiently largeμ∼log1/λ, so that for distancel≥l_α/4andl_α>γL(withγconstant), one can find an operator𝒱^(<l)involving at mostl-body interactions, satisfying

    |Φ_α,i(λ)⟩-
     𝒱^(<l)|Φ_α,i(0)⟩<C e^-μ l,

where0<C<1is a constant. We also give an intuitive argument for its validness in disordered DTC, in Appendix <ref>. We remark that the conjecture is a weaker version of the “local perturbations perturb locally” (LPPL) <cit.>,  an analog of a similar property in gapped systems <cit.>. The LPPL can be equivalently written in a similar form as Eq. (<ref>), with𝒱^(<l)limited to short-range interactions, whereas our conjecture has no restriction on the interaction length. Moreover, LPPL requires that Eq. (<ref>) holds forlgreater than a finite number, whereas we only require that Eq. (<ref>) holds forlgreater than a portion of system sizeL. The LPPL is usually satisfied and well accepted in MBL-DTC systems <cit.> with local spectral gaps and served as a stepping stone to prove their stability. 


When applying our unitary perturbation theory to disordered DTC systems, eigenstates (<ref>) ofU_F^0have a special property against local perturbations: the magnitudes of terms in the quasi-energy perturbation (<ref>) for different eigenstates in the same subspace ({|Φ_α,j⟩}with differentj's) are almost identical, with deviations at the order ofλ^O(L). The terms can have relative phases of the forme^i c 2π/n, withcbeing an integer, but they can be eliminated if states involved are not too far away from{|Φ_α,j⟩}. In addition, using the above conjecture, contributions from globally distinct states to the quasi-energy changes can be shown to be exponentially small. Furthermore, the relation holds up toO(L)-order perturbation theory. Thus, the overall quasi-energy changes for states in the same sector are almost identical, with deviations at the order ofλ^O(L). Thus, the difference in the perturbed quasi-energiesE_α,j+1-E_α,jremains almost2π/n, withλ^O(L)corrections, despite that the change from original quasi-energyE_α,j-ε_α,jcan be still as large asO(λ). We leave the rigorous statements and details in Appendix <ref>.

Moreover, the analytical analysis also works for many existing model (see discussions in Appendices <ref> and <ref>), and offers a systematic way to construct newn-tuple DTC models. To do this, one can firstly consider a solvable model, where the Floquet operatorU_Fexactly transforms one state|χ_m⟩to another globally different orthogonal state|χ_m+1⟩, in a cyclical way with periodn>1. Then one can add disorders that do not mix|χ_m⟩'s, to break the quasi-spectral degeneracies. Note that such disorder may not always lead to MBL, but is usually sufficient to construct a DTC model. According the proof, the2π/ngap for such model is exponentially robust against local perturbations, which leads to a DTC phase.






§ N-TUPLE DTC MODEL FROM PERMUTATION

We define the period-nTdiscrete time crystal on a 1-D length-Lspin-1/2chain with open boundary condition [different boundary conditions will not bring qualitative differences in this model], whose evolution is governed by a period-THamiltonianH(t+T)=H(t). Within one period (T=t_1+t_2+t_3), the Hamiltonian is (see Fig. <ref>)

    H^(n)(t)= 
        H_1^(n),     for  0 ≤ t<t_1 ,
     
        H_2^(n),     for  t_1 ≤ t<t_1+t_2,
    
        H_ int,     for  t_1+t_2 ≤ t<t_1+t_2+t_3,

whereH_1^(n)andH_2^(n)both consist of swap gates acting as permutations and the interaction HamiltonianH_intincludes disordered couplings and fields 

    H_ int=∑_i<j J_ijσ_i^z σ_j^z+∑_i=1^L h_i^z σ_i^z+∑_i=1^L ϵ_i^x σ_i^x,

where the parametersh_i^z∈[0,2h̅^̅z̅]are randomly chosen from uniform distributions, andϵ_i^x∈[-ϵ^x,ϵ^x]'s are small uniformly random perturbations. We choose the couplingJ_ijto be a disordered power-law interaction obtained by

    J_ij=1/ℒ_L, κJ̃_i, j/|i-j|^κ,

whereJ̃_ij∈[1/2J̅,3/2J̅]is uniformly and randomly chosen, andℒ_L, κis a coefficient to make the energy extensive <cit.>
    ℒ_L, κ= 1     if κ⩾ 1, 
    ln L     if κ=1, 
     L^1-κ    if κ<1.

We remark that power-law interaction is not the only choice [For example, one can randomly choose J_ij uniformly for |i-j|≤ 2n, and set J_ij=0 for |i-j| > 2n, where n is the length of the spatial permutation units.] for our model being DTC; however, nearest-neighbor interactions are not sufficient <cit.> to break all the degeneracies in the quasi-spectrum. For the even periodn≥2, theH_1^(n)andH_2^(n)are






    H_1^(n)   =   π/2 t_1(1-ϵ_1)∑_i=1^L/2-11/2(σ̂_2i-1·σ̂_2i-1),
    
        H_2^(n)   =   π/2 t_2(1-ϵ_2)( ∑_i=1^L/2-11/2 (σ̂_2i·σ̂_2i+1-1).
          +.∑_i=1^L/n-11/2(1-σ̂_ni·σ̂_ni+1) ),

whereσ̂_i·σ̂_j=σ_i^x σ_j^x+σ_i^y σ_j^y+σ_i^z σ_j^zis the Heisenberg type of interaction. For odd periodn≥3, theH_1^(n)andH_2^(n)are

    H_1^(n)   =   π/2 t_1(1-ϵ_1)( ∑_i=1^L/2-11/2(σ̂_2i-1·σ̂_2i-1)+∑_i=1^L/(2n)-11/2(1-σ̂_(2i-1)n·σ̂_(2i-1)n+1) ),
    
        H_2^(n)   =   π/2 t_2(1-ϵ_2)( ∑_i=1^L/2-11/2( σ̂_2i·σ̂_2i+1-1)+∑_i=1^L/(2n)-11/2(1- σ̂_2ni·σ̂_2ni+1) ),

whereϵ_1andϵ_2are close to 0, implying small imperfections. For convenience, in all the following discussions, we assume the system lengthLis divisible [This assumption is for convenience and is not essential for DTC properties. If L is not divisible by n or 2n, one can always set the H_1^(n) and H_2^(n) to be zero at the remainder sites.] bynfor evenLor2nfor oddL, respectively.  The Floquet evolution operator for a single period is thenU_F^(n)(ϵ̂)=U_intU_2^(n)U_1^(n)≡e^-i H_int t_3e^-i H_2^(n)t_2e^-i H_1^(n)t_1.


To grasp thenTperiods of such evolution, it is helpful to look into the solvable, unperturbed case (ϵ_i^x=ϵ_1=ϵ_2=0). We define the unperturbed Floquet operator asU_F^(n)≡U_F^(n)(ϵ̂=0)for simplicity. Firstly, we notice that the unitariesU_1^(n), U_2^(n)are groups of swap gate, asSWAP_i,j≡1/2(I+σ̂_i·σ̂_j)=e^-i π/4 (σ̂_i·σ̂_j-1) . Thus,U_1^(n)andU_2^(n)effectively perform swap gates on odd and even bonds, respectively (see Fig. <ref>). The gate structure is periodic in space, with each smallest spatial permutation unit consisting ofnspins [For odd period n, the minimum permutation unit consisting of 2n spins, but the left half and the right half are symmetrical (see Fig. <ref>) and do not have qualitative differences in analytical analysis.]. Within each unit, it can be proved (see Appendix <ref>) that such spin permutation has a period ofn. In addition, theH_intonly adds a phase forZ-basis states. Thus, the whole system follows annT-period evolution forZ-basis initial states. One may notice that someZ-basis states can oscillate with a smaller periodk<n, wherekdividesn, i.e.,k|n. Taken=L=4, for example. While all of the states satisfy(U_F^(4))^4 |ψ^(4)⟩=|ψ^(4)⟩, some of the states (e.g.,|↑↓↓↑⟩) also satisfy(U_F^(4))^2|ψ^(2)⟩=|ψ^(2)⟩and the ferromagnetic states (e.g.,|↑↑↑↑⟩) satisfyU_F^(4)|ψ^(1)⟩=|ψ^(1)⟩, where we label thek-period state as|ψ^(k)⟩. We prove that such lower-period states are exponentially rare, and most of them also have robustk-period subharmonic oscillation and will not be mixed withn-period states even if perturbations are added (see discussions in Appendices <ref> and <ref>). Thus, we will mostly considern-period states in the rest of the discussion in the main text.

When local perturbations are added, we expect that the subharmonic response will be preserved. However, different initial states may have different robustness. As shown in Sec. <ref>, the robustness is determined by the minimum distance among states in the same subspace. If ann-period state has ak-period unit (withk<n), then this unit does not contribute to the distance. Thus, the robustness only relies on the number ofn-period units for states inα-th subspace, which is denoted byl_αin our model. Note that thel_αhere may differ from the definition for general DTC models in Sec. <ref> by a factor, but will not change the overall scaling.  In Appendix <ref>, we prove that for the vast majority of the Hilbert space,l_αisO(L/n). Thus, for most states, the subharmonic oscillation will maintain up to ane^O(L/n)time scale. In the thermodynamic limit, the time scale approaches infinity.


In the following discussions on numerics, we set, for convenience,J̅=4, h̅^̅z̅=12, κ=0.5, and fix the evolution time for each Hamiltoniant_1=t_2=1/2andt_3=1. For simplicity, we set allϵproportional to a parameterλ, controlling the overall perturbation strength, whereϵ^x=λ, ϵ_1=0.9λand ϵ_2=1.1λ.








 §.§ Symmetries and DTC-charges in solvable case


In the unperturbed case (allϵ=0), the Floquet operatorU_F^(n)exactly permutes spins to different sites. For convenience, we relabel them asσ_i,1≡σ_(i-1)n+1and

    σ^z_i,j≡(U_F^(n))^j-1∘ (σ^z_i,1),

whereU∘σ≡U^†σUand we useU_F^(n)to denote the unperturbed Floquet operator forn-DTC throughout the paper. The actions ofU_F^(n)on them are 

    U_F^(n)∘ (σ^z_i,j)=σ^z_i,j+1,
    
        [U_F^(n)]^n ∘ (σ^z_i,j)=σ^z_i,j.





A slightly different picture arises from the perspective of states (see Fig. <ref>b). ForZ-basis states , the Floquet operator cyclically transforms oneZ-basis state|z⟩to anotherU_F^(n)|z_m⟩=e^-i E_int(z_m+1)|z_m+1⟩, whereE_int(z)is the energy of theH_intfor the state|z⟩. Thus, the whole Hilbert space is divided into exponentially many dynamically disjoint fragments, with each subspace consisting ofkZ-basis states, wherek(with k| n) is the period of them. Forn-period states in theα-th subspace, the corresponding eigenstates|ϕ_α,j^(n)⟩and quasi-energiesε_α,j^(n)are (see Eq. (<ref>) and Appendix <ref> for more details)

    |ϕ_α,j^(n)⟩   =1/√(n)∑_m=1^n e^i (θ_α,m^(n)+m j2π/n )|z^(n)_α,m⟩,
    ε_α,j^(n)   =1/n∑_m=1^n E_ int(z^(n)_α,m) + j 2π/n,

whereθ_α,m^(n)=1/n∑_i=1^n (m-in) E_int(z_α, i^(n)). We thus see that the quasi-energies within each subspace are separated by exactly2π/nspacing. Similarly, one can obtain2π/kspacing fork-period states.

With the exact2π/k(wherek|n) quasi-spectral separations, one can identify a symmetry groupℤ_nfor the system using the group generator

    S   ≡   ∑_α∑_j=1^k_α e^-i j2π/k_α|ϕ^(k_α)_α,j⟩⟨ϕ^(k_α)_α,j|,

wherek_α(which dividesn) is the dimension of the subspaceα.
By construction, it is easy to see thatS^n=1and thatScommutes withU_F^(n). The non-trivial relation of theU_F^(n)andSis 

    .U_F^(n)|_α= .e^-i E_α/k_α S|_α,

whereE_α≡∑_m=1^k_α E_int(z^(k_α)_α,m)andU|_αis the operatorUrestricted in theα-th subspace. This indicates that, when acting within theα-th subspace, the Floquet operatorU_F^(n)has exactly the same effects as theℤ_nsymmetry generatorS, up to a global phasee^-i E_α/k_α. This immediately leads to[U_F^(n)]^n |_α= e^-i E_αn/k_αbeing a pure phase in theα-th sector. 



Thus, allZ-basis states within the same sector have the same quasi-energies w.r.t. then-th power ofU_F^(n), i.e.,[U_F^(n)]^n, which gives rise to more symmetries thanU_F^(n); see Eq. (<ref>). We are interested in symmetries that only appear in[U_F^(n)]^nbut not inU_F^(n)so that they naturally exhibit subharmonic oscillations underU_F^(n). As expectation values of symmetry operators are conserved charges for[U_F^(n)]^n, but they oscillate underU_F^(n)within eachn-period evolution,  so we also refer to them as  DTC-charges. One can verify that allσ^z's , along with most their linear combinations, are DTC-charges ofU_F^(n). By measuring those DTC-charges, we can observe the subharmonic oscillations, which can be seen by Eq. (<ref>).



 §.§ Emergent symmetries and robust subharmonic oscillation against local perturbations


We now show that the subharmonic oscillations remain robust against local perturbations. In the unperturbed case, forZ-basis states in theα-th subspace, their components within one spatial unit are distinct if they are oscillating with period-nin the unit. We usel_αto count the number of period-nunits in the whole chain. From the discussion in Sec. <ref>, the correction for the2π/nquasi-spectral separation isΔ^(n) ∼λ^O(l_α). Thus, the time-crystalline structure is exponentially robust to local perturbations and becomes exact in the thermodynamic limit. We have also numerically simulated the quasi-spectral gap dependence withλfor 4-DTC (see Fig. <ref>). With sorted quasi-energies in full space{ε_i}andα-th subspace{ε_α,i}, we define the quasi-spectral spacing and the deviation from2π/nquasi-spectral separation, respectively

    Δ^(0)_i   =ε_i+1-ε_i, 
    Δ^(n)_α,i   =|ε_α,i+1-ε_α,i-2 π / n|.

We then average them as 

    log_10Δ^(0)   ≡𝐚𝐯𝐠.[1/𝒩∑_ilog_10Δ^(0)_i], 
    log_10Δ^(n)   ≡𝐚𝐯𝐠.[max_α,i(log_10Δ^(n)_α,i)],

where the symbol𝐚𝐯𝐠.[⋯]means averaging over disorders,𝒩counts the number of gaps, and the maximum is only picked fromn-period sectors. Sincel_αaffects the scaling behavior of2π/ngap, in numerical simulations, we average theΔ^(0)_iover the whole Hilbert space, but only pick the maximum ofΔ^(n)_α,iin a subspace with the samel_α's for the consistency of the scaling. In the unperturbed caseλ=0, the subspace is chosen in theZ-basis, where each spatial unit contains exactly one spin down. Forλ>0, the subspace is induced by applying𝒱(λ)to the unperturbed subspace, where𝒱(λ)relates the unperturbed eigenstates to perturbed ones. Forn=4, we see from Fig. <ref> that the average quasi-spectral spacing is almost fixedΔ^(0) ∼O(1/2^L),  whereas the average changes of the2π/n=π/2(forn=4) quasi-spectral separation show a clear exponential scalingΔ^(4)∼O(λ^L/4). The scaling can also be analytically obtained (see Appendix <ref>).



With robust2π/ngaps, the symmetryℤ_nand the DTC-chargesσ^z's are perturbed into approximate symmetries and DTC-charges, where the deviations from their exact definitions are exponentially small. Suppose that the original eigenstates ofU_F^(n)(0)are{|ϕ_α,j⟩}and the perturbed eigenstates ofU_F^(n)(λ)are{|ψ_α,j(λ)⟩}, respectively. We define a unitary operator𝒱that connects the original eigenstates to the perturbed ones:𝒱(λ)|ϕ_α,j⟩=|ψ_α,j(λ)⟩. Defining|z̃(λ)⟩≡𝒱(λ)|z⟩, one can verify that (see Appendix <ref> for details)

    U_F^(n)(λ)|z̃_α,m(λ)⟩=e^-i Ẽ_α,m+1|z̃_α,m+1(λ)⟩+λ^O(l_α).

Theℤ_nsymmetry is then perturbed toS̃(λ)≡𝒱(λ)S𝒱^†(λ).
One can verify thatS̃(λ)is still an exactℤ_ngenerator satisfyingS̃^n=1, and remains an exact symmetry ofU_F^(n)(λ). However, the relation Eq. (<ref>) no longer holds exactly; instead, it becomes an approximate relation

    .U_F^(n)(λ)|_α= .e^-i Ẽ_α/k_αS̃|_α + λ^O(l_α),

where theα-th subspace is induced by the perturbed eigenstates|ψ_α⟩'s. These results are a direct consequence of the fact that the2π/k_αquasi-energy separation is perturbed at the order ofλ^O(l_α)(see Appendix <ref> for details). Therefore,U_F^(n)(λ)approximately acts like aℤ_nsymmetry generator in the new basis, up to aλ^O(l_α)correction.

This leads to the perturbed DTC-chargesτ_i^z's ofU_F^(n)(λ), which can be obtained from similar transformations onσ^z's. We can defineτ_i^a≡𝒱(λ)σ_i^a𝒱^†(λ)witha ∈{x,y,z}. It is easy to check that |z̃⟩'s are eigenstates ofτ_i^z's. From Eq. (<ref>), we see thatτ^z's become approximate emergent symmetries of[U_F^(n)(λ)]^n, and act like exact DTC-charges in the thermodynamic limit, i.e., 

    [(U_F^(n)(λ))^n,τ^z]|_α   ∼λ^O(l_α),
    
        U_F^(n)(λ) ∘ (τ^z_i,j)|_α   =τ^z_i,j+1|_α+λ^O(l_α).

In addition, we expect𝒱(λ)to be close to identity whenλis small; we can then expandσ^z's inτ-basis

    σ^z_i,j=𝒱(λ)^†τ^z_i,j𝒱(λ)=[1-O(λ)]τ^z_i,j+O(λ).

Thus, by measuringσ^z's, we will observe that their major components (dominated byτ^z's) have robust subharmonic oscillations, up to a time scalee^O(L/n), which tends to infinity in the thermodynamic limit.