Emergent devil's staircase without particle-hole symmetry in Rydberg quantum gases with competing attractive and repulsive interactions

The devil's staircase is a fractal structure that characterizes the ground state of one-dimensional classical lattice gases with long-range repulsive convex interactions. Its plateaus mark regions of stability for specific filling fractions which are controlled by a chemical potential. Typically such staircase has an explicit particle-hole symmetry, i.e., the staircase at more than half-filling can be trivially extracted from the one at less than half filling by exchanging the roles of holes and particles. Here we introduce a quantum spin chain with competing short-range attractive and long-range repulsive interactions, i.e. a non-convex potential. In the classical limit the ground state features generalized Wigner crystals that --- depending on the filling fraction --- are either composed of dimer particles or dimer holes which results in an emergent complete devil's staircase without explicit particle-hole symmetry of the underlying microscopic model. In our system the particle-hole symmetry is lifted due to the fact that the staircase is controlled through a two-body interaction rather than a one-body chemical potential. The introduction of quantum fluctuations through a transverse field melts the staircase and ultimately makes the system enter a paramagnetic phase. For intermediate transverse field strengths, however, we identify a region, where the density-density correlations suggest the emergence of quasi long-range order. We discuss how this physics can be explored with Rydberg-dressed atoms held in a lattice.

