Activity-induced phase transition in a quantum many-body system
AActivity-induced phase transition in a quantum many-body system
Kyosuke Adachi,
1, 2
Kazuaki Takasan, and Kyogo Kawaguchi
1, 4, 5 Nonequilibrium Physics of Living Matter RIKEN Hakubi Research Team,RIKEN Center for Biosystems Dynamics Research, 2-2-3 Minatojima-minamimachi, Chuo-ku, Kobe 650-0047, Japan. RIKEN Interdisciplinary Theoretical and Mathematical Sciences Program, 2-1 Hirosawa, Wako 351-0198, Japan. Department of Physics, University of California, Berkeley, California 94720, USA RIKEN Cluster for Pioneering Research, 2-2-3 Minatojima-minamimachi, Chuo-ku, Kobe 650-0047, Japan. Universal Biology Institute, The University of Tokyo, Bunkyo-ku, Tokyo 113-0033, Japan. (Dated: August 4, 2020)A crowd of nonequilibrium entities can show phase transition behaviours that are prohibited in conventionalequilibrium setups. One of the simplest and well-studied phenomena observed in such active matter is themotility-induced phase separation (MIPS), where self-propelled particles spontaneously aggregate. An inter-esting question is whether a similar activity-driven phase transition occurs in a pure quantum system. Herewe introduce a non-Hermitian quantum many-body model that undergoes quantum phase transitions, whichincludes the analogue of the classical MIPS within its parameter space. The phase diagram, which can be exper-imentally tested in an open quantum system, indicates that the addition of spin-dependent asymmetric hoppingto a simple model of hard-core bosons is su ffi cient to induce flocking as well as phase separation. Moreover,we find that the quantum phase transitions in the model are equivalent to the transitions of dynamical pathsin classical kinetics upon the application of biasing fields. This example demonstrates that quantum systemscan indeed show activity-induced phase transitions, and sheds light on the rich connection between classicalnonequilibrium kinetics and non-Hermitian quantum physics. The collective dynamics of active or driven components canlead to phase transitions and pattern formations that are pro-hibited in equilibrium systems [1]. Recent works have shownthe properties of materials such as surface flow [2, 3], odd vis-cosity / elasticity [4, 5], and anomalous topological defect dy-namics [6–8] that can be realized by introducing activity intothe design. Additional to the application in biophysical stud-ies [9, 10], combining the understandings of active systemswith a broader range of models in condensed matter shouldbring progress not only in nonequilibrium physics but also inmaterial science [11].Although the scope of active matter has greatly widenedin the past years [12], its quantum analogue has so far notbeen explicitly proposed. Part of the reason is that the corre-sponding quantum system must be open (i.e., non-Hermitian),which is typically more di ffi cult to control both in experimentand theory compared with closed (i.e., Hermitian) systems.In recent years, however, there is rapidly growing interest instudying non-Hermitian quantum systems [13, 14] due to theadvances in atomic-molecular-optical experiments allowingcontrol over open quantum systems [15, 16] and theoreticaldevelopments exploring quantum critical phenomena [17–21]and topological phases [22–28]. We are therefore in positionto ask whether there exist new phases of matter induced byactivity (i.e., non-Hermitian terms that can be interpreted asself-driving) in quantum many-body systems, and if so, howthey can be realized in experiments.In this work, we start by considering the motility-inducedphase separation (MIPS) in a quantum system. MIPS is atypical example of activity-driven phenomena studied in clas-sical stochastic systems [29–37]. We here introduce a non-Hermitian model of hard-core bosons to demonstrate that aquantum system can undergo activity-induced phase transitionin the ground state that is equivalent to the classical MIPS.We further explore the richness of the quantum phase dia- gram, and discuss its relation to dynamical phase transitionsobserved in classical stochastic kinetics. Non-Hermitian hard-core bosons and classical MIPS
The model we study here involves quantum hard-core bosonswith “spin” s ( = ±
1) in a L x × L y rectangular lattice with peri-odic boundary condition (PBC): H = − J (cid:88) (cid:104) i , j (cid:105) , s ( a † i , s a j , s + a † j , s a i , s ) − ε J (cid:88) i , s s ( a † i , s a i − ˆ x , s − a † i , s a i + ˆ x , s ) − h (cid:88) i , s a † i , s a i , − s − U (cid:88) (cid:104) i , j (cid:105) ˆ n i ˆ n j − U (cid:88) i ˆ m i (ˆ n i + ˆ x − ˆ n i − ˆ x ) + (4 J + h ) N , (1)where ˆ n i , s : = a † i , s a i , s is the local density of particles with spin s , ˆ n i : = ˆ n i , + + ˆ n i , − , and ˆ m i : = ˆ n i , + − ˆ n i , − . ˆ x is the unit hor-izontal translation, and N is the fixed total number of parti-cles. The second term in Eq. (1) describes the spin-dependentasymmetric hopping ( J > − ≤ ε ≤ ε (cid:44)
0. The fourth and fifth terms representthe spin-independent and dependent nearest-neighbor interac-tions, respectively, with its general form discussed in Supple-mentary Section I.A. We take h > E ), and thecorresponding eigenstate can be taken to have all its elementsreal and positive (which we denote as | ψ (cid:105) ). This is due tothe Perron-Frobenius theorem, which can be applied since theo ff -diagonal elements of H in the Fock-space representationare all real and negative. In this work, we focus on how theground state | ψ (cid:105) (with ground state energy E ) changes ac- a r X i v : . [ c ond - m a t . s t a t - m ec h ] A ug a b cd + + –– FIG. 1.
The active lattice gas shows MIPS with the Ising criticality. a , In the ALG, particles with spin + (red) or − (blue) can hop tothe nearest-neighbor sites with spin-dependent rates and can flip their spins. b , ρ - ε phase diagram (heatmap of φ ) of the ALG with typicalconfigurations for h = . J in 100 ×
10 systems. In the top configuration, we show an example of L × L sub-boxes used in the finite-sizescaling analysis. c , The Binder cumulant Q when varying ε and L ( = , , , , ,
14) for ρ = .
55 and h = . J . d , Size scalingof ∂ Q /∂ε | ε = ε c , χ L | ε = ε c , and (cid:104) ( ∆ ρ L ) (cid:105) | ε = ε c for ρ = . , . , . h = . J . The black dashed lines correspond to the 2D Ising criticalexponents ( β = / γ = /
4, and ν = cording to the change of parameters in H .Within the parameter space of (1), there is a special sub-space defined by U = J and U = ε J , where the Hamilto-nian can be mapped to the transition rate matrix of an activelattice gas model (ALG) [38–40] (see Methods). The ALGhere is an N -particle model where the particles are exclusivelyhopping within the L x × L y rectangular lattice with the PBC(Fig. 1a). Each particle has a spin s ( = ±
1) as its internalvariable, which sets the rate of asymmetric hopping in the x -direction as (1 + ε s ) J and (1 − ε s ) J for the positive and neg-ative directions, respectively. The y -directional hopping rateis J , the spin flipping rate is h , and we define the density as ρ : = N / L x L y (0 < ρ < L x = L and L y = L ) to minimize the e ff ects of phaseboundaries, and to apply the sub-box method [41, 42]. InFig. 1b, we show the steady-state phase diagram of the ALGin the ρ - ε plane for h = . J obtained by Monte Carlo sim-ulations. We find that MIPS appears when the self-propulsionstrength ε is high, which is confirmed using φ : = √− min C ( x )as an order parameter, where C ( x ) is the density correlationfunction in the x -direction (see Methods). We expect the ex-tremum of the coexistence curve to be a critical point, ( ρ c , ε c ),in the thermodynamic limit ( L x , L y → ∞ with fixed ρ ). Todetermine the critical point, we apply a modified version of arecently proposed method [41] based on the finite-size scal-ing and consider expectation values of physical quantities in L × L sub-boxes (Fig. 1b, Methods). Since ρ c is expected to bearound 0 .
55 from the phase diagram (Fig. 1b), we examine ε c using the finite-size scaling for ρ = .
5, 0 .
55, and 0 .
6. Defin-ing ∆ ρ L : = ρ L − ρ , where ρ L is the density in the sub-box,we find that the Binder cumulant Q : = (cid:104) ( ∆ ρ L ) (cid:105) / (cid:104) ( ∆ ρ L ) (cid:105) approximately crosses at a single point with varying ε and thesystem size, which indicates ε c (cid:39) .
39 for ρ = .
55 (Fig. 1c).We note that the same kinetics in one-dimension (i.e., no hop-ping to the y -direction) does not induce macroscopic phaseseparation (Supplementary Section I.E).Previous works have studied similar models to the ALG, but with di ff erent setups [30, 42–46]. According to [45], ahydrodynamic description for the ALG can be derived exactlywhen taking L → ∞ with ε = O ( L − ) and h / J = O ( L − ),suggesting that the MIPS transition will show the mean-fieldcritical exponents in this limit. On the other hand, the MIPSin our setup of ALG may belong to the 2D Ising universalityclass, as we have assumed ε and h to be size-independent [42].To check this, we estimate the critical exponents using the fol-lowing expressions based on the scaling hypothesis [41, 42]: ∂ Q /∂ε | ε = ε c ∼ L /ν , χ L | ε = ε c ∼ L γ/ν , and (cid:104) ( ∆ ρ L ) (cid:105) | ε = ε c ∼ L − β/ν , where χ L [: = ( (cid:104) N L (cid:105) − (cid:104) N L (cid:105) ) / (cid:104) N L (cid:105) ] is the particle-number fluctuation in the sub-box (Fig. 1d). We find that thecritical exponents ( β , γ , and ν ) are consistent with the 2D Isinguniversality as in [42].In the quantum model, the classical condition ( U = J and U = ε J ) induces E =
0. The corresponding righteigenstate | ψ (cid:105) is equivalent to the steady-state distributionof the ALG, and the left eigenstate is the coherent state, (cid:104) ψ (cid:48) | = (cid:104) P | : = (cid:104) | exp( (cid:80) i , s a i , s ). For the case of ε = U = H is Hermitian and equivalent to the ferromagneticXXZ model with fixed magnetization [47, 48] (Supplemen-tary Section I.B), where a first-order transition between thesuperfluid and phase-separated states occurs at the Heisenbergpoint ( U = J ) [49]. The Heisenberg point is also special inthat the right and left ground states are both the coherent state. Quantum phase diagram and dynamical phase transition
To explore the phase diagram outside of the classical condi-tion, we conducted the di ff usion Monte Carlo (DMC) sim-ulation [50, 51]. In short, we run the Monte Carlo simula-tion for the ALG but with the additional steps of re-samplingthe states based on the calculated weights of the paths. Thisworks since the Hamiltonian can be divided into two parts H = − W − D , where W : = − H ( U = J , U = ε J ) cor-responds to the classical dynamics and D , being a diagonalmatrix, can be interpreted as the re-sampling weights (seeMethods). To discuss the phases, we focus on physical quan-tities which are functions of the configuration of the parti- ab -0.2 -0.1 0 0.1 0.200.10.20.30.50.4 -0.2 -0.1 0 0.1 0.2 -0.2 -0.1 0 0.1 0.200.050.10.15 -0.2-0.100.1-0.2 -0.1 0 0.1 0.2 -0.2 -0.1 0 0.1 0.2 -0.2 -0.1 0 0.1 0.200.10.20.30.50.4 00.050.10.15 -0.04-0.0200.02 FIG. 2.
Breaking the classical condition leads to quantum phasetransitions. a , b , U dependence with U = ε J ( a ) or U depen-dence with U = J ( b ) of the order parameters, φ PS and φ mPS ,and the ground-state energy, E , for ρ = . h = . J , and ε = , . , . × E , we also plotted the analytical results of (cid:104) H (cid:105) C for a dis-ordered state (dashed) and a PS state (dotted) ( a ) or for a mPS statewith one (dashed) or four (dotted) clusters ( b ) (see SupplementarySection I.C). cles, A ( { ˆ n i , s } ), and calculate (cid:104) A (cid:105) C : = (cid:104) P | A ( { ˆ n i , s } ) | ψ (cid:105) / (cid:104) P | ψ (cid:105) .Phase-separated (PS) states are characterized by φ PS : = ( L x L y ) − (cid:80) (cid:104) i , j (cid:105) ( (cid:104) ˆ n i − ρ )(ˆ n j − ρ ) (cid:105) C . For microphase-separated(mPS) states, in which the number of clusters is O ( L x ) (seeconfigurations in Fig. 2b), we utilize φ mPS as the order pa-rameter, which is the density of clusters with oppositely po-larized edges: φ mPS : = L x − (cid:80) L x i = (cid:104) ˆ m X i (ˆ n X i + − ˆ n X i − ) (cid:105) C , whereˆ n X i : = L y − (cid:80) L y j = ˆ n i ˆ x + j ˆ y and ˆ m X i : = L y − (cid:80) L y j = ˆ m i ˆ x + j ˆ y . In thelarge-size limit ( L x , L y → ∞ ), φ PS > φ mPS = φ PS > | φ mPS | > U from 2 J . As shown inFig. 2a, φ PS increases rapidly as a function of U at around U = J for a broad range of ε ( = , . , .
6) and U ( = ε J ),with the ground-state energy E having a kink at U = J .This line of phase separation transition extends from the first-order transition in the XXZ model ( ε =
0) [47, 49]. Second,for high enough ε ( = . φ PS and an increase in φ mPS occur simultaneously as U crosses ε J (Fig. 2b). As alsoindicated from the typical configuration and the kink in E (Fig. 2b), this is expected to be a discontinuous transition be-tween the PS and mPS states. For low ε ( = , . ε while fixing U = J and U =
0. Intriguingly, we find that a ferromagnetic order ap-pears without phase separation for high ε ( (cid:38) . M : = N − (cid:104) (cid:0)(cid:80) i ˆ m i (cid:1) (cid:105) C (Fig. 3a). Such polar order, which a b FIG. 3.
Spin-dependent asymmetric hopping drives phase sep-aration and polar order. a , ε dependence of the squared mag-netization M and typical configurations in 50 × b , ε dependence of M and φ PS in one-dimensional (1D) systems with L x = , , , , a and b , we set ρ = . h = . J , U = J , and U = should be accompanied by flow due to the asymmetric hop-ping, is reminiscent of the flocking of self-propelled particlesobserved, e.g., in the Vicsek model [52] and the active Isingmodel [53–55], although our model (1) does not include ex-plicit polar interactions.To investigate whether the polar order remains in larger sys-tems, we further performed simulations in one-dimensionalsystems. The size dependence of M and φ PS in one-dimensional systems (Fig. 3b) shows that the polar state isdestabilized and instead the PS state appears as the systemsize becomes larger. In addition, the discontinuous changesin M and φ PS indicate bi-stability of the polar and PS statesin finite systems. Similarly in large two-dimensional systems,the PS state can replace the polar state, as observed in the U -dependence of M and φ PS for the system with size 50 × U = ×
3. First, Fig. 4a is the U - ε phase diagram aroundthe classical line ( U = J and U = ε J ) indicated in red. Inaddition to the classical MIPS, the PS-mPS transition occursin crossing the classical line at high ε (see Fig. 2b). Next,Fig. 4b,c display the U - U phase diagrams around the clas-sical line. For low ε ( = .
2) (Fig. 4b), we find that the U -induced phase separation transition (Fig. 2a) occurs robustlyagainst U -perturbation from the classical line. In contrast, forhigh ε ( = .
6) (Fig. 4c), slight changes in U and U aroundthe classical line can lead to the mPS and polar states.The DMC simulation becomes less reliable for the param-eter regions far away from the classical line. Nevertheless,there are symmetries in this system that hint the positions ofthe phase boundaries in a wider parameter region (Fig. 4d).First, we have E ( − ε, − U ) = E ( ε, U ) which is due to H ( ε, U ) = ˆ U † H ( − ε, − U ) ˆ U , where ˆ U is the unitary opera-tor of spin reversal. We also have E ( ε, U ) = E ( − ε, U ),which follows from H ( ε, U ) † = H ( − ε, U ). Since the analyt-ical property of E indicates the positions of the phase bound-aries, we expect that the boundaries calculated in Fig. 4a-cmay have corresponding phase boundaries in ε < / or U < | ε | in crossing the dual classical line defined by U = J and U = − ε J , which is where E = | ψ (cid:105) = | P (cid:105) (Fig. 4d).The scheme of the DMC implies an interesting connectionbetween the quantum model and the classical kinetics. Forthe ALG with the transition rate matrix W , we denote the con-figuration of the particles at time t as C t = { n i , s ( t ) } , and itsstochastic trajectory as C t = C k ( t k ≤ t < t k + ) with t k beingthe time point of the k -th jump. For a path-dependent quantity¯ B τ : = (cid:82) τ dtB C t , C t + (cid:80) k B C k , C k + defined using some real matrix B that acts on the Fock space, we introduce λ W ( B ) : = lim τ →∞ τ ln (cid:104) exp( ¯ B τ ) (cid:105) W , (2)where the ensemble average (cid:104)· · ·(cid:105) W is taken over the tra-jectories in the ALG. λ W ( B ) is equivalent to the dominanteigenvalue of a biased transition rate matrix W B C , C (cid:48) : = (1 − δ C , C (cid:48) ) W C , C (cid:48) e B C , C(cid:48) + δ C , C (cid:48) ( W C , C (cid:48) + B C , C (cid:48) ) [56, 57]. Typical pathsthat appear in the biased dynamics can become dramaticallydi ff erent from the original dynamics, which is the hallmark ofdynamical phase transition that can be captured by the (non-)analytical behaviour of λ W ( B ) [57]. Biased kinetics and dy-namical phase transition have been studied with interests inexploring glassy systems [57–59] and in characterizing phasesin models of active matter [46, 60–62].The quantum Hamiltonian (1) can be interpreted as the tran-sition rate matrix with bias by writing H = − W B , where thebias is B = u D + u D with D C , C (cid:48) : = (cid:104)C| (cid:80) (cid:104) i , j (cid:105) ˆ n i ˆ n j |C (cid:48) (cid:105) , D C , C (cid:48) : = (cid:104)C| (cid:80) i ˆ m i (ˆ n i + ˆ x − ˆ n i − ˆ x ) |C (cid:48) (cid:105) , where |C(cid:105) is the Fock-space basis corresponding to the configuration C (: = { n i , s } ),and u : = U − J , u : = U − ε J . We then arrive at E ( ε, U , U ) = − λ W ( u D + u D ) , (3)which means that the quantum phase transitions, captured bythe property of E , are equivalent to the dynamical phase tran-sitions induced by the bias u D + u D . The bias here has aclear interpretation: increasing u and u favors larger φ PS and φ mPS , respectively.More generally, we can bring in an arbitrary transitionrate matrix W (cid:48) and define the bias B (cid:48) appropriately so that H = − W (cid:48) B (cid:48) . One interesting choice is W (cid:48) = ˜ W and B (cid:48) = ˜ B ,where ˜ W is the transition rate matrix that has the same diago-nal elements as H , and ˜ B is non-diagonal and non-Hermitian(Fig. 4d). From this point of view, the ε -dependent transi-tion toward the polar (flocking) phase (Fig. 3) can be un-derstood as the consequence of biasing the kinetics towardlarger ˜ B , which encourages more spin-dependent asymmet-ric hopping and therefore dissipation (see Supplementary Sec-tion I.F). Consistent with this, dynamical phase transition in-duced by biasing toward higher dissipation has been reportedin the studies of active Brownian particles [60–62]. Relevance to experiments
Lastly, we discuss the possibility of observing the activity-induced phase transitions in quantum experiments. In ultra-cold atom experiments, the Bose-Hubbard model has already a b c -0.25 0 0.25 0.5 0.75 1 1.2500.20.40.60.81 1.8 1.9 2 2.1 2.200.10.20.30.4 1.8 1.9 2 2.1 2.20.40.50.60.70.8
PSmPSDPD PSD P mPSPS e (ii) OBC, DPmPS PSSF+PP mPS PS DP PSmPSmPSP PS (i) PBC,(iii) PBC, (iv) OBC, d FIG. 4.
Ground-state phase diagrams of the quantum model. a , U - ε phase diagram for U = J around the classical line (red box),with PS ( φ PS > . φ mPS ≤ . φ mPS > . M > . b , c , U - U phasediagrams for ε = . .
6, respectively, around the cross sectionof the classical line (red box). d , Schematic of the U - ε plane at U = J . The Hamiltonian has two symmetries (see main text), meaningthat the points indicated by squares all have the same value of E .The classical line ( U = ε J ) and the dual classical line ( U = − ε J )have E =
0. The same Hamiltonian (e.g., black square) can bedescribed in multiple ways of classical stochastic dynamics (e.g., W and ˜ W ) with bias (e.g., green and magenta arrows). e , U - U phasediagrams for ε = . L x = φ PS > .
05 and φ mPS ≤ . φ mPS > . M > . φ SF > . (cid:104)· · ·(cid:105) C (i,ii) or (cid:104)· · ·(cid:105) Q (iii, iv), for the PBC (i, iii) or OBC (ii, iv). Superfluidstates cannot be identified in DMC calculations or by using (cid:104)· · ·(cid:105) C (see Supplementary Section II.C). In all figures, we set ρ = . h = . J . been well-simulated [16, 63]. Although the non-Hermitianterms are di ffi cult to implement in general, the asymmetrichopping terms of Eq. (1) can be realized by introducing acoherent coupling between the original square lattice and adissipative auxiliary lattice, as proposed in Ref. [25]. Suchopen system is described by the quantum Master equation,which can be shown to reduce to a non-Hermitian quantummechanical system by post-selected quantum trajectories (seeSupplementary Section II.A). To observe the ground state forthe non-Hermitian quantum system, we can first prepare theground state of the Hermitian model and adiabatically intro-duce ε . Since the ground state is unique within the parameterspace of Eq. (1), the system should stay in the ground state fora finite time [64] (see Supplementary Section II.B).The measurable quantity in quantum experiments is (cid:104)· · ·(cid:105) Q : = (cid:104) ψ | · · · | ψ (cid:105) rather than (cid:104)· · ·(cid:105) C [65–67]. Furthermore,typical cold atom experiments are in open boundary condition(OBC) [16], in which case the exact mapping to a classicalsystem does not exist (see Supplementary Section I.A). To ad-dress these points, we conducted exact diagonalization for asmall one-dimensional system to check how the order param-eters redefined using (cid:104)· · ·(cid:105) Q and the di ff erent boundary con-ditions will change the result of the phase diagram (Fig. 4e).We found that all of the phases exist in the various setups, withan additional polar-superfluid phase which can be captured byan o ff -diagonal order parameter, indicating that experimentswith small systems can already lead to interesting results (SeeSupplementary Section II.C).Here we have shown that a quantum many-body systemcan undergo activity-induced phase transition in the similarmanner as in the classical MIPS but with a richer phase dia-gram. The fact that the addition of a simple spin-dependenthopping can lead to non-trivial phases indicates the poten-tial of open quantum systems. Models with asymmetric hop-ping have been studied extensively in the recent context ofnon-Hermitian topological phases [25, 27, 28, 68–70]. It willbe interesting to consider the topological characterization ofphases in strongly interacting systems such as in the modelstudied here. Furthermore, the correspondence between thequantum Hamiltonian and the classical transition rate matrixwith bias indicates that dynamical phase transitions in gen-eral classical kinetics can in principle be probed by zero-temperature phase transitions in quantum experiments. Thisconnection is so far restricted to a stoquastic Hamiltonian (i.e., matrix with all its o ff -diagonal terms being real negative [71]);exploring other models of quantum active matter, especiallynon-stoquastic models that have no classical analogues, willbe an interesting next step. Acknowledgements
We thank Shin-ichi Sasa, MasatoItami, Hiroyoshi Nakano, Tomohiro Soejima, Masaya Naka-gawa, Yuto Ashida, and Hosho Katsura for the scientific dis-cussions. We are also thankful to Zongping Gong, TakahiroNemoto, Takaki Yamamoto, and Yoshihiro Michishita forhelpful comments. K.A. is supported by JSPS KAKENHIGrant No. JP20K14435, and the Interdisciplinary Theoreticaland Mathematical Sciences Program (iTHEMS) at RIKEN.K.T. thanks JSPS for support from Overseas Research Fel-lowship. K.K is supported by JSPS KAKENHI Grants No.JP18H04760, No. JP18K13515, No. JP19H05275, and No.JP19H05795. The numerical calculations have been per-formed on cluster computers at RIKEN iTHEMS.
Author contributions
K.A., K.T., and K.K. conceivedthe project. K.A. performed the Monte Carlo simulations andmade all the plots. K.T. and K.K. conducted the exact diag-onalization calculations. K.A., K.T., and K.K. discussed theresults and wrote the manuscript and the supplementary. [1] Marchetti, M. C. et al.
Hydrodynamics of soft active matter.
Rev. Mod. Phys. , 1143–1189 (2013).[2] Nash, L. M. et al. Topological mechanics of gyroscopic meta-materials.
Proc. Natl Acad. Sci. USA , 14495–14500 (2015).[3] Soni, V. et al.
The odd free surface flows of a colloidal chiralfluid.
Nat. Phys. , 1188–1194 (2019).[4] Banerjee, D., Souslov, A., Abanov, A. G. & Vitelli, V. Oddviscosity in chiral active fluids. Nat. Commun. , 1–12 (2017).[5] Scheibner, C. et al. Odd elasticity.
Nat, Phys. , 475–480(2020).[6] Sanchez, T., Chen, D. T., DeCamp, S. J., Heymann, M. &Dogic, Z. Spontaneous motion in hierarchically assembled ac-tive matter. Nature , 431–434 (2012).[7] Saw, T. B. et al.
Topological defects in epithelia govern celldeath and extrusion.
Nature , 212–216 (2017).[8] Kawaguchi, K., Kageyama, R. & Sano, M. Topological defectscontrol collective dynamics in neural progenitor cell cultures.
Nature , 327–331 (2017).[9] Prost, J., J¨ulicher, F. & Joanny, J.-F. Active gel physics.
Nat.Phys. , 111–117 (2015).[10] Trepat, X. & Sahai, E. Mesoscale physical principles of collec-tive cell organization. Nat. Phys. , 671–682 (2018).[11] Needleman, D. & Dogic, Z. Active matter at the interface be-tween materials science and cell biology. Nat. Rev. Mater. ,1–14 (2017).[12] Gompper, G. et al. The 2020 motile active matter roadmap.
J.Phys. Condens. Matter , 193001 (2020).[13] El-Ganainy, R. et al. Non-hermitian physics and PT symmetry.
Nat. Phys. , 11–19 (2018).[14] Ashida, Y., Gong, Z. & Ueda, M. Non-hermitian physics.Preprint at https: // arxiv.org / abs / Adv. Atom. Mol. Opt. Phys. , vol. 61 of
Advances In Atomic, Molecu-lar, and Optical Physics , 1 – 80 (Academic Press, 2012).[16] Sch¨afer, F., Fukuhara, T., Sugawa, S., Takasu, Y. & Takahashi,Y. Tools for quantum simulation with ultracold atoms in opticallattices.
Nat. Rev. Phys. (2020).[17] Ashida, Y., Furukawa, S. & Ueda, M. Parity-time-symmetricquantum critical phenomena.
Nat. Commun. , 15791 (2017).[18] Nakagawa, M., Kawakami, N. & Ueda, M. Non-hermitiankondo e ff ect in ultracold alkaline-earth atoms. Phys. Rev. Lett. , 203001 (2018).[19] Hamazaki, R., Kawabata, K. & Ueda, M. Non-hermitian many-body localization.
Phys. Rev. Lett. , 090603 (2019).[20] Kawabata, K. & Ryu, S. Nonunitary scalingtheory of non-hermitian localization. Preprint athttps: // arxiv.org / abs / arXiv e-prints (2020). arXiv:2002.08478.[22] Hu, Y. C. & Hughes, T. L. Absence of topological insulatorphases in non-hermitian pt -symmetric hamiltonians. Phys. Rev.B , 153101 (2011).[23] Esaki, K., Sato, M., Hasebe, K. & Kohmoto, M. Edge statesand topological phases in non-hermitian systems. Phys. Rev. B , 205128 (2011).[24] Shen, H., Zhen, B. & Fu, L. Topological band theory for non-hermitian hamiltonians. Phys. Rev. Lett. , 146402 (2018).[25] Gong, Z. et al.
Topological phases of non-hermitian systems.
Phys. Rev. X , 031079 (2018).[26] Yao, S. & Wang, Z. Edge states and topological invariants ofnon-hermitian systems. Phys. Rev. Lett. , 086803 (2018).[27] Kawabata, K., Shiozaki, K., Ueda, M. & Sato, M. Symmetryand topology in non-hermitian physics.
Phys. Rev. X , 041015(2019). [28] Okuma, N., Kawabata, K., Shiozaki, K. & Sato, M. Topolog-ical origin of non-hermitian skin e ff ects. Phys. Rev. Lett. ,086801 (2020).[29] Tailleur, J. & Cates, M. Statistical mechanics of interactingrun-and-tumble bacteria.
Phys, Rev. Lett. , 218103 (2008).[30] Thompson, A. G., Tailleur, J., Cates, M. E. & Blythe, R. A.Lattice Models of Nonequilibrium Bacterial Dynamics.
J. Stat.Mech.
P02029 (2011).[31] Fily, Y. & Marchetti, M. C. Athermal phase separation of self-propelled particles with no alignment.
Phys. Rev. Lett. ,235702 (2012).[32] Bialk, J., Lwen, H. & Speck, T. Microscopic theory for thephase separation of self-propelled repulsive disks.
EPL ,30008 (2013).[33] Redner, G. S., Hagan, M. F. & Baskaran, A. Structure anddynamics of a phase-separating active colloidal fluid.
Phys. Rev.Lett. , 055701 (2013).[34] Cates, M. E. & Tailleur, J. Motility-Induced Phase Separation.
Ann. Rev. Condens. Matter Phys. , 219–244 (2015).[35] Fodor, E. et al. How far from equilibrium is active matter?
Phys. Rev. Lett. , 038103 (2016).[36] Buttinoni, I. et al.
Dynamical Clustering and Phase Separationin Suspensions of Self-Propelled Colloidal Particles.
Phys. Rev.Lett. , 238301 (2013).[37] Liu, G. et al.
Self-driven phase transitions drive myxococcusxanthus fruiting body formation.
Phys. Rev. Lett. , 248102(2019).[38] Doi, M. Second quantization representation for classical many-particle system.
J. Phys. A , 1465–1477 (1976).[39] Doi, M. Stochastic theory of di ff usion-controlled reaction. J.Phys. A , 1479–1495 (1976).[40] Peliti, L. Path integral approach to birth-death processes on alattice. J. Phys. France , 1469–1483 (1985).[41] Siebert, J. T. et al. Critical behavior of active brownian parti-cles.
Phys. Rev. E , 030601 (2018).[42] Partridge, B. & Lee, C. F. Critical motility-induced phase sep-aration belongs to the ising universality class. Phys. Rev. Lett. , 068002 (2019).[43] Peruani, F., Klauss, T., Deutsch, A. & Voss-Boehme, A. Tra ffi cjams, gliders, and bands in the quest for collective motion ofself-propelled particles. Phys. Rev. Lett. , 128101 (2011).[44] Soto, R. & Golestanian, R. Run-and-tumble dynamics in acrowded environment: Persistent exclusion process for swim-mers.
Phys. Rev. E , 012706 (2014).[45] Kourbane-Houssene, M., Erignoux, C., Bodineau, T. &Tailleur, J. Exact hydrodynamic description of active latticegases. Phys. Rev. Lett. , 268003 (2018).[46] Whitelam, S., Klymko, K. & Mandal, D. Phase separation andlarge deviations of lattice active matter.
J. Chem. Phys. ,154902 (2018).[47] Matsubara, T. & Matsuda, H. A Lattice Model of Liquid He-lium, I.
Prog. Theor. Phys. , 569–582 (1956).[48] Whitlock, R. T. & Zilsel, P. R. Pseudospin model for hard-corebosons with attractive interaction. zero temperature. Phys. Rev. , 2409–2420 (1963).[49] Saryer, O. S. Two-dimensional quantum-spin-1 / Philos. Mag. , 1787–1824 (2019).[50] Giardin`a, C., Kurchan, J. & Peliti, L. Direct evaluation of large-deviation functions. Phys. Rev. Lett. , 120603 (2006).[51] Giardin`a, C., Kurchan, J., Lecomte, V. & Tailleur, J. Simulatingrare events in dynamical processes. J. Stat. Phys. , 787–811(2011).[52] Vicsek, T., Czir´ok, A., Ben-Jacob, E., Cohen, I. & Shochet, O. Novel type of phase transition in a system of self-drivenparticles.
Phys. Rev. Lett. , 1226–1229 (1995).[53] Solon, A. P. & Tailleur, J. Revisiting the flocking transitionusing active spins. Phys. Rev. Lett. , 078101 (2013).[54] Solon, A. P. & Tailleur, J. Flocking with discrete symme-try: The two-dimensional active ising model.
Phys. Rev. E ,042119 (2015).[55] Solon, A. P., Chat´e, H. & Tailleur, J. From phase to microphaseseparation in flocking models: The essential role of nonequilib-rium fluctuations. Phys. Rev. Lett. , 068101 (2015).[56] Lebowitz, J. L. & Spohn, H. A gallavotti–cohen-type symmetryin the large deviation functional for stochastic dynamics.
J. Stat.Phys. , 333–365 (1999).[57] Garrahan, J. P. et al. Dynamical first-order phase transition inkinetically constrained models of glasses.
Phys. Rev. Lett. ,195702 (2007).[58] Hedges, L. O., Jack, R. L., Garrahan, J. P. & Chandler, D. Dy-namic order-disorder in atomistic models of structural glass for-mers. Science , 1309–1313 (2009).[59] Merolle, M., Garrahan, J. P. & Chandler, D. Space–time ther-modynamics of the glass transition.
Proc. Natl Acad. Sci. USA , 10837–10840 (2005).[60] Nemoto, T., Fodor, E., Cates, M. E., Jack, R. L. & Tailleur, J.Optimizing active work: Dynamical phase transitions, collec-tive motion, and jamming.
Phys. Rev. E , 022605 (2019).[61] Tociu, L., Fodor, E., Nemoto, T. & Vaikuntanathan, S. Howdissipation constrains fluctuations in nonequilibrium liquids:Di ff usion, structure, and biased interactions. Phys. Rev. X ,041026 (2019).[62] Fodor, ´E., Nemoto, T. & Vaikuntanathan, S. Dissipation con-trols transport and phase transitions in active fluids: mobil-ity, di ff usion and biased ensembles. New J. Phys. , 013052(2020).[63] Gross, C. & Bloch, I. Quantum simulations with ultracoldatoms in optical lattices. Science , 995–1001 (2017).[64] Wang, H., Lang, L.-J. & Chong, Y. D. Non-hermitian dynam-ics of slowly varying hamiltonians.
Phys. Rev. A , 012119(2018).[65] Dalibard, J., Castin, Y. & Mølmer, K. Wave-function approachto dissipative processes in quantum optics. Phys. Rev. Lett. ,580–583 (1992).[66] Dum, R., Zoller, P. & Ritsch, H. Monte carlo simulation of theatomic master equation for spontaneous emission. Phys. Rev. A , 4879–4887 (1992).[67] Carmichael, H. An Open Systems Approach to Quantum Optics (Springer, Berlin, 1993).[68] Hatano, N. & Nelson, D. R. Localization transitions in non-hermitian quantum mechanics.
Phys. Rev. Lett. , 570–573(1996).[69] Hatano, N. & Nelson, D. R. Vortex pinning and non-hermitianquantum mechanics. Phys. Rev. B , 8651–8673 (1997).[70] Murugan, A. & Vaikuntanathan, S. Topologically protectedmodes in non-equilibrium stochastic systems. Nat. Commun. , 1–6 (2017).[71] Bravyi, S., Divincenzo, D. P., Oliveira, R. I. & Terhal,B. M. The complexity of stoquastic local hamiltonian prob-lems. Quantum Inf. Comput , 361385 (2008).[72] Martino, L., Elvira, V. & Louzada, F. E ff ective sample size forimportance sampling based on discrepancy measures. SignalProcess. , 386 – 401 (2017).[73] Cha, M.-C., Shin, J.-G. & Lee, J.-W. Momentum distribution ofthe hard-core extended boson hubbard model in one dimension.
Phys. Rev. B , 193104 (2009).[74] Nemoto, T., Bouchet, F., Jack, R. L. & Lecomte, V. Population- dynamics method with a multicanonical feedback control. Phys. Rev. E , 062123 (2016).[75] Lecomte, V., Appert-Rolland, C. & van Wijland, F. Chaoticproperties of systems with markov dynamics. Phys. Rev. Lett. , 010601 (2005).[76] Touchette, H. The large deviation approach to statistical me-chanics. Phys. Rep. , 1–69 (2009).[77] Kurchan, J. Fluctuation theorem for stochastic dynamics.
J.Phys. A , 3719 (1998).[78] Matthews, M. R. et al. Dynamical response of a bose-einsteincondensate to a discontinuous change in internal state.
Phys.Rev. Lett. , 243–247 (1998).[79] Bornheimer, U., Vasi´c, I. & Hofstetter, W. Phase transitionsof the coherently coupled two-component bose gas in a squareoptical lattice. Phys. Rev. A , 063623 (2017).[80] Landig, R. et al. Quantum phases from competing short- andlong-range interactions in an optical lattice.
Nature , 476–479 (2016).[81] Browaeys, A. & Lahaye, T. Many-body physics with individu-ally controlled rydberg atoms.
Nat. Phys. , 132–142 (2020).[82] Trefzger, C., Menotti, C., Capogrosso-Sansone, B. & Lewen-stein, M. Ultracold dipolar gases in optical lattices. J. Phys. B , 193001 (2011).[83] Yan, B. et al. Observation of dipolar spin-exchange interactionswith lattice-confined polar molecules.
Nature , 521–525(2013).[84] Zhao, H., Knolle, J. & Mintert, F. Engineered nearest-neighborinteractions with doubly modulated optical lattices.
Phys. Rev.A , 053610 (2019).[85] Syassen, N. et al.
Strong dissipation inhibits losses and inducescorrelations in cold molecular gases.
Science , 1329–1331(2008).[86] Tomita, T., Nakajima, S., Danshita, I., Takasu, Y. & Takahashi,Y. Observation of the mott insulator to superfluid crossover of adriven-dissipative bose-hubbard system.
Sci. Adv. , e1701513(2017).[87] Mark, M. J. et al. Preparation and spectroscopy of a metastablemott-insulator state with attractive interactions.
Phys. Rev. Lett. , 215302 (2012).[88] Barontini, G. et al.
Controlling the dynamics of an open many-body quantum system with localized dissipation.
Phys. Rev.Lett. , 035302 (2013).[89] Bloch, I., Dalibard, J. & Zwerger, W. Many-body physics withultracold gases.
Rev. Mod. Phys. , 885–964 (2008).[90] Gong, Z., Higashikawa, S. & Ueda, M. Zeno hall e ff ect. Phys.Rev. Lett. , 200401 (2017).[91] Gerbier, F. & Dalibard, J. Gauge fields for ultracold atoms inoptical superlattices.
New J. Phys. , 033007 (2010).[92] Reiter, F. & Sørensen, A. S. E ff ective operator formalism foropen quantum systems. Phys. Rev. A , 032111 (2012).[93] Yamamoto, K. et al. Theory of non-hermitian fermionic super-fluidity with a complex-valued interaction.
Phys. Rev. Lett. ,123601 (2019).[94] Bakr, W. S. et al.
Probing the superfluid–to–mott insulator tran-sition at the single-atom level.
Science , 547–550 (2010).[95] Sherson, J. F. et al.
Single-atom-resolved fluorescence imagingof an atomic mott insulator.
Nature , 68–72 (2010).
Methods
Mapping to the classical model.
We will show that theHamiltonian (1) is mapped to the active lattice gas model(ALG) under the classical condition ( U = J and U = ε J ).First, defining W : = − H ( U = J , U = ε J ), we can obtain W = ˆ P (cid:26) J (cid:88) (cid:104) i , j (cid:105) , s ( a † i , s a j , s + a † j , s a i , s ) + ε J (cid:88) i , s s ( a † i , s a i − ˆ x , s − a † i , s a i + ˆ x , s ) + h (cid:88) i , s a † i , s a i , − s − J (cid:88) (cid:104) i , j (cid:105) , s [ˆ n i , s (1 − ˆ n j ) + ˆ n j , s (1 − ˆ n i )] − ε J (cid:88) i , s s [ˆ n i , s (1 − ˆ n i + ˆ x ) − ˆ n i , s (1 − ˆ n i − ˆ x )] − h (cid:88) i , s ˆ n i , s (cid:27) ˆ P . (4)Here, we explicitly introduce the projection operator ˆ P to apartial Fock space where the total particle number is N withno multiple occupancy.Using W C , C (cid:48) : = (cid:104)C| W |C (cid:48) (cid:105) , where |C(cid:105) is the Fock-space basiscorresponding to a N -particle configuration C (: = { n i , s } ), wecan show that (i) (cid:80) C W C , C (cid:48) = W C , C (cid:48) > C (cid:44) C (cid:48) .Thus, we can think of W C , C (cid:48) as a transition rate matrix of aclassical Markov process which yields the master equation:d P ( C , t )d t = (cid:88) C (cid:48) W C , C (cid:48) P ( C (cid:48) , t ) . (5)where P ( C , t ) is the probability of configuration C at time t .The first three terms of (4) (non-diagonal elements of W C , C (cid:48) )represent a symmetric hopping rate, a spin-dependent changein the hopping rate, and a spin flipping rate; the last threeterms of (4) (diagonal elements of W C , C (cid:48) ) represent the cor-responding escape rates.Using a state vector | ψ ( t ) (cid:105) = (cid:80) C P ( C , t ) |C(cid:105) according tothe Doi-Peliti method [38–40], we can find that (5) is nothingbut the imaginary-time Schr¨odinger equation, d | ψ ( t ) (cid:105) / d t = − H ( U = J , U = ε J ) | ψ ( t ) (cid:105) . Thus, the steady state of theALG represented by | ψ ( t → ∞ ) (cid:105) is equivalent to the groundstate of the Hamiltonian, | ψ (cid:105) . Also, using the coherentstate (cid:104) P | = (cid:104) | exp( (cid:80) i , s a i , s ), we can express the expectationvalue of a classical physical quantity A ( { n i , s } ) as (cid:104) A (cid:105) ( t ) = (cid:80) C A ( C ) P ( C , t ) = (cid:104) P | A ( { ˆ n i , s } ) | ψ ( t ) (cid:105) . Especially for the steadystate ( t → ∞ ), we obtain (cid:104) A (cid:105) ( t → ∞ ) = (cid:104) P | A ( { ˆ n i , s } ) | ψ (cid:105) = (cid:104) A (cid:105) C . Indicator of phase separation in the ALG.
Using n X i : = L y − (cid:80) L y j = (cid:80) s n i ˆ x + j ˆ y , s , we define the density correlation func-tion in the x -direction as C ( x ) : = L x − (cid:80) i (cid:104) ( n X i − ρ )( n X i + x − ρ ) (cid:105) .Here, we introduce an indicator of phase separation, φ : = √− min x C ( x ). In a disordered state with no spatial correla-tion, we obtain C ( x ) = φ =
0. On the other hand,in a PS state, we can show φ = min { f h , f l } ( ρ h − ρ l ), wherethe densities (fractions) of the high-density and low-densityphases are ρ h and ρ l ( f h and f l ), respectively. Monte Carlo simulation of the ALG.
Setting a time step ∆ t [ = O ( N − )], we first randomly choose a particle from N parti-cles. Then, we flip the particle’s spin from s to − s with prob-ability hN ∆ t or move the particle to a neighboring empty sitewith probability s (1 + ε ) JN ∆ t , s (1 − ε ) JN ∆ t , or JN ∆ t de-pending on the hopping direction. We repeat this procedure M times until the total time T (: = M ∆ t ) is reached.In simulations, we took ∆ t = / N (4 J + h ). For Fig. 1b,we used M = N with 1000 samples. For Fig. 1c,d, wetypically took 5000 samples and used M = N for L = , , , M = × N for L =
12; and M = × N for L =
14. In all simulations, we set the disordered state with nospatial correlation as the initial state.
Sub-box method in the finite-size scaling analysis.
We ap-ply the sub-box method [41, 42], using six L × L sub-boxes(see green boxes in Fig. 1b). To sample configurations fromboth the high-density and low-density phases, we set the x -coordinate of the center of each sub-box as x c − L , x c , x c + L , x c + L , x c + L , and x c + L (mod L ), where x c is the cen-ter of mass of the system in the x -direction. In the finite-sizescaling analysis, we take the average of a quantity such as N L and ( ∆ ρ L ) over all the samples and sub-boxes. Di ff usion Monte Carlo simulation. For the quantum model[Eq. (1)], we first divide the Hamiltonian into two parts H = − W − D , where W is given by (4) and D is diago-nal in the Fock space. To numerically calculate the quantity (cid:104) A (cid:105) C = (cid:104) P | A ( { ˆ n i , s } ) | ψ (cid:105) / (cid:104) P | ψ (cid:105) for the ground state | ψ (cid:105) , wetransform (cid:104) A (cid:105) C as (cid:104) A (cid:105) C = lim T →∞ (cid:104) P | A ( { ˆ n i , s } )e ( W + D ) T | ψ ini (cid:105)(cid:104) P | e ( W + D ) T | ψ ini (cid:105) = lim T →∞ (cid:80) C , C A ( C ) (cid:104)C| e ( W + D ) T |C (cid:105) P ini ( C ) (cid:80) C , C (cid:104)C| e ( W + D ) T |C (cid:105) P ini ( C ) (6)where | ψ ini (cid:105) : = (cid:80) C P ini ( C ) |C(cid:105) with P ini ( C ) ≥ T → ∞ , we consider a finite butlarge enough T for the initial state to relax to the ground state.Splitting the total time T as T = M ∆ t with a time step ∆ t [ = O ( N − )] and writing C = C M for convenience, we candivide the time evolution into small steps: (cid:104)C| e ( W + D ) T |C (cid:105) = (cid:88) C , ··· , C M − M (cid:89) m = (cid:104)C m | e ( W + D ) ∆ t |C m − (cid:105)≈ (cid:88) C , ··· , C M − M (cid:89) m = ( δ C m , C m − + W C m , C m − ∆ t )(1 + D C m − ∆ t ) , (7)where D C : = (cid:104)C| D |C(cid:105) and the approximation in the third line iscorrect up to O ( ∆ t ). Since δ C m , C m − + W C m , C m − ∆ t is a stochasticmatrix for the ALG, we can approximately calculate Eq. (7)by assigning the weight (cid:81) Mm = (1 + D C m − ∆ t ) to the sampled path C → C → · · · → C M in the Monte Carlo (MC) simulationsof the ALG.To e ffi ciently sample the configurations that have highprobability weights but rarely appear in the MC simula-tion, we use the re-sampling technique [50, 51]. We con-sider a set of configurations, {C ( i ) m } N c i = , which evolve indepen-dently through the MC dynamics. Correspondingly, we intro-duce a set of cumulative weights, { w ( i ) m } N c i = , according to the paths {C ( i )0 → · · · → C ( i ) m } N c i = . Whenever the e ff ective sam-ple size [72], ( (cid:80) N c i = w ( i ) m ) / (cid:80) N c i = ( w ( i ) m ) , becomes smaller than0 . N c during the MC dynamics, we perform re-sampling ofconfigurations from the distribution of {C ( i ) m } N c i = weighted by { w ( i ) m } N c i = and then reset the weights as w ( i ) m = i . Us-ing the final-time configurations and weights, {C ( i ) M } N c i = and { w ( i ) M } N c i = , we estimate (cid:104) A (cid:105) C as (cid:104) A (cid:105) C ≈ (cid:80) N c i = w ( i ) M A ( C ( i ) M ) (cid:80) N c i = w ( i ) M . (8)In 2D simulations, we typically took ∆ t = / N (4 J + h )and used ( N c , M ) = (5 × , N ) for Fig. 2a; ( N c , M ) = (10 , N ) for Figs. 2b and 3a; and ( N c , M ) = (2 × , × N ) for Fig. 4a-c. In 1D simulations for Fig. 3b, we took ∆ t = / N (2 J + h ) and ( N c , M ) = (10 , N ). In all simula-tions, we set the disordered state with no spatial correlation asthe initial state, while we confirmed that there is no qualitativedependence on the initial state (Supplementary Section I.D.). Supplementary Information:Activity-induced phase transition in a quantum many-body system
Kyosuke Adachi, Kazuaki Takasan, and Kyogo Kawaguchi(Dated: August 4, 2020)
I. PROPERTIES OF THE MODEL AND DETAILS OF THE ANALYSISA. Generalized quantum model and classical condition
We consider a generalized version of the two-component hard-core boson model (1) in the main text: H gen = ˆ P − (cid:88) i (cid:88) l = x , y (cid:88) s , r = ± J ( l ) s , r a † i + r ˆ l , s a i , s − (cid:88) i (cid:88) a = , , , (cid:88) s , s (cid:48) = ± h a σ as , s (cid:48) a † i , s a i , s (cid:48) − (cid:88) i (cid:88) l = x , y (cid:88) s , r = ± U ( l ) s , r ˆ n i , s ˆ n i + r ˆ l ˆ P , (S1)where σ is the 2 × σ a ( a = , ,
3) are the Pauli matrices, and ˆ P is the projection to a partial Fock space wherethe total particle number is N with no multiple occupancy. We assume [ a i , s , a † j , s (cid:48) ] = [ a i , s , a j , s (cid:48) ] = [ a † i , s , a † j , s (cid:48) ] = i , s ) (cid:44) ( j , s (cid:48) ); { a i , s , a † i , s } = a i , s = ( a † i , s ) =
0. The first term of (S1) represents hopping, which is in general non-Hermitian and dependenton the spin or the hopping direction. The second and third terms represent the e ff ect of external fields and the generalizednearest-neighbor interactions, respectively. For J ( l ) s , r = (1 + sr εδ l , x ) J , h a = − (4 J + h ) δ a , + h δ a , , and U ( l ) s , r = U / + srU δ l , x , wecan reproduce the model (1) in the main text.Here, we take U ( l ) s , r = J ( l ) s , r , h = − (cid:80) l , s , r J ( l ) s , r / − h , h =
0, and h = − (cid:80) l , s , r sJ ( l ) s , r / J ( l ) s , r > h >
0, which isthe generalized classical condition (see the main text and Methods). Defining W : = − H gen under this classical condition, we canobtain W = ˆ P (cid:88) i (cid:88) l = x , y (cid:88) s , r = ± J ( l ) s , r (cid:20) a † i + r ˆ l , s a i , s − ˆ n i , s (cid:0) − ˆ n i + r ˆ l (cid:1)(cid:21) + (cid:88) i (cid:88) s = ± h (cid:16) a † i , s a i , − s − ˆ n i , s (cid:17) ˆ P . (S2)Defining W C , C (cid:48) : = (cid:104)C| W |C (cid:48) (cid:105) , where |C(cid:105) is the Fock-space basis, we can show that (i) (cid:80) C W C , C (cid:48) = W C , C (cid:48) > C (cid:44) C (cid:48) ,and thus we can interpret W C , C (cid:48) as a transition rate matrix of a classical Markov process. Under this interpretation, J ( l ) s , r is thehopping rate of a particle with spin s from a site i to the adjacent site i + r ˆ l , and h is the spin flipping rate.Lastly, we briefly discuss the quantum model (1) in the main text for the open boundary condition (OBC). OBC in a quantumsystem is when the hopping to the outside of the L x × L y region ( Ω ) is prohibited and there are no interactions between theparticles inside and the outside of Ω . This is di ff erent to the OBC in the classical system such as in ALG, meaning that thereis no classical line in the case of OBC. We conducted exact diagonalization calculations for a small one-dimensional quantumsystem to check the e ff ect this open boundary condition to the phase diagram (see II C and Fig. S8 for more details). On theother hand, we can think of a quantum system that corresponds to the ALG with OBC by setting U = J and U = ε J andadding a boundary term: W C , C (cid:48) = − (cid:104)C| H + H bd |C (cid:48) (cid:105) with H bd : = − J ˆ P [ (cid:80) i ∈ ∂ Ω \ ∂∂ Ω ˆ n i + (cid:80) i ∈ ∂∂ Ω ˆ n i + ε (cid:80) L y j = ( ˆ m L x ˆ x + j ˆ y − ˆ m x + j ˆ y )] ˆ P .Here we denoted the boundary points of Ω as ∂ Ω and the four corner points as ∂∂ Ω . B. Correspondence to the ferromagnetic XXZ model
We consider the case where ε = U =
0. Since there is no spin-dependence in this model, it is equivalent to thesingle-component hard-core boson model ( J > U > H HCB = − J (cid:88) (cid:104) i , j (cid:105) (cid:16) a † i a j + a † j a i (cid:17) − U (cid:88) (cid:104) i , j (cid:105) ˆ n i ˆ n j + const . (S3)Mapping the Fock bases to spin-1 / | n i = (cid:105) → | s zi = − / (cid:105) and | n i = (cid:105) → | s zi = + / (cid:105) , or equivalently, a i → ˆ S − i and a † i → ˆ S + i with ˆ S ± i : = ˆ S xi ± i ˆ S yi , we obtain H HCB → H XXZ = − (cid:88) (cid:104) i , j (cid:105) (cid:104) J (cid:16) ˆ S xi ˆ S xj + ˆ S yi ˆ S yj (cid:17) + U ˆ S zi ˆ S zj (cid:105) + const . (S4)For U > H XXZ represents the ferromagnetic XXZ model. Here, the total particle number N and the system size L x L y in thehard-core boson model are related to the total magnetization M z tot in the XXZ model as M z tot = N − L x L y /
2. In particular, when U = J , H XXZ is nothing but the ferromagnetic Heisenberg Hamiltonian [47, 48, 73].0 a + –+++ ––– + –+++ ––– + –+++ ––– + –+++ –––+ –+++ – –– –++ +– ––– ++++ +–+ +– ––+ –––––++ – ––++ + – +++ +–+ +– ––+ ––– +– ++ – ++– b c Fully phase-separated Fully polarized microphase-separatedFully polarized microphase-separated
FIG. S1. Schematic figures of representative states. a , In the fPS state, a single cluster with random spins is formed and its circumference isminimized. b , In the fpmPS state, there are N cl clusters with oppositely polarized edges. c , For large enough U ( (cid:29) J , h , U ), the fpmPS stateis stable with a macroscopic number of clusters, N cl = ρ L x / a b FIG. S2. Illustration of the convergence of the DMC simulations. a , U -dependence of the order parameters, φ PS and φ mPS , obtained with thedisordered (solid line with circles) or the PS (dashed line with triangles) initial state for ε = , . , . b , Time evolution of φ PS and E insimulations, obtained with the disordered initial state for ε = . U / J = . , , .
2. In both a and b , we considered 50 × ρ = . h = . J , and U = ε J . Simulation parameters are ∆ t = / N (4 J + h ), N c = × , and M = N as used in Fig. 2a of themain text (see Methods). Note that, since we set (cid:126) =
1, time and inverse of energy have the same dimension.
C. Energy of di ff erent states For an arbitrary state | ψ (cid:105) = (cid:80) C P ( C ) |C(cid:105) , where |C(cid:105) is the Fock-space basis, we can calculate (cid:104) H (cid:105) C as (cid:104) H (cid:105) C = (cid:88) C (2 J − U ) (cid:88) (cid:104) i , j (cid:105) n i n j + ( ε J − U ) (cid:88) i m i ( n i + ˆ x − n i − ˆ x ) P ( C ) (cid:44) (cid:88) C P ( C ) . (S5)Here, n i and m i are the local density and magnetization for the configuration C , respectively.In Fig. 2a,b in the main text, we plotted (cid:104) H (cid:105) C calculated for the disordered state with no spatial correlation, the fully phase-separated (fPS) state, and the fully polarized microphase-separated (fpmPS) state. First, the disordered state with no spatialcorrelation is defined as | ψ (cid:105) = ( (cid:80) i , s a † i , s ) N | (cid:105) , and the corresponding energy is (cid:104) H (cid:105) C = J − U ) ρ L x L y by neglecting o ( L x L y ),which we plot in Fig. 2a (dashed line). Second, we define a fPS state as | ψ (cid:105) = (cid:81) i ∈ Ω ( a † i , + + a † i , − ) | (cid:105) , where Ω is an area containing N sites and minimizing the circumference (Fig. S1a). The corresponding energy is (cid:104) H (cid:105) C = J − U ) ρ L x L y by neglecting o ( L x L y ), which we plot in Fig. 2a (dotted line). Lastly, we define a fpmPS state with N cl clusters, assuming commensurability,as | ψ (cid:105) = (cid:81) N cl n = [ (cid:81) i ∈ Ω n \ ( ∂ Ω L n ∪ ∂ Ω R n ) ( a † i , + + a † i , − ) (cid:81) i ∈ ∂ Ω L n a † i , + (cid:81) i ∈ ∂ Ω R n a † i , − ] | (cid:105) , where Ω n is the n -th rectangular area and ∂ Ω L(R) n is its left(right) boundary (Fig. S1b). The corresponding energy is (cid:104) H (cid:105) C = (2 J − U )(2 ρ L x L y − N cl L y ) + ε J − U ) N cl L y , which we plotwith N cl = N cl = U (cid:29) J , h , U ( > N cl = ρ L x / h / U , J / U , and U / U . D. Convergence of simulations and asymmetric-hopping-induced phase separation
In the di ff usion Monte Carlo (DMC) simulations, we checked the convergence to the steady-state by examining the initial-state dependence of the results and the relaxation of the order parameters and the ground-state energy. As an illustration, weshow the U -dependence of φ PS and φ mPS obtained with the fPS initial state, compared with that obtained with the disorderedinitial state (Fig. S2a and also see Fig. 2a in the main text). Apart from statistical errors, we do not see di ff erences due to initialconditions for the case of system size 30 ×
3, but there is discrepancy in the case of 50 × ffi cient for the large system size simulation [74]. Further, we show an example of the time1 FIG. S3. U -dependence of φ PS , φ mPS , and M for ρ = . ε = . h = . J , and U = J in 30 × × ∆ t = / N (4 J + h ), N c = , and M = N . a b FIG. S4. Macroscopic MIPS is not stable in the 1D ALG. a , ρ -dependence of φ for di ff erent values of ε (in increments of 0 .
1) in the 1D ALGwith L x = , , b , L x -dependence of φ for ρ = . ε =
1, which indicates φ ∼ L x − . for large L x . In both a and b , we used h = . J . Also, we used ∆ t = / N (2 J + h ) and took 10 samples with M = N for L x = , , × samples with M = × N for L x = × samples with M = × N for L x = dependence of φ PS and E evolving from the disordered initial state in the DMC simulations (Fig. S2b), which indicates that thesteady-state is achieved in the final state. Note that, for U = J and U = ε J (classical condition), E is trivially zero accordingto the probability conservation.We show the U -dependence of the order parameters for 30 × × M is destabilized and instead the PS state with finite φ PS dominates broader parameter regions as the system becomeslarger, though the dependence on the initial state remains around the phase boundary in the 50 × U = E. One-dimensional model
We show the ρ and ε -dependence (Fig. S4a) and the size-dependence (Fig. S4b) of φ [: = − min C ( r )] for the 1D counterpartof the ALG, suggesting that the macroscopic motility-induced phase separation (MIPS) is not stable in the thermodynamiclimit ( φ → L x → ∞ ). This result is consistent with preceding studies of similar 1D models [29, 30, 34, 44], where themacroscopic MIPS does not occur due to the spontaneous formation of domain boundaries.For the quantum model, Figs. S5 and S6 show the 1D counterparts of Figs. 2 and 4a-c in the main text, respectively. We cansee that the discontinuous transition occurs in crossing the classical line (Fig. S5) as observed in 2D systems, and the topologyof the phase diagrams (Fig. S6) is also similar. Note that in 1D systems with finite ε or U , the mPS order parameter φ mPS [: = L x − (cid:80) L x i = (cid:104) ˆ m i (ˆ n i + − ˆ n i − ) (cid:105) C ] is generically non-zero even for L x → ∞ , and consequently the disordered and mPS states areindistinguishable from the symmetry perspective.2 ab FIG. S5. Quantum phase transitions in 1D systems. a , b , U -dependence with U = ε J ( a ) or U -dependence with U = J ( b ) of φ PS , φ mPS ,and E for ρ = . h = . J , and ε = , . , . L x = E , we also plotted (cid:104) H (cid:105) C for the disordered state with no spatial correlation (dashed) and the fPS state (dotted)( a ) or for the fpmPS states with N cl = N cl = N cl =
13 (dash-dotted) ( b ). We set ∆ t = / N (2 J + h ) and took N c = and M = × N for a , while N c = and M = N for b . a b -0.25 0 0.25 0.5 0.75 1 1.2500.20.40.60.81 1.8 1.9 2 2.1 2.2 1.8 1.9 2 2.1 2.2 1.8 1.9 2 2.1 2.2-0.2-0.100.10.2 00.10.20.30.4 0.40.50.60.70.8 D DD mPSPS PS mPSP PmPSPS PS
FIG. S6. Ground-state phase diagrams of the 1D quantum model. a , U - ε phase diagram for U = J around the classical line (red box) withPS ( φ PS > . φ mPS ≤ . φ mPS > . M > . b , U - U phase diagrams for ε = , . , . ρ = . h = . J . Simulation parameters are ∆ t = / N (2 J + h ), N c = , and M = N . F. Classical system with dynamic bias
Let us introduce a real matrix B that acts on the Fock space (i.e., its indices are the configurations C ). We define the path-dependent observable between time t ∈ [0 , τ ] constructed from B as¯ B τ : = (cid:90) τ dtB C t , C t + (cid:88) k B C k , C k + , (S6)where the path is denoted as C t = C k ( t k ≤ t < t k + ) with t k being the time point of the k -th jump. We introduce the followingquantity: λ W ( B ) : = lim τ →∞ τ ln (cid:104) exp( ¯ B τ ) (cid:105) W , (S7)3defined for a matrix B , where (cid:104)· · ·(cid:105) W indicates the ensemble average taken with the path probability generated by W . It is knownthat λ W ( B ) is the dominant eigenvalue of W B , which is defined as [56, 75] W B C , C (cid:48) = W C , C (cid:48) e B C , C(cid:48) ( C (cid:44) C (cid:48) ) W C , C + B C , C ( C = C (cid:48) ) . (S8)We may think of a biased Markovian dynamics defined byd P ( C , t )d t = (cid:88) C (cid:48) W B C , C (cid:48) P ( C (cid:48) , t ) , (S9)which does not conserve the sum of the probability unless λ W ( B ) = (cid:88) C P ( C , t ) ∝ e λ W ( B ) t ( t → ∞ ) . (S10)For the special case of B = uA , where u is a real parameter, Eq. (S7) is the scaled cumulant generating function (SCGF) of A : λ WA ( u ) : = lim τ →∞ τ ln (cid:104) exp( u ¯ A τ ) (cid:105) W (S11)An interesting property of the SCGF is that it is related to the rate function, I WA ( a ) : = − lim τ →∞ τ ln Prob W ( ¯ A τ (cid:39) a τ ) , (S12)which characterizes the rare events of the paths taking atypical values of ¯ A τ in the stochastic dynamics following W . The relationis called the G¨artner-Ellis theorem [76]: I WA ( a ) = sup u ∈ R [ ua − λ WA ( u )] . (S13)The biased rate, such as W uA , can therefore be thought of as conditioning on the paths so that the average ¯ A τ will become acertain value a τ : d λ WA ( u ) du = lim τ →∞ (cid:104) ¯ A τ (cid:105) W uA τ = a . (S14)This equation determines u (i.e., the strength of the bias) required to achieve lim τ →∞ (cid:104) ¯ A τ (cid:105) W uA /τ = a . Here, (cid:104)· · ·(cid:105) W uA = (cid:104)· · · exp( u ¯ A τ ) (cid:105) W / (cid:104) exp( u ¯ A τ ) (cid:105) W denotes the ensemble in the biased kinetics.Let us write H = H ( J , ε, U , U , h ) to describe the parameter dependence of the Hamiltonian [(1) in main text]. Then wecan take ˜ W = − H ( J , ε , U , U , h ) with J = U / ε = U / U , which is a transition rate matrix that is distinct from W [ = − H ( J , ε, J , ε J , h )] (see Fig. 4d in the main text). According to Eq. (S8), we should take the bias ˜ B as˜ B C , C (cid:48) = | V C , C (cid:48) | ln JJ + ln 1 + ε V C , C (cid:48) + ε V C , C (cid:48) (S15)in order to obtain − H ( J , ε, U , U , h ) = ˜ W ˜ B . Here, V is a skew-Hermitian matrix given by V C , C (cid:48) = (cid:88) i , s s (cid:104)C| ( a † i , s a i − ˆ x , s − a † i , s a i + ˆ x , s ) |C (cid:48) (cid:105) (S16)Biasing the system toward larger ε induces the absolute value of ˜ B to become larger. To see its relation to the dissipation, let usintroduce the entropy production σ ( W (cid:48) ) defined for general transition rate matrices [56]: σ C , C (cid:48) ( W (cid:48) ) = ln W (cid:48)C , C (cid:48) W (cid:48)C (cid:48) , C . (S17)Then we find, ˜ B − ˜ B † = σ ( ˜ W ˜ B ) − σ ( ˜ W ) , (S18)which indicates that the di ff erence of entropy production in the biased and unbiased kinetics is exactly the non-Hermiticity ofthe bias. We also note that there is a fluctuation theorem-like relation [56, 77]: λ ˜ W ( ˜ B ) = λ ˜ W ( ˜ B † − σ ( ˜ W )) , (S19)which follows from ( ˜ W ˜ B ) † = ˜ W − σ ( ˜ W ) + ˜ B † . This symmetry, which is nothing but the E ( − ε, U ) = E ( ε, U ) symmetry noted inthe main text, is depicted as magenta arrows in Fig. 4d.4 II. EXPERIMENTAL SETUPA. Implementation of the quantum model
We describe how to implement our quantum model [Eq. (1) in the main text] in ultracold atomic systems. The basic model isa two-component Bose-Hubbard model on a square lattice, which is realized with bosonic ultracold atoms (e.g., Rb) in opticallattices [16, 63]. The two components correspond to the two internal states of atoms, and we call them spins relating to the modelsetup. The Hubbard interaction is controllable via the Feschbach resonance and we assume strong repulsive interaction to reachthe hard-core limit. Other ingredients to be implemented are following: (i) transverse magnetic field, (ii) nearest-neighbourinteraction, (iii) spin-dependent asymmetric hopping.For (i), the transverse magnetic field is implemented by a coherent coupling between the two internal states. Such coherentcoupling is well-studied and widely used in two-component bosonic atomic gases [78, 79]. The ingredient (ii) is a little morechallenging. The fourth (fifth) term of Eq. (1) in the main text are written as the (spin-dependent) density-density nearest-neighbour interaction. In optical lattice systems, on-site interaction is introduced naturally because most scattering occurslocally. In contrast, long-range interaction (including the nearest-neighbor interaction in a lattice system) is di ffi cult to berealized in general. However, there have been various proposals to overcome this di ffi culty, such as the use of optical cavity [80],Rydberg states [81], dipolar interaction [82, 83], or Floquet engineering [84]. Another idea to realize attractive interaction underhard-core condition using dissipation is to consider an attractive Bose-Hubbard model with strong two-body loss, which can beintroduced both intrinsically [83, 85] and artificially, e.g., via photoassociation [86]. In the Zeno limit, the quantum Zeno e ff ectsuppresses the double occupancy and the hard-core condition should be e ff ectively satisfied. A similar phenomenon in the caseof three-body loss has been observed in experiment [87].The ingredient (iii) is qualitatively di ff erent from the others since this is non-Hermitian. However, as detailed in Ref. [25], it ispossible to implement the non-Hermitian e ff ect by using the dissipative optical lattice setup. Here we briefly outline the method,with remarks on how to consider spin-dependent non-Hermitian hopping.The typical description of dissipative cold atomic systems is by a Lindblad-type quantum master equation which is d ρ ( t ) dt = − i [ H , ρ ( t )] + (cid:88) k D [ L k ] ρ ( t ) , (S20)where D [ L ] ρ ( t ) = L ρ ( t ) L † − { L † L , ρ ( t ) } /
2. Under postselection, where we only leave the experimental data that the loss processdid not happen, this master equation is simplified as d ρ ( t ) dt = − i ( H e ff ρ ( t ) − ρ ( t ) H † e ff ), where the e ff ective Hamiltonian is defined as H e ff = H − i (cid:88) k L † k L k . (S21)Intuitively, this e ff ective Hamiltonian contains the back action of the dissipation since we only see the Hilbert space where theloss event does not happen, which becomes the the origin of non-Herminicity. This non-Hermitian Hamiltonian has been well-examined in the studies of quantum trajectory method, which is an e ffi cient approach to simulate the dynamics of open quantumsystems [65–67]. The form of the e ff ective Hamiltonian (S21) suggests that we can engineer the non-Hermitian Hamiltonian bychoosing the adequate dissipators { L k } .To engineer the spin-dependent asymmetric hopping, we can use the following dissipators: L k = √ ε J ( a j , s + isa j + ˆ x , s ) , (S22)where k denotes a pair of the indices ( j , s ). Assuming these dissipators, we obtain − i (cid:88) k L † k L k = − ε J (cid:88) j , s s ( a † j , s a j − ˆ x , s − a † j , s a j + ˆ x , s ) − i ε JN (S23)where N denotes the total number operator N = (cid:80) j , s ˆ n j , s . The first term is nothing but the spin-dependent asymmetric hopping.We consider a subspace with a fixed total particle number and thus the second term only gives a constant energy shift.Finally, we explain how to implement these dissipators in experiments. It is not straightforward to realize a nonlocal one-bodyloss like Eq. (S22) since the usual loss occurs locally [88]. However, an experimental scheme for the spin- independent casehas been already proposed in Ref. [25] and the spin- dependent case [Eq. (S22)] is also possible by combining this scheme withspin-dependent lattice and coupling [89].First, we consider the spin-independent case. Following Ref. [90], the basic idea is introduction of a nonlocal coherentcoupling to an auxiliary dissipative lattice, schematically shown in Fig. S7. The coherent coupling to the dissipative latticedisplaced by a half of the lattice constant naturally induces the hopping from the j - and ( j + ˆ x )-sites, which becomes the origin5 FIG. S7. Experimental implementation of the asymmetric hopping in cold atomic systems. The original optical lattice (dark blue), thedissipative optical lattice (light blue), the coherent coupling between two lattices (red and green arrows) and the running wave in the x -direction(orange line) are introduced and bosonic atoms (blue balls) are loaded in the optical lattices. of the nonlocal loss. This setting of half-lattice is possible by using the internal atomic states with opposite Stark shifts. Forinstance, S and P of Yb atoms have the opposite Stark shift [91] (for detail, see Appendix F of [25]). Writing down theMaster equation within the tight-binding approximation and eliminating the fast decay mode [92], we can obtain the nonlocalone-body loss. In order to obtain Eq. (S22), we need to additionally introduce a running wave laser whose wavelength is equalto that of the lattice constant. This running wave provides the phase di ff erence between the couplings at the j - and ( j + ˆ x )-sites.Taking this e ff ect into account, we obtain Eq. (S22) except for the spin-dependency.To reach Eq. (S22), we need to change the asymmetric direction depending on the spins. The asymmetric direction is deter-mined by the propagation direction of the running wave laser [25]. Therefore, the oppositely-directed running wave lasers coupleto the each spin component respectively as in the case of spin-selective optical lattice [89], which will lead to spin-dependentasymmetric hopping. We remark that we need to select atom species having adequate internal states and specific laser settingsfor the real implementation of the oppositely directed running wave laser. B. Preparation of the ground state
In our study, we focus on the ground state | ψ (cid:105) , which is defined here as an eigenstate having the smallest real part of theenergy eigenvalue. In Hermitian systems, the ground state is e ff ectively realized at low temperature, but it is not clear how torealize the ground state in non-Hermitian systems. This is because the non-Hermitian system is expected to reach a state withthe largest imaginary part of the energy eigenvalues in the long-time limit. Thus, we need a some special protocol to realize theground state in non-Hermitian systems.Here, we propose an approach of adiabatic preparation. First, we prepare the Hermitian system ( ε =
0) and realize the lowtemperature state via thermalization in a closed quantum system. Then, we introduce dissipation adiabatically; i.e., turn on theasymmetric hopping term very slowly. Thanks to the Perron-Frobenius theorem, the uniqueness and the realness of the groundstate energy is guaranteed, and thus the energy gap should remain open through this process at least for a finite time (in a finitesystem). We note that there are various types of energy gaps in non-Hermitian systems [27], and here we have in mind thegap ∆ = | E − E | between the two energy eigenvalues E and E . It is not clear that the adiabatic theorem is valid even innon-Hermitian systems, but it has been shown by a numerical calculation that, when there is a finite gap in the sense of ∆ , thestate keeps sitting on the same state for a finite time under varying the parameters slowly [64]. Therefore, the ground state isrealized for a desired value of ε by this protocol. We remark that a similar approach has also been used in previous works onnon-Hermitian quantum many-body systems [17, 93]. C. Measurable quantities and its relation to the results from the Monte Carlo simulation
The most promising method to detect activity-induced phase transition such as MIPS is a quantum gas microscope (QGM) [16,63, 94, 95]. This enables us to measure the observable in a spatially-resolved way. Using the observed quantities, we cancalculate the order parameters of each phase transition. For instance, the indicator of MIPS, φ = √− min C ( r ) (for the definition,see the main text), is calculated from the local density data. The technique of QGM is growing rapidly and the measurements inthe Bose-Hubbard systems have already been conducted [94, 95].6In real experiments in open quantum systems, the measurable expectation value for the ground state | ψ (cid:105) is (cid:104) A (cid:105) Q = (cid:104) ψ | A ( { ˆ n i , s } ) | ψ (cid:105) [65–67] for a physical quantity A ( { ˆ n i , s } ), which depends on the particle configuration. On the other hand,the corresponding quantity calculated in the DMC simulation is (cid:104) A (cid:105) C = (cid:104) P | A ( { ˆ n i , s } ) | ψ (cid:105) / (cid:104) P | ψ (cid:105) with (cid:104) P | = (cid:104) | exp( (cid:80) i , s a i , s ).Another choice, which has been utilized in theoretical studies (e.g. Ref. [93]), is (cid:104) A (cid:105) LR = (cid:104) ψ (cid:48) | A ( { ˆ n i , s } ) | ψ (cid:105) / (cid:104) ψ (cid:48) | ψ (cid:105) with (cid:104) ψ (cid:48) | being the left ground state.For each of (cid:104)· · ·(cid:105) Q , (cid:104)· · ·(cid:105) C , and (cid:104)· · ·(cid:105) LR , we define the order parameters, φ PS , φ mPS , and M (see the main text). For (cid:104)· · ·(cid:105) Q(LR) ,we also define the order parameter for the superfluid (SF) state, which is characterized by the o ff -diagonal long-range order,as φ SF = L x − (cid:80) s (cid:80) | i − j | = L x / (cid:104) a † i , s a j , s (cid:105) Q(LR) . Note that φ SF for (cid:104)· · ·(cid:105) C is meaningless since the SF order and the density order areequivalent ( (cid:104) a † i , s a j , s (cid:105) C = (cid:104) ˆ n i , s ˆ n j , s (cid:105) C ). In addition to the periodic boundary condition (PBC), we here consider the open boundarycondition (OBC), which is relevant to typical cold atom experiments [16, 63].To clarify how the phase diagrams depend on the definition of order parameters and the boundary condition, we calculatedthe order parameters using exact diagonalization (ED) in small 1D systems (Fig. S8). First, we find that all the states predictedusing (cid:104)· · ·(cid:105) C with the PBC [(i) in Fig. S8a-c] appear, regardless of the definition of order parameters or the boundary condition.Thus, the DMC simulation, which is applicable to larger systems as demonstrated in the main text, is useful in qualitativelypredicting the phase diagram (apart from the SF order) in the experimentally relevant case, where we use (cid:104)· · ·(cid:105) Q with the OBC.Next, focusing on the cases with the PBC, we see that the SF state appears for ε = ε [(iii, v) in Fig. S8b,c]. Lastly, since the OBC prevents the particles from flowing, the polar order is suppressed [(ii, iv, vi) inFig. S8a,b] unless ε is large enough [(ii, iv, vi) in Fig. S8c].7 a b (i) PBC,(iii) PBC, (ii) OBC,(iv) OBC,(v) PBC, (vi) OBC, (i) PBC,(iii) PBC, (ii) OBC,(iv) OBC,(v) PBC, (vi) OBC, c (i) PBC,(iii) PBC, (ii) OBC,(iv) OBC,(v) PBC, (vi) OBC, D SF+PPmPS PS SF+PSFSF PSD D PS PSDmPS mPS mPSmPSmPSmPSmPS PSPSPSPSPSPSPSPS DDDD PP DD DDSF+PSF+P PPPP P DmPSmPSD PSPSPS PSPSPSDSF+PmPSmPSmPSmPS
FIG. S8. a , b , c , U - U phase diagrams in small 1D systems ( L x =
12) for (a) ε =
0, (b) 0 .
2, and (c) 0 .
6, respectively, with PS ( φ PS > .