Quasiparticles of Decoherence Processes in Open Quantum Many-Body Systems: Incoherentons

The relaxation dynamics of an open quantum system is determined by the competition between the coherent Hamiltonian dynamics of a system and the dissipative dynamics due to interactions with environments. It is therefore of fundamental interest to understand the transition from the coherent to incoherent regimes. We find that hitherto unrecognized quasiparticles -- incoherentons -- describe this coherent-to-incoherent transition in eigenmodes of a Liouvillian superoperator that governs the dynamics of an open quantum many-body system. Here, an incoherenton is defined as an interchain bound state in an auxiliary ladder system that represents the density matrix of a system. The Liouvillian eigenmodes are classified into groups with different decay rates that reflect the number of incoherentons involved therein. We also introduce a spectral gap -- quantum coherence gap -- that separates the different groups of eigenmodes. We demonstrate the existence of incoherentons in a lattice boson model subject to dephasing, and show that the quantum coherence gap closes when incoherentons are deconfined, which signals a dynamical transition from incoherent relaxation with exponential decay to coherent oscillatory relaxation. Furthermore, we discuss how the decoherence dynamics of quantum many-body systems can be understood in terms of the generation, localization, and diffusion of incoherentons.


I. INTRODUCTION
Understanding the role of environments on quantum coherence presents a key challenge in quantum physics [1][2][3].The concomitant decoherence of quantum superposition of a system places a major obstacle in the development of quantum technologies [4][5][6][7].Moreover, there has been a surge of interest in nonequilibrium dynamics of open quantum many-body systems owing to experimental progress in atomic, molecular, and optical (AMO) systems, enabling one to control not only the Hamiltonian of a quantum system but also its coupling to an environment [8][9][10][11][12][13][14][15][16][17][18].
The dynamics of an open quantum system can, in general, be described by a quantum master equation for its density matrix.In particular, in a typical AMO system, the weak coupling and the separation of time scales between the system and an environment allow the dynamics of the density matrix to be described by a Markovian quantum master equation [3].The superoperator that generates the time evolution of the density matrix is referred to as the Liouvillian L. The relaxation dynamics of an open quantum system is fully characterized by the complex spectrum and eigenmodes of L. In general, the Liouvillian L consists of a coherent part describing the unitary time evolution governed by the Hamiltonian of the system and an incoherent part due to the coupling with the environment.The competition between these contributions causes a transition from a coherent regime to an incoherent one.Such a coherent-to-incoherent transition, a phenomenon found in many quantum systems [19][20][21][22][23][24][25][26][27][28], is detrimental to quantum technologies, including quantum computation.However, it is a formidable task to understand how decoherence proceeds in * taiki.haga@omu.ac.jp open quantum many-body systems because of exponentially large Hilbert-space dimensions.It is highly desirable to establish an effective description of the competition between coherent and incoherent dynamics in many-body systems.In this regard, it should be recalled that the concepts of spectral gaps and quasiparticles play a pivotal role in quantum many-body physics.In isolated systems, quantum phase transitions in the ground state are characterized by the closing of the spectral gap [29], and the low-energy behavior is governed by quasiparticle excitations, which allow an effective description of complex many-body systems [30].
In the present paper, we investigate spectral gaps and quasiparticles that characterize the physics of Markovian open quantum many-body systems described by a Liouvillian superoperator.Note, however, that the incoherent-coherent transition often occurs far from the steady state, in which the conventional low-energy description in terms of quasiparticles is inapplicable.We here discover quasiparticles, "incoherentons," that naturally describe the incoherent-coherent transition in Liouvillian eigenmodes of open quantum many-body systems.As opposed to the conventional notion of quasiparticles, incoherentons are applicable to far-from-equilibrium regions.Incoherentons are defined on a space of operators due to the matrix nature of the density matrix.To show this, we use the fact that any density matrix of a system can be mapped to a vector in the tensor product space of bra and ket spaces.Since this product space can be interpreted as the Hilbert space of a ladder system consisting of two chains of bra and ket spaces [see Fig. 1(a)], we call such a mapping the ladder representation of the density matrix.In the ladder representation, the Liouvillian L is mapped to a non-Hermitian Hamiltonian L of the ladder system.The coherent part of L governs the independent dynamics of particles in each chain, while the incoherent part acts as a non-Hermitian interaction between different chains.Thus, depending on which contribution is dominant, the Liouvillian eigenmode either forms a scattering state extended over the entire ladder or an interchain bound state in which the degrees of freedom of the two chains are strongly correlated.Since the existence of such a bound state implies the localization of matrix elements near diagonal components in the original matrix representation, we will refer to it as an incoherenton [see Fig. 1(a)].
The concept of incoherentons provides several insights into the dynamics of open quantum many-body systems and allows us to discover a universal mechanism for incoherentcoherent transitions.They are summarized below as well as in Figs.1(b)-(d), and will be discussed in the following sections.

Deconfinement of incoherentons (Sec. III):
The incoherent-coherent transition of Liouvillian eigenmodes can be understood by the deconfinement of incoherentons.Since the dissipation corresponds to chain-to-chain interactions in the ladder representation, the confinement length of an incoherenton increases with decreasing dissipation, and eventually a transition from bound to scattering states occurs at some critical strength of dissipation [see the left panel of Fig. 1(b)].

Quantum coherence gap closing (Sec. III):
When the dissipation is sufficiently strong, a gap exists between groups of Liouvillian eigenvalues with different numbers of incoherentons.We call the gap between such groups quantum coherence (QC) gap because it separates groups of eigenmodes with different degrees of quantum coherence.The QC gap ∆ QC closes at the deconfinement transition (see the middle panel of Fig. 1(b), where λ denotes the eigenvalues of L).
Incoherent-coherent dynamical transition (Sec.IV): The QC gap closing signals the onset of a dynamical transition from overdamped relaxation, where an inhomogeneous initial state relaxes exponentially to a uniform steady state, to underdamped relaxation, where the local density exhibits oscillatory behavior (see the right panel of Fig. 1(b), where O is an appropriate observable).We argue that this provides a hitherto unknown type of incoherent-coherent transitions in an extended lattice system.

Hierarchy of eigenmodes (Secs. V, VI, and VII):
The many-body eigenmodes are classified into groups with different decay rates characterized by the number of incoherentons involved therein [see Fig. 1(c), where the small gray boxes represent incoherentons].Each group of eigenvalues is separated from the others by the QC gaps.The more incoherentons an eigenmode involves, the slower it decays.

Many-body decoherence (Sec. VIII):
The number of incoherentons in the density matrix increases as relaxation proceeds, which means that the relaxation of a many-body state is accompanied by the production of incoherentons.Furthermore, the late decoherence process is characterized by the localization and diffusion of incoherentons [see Fig. 1(d)].Here, we highlight the distinction between the incoherentcoherent transition described in this work and conventional phase transitions in isolated and open quantum many-body systems.Table I summarizes different types of transitions, spectral gaps, and characteristic length scales.In isolated quantum systems (see the left column of Table I), the energy gap ∆ E of a Hamiltonian is defined as the energy difference between the ground state and the first excited state.The correlation length ξ g of the ground state and ∆ E are related to each other by ξ g ∼ v/∆ E , where v is the propagation velocity of low-energy excitations with wavelengths comparable to ξ g .Here and henceforth, the Planck constant ℏ is set to unity.At a quantum phase transition of the ground state, ξ g diverges, accompanied by the closing of ∆ E and the divergence of characteristic time scales of low-energy excitations [29].
In open quantum systems, a phase transition of the steady state, known as the dissipative phase transition [31][32][33][34][35][36][37][38][39], is characterized by the Liouvillian gap ∆ L , which is defined as the smallest absolute value of the real parts of nonzero Liouvillian eigenvalues (see the middle column of Table I).The relation between the correlation length ξ s of the steady state and ∆ L is given by ξ s ∼ v/∆ L , where v is the propagation velocity of excitations near the steady state.The dissipative phase transition is characterized by the divergence of ξ s and the closing of ∆ L [40][41][42][43][44][45][46].The longest timescale for the system to reach the steady state is expected to be inversely proportional to ∆ L [47][48][49] (see, however, Refs.[50,51] for exceptions).Thus, the closing of ∆ L leads to the divergence of the relaxation time.
The deconfinement of incoherentons together with the QC gap closing constitutes the third type of transition in quantum many-body systems (see the right column of Table I).The relation between the confinement length ξ con of incoherentons and the QC gap ∆ QC is given by where Γ is the decay rate of relevant eigenmodes and ξ con is measured in units of the lattice spacing.An important distinction of the deconfinement of incoherentons from other wellknown transitions is that it is a transition of non-steady eigenmodes having finite lifetimes.Thus, the deconfinement of in- This paper is organized as follows.Section II details the ladder representation of the Liouvillian and introduces a system of hard-core bosons subjected to on-site dephasing, serving as a representative model for open quantum many-body systems.In Sec.III, the concept of incoherenton is introduced for the one-particle case.We describe the deconfinement of incoherentons and the QC gap closing in terms of the prototypical model.Section IV demonstrates that the relaxation dynamics of particle density display an incoherent-coherent transition corresponding to the parameter at which the QC gap closes.In Sec.V, the concept of incoherentons is generalized to many-body systems.By numerically diagonalizing the Liouvillian of the dephasing hard-core bosons, we demonstrate the deconfinement of incoherentons and the closing of the QC gap for the many-body case.In Sec.VI, we obtain an exact many-body solution of the dephasing hard-core boson model with the Bethe ansatz method, which analytically confirms the existence of incoherentons and their deconfinement transitions.In Sec.VII, we discuss how the incoherenton framework can be applied in the presence of particle exchange with the environment, and demonstrate that the phenomenology of incoherentons remains intact for small loss and gain rates of particles.Section VIII introduces a simple description of many-body decoherence via incoherentons, identifying three distinct decoherence regimes related to the production, localization, and diffusion of incoherentons.In Sec.IX, we summarize our results and discuss prospects for future work.In Appendix A, general properties of the Liouvillian spectrum and eigenmodes are summarized.In Appendix B, we present a thorough analysis of the Liouvillian spectrum and eigenmodes for the one-particle case without resorting to the Bethe ansatz.In Appendix C, we show that incoherentons do not exist in continuous systems.This fact implies that the spatial discreteness of lattice systems is crucial for the existence of incoherentons.In Appendix D, we discuss measuring incoherenton correlation functions in ultracold atomic systems.In Appendix E, we explore the Liouvillian spectra of a dephasing Bose-Hubbard model through numerical diagonalization and shows the deconfinement of incoherentons within this model.
In Appendix F, we present the results for dephasing hard-core bosons with next-nearest-neighbor hopping.In Appendices E and F, we provide evidence supporting the universality of the incoherenton framework.