Introduction.-Recent years have seen a great success in emulating and studying many-body Hamiltonians with cold atomic systems (see e.g.[1]).One feature which is currently receiving much attention is the formation and exploration of crystal structures in systems composed of cold trapped ions [2], polar molecules [3][4][5] and also gases of Rydberg atoms [6,7] in all of which constituents interact with a 1/r α potential (e.g.α = 3 for dipolar interactions and α = 6 for van der Waals (vdW) interactions).In particular Rydberg gases have witnessed a recent experimental breakthrough which enabled the adiabatic preparation of a crystalline state [8].
The above-mentioned systems have in common that in the case of repulsive interaction the corresponding potential is strictly convex, i.e. particles lower their interaction energy by assuming an as uniform as possible distribution in space.Remarkably this property makes the ground state configuration independent of the actual details of the interactions, e.g. the specific power-law, and leads to the formation of a so-called generalized Wigner crystal [9,10].Furthermore, when the particles are held in a lattice, the permitted filling fractions of the ground state configuration form a fractal curve known as the complete devil's staircase [11].
The case of non-convex interactions, however, is far less understood.Although ground states can in principle feature crystalline structures [12], there is currently a quite limited understanding of the mechanism behind their formation and their melting in the presence of (quantum) fluctuations.Systems featuring non-convex interactions are know to exhibit a quite intricate behavior.A wellknown example is competing short-range attractive and long-range repulsive interactions.This kind of interactions was considered e.g. in classical fluids and leads to a much more variate phase diagram than the typical Lennard-Jones potential [13,14].The competition between the attractive and the repulsive terms in the po-tential leads to the presence of cluster phases, as found e.g. for mixtures of colloids and polymers [15,16].g In this work we study the crystal structures in the ground state of a 1D spin model in a regime where spins are coupled by competing attractive and repulsive interactions.Using Rydberg dressing of electronic ground state atoms [17,18], we engineer the potential in a way arXiv:1504.07934v1[cond-mat.quant-gas]29 Apr 2015 that the nearest-neighbor interaction W can be tuned from repulsive to attractive while the longer range interactions follow the usual vdW form V (k) = C 6 /(ak) 6 for each k ≥ 2 with C 6 and a being the dispersion constant and the lattice spacing.For negative W we show that in the classical limit the filling fraction of excitations follows a staircase that can be explained in terms of the union of two staircases, the first of dimers of excitations (for filling fractions f ≤ 1/2), and the second of dimers of holes (for f ≥ 1/2).In the region f > 1/2 we observe the presence of larger aggregates of excitations, and we pinpoint the region of stability of these polymers of different lengths.We eventually perform a qualitative study of the melting of the crystal state through quantum fluctuations introduced by a transverse field finding an extent of the density-density correlation function which hints towards quasi long-range order.
The model.-We consider the one-dimensional spin model described by the following Hamiltonian where ni = (ŝ z i + 1)/2 = |1 i 1| and ŝx i = (|1 i 0| + |0 i 1|) with ŝx and ŝz being Pauli matrices.A site i can either be occupied |1 i , or empty |0 i .The two states are coupled by a transverse magnetic field with strength g, and interactions occur only among occupied sites.In what follows we set a as the length unit without loss of generality.
The key feature of the spin model Eq. ( 1) is a controllable non-convex interaction, which can be realized by Rydberg dressing of electronic ground state atoms (more details on the dressed potential are discussed at the end of this work).The relevant levels are depicted in Fig. 1(a), where the two electronic ground states |0 and |1 (typically Zeeman or hyperfine states) form the pseudo-spin and are Raman coupled with coupling g.Rydberg dressing of the atoms in state |1 induces a potential between these atoms, whose shape is illustrated in Fig. 1(b).In an optical lattice, this allows us to establish a long-range non-convex interaction which is attractive for nearestneighbor excitations while it follows a vdW form for more distant excitations.Note that by controlling the lattice spacing and the distance R res of the two-atom resonant excitation, it is possible to put more lattice sites on the attractive branch (see Fig. 1(b)) which realizes a whole class of models with non-convex potentials.However, in the present work we focus on the case where only the nearest-neighbor interaction is attractive.
The generalized Wigner crystal.-Letus firstly determine the ground state configuration of the model Eq.
(1) in the classical limit (g = 0).One emerging feature is that the crystal formed in the ground state depends on the sign of the nearest neighbor interaction W .In the case of positive W the introduction of an excitation is energetically unfavorable and the ground state of Eq. ( 1) is the empty state ( i ⊗|0 i ).On the other hand in the limit of infinite negative W the ground state is the ferromagnetic state ( i ⊗|1 i ).We want to study how the ground state interpolates between these two cases by concentrating on negative W .Note that the nearest neighbour interaction W induced staircase studied in this work is in contrast with all the previous studies (for both convex [9,11] and non-convex interactions [12]) where the ground states are driven by a chemical potential term.This difference has a crucial effect on the physics involved, e.g., for the model studied in [12], because the chemical potential term can absorb the linear terms of the particle-hole transformations and the model has an explicit particle-hole symmetry.However, since our system has no manifest particle-hole symmetry, the emergence of a dimer hole staircase for f > 1/2 is rather unusual and related to the convex vdW tail itself [19].We first perform a numerical analysis of the infinite chain as follows.We consider configurations which are periodic with period q, i.e. considering rational filling fractions f = p/q, with p ≤ q, and we set the maximal q = 20.For each given W we compute numerically the energy density associated with each configuration, that is we evaluate the energy of a single period and then we divide it by the period length.We eventually identify the ground state as the minimum energy density configuration.In the simulation we also set V (2) as the energy unit, which is used throughout the paper.Using this method, we obtain a staircase pattern connecting the empty state (f = 0) at W = 0 with the completely filled state (f = 1) for large negative W , which is shown in Fig. 2(a).The data show that we never observe single excitations, but rather polymers, as aggregating at least two excitations is always energetically favorable for negative W .We find that the excitations arrange in dimers for filling fractions f < 1/2, which is in agreement with the results of [12].Once the threshold of f = 1/2 is reached the ground state configuration becomes

