Z2 parton phases in the mixed-dimensional $t-J_z$ model

We study the interplay of spin- and charge degrees of freedom in a doped Ising antiferromagnet, where the motion of charges is restricted to one dimension. The phase diagram of this mixed-dimensional $t - J_z$ model can be understood in terms of spin-less chargons coupled to a Z2 lattice gauge field. The antiferromagnetic couplings give rise to interactions between Z2 electric field lines which, in turn, lead to a robust stripe phase at low temperatures. At higher temperatures, a confined meson-gas phase is found for low doping whereas at higher doping values, a robust deconfined chargon-gas phase is seen which features hidden antiferromagnetic order. We confirm these phases in quantum Monte Carlo simulations. Our model can be implemented and its phases detected with existing technology in ultracold atom experiments. The critical temperature for stripe formation with a sufficiently high hole concentration is around the spin-exchange energy $J_z$, i.e., well within reach of current experiments.

Introduction.-Ultracold atoms in optical lattices provide an excellent platform to perform analogue quantum simulations: they can mimic the behavior of tunable model Hamiltonians that are difficult or impossible to solve with current numerics. Since the advent of quantum simulators, an application to the 2D Fermi-Hubbard model has been a central goal: This model is believed to describe some of the most essential but theoretically poorly understood properties of strongly correlated electrons in the context of high-temperature superconductors. In the past years, significant steps have been taken towards simulating the Hubbard model, including the observation of long-range [1] and canted [2] antiferromagnetism (AFM), bad metallic- [3] and spin- [4] transport, magnetic polarons [5,6], string patterns [7,8], and in 1D spin-charge separation [9,10] and incommensurate magnetism [11]. Nevertheless, the critical temperatures of the expected ordered phase (stripes [12], superconductivity [13]) are too low and have not yet been reached in ultracold fermion experiments.
In this Letter we make use of the versatility of ultracold atoms to study a closely related cousin of the 2D Hubbard model. Its two main advantages are: (i) significantly enhanced critical temperatures for the formation of stripe order amenable to quantum simulation; (ii) thorough theoretical understanding and numerical control of the underlying physics. Both (i) and (ii) provide a promising starting point, in experiment and theory, for a systematic exploration of the 2D Hubbard model.
Specifically, we consider a t − J z model with mixed dimensionality [14] as elucidated in Fig. 1 (a), The dopants (holes) are free to move only along the xdirection, with tunneling rate t, while nearest-neighbor (NN) AFM Ising interactions between the spins, of strength J z , are present along all dimensions of the lattice. In Eq. (1) i, j denotes a pair of NN sites in a twodimensional square lattice (every bond is counted once in the sum). Similarly, i, j x denotes a nearest neighbor bond oriented along the x-axis. We consider spin-1/2 particlesĉ j,σ with a hard-core constraint imposed by the projectorP onto the subspace without double occupancies. The statistics of the particlesĉ j,σ plays no role: By introducing Jordan-Wigner strings along the chains in x-direction one can switch between fermions and bosons. arXiv:2007.13759v1 [cond-mat.str-el] 27 Jul 2020 Symmetries and mapping to Z 2 lattice gauge theory.-Since holes cannot tunnel along y, their number N h y is conserved in each chain y. In the following we restrict ourselves to equal doping n h in every chain, N h y = n h L x with the system size L x(y) along x (y). In addition to the global U (1) ⊗Ly charge-conservation symmetries, and the conservation of total spin j S z j , the system exhibits hidden symmetries. Namely, when the holes move they only change the positions of the spins in the 2D lattice, while it is impossible to permute their configurations within any given chain. This is formalized by the notion of squeezed space, introduced to describe 1D doped spin chains [15,16]: To this end Fock states ⊗ y |σ (1,y) , ..., σ (Lx−1,y) , σ (Lx,y) , with σ j =↑, h, ↓ denoting local spin-and charge configurations, are re-labeled by j creates a hard-core chargon with the same statistics asĉ j,σ on the sites occupied by holes. The spin states in squeezed space are related to spins in the lattice bỹ where n h j denotes the chargon occupation numbers. After this re-labeling, the eigenfunctions of (1) become |Ψ = |Ψ ⊗ |Ψ c , where |Ψ = |{σj}j denotes a Fock configuration of spins in squeezed space and |Ψ c is a (generally correlated) chargon wavefunction [15]. Since we consider classical Ising interactions, every Fock configuration |Ψ defines a separate hidden-symmetry sector ofĤ. In the following we restrict ourselves to Néel states in squeezed space: |Ψ = |N ≡ |... ↑↓↑ ... .
If projected to the subspace |Ψ = |N , the Hamiltonian for the chargons (with densityn h j =ĥ † jĥ j ) becomeŝ where the sign of the tunnelling term is irrelevant.
To express the non-local (but instantaneous) interaction energyĤ int [{n h j }] in a compact form, we introduce the following string operators, By definition, each pair of holes is connected by a string of link variables τ x i,j = −1 [see Fig. 1 (a)] and the following Z 2 Gauss law is satisfied for all sites j: Owing to this Gauss law, the two link variables including a site j occupied by a spin σ j are equal, τ x j−ex,j x = τ x j,j+ex x = (−1) π j . Their value is given by the sublattice parity π j = 0, 1 of this spin, i.e., the number of times mod 2 the spin has switched sub-lattice (starting from a Néel state with all holes located on the right edge).
The Ising interaction between neighboring spins i, j y along y can be expressed in terms of the sub-lattice parities, J zŜ z iŜ z j = −J z (−1) π i +π j /4, since we use |Ψ = |N . Along the chains each bond i, j x gives J zŜ z iŜ z j = −J z /4 unless one of the sites is occupied by a chargon.
We proceed by promoting the link variables to a Z 2 lattice gauge theory (LGT) subject to the Z 2 Gauss law (5). This requires adding a termτ z i,j x in the tunneling term in Eq. (3) which correctly flips the sign of τ x has a concrete physical meaning as it can be measured from the local spin configuration.
Finally we arrive at the exact representation of Eq. (1), where A = L x L y is the total area. We introduced the dimensionless inter-chain coupling parameter α, which is α = 1 for our model in Eq. (1). Many-body phase diagram.- Fig. 1 (b) shows the phase diagram of the model in Eq. (6) as a function of temperature k B T and doping n h . The phase boundaries are estimated using a parton-based mean-field description (see supplements). Each phase is also found in our quantum Monte Carlo (QMC) simulations.
For the ground state (T = 0) we predict a vertical stripe phase, characterized by charge modulations with wavelength λ = 1/n h . The Z 2 electric field changes sign across each stripe, respecting the Z 2 Gauss law.
As a result, incommensurate long-range spin correlations are found along x, see Fig. 2 (a): The binding mechanism into stripes can be readily understood from Eq. (6): The interactions of the Z 2 electric field lines favor alignment of the latter along y, which is achieved by creating strong charge correlations along the y-direction. Such localization along y is cheap due to the absence of chargon tunneling in this direction. On the other hand, strong anti-bunching along x allows each chargon to delocalize as much as possible, in direct competition with the attraction of Z 2 electric field lines.
As shown in Fig. 2 (b), stripes are indeed characterized by long-range AFM order in the y-direction (correspond- ing to aligned Z 2 electric field lines): Numerically we find that long-ranged correlations C Y (d) develop below a non-zero critical temperature T S > 0. Our QMC simulations in Fig. 2 (c) and (d) show that k B T S ≈ 1.0(5)J for the chosen value of t/J z = 3 and 20% hole doping for linear system size L = 15. Within each chain our system has a conserved number of holes, associated with separate U (1) symmetries. In the long-wavelength limit, the corresponding effective field theory describes a U (1) symmetric field without quantum fluctuations of the charge along y. Integrating out thermal fluctuations at temperatures k B T > 0 yields an effective action of a 1 + 1 dimensional quantum system. With the global U (1) symmetry along y, we thus expect power-law correlations along x and y in the stripe phase: Below the critical temperature for stripe formation, T S > 0, these replace the infinite-range correlations Eqs. (7), (8) expected in the true ground state.
We find that our finite-size simulations with open boundaries are consistent with very weak power-law correlations C Y (d) when 0 < T T S . The detailed nature of the transition at T S remains a subject of future investigation, but we expect it to be in BKT class.
At higher temperatures and beyond a rather small critical doping value n h ≥ n c h (T ) we predict a chargon gas. It has no long-range AFM order in either direction, The loss of antiferromagnetism is entirely due to chargon dynamics, however: in squeezed space the spin wavefunction is still the classical Néel state. Hence the chargon gas is characterized by its hidden AFM order, which manifests itself in the non-local string correlations defined by the Z 2 Gauss law (5). Related string correlations have been observed in 1D Hubbard models [9,17] and are commonly used to characterize topological order in 1D systems [18,19].
In contrast to the stripe phase, the chargon gas is characterized by a disordered Z 2 electric field: Chargons are hence deconfined and form a gapless phase [20], corresponding to free fermionic holes at the meanfield level. Finally, at very low doping n h < n c h (T ) but above the critical temperature T > T S (n h ) for stripe formation, we predict a meson gas. It is characterized by a uniform Z 2 electric order parameter This should be contrasted to the T = 0 stripe phase with incommensurate magnetism, where e i,j x = 0 is modulated in space with a wavelength λ = 2/n h , such that i,j e i,j x = 0; in the stripe phase at 0 < T < T S the thermal average e i,j x = 0 is expected to be strictly zero in the thermodynamic limit. As a direct consequence of ν cc = 0, the meson gas has commensurate long-range AFM order along both directions, Physically, the meson gas can be understood as a paired phase of chargons. The Z 2 electric string connecting two chargons is associated with a linear string tension ∝ ν cc , which precludes one-chargon excitations in the thermodynamic limit; i.e. the meson gas corresponds to a confined phase which has, even in the zero-temperature limit, ĥ † Due to the restriction of chargon dynamics along one direction, the meson gas also corresponds to a Luttinger liquid with fractionalized excitations in the zero-temperature limit [20].
To identify the meson gas phase in our QMC numerics, we calculate histograms of chargon separations in Fig. 3. The hallmark of chargon-chargon meson formation is a narrow distribution p 2n−1,2n (r) of separations r between chargon numbers 2n − 1 and 2n, with n = 1, 2, ... and counting from the left edge, and a broader and different distribution p 2n,2n+1 (r) between chargons 2n and 2n + 1. This feature is clearly visible in Fig. 3 (a) in the expected low-doping regime, where we also find a non-vanishing Z 2 electric order parameter, τ x i,j x = 0.8842 (8). In the other two phases, the histograms show significantly different features, see Fig. 3 (b), (c). We also confirmed at J z = t and for fixed doping n h = 5% that the meson gas phase only exists at intermediate temperatures (see supplements). In the phase diagram, the meson gas is associated with an unusual re-entrant behavior as one increases temperature along a line of constant, but small, doping: at T = 0 the system has incommensurate long-range AFM correlations, which we predict to be destroyed by thermal fluctuations of the stripes when 0 < T < T S . When the meson gas phase is entered for T > T S , true long-range AFM correlations are restored. This counter-intuitive behavior is possible since AFM correlations are merely hidden in the intermediate fluctuating stripe regime.
Finally, we want to make a connection with Ref. [14], where a single mobile dopant but with SU (2) invariant Heisenberg interactions has been studied. It was found that the hole forms a magnetic polaron [5,[21][22][23] which can be understood as a meson-like bound state of a spinon and a chargon [24] connected by a geometric string of displaced spins [14]. Our meson phase is an analog of this finding but at finite hole concentration and for Ising-type interactions.
Methods.-Our calculations are based on a number of different but standard techniques such as bosonization, mean-field parton theory, and QMC simulations.
In order to work with a 1D field theory amenable to bosonization, our crucial approximation is the decoupling ansatẑ for i, j y NN along y, i.e. i x = j x . The different phases correspond to different solutions for V MF (i x ). These approximations are justified because we find the same phases in the quantum Monte Carlo simulations. We find the critical Luttinger parameter below which the ground state forms stripes to be large, K c = 8. We refer to the supplementary material for details.
Discussion and outlook.-In summary, we showed that the mix-D t − J z model can be directly mapped onto a Z 2 LGT. The many-body phase diagram of our model features in the ground state a stripe phase where the holes form vertical walls. Above a critical temperature T S , but within the Néel Z 2 gauge sector (which has the lowest energy at zero doping), we find two gaseous phases: a confined meson gas, with long-range AFM order, and a deconfined chargon gas with hidden AFM correlations. Experimentally, the model Eq. (1) can be realized in the large U/t-limit of a bosonic Hubbard model with a strong tilt ∆ t along the y-direction: The strong tilt suppresses resonant tunnelling of dopants along y, whereas the super-exchange mechanism remains intact in both directions [14,25]; to obtain AFM Ising interactions, one can use spin-dependent scattering lengths [25,26]. Rydberg atoms, which have already demonstrated Ising spin systems [27][28][29][30][31][32], are an alternative option: By using multiple hyperfine levels to encode both spin and charge degrees of freedom, our mix-D t − J z model should also be realizable; see also Ref. [33] for a discussion of generic t−XY Z models in polar molecules, and Refs. [34][35][36][37] for direct implementations of Z 2 LGTs. For all systems, we propose to start from a classical Néel state without holes, which can be doped with mobile charges e.g. by adiabatic deformations of the trapping potentials. This should guarantee that thermal fluctuations outside the gauge sector of our Z 2 LGT are negligible.
In spite of the overwhelming simplifications of our model, the presence of a stripe phase and a confinementto-deconfinement transition at elevated temperatures draws one's attention to the cuprates. A goal for future investigations is to study related models which are more closely related to the 2D t−J model: as a first step, other gauge sectors with domain walls in squeezed space -corresponding to spinons -can be considered. By replacing Ising interactions with SU (2) invariant Heisenberg couplings, a much richer model is expected and it remains to be seen if any connections to Z 2 LGTs can be drawn. Finally, the goal is to include charge dynamics along the second direction: this may provide an adiabatic route to the stripe phase observed in cuprates. To understand the many-body phase diagram, we bosonize the Z 2 LGT chains along x. As a starting point we consider decoupled chains (α = 0 in Eq. (6)) and study when interactions between the chains (α = 0) become relevant. For α = 0 the chargons form a Luttinger liquid with K ≈ 1 [38] which corresponds to the chargon gas phase; the weak NN attraction −J zninj /4 is expected to lead to slight deviations from a free Fermi gas with K 1 [39].
To continue working with a 1D field theory amenable to bosonization, we make a mean-field ansatz as written already in the main text, and writê for i, j y NN along y, i.e. i x = j x . The different phases we predict correspond to qualitatively different meanfield solutions V MF (i x ).
We start with stripes at zero temperature, characterized by periodic charge modulations with wavenumber k = 2πn h . Since the Z 2 electric field changes sign across every stripe, V MF (x + 1/n h ) = −V MF (x). Hence we obtain the following Fourier expansion, Using the Z 2 Gauss law (5) we can expressτ x j,j+ex x = ix≤jx cos(πn h (ix,jy) ) in terms of the chargon density. From the bosonization formulan h j n h − ∂ xφ (x j )/π we find that a non-oscillating term in x survives in the resulting effective action per chain, where S 0 is the free action with Luttinger parameter K.
The interaction in Eq. (15) is relevant for K < 8, showing that our model with K ≈ 1 for t J z should form stripes. The comparatively large critical value of K c = 8 indicates that long-ranged attractive interactions would be required to de-stabilize the stripe phase. The effective action in Eq. (15) highlights another peculiarity of the Z 2 LGT: the gapped elementary excitations correspond to 2π phase slips of the fieldφ(x), which in turn correspond to pairs of two chargons. This demonstrates that chargon excitations on top of the stripe ground state are confined and form mesons.
The meson gas phase corresponds to a uniform meanfield solution V MF (x) = ν cc = 0. In this case only V 0 = 0 in Eq. (14) and the resulting effective action S α has a term oscillating with cos(πn h x), making the inter-chain coupling α irrelevant in this phase -see Ref. [20].

B. QMC simulations
The quantum Monte Carlo simulations are based on a sampling of the path integral representation of the partition function with the help of worm-type updates [40], here in the variant of Ref [41]. The structure of our Hamiltonian in Eq. (1) is similar to a twocomponent Bose-Hubbard system with nearest-neighbor density-density interactions, which is standard. The Néel gauge sector is implemented via the initial configuration as well as making sure that the imaginary time difference between the worm ends never exceeds β/2 in absolute value. This also makes the simulations canonical.
We observed rather large autocorrelation and thermalization times in the simulations. We attribute this to the physics of the phases we found: In the stripe phase, effectively shifting the position of the hole wall has a large barrier, whereas in the meson phase the binding energy of two holes makes it likewise difficult to move such pairs around.
In Fig. 4 we show more QMC results for the chargon distance histograms. We consider t/J z = 1 and fix the doping to n h = 5%. For the lowest temperature (a) we find no significant differences between the distributions p 1,2 (r) and p 2,3 (r) of separations r between particles 1 and 2 (respectively, 2 and 3). For intermediate temperatures k B T = J z (b) we find that particles 2 and 3 are significantly further apart than 1 and 2 or 3 and 4. This is a signature of the meson gas phase. At higher temperatures k B T = 20J z in (c) we find a significantly broader distribution where the distances between particles 1 and 2, 2 and 3, and 3 and 4 behave very similarly. These features indicate a chargon gas phase.
In Fig. 5 we show static spin and charge structure factors in the stripe phase. As expected we observe sharp delta-peaks at ±2πn h in the charge sector, and ±πn h in the spin sector relative to the (0, 0) and (π, π) values expected in the hole-free background, respectively -another hallmark of stripe formation. We also checked (not shown) that in the meson phase the static charge and spin structure factors show peaks at (0, 0) and (π, π), re- spectively -as expected.

C. Microscopic parton theory
In this section we introduce the microscopic mean-field descriptions of parton phases which we used to determine the phase diagram in Fig. 1 (b) of the main text. We focus on the low doping regime n h 1, except for the chargon gas phase. Our starting point is the Hamiltonian Eq. (6) in the gauge sectorĜ j |ψ = |ψ for all j, and we ignore finite-size effects in the following. The details of our calculation in Fig. 1 (b) of the main text are described in I C 4.

Mean-field theory of stripes
Our strategy in this section is to describe how individual stripes along y can form. These stripes contain exactly one hole per chain and stripe; they are expected to have long-range charge correlations along y. The basic binding mechanism is provided by the attractive interac-tion of parallel Z 2 electric field lines. For simplicity we neglect interactions between stripes at different positions x 1 and x 2 . To describe the individual stripes, we focus on the simplest case with exactly one hole per chain.
Zero temperature.-We describe the one-stripe ground state by a variational Gutzwiller wavefunction, It represents a product of identical chargon wavefunctions |φ (0) y per chain y. We can express the latter as where |0 denotes the undoped Néel state. Assuming that the stripe is centered around x = 0, we denote the basis states by | y in the second line of Eq. (17): here can be understood as the length of the Z 2 electric string measured relative to center of the stripe.
To approximate the ground state of the stripe, we optimize the wavefunction φ (0) in order to minimize the variational energy: where the inter-chain potential is given by The optimization problem is easy to solve numerically, by introducing a maximum string length | | ≤ max in the one-particle Hilbertspace | . In Fig. 6 we show the probability amplitude of the optimized wavefunction for different values of t/J z . From the variational ground state energy Ĥ 0 we obtain the energy per hole in one stripe, which we measure relative to the undoped Néel state: Instead of variationally optimizing φ (0) to minimize the overall energy, we can derive a self-consistency mean-field equation: Averaging over neighboring chains, the effective mean-field Hamiltonian per chain becomeŝ where the self-consistent mean-field potential is We find that the variationally optimized solution is an eigenstate of the mean-field Hamiltonian. Phonon spectrum.-To describe the effects of thermal fluctuations later on, the excitation spectrum of stripes is required. To this end we make a variational ansatz with momentum k y along the y-direction. Here |φ (0) y is the optimized ground state wavefunction obtained by minimizing Ĥ in Eq. (18). The resulting phonon dispersion is given by Since components in |φ (1) which are not orthogonal to |φ (0) vanish in Eq. (23), we may assume Using this relation, Eq. (24) yields where (lower indices denote chain numbers): , and we defined Note that the k y -dependence of ∆ S (k y ) and J ex (k y ) enters because the optimized wavefunction |φ (1) (k y ) depends on k y (we dropped the explicit dependence in the above expressions for clarity).
Lattice gas model.-When more than one hole per chain is considered we expect stripes to repel each other weakly along x-direction: Chargons are mutually hard-core particles, which leads to an additional localization energy cost when two chargons come close to each other. To estimate the strength of the repulsive stripe-potential, we calculate the ground state energy per hole for a stripe in a box from x = − max ... max , see Fig. 7. On one hand this gives an estimate of the energy cost to localize stripes within a region ± max . On the other hand our results show for the dilute regime (no stripe interactions) that our calculations quickly converge with max .
The classical energy of multiple stripes at positions x j , where j = 1...L x n h , can be modeled by this corresponds to the Hamiltonian of a classical 1D Ising model (lattice gas model). We estimate the effective potential V S (d) ≈ 0 (d) − 0 (∞) from the energy per hole for one string in a box, see Fig. 6. For the smallest distances V S (1) ≈ 0.2t is a small fraction of t. Since the solution φ (0) asymptotically behaves like the Airy function Ai( ) (due to the linear string potential Eq. (19) at large distances), this potential decays super-exponentially: The ground state of the model Eq. (27) has long-range order (charge density wave, CDW) along x.
Non-zero temperatures: phonons.-Next we discuss the stability of the stripes forming along y when thermal fluctuations are included. The phonon excitation spectrum of a single stripe, see Eq. (24), is gapped: ∆ = ∆ S (0) − 2J ex (0). Therefore we expect that stripes are robust as long as k B T ∆. In the limit of small doping n h → 0 we show below, however, that the meson gas has a smaller free energy.
To estimate the free energy per hole in the stripe phase, F S /N h , we neglect phonon (b ky ) non-linearities and describe fluctuations on top of the stripe phase by the free phonon Hamiltonian This leads to an ideal-Bose gas contribution to the free energy of There is a second contribution F S,cl to the free energy, related to classical fluctuations of the locations x j of the stripes, j = 1...n h L x . This contribution to the free energy vanishes in the thermodynamic limit L x , L y → ∞. To see this, we give a simple estimate by focusing on the entropic contribution: F S,cl ≈ −k B T log(M S ) where M S denotes the number of microscopic stripe configurations. In the dilute regime when n h 1 we expect that the energetic contribution beyond the ground state is negligible since the stripe interaction potential decays super-exponentially, see Eq. (28). Assuming further that the width of each stripe, σ S = [ φ (0) |ˆ 2 |φ (0) ] 1/2 , is small compared to the average distance, σ S 1/n h , we obtain: This result is independent of L y , and thus sub-dominant compared to F S,ph . Combining our results, we obtain where S 0 is the ground state energy per hole relative to the undoped Néel state, see Eq. (20). Eq. (32) is valid when the average distance of stripes is well beyond their width, n h σ S : this justifies neglecting interactions between stripes.
Non-zero temperatures: Mean-field.-At higher temperatures k B T ∆ S − 2J ex , the inclusion of only the lowest phonon band (29) is no longer sufficient. To describe thermal fluctuations in this regime, we generalize the mean-field ansatz from Eq. (16) to a product of density matrices: whereĤ MF is the mean-field Hamiltonian from Eq. (21) with a temperature dependent potential V eff ( ; T ). As for the ground state, one can either solve for the mean-field stateρ S,MF by determining the mean-field potential self-consistently: Alternatively, we choose to treat V eff ( ) as a set of free parameters which determineρ (0) and can be optimized variationally. To this end we minimize the free energy per hole (relative to the energy of the undoped Néel state): Here E MF n denote the eigenenergies of the mean-field HamiltonianĤ MF + L x J z /2 relative to the undoped Néel state and the mean-field string-length distribution is given by This string length distribution is shown for the variationally optimized mean-field state in Fig. 8.

Mean-field theory of the meson gas
Next we describe individual mesons, composed of two chargons, in the meson gas regime. We start by considering a single pair of two holes, and extend our results to finite densities afterwards.
Two individual chargons.-Our starting point is the undoped Néel state, for which τ x i,j ≡ 1 everywhere. Next we add two chargons in the same chain, which are connected by a Z 2 electric string along which τ x i,j = −1. This leads to a linear potential depending on the distance of the two chargons, where dE/d = J z is the linear string tension and m = 2J z a zero-point energy shift. Note that the string length = 0 cannot become zero since chargons are mutually hard-core.
Since the motion of chargons is restricted to one dimension, along x, we can treat them as distinguishable. Transforming into the co-moving frame with the left chargon yields the following effective Hamiltonian in the relative coordinate > 0: (38) Here energies are measured relative to the undoped Néel state; k x denotes the total conserved momentum of the meson and the effective tunnel coupling is given by Using exact diagonalization, the spectrum n cc (k x ) of (38) can be easily computed. This yields the ground state, n = 0, and allows to include thermal fluctuations by including higher vibrational excitations n ≥ 1 with their respective Boltzmann weights.
Finite doping.-To calculate the free energy per hole at non-zero doping, n h > 0, we treat the chargon-chargon mesons as point-like hard-core bosonsd ( †) kx,n . This is justified if the average string length cc = | | is small compared to the average meson distance: cc 2/n h . By introducing Jordan-Wigner strings, the mesonŝ d n denotes the internal vibrational states of the mesons. Neglecting interactions between mesons in different chains for now, and working in the grand-canonical ensemble, we obtain the following free energy per hole: As usual, we define F Mean-field interaction effects.-On a mean-field level, interactions between mesons from different chains can be easily taken into account: They lead to a reduction of the linear string tension dE/d and the zero-point energy shift m. To calculate the renormalized meson potential parameters we start from the mean-field Hamiltonian: whereĤ t denotes the hopping part of chargons; note a factor of 2 in the second term which accounts for interactions with both neighboring chains in y-direction.
In the meson mean-field phase, the average Z 2 electric field is related to the average meson string length: Using this result we read off from the mean-field chargon Hamiltonian, Eq. (43): Note that the resulting chargon potential V MF cc ( ) from Eq. (37) counts the mean-field energy Eq. (43) of mesons relative to the energy E MF 0 of an undoped chain described by the same mean-field Hamiltonian (the latter including doping in the neighboring chains).
To evaluate the resulting free energy, F cc = Ĥ MF − T S MF , we note that the mean-field Hamiltonian Eq. (43) yields F MF cc = Ĥ cc MF MF − T S MF . This gives Since the true and mean-field Hamiltonians only differ by their potential energy terms, we can write the free energy of the mean-field state -now measured relative to the energy E N 0 = −L x L y J z /2 of the undoped Néel state -as: Here f MF,0 cc is the mean-field meson free energy per chain (described by the chargon potential Eq. (37) with meanfield parameters) measured relative to the energy of an undoped chain in the mean-field Hamiltonian. Simplifying all terms in Eq. (49) yields The mean-field free energy obtained per meson, , can be calculated as in Eqs. (41), (42) but using the renormalized mean-field chargon spectrum n cc (k x ) → MF,n cc (k x ): and:

Mean-field theory of the chargon gas
For the mean-field theory of the chargon gas we assume that chargons form independent free Fermi gases in each chain: with the free chargon dispersion Ignoring the weak nearest-neighbor attraction −J z /4 leads us to the variational ansatz of a free Fermi gas: with k F = πn h denoting the Fermi momentum, and similarly for density matrices at non-zero temperatures. This leads to the following free energy per hole, measured relative to the undoped Néel state. Here we introduced the nearest-neighbor chargon correlator, and the chargon chemical potential µ c is determined by

Results
Now we compare the results collected in the last few subsections. We start at zero temperature, where the free energy F = E reduces to the (variational) ground state energy. In the following we calculate energies per hole, where N h = L x L y n h denotes the total number of holes in the system. All energies are measured relative to the undoped 2D Néel state.
Mean-field phase diagram.-The mean-field phase diagram in Fig. 1 (b) of the main text indicates which of the three considered phases has the lowest free energy. In the calculation of the meson gas free energy, we included meson-meson interactions on a mean-field level, using Eq. (52). We used the finite-temperature Gutzwiller ansatz, Eq. (33), for the stripes, but found essentially identical results by using the one-phonon approximation Eq. (30).
We only calculated the meson gas free energy in a regime where the average string length cc < 1/n h , which is required to assume point-like mesons in our calculation. Likewise, we only calculated the stripe free energy in regimes where ˆ < 0.5/n h , beyond which stripestripe interactions must be included. This leads to the cusps and re-entrant behaviors we find in the mean-field phase diagram. We do not expect these to be physically meaningful.
Zero temperature.-Neglecting interactions between stripes along x, we find for the ground state energy per hole in the stripe phase, see Eq. (32): For the meson gas we obtain a similar result, see Eq. (51): where we ignored mean-field effects yielding O(n h ) corrections. For the chargon gas Eq. (60) gives The last term captures finite-density corrections O(n h ). Comparing Eqs. (63) -(65) we note that the energy per hole of the chargon gas diverges when n h → 0. This is understood by noting that the chargon gas has vanishing AFM correlations along y at any non-zero doping n h > 0, and thus does not approach the 2D Néel state when n h → 0. This leads to a total energy cost O(L x L y ) at n h = 0, which translates into a diverging energy per hole in this regime. We conclude that the chargon gas can only exist beyond a critical doping n ch h > 0; below n ch h only the stripe and meson gas compete.
To compare stripes and the meson gas, we first compare their asymptotic ground state energies per hole at infinitesimal doping: S 0 and 1 /2 n=0 cc (k x = 0). We find that the stripe consistently has a lower energy per hole than the two-chargon meson.
Low-temperature, low doping: meson gas vs. stripes.-To understand the behavior at small but non-zero temperatures, we expand the stripe free energy in k B T . Further, since we also focus on small doping n h 1, we only compare stripe and meson gas phases relevant to this regime in the following.
For the stripe phase we can capture thermal fluctuations using low-energy phonon excitations, see Eq. (32). The non-zero phonon gap, ∆ = ∆ S (0) − 2J ex (0), leads to a finite activation energy and corresponding exponential suppression of thermal fluctuations for stripes: For the meson gas in the same regime, we can neglect vibrational excitations with n>0 cc > n=0 cc and consider a free Fermi gas (note that mesons are hard-core particles moving in 1D chains) in the lowest band 0 cc (k x ). Expanding 0 cc (k x ) = 0 cc + k 2 x /2m cc and working in the low doping regime, we obtain a classical gas with free energy per hole: Comparing Eqs. (66), (67) yields the transition temperature from stripes to the meson gas when n h 1: High-temperature expansion: chargon vs. meson gas.-We have already seen in the previous paragraph that stripes have too low entropy at high temperatures. Therefore we only consider the gaseous phases in the following.
At high temperatures, the free energy per hole of the chargon gas, Eq. (60), becomes This term is still dominated by the first divergent contribution at low loping n h → 0.
To estimate the free energy per hole at high temperatures k B T t and low doping n h 1/ cc for the meson gas phase, we note that the meson dispersion can be approximated by n cc (k x ) ≈ nJ z (70) when n cc t. In this regime, the eigenstates correspond to Wannier-Stark states in the linear potential V cc ( ) ≈ J z . The corresponding free energy per hole is purely classical and becomes: The average string length can be estimated in a similar manner, when k B T t. Comparison of Eqs. (69), (69) shows that the meson gas is thermodynamically favored in the low doping, high temperature regime.