A. Liouvillian superoperator
We focus on Markovian open quantum lattice systems with bulk dissipation, in which the dissipation acts uniformly on every site.Within the Born-Markov approximation [3], the time evolution of the density matrix ρ is described by a quantum master equation, which is generated by a Liouvillian superoperator L [52,53]: where [A, B] := AB − BA, {A, B} := AB + BA, and L ν is a Lindblad operator.The quantum master equation ( 2) is justified when the time scale of dynamics induced by the systemenvironment coupling is much longer than the characteristic time scale of the environment.This condition is well satisfied for typical AMO systems such as trapped two-level atoms with spontaneous emission and an optical cavity with photon loss [16,18,38].The index ν for the Lindblad operator L ν denotes the lattice sites and the types of dissipation.We assume that each L ν has support on a finite number of sites.The master equation ( 2) can be rewritten as where the non-Hermitian effective Hamiltonian H eff reads It is convenient to define and In the quantum trajectory description [18], where the dynamics of an open quantum system is described by stochastic trajectories of pure states, L H describes a deterministic time evolution generated by the effective Hamiltonian H eff , and L jump describes quantum jump processes.If the Liouvillian is diagonalizable, its eigenmodes ρ α can be defined by where λ α is the αth eigenvalue and D is the dimension of the Hilbert space H of the system.A steady state ρ ss corresponds i, j ρ i j |i⟩ ⟨ j| is mapped to a vector |ρ) = i, j ρ i j |i⟩ ⊗ | j⟩.In a onedimensional tight-binding model with L sites, a basis vector |i⟩ of the Hilbert space can be written as |n 1 , ..., n L ⟩, where n l = 0, 1, ... is the occupation number of particles at site l = 1, ..., L. Similarly, ⟨ j| is represented as ⟨n ′ 1 , ..., n ′ L |.In the ladder representation, a basis vector |i⟩ ⊗ | j⟩ describes a state of a two-leg ladder.
to an eigenmode with zero eigenvalue.We arrange the eigenvalues General properties of the Liouvillian spectrum and eigenmodes are summarized in Appendix A. In terms of the Liouvillian eigenmodes, the time evolution of the density matrix is given by where c α is the coefficient of eigenmode expansion of the initial density matrix.We have assumed that the steady state ρ 0 = ρ ss is unique.Equation ( 8) implies that the relaxation dynamics of an open quantum system is fully characterized by the spectrum and eigenmodes of the Liouvillian.Let {|i⟩} i=1,...,D be an orthonormal basis set of H that specifies real-space configurations of particles or spins.For example, one can consider the real-space Fock basis |n 1 , ..., n L ⟩, where n l = 0, 1, ... denotes the occupation number of particles at site l, and L is the system size.In terms of this orthonormal basis, the density matrix ρ is written as where ρ i j := ⟨i|ρ| j⟩.Let us identify an operator |i⟩ ⟨ j| on H with a vector |i⟩ ⊗ | j⟩ in the tensor product space H ⊗ H, the first (second) space of which will be referred to as the ket (bra) space [49,[54][55][56][57].Then, the density matrix ( 9) is mapped onto the following vector: where we have used a round ket symbol |...) to emphasize that it belongs to H ⊗ H rather than H.It should be noted that, for one-dimensional cases, H ⊗ H can be considered as the Hilbert space of a ladder system composed of two chains (see Fig. 2).Thus, in the following, we refer to Eq. ( 10) as the ladder representation of the density matrix.
In the ladder representation, the Liouvillian L is mapped to a non-Hermitian Hamiltonian L of the ladder system.The ladder representations of L H and L jump are given by and where where ρ α in Eq. ( 7) and |ρ α ) in Eq. ( 13) are related to each other by Eqs. ( 9) and (10).We comment on the diagonalizability of the Liouvillian.Contrary to Hermitian operators, a non-Hermitian operator is not diagonalizable at exceptional points (EPs) [58][59][60][61][62][63][64][65][66][67].While the set of EPs has zero measure in the parameter space (see, e.g., Sec.2.6.1 in Ref. [66]), the system can encounter an exceptional point when a certain parameter is continuously adjusted while keeping others fixed.We note, however, that the diagonalizability of a Liouvillian is unimportant for our argument in this work.An EP only indicates that the Liouvillian contains a Jordan block with size larger than one.In the most typical case of the lowest-order EP, two eigenvectors coalesce, and thus the Liouvillian involves a two-by-two Jordan block.Nevertheless, the remaining eigenvectors, corresponding to one-by-one Jordan blocks, are unaffected, and it is worth studying their structure.Consequently, even if the Liouvillian is not diagonalizable, our argument based on Liouvillian eigenmodes is applicable due to the predominance of one-by-one Jordan blocks in all eigenmodes.

B. Example: hard-core bosons under dephasing
We introduce a prototypical model of open quantum manybody systems, which will be analyzed in the following sections to demonstrate the concept of incoherentons.The system is defined on a one-dimensional lattice with size L under the periodic boundary condition.The Hamiltonian of the system is given by where b † l and b l are the creation and annihilation operators of a boson at site l, and J represents the tunneling amplitude.We assume the hard-core condition (b † l ) 2 = 0, which prohibits more than two particles from occupying a single site.The Lindblad operators for on-site dephasing are given by where γ denotes the strength of dephasing.Note that the total particle number N = L l=1 b † l b l is conserved, i.e., Tr[NL(ρ)] = 0 for any density matrix ρ.The steady state of the corresponding master equation is the infinite-temperature state ρ ss = D −1 I, where I is the identity operator, which is a consequence of the Hermiticity of the Lindblad operator L l .Figure 3(a) shows a schematic illustration of the model.
Let n l = 0, 1 be the occupation number of the hard-core bosons at site l.We define 2 L orthonormal basis vectors where |v⟩ is the vacuum state of the system.In the ladder representation, an operator |{n l }⟩ ⟨{m l }| on the Hilbert space of the system is interpreted as a state |{n l }⟩ ⊗ |{m l }⟩ of the ladder.The Liouvillian is then rewritten as where b l,+(−) represents the annihilation operator on the first (second) chain of the ladder.The on-site dephasing can be considered as an interchain interaction acting on a particle pair occupying the same rung.Figure 3(b) shows a schematic illustration of the ladder representation.This model can be realized with ultracold atomic gases in an optical lattice.The on-site dephasing can be induced by the combined effect of coherent laser fields coupled to two internal atomic levels and spontaneous emission [68][69][70][71].Suppose that a ground-state atom is excited by a laser field with frequency ω L and subsequently returns to its ground state through spontaneous emission with rate Γ s .Figure 3(c) shows a level diagram of the atom excited by the laser.The transition rate between the ground and excited states is characterized by the Rabi coupling Ω, which is proportional to the intensity of the laser.The detuning of the laser is given by ∆ = ω L − ω eg , where ω eg is the excitation energy of the atom.When |∆| ≫ Ω, Γ s , the excited state can be adiabatically eliminated and one obtains the Lindblad master equation with an on-site dephasing γ = Γ s Ω 2 /∆ 2 [38,68,69].The dephasingtype Lindblad operator given by Eq. ( 15) also appears in a master equation of ultracold atoms in an optical lattice driven by a stochastically fluctuating on-site potential [72].

III. INCOHERENT-COHERENT TRANSITION AS DECONFINEMENT OF INCOHERENTONS
A. Incoherenton: an interchain bound state The non-Hermitian Hamiltonian LH defined by Eq. ( 11) independently acts on each chain of the ladder and does not create correlations between the bra space and the ket space.If the Hamiltonian only contains kinetic energy terms which cause hopping of particles along each chain, LH prefers planewave eigenmodes extended over each chain of the ladder.On the other hand, Ljump defined by Eq. ( 12) plays a role of an interchain interaction, e.g., see Eq. ( 17).Since each Lindblad operator L ν has its support on a finite number of sites, Ljump describes a local interaction between chains.The interchain Hamiltonian Ljump leads to the formation of an interchain bound state, in which the degrees of freedom in each chain are strongly correlated.As a consequence, in the case of a one-particle system, the eigenmodes of L can be classified into the following two groups depending on which of the contributions from LH and Ljump is dominant: 1. Deconfined eigenmode, where the intrachain kinetic energy dominates the interchain interaction and the eigenmodes are extended over the entire ladder.
2. Confined eigenmode, where the interchain interaction dominates the intrachain kinetic energy and an interchain bound state is formed.
In terms of this classification of eigenmodes, the interplay between the coherent and incoherent dynamics in open quantum systems is understood as a competition between the intrachain kinetic energy and the interchain interaction in Liouvillian eigenmodes.Here, the existence of an interchain bound state in the Lindblad ladder is nontrivial since the interchain interaction Ljump has no clear notion of repulsiveness or attractiveness due to non-Hermiticity of L. Deconfined and confined eigenmodes are schematically illustrated in the right column of Fig. 4. We call the interchain bound state in a confined eigenmode as an "incoherenton."For the one-particle case, an incoherenton is defined as follows.In the ladder representation, an eigenmode can be written as where b † l,+(−) is the creation operator of a particle at site l on the first (second) chain of the ladder and |v) is the vacuum herent eigenmodes for a one-particle system.We call an interchain bound state on the ladder as an incoherenton, for which the confinement length is denoted by ξ con .The gray scale shows the magnitude of each matrix elements, where the darker one shows the larger magnitude.The coherent eigenmode has both diagonal and off-diagonal matrix elements that are comparable in magnitude, whereas the incoherent eigenmode has the predominant diagonal matrix elements.state of the ladder.An incoherenton is represented by matrix elements ρ α,lm that decay exponentially with respect to the relative coordinate, where ξ con is the confinement length of the incoherenton (see the right-bottom panel of Fig. 4).The divergence of ξ con signals deconfinement of an incoherenton.It should be noted that the critical values of control parameters at which the deconfinement transition occurs depend on the eigenmode under consideration.
In the matrix representation of the density matrix, the presence of an incoherenton implies the localization of the eigenmodes near diagonal matrix elements, and its deconfinement implies the delocalization over off-diagonal matrix elements.Since the off-diagonal elements of the density matrix measure the degree of quantum coherence, we refer to the deconfined (confined) eigenmodes in the ladder representation as coherent (incoherent) eigenmodes in the matrix representation.The left column of Fig. 4 illustrates the coherent and incoherent eigenmodes.The confinement length ξ con of an incoherenton in a confined (incoherent) eigenmode quantifies the characteristic length scale in which the quantum coherence in the eigenmode is retained.We also call eigenvalues associated with these eigenmodes as coherent-mode eigenvalues or incoherent-mode eigenvalues.