R-mer
Inspired by the fact that two excitations are always on neighbouring sites when filling fraction f < 1/2, we develop an effective dimer theory to describe the resulting crystal structure.
Let us consider the presence of a dimer particle (dp) on sites (i, i + 1) as an effective single particle with associated dimer particle number operator pi = ni ni+1 .The effective Hamiltonian of the dimers reads where W assumes the role of a chemical potential for the dimers.Here V p (r) = 2V (r) + V (r + 1) + V (r − 1) is the interaction potential between two dimers separated by r sites, which is convex.Note, that for f < 1/2 the convexity of the potential ensures that no neighbouring dimers are allowed.We can now employ a method due to Hubbard [11] to identify stability regions of different densities of dimers with respect to W .The result is plotted in Fig. 2(b).One observes a remarkable agreement between the dimer model and the numerical result by solving the full Hamiltonian.
In the region f > 1/2 larger polymers form and the numerical results suggest that the ground state exhibits a rather intricate crystal structure.In fact we find regions in which only one kind of polymer is stable, and regions in which different kinds of polymers coexist.Here the cases in which only polymers of given R neighboring excitations (so-called R-mers) are present correspond to the wider plateaus in the staircase (black data points in Fig. 2(c)).The "mixed" cases correspond to the narrow regions interpolating between two plateaus (cyan data points in Fig. 2(c)).
In order to simplify the description, we neglect the presence of mixtures as the dominant plateaus correspond to a single kind of polymers.Our aim is to determine when the crystal changes from R-mer to R + 1-mer as we change W , which is achieved by examining the stability of R-mers against W . Considering an R-mer as a single particle, we define the effective filling fraction of R-mers as f eff = P/Q, which is related to the real filling fraction f = p/q through f eff = p/[R(q − p) + p].In the following, we consider f eff = 1/Q for reasons that become clear later on.The energy of the R-mer in the thermodynamic limit reads where is the bond energy of the R-mer, which is associated to the mutual interactions of its internal excitations.
is the interaction energy between two R-mers separated by x (x = lu), where u is the length of the unit cell.
When the crystal changes from R-mer to R+1-mer, the corresponding unit cell length changes from u = Q+R−1 (for R-mer) to u = Q+R (for R+1-mer) at the same f eff [20].The infinite chain can be divided in equal periods on which either u R-mers or u R + 1-mers are disposed, so that we can compute the transition point by solving In order to carry out this computation, one however needs to know the filling fraction at which the transition takes place.Both theoretical and numerical arguments suggest f eff = 1/3 for all transitions [20].The comparison of W obtained from the simple theoretical model Eq. ( 4) with the numerical data is shown in Fig. 2(c).One can see that the transition values of W are in excellent agreement.
Here f eff = 1/3 also means that there are one R-mer and two neighboring holes in each unit cell.This stimulates us to search for an effective theory of dimers of holes (dh).In analogy to the dimers of particles Eq. ( 2), the effective Hamiltonian reads where V h (r) = V p (r) is the interaction energy between dimer holes, ĥi = (1 − ni )(1 − ni+1 ) is the dimer hole number operator and µ h = 3W + 4 ∞ j=2 V (j) can be regarded as the chemical potential.We again study the phase diagram with the Hubbard method as the interaction is convex.In Fig. 2(d) we report the filling fraction of excitations found in this way, and we compare it with the numerical results of the full Hamiltonian, finding a perfect agreement in the region f > 1/2.Especially, from the expression of the chemical potential µ h , we can find the transition value of W to the ferromagnetic state exactly, i.e., W c = −4/3 ∞ j=2 V (j) = −1.47994.Melting of the crystal.-In the case of a convex interaction, it was shown in [6] that the introduction of a quantum term leads to a two-stage melting process when the excitations are dilute.For weak transverse fields the crystal structure melts into a floating solid which will eventually melt into a paramagnet for higher values of g.Here we explore the qualitative features of the melting process in the case of a dense crystal, namely at filling fraction f = 1/2.We perform numerical studies by means of Matrix Product States implementations of the Density Matrix Renormalization Group (DMRG) [21] and its infinite version iDMRG [22].The results are presented in Fig. 3.The left panel shows the transverse magnetization i ŝx i /L as function of W and g, where L is the length of the chain.The lobe structure is similar to the one found in [6,23] and displays two main lobes corresponding to f = 1/2 and f = 1 (Note, the fine structure of the staircase is too narrow to be resolved at the precision level of the plot).
We use the density-density correlation function g (2) (k) = ni ni+k − n 2 (which is independent of the initial site i for an infinite system) to characterize the melting process, and study its behavior when exiting the f = 1/2 lobe along the black line in Fig. 3. Performing a numerical simulation for very small values of g (close to the classical limit) is challenging due to the presence of quasi translational invariance which leads to a highly degenerate ground state.We consider g = 0.05 (see point A in Fig. 3), which allows the degeneracy to be slightly lifted while still inside the lobe, finding a periodic g (2) (k).In our case (f = 1/2) the ground state is quasi-four-fold degenerate and for small values of g has the structure |ψ g ∝ (|..1100.. + |..0110.. + |..0011.. + |..1001.. )/2 + O(g 2 ).This leads to the behavior of g (2) (k) = cos (kπ/2) /4 + O(g 2 ) which reproduces remarkably well the results shown for point A in Fig. 3.
When g is sufficiently large, the ground state of the system is in a paramagnetic phase and g (2) (k) shows an exponential decay (point C in Fig. 3).By considering intermediate values of g it is possible to study the melting of the crystal phase.We find that in this regime (point B in Fig. 3) the correlations are significantly longer ranged, and in the right panel of Fig. 3, we show a typical result with bond dimensions of 100.The simulated correlations in the vicinity of the point B show similar behaviour (not shown), which hints towards the presence of quasi longrange order in a region extending from g 0.1 to g 0.2.The detailed numerical study of this melting process is left for future work.
Summary.-In conclusion, we have studied a quantum spin model with long-range non-convex interactions which is exactly solvable in the classical limit.The ground state of the system has an intriguing structure, in the sense that it can be explained by two complete devil's staircases of dimer particles and dimer holes which merge at filling fraction f = 1/2.We also studied some aspects of the melting of the staircase induced by a transverse quantum term.Our numerical results suggest the emergence of quasi long-range order as the dimer crystal melts into a paramagnet.The exact extent of this phase, and the precise nature of the melting mechanism will be matter of future research.Moreover, we have also shown that the theoretical predictions in this work can be veri-fied experimentally by using dressed interactions between Rydberg atoms.Finally we remark that in the Rydberg dressing scheme we propose the attractive range can be tuned freely, and extended over more lattice spacings.4) with f eff = 1/3.

