Steady-state Phase Diagram of a Weakly Driven Chiral-coupled Atomic Chain

A chiral-coupled atomic chain of two-level quantum emitters allows strong resonant dipole-dipole interactions, which enables significant collective couplings between every other emitters. We numerically obtain the steady-state phase diagram of such system under weak excitations, where interaction-driven states of crystalline orders, edge or hole excitations, and dichotomy of chiral flow are identified. We distinguish these phases by participation ratios and structure factors, and find two critical points which relate to decoherence-free subradiant sectors of the system. We further investigate the transport of excitations and emergence of crystalline orders under spatially-varying excitation detunings, and present non-ergodic butterfly-like system dynamics in the phase of extended hole excitations with a signature of persistent subharmonic oscillations. Our results demonstrate the interaction-induced quantum phases of matter with chiral couplings, and pave the way toward simulations of many-body states in nonreciprocal quantum optical systems.

A chiral-coupled atomic chain of two-level quantum emitters allows strong resonant dipole-dipole interactions, which enables significant collective couplings between every other emitters. We numerically obtain the steady-state phase diagram of such system under weak excitations, where interaction-driven states of crystalline orders, edge or hole excitations, and dichotomy of chiral flow are identified. We distinguish these phases by participation ratios and structure factors, and find two critical points which relate to decoherence-free subradiant sectors of the system. We further investigate the transport of excitations and emergence of crystalline orders under spatially-varying excitation detunings, and present non-ergodic butterfly-like system dynamics in the phase of extended hole excitations with a signature of persistent subharmonic oscillations. Our results demonstrate the interaction-induced quantum phases of matter with chiral couplings, and pave the way toward simulations of many-body states in nonreciprocal quantum optical systems.
Introduction.-A chiral-coupled atomic system [1,2] from an atom-fiber [3] or an atom-waveguide [4,5] interface presents the capability to engineer the directionality of light transmissions. This leads to a broken time-reversal symmetry of light-matter couplings, and results in nonreciprocal decay channels. Unidirectional coupling [6] can therefore be enabled by spin-momentum locking [7,8], where light propagation highly correlates to its transverse spin angular momentum. Such one-dimensional (1D) nanophotonics system has been studied to create mesoscopic quantum correlations [9] and quantum spin dimers [10][11][12], to enable simulations of long-range quantum magnetism [13], and to manifest emerging universal dynamics [14] and strong photonphoton correlations [15]. These compelling predictions rely on the emergence of infinite-range resonant dipole-dipole interactions (RDDI) [16] in 1D system-reservoir interactions, in huge contrast to the ones in free-space, which decrease with an inter-atomic distance in the long range [17].
The nonreciprocal decay channels of chiral-coupled systems can be tuned by external magnetic fields [3][4][5], such that the amount of light transmissions in the allowed direction can be controlled by the internal states of quantum emitters [3]. This can be attributed to reservoir engineering, which has spurred many interesting studies of non-equilibrium phase transitions in driven-dissipative quantum systems [18][19][20]. In such open systems, self-organized supersolid phase [21] of Bose-Einstein condensate can be realized by coupling to an optical cavity, and exotic spin phases and multipartite entangled states respectively can be dissipatively prepared in laserexcited Rydberg atoms [22] and ions [23]. Under engineered dissipations, Majorana edge modes as topological states of matter [24,25] and critical phenomena at steady-state phase transitions [26] are also predicted in lattice fermions. These quantum phases of matter under non-equilibrium phase transitions show the potential to explore dynamical phases driven by competing dissipation and interaction strengths [27].
In this Letter, we consider two-level quantum emitters coupled to a nanofiber or waveguide with equal inter-atomic distances. We obtain the steady-state phase diagram of such chiral-coupled atomic chain in the low saturation limit, de-termined by two competing parameters of directionality and dipole-dipole interaction strengths. The interaction-driven phases include the states with extended distributions (ETD), crystalline orders (CO), bi-edge/hole excitations (BEE/BHE), and of chiral-flow dichotomy (CFD), which we classify by participation ratios and structure factors. Two critical points are also located, where time to reach steady states is longer than a power-law dependence of system sizes. This critically slow equilibrium of system dynamics relates to the decoherence-free sectors of the eigen-spectrum. We further explore the possibility to relocate the atomic excitations via spatially-varying field detunings or excitation directions. Finally, in the phases of BEE and BHE, non-ergodic signatures of subharmonic oscillations emerge, and we specifically present a butterfly-like system dynamics as an example. The steady-state phases and their dynamical evolutions investigated here present a distinct interplay between RDDI and directionality of a chiral-coupled atomic chain, which give insights to preparations and simulations of many-body states in nonreciprocal quantum optical systems.
Model.-We consider a generic driven-dissipative model in Lindblad forms for a 1D chiral-coupled atomic chain [12], The external and coherent field excitation is denoted as which drives N two-level quantum emitters (|g and |e being the ground and excited state respectively) with spatially dependent Rabi frequencies Ω µ and detunings δ µ , and dipole operators are σ † µ ≡ |e µ g| with σ µ = (σ † µ ) † . The coherent parts of RDDI are 1. Steady-state phase diagram of a weakly driven chiral-coupled atomic chain. The chiral-coupled system consists of two-level quantum emitters with nonreciprocal decay rates of γL = γR, mediated by a nanofiber or waveguide. (a) In the steady state, the phases of extended distributions (ETD), finite crystalline orders (CO), bi-edge/hole excitations (BEE/BHE), and chiral-flow dichotomy (CFD, wavy dashes) are identified under parameter spaces of directionality factor D and dimensionless inter-atomic distance ξ. For different scalings of participation ratios on the number of atoms (N α ), we further separate the BHE phase into strong (s), moderate (m), and weak (w) regimes with α > 0.5 (below D ≈ 0.37), 0.1 < α < 0.5, and α < 0.1 (above D ≈ 0.67), respectively. Shaded regions at two lower corners of the diagram represent the critical regimes. (b) Participation ratios (PR) at different horizontal cuts of the diagram (increasing D from the lower to upper panels), showing abrupt changes crossing the CO phases, with N = 50 (solid line), 100 (dashes), and 150 (dots). and the dissipative ones in Lindblad forms are They determine the collective energy shifts and decay rates, respectively, which mediate the whole system in infinite-range of sinusoidal forms. The subscripts L(R) label the left(right)propagating decay channels, and k s = 2π/λ denotes the wave vector for the transition wavelength λ. We proceed to solve for the steady-state solutions of Eq. (1) in the low saturation limit, which truncates the hierarchycoupled equations [28,29] to self-consistently coupled dipole operators. We obtain steady-state σ (s) satisfyingσ µ = 0, where a uniform θ s with small Rabi frequency Ω represents the excitation angle to the alignment of the chain, and r denotes the atomic distributions. In the coupling matrix M , asymmetry of off-diagonal matrix elements arises due to unequal γ L(R) , and they are where r µ,ν = r µ − r ν . From Eq. (6), we determine the interaction-driven quantum phases of matter which mainly rely on the interplay between the chirality of chiral-coupled systems and RDDI determined by inter-atomic distances. The external spatially-varying excitation detunings and excitation angles play extra roles in relocating the atomic excitations and imprinting extra phases on the atoms, which we investigate below.
Steady-state phase diagram.-We use the directionality factor D ≡ (γ R − γ L )/γ [3] to quantify the degree of light transmissions of the system with a normalized decay channel γ R + γ L = γ. D = ±1 and 0 present the unidirectional [10,33,34] and reciprocal couplings respectively. The RDDI strength can be quantified by ξ ≡ k s r µ+1,µ , where ξ = 0 or 2π represents the strong coupling regime with vanishing dispersions that cannot be achieved in free space [17]. In Fig.  1, we numerically obtain the steady-state phase diagram of a weakly driven chiral-coupled atomic chain at θ s = π/2. The states with CO are determined by finite structure factors S(k) = N j,m e ik(j−m)P jPm [28] with normalized ex- . This phase boundary extends symmetrically from ξ = π/2 toward two critical points at the lower corners of the diagram, which also coincides with abrupt changes of participation ra- with the Heaviside step function Θ evaluates the state variations from an ETD phase (narrow region close to ξ = 0) of uniform distributions N −1 . As shown in Fig. 1(b), sharp rises of P R emerge when ξ crosses extended CO or BHE (P 1(N ) < N −1 ) phases. This is in huge contrast to the localized BEE phase (P 1(N ) > N −1 ), where P R stays constantly small as system size increases. At D = 0, the phase boundary between BEE and BHE starts from ξ = π/3 [28] and collapses to π/2 as D increases. For an increasing D, P R of BHE phase suppresses, and we further distinguish it by different scalings of system sizes [28]. The strong (s) BHE region suggests a significant extended state distributions with hole excitations at both ends of the chain, which endows the emergence of persistent subharmonic dynamics we will present below.
We further locate two critical points at D = 0, ξ = 0 and π respectively, and a phase of CFD which shows very different P R for an even or odd number of atoms. The critically slow behavior can be identified by the time of atomic excitations to reach the other end, which has an algebraic dependence of directionality D as shown in Fig. 2(a). This slow equilibrium of system dynamics also manifests in other extended phases. Criticality arises due to (N − 1) decoherencefree modes in the eigenvalues of (−N/2, 0, ..., 0)γ, obtained from M , where zero decay modes are not allowed at other ξ's. These decoherence-free modes mark a distinct region separating from other noncritical ones at D = 0 with subradiant sectors of small but finite decay rates.
We demonstrate two examples of CO and BEE phases in Fig. 2(b), where a period of eight sites and significant excitations at both edges emerge, respectively. CFD phase in Fig.  2(c) can be characterized by a dichotomy of steady states for even and odd number of atoms. For the first impression with a finite D, we expect of a flow of atomic excitations toward the preferred direction with its minimum in the opposite side. By contrast another configuration with its minimum moving toward the center of the chain as D decreases shows up in odd number of atoms. It is the extra atom in an odd chain that distinguishes the two configurations near D = 0 as shown in the inset of Fig. 2(c), which should satisfy the inversion symmetry as D → 0. Even(odd) number of atoms in this particular phase presents balanced(unbalanced) RDDI of alternating ±γ R/L , and this dichotomy can be further classified by two distinguishing P R's. The more delocalized state for an even chain reaches its maximal P R close to the critical point, whereas the P R for an odd case never exceeds the value at D = 1, and instead shows a minimum at a finite and universal value of D = 0.28. The population accumulation toward both edges of the chain, representing a more localized phase, resembles an unpaired spin localization in a tight-binding Su-Schrieffer-Heeger model [36] or an edge state of bosons in a superlattice Bose-Hubbard model [37]. We note that the cross of P R near D = 1 is due to finite size effect, which collapses in thermodynamic limit.
Transport of atomic excitations and emergence of crystalline order.-Next we investigate the effect of spatiallydependent δ µ and θ s , which lead to redistribution of atomic excitations. We quantify the transport of atomic excitations by a difference of them between the left and right parts of the chain [38], where we have excluded the central one for an odd N . Positive or negative T p represents that the left or right parts of the chain are more occupied, and right(left) linearly-increasing δ µ should favor positive(negative) T p since atoms are less excited under off-resonant driving fields, in a sense of noninteracting regime with negligible RDDI. By contrast in Fig. 3(a), for an onset of small slope of linearly-increasing detunings to the right, the interaction-driven transport changes a positive T p to negative one when ξ π/8, comparing T p = 0 for both BEE and sBHE phases at D = 0 with a vanishing slope (symmetric distributions with inversion symmetry). This is more evident for moderate slopes, where negative T p shows up in different ranges, presenting the competition between RDDI strength and excitation detunings, that is, the case with a larger ξ covers less ranges of slopes for negative T p , indicating less effect of RDDI on transport properties. Eventually all T p turn to a positive side under large slopes as in noninteracting regimes. For transport of the excitations in BEE or sBHE phases, in Fig. 3(b) we show the atom-atom correlations which highlight the central localization of excitations as ξ increases under a harmonic-potential-like detuning. This presents the repulsion of atomic excitations in BEE phase and the dominance of external potential over the hole excitations in sBHE phase. An example of edge accumulations induced from CO phase and finite CO emerged from CFD phase in Fig. 3(c) further shows the effect of δ µ , which raises broadened CO. On the other . Arrows indicate the cases of ξ = π/32, π/16, π/8, π/4, π/2, 3π/4, and 7π/8, at D = 0. (b) A localization of atomic excitations appears in atom-atom correlationsPiPj under a harmonic-potential-like δµ = h[µ − (N + 1)/2] 2 with h = N −1 as ξ increases from π/4 (left), π/2 (middle), to 3.9π/4 (right), at D = 0. (c) Edge accumulation is evident inPiPj at D = 1 with ξ = 3π/4 in the left, and in the right broad CO emerge at D = 0.01 and ξ = π, both under a smaller h = 0.01N −1 of harmonic-potential-like detunings with the same color bar in (b). (d) Finite structure factors S(k) move from low to high k from θs = 7π/8, 6π/8, ..., to 2π/8, at D = 1, ξ = π/4, and N = 101. In the right panel, for D = 0.3, ξ = π/4, the localized edge excitation moves from the right end at θs = 7π/8 (solid), to the left end at θs = π/8 (dash-dots), comparing BEE phase at θs = π/2 (dots). N = 51 for (a-c) and right panel of (d).
hand, in Fig. 3(d) the effect of excitation angles θ s manifests in moving the locations of finite structure factors by imprinting spatially-dependent phases on the chain, and in relocating the edge excitations to the left-or right-most of the chain as if a strong linear potential is induced for single-edge excitations. A joint manipulation of excitation detunings and angles thus enables a controllable transport of localized atomic excitations.
Persistent subharmonic oscillations.-Here we present the long time dynamics of equilibration in steady states, specifically in the sBHE phase which shows strong extended features of hole excitations. There is no definite region for this collectively coupled dynamics, but can be approximately confined in mBHE and sBHE for D 0.5 and BEE for smaller D.
In Fig. 4, we show two examples in non-cascaded (D = 0) and reciprocal (D = 0) coupling regimes. The non-cascaded time dynamics shows oscillating spread of populations from left to right initially since γ R > γ L , and then back and forth until equilibrium. We present an oscillation of total population P t with two time scales in the lower panel of Fig. 4(a), a signature of persistent subharmonic evolutions. This is further enhanced in the reciprocal coupling regime where much longer time is required to reach steady states. We attribute this phenomena as non-ergodic, associated to discrete timecrystalline order in subharmonic temporal responses [39] and nonequilibrium many-body scars [40] with a constrained local Hilbert space. The subradiant sectors of the eigen-spectrum under reciprocal couplings can provide some clues on this long time behaviors, where large energy splittings due to collective frequency shifts emerge toward the lowest subradiant decay rates as we go from BEE to sBHE phases along the line of D = 0 [28]. This leads to highly dispersive couplings arising from strong RDDI and thus system dynamics that is far from equilibrium. The reoccurring patterns in Fig. 4(b) re-semble a butterfly, where agglomeration of atomic excitations resides around the center of the chain before expanding to the edges, presenting another feature of non-ergodicity. As a final remark, we note of the scale of P t by Ω/γ instead of (Ω/γ) 2 , showing the excitation enhancement due to the collective couplings of RDDI, compared to the noninteracting case which should be order of (Ω/γ) 2 .
In conclusion, the chiral-coupled atomic chain presents fruitful interaction-driven quantum phases of matter under driven-dissipative settings, with competing interactions between intrinsic 1D RDDI, directionality, and external excitation parameters. Under a weak excitation condition, the system shows critically slow behaviors, crystalline orders, and localized edge or extended hole excitations. The nonequilibrating states of extended hole excitations specifically manifest long time oscillations, in huge contrast to the states with CO and phases of BEE and wBHE close to the unidirectional coupling regime. Future directions can lead to unraveling clear mechanism for initiation of persistent subharmonic time evolutions, its relation to ergodicity of the chiral-coupled system, or many-body simulations of exotic or topological states in nonreciprocal quantum optical systems.
We acknowledge the support from the Ministry of Science and Technology (MOST), Taiwan we are able to derive the coupled equations with a hierarchy of multiple atomic correlations. We will see below that the hierarchy arises due to the resonant dipole-dipole interactions of H L(R) and L L(R) [ρ], where their explicit forms can be found in the main paper.
As a demonstration of hierarchy-coupled equations, we consider a uniform excitation field with Rabi frequency Ω at a right angle θ s = π/2. First, the time evolutions of coherence operators σ µ readṡ where γ L(R) respectively quantifies the couplings between the atom to the rest of left(right) of the chain. Next, the evolution of excitation population readṡ where H.c is Hermitian conjugate. We can see that the above two single atomic operators couple with each other via two-body correlations of σ gg µ σ ν , σ ee µ σ ν , and σ † µ σ ν . This indicates the hierarchy relations between nth and (n + 1)th moments of operators up to N th ones.
To have another taste of hierarchy-coupled equation of two-body correlations, we show the time evolutions of σ † µ σ ν as an example, which reads where three-body operators emerge. In treating the effect of quantum fluctuations in dissipative quantum optical systems [31], we can usually truncate the hierarchy up to two-body noise correlations due to the stochastic nature of quantum noises. For a general driven-dissipative system we consider here, many-body states can be all explored (a total of 2 N configurations for N two-level quantum registers), and this is exactly the merit of quantum resource for quantum information processing and quantum computation. However, this also prevents appropriate numerical simulations of general quantum dynamics by classical computers, except for systems with smaller higher order moments or with a technique of matrix product states which are efficient in finding the many-body ground states. Next we simplify the hierarchy-coupled equations by taking a low saturation regime, where σ gg µ ≈ 1 ≫ σ ee µ . Under this condition, Eq. (9) reduces tȯ and the evolutions of all other higher order moments can be decomposed in terms of the above one. For example, σ † µ σ ν = σ † µ σ ν , where · denotes the expectation values, such that d σ † µ σ ν /dt = σ † µ σ ν + σ † µ σ ν . In the main paper, we solve the steady-state solutions of Eq. (12), which manifest many interesting quantum phases in the steady state.

