Frustrated plane-polarized dipoles in one dimension
FFrustrated plane-polarized dipoles in one dimension
Niraj R. Ghimire and Susanne F. Yelin
1, 2 Department of Physics, University of Connecticut, Storrs, CT 06269, USA Department of Physics, Harvard University, Cambridge, MA 02138, USA (Dated: July 29, 2019)We investigate the zero-temperature quantum phases of a quasi-one-dimensional zigzag chain ofdipoles that are polarized in a plane by an external electric field. Since the Hamiltonian containsnearest-neighbor (NN) and next-nearest-neighbor (NNN) hopping and interaction terms, this modelallows frustration which induces phases that can be interesting and unusual. By using the densitymatrix renormalization group (DMRG) algorithm, we produce a complex phase diagram. This isan extension of an earlier work by Wang et. al. [Phys. Rev. A , 043615 (2017)]. I. INTRODUCTION
Ultracold atoms in optical lattices serve as ideal plat-form for quantum simulation, which is known to be adifficult problem even for the most advanced supercom-puters of today, especially when the system size is large .Because the geometry, dimension, and depth of an opti-cal lattice can be controlled to a high degree, ultracoldatom-based simulators have already been used to inves-tigate quantum many-body problems applicable to fieldsranging from condensed matter physics to high energyphysics . Although atoms interact via short-range con-tact interactions in most cold atom experiments, many-body systems with longer-range interactions are pre-dicted to exhibit intriguing quantum phases .In the presence of geometrical frustration, a situationwhere not all the interactions are satisfied, the systemis expected to exhibit even more interesting features.For instance, quantum spin liquid phases have beenfound in frustrated spin − and in frustrated spin − / . Similarly, Haldane phases havebeen shown in a spin − / and in a frustrated zigzag optical lattice of ultra-cold bosons . One of the questions that therefore arises iswhether frustration in a zigzag lattice of plane-polarizeddipoles leads to phases with non-trivial correlations be-tween lattice points.Wang et. al. have shown a rich phase diagram forthis system with the chain opening angle γ ≥ π/ γ . With the introduction of the NNN hopping, itbecomes impossible to do exact calculations for a sys-tem size large enough to exhibit many-body effects, wetherefore need a numerical approximation method. Weuse the Density Matrix Renormalization Group (DMRG)method because it is the most powerful numericalmethod to simulate one-dimensional systems . II. THE MODEL (a) Dipoles polarized at an angle θ in the plane of the zigzagchain(b) Spin − / Figure 1: (Color online) A zigzag chain of dipolesmapped to one of spin − / N = 100 sites but thefigure shows only seven sites labeled 1 through 7. Thehopping is allowed in a leg/direction (odd, even orNNN) of the chain only if the ends of the leg containopposite spins.Fig. 1 shows the spin − / | (cid:105) ≡ |↑(cid:105) , while an empty site is represented bya spin down, | (cid:105) ≡ |↓(cid:105) . With the constraint that doubleoccupancy is not allowed on any lattice sites, we map thisquasi-one-dimensional model of dipoles to a spin − / .Over the years, there has been a lot of work to studythe phase diagram of frustrated two-leg spin ladders usingvarious models, for instance, Refs. . As compared tothose, our model is simple because it is one-dimensional,has fewer degrees of freedom, and still exhibits frustra-tion. a r X i v : . [ c ond - m a t . s t r- e l ] J u l The Hamiltonian of the system is written as H = − J (cid:88) j ( S + j S − j +1 + h.c. ) − J (cid:88) j ( S + j S − j +2 + h.c. )+ V odd (cid:88) j = odd S zj S zj +1 + V even (cid:88) j = even S zj S zj +1 + V (cid:88) j S zj S zj +2 + h (cid:88) j S zj (1)where J > J > h is the magnetic field. The systemis half-filled, therefore the field term can be neglected.The spin operator S z is defined such that S z |↑(cid:105) = |↑(cid:105) and S z |↓(cid:105) = − |↓(cid:105) . V even and V odd are NN dipolarinteractions along even and odd legs of the chain re-spectively and V is the NNN dipolar interaction. Theinteractions are related to the dipole coupling strength (cid:15) dd = µ e / (4 π(cid:15) o | (cid:126)r − (cid:126)r | ), chain opening angle γ and po-larization angle θ as : V even = (cid:15) dd (cid:20) − (cid:18) π − γ − θ (cid:19)(cid:21) (2) V odd = (cid:15) dd (cid:20) − (cid:18) γ − θ (cid:19)(cid:21) (3) V = (cid:15) dd [2(1 − cos( γ ))] / (cid:20) − (cid:18) π − θ (cid:19)(cid:21) (4)where (cid:15) o and µ e are the vacuum permittivity and electricdipole moment, and (cid:126)r and (cid:126)r are the position of the twointeracting molecules.Before running any numerical simulations, we want toget an intuitive understanding of the model. We startwith some fundamental questions: Is there any regimewhere we can predict the ground state of the system andthen use numerics to validate our prediction? Can weidentify the frustrated and non-frustrated regimes andmap them to the physical parameter regime of γ and θ ?How are the NN and NNN hopping amplitudes relatedto one another and to γ and lattice depth? How differentdo the ground state phase diagrams look like for differ-ent lattice depths? As shown in Fig. 1, there are pairwiseinteractions in odd, even and NNN directions, each ofwhich can be attractive or repulsive. We will study theeffect of each interaction separately and put them to-gether afterwards to analyze their collective effect on thesystem.We write the Hamiltonian for any two interacting sites i and j , where j = i + 1 or i + 2, as H two-site-term = β (cid:32) −
12 ( S + i S − j + h.c. ) + αS zi S zj (cid:33) (5)where β = 2 J and α = V / J , and we refer to themas “relative” hopping and interaction strengths respec-tively. If we exactly solve this “two site term” in thebasis {|↑↑(cid:105) , |↑↓(cid:105) , |↓↑(cid:105) , |↓↓(cid:105)} , we will obtain the follow-ing result: Regardless of the value of β , the two sites prefer parallel alignment, ↑↑ or ↓↓ , represented by theletter “F” (for “ferromagnetic”) if the pairwise interac-tion α < − /
4, and antiparallel alignment, ↑↓ or ↓↑ ,represented by the letter “A” (for “antiferromagnetic”)if α > − /
4. It is worth noting that the critical value α c = − / H = (cid:88) j = odd β (cid:32) −
12 ( S + j S − j +1 + h.c. ) + α o S zj S zj +1 (cid:33) + (cid:88) j = even β (cid:32) −
12 ( S + j S − j +1 + h.c. ) + α e S zj S zj +1 (cid:33) + (cid:88) j β (cid:32) −
12 ( S + j S − j +2 + h.c. ) + α S zj S zj +2 (cid:33) (6)which is the sum of all the two-site terms in the threedirections, where β = 2 J , β = 2 J ,α o = V odd J , α e = V even J , α = V J . (7)The Hamiltonian written in this form helps us identifythe frustrated and non-frustrated regimes and predict theground state of the system prior to any simulations as wewill discuss in the next section.The relative hopping amplitudes β and β dependon the distance between interacting sites, chain open-ing angle, and lattice depth. If d and d are the lengthsof the odd (or even) and NNN legs respectively, then d = 2 d sin( γ/ β and β decrease exponentially with distance, we canshow that β β = exp (cid:20) − d λ (cid:16) γ/ − (cid:17)(cid:21) (8)where λ is a function of the lattice depth, and has theunits of length. Although d and λ can change when γ is varied, we can always set the ratio d /λ to a desiredvalue by tuning the lattice depth and thereby fixing λ independent of d or γ . The larger the value of the ratio d /λ , the deeper the lattice. Since γ, θ and d /λ can bevaried independently in real experiments, our model andall the results associated with it depend on these threeparameters.Throughout this paper, we use zero temperature, openboundary conditions, and (cid:15) dd = 1, and unless otherwisestated, d /λ = 0 .
1. In addition, we set β = 1, and withthis choice of β we allow the interactions to be muchstronger than the hopping.Fig. 2 shows how α o , α e and α depend on γ and θ while Fig. 3 illustrates how β varies with γ for differentlattice depths.Figure 2: (Color online) Mapping of the relative interaction strengths α o , α e and α to the physical parameterregime of the lattice, chain opening angle γ and polarization angle θ : Since β and α diverge as γ →
0, we take π/ γ . With π/ ≤ γ ≤ π and − π/ ≤ θ ≤ π/
2, we observe that both α o and α e varybetween − .
00 and 1 .
00, while α varies between − .
74 and 6 . π π π π π π γβ π π π π π π γβ d / λ = d / λ = d / λ = Figure 3: (Color online) Plot of β against γ for threedifferent values of d /λ . Since d /λ = 10 corresponds toa deep lattice, β increases much more exponentiallywith decreasing γ as compared to the other two valuesof d /λ . The inset shows a zoomed in plot for d /λ = 0 . d /λ = 1 which corresponds to a lattice ofintermediate depth as compared to the other two ratios.Before we proceed to the next section, we want to clar-ify that by setting the temperature to absolute zero wenullify thermal fluctuations. However, the experimentalrealization of this model would be a system at nanokelvintemperature with small but negligible thermal fluctua-tions. An example of such a system would be an ultracoldbosonic gas of Na Rb molecules that are stable againstchemical reaction in their absolute ground state , havea large permanent electric dipole moment (for instance,as large as 3.3 Debye ) which can lead to strong dipolarinteractions, and can be easily polarized by a moderateelectric field. For instance, a 5 kV cm − electric field caninduce a dipole moment larger than 2 Debye . As forthe zigzag optical lattice, which can be produced by usingthree laser beams as explained in Ref. , it would be nat-ural to set d ∼ − ), the dipolar coupling strength (cid:15) dd ≈ µ e / π(cid:15) o d ≈ . × − Joules.A natural energy scale for molecules in optical latticepotentials is the molecular recoil energy E r = (cid:126) k / m where m is the molecular mass. Since recoil energies (di-vided by the Plank constant h ) are of the order of sev-eral kilohertz , we estimate that E r /h ∼
10 kilohertz formolecular dipoles which means E r ≈ . × − Joules.With this estimate, we obtain (cid:15) dd ≈ . E r . By setting β = 1 and (cid:15) dd = 1, we are using (cid:15) dd as our energy scale sothat J = 0 . (cid:15) dd , a value that might be too small to probeexperimentally but could be increased by using smallerlattice constant (i.e., < > J , we canreadily see how the interaction strength in each of thethree directions scales with the corresponding hoppingstrength. For instance, when ( γ, θ ) = ( π/ , π/ | J /V even1 | = 0 . , | J /V odd1 | = 0 . | J /V | = 0 . III. FRUSTRATED AND NON-FRUSTRATEDREGIMES
As mentioned in the previous section, the pairwise in-teraction α in any direction is ferromagnetic or attrac-tive if α < − /
4, and antiferromagnetic or repulsive if α > − /
4. If we arrange the interactions in all the direc-tions based on whether they are attractive or repulsive,we find eight different combinations/regions as shown inFig. 4. Although this figure corresponds to the value of d /λ equal to 0 .
1, we get qualitatively similar plots forany other value of d /λ (see the Appendix); this impliesthat the phase diagrams should also be similar regardlessof the value of d /λ . Of the eight regions, four (AAA,AFF, FAF and FFA) are in the frustrated regime whilethe other four (FFF, AAF, AFA and FAA) are in thenon-frustrated regime.We will first explain and analyze non-frustrated regionsin the absence of hopping and then discuss the potentialscenario when the hopping is allowed. The simplest caseof a non-frustrated regime is the region FFF where thepairwise interactions in all the directions are ferromag-Figure 4: (Color online) Mapping of the frustrated andnon-frustrated regimes to the physical parameter regimeof chain opening angle γ and polarization angle θ .There are eight regions each with a unique color andlabeled with three letters which correspond, from left toright, to the odd, even and NNN directions respectively(Frustrated: AAA, AFF, FAF and FFA; non-frustrated:FFF, AAF, AFA and FAA). The black solid, bluedashed and red solid lines represent the contours for α o , α e and α respectively, each of which is equal to − / {| . . . ↑↑↑↓↓↓ . . . (cid:105) , | . . . ↓↓↓↑↑↑ . . . (cid:105)} would bethe exact ground states (from now on, the curly braces {} will represent states with the same energy). Anothernon-frustrated region is AAF where the pairwise interac-tions in the odd and even directions prefer antiferromag-netic (AFM) alignment while that in the NNN directionprefers FM alignment. In the absence of hopping, thetwo Neel states {|↑↓↑↓↑↓ . . . (cid:105) , |↓↑↓↑↓↑ . . . (cid:105)} are equallylikely configurations to have the lowest energy and there-fore, we expect the ground state to be AFM. Similarly,the ground state is expected to be a dimer of the type {|↑↑↓↓↑↑↓↓ . . . (cid:105) , |↓↓↑↑↓↓↑↑ . . . (cid:105)} in the non-frustrated re-gion FAA, and of the type {|↑↓↓↑↑↓↓ . . . (cid:105) , |↓↑↑↓↓↑↑ . . . (cid:105)} in the non-frustrated region AFA. In the presence of hop-ping, however, the four non-frustrated regions could fea-ture phases that become superfluid instead of solid, par-ticularly when the hopping dominates over the interac-tions.The four regions in the frustrated regime are poten-tially more interesting. The first such region is AFFwhere the pairwise interaction in the odd leg prefers AFMalignment while those in the even and NNN legs preferFM alignment. It is impossible for the spins to satisfy theinteractions in all directions simultaneously, and hence the system is frustrated. We can make similar argumentsto conclude that the other three regions FAF, FFA andAAA are also frustrated. As we will see later, there areregions in the frustrated regime where the pairwise inter-actions in the three directions are of similar strength andthus compete against one another. These regions requireparticular attention. IV. PHASE DIAGRAM
Fig. 5 shows the zero-temperature ground state phasediagram of the system for different values of γ and θ .This diagram has been produced with several DMRG tri-als each with a different initial state/condition, and themost appropriate ground state (the one with the low-est energy possible) has been considered. The differentphases, the order parameters and correlation functionsused to identify them, and the crossover between thosephases will be discussed in the subsequent paragraphs(see the Appendix for additional correlations). We la-bel the initial state as | init (cid:105) . We name the initial statewith spins randomly distributed in the lattice as “ran-dom initial state” and label it as | random (cid:105) . The letter“ E ” with a value attached to it will represent the energyof the ground state returned by a simulation. We will of-ten show ground states for two different initial states todemonstrate how the initial conditions affect the final re-sults obtained from DMRG simulations. When we showthe results for only one initial state, it means that thestate has led to the most appropriate ground state. Thecolor brightness for each phase represents the value of itsorder parameter while the black color represents the re-gion where all the order parameters vanish. We producethis phase diagram for the finite system size N = 100and we extrapolate the boundary between phases in thethermodynamic limit N → ∞ using finite-size scalinganalysis which we will discuss later. We find a sharp tran-sition between FM and AFM phases, and hence DRMGpinpoints the boundary between these two phases, whilewe find a smooth transition everywhere else as we willdiscuss later.It should be noted that the Hamiltonian Eq.(6) re-mains unchanged under the transformation θ → − θ (where α o and α e swap their values while α stays thesame). This implies that the phase diagram gives sim-ilar results in the range θ ∈ [ − π/ ,
0] as in the range θ ∈ [0 , π/ A. Dimerized phases
In the earlier section, we mentioned two distinct sets ofexpected ground states: {|↑↑↓↓↑↑↓↓ ... (cid:105) , |↓↓↑↑↓↓↑↑ ... (cid:105)} and {|↑↓↓↑↑↓↓ ... (cid:105) , |↓↑↑↓↓↑↑ ... (cid:105)} . We call this type ofdimer a “z-dimer” and although the non-frustrated re-gions FAA and AFA are the natural candidates for thisFigure 5: (Color online) Ground state phase diagram.These results depend on three independent parameters:chain opening angle γ , polarization angle θ , and theratio d /λ which we have set equal to 0 .
1. Each color isassociated with a different phase; the brighter a color,the deeper the system in that phase. The black colorcorresponds to the region where the order parametersvanish for all phases. AFM1 and AFM2 are bothantiferromagnetic phases labeled differently because ofthe nature of the ground state returned by DMRG. Thewhite curve labeled as “ α o + α e = − /
2” represents thephysical parameter regime where one of the pairwiseinteractions in NN directions is attractive while theother repulsive, and they both are the same distanceaway from their critical values α o,c = α e,c = − /
4. Thesuperfluid phase has been drawn using the values of thecorrelation function for the finite system size N = 100.All the other phases and their boundaries have beendrawn using the values of order parameters for theaforementioned system size. The white dots with verysmall error bars, obtained using finite-size scalinganalysis, represent the phase boundaries in thethermodynamic limit N → ∞ .phase, a frustrated region can also exhibit this type ofphase as shown in Fig. 6.Before discussing the other type of dimer that appearsin the phase diagram, let us define | + (cid:105) ≡ (1 / √ |↑↓(cid:105) + |↓↑(cid:105) ). Then a “xy-dimer” is simply the triplet boundstate | + (cid:105) ⊗ ... ⊗ | + (cid:105) or the one with free spins at theedges (often referred to as “dangling spins”) {|↑(cid:105) ⊗ | + (cid:105) ⊗ ... ⊗| + (cid:105)⊗|↓(cid:105) , |↓(cid:105)⊗| + (cid:105)⊗ ... ⊗| + (cid:105)⊗|↑(cid:105)} . The xy-dimer withdangling spins (or bound spins at the edges) is plausiblewhen the interaction in the even (or odd) direction ishighly repulsive while that in the other two directionsis weak as shown in Fig. 7. If the hopping amplitudeswere positive (i.e., J < J < .Spin liquid phases, which are phases with no magneticlong-range Neel order, are expected to be stable in sys-tems where quantum fluctuations can strongly suppress -1 0 1 0 50 100 〈 S j z S j + z 〉 j (a) ( γ, θ ) = ( π/ , π/ -1 0 1 0 50 100 〈 S j z S j + z 〉 j (b) ( γ, θ ) = ( π/ , π/ Figure 6: (Color online) Z-dimer phase. | init (cid:105) = |↓↑↓↑↓↑ . . . (cid:105) . The left plot shows a z-dimer inthe non-frustrated region FAA as expected. The rightplot shows a similar phase in the frustrated region FFAwhich clearly indicates that the attractive interaction inthe odd (or even) direction and the repulsive interactionin the NNN direction dominate over the attractiveinteraction in the third direction. -1 0 1 0 50 100 〈 S j z S j + z 〉 j (a) -1 0 1 0 50 100 〈 S j + S j + - 〉 j (b) Figure 7: (Color online) XY-dimer phase.( γ, θ ) = (5 π/ , . π ). Region: AAA. | init (cid:105) = | random (cid:105) . These two plots have been producedwith exactly the same initial condition. What we see isan example of a xy-dimer with dangling spins, whichmeans the repulsive interaction in the odd direction hasa dominating effect over that in the even and NNNdirections.magnetism, and these situations are found in low dimen-sions and in frustrated systems . Our model is com-prised of both. In the following paragraph, we explorethe possibility of such a phase.For a finite lattice, a xy-dimer phase with bound spinsat the edges is lower in energy than the one with danglingspins at the edges, and the system chooses as its groundstate the former or the latter depending on the values ofthe pairwise interactions. In the thermodynamic limit,however, the two phases would have the same energy.Therefore, one would expect the frustrated region thatresults in the xy-dimer phase to be an ideal candidate fora spin liquid phase when the interactions in the odd andeven directions are equally repulsive; this would allowthe ground state to be in the superposition of the two xy-dimer phases, a state similar to a resonating valence bond (see Ref. for a nice review of this state) but with thexy-singlets replaced with xy-triplets. In other words, aspin liquid phase may occur if the triplet bond connectingtwo adjacent sites can freely switch between odd and evendirections. The fact that the pairwise interactions in thetwo NN directions are always unequal in the xy-dimerregime of our model eliminates the possibility of a spinliquid phase.Similarly, because of the existence of triplet bonds, theregion in the phase diagram where a xy-dimer is observedis the only one where there could potentially be a Haldanephase. The existence of such a phase can be numericallyinvestigated using a string correlation function . Weconsider the one employed by Furukawa et. al. : O z str ( l, l + 2 r ) = − (cid:42) ( S zl + S zl +1 ) exp (cid:32) iπ l +2 r − (cid:88) m = l +2 S zm (cid:33) × ( S zl +2 r + S zl +2 r +1 ) (cid:43) (9)To explain how this correlation function is associatedwith a Haldane phase, we consider a pair of spins at adja-cent sites l + 2 j and l + 2 j + 1. If there were such a phase,the sum of the spins S zl +2 j + S zl +2 j +1 measured along thezigzag chain would alternate between +1 and − O z str ( l, l + 2 r )would detect this hidden order and take non-zero valuesas r becomes large. We calculate this correlation functionfor all j and r but we do not see a pattern as explainedbefore, and therefore we claim that we do not find a Hal-dane phase. And although we are unable to find one, wenote that Xu et. al. have shown the existence of sucha phase in an experimentally realizable spin-1 model ofbosons in a zigzag optical lattice. B. Superfluid phase
The reason that there are only small regions of su-perfluid (SF) phase in our phase diagram is that wechoose our parameters such that the interactions aremuch stronger than the hopping. Depending on the val-ues of β and β , there can be various regions of SF phase.The existence of this phase is confirmed by the polynomi-ally decaying long-range correlation (cid:104) S +1 S − j (cid:105) , knownas the “superfluid correlation”, as shown in Fig. 8 (see theAppendix for additional correlations).These two plots also show that that the two differ-ent frustrated regions AAA and FFA can feature thesame phase (SF in this case). It is worthwhile to lookat the values of the pairwise interactions for the leftplot: ( α o , α e , α ) = (0 . , . , . -5-2.5 0 0 2.5 5 l og e | 〈 S + S j - 〉 | log e j (a) ( γ, θ ) = (2 π/ , -3-1.5 0 0 2.5 5 l og e | 〈 S + S j - 〉 | log e j (b) ( γ, θ ) = ( π, . π ).Region: FFA. Figure 8: (Color online) SF phase. | init (cid:105) = | random (cid:105) .The two plots show the polynomially decayingsuperfluid correlation; the non-polynomial decay nearthe open ends of the chain is due to the edge effect. C. Ferromagnetic phase
Fig. 9 shows the ferromagnetic (FM) phase in this sys-tem. We show results subject to two different initial con-ditions in order to highlight the nature of the phase re-turned by DMRG. When the system is in the FM regime,the FM state with a single domain wall is the true groundstate because it has the lowest energy as compared to thestates produced with any other initial conditions. -1 0 1 0 50 100 〈 S j z S j + z 〉 j (a) | init (cid:105) = | random (cid:105) . E = − . -1 0 1 0 50 100 〈 S j z S j + z 〉 j (b) | init (cid:105) = | . . . ↓↓↓↑↑↑ . . . (cid:105) . E = − . Figure 9: (Color online) FM phase. ( γ, θ ) = ( π, π/ α o + α e = − /
2” represents the points where α o and α e are equally far away from their critical values α o,c = α e,c = − /
4, one being attractive while the otherrepulsive. So one would expect a FM phase on one sideof this line and an AFM phase on the other. Our results,however, show that the attractive interaction in the odd(or even) direction of the spin chain dominates over therepulsive interaction in the even (or odd) direction to acertain threshold, thus resulting in a FM phase on bothsides of this line. It should be noted that this line disap-pears when γ → . π because above this value of γ ,the system would be deep in the FM regime and there-fore, we do not obtain an AFM phase regardless of thevalue of θ . D. Antiferromagnetic phase
In Fig. 10 and Fig. 11, it can be seen that the accuracyof DMRG depends on the choice of initial state. Thereare obviously two different AFM regimes. We label thephase as “AFM1” when the NN correlations (cid:104) S zj S zj +1 (cid:105) are negative but greater than -1 for each site index j as shown in Fig. 10. A look at the values of the long-range correlation (cid:104) S z S zj (cid:105) (see the Appendix) confirmsthat this is an AFM phase. Simlarly, we label the phaseas “AFM2” when the system is deep in the AFM regimeso that (cid:104) S zj S zj +1 (cid:105) ≈ −
1. It is worth noting that althougha pure AFM phase is expected in the non-frustrated re-gion AAF, a simulation with a random initial state re-sults in a phase that has mostly AFM correlations butwith one or more clusters of identical spins, which we call“trapped regions”. It is clearly not a true phase but stillmakes sense from an experimental point of view, whichwe will explain later. -1 0 1 0 50 100 〈 S j z S j + z 〉 j (a) | init (cid:105) = | random (cid:105) . E = − . -1 0 1 0 50 100 〈 S j z S j + z 〉 j (b) | init (cid:105) = |↓↑↓↑↓↑ . . . (cid:105) . E = − . Figure 10: (Color online) AFM1 phase. ( γ, θ ) = ( π, E. Phase transitions and DMRG
Fig. 12 shows how initial states affect the ground stateenergy in DMRG simulations and why it is importantto perform multiple trials with various initial conditions.If we look at these results with reference to the phasediagram (Fig. 5), we can see that in the regime wherethe ground state is expected to be dimerized or AFM,the best choice for the initial state would be a z-dimer,a Neel state or a xy-dimer because these three statesresult in exactly the same ground state. Similarly, inthe regime where the ground state is expected to be FM, -1 0 1 0 50 100 〈 S j z S j + z 〉 j (a) | init (cid:105) = | random (cid:105) . E = − . -1 0 1 0 50 100 〈 S j z S j + z 〉 j (b) | init (cid:105) = |↓↑↓↑↓↑ . . . (cid:105) . E = − . Figure 11: (Color online) AFM2 phase.( γ, θ ) = ( π/ , π/ -250-150-50 0 π /4 π /2 E θ | init 〉 = ... ↓↓↓↑↑↑ ... 〉| init 〉 = ↓↑↓↑↓↑ ... 〉| init 〉 = ↑↓↓↑↑↓↓ ... 〉| init 〉 = ↑↑↓↓↑↑↓↓ ... 〉| init 〉 = random 〉| init 〉 = xydimer 〉 Figure 12: (Color online) Ground state energy of thesystem, E , plotted as a function of polarization angle, θ ,for γ = π/
3. The state | init (cid:105) has been used to denotethe “initial state” for a DMRG simulation, | random (cid:105) denotes the “random initial state” and | xydimer (cid:105) denotes the triplet bound state | + (cid:105) ⊗ . . . ⊗ | + (cid:105) . Thisfigure clearly shows that in a regime where a FM phaseis expected, only a simulation with a FM initial stateresults in a true ground state. It also shows that severalcurves meet at two points: θ = 0 . π , which belongsto a smooth crossover between z-dimer and FM phases(see Fig. 13), and θ = 0 . π , which lies at a sharpcrossover between FM and AFM phases (see Fig. 14).a simulation must start with a single domain wall FMstate.Simulations with various initial conditions clearly showthat there is a sharp transition between FM and AFMphases, and a smooth transition between z-dimer andFM phases and between SF and other phases (see theAppendix for detailed explanation of transition betweenSF and AFM phases). Experiments, however, can beexpected to confirm the unclear DMRG results in thefollowing way: Suppose we build a system from a sampleof randomly distributed spins and slowly cool it down sothat the spins restribute in the lattice to minimize theirenergy. If the sample consists of one or more trappedregions, the system must overcome an enormous energyhurdle to flip the spins in these regions, therefore thespin configuration would be expected to show signaturesof these trapped regions (as we saw earlier in Fig. 11)although it is not the lowest energy configuration. -1 0 1 0 50 100 〈 S j z S j + z 〉 j (a) | init (cid:105) = | . . . ↓↓↓↑↑↑ . . . (cid:105) . E = − . -1 0 1 0 50 100 〈 S j z S j + z 〉 j (b) | init (cid:105) = |↓↑↓↑↓↑ . . . (cid:105) . E = − . Figure 13: (Color online) Ground states and theirenergies subject to two initial conditions.( γ, θ ) = ( π/ , . π ). Region: FAF. -1 0 1 0 50 100 〈 S j z S j + z 〉 j (a) | init (cid:105) = | . . . ↓↓↓↑↑↑ . . . (cid:105) . E = − . -1 0 1 0 50 100 〈 S j z S j + z 〉 j (b) | init (cid:105) = |↑↑↓↓↑↑↓↓ . . . (cid:105) . E = − . Figure 14: (Color online) Ground states and theirenergies subject to two initial conditions.( γ, θ ) = ( π/ , . π ). Region: FAF. F. Order parameters
We define the order parameters for ferromagnetic, an-tiferromagnetic, z-dimer and xy-dimer phases as follows: O ferro = 4 N N (cid:88) i = N +1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) N +3 (cid:88) j = N (cid:104) S zi S zj (cid:105) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (10) O neel = 4 N N (cid:88) i = N +1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) N +3 (cid:88) j = N ( − j (cid:104) S zi S zj (cid:105) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (11) O zdimer = 2 N N (cid:88) i = N +1 (cid:12)(cid:12)(cid:12) (cid:104) S zi S zi +1 − S zi +1 S zi +2 (cid:105) (cid:12)(cid:12)(cid:12) (12) O xydimer = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) N +4 (cid:88) i = N +1 ( − i (cid:12)(cid:12)(cid:12) (cid:104) S zi S zi +1 (cid:105) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (13) π /4 π /2 o r de r pa r a m e t e r θ O ferro O neel O zdimer O xydimer (a) γ = π/ π /4 π /2 o r de r pa r a m e t e r θ O ferro O neel O zdimer O xydimer (b) γ = 7 π/ Figure 15: Order parameter for various phases as afunction of polarization angle θ Although we use correlation functions to explain how weidentify each phase, we use order parameters to find howdeep the system is in a given phase and also to find thecrossover between the phases. For the dimerized phases,we use the definitions given by Furukawa et. al. To min-imize the edge effects due to open boundaries, we use themethod employed by Rossini et. al. − we define the or-der parameters for ferromagnetic, antiferromagnetic andz-dimer phases as average expectation values of the cor-relators between spins in the middle part of the chain.For the xy-dimer phase, however, we only consider thethe correlations N/ (cid:104) S +1 S − j (cid:105) which, as mentioned earlier, de-cays polynomially in this phase. It should be noted thatwe have defined the order parameters such that they arealways non-negative.Fig. 15 shows how the order parameters for differentphases vary with polarization angle θ for a given value of γ . By definition, the order parameter for a given phaseshould vanish in all other phases and our results for FMand AFM phases are consistent with this. However, thedimerized phases consist of two flavors, xy and z, whichpair neighboring spins in different directions. Therefore,their order parameters overlap. The finite size scaling,which we will discuss later, along with the values of cor-relation functions allows us to find the boundary betweenthese two phases. G. Finite-size scaling and extrapolation
As mentioned earlier, the phase diagram (Fig. 5) hasbeen drawn using the values of order parameters and cor-relation functions for the finite system size N = 100. Weextrapolate the phase boundaries in the thermodynamiclimit N → ∞ using the finite-size scaling method ex-plored by Rossini et. al. We calculate the energy gap fordifferent system sizes N and find the value of θ for whichthe gap is minimum for each N , as shown in Fig. 16a.We call this value θ min . We then plot these θ min against1 /N and extrapolate the value of θ min when 1 /N → N . Our analysisshows that the energy gap scales polynomially with thesystem size near the boundary between z-dimer and fer-romagnetic phases as shown in Fig. 16b. We are unableto find the boundary between xy-dimer and superfluidphases. π /180 54 π /180 56 π /180 ∆ θ N = 100N = 200N = 300N = 400N = 500 (a) ∆ N ∆ = 32.4072 N -1.00003 (b) π /18055 π /18057 π /180 0 0.005 0.01 θ m i n (c) Figure 16: Finite-size scaling of the energy gap: γ = π/
6. (a) The energy gap ∆ is plotted as a functionof the polarization angle θ and different system sizes N .The gap is minimum at the phase transition point, wedenote the corresponding value of θ by θ min . (b) Theenergy gap is plotted as a function of the system size at θ = 0 . π which is near the phase transition point.The line of best fit is ∆ = 32 . N − . , whichimplies that the energy gap scales polynomially withthe system size. (c) By plotting θ min against 1 /N , weextrapolate the phase transition point in thethermodynamic limit N → ∞ as θ = (0 . ± . × − ) π .0 V. CONCLUSION
In conclusion, we have numerically studied the ground-state properties of a quasi-one-dimensional model thatcontains hopping and interactions up to second neigh-bors. Even though this is a rather simple model, it com-prises of frustrated regimes that lead to a rich phasediagram. We have used a novel approach to write theHamiltonian that gives an intuitive understanding of themodel, makes it convenient to identify frustrated andnon-frustrated regimes, and helps predict the groundstates beforehand so that the results obtained from nu-merical simulations can be verified. We have observedall the phases that Wang et. al. investigated. Nev-ertheless, in contrast to what was shown in their phasediagrams, we have observed a sharp transition betweenFM and AFM phases. We are, however, unable to findany spin liquid, Haldane or topological phase in this sys-tem. ACKNOWLEDGMENTS
We thank Q. Wang, H. Pichler, A. V. Balatsky and J.Javanainen for many helpful discussions. DMRG simu-lations were performed using the ITensor library . Weare extremely grateful to E. M. Stoudenmire, the leaddeveloper of ITensor, for helping us write the codes forour model and checking them for errors. This researchproject is supported by National Science Foundation. APPENDIXA. Frustrated and non-frustrated regimes
In Fig. 4, we saw how the eight regions - four frustrated(AFF, FAF, FFA and AAA) and four non-frustrated(FFF, AAF, AFA and FAA) - were related to the chainopening angle γ and polarization angle θ given the ratio d /λ = 0 .
1. Fig. 17 illustrates how these regions dependon the angles γ and θ for other lattice depths. We findthat all the eigght regions exist in our system, althoughtheir shape and size vary, regardless of the value of d /λ . (a) d /λ = 1(b) d /λ = 10 Figure 17: (Color online) Frustrated and non-frustratedregions for other lattice depths.
B. Correlation functions for various phases
In the body of this paper, we have shown the values ofonly one or two correlation functions to confirm a givenphase. In this section, we will show additional plots tosupport our claim. We will also include the values ofthe interactions to show which frustrated/non-frustratedregion the example point under consideration belongs to.1
1. Z-dimer phase
Fig. 18 shows additional plots for the z-dimer phaseshown in Fig. 6a, which belongs to the non-frustratedregion FAA. In principle, one should obtain (cid:104) S zj (cid:105) =0 for each site index j because the ground stateis expected to be a superposition of the two states {|↓↓↑↑↓↓↑↑ . . . (cid:105) , |↑↑↓↓↑↑↓↓ . . . (cid:105)} . However, DMRG re-turns one of these two states rather than a superposition.A similar argument is valid for all other phases.The other three plots are straightforward. We wouldexpect the same results regardless of whether the groundstate is a single z-dimer state, as is the result fromDMRG, or a superposition of two degenerate z-dimerstates, as is the result from ab-intio calculations. A sim-ilar argument is valid for all other phases. -1 0 1 0 50 100 〈 S j z 〉 j (a) -1 0 1 0 50 100 〈 S z S j z 〉 j (b) -1 0 1 0 50 100 〈 S j + S j + - 〉 j (c) -1 0 1 0 50 100 〈 S + S j - 〉 j (d) Figure 18: Additional correlations for the z-dimer phase.( γ, θ ) = ( π/ , π/ α o , α e , α ) = ( − . , . , . | init (cid:105) = |↓↑↓↑↓↑ . . . (cid:105) .
2. XY-dimer phase
Fig. 19 shows additional plots for the xy-dimer phaseshown in Fig. 7, which belongs to the frustrated regionAAA. -1 0 1 0 50 100 〈 S j z 〉 j (a) -1 0 1 0 50 100 〈 S z S j z 〉 j (b) -1 0 1 0 50 100 〈 S + S j - 〉 j (c) Figure 19: Additional correlations for the xy-dimerphase. ( γ, θ ) = (5 π/ , . π ).( α o , α e , α ) = (0 . , . , . | init (cid:105) = | random (cid:105) .
3. Superfluid phase
Fig. 20 shows additional plots for the SF phase shownin Fig. 8b, which belongs to the frustrated region FFA. -1 0 1 0 50 100 〈 S j z 〉 j (a) -1 0 1 0 50 100 〈 S z S j z 〉 j (b) -1 0 1 0 50 100 〈 S j z S j + z 〉 j (c) -1 0 1 0 50 100 〈 S j + S j + - 〉 j (d) Figure 20: Additional correlations for the SF phase.( γ, θ ) = ( π, . π ).( α o , α e , α ) = ( − . , − . , − . | init (cid:105) = | random (cid:105) .2
4. Ferromagnetic phase
Fig. 21 shows additional plots for the ferromagneticphase shown in Fig. 9, which belongs to the non-frustrated region FFF. -1 0 1 0 50 100 〈 S j z 〉 j (a) -1 0 1 0 50 100 〈 S z S j z 〉 j (b) -1 0 1 0 50 100 〈 S j + S j + - 〉 j (c) -1 0 1 0 50 100 〈 S + S j - 〉 j (d) Figure 21: Additional correlations for the FM phase.( γ, θ ) = ( π, π/ α o , α e , α ) = ( − . , − . , − . | init (cid:105) = | . . . ↓↓↓↑↑↑ . . . (cid:105) .
5. Antiferromagnetic phase: AFM1
Fig. 22 shows additional plots for the AFM1 phaseshown in Fig. 10, which belongs to the frustrated regionAAA. -1 0 1 0 50 100 〈 S j z 〉 j (a) -1 0 1 0 50 100 〈 S z S j z 〉 j (b) -1 0 1 0 50 100 〈 S j + S j + - 〉 j (c) -1 0 1 0 50 100 〈 S + S j - 〉 j (d) Figure 22: Additional correlations for the AFM1 phase.( γ, θ ) = ( π, α o , α e , α ) = (1 . , . , . | init (cid:105) = |↓↑↓↑↓↑ . . . (cid:105) .
6. Antiferromagnetic phase: AFM2
Fig. 23 shows additional plots for the AFM2 phaseshown in Fig. 11, which belongs to the non-frustrated re-gion AAF. -1 0 1 0 50 100 〈 S j z 〉 j (a) -1 0 1 0 50 100 〈 S z S j z 〉 j (b) -1 0 1 0 50 100 〈 S j + S j + - 〉 j (c) -1 0 1 0 50 100 〈 S + S j - 〉 j (d) Figure 23: Additional correlations for the AFM2 phase.( γ, θ ) = ( π/ , π/ α o , α e , α ) = (0 . , . , − . | init (cid:105) = |↓↑↓↑↓↑ . . . (cid:105) .3 C. Transition between antiferromagnetic andsuperfluid phases
In the phase diagram, it is hard to locate the exactboundary between AFM and SF phases for the finite sys-tem size N = 100. To understand the transition betweenthese two phases, we neglect the hopping and interactionin the NNN direction (i.e., we set β = 0 and α = 0.).We are interested in the situation where α o,e > − / β = 1 and forconvenience, we consider α o = α e .Fig. 24 and Fig. 25 show the various correlations for thecases α o,e = 0 . α o,e = 0 .
4. It is interesting to notethat the nature of the correlations (cid:104) S zj S zj +1 (cid:105) and (cid:104) S z S zj (cid:105) is not very different for the two cases; in fact, these corre-lations suggest the likelihood of an AFM phase. However,a SF phase in the former case is confirmed by the poly-nomial decay of the correlation (cid:104) S +1 S − j (cid:105) while an AFMphase in the latter is confirmed by the tendency of thespins to localize in lattice sites as indicated by the al-ternating sign for the values of the correlation (cid:104) S zj (cid:105) andthe exponential decay of the correlation (cid:104) S +1 S − j (cid:105) whichclearly indicates an insulating phase.Therefore, depending on the strength of β (with hop-ping and interaction between nearest-neighbors only),the system can be in a SF or AFM phase when thepairwise interactions in the odd and even directions pre-fer antiparallel alignment. We also notice that there isa smooth crossover somewhere between α o,e = 0 . α o,e = 0 .
4. Based on these results, it is safe to concludethat the nature of the transition between SF and AFMphases in the phase diagram (Fig. 5) is qualitatively thesame. -1 0 1 0 50 100 〈 S j z 〉 j (a) -1 0 1 0 50 100 〈 S z S j z 〉 j (b) -1 0 1 0 50 100 〈 S j z S j + z 〉 j (c) -10-5 0 0 2.5 5 l og e | 〈 S + S j - 〉 | log e j (d) Figure 24: Correlations for the SF phase with β = 1 , α o = α e = 0 . | init (cid:105) = | random (cid:105) . -1 0 1 0 50 100 〈 S j z 〉 j (a) -1 0 1 0 50 100 〈 S z S j z 〉 j (b) -1 0 1 0 50 100 〈 S j z S j + z 〉 j (c) -12-6 0 0 50 100 l og e | 〈 S + S j - 〉 | j (d) Figure 25: Correlations for the AFM phase with β = 1 , α o = α e = 0 . | init (cid:105) = | random (cid:105) . I. M. Georgescu, S. Ashhab, and F. Nori, Rev. Mod. Phys. , 153 (2014). C. Gross and I. Bloch, Science , 995 (2017). T. Lahaye, C. Menotti, L. Santos, M. Lewenstein, andT. Pfau, Reports on Progress in Physics , 126401 (2009). M. A. Baranov, M. Dalmonte, G. Pupillo, and P. Zoller,Chemical Reviews , 5012 (2012). K. R. A. Hazzard, M. van den Worm, M. Foss-Feig, S. R. Manmana, E. G. Dalla Torre, T. Pfau, M. Kastner, andA. M. Rey, Phys. Rev. A , 063622 (2014). F. L. Buessen, M. Hering, J. Reuther, and S. Trebst, Phys.Rev. Lett. , 057201 (2018). S. Yan, D. A. Huse, and S. R. White, Science , 1173(2011). S. Furukawa, M. Sato, S. Onoda, and A. Furusaki, Phys.Rev. B , 094417 (2012). S. Greschner, L. Santos, and T. Vekua, Phys. Rev. A ,033609 (2013). Q. Wang, J. Otterbach, and S. F. Yelin, Phys. Rev. A ,043615 (2017). S. R. White, Phys. Rev. Lett. , 2863 (1992). S. R. White, Phys. Rev. B , 10345 (1993). U. Schollw¨ock, Rev. Mod. Phys. , 259 (2005). K. A. Hallberg, Advances in Physics , 477 (2006). U. Schollw¨ock, Annals of Physics , 96 (2011). M. Azuma, Z. Hiroi, M. Takano, K. Ishida, and Y. Ki-taoka, Phys. Rev. Lett. , 3463 (1994). A. L¨auchli, G. Schmid, and M. Troyer, Phys. Rev. B ,100409 (2003). T. Tonegawa, K. Okamoto, T. Hikihara, and T. Sakai,Journal of Physics: Conference Series , 012003 (2017). G. Giri, D. Dey, M. Kumar, S. Ramasesha, and Z. G.Soos, Phys. Rev. B , 224408 (2017). S. Wessel, B. Normand, F. Mila, and A. Honecker, SciPostPhys. , 005 (2017). P. S. ˙Zuchowski and J. M. Hutson, Phys. Rev. A , 060703(2010). M. Aymar and O. Dulieu, The Journal of Chemical Physics , 204302 (2005). F. Wang, X. He, X. Li, B. Zhu, J. Chen, and D. Wang,New Journal of Physics , 035003 (2015). C. Becker, P. Soltan-Panahi, J. Kronj¨ager, S. D¨orscher,K. Bongs, and K. Sengstock, New Journal of Physics ,065025 (2010). M. Lemeshko, R. V. Krems, and H. Weimer, Phys. Rev.Lett. , 035301 (2012). A. Eckardt, Rev. Mod. Phys. , 011004 (2017). It is worth noting that the xy-dimer with dangling spinswill look similar to the valence bond solid state of theAKLT model if we replace the xy-triplets with xy-singlets;however, the absence of non-local correlation in the formermakes it strikingly different from the latter. H. T. Diep,
Frustrated Spin Systems (World Scientific Pub-lishing Company, 2005). Y. Zhou, K. Kanoda, and T. K. Ng, Rev. Mod. Phys. ,025003 (2017). M. den Nijs and K. Rommelse, Phys. Rev. B , 4709(1989). H. Tasaki, Phys. Rev. Lett. , 798 (1991). H. Watanabe, K. Nomura, and S. Takada, Journal of thePhysical Society of Japan , 2845 (1993). Y. Nishiyama, N. Hatano, and M. Suzuki, Journal of thePhysical Society of Japan , 1967 (1995). S. R. White, Phys. Rev. B , 52 (1996). E. H. Kim, G. F´ath, J. S´olyom, and D. J. Scalapino, Phys.Rev. B , 14965 (2000). J. Xu, Q. Gu, and E. J. Mueller, Phys. Rev. Lett. ,085301 (2018). D. Rossini and R. Fazio, New Journal of Physics , 065012(2012). B. Pandey, S. Sinha, and S. K. Pati, Phys. Rev. B ,214432 (2015).39