B. Deconfinement transition and quantum coherence gap
We demonstrate the coexistence of the confined and deconfined eigenmodes for the one-particle case of the model in-troduced in Sec.II B. Let |l⟩ = b † l |v⟩ be the state in which the particle is located at site l.Then, {|l⟩} l=1,...,L provides an orthonormal basis set of the Hilbert space of the one-particle sector.In terms of this basis, an eigenmode of L is written as where we assume the normalization L l,m=1 |ρ α,lm | 2 = 1.In the absence of coherent hopping (J = 0), the matrix elements of L are given by Thus, the action of L is decoupled into a "diagonal" subspace spanned by {|l⟩ ⊗ |l⟩} l=1,...,L and an "off-diagonal" subspace spanned by {|l⟩ ⊗ |m⟩} l,m=1,...,L;l m .In the diagonal subspace, there is an L-fold degenerate eigenvalue λ = 0, and in the offdiagonal subspace, there is an (L 2 − L)-fold degenerate eigenvalue λ = −γ.
We next consider the cases of J 0. A detailed analysis of the one-particle eigenmodes is presented in Appendix B. In the presence of a nonzero coherent hopping, since LH mixes the diagonal subspace with the off-diagonal one, the diagonal eigenmodes with eigenvalues near 0 are no longer exactly diagonal.However, when J ≪ γ, the matrix elements ρ α,lm of the eigenmodes are still localized near the diagonal elements as in Eq. ( 19). ( The smaller S α,off /S α,diag is, the stronger the localization of the eigenmode is.In Fig. 5(a), the incoherent-mode eigenvalues whose eigenmodes satisfy S α,off /S α,diag < 0.1 are shown by red squares, and the other eigenvalues by blue circles.For a weak coherent hopping (J = 0.15 or 0.2), the Liouvillian spectrum consists of the incoherent-mode eigenvalues on the real axis and the coherent-mode eigenvalues accumulated around Re[λ] = −γ = −1.For J = 0, these two types of eigenvalues are highly degenerate at λ = 0 and −γ.The presence of a nonzero J lifts such degeneracy and leads to two elongated bands parallel to the real and imaginary axes.Let us define the quantum coherence (QC) gap ∆ QC as where {λ (c) α } and {λ (i) β } are the coherent-mode eigenvalues (blue circles) and the incoherent-mode eigenvalues (red squares), respectively.The QC gap ∆ QC should not be confused with the Liouvillian gap ∆ L .While ∆ L = |Re[λ 1 ]| is the gap between the steady state and the slowest decaying eigenmode, ∆ QC is the gap between spectral bands of non-steady eigenmodes.
As J increases, ∆ QC decreases, and for J ≥ J c = 0.25, the bands of the coherent-mode and incoherent-mode eigenvalues touch one another.The arrows ( i )-(iv) in Fig. 5(a) track an evolution of one eigenvalue that has the smallest real part in the incoherent-mode spectrum for J < J c . Figure 5(b) shows the color plots of |ρ α,lm | corresponding to these eigenvalues.For ( i ) and (ii), ρ α,lm is well localized near the diagonal elements.In contrast, for (iv), ρ α,lm is delocalized over the offdiagonal elements.Thus, in the ladder representation, the deconfinement transition of an incoherenton occurs at J c = 0.25.From Eq. (B12) in Appendix B, the incoherent-mode eigenvalue with the maximal |Re[λ]| is given in the limit of L → ∞ by Since the real parts of the coherent-mode eigenvalues are identical to −γ in this limit (see Appendix B), we have Thus, the critical value of J at which ∆ QC closes is given by While Eq. ( 24) implies that a real-complex transition occurs at J = J c for the infinite system, the incoherent-mode eigenvalues for a finite system indicated by the arrows ( i )-(iv) in Fig. 5(a) remain real for J > J c .
For each eigenmode |ρ α ), we define the fraction of an onsite bound pair in |ρ α ) as where n l,± = b † l,± b l,± is the number-density operator, and N b,α ∈ [0, 1].For J = 0, N b,α = 1 for the incoherent eigenmodes and N b,α = 0 for the coherent eigenmodes.Figure 5(c) shows {N b,α } α=0,...,L 2 −1 for different values of J and the system sizes L = 10 and L = 20.For J < J c = 0.25, there exists a gap between clusters of N b,α around 0 and 1, and it closes at J = J c .The width of the cluster around N b,α = 0 decreases in inverse proportion to L, which implies that the eigenmodes in this cluster are scattering states that extend over the entire system.In contrast, the width of the cluster around N b,α = 1 is independent of L because the eigenmodes in this cluster are localized with a confinement length ξ con , which is independent of L. Thus, in the limit of L → ∞, N b,α can be considered as an order parameter, which has a nonzero value for incoherent eigenmodes but vanishes for coherent eigenmodes.
The relation between ∆ QC and the confinement length ξ con that is maximized over all incoherentons is given by Eq. ( 1).This relation may be interpreted as follows.Let us denote the typical decay rate of coherent eigenmodes (without incoherenton) as Γ coh and that of incoherent eigenmodes (with an incoherenton) as Γ inc .For on-site dissipation, in which each L ν acts on a single lattice site, the decay rate of extended coherent eigenmodes is larger than that of localized incoherent eigenmodes because the dissipation suppresses offdiagonal elements of the density matrix.Note that Γ inc approaches Γ coh as ξ con diverges to infinity.Thus, we assume (Γ coh − Γ inc )/Γ coh ∼ 1/ξ con , where ξ con is measured in units of the lattice spacing.While we have no general proof for this assumption, we can prove it for the one-particle model under on-site dephasing (see Appendix B).Since ∆ QC ≃ Γ coh − Γ inc , we obtain Eq. ( 1).
It should be noted that the examination of the one-particle spectra and eigenmodes provided above is essentially the same as the analysis in Sec.III of Ref. [49].The author of this reference considered the XX spin chain with bulk dephasing, which is equivalent to dephasing hard-core bosons discussed in our study.We note that Eq. ( 8) and Figs. 3 and 4 in Ref. [49] correspond to Eq. (B12) and Fig. 5 in our study, respectively.The primary distinction between Ref. [49] and our study is that the former emphasizes the most slowly decaying mode which determines the Liouvillian gap ∆ L , whereas the latter focuses on the most rapidly decaying mode in the incoherent eigenmodes which governs the QC gap ∆ QC [see the arrows in Fig. 5(a)].
We draw attention to the connection between our findings and exceptional points (EPs) in non-Hermitian physics [58][59][60][61][62][63][64][65][66][67].The eigenvalues typically exhibit a square-root dependence on the parameter near an EP.As indicated by Eq. ( 24), the critical value J c , where the deconfinement of an incoherenton takes place, signifies an EP of the Liouvillian in the limit of infinite system size.Thus, we have identified a novel class of EPs associated with the transition between coherent and incoherent attributes of eigenmodes.In addition, as detailed in Secs.V and VI, the deconfinement of incoherentons can be generalized to many-body cases.The deconfinement of incoherentons offers a generic mechanism of producing EPs that have a significant consequences on the dynamics in open quantum many-body systems.
The discrete nature of a lattice system is essential for the formation of incoherentons.In fact, for a free particle under a dephasing-type dissipation in continuous space, we can show the absence of such an interchain bound state (see Appendix C).The creation of a bound state due to spatial discreteness in a lattice system has also been known in the conventional two-body problem with a repulsive interaction [73]; while in continuous space no bound state is allowed between particles with a repulsive interaction, on a lattice a bound state exists for an arbitrarily strong repulsive interaction.

IV. INCOHERENT-COHERENT DYNAMICAL TRANSITION
We have shown that the QC gap ∆ QC closes at a certain critical point.Since the dynamics of open quantum systems are intimately related to the Liouvillian spectrum, it is natural to expect that such a change in the structure of the spectrum would significantly alter the transient dynamics to the steady state.In this section, we demonstrate that the QC gap closing is accompanied by an incoherent-coherent dynamical transition from overdamped relaxation dominated by dissipation to underdamped relaxation dominated by unitary time evolution.
In the one-particle sector, {|l⟩ = b † l |v⟩} l=1,...,L provides an orthonormal basis set of the Hilbert space.We denote the corresponding matrix elements of the density matrix as ρ lm , which satisfies a normalization condition l ρ ll = 1.We note that the steady state ρ ss is given by the infinite-temperature state, (ρ ss ) lm = L −1 δ lm .We consider the following initial state, whose particle density is modulated with wavenumber k, where k = 2πs/L (s = −L/2 + 1, ..., L/2), and ∆n represents the amplitude of the density modulation.For J = 0, the initial state given by Eq. (28) shows no time evolution because the action of the Liouvillian to any diagonal density matrix vanishes.In the presence of a nonzero J, the density matrix relaxes toward the uniform steady state ρ ss .The perturbation with wavenumber k to the steady state can selectively excite the incoherent eigenmodes with the same wavenumber k.Thus, we expect that the decay rate of the particle density starting from the initial state (28) with each k reflects the structure of the incoherent-mode spectrum.From the density profile n l (t) = ρ ll (t) at time t, we define which relaxes to zero in the limit of t → ∞ for k 0. Figures 6 (a) and (b) show n(π, t) and n(π/2, t) for different values of J, respectively.The decay rate is a decreasing function of J and it vanishes at J = 0. From these figures, one finds that there exists a k-dependent critical value J rel c (k) below which n(k, t) exhibits an exponential decay e −Γt and above which it shows a damped oscillation e −γt cos ωt.The critical value is estimated as J rel c (π) ≃ 0.25 and J rel c (π/2) ≃ 0.35.An important observation is that J rel c (k) is close to the value of J at which the incoherent eigenmode with wavenumber k exhibits the deconfinement transition.From Eq. (B12), the incoherentmode eigenvalue is written as in terms of the wavenumber k.This eigenvalue becomes complex at a critical value: We then obtain J c (π) = 0.25γ and J c (π/2) ≃ 0.354γ, which are close to J rel c (π) and J rel c (π/2), respectively.It is reasonable that the real-complex transition of λ inc (k) is accompanied by an incoherent-coherent dynamical transition from overdamped to underdamped relaxations starting from an incoherent initial state.Recall that the QC closing occurs at , where the minimum is attained at k = π.
We close this section by stressing that our work provides a new type of incoherent-coherent dynamical transition in open quantum lattice systems.Previous incoherent-coherent dynamical transition has mostly been studied with respect to the spin-boson model, where the expectation value of the spin variable shows a transition from overdamped to underdamped relaxation [19][20][21][22][23][24][25][26].An important distinction from the spinboson model is that our model has spatially extended degrees of freedom, which play an essential role in the deconfinement of incoherentons.Furthermore, it should be noted that the transition discussed here becomes sharp only in the limit of infinite system size.