FINITE STRUCTURE FACTORS IN THERMODYNAMIC LIMIT
To identify crystalline orders (CO), we use structure factors defined as to characterize the phase with CO.P j denotes the normalized steady-state population,r The phase boundaries of CO-BEE (bi-edge excitations) and CO-BHE (bi-hole excitations) in the main paper are determined respectively by abrupt occurrences of edge or hole excitations at both ends of the chain. This coincides with clear differences of participation ratios between extended and localized states, and with a distinction by a finite CO. The ETD (extended distributions)-CO phase boundary, however, is not obvious since ETD and CO phases are both extended. Close to this phase boundary, the structure factors are decreasing as N increases. Therefore, we fit the maximum (Max.) of S(k = 0) in terms of N by (aN −1 + b) with fitting parameters a and b, where a finite b > 0 represents a finite S(k) in thermodynamic limit of N → ∞, and we characterize it as CO phase.
We  Here we show the analytical results for the phase boundary at ξ = π/3 from BEE to BHE at D = 0, and the numerically obtained boundaries between sBHE, mBHE, and wBHE. For δ µ = 0 and θ s = π/2 considered in the phase diagram in the main paper, we have the coupling matrix at D = 0, which is symmetric, and so is its inverse M −1 . At the BEE-BHE phase boundary,P 1(N ) should be equal to N −1 to distinguish from edge (P 1(N ) > N −1 ) and hole (P 1(N ) < N −1 ) excitations. We can obtainP j from M −1 . For N = 3, and takeP 1 (P N =P 1 ) as an example, which gives ξ = π/3 for A = 1 whenP 1 = 1/3. Another solution of ξ ≈ 0 from A = 0 is ignored due to the divergence ofP 1 . This relates to critical regimes, which we will discuss later. For a general N , we obtain the solution of ξ from which gives again ξ = π/3 for A = 1. When ξ ≈ 0, we obtainP 1(N ) = 1/2, which indicates that the atomic excitations are populated equally to the edges of the chain. When ξ = π, we haveP 1(N ) = 1/(4N ) for N ≫ 1, a signature of hole or null excitations at the edges.
Phase boundaries between sBHE, mBHE, and wBHE In the BHE phase, we note of the suppression of the participation ratios (PR) as D increases, which indicates that BHE becomes less delocalized. To quantify various regions of PR, we again fit the maximum of PR in terms of N by βN α with fitting parameters β and α, where α denotes the degree of scalings for strong (s), moderate (m), and weak (w) regimes with α > 0.5 (below D ≈ 0.37), 0.1 < α < 0.5, and α < 0.1 (above D ≈ 0.67), respectively. This algebraic dependence of N is suggested at D = 0, where PR∝ N −1 .
Below we show the PR and their fitting parameter of α around sBHE, mBHE, and wBHE regimes by three data points, Here we take N = 100 as an example, and show below the energy splittings of the eigen-spectrum from BEE to sBHE in the phase diagram at D = 0. The eigen-spectrum is directly obtained from the coupling matrix M defined in the main paper. FIG. 5. Frequency shifts of the eigen-spectrum for N = 100. Here we sort the eigenvalues Im(λn) according to ascendant decay rates of Re(λn). In the upper panel, the inset shows the case of ξ/π = 0.1 in logarithmic scales. We specifically focus on the subradiant sectors of n 80, which is way below Re(−λn) = γ. In Fig. 5, we sort the frequency shifts of λ n according to the sorting from low to high decay rates as n increases. The trend of the logarithmic scales of eigen decay rates is similar for all cases of ξ, except some larger superradiant or lower subradiant eigenvalues appear for different parameters of ξ. We scan ξ from BEE (ξ < π/3) to sBHE (ξ > π/3), and find the emergence of energy splittings appearing toward the lowest subradiant sectors. Since the long time behaviors of subharmonic oscillations observed in the main paper are more significant in the phase of sBHE, we attribute the non-ergodic oscillations to the highly dispersive subradiant sectors of the spectrum. This can be further clarified in Fig. 6, where we plot the frequency shifts in an ascendant direction. Again as we go from BEE to sBHE phases, the eigen frequency shifts close to Im(λ n ) = 0 start to split, and an energy gap-like jump appears from red-detuned to blue-detuned shifts. This opening of energy gap corresponds to the subradiant sectors which allow relatively lower decay rates with finite energy shifts. In the perspective of steady states under driven-dissipative settings, these constrained subradiant states are responsible for the persistent subharmonic oscillations or butterfly-like time dynamics we present in the main paper.