R-mer staircase
In the main text, we have pointed out, that in order to determine the transition points between different R-mers using Eq.( 4), one needs to know the filling fraction at which the transition occurs.Apart from determining the transition filling fractions from numerical simulation, one can construct a purely theoretical self consistent argument.To show how this works, we first determine the stability regions of a gas of R-mers.The boundaries of a stability region are determined from the energy condition where the energies of a certain configuration of R-mers and the energy of a new configuration where a single R-mer was subtracted/added to the chain are equal (we refer the reader to the original work [11] or to [34,35] for more details; the generalization to R-mers is straightforward).The expressions for the stability regions of a chain of R-mers read r l E int R (r l − 1) − (r l − 1)E int R (r l ) (6) where E b,int R = R−1 l=1 (R − l)V l is the repulsive interaction energy between the constituents of a single R-mer and r l is the distance between first spins of l-th nearest R-mers.The expressions Eqs.(6,7) are derived under the assumption that the interaction between the R-mers is repulsive, i.e. in our case it breaks down when the configuration of the gas can allow two R-mers sitting next to each other, which happens for high filling fractions.We thus need to restrict to filling fractions, where Eqs.(6,7) are valid (no adjacent R-mers).This is ensured by condition r 1 > R + 1.We plot the devils staircase of different R-mers in Fig. 4.
With the stability regions at hand, one should now compare the ground state energies of all R-mers configurations which are stable at a given value of W in order to decide which R-mer provides the minimal energy.Similarly, one can use the condition Eq.( 4) (when restricting to filling fractions of the form f eff = 1/Q).To illustrate this on an example, lets consider e.g.filling fraction f = 0.5 which admits both dimers (see Fig. 4, black line, f eff = 1/3) and trimers (red line, f eff = 1/4).The condition Eq.( 4) then yields W = −1.6858,which is inconsistent with the given stability region W ∈ [−0.8319, −0.5759] of trimers.In this way one can show that as W is decreased, the chain is filled with dimers only down to W = −1.3066,where the first transition to trimers occurs (orange cross at f = 0.5 in Fig. 4).Following this logic it is then straightforward to show, that all the remaining transitions between R-mers and R + 1-mers occur at f eff = 1/3 and quickly approach the asymptotic value W (R → ∞) ∼ = −1.4799.The values of W for some of the transitions are shown in Table I and plotted in Fig. 2.
To summarize, the simple theoretical model where one considers a chain of R-mers only yields results fully consistent with numerical simulations as well as with the description in terms of dimers of holes, which actually appears naturally when the condition Eq.( 4) is combined with the devils staircase picture -all the transitions happen well inside the large stability regions, which all correspond to f eff = 1/3, i.e. dimers of holes.