V. HIERARCHY OF EIGENMODES
In the following sections, we discuss the generalization of incoherentons and QC gap to many-body systems.In these cases, incoherentons and deconfined particles, in general, coexist.Furthermore, two or more particles can form a single bound state.We refer to such a 2m-particle composite incoherenton as an mth-order incoherenton [see Fig. 7(a)].An mth-order incoherenton can be represented by the m-particle reduced eigenmode, If all particles form a single mth-order incoherenton, the mparticle reduced eigenmode is expected to behave as where ξ con is the confinement length.Figure 7 (b) shows the ladder representation of a typical Liouvillian eigenmode.The eigenmodes of the N-body system can be classified according to how many mth-order incoherentons (1 ≤ m ≤ N) they contain.In this section, we describe a typical scenario for the structure of eigenmodes and spectra, which is expected to be applicable to a broad class of dephasing-type dissipation.
Let J be a parameter of the Hamiltonian that creates quantum coherence between lattice sites, such as the tunneling amplitude between adjacent sites, and let γ be a parameter that describes the strength of dissipation.The balance between J and γ characterizes the competition between the intrachain kinetic energy LH and the interchain interaction Ljump .Figure 8 shows a schematic diagram of the Liouvillian eigenmodes.The structure of the eigenmodes in the ladder representation is depicted in boxes, with circles and squares representing particles and incoherentons, respectively.The eigenmodes are ordered from top to bottom in decreasing order of their decay rates |Re[λ]|.Note that eigenmodes with a larger number of incoherentons decay more slowly.This is because the dissipation suppresses off-diagonal elements of the density matrix and thus coherent eigenmodes with less incoherentons decay faster.
We first consider the case where the dissipation is dominant over the coherent tunneling (J ≪ γ).In general, strong dissipation has the effect of suppressing the coherent evolution of quantum systems, known as the quantum Zeno effect [74][75][76], which has recently attracted considerable attention in the context of AMO systems [77][78][79][80][81][82].As shown in Fig. 8(a), the eigenmodes of an N-particle system are divided into N + 1 groups with different decay rates.Each group of eigenmodes is characterized by the number of particles that do not form incoherentons.The eigenmodes in the group with the largest decay rate have no incoherenton, while in the eigenmodes belonging to the group with the smallest decay rate, all particles form incoherentons. Furthermore, each group contains eigenmodes with various types of incoherenton.For N = 3, the group with the smallest decay rate consists of eigenmodes in which particles form (i) three first-order incoherentons or (ii) one first-order incoherenton and one second-order incoherenton or (iii) a single third-order incoherenton [see the bottom box of Fig. 8(a)].The nth QC gap ∆ (n)  QC is defined by where {λ (n)  α } are the eigenvalues with n unbound pairs [see Fig. 8(a)].
Let us consider what would happen if the coherent tunneling amplitude J is gradually increased while keeping γ fixed (or if the dissipation strength γ is decreased while keeping J fixed).Since J represents the amplitude of the intrachain tunneling LH , the confinement length ξ con (the QC gap ∆ QC ) of an incoherenton increases (decreases) with increasing J.At some critical point J = J c , the confinement length ξ con diverges, and a deconfinement transition of an incoherenton takes place.The QC gap ∆ QC closes at J = J c . Figure 8(b) shows the hierarchy of eigenmodes for J > J c , which forms a continuum where the number of incoherentons can vary continuously with respect to the decay rate.As in the one-particle case, the relationship between ξ con and ∆ QC is expected to be given by Eq. (1).We highlight the uniqueness of our findings in the context of the segment structure of the Liouvillian spectrum for strong dissipation, demonstrated in Fig. 8(a), which has been reported for several open quantum many-body systems in recent literature [81][82][83].Firstly, Ref. [81] primarily relies on a perturbation theory relevant only to strong dissipation.However, the closing of the QC gap or merging of spectral bands is a nonperturbative phenomenon that cannot be adequately captured by perturbation theory.Secondly, Ref. [83] is based on a general concept of locality and employs a randomly constructed Liouvillian.Consequently, it lacks an intuitive picture of the hierarchical structure of the Liouvillian spectrum.In contrast, our incoherenton framework offers a clear physical picture in terms of the number or order of incoherentons, which elucidates the relationship between the hierarchical structure of the spectrum and the eigenmodes associated with each spectral band.
Let us now demonstrate the scenario of Fig. 8 for the dissipative hard-core boson model introduced in Sec.II B. We consider the case of particle number N = 3.For each manybody eigenmode |ρ α ), the number of bound pairs N b,α is defined by Eq. (27).For J = 0, there are four highly degenerate eigenvalues 0, −γ, −2γ, and −3γ.The number of bound pairs Figure 9(a) shows the Liouvillian spectra with J = 0.1, 0.15, 0.2, and 0.25.The colors of the dots represent N b,α for the corresponding eigenmodes |ρ α ).In the presence of a nonzero but small J, there are four bands around 0, −1, −2, and −3.As J increases, the widths of these bands increase, and at J = J c ≃ 0.2, they merge almost simultaneously.27) with L = 6 and 8 for different values of J from 0 to 0.25.For J = 0, N b,α is degenerate at 0, 1, 2, and 3.In the presence of nonzero J, since the coherent hopping mixes the eigenmodes with different N b,α , the values of N b,α distribute over a finite width around 0, 1, 2, and 3.The gaps between these clusters close at J ≃ 0.2, which is identical to the value of J at which the QC gap in the Liouvillian spectrum closes.It should be noted that the critical hopping amplitude J c ≃ 0.2 is slightly shifted from that for the one-particle case J c = γ/4 = 0.25 owing to the interactions among incoherentons and deconfined particles.
It is a nontrivial issue whether the critical hopping amplitude J c , at which the QC gap closes, remains nonzero when the limit of infinite system size is taken at a constant density N/L.In Sec.VI, we will show that a certain class of incoherentons exhibits deconfinement at a value of J independent of the system size.However, this does not mean that J c is generally independent of the system size because J c could depend on spectral bands.Specifically, it is widely believed that an infinitesimally small integrability-breaking perturbation leads to random matrix statistics at the center of the spectrum [84][85][86], implying that J c becomes zero for bands located at the center of the spectrum in the thermodynamic limit.However, even if this is the case, J c for spectral bands near the steady state may remain nonzero in the thermodynamic limit.
Note that N b,α counts the number of the particle pairs forming incoherentons, but does not indicate the order of incoherentons.As illustrated in Fig. 8, the spectral band with the smallest decay rate comprises three types of eigenmodes, which N b,α fails to distinguish.To quantify the number of incoherentons of each order, we introduce incoherenton correlation functions: and so on.The qualitative features of C (s) α (m 1 , ..., m s−1 ) (s = 2, 3, ...) allow the identification of the type of eigenmodes.For instance, let us consider the spectral band with the second smallest decay rate in Fig. 8.It includes two types: (i) one with two first-order incoherentons and (ii) the other with a single second-order incoherenton.While C (2)  α (m) takes an almost constant value for type (i) eigenmodes, it exponentially α (l, m) at l = 0, m = 0, and l = m are set to zero.For eigenmodes with three firstorder incoherentons, C (3)  α (l, m) is delocalized across the full range of l and m.For eigenmodes with one first-order incoherenton and one second-order incoherenton, C (3)  α (l, m) is localized at the edge of the (l, m) space.For eigenmodes with one third-order incoherenton, C (3)  α (l, m) is localized at the corner of the (l, m) space.In panel (b), the eigenmodes in the band with the smallest decay rate are classified on the basis of whether the point (l * , m * ) that maximizes C (3)  α (l, m) is located in the bulk, edge, or corners of the (l, m) space.
decreases with respect to m for type (ii) eigenmodes.Similarly, the eigenmodes within the band with the smallest decay rate can be classified based on the behavior of C (3)  α (l, m). Figure 10(b) shows a spectrum where color differentiates the type of eigenmodes.In particular, color maps of C (3)  α (l, m) for three representative eigenmodes in the band with the smallest decay rate are displayed in Fig. 10(c).These three types of eigenmodes can be distinguished on the basis of whether C (3)  α (l, m) is delocalized, localized on the edge, or localized at the corners.It is worth emphasizing that quantities like N b,α and C (2)  α (m) can be experimentally measured in ultracold atoms on optical lattices through a process involving the interference of two system copies and the enumeration of atoms within each copy [87,88].Detailed experimental protocols are explained in Appendix D.
We here remark on the universal validity of the scenario illustrated in Fig. 8.We expect that it is applicable to systems with local dissipation that satisfies the detailed balance condition.In the context of Markovian open quantum systems, the detailed balance condition is expressed as ρ G L † (A) * = L(ρ G A * ) for any operator A, where ρ G = e −βH /Tr[e −βH ] is the Gibbs state with an inverse temperature β and * denotes the complex conjugation [89].If the system satisfies the detailed balance condition, it relaxes to the equilibrium state described by ρ G .The on-site dephasing L l = b † l b l satisfies this condition with β = 0.In Appendices E and F, we present numerical results for other models, i.e., dephasing Bose-Hubbard model and dephasing hard-core bosons with next-nearest-neighbor hopping.In the Bose-Hubbard model, unlike the hard-core model discussed above, a single lattice site can be occupied by multiple particles.The Liouvillian eigenmodes of this model also show a hierarchical structure depending on the numbers of interchain bound states (incoherentons) and intrachain bound states, and the QC gap in the spectrum closes at a certain value of the hopping amplitude.The results presented in Appendices E and F support the universality of our scenario.
The impact of dissipation violating the detailed balance condition on the hierarchical structure of spectra is nontrivial.In Sec.VII, we explore the scenario involving particle loss and gain and demonstrate that the hierarchical structure is preserved.However, when dissipation induces a current, the hierarchical structure depicted in Fig. 8 can undergo substantial changes.A simple example of such dissipation is realized by the Lindblad operators L l = b † l+1 b l , which describe stochastic hopping induced by external driving [51,90,91].In this case, we observe that the spectrum possesses a topologically distinct structure from the striped band structure shown in Fig. 9, which will be elaborated in a future publication.

VI. EXACT MANY-BODY SOLUTION
While we have shown the numerical results for the manybody system in the previous section, analytical solutions are also possible since the dissipative hard-core boson model can be exactly solved by the Bethe ansatz [92].In the following, we derive an exact many-body solution of the model to show that some Liouvillian eigenmodes have higher-order incoherentons and that they exhibit the deconfinement transition at a critical hopping amplitude.First, it should be noted that the hard-core bosons can be mapped to the XX chain, by using the following identification between the bosonic operators b l and the Pauli matrices σ µ l : Thus, our model can be regarded as a dissipative spin chain with Lindblad operators The boundary condition is set to be periodic: σ µ L+1 = σ µ 1 (µ = x, y, z).In this section, we assume that the system size L is even.
Employing the Jordan-Wigner transformation, we introduce a fermion annihilation operator which satisfies the anticommutation relations {c l , c † m } = δ lm and {c l , c m } = 0, and rewrite the Hamiltonian and the Lindblad operator as and where N = L l=1 c † l c l is the number of fermions.We here use a pseudospin index σ =↑, ↓ to express each chain in the ladder representation of the Liouvillian, and denote the annihilation operator of a fermion on each chain as c l,σ .Then, the ladder representation of the Liouvillian (multiplied by i) where J ↑ = −J ↓ = J, and N σ := L l=1 c † l,σ c l,σ is the number of fermions on each chain.The dependence of the hopping amplitude J σ on the pseudospin can be removed by a unitary transformation and the transformed Liouvillian iU † LU is equivalent to the Hubbard model with an imaginary interaction strength and an imaginary chemical potential [92,93].Since we focus on the Liouvillian L of an N-particle system, we assume N ↑ = N ↓ = N.Note that if N is odd (even), the periodic (antiperiodic) boundary condition is imposed on fermions on each chain.For this reason, we hereafter consider the non-Hermitian Hubbard model under a flux with the periodic boundary condition, c L+1,σ = c 1,σ .If N is odd or even, ϕ is set to 0 or π, respectively.Instead of twisting the particular bond between the sites l = 1 and L, we have introduced a uniform complex hopping amplitude to ensure the translation invariance of the model.The original model with the twist at the particular bond is obtained by performing the gauge transformation c l,σ → exp(iϕl/L)c l,σ on H ϕ .The Hubbard model (46) has the spin SU(2) symmetry [H ϕ , S µ ] = 0 (µ = +, −, z), where In addition, for ϕ = 0 or π, the model ( 46) possesses the η-SU(2) symmetry [94,95] [H ϕ , η µ ϕ ] = 0 (µ = +, −) and [H ϕ , η z ] = 0, where which are generalized in order to incorporate the case of antiperiodic boundary condition [96].These symmetries are called the weak symmetry of the Lindblad equation and lead to a block-diagonal structure of the Liouvillian [97,98].
The one-dimensional Hubbard model ( 46) is exactly solvable by using the Bethe ansatz method [99,100].The Yang-Baxter integrability of the Hubbard model is preserved even when the interaction strength is complex valued [92,93,101].The Bethe equations for the Hubbard model (46) are given by [99,100,102] where k a (a = 1, ..., N ↑ + N ↓ ) is a quasimomentum, Λ α (α = 1, ..., N ↓ ) is a spin rapidity, and u = iγ/(4J) is the pureimaginary dimensionless interaction strength.An eigenvalue λ of L is obtained from a solution of the Bethe equations as A Bethe wave function constructed from a solution of the Bethe equations ( 51) and ( 52) provides a Bethe eigenstate of the Hubbard model (46), which can be interpreted as a Liouvillian eigenmode in the ladder representation.Since Bethe eigenstates satisfy the highest-weight (lowest-weight) condition of the spin SU(2) (η-SU(2)) symmetry, a general eigenstate can be obtained by acting S − or η + ϕ on a Bethe eigenstate [96,100,[103][104][105].Noting the commutation relation [H ϕ , η + ϕ ] = 0, the steady state, H ϕ |Ψ 0 ) = 0, is given by where |v) is the vacuum state of fermions [92].Since η + ϕ creates a bound pair of particles of the two chains, the steady state is composed of N first-order incoherentons.With the unitary transformation (45), the state |Ψ 0 ) is equivalent to the infinite-temperature state of the N-particle sector in the original problem.We note that an incoherenton created by η + ϕ is localized on a rung of the ladder, while an incoherenton in an excited eigenmode can have a nonzero confinement length.
The Bethe equations ( 51) and ( 52) for sufficiently large L allow k-Λ string solutions, in which a part of quasimomenta and spin rapidities forms a string pattern [100,106].Since a k-Λ string solution of length 2m describes a bound state made of m spin-up particles and m spin-down ones [105], it offers an mth-order incoherenton.A k-Λ string of length 2m is composed of 2m quasimomenta k 1 , ..., k 2m and m spin rapidities Λ 1 , ..., Λ m that satisfy and where µ ∈ R is the center of the k-Λ string (see Fig. 11), and we set the branch so that −π/2 < Re[arcsin x] ≤ π/2 [92,93].
The deconfinement of incoherentons is diagnosed from the disappearance of a k-Λ string solution.Let us consider a situation in which all quasimomenta and spin rapidities form a single length-2m k-Λ string solution of the Bethe equations for N ↑ = N ↓ = m.By multiplying the Bethe equations ( 51) for a = 1, ..., N, we obtain exp where we have used Eq. ( 52) and the fact that ϕ is set to 0 or π.Since Eqs. ( 55) and ( 57) imply that k 1 + k 2m is real, we can set without loss of generality.From Eq. ( 55), k 1 and k 2m satisfy which leads to For our case with u = iγ/(4J), the solution of Eq. ( 60) is given by and For a given p = Re[k 1 ], the imaginary part κ is obtained from Eq. ( 62), and by substituting k 1 = p − iκ into sin k 1 = iµ + miu, the center of the string can be calculated as Since cosh κ > 1, the solution of Eq. ( 62) exists only for For mγ/(4J) > 1, arbitrary −π < p < 0 satisfies the above condition.However, for mγ/(4J) ≤ 1, the string solution for some p around −π/2 disappears, indicating the deconfinement of the string solution at The eigenvalue given by Eq. ( 53) can be calculated from Eqs. ( 55) and ( 62) as Thus, the deconfinement of the mth-order incoherenton occurs at Liouvillian eigenmodes with eigenvalues λ = −mγ.
We here define the total momentum K of the mth-order incoherenton by where we have introduced the phase shift mπ to compensate the unitary transformation (45).From Eqs. ( 55), (58), and (61), the total momentum reads K = 2p.Thus, in terms of K, Eq. ( 66) is rewritten as The η-SU( 2) symmetry [H ϕ , η + ϕ ] = 0 yields an eigenstate of H ϕ in the sector of where |ψ 2m ) is a length-2m k-Λ string solution of the Bethe equations.Since η + ϕ generates an on-site pair of particles with opposite spins, |Ψ) is interpreted as a state that involves an mth-order incoherenton and N − m first-order incoherentons.As the action of η + ϕ does not change the eigenvalue, the deconfinement transition of |Ψ) occurs in Liouvillian eigenmodes with eigenvalues near λ = −mγ.Thus, we conclude that the N-body dissipative dynamics governed by the Liouvillian L shows the deconfinement transition of mth-order incoherentons for m = 1, 2, ..., N.
It is worth noting that the deconfinement of bound states does not occur in the ordinary Hermitian Hubbard model with real u.In this case, the solution of Eq. ( 60) is given by and which can be satisfied for any value of p because the range of sinh κ is (−∞, ∞).Thus, the deconfinement transition in the string solution is unique to the dissipative system which can be mapped to the non-Hermitian Hubbard model with an imaginary interaction strength.

VII. EFFECTS OF PARTICLE LOSS AND GAIN
In the system of hard-core bosons subject to on-site dephasing, the total number of particles is conserved during time evolution.A question arises as to whether the incoherenton framework is applicable to situations where particle exchange with the environment occurs.Such situations appear in cases like driven optical cavities [39] and exciton-polariton systems [38].In the following, we confirm that the incoherenton framework essentially holds in the presence of particle loss and gain, at least for small loss and gain rates.
We incorporate particle loss and gain by considering the following Liouvillian: where H is the Hamiltonian ( 14) of hard-core bosons, and κ 1 and κ 2 represent the rates of particle loss and gain, respectively.When the hard-core boson model is mapped to a spin model, the loss and gain terms in Eq. ( 72) correspond to dissipative processes that flip spins down and up at rates κ 1 and κ 2 , respectively.Note that the loss and gain terms in Eq. ( 72) mix sectors with different particle numbers.The ladder representation of Eq. ( 72) can be expressed as where Ld is the Liouvillian (17) of dephasing hard-core bosons, and N = l (n l,+ + n l,− )/2 is the total particle number.
Firstly, let us consider the scenario with particle loss but without gain, i.e., κ 2 = 0. Importantly, L can be expressed in a block-upper-triangular form because L1 reduces the particle number but does not increase it.The eigenvalues of a blocktriangular matrix are given by those of its diagonal blocks.Consequently, the spectrum of L is simply the union of spectra of each particle sector: where L(N) represents Ld in the N-particle sector.This is a general property of Liouvillian with loss but without gain [101,107,108].Equation (75) implies that the singular behavior of the spectra of Ld , linked to the deconfinement of incoherentons, is directly transferred to the spectra with loss.Thus, the presence of particle loss does not affect the incoherenton picture.
Next, let us consider the situation with both particle loss and gain.We numerically diagonalize L in the subspace with l n l,+ = l n l,− .We calculate the total particle number N and the number N b of interchain bound states, as defined by Eq. ( 27), for each eigenmode.Figure 12 shows the scatter The QC gap closing in the spectrum can also be observed in Fig. 13.Note that the data in Figs. 12 and 13 primarily deal with situations where loss and gain rates are sufficiently smaller than the dephasing rate.Extending the incoherenton framework to the case with strong loss and gain deserves further study.

VIII. EFFECTIVE DESCRIPTION OF MANY-BODY DECOHERENCE
In Fig. 8 of Sec.V, we have shown that the Liouvillian eigenmodes are arranged in a hierarchy characterized by incoherentons.Each group in this hierarchy has a different decay rate and is separated from each other by the QC gap when the dissipation dominates.In this section, we discuss the consequences of this hierarchy of Liouvillian eigenmodes for the process of quantum decoherence.First, we introduce multiorder quantum coherence as a quantitative measure of how many incoherentons a given density matrix contains.The time evolution of the quantum coherence is investigated for the dissipative hard-core boson model by numerically solving the master equation.We argue that the decay process of the quantum coherence can be understood in terms of the production, diffusion, and localization of incoherentons.

A. Multi-order quantum coherence
First, the general concept of "quantum coherence" is outlined according to Refs.[109] and [110].We consider the dissipative hard-core boson model introduced in Sec.II B. Similarly to Eq. ( 33), we define the one-particle reduced density matrix G (1) by where b † l and b l are the creation and annihilation operators of a boson at site l.Then, the state ρ is said to have the first-order coherence if G (1) satisfies the following relation: for any l 1 and l 2 , where n l = b † l b l is the number operator at site l.In other words, a strong correlation between distantly separated points exists in a coherent state.Similarly, the twoparticle reduced density matrix G (2) is defined by The state ρ is said to have the second-order coherence if G (1)  and G (2) satisfy Eq. ( 77) and for any l 1 , l 2 , l 3 , and l 4 .The notion of the sth-order coherence (s ≥ 3) can also be defined from the s-particle reduced density matrix G (s) in a similar manner.We define the amount of the first-order coherence by which simply measures the amount of off-diagonal components of G (1) .Since the steady state of the model is the infinite-temperature state ρ ss = D −1 I due to dephasing, all off-diagonal components of ρ vanish in the long-time limit, and then we have lim t→∞ χ 1 (t) = 0.If the state ρ has the first-order coherence, e.g., a Bose-condensed pure state, the amount of the first-order coherence is given by χ 1 ∝ NL because |G (1)  l 1 ,l 2 | ∝ N/L from Eq. ( 77).It should be noted that the expectation value of an arbitrary one-body observable, can be written as In particular, G (1) is related to the momentum distribution of particles, which is accessible by time-of-flight experiments in ultracold atomic gases.We also define the amount of the second-order coherence χ 2 by taking the summation of G (2)  l 1 ,l 2 ;l 3 ,l 4 over all off-diagonal indices: Note that G (2) l,l;l 3 ,l 4 = G (2) l 1 ,l 2 ;l,l = 0 due to the hard-core condition.If the state ρ has the second-order coherence, the amount of the second-order coherence is given by χ 2 ∝ N 2 L 2 .The amount of the higher-order coherence χ s (s = 3, 4, ...) can also be defined from the s-particle reduced density matrix G (s) in a similar manner.
When dissipation is dominant (J ≪ γ), the Liouvillian eigenmodes are arranged in N + 1 bands with different decay rates sγ (s = 0, ..., N), as shown in Fig. 8 (a).We denote the set of eigenmodes that belong to each band as {ρ 0,α }, {ρ 1,α }, ..., {ρ N,α }.Each ρ s,α involves s deconfined pairs in the ladder representation.Then, the eigenmode expansion of the density matrix can be rearranged as where S r denotes the set of indices for {ρ r,α }.Note that the steady state ρ ss = D −1 I belongs to {ρ 0,α }.We refer to Eq. ( 84) as the hierarchical expansion of the density matrix.The decay rate of each ρ s,α is given by sγ + O((J/γ) 2 ).From the definition, the dominant contribution to the amount of the sth-order coherence χ s comes from ρ s,α , because the s-particle reduced density matrix G (s) of eigenmodes with s deconfined pairs has large off-diagonal components.Thus, when J ≪ γ, χ s initially decays as In the next subsection, it is argued that the initial decay of χ s is due to the generation of incoherentons and that the relaxation of χ s at long times is characterized by the localization and diffusion of incoherentons.