R res 1 2 FIG. 1 .
FIG. 1. (color online) (a) Level scheme.Two ground states |0 and |1 are coupled coherently through a Raman process.The |1 state is weakly dressed to two Rydberg states |R1 and |R2 by two blue-and red-detuned lasers.(b) Interaction potential between the Rydberg dressed |1 states.Varying the laser parameters and the lattice spacing a, we achieve a non-convex interaction where the nearest neighbor interaction is attractive and the longer range interactions are repulsive.The blue-detuned laser induces the two-atom resonant excitation at a distance Rres (dashed line in the figure).See text for details.(c) Tuning the non-convex interaction by changing the laser detuning ∆2.In the example, two Rydberg states |R1 = |60S and |R2 = |70S of Rubidium atoms are considered.See text for more details.In (b), we assume ∆2 = 250 MHz.Other parameters used in plotting (b) and (c) are Ω1 = 30 MHz, ∆1 = −300 MHz, Ω2 = 22 MHz and a = 2 µm.

FIG. 2 .
FIG. 2. (color online)The black curves in all panels are the numerically extracted staircases taking into account any rational filling fractions f = p/q with q ≤ 20 where the range of the interaction has been cut off at k ≤ 1000 lattice sites.(a) Staircase pattern of the filling fraction f as function of W .The interactions from next-nearest neighbors onwards follow the form V (r) ∝ 1/r 6 .The ground state can be thought of as arrangements of dimer particles (f ≤ 1/2) and dimer holes (f ≥ 1/2) (as indicated by the shaded areas), which merge at filling fraction of f = 1/2 where the spin configuration is given by • • • 11001100 • • • .The green curves are the results from the effective models based on respectively dimer particles, (b), and dimer holes, (d).(c) The red stars are the transition points between different kinds of R-mers as given by Eq. (4).Plateaus shown in cyan correspond to mixtures of different R-mers.

FIG. 3 .
FIG. 3. (color online) Quantum melting of the dimer crystal with increasing transverse field strength g.Left panel shows the transverse magnetization i ŝx i /L as function of W and g, where the simulation was done for L = 60 lattice sites with open boundary conditions and the range of the interactions has been truncated to 10 lattice sites.Right panel shows the iDMRG calculations of the density-density correlation function g (2) (k) vs. k for three different values of g = (0.05, 0.15, 0.3) at W = −0.91 corresponding to the points A,B,C in the left panel (see text for details).

FIG. 4 .
FIG.4.(color online) The devil's staircase for a gas of R-mers.The black, red, green and blue staircases correspond to the gas of dimers, trimers, pentamers and octamers respectively.The orange crosses indicate the transition points between different polymers (e.g. the transition point between dimers and trimers is represented by the orange cross at filling fraction 0.5).As R is increasing the transition value W quickly approaches the asymptotic value W (R = ∞) ≈ −1.48, which is indicated by the orange dashed line.The parts of the staircase for higher filling fractions (after the large plateaus as W is decreased) are not shown as they are not captured by our simple theory (see text for details).

TABLE I .
Whether this introduces new physics stands as an open point.Values of W where a transition between R/(R + 1)mers occurs, determined from Eq.(