B. Numerical results for relaxation of quantum coherence
By solving the quantum master equation numerically, we investigate the time evolution of χ s for the dissipative hardcore boson model.We take a random pure state as an initial state, i.e., ρ ini = |ψ r ⟩ ⟨ψ r | where |ψ r ⟩ is a normalized vector uniformly sampled from the set of unit vectors in the Hilbert space.Figure 14(a) shows the time evolution of the absolute values of the one-particle reduced density matrix G (1) .Two regimes can be clearly distinguished.In the first regime, the off-diagonal components of G (1) decay rapidly, which implies the production of incoherentons.In the second regime, a slow diffusion of the diagonal components is observed, and at long  (1) with γ = 1 and J = 0.1.The system size is L = 10 and the particle number is N = 3.The initial state is a random pure state.In the regime of production of incoherentons, a fast decay of the off-diagonal components can be seen.In the diffusion regime of incoherentons, on the other hand, a slow diffusion of the diagonal components occurs, leading to an infinite-temperature state in the long-time limit.To highlight the variation of the off-diagonal components, they are multiplied by a factor of 4, i.e., these plots represent (4 − 3δ l,m )|G (1)  l,m |.(b) Time evolution of Γ 1 and Γ 2 with γ = 1 and J = 0.05, 0.1, 0.15, 0.2, and 0.25 from bottom to top at t = 10.The horizontal and vertical axes are plotted on a logarithmic scale.Γ s is calculated from χ s averaged over 100 initial random states.The system sizes are L = 10 and 12, and the particle number is N = 3.The dotted lines show 1/t for Γ 1 and 2/t for Γ 2 .The arrows indicate the beginning and end of the localization regime (τ 1 and τ 2 ) of incoherentons for J = 0.1.
times they converge to N/L.We refer to the first (second) regime as the incoherenton production (diffusion) regime.
It is convenient to define the time-dependent decay rate Γ s of the sth-order coherence as For J = 0, we have Γ s (t) = sγ for all t. Figure 14(b) shows the time evolution of Γ 1 and Γ 2 with dephasing γ = 1.For a small hopping amplitude such as J = 0.05 or 0.1, two plateaus and subsequent algebraic decay Γ s ∼ 1/t are observed.(For J = 0.1, the beginning and end of the second plateau are indicated by arrows.)The height of the first plateau is sγ, which is the initial decay rate of χ s .These numerical data suggest the existence of three regimes for J ≪ γ : with some intermediate decay rate κ s and exponent η s .The algebraic decay is fitted roughly by 1/t for Γ 1 and 2/t for Γ 2 [the dotted lines in Fig. 14(b)], which implies η 1 ≃ 1 and η 2 ≃ 2.
The first and third regimes of Eq. ( 87) correspond to the production and diffusion regimes, respectively.For reasons that will be explained in the next subsection, we refer to the second regime of Eq. ( 87) as the localization regime of incoherentons.As J approaches the transition point J c ≃ 0.2 where the QC gap closes, the second plateau shrinks and eventually disappears at J = J c , leaving a small bump.
It should be noted that the power-law decay of χ s in the diffusion regime does not last forever in a finite system.For t ≫ ∆ −1 L , where ∆ L is the Liouvillian gap, the relaxation of χ s is determined by the slowest eigenmode, and thus χ s decays as e −∆ L t .It is known that the Liouvillian gap closes as ∆ L ∼ L −2 for our model [92].In Fig. 14(b), a third plateau of Γ s appears in the long-time regime (see, e.g., the data of J = 0.25 at t > 30).We can confirm that the height of this plateau scales as L −2 .The difference in the height of the third plateau for Γ 1 and Γ 2 is due to the difference in the slowest decaying eigenmodes which contribute to χ 1 and χ 2 .
We present a simple theoretical argument showing the existence of the three relaxation regimes, i.e., regimes described by incoherenton production, localization, and diffusion.We assume that the relaxation of χ s is given by a superposition of exponential functions: where D s (µ) is a weighted density of states, which expresses how many eigenmodes with decay rate µ contribute to χ s .More precisely, Eq. ( 88) can be obtained by substituting the hierarchical expansion (84) into the definition of χ s and replacing the sum over eigenmodes with an integral over the decay rate.For simplicity, we focus on the case of s = 1.Let us consider the following D 1 (µ): Figure 15(a) shows a schematic illustration of D 1 (µ).The support of D 1 (µ) near µ = 0 represents the contribution from incoherent eigenmodes where all particles form incoherentons.

C. Characterization of many-body decoherence by incoherentons
The numerical results on the relaxation processes in the previous subsection can be explained in terms of the dynamics of incoherentons, as summarized in Fig. 16.When the dissipation dominates (J < J c ), three distinct relaxation processes emerge: Incoherenton production: (t < τ 1 ∼ γ −1 ln(Lγ/J)) In this regime, the coherence decays exponentially as χ s (t) ∼ e −sγt .Thus, Γ s has a plateau of height sγ [see Fig. 16(c-2)].Since eigenmodes with a smaller number of incoherentons decay faster, the number of incoherentons increases with time in this regime.Let us now estimate the timescale τ 1 at which the deviation from χ s (t) ∼ e −sγt begins.The magnitude of the contribution of the incoherent eigenmodes ρ 0,α to χ s can be estimated as L s (J/γ) s from perturbation theory [111].
When the contributions of ρ s,α and ρ 0,α to χ s are comparable, the incoherenton production regime ends.This condition is expressed as L 2s e −sγτ 1 ∼ L s (J/γ) s , where the factor L 2s comes from the sum over off-diagonal components in, e.g., Eqs. ( 80) and (83).
We denote the eigenvalue with the smallest real part among the incoherent eigenmodes {ρ 0,α } as λ * inc .In this regime, the relaxation of χ s is determined by eigenmodes with decay rates of O(|λ * inc |), and thus, it decays as χ s ∼ e −κ s t with κ s = O(|λ * inc |), which, in general, depends on s.That is, Γ s has a plateau of height κ s [see Fig. 16(c-2)].In general, eigenmodes that contain incoherentons with larger confinement length ξ con decay faster.Thus, the decay of eigenmodes with eigenvalues of O(λ * inc ) leads to a reduction of ξ con , i.e., the localization of incoherentons.Incoherenton diffusion: (τ 2 < t) In this regime, the relaxation of χ s is determined by the incoherent eigenmodes {ρ 0,α } with small decay rates (≪ |λ * inc |).In the hard-core boson model under on-site dephasing, since there exist eigenvalues arbitrarily close to 0 for an infinitely large system, the coherence exhibits a power-law decay χ s ∼ t −η s .Note that algebraic decay is a general feature of open quantum manybody systems in which the Liouvillian gap ∆ L vanishes in the thermodynamic limit [47].This relaxation process proceeds by rearrangement of the positions of welllocalized incoherentons, i.e., the diffusion of incoherentons.In a system with a nonzero ∆ L in the thermodynamic limit, this regime is expected to be absent .
When the QC gap closes (J > J c ), we have τ 1 ∼ γ −1 ln L and τ 2 ∼ γ −1 , and thus, the incoherenton localization regime cannot be observed.Instead, the incoherenton diffusion regime directly follows the incoherenton production region [see Figs.16(b-1) and (c-1)].In this case, incoherenton production and localization occur simultaneously, and these regimes cannot be clearly separated.
As mentioned above, the crossover timescale τ 1 depends logarithmically with respect to the system size L.This is because the sum in Eq. ( 80) or ( 83) is taken over all off-diagonal components.It may be reasonable to restrict this sum to the off-diagonal components that are close to the diagonal components.For a typical one-body observable (81), O (1) l 1 ,l 2 rapidly decays to zero as |l 1 − l 2 | increases, so only off-diagonal components satisfying l 1 ≃ l 2 contribute to Eq. ( 82).Then, we can also define χ1 := where the exponential factor with c = O(1) suppresses the contribution of off-diagonal components with large |l 1 − l 2 |.
If we focus on the relaxation of χ1 , the crossover timescale τ 1 can be given by γ −1 ln(γ/J), which is independent of the system size.

IX. CONCLUSION
We have proposed the notion of incoherenton in open quantum many-body systems, which characterizes the hierarchical structure of Liouvillian eigenmodes and their incoherentcoherent transitions.Under the mapping of the Liouvillian to a non-Hermitian ladder Hamiltonian, incoherentons are defined as interchain bound states.The decay rate of each eigenmode is determined approximately by the number of incoherentons that the relevant eigenmode involves.The quantum coherence (QC) gap is defined as the minimum difference in decay rates between eigenmodes with different numbers of incoherentons.As the coherence parameter of the system increases, the deconfinement of an incoherenton occurs at a certain critical point, causing the QC gap closing.For a hard-core boson system under on-site dephasing, we have demonstrated numerically and analytically the deconfinement of incoherentons.Furthermore, the process of many-body decoherence is discussed in terms of incoherentons.Three relaxation regimes corresponding to the production, localization, and diffusion of incoherentons are identified.
Note that our framework of incoherentons may not suffice to capture every intricate detail of the Liouvillian spectra and eigenmodes.As highlighted in Appendix E, the dephasing Bose-Hubbard model incorporates intrachain bound states, complementing the role of interchain bound states, i.e., incoherentons.However, our primary objective is to describe and explore universal characteristics found within the complex behaviors of the Liouvillian spectra and eigenmodes, which are summarized as follows: 1.Under strong dissipation, the spectrum displays multiple bands, with eigenmodes distinctly marked by both inter-and intra-chain bound states.
2. Interchain and intrachain bound states affect decay rates and frequencies, respectively, thereby governing the temporal dynamics of eigenmodes.
3. Decreasing dissipation leads to the merging of specific bands, signaling the deconfinement of interchain bound states.
While we believe in the robustness of these observations, these characteristics may require refinement or further extension in more complicated situations.The presence of additional degrees of freedom could introduce other types of bound states.We expect that our current work provides a solid foundation for subsequent research that further refines the quasiparticle descriptions of Liouvillian eigenmodes.
In this study, we focused on systems with local bulk dissipation.It is natural to ask whether the incoherenton picture summarized above holds for other types of dissipation as well.When there is dissipation only at the boundary of the system [49], the localization of incoherentons near the boundary is expected.The effect of nonlocal dissipation is also nontrivial, which leads to long-range interactions between the chains in the ladder representation of a Liouvillian.The influence of long-range interactions on the formation of bound states is well studied in the context of isolated quantum many-body systems, and it has been pointed out that a new type of bound state can be realized in quantum spin chains with long-range interactions [112][113][114].Understanding the impact of different types of dissipation on the deconfinement transition of eigenmodes and the QC gap closing in spectra deserves further study.
Quasiparticles are a key concept in many-body physics.In isolated quantum many-body systems, the existence of welldefined quasiparticles ensures the validity of low-energy effective field theories, which describe the thermodynamic and transport properties of the system through statistical mechanics of weakly interacting quasiparticles.Identifying the quasiparticles is to distinguish a set of relevant variables for characterizing the low-energy behavior of the system from many irrelevant variables.It is expected that complex relaxation processes in open quantum many-body systems are described by a simple kinetic theory of various incoherentons, which should be studied in detail in future works.A better understanding of incoherentons can provide an efficient way to predict decoherence effects in the control of large-scale quantum devices.
The formation of bound states between interacting particles is a universal phenomenon from particle physics to condensed matter physics.Phenomena that arise from the formation of specific types of bound states include, for example, BCS-BEC crossover in interacting Fermi gases [115][116][117] and Efimov resonances in three-body bound states of atoms with large scattering lengths [118][119][120][121].It is an important task to investigate the effect of the formation of various types of bound states in Lindblad ladder systems on the structure of Liouvillian spectra and dynamical features in open quantum systems.

If
α ).This implies that (i) the Liouvillian spectrum on the complex plane is symmetric with respect to the real axis, and (ii) if ρ α is Hermitian, the corresponding eigenvalue λ α is real.
Appendix B: One-particle solution of Liouvillian eigenmodes In this Appendix, we present a detailed analysis of the Liouvillian spectrum and eigenmodes in the one-particle case without resorting to the Bethe ansatz.In the ladder representation, the Liouvillian L is mapped to a non-Hermitian Hamiltonian L of a two-particle system on the ladder.We denote by |l⟩ ⊗ |m⟩ the state in which each particle is located at sites l and m of each chain of the ladder.Since L is translationally invariant, it is convenient to introduce a basis where k = 2πs/L (s = −L/2 + 1, ..., L/2) is the momentum of the center of mass and l = −L/2 + 1, ..., L/2 is the relative coordinate.We write matrix elements in this basis as where the matrix elements between different momenta k vanish owing to the translational symmetry of L. The matrix elements of L(k) are given by Note that the indices l and m satisfy the periodic boundary condition.Equation (B3) defines an effective tightbinding model for the relative coordinate, which has an imaginary hopping amplitude between neighboring sites.Let {λ j (k)} j=1,...,L be the eigenvalues of L(k).For k = 0, L(k) becomes diagonal for arbitrary J, and it has a single zero eigenvalue and an (L − 1)-fold degenerate eigenvalue λ = −γ.In particular, the zero mode is given by ρ 0,lm = L −1/2 δ lm .
In terms of the effective tight-binding model (B3), the coherent and incoherent eigenmodes correspond to scattering and bound states, respectively.We write an eigenmode as From Eqs. (B6)-(B9), α, β, and λ can be calculated as functions of k.
(B10) Substitution of Eq. (B10) into Eq.(B9) yields for L → ∞.From Eqs. (B7), (B8), and (B11), the eigenvalue λ inc (k) associated with the bound state can be calculated as which coincides with Eq. ( 68) with m = 1.Figures 17(a) and (b) show Re[α] and Im[α] as functions of k, respectively.Note that Re[α] is the inverse of the confinement length ξ con .For k = 0, since the zero mode of L(k) is completely localized at l = 0, we have Re[α] = ∞.Figure 17(a) shows that there is a critical hopping amplitude J c below which a bound state exists (Re[α] > 0) for all k.On the other hand, for J > J c , there is a critical wavenumber k c (J) such that the bound state disappears (Re[α] = 0) for k ≥ k c (J).From Fig. 17 (b), one finds that Im[α] = −k/2 for J < J c , which is consistent with the fact that λ inc (k) is real.For J > J c , Im [α] shows a cusp at k = k c (J) owing to the disappearance of the bound state.Figure 17(c) shows trajectories of λ inc (k) for different values of J.As k increases from −π, λ inc (k) initially moves from the left end point of the incoherent-mode spectrum to the right, reaches the origin for k = 0, and returns to the left end point for k = π.The quantum coherence gap ∆ QC decreases with increasing J, and eventually, it closes at J = 0.25.Figure 17(c) is consistent with the numerical results shown in Fig. 5(a).
Let us determine the critical hopping amplitude J c at which the confinement length ξ con diverges.As J increases, the disappearance of the bound state firstly occurs at k = π.Thus, from Eqs. (B7), (B8), and (B11), we have α = β and The hyperbolic sine function has the following property: if sinh α is purely imaginary, | sinh α| ≤ 1 and | sinh α| > 1 imply Re[α] = 0 and Re[α] 0, respectively.Thus, we have which agrees with Eq. ( 65) with m = 1.It should be noted that the incoherent-mode eigenvalue λ inc (π) given by Eq. (B12) becomes complex for J > J c .In Fig. 17, J c = 0.25 since γ = 1.From Eq. (B13), the confinement length ξ con of the eigenmode with k = π reads for J ≃ J c .Since the QC gap ∆ QC is given by Eq. ( 25), we have which is consistent with Eq. (1). Figure 17(d) shows the Liouvillian spectrum in the limit of L → ∞.It consists of the coherent-mode spectrum at Re[λ] = −γ parallel to the imaginary axis and the incoherentmode spectrum on the real axis.The reason why the real parts of the coherent-mode eigenvalues are equal to −γ in the limit of L → ∞ can be understood as follows.Since the coherent eigenmodes are extended over the relative coordinate, Re[α] should vanish.By substituting Re[α] = 0 into Eq.(B7), we obtain Re[λ(k)] = −γ.By numerical diagonalization of the Liouvillian, one can also verify that the width of the coherentmode spectrum along the real axis decreases as L increases.

Appendix C: Absence of incoherenton in continuous systems
As mentioned at the end of Sec.III, the spatial discreteness of the lattice system is crucial for the formation of incoherentons.In this Appendix, we show that incoherenton does not exist in systems where a free particle in continuous space undergoes dephasing.
The Hamiltonian of a free particle in one-dimensional continuous space is given by where ψ † (x) and ψ(x) are the creation and annihilation operators of a boson at position x, which satisfy the canonical commutation relations: [ψ(x), ψ(x ′ )] = [ψ † (x), ψ † (x ′ )] = 0, [ψ(x), ψ † (x ′ )] = δ(x − x ′ ).(C2) The Liouvillian L that governs the time evolution of the density matrix ρ is given by We consider the following Lindblad operator where g(x) is a short-ranged function that rapidly decays for large |x|.The Lindblad operator given by Eq. (C4) describes a dephasing process of a particle near position x.
We focus on the one-particle sector of the Hilbert space.Let |x⟩ = ψ † (x) |v⟩ (|v⟩: vacuum state) be the state in which a particle is located at position x.Then, {|x⟩} x∈(−∞,∞) is an orthonormal basis in the one-particle sector.In terms of this basis, the density matrix ρ is written as where ρ(x, y) is the matrix element of ρ.In the ladder representation, |x⟩ ⟨y| is mapped to a tensor-product state |x⟩ ⊗ |y⟩, which specifies a two-particle state in a ladder.The Liouvillian L is also mapped to a non-Hermitian operator The matrix element of L is calculated as axis (frequencies) specified by intrachain bound states.While the correspondence between decay rates (frequencies) and the number of incoherentons (intrachain bound states) is not perfectly one-to-one, such a hierarchical picture provides a qualitative understanding of the spectral structure.
Appendix F: Dephasing hard-core bosons with next-nearest-neighbor hopping In this Appendix, we present the results for dephasing hardcore bosons with next-nearest-neighbor hopping to test the robustness of the incoherenton picture to integrability-breaking perturbations.The Hamiltonian of this system is given by where J 1 and J 2 denote the hopping amplitudes between nearest-neighbor and next-nearest-neighbor sites, respectively.The next-nearest-neighbor hopping introduces the simplest perturbation that breaks the integrability of the original model without altering the number of spectral bands.In the following calculation, we assume J := J 1 = J 2 .Figure 19(a) presents the Liouvillian spectra for J = 0.1, 0.15, 0.2, and 0.25.The colors of the dots represent N b,α , defined by Eq. (27).With a small, nonzero J, four spectral bands   3)  α (l, m) at l = 0, m = 0, and l = m are set to zero.For eigenmodes with three first-order incoherentons, C (3)  α (l, m) is delocalized across the full range of l and m.For eigenmodes with one first-order incoherenton and one second-order incoherenton, C (3)  α (l, m) is localized at the edge of the (l, m) space.For eigenmodes with one third-order incoherenton, C (3)  α (l, m) is localized at the corner of the (l, m) space.In panel (b), the eigenmodes in the band with the smallest decay rate are classified based on whether the point (l * , m * ) that maximizes C (3)  α (l, m) is located in the bulk, edge, or corners of the (l, m) space.
emerge around 0, −1, −2, and −3.As J increases, the widths of bands increase, and they merge almost simultaneously at J = J c ≃ 0.2.This behavior is the same as the case without next-nearest-neighbor hopping [see Fig. 9(a)].
Each eigenmode can be classified according to the number N b,α of incoherentons [Eq.(27)] and the qualitative behavior of incoherenton correlations C (s)  α (m 1 , ..., m s−1 ) (s = 2, 3, ...) defined by Eqs.(36) and (37).Figure 19(b) depicts a spec-trum where color variations differentiate eigenmode types.Additionally, color maps of C (3)  α (l, m) for three representative eigenmodes in the band with the smallest decay rate are shown in Fig. 19(c).The three types of eigenmodes can be distinguished by whether C (3)  α (l, m) delocalized, on the edge, or localized at the corners.All features presented in Fig. 19 are consistent with those in Fig. 10, presenting further evidence that our incoherenton picture is applicable to nonintegrable systems.

FIG. 1 .
FIG. 1. Quasiparticle description of relaxation processes in terms of incoherentons.(a) An incoherenton is a bound state between the degrees of freedom of ket and bra spaces of a density matrix.(b) The incoherent-coherent transition of eigenmodes can be characterized by the deconfinement of incoherentons (left panel) and the closing of the quantum coherence (QC) gap ∆ QC (middle panel), where λ is the eigenvalue of the Liouvillian.The QC gap closing causes a dynamical transition from incoherent exponential relaxation to coherent oscillatory relaxation (right panel).(c) The many-body eigenmodes are classified into groups according to the number of incoherentons involved.Each group of eigenvalues is separated from the others by the QC gaps.Since the system loses coherence over time, the more coherent the mode, the larger the decay rate.(d) The relaxation dynamics of open quantum many-body systems is effectively described by the production, localization, and diffusion of incoherentons.
coherentons significantly alters the transient dynamics of open quantum systems, where incoherent-coherent transitions are expected to take place [see the right panel of Fig. 1(b)].

FIG. 3 .
FIG. 3. (a) Schematic illustration of a system of hard-core bosons with on-site dephasing.(b) Ladder representation of the model.While coherent hopping acts on individual particles, the on-site dephasing acts on a particle pair occupying the same rung (vertical dashed line).(c) Physical implementation of on-site dephasing for an atom in an optical lattice.The double arrow shows the Rabi coupling Ω induced by a laser with frequency ω L .The wavy arrow shows spontaneous decay with rate Γ s .∆ = ω L − ω eg is the detuning of a laser, where ω eg is the excitation energy of the atom.

FIG. 4 .
FIG.4.Schematic illustrations of coherent eigenmodes and incoherent eigenmodes for a one-particle system.We call an interchain bound state on the ladder as an incoherenton, for which the confinement length is denoted by ξ con .The gray scale shows the magnitude of each matrix elements, where the darker one shows the larger magnitude.The coherent eigenmode has both diagonal and off-diagonal matrix elements that are comparable in magnitude, whereas the incoherent eigenmode has the predominant diagonal matrix elements.

FIG. 8 .
FIG.8.Hierarchy of eigenmodes for a three-particle case.The eigenmodes are ordered from top to bottom in decreasing order of their decay rates |Re[λ]|.The lower end of the most incoherent band represents the steady state.For (a) with J < J c , the groups of eigenmodes with different decay rates are separated by the QC gap ∆ QC .For (b) with J > J c , the QC gap closes and all groups are continuously connected.

FIG. 10 .
FIG. 10.Classification of eigenmodes of the dephasing hard-core bosons.(a) Complete spectrum with γ = 1 and J = 0.15.The system size is L = 10 and the particle number is N = 3.(b) Part of the spectrum corresponding to the red square in panel (a).The color differentiates between the types of eigenmodes, which are illustrated in the insets.(c) Color maps of C (3)α (l, m) for three representative eigenmodes with eigenvalues λ = −0.314,−0.194, and −0.03.To emphasize off-diagonal components, the values of C(3)  α (l, m) at l = 0, m = 0, and l = m are set to zero.For eigenmodes with three firstorder incoherentons, C(3)  α (l, m) is delocalized across the full range of l and m.For eigenmodes with one first-order incoherenton and one second-order incoherenton, C(3)  α (l, m) is localized at the edge of the (l, m) space.For eigenmodes with one third-order incoherenton, C(3)  α (l, m) is localized at the corner of the (l, m) space.In panel (b), the eigenmodes in the band with the smallest decay rate are classified on the basis of whether the point (l * , m * ) that maximizes C(3)  α (l, m) is located in the bulk, edge, or corners of the (l, m) space.

FIG. 11 .
FIG. 11.Schematic illustration of a k-Λ string solution in the complex plane of sin k.

FIG. 12 .
FIG. 12. Scatter plots of the particle number N and the number N b of interchain bound states for eigenmodes of dephasing hard-core bosons with loss and gain.The system size is L = 8.Panels (a) through (d) represent varying hopping amplitudes: (a) J = 0.1, (b) J = 0.15, (c) J = 0.2, and (d) J = 0.25.The rates of dephasing, loss, and gain are set to γ = 1, κ 1 = 0.02, and κ 2 = 0.01, respectively.Histograms of N and N b are also included.For smaller values of J, clusters are observed at integer values of N and N b .As J increases, these clusters stretch along the N b axis, eventually merging around J ≃ 0.2.

FIG. 14 .
FIG.14.Relaxation of the reduced density matrix and the amount of the quantum coherence.(a) Time evolution of the absolute values of G(1) with γ = 1 and J = 0.1.The system size is L = 10 and the particle number is N = 3.The initial state is a random pure state.In the regime of production of incoherentons, a fast decay of the off-diagonal components can be seen.In the diffusion regime of incoherentons, on the other hand, a slow diffusion of the diagonal components occurs, leading to an infinite-temperature state in the long-time limit.To highlight the variation of the off-diagonal components, they are multiplied by a factor of 4, i.e., these plots represent (4 − 3δ l,m )|G(1)  l,m |.(b) Time evolution of Γ 1 and Γ 2 with γ = 1 and J = 0.05, 0.1, 0.15, 0.2, and 0.25 from bottom to top at t = 10.The horizontal and vertical axes are plotted on a logarithmic scale.Γ s is calculated from χ s averaged over 100 initial random states.The system sizes are L = 10 and 12, and the particle number is N = 3.The dotted lines show 1/t for Γ 1 and 2/t for Γ 2 .The arrows indicate the beginning and end of the localization regime (τ 1 and τ 2 ) of incoherentons for J = 0.1.

FIG. 16 .
FIG. 16.Summary of the effective description of many-body decoherence.The left panels show the case of J > J c , while the right panels show the case of J < J c .(a) Schematic illustrations of the Liouvillian spectra, where λ * inc denotes the incoherent-mode eigenvalue with the smallest real part.(b) ln χ s as a function of time t.For J > J c , χ s initially decays as e −sγt , and at long times t ≫ γ −1 , it exhibits a power-law behavior t −η .For J < J c , there is an intermediate regime where χ s decays as e −|λ * inc |t .(c) ln Γ s as a function of ln t.The early regime where χ s decays exponentially and the late regime where it decays with a power-law correspond to the incoherenton-production and incoherenton-diffusion regimes, respectively.The intermediate regime is characterized by the localization of incoherentons.

FIG. 17 .
FIG. 17. (a), (b) Real and imaginary parts of α as functions of k with γ = 1 for J = 0.15, 0.2, 0.25, 0.3, and 0.35.(c) Trajectories of λ inc (k) on the real axis as k increases from −π to π. J increases from 0 to 0.25 in increments of 0.025.(d) Schematic illustration of the Liouvillian spectrum in the limit of L → ∞.The purple and blue bands represent the incoherent-mode and coherent-mode spectra, respectively.

(∂y 2 + 25 FIG. 18 .
FIG. 18. Three-particle Liouvillian spectra and eigenmodes of the Bose-Hubbard model under dephasing.Schematic illustrations of the spectra and scatter plots of (N inter , N intra ) are shown for the gapped case (a) and the gapless case (b).Panel (c) illustrates the eigenmodes indicated in panels (a) and (b), where the red (blue) squares represent an interchain (intrachain) bound state.In panels (a) and (b), series of spectral bands connected by breaking or creating a single first-order incoherenton are highlighted by the green, purple, orange, and black squares.The QC gaps are defined as the gaps between spectral bands belonging to the same series.In the gapless case, the QC gaps close and all bands belonging to each series merge into a single band.Panels (d) and (e) show the spectra and scatter plots of (N inter , N intra ) with J = 0.1.0.15, 0.2, and 0.25.The system size is L = 8 and the particle number is N = 3.The other parameters are set to U = 4 and γ = 1.The double arrows in panel (d-1) represents the QC gaps corresponding to each series of bands.

FIG. 19 .
FIG. 19.Liouvillian spectra for the dephasing hard-core bosons with next-nearest-neighbor hopping.(a) Spectra with γ = 1 and J = 0.1, 0.15, 0.2, and 0.25, where J := J 1 = J 2 .The systems size is L = 10 and the particle number is N = 3.The eigenvalues satisfying 0 ≤ N b,α < 3/4, 3/4 ≤ N b,α < 3/2, 3/2 ≤ N b,α < 9/4, and 9/4 ≤ N b,α ≤ 3 are colored by blue, light blue, green, and red, respectively.(b) Part of the spectrum corresponding to the square in panel (a-2).The color differentiates between the types of eigenmodes, which are illustrated in the insets.(c) Color maps of C (3)α (l, m) for three representative eigenmodes with eigenvalues λ = −0.543,−0.410, and −0.131.To emphasize off-diagonal components, the values of C(3)  α (l, m) at l = 0, m = 0, and l = m are set to zero.For eigenmodes with three first-order incoherentons, C(3)  α (l, m) is delocalized across the full range of l and m.For eigenmodes with one first-order incoherenton and one second-order incoherenton, C(3)  α (l, m) is localized at the edge of the (l, m) space.For eigenmodes with one third-order incoherenton, C(3)  α (l, m) is localized at the corner of the (l, m) space.In panel (b), the eigenmodes in the band with the smallest decay rate are classified based on whether the point (l * , m * ) that maximizes C(3)  α (l, m) is located in the bulk, edge, or corners of the (l, m) space.
FIG. 19.Liouvillian spectra for the dephasing hard-core bosons with next-nearest-neighbor hopping.(a) Spectra with γ = 1 and J = 0.1, 0.15, 0.2, and 0.25, where J := J 1 = J 2 .The systems size is L = 10 and the particle number is N = 3.The eigenvalues satisfying 0 ≤ N b,α < 3/4, 3/4 ≤ N b,α < 3/2, 3/2 ≤ N b,α < 9/4, and 9/4 ≤ N b,α ≤ 3 are colored by blue, light blue, green, and red, respectively.(b) Part of the spectrum corresponding to the square in panel (a-2).The color differentiates between the types of eigenmodes, which are illustrated in the insets.(c) Color maps of C (3)α (l, m) for three representative eigenmodes with eigenvalues λ = −0.543,−0.410, and −0.131.To emphasize off-diagonal components, the values of C(3)  α (l, m) at l = 0, m = 0, and l = m are set to zero.For eigenmodes with three first-order incoherentons, C(3)  α (l, m) is delocalized across the full range of l and m.For eigenmodes with one first-order incoherenton and one second-order incoherenton, C(3)  α (l, m) is localized at the edge of the (l, m) space.For eigenmodes with one third-order incoherenton, C(3)  α (l, m) is localized at the corner of the (l, m) space.In panel (b), the eigenmodes in the band with the smallest decay rate are classified based on whether the point (l * , m * ) that maximizes C(3)  α (l, m) is located in the bulk, edge, or corners of the (l, m) space.

TABLE I .
Quantum phase transition of the ground state (GS), dissipative phase transition of the steady state (SS), and deconfinement of incoherentons, which are characterized by the energy gap ∆ E , the Liouvillian gap ∆ L , and the QC gap ∆ QC , respectively.The closing of the QC gap is accompanied by a divergence of the confinement length of incoherentons.The critical slowing down of relaxation dynamics is characteristic of both quantum and dissipative phase transitions.The deconfinement of incoherentons signals a dynamical transition from incoherent exponential relaxation to coherent oscillatory one.