Phase diagram of the hexagonal lattice quantum dimer model: order parameters, ground-state energy, and gaps

1.2 MB·English

by  Thiago M. Schlittler

Phase diagram of the hexagonal lattice quantum dimer model:
order parameters, ground-state energy, and gaps
Thiago M. Schlittler,1 R´emy Mosseri,1 and Thomas Barthel2
1Laboratoire de Physique Th´eorique de la Mati`ere Condens´ee, UMR CNRS/Universit´e Pierre
et Marie Curie/Sorbonne Universit´es, 4 Place Jussieu, 75252 Paris Cedex 05, France
2Laboratoire de Physique Th´eorique et Mod`eles Statistiques,
Universit´e Paris-Sud, CNRS UMR 8626, 91405 Orsay Cedex, France
(Dated: December 24, 2014)
Thephasediagramofthequantumdimermodelonthehexagonal(honeycomb)latticeiscomputed
numerically, extending on earlier work by Moessner et al. The different ground state phases are
studied in detail using several local and global observables. In addition, we analyze imaginary-
time correlation functions to determine ground state energies as well as gaps to the first excited
states. This leads in particular to a confirmation that the intermediary so-called plaquette phase
is gapped – a point which was previously advocated with general arguments and some data for an
order parameter, but required a more direct proof. On the technical side, we describe an efficient
world-lineQuantumMonteCarloalgorithmwithimprovedclusterupdatesthatincreaseacceptance
5 probabilitiesbytakingaccountofpotentialtermsoftheHamiltonianduringtheclusterconstruction.
1 The Monte Carlo simulations are supplemented with variational computations.
0
2
PACSnumbers: 75.10.Jm,05.50.+q,05.30.-d,05.30.Rt

n
a
J I. INTRODUCTION the hexagonal lattice, these are ground-state configura-
9 tions of a classical Ising-spin model with antiferromag-
Interacting spin systems in two dimensions have been netic interactions on the (dual) triangular lattice, planar
] rhombus tilings, and height models [10, 11]. Topological
l widely studied over the last decades, both from experi-
e mental and theoretical points of view. Of importance in sectors can be characterized by so-called fluxes (see be-
- low). These sectors are invariant under the flips. The
r this context is the so-called resonating valence bond ap-
t topologicalpropertiesdependforexampleonthebound-
s proachputforwardbyP.W.Andersonin1973[1]inorder
. to analyze the physics of spin 1/2 Heisenberg antiferro- ary conditions and have consequences on the physics of
t
a magnets. Thishaslaterbeenadvocatedasawaytostudy the quantum dimer model.
m
the yet unsolved problem of high-temperature supercon- The quantum version, as proposed by Rokhsar and
- ductivity. Following Rokhsar and Kivelson [2], it proves Kivelson, corresponds to considering the set of all dimer
d
interesting, when studying the low energy properties of coveringsoftheclassicalproblemasanorthonormalbasis
n
o these phases, to consider a simpler model, called the spanning the Hilbert space. The Hamiltonian contains
c quantum dimer model (QDM). In the latter, the SU(2) kinetic terms that correspond precisely to the elemen-
[ singlet bonds are replaced by hard core dimers defined tary flips described above and an additional potential

on the edges of the lattice. Quantum dimer models have term, proportional to the number of flippable plaque-
1
v been employed to study for example superconductivity ttes. The competition between these kinetic and poten-
2 [2, 3], frustrated magnets [4–8], or hardcore bosons [9]. tial terms leads to a non-trivial phase diagram. For ex-
4 They can feature topological order, spin liquid phases, ample, when the potential term dominates in amplitude
2 and deconfined fractional excitations [7]. andisofnegativesign,thegroundstateisexpectedtobe
2 Before enlarging on the quantum systems, let us say a dominated by configurations which maximize the num-
0
few words about the classical case. Lattice dimer cover- ber of flippable plaquettes; for the opposite sign, one ex-
.
1 ings–thebasisstatesoftheHilbertspaceinthequantum pects a ground state dominated by dimer configurations
0 case – represent already a rich mathematical problem without flippable plaquettes. As will be discussed, such
5 with many connections to statistical physics problems. configurations exist and correspond to the so-called star
1
For a graph defined by its vertices and edges (defining andstaggeredphases,respectively. In-betweenthesetwo
:
v faces, often called plaquettes in the present context), a extremes, the phase diagram can display intermediary
Xi dimer covering is a decoration of the bonds, such that phases. The ground state is known exactly for the point
every vertex is reached by exactly one dimer. The sim- where kinetic and potential terms are of equal strength.
r
a plest rearrangement mechanism for dimer coverings is Thephysicsaroundthisso-calledRokhsar-Kivelson(RK)
provided by so-called plaquette flips. These are appli- point is expected to be different for bipartite and non-
cableforplaquettesaroundwhicheverysecondbondhas bipartite lattices [7].
a dimer and the flip amounts to exchanging covered and In this paper, we provide an extensive study of the
uncoveredbonds,yieldingadifferentvaliddimercovering quantumdimermodelonthebipartitehexagonal(honey-
(e.g., ←→ forahexagonallattice). Dimercoverings comb) lattice along the lines already followed by Moess-
arecloselyrelatedtootherconfigurationalproblems. For ner et al. [12]. In their seminal work, these authors nu-
2
A A A A
C C C C
B B B B
A A
C C
B B
(a) (b) (c)
Figure 1: Prototypes of quantum dimer ground states on a honeycomb lattice: (a) star phase, (b) plaquette phase, (c)
staggeredphase. Edgeswithahighprobabilityofcarryingadimerareindicatedinblack,andedgeswitha∼50%probability
areindicatedingray. The(dual)latticecanbedecomposedintothreetriangularsublatticesA,B,andC asshown. Inthestar
state, flippable plaquettes occupy two of the sublattices, while, in the plaquette phase, all plaquettes of one of the sublattices
are in a benzene-like resonance state (gray hexagons).
merically investigated the phase diagram by studying a II. QUANTUM DIMER MODEL
local order parameter which, in addition to the generic
RK transition point, shows a first order transition which A. Hilbert space and Hamiltonian
separatesthestarphasefromanintermediaryphase,the
so-calledplaquettephase. SeeFig.1forasketchofthese
Weconsiderthe2Dhexagonallatticeofspins-1/2with
twophases. Basedonthedataforthreedifferenttemper-
periodic boundary conditions. As described in the intro-
atures, Moessner et al. argued that the plaquette phase
duction, the quantum dimer models are defined on the
shouldbegapped–apointwhichconflictswithanearlier
subspace spanned by dimer configurations where every
√
analyticalanalysis[13]. InthepresentworkweuseQuan-
spin forms a singlet (|↑,↓(cid:105)−|↓,↑(cid:105))/ 2 with one of its
tum Monte Carlo simulations to extend the numerical
three nearest neighbors. These different dimer configu-
work by studying order parameters for different system
rations are used as an orthonormal Hilbert space basis.
sizes and temperatures as well as ground-state energies
Modelsofthistypeareforexampleimportantinthecon-
andexcitationgapswhichweobtainfromimaginary-time
text of resonating valence bond states and superconduc-
correlation functions. This leads to a clear confirmation
tivity[2,4,6]. Notethatdifferentdimercoveringsofthe
of the gapped nature of the plaquette phase. We shortly
lattice(dimerproductstates)arenotorthogonalwithre-
explain the reason for conflicting results of Ref. [13] and
spect to the conventional inner product for spin-1/2 sys-
supplement the Monte Carlo results with a variational
tems ((cid:104)σ |σ(cid:48)(cid:105)=δ ). However, as explained in Ref. [2],
σσ(cid:48)
treatment.
the two inner products can be related to one another
through additional longer-ranged terms in the Hamilto-
nian that turn out to be not essential. The Hamiltonian
Hˆ =−t(cid:88)(| (cid:105)(cid:104) |+h.c.)
QDM i i
i
The outline of this paper is the following. In sec- (cid:88)
+V (| (cid:105)(cid:104) |+| (cid:105)(cid:104) |) (1)
i i i i
tion II, the quantum dimer Hamiltonian is detailed and
i
the nature of the different phases is explained. In sec-
tion III, we describe the employed world-line Quantum containsakineticterm∝tthatflipsflippableplaquettes
Monte Carlo algorithm which is based on a mapping of (those with three dimers along the six plaquette edges)
the two-dimensional (2D) quantum model to a 3D clas- and a potential term ∝ V that counts the number of
sical problem, and which we accelerate through suitable flippable plaquettes. The sums in Eq. (1) run over all
cluster updates. Section IV introduces the employed ob- plaquettes i of the hexagonal lattice on a torus. The
servables. In section V, we present the results of the nu- potential term favors (V < 0) or disfavors (V > 0) flip-
mericalsimulations,andcharacterizethedifferentphases pable plaquettes. The only free parameter of this model
and phase transitions on the basis of different observ- ishencetheratioV/t. Inthefollowing,aplaquettecarry-
ables, ground-state energies, and energy gaps. Supple- ingj dimersiscalledaj-plaquettesuchthat3-plaquettes
mentary variational computations are described in sec- are the flippable ones.
tion VI. Section VII gives a summary of the results. De- The configuration space of the system is not simply
tailed discussions of some technical issues are delegated connected but consists of different topological sectors
to the appendices. whicharenotflip-connected. Eachsectorischaracterized
3
by two (flux) quantum numbers, also known as winding the Ising-type quantum model
numbers: Call A and B the two triangular sublattices
of the hexagonal lattice such that all nearest neighbors HˆQIM =Jz(cid:88)σˆizσˆjz−t(cid:88)σˆix+V (cid:88)δBˆi,0 (2)
of any site from A are in B. To compute the flux W (cid:104)i,j(cid:105) i i
through a cut C of the lattice, first orient all cut edges, on the triangular lattice, where {σˆx,σˆy,σˆz} denote the
say,fromAtoB,weightthemby+2or−1,dependingon i i i
Pauli spin matrices for lattice site i. The operator
whether they are covered by a dimer or not, and multi- Bˆ := (cid:80) σˆz, with N being the set of the six near-
plyeachweightby±1accordingtotheorientationofthe i j∈Ni j i
est neighbors of site i, yields for an {σˆz}-eigenstate the
edge with respect to C. The flux W is then computed by i
value zero, if exactly three of the six bonds starting at
summing the contributions of all cut edges. Such fluxes
site i are frustrated. A bond is called frustrated if the
W are invariant under plaquette flips. As fluxes through
corresponding two spins are parallel.
closed contractible curves C are zero, one has two flux
Atthecenterofeachtriangleliesavertexofthehexag-
quantumnumbersW andW ,correspondingtothetwo
x y onal lattice. For a given dimer covering, one dimer is
topologically distinct closed non-contractible curves on
shared by this vertex and the dimer crosses exactly one
the torus. Notice that these two fluxes characterize an
of the three edges of the triangle at an angle of 90°. See
averageslopeintheheightrepresentation[10]ofthesys-
Fig. 2. For sufficiently strong J , the physics of the
tem. z
quantum Ising model (2) is restricted to the subspace
Let us briefly recall the phase diagram obtained in
spanned by the classical ground states. Those have ex-
Ref. [12]. Three phases belonging to two different topo-
actly one frustrated bond per triangle (all other configu-
logical sectors have been described. The ground states
rations have higher energy). The identification of dimer
for the so-called star phase (−∞ < V/t < (V/t) ) and
C basisstatesandIsingbasisstatesisthenstraightforward.
the plaquette phase ((V/t) < V/t < 1) are found in
C Given a certaindimer configuration, put aspin up onan
the zero flux sector, while the staggered phase ground
arbitrary site. Associating frustrated Ising bonds with
states (1<V/t<∞) are in the highest flux sector. See
those that are crossed by a dimer in the given state, we
Fig. 1. The ground states in the zero flux sector can be
canworkinward-out,assigningfurtherIsingspinstillthe
distinguished using sublattice dimer densities. For that
triangular lattice is filled. The state, up or down, for a
purpose, we recall that the plaquettes of the hexagonal
newsitedependsonthespinstateofanalreadyassigned
lattice can be separated into three subsets – triangular
neighboring site and on whether the corresponding bond
sublatticesA,B andC ofdisjointplaquettes,asdepicted
is frustrated or not.
in Fig. 1, such that every hexagon of a set shares bonds
This mapping of dimer configurations on the hexago-
with three hexagons of the two other sets each.
nal lattice to spin-1/2 configurations on the triangular
lattice implies that, for the quantum Ising model, we
employ the conventional inner product (cid:104)σ |σ(cid:48)(cid:105) = δ
σσ(cid:48)
III. QUANTUM-CLASSICAL MAPPING AND that makes different {σˆz}-eigenstates orthonormal. In
MONTE CARLO SIMULATION i
the Hamiltonian (2), the spin-flip terms ∝ t correspond
tothekineticterminthequantumdimermodel(1). Due
AsdonebyMoessner,Sondhi,andChandra[5,12,14], totheenergeticconstraintimposedbyJ →∞,theyare
z
the 2D quantum dimer model on a hexagonal lattice can onlyeffectiveforsiteswherethespinflipdoesnotchange
be studied by mapping it first to a 2D quantum Ising the number of of frustrated bonds, corresponding to the
model on the (dual) triangular lattice. The resulting flippable plaquettes in the dimer model. The term ∝ V
Ising-typequantummodelcanbestudiedefficientlyusing corresponds exactly to the potential term in the dimer
world-line Quantum Monte Carlo [15] by approximating model.
its partition function and observables by those of a clas-
sical 3D Ising-type model (CIM) on a stack of triangu-
lar2Dlattices(quantum-classicalmapping)asdescribed
in the following sub-sections. We accelerate the Monte
Carlo simulation of the classical 3D model through suit-
Flip
able cluster updates.
A. Equivalence to a quantum Ising model on the
dual lattice
As shown in Fig. 2, the dual of the hexagonal lattice Figure 2: Equivalence of dimer coverings of the hexagonal
is the triangular lattice whose vertices are located at the lattice and Ising-spin configurations on the (dual) triangular
hexagoncenters. Weassignaspin-1/2(σ =±1)toeach lattice. Every dimer corresponds to a frustrated bond (↑−↑
i
of the vertices and, as explained in the following, the or ↓ − ↓). Flipping a plaquette in the hexagonal lattice is
quantum dimer model (1) maps for the limit J →∞ to equivalent to flipping a spin in the dual lattice.
z
4
The equivalence of the two models is slightly compli- As shown in appendix B, the parameters Kz and Kτ
cated by two issues. First, as we are free to choose the of the classical Ising model (3) are given by
orientationofthefirstassignedspin,agivendimerconfig-
urationcorrespondstotwospinconfigurationsthatdiffer Kz =∆β and e−2Kτ =tanh(∆βt). (6)
by a global spin-flip. Second, periodic boundary condi-
tions correspond, for certain topological sectors of dimer
configurations, to anti-periodic boundary conditions in C. Monte Carlo algorithm with cluster updates
the Ising model.
The representation (5c) of expectation values (5b) of
quantum observables as expectation values of classical
B. Approximation by a classical 3D Ising model observables is of great value, as it can be evaluated effi-
cientlywithaMonteCarloalgorithmbysamplingclassi-
To apply world-line Quantum Monte Carlo [15], we cal states σ. Specifically, one generates a Markov chain
can approximate the partition function and observables of classical states σ with probabilities e−ECIM(σ)/ZCIM
of the quantum Ising model (2) on the 2D triangular and averages O(σ) over these states.
latticebythoseofa3DclassicalIsingmodelonastackof The most simple scheme would be to choose in every
2DtriangularlatticesbyaTrotter-Suzukidecomposition iteration of the algorithm one of the flippable spins (a
[16, 17]. To this purpose, we separate the Hamiltonian spin on site j of time slice n is flippable, iff (cid:80) σn =
(2) into two parts 0), compute the energy difference E (σ(cid:48))−i∈ENj i(σ)
CIM CIM
thattheflippingofthespinwouldcause, andflipitwith
HˆQIM =Hˆz+Hˆx with Hˆx :=−t(cid:88)σˆix and a probability that is given by the so-called Metropolis
i rule as detailed in appendix C.
Hˆz :=Hz({σˆz}):=J (cid:88)σˆzσˆz+V (cid:88)δ . However, as one increases the accuracy by reducing
i z i j Bˆi,0 ∆β (for a fixed inverse temperature β =N∆β), the cou-
(cid:104)i,j(cid:105) i
pling Kτ between the time slices increases with Kτ ∝
AsdetailedinappendixB,onecanusetheTrotter-Suzuki log(1/∆β) and the classical Ising model, hence, becomes
decomposition stiff with respect to the time direction. In the generated
states σ, there will occur larger and larger 1D clusters
(cid:16) (cid:17)N
e−βHˆQIM = e−∆2βHˆze−∆βHˆxe−∆2βHˆz +O(∆β3) of spins along the time-direction that have the same ori-
entation, σm = σm+1 = ··· = σm+n. Flipping one of
i i i
of the density operator with imaginary-time step ∆β ≡ the spins inside such a cluster becomes less and less fre-
β/N to determine the parameters Kz and Kτ for the quent as the associated energy change increases with the
classical Ising model increasing coupling Kτ. This would result in an ineffi-
cient Monte Carlo sampling with high rejection rates for
(cid:88) (cid:88)
E (σ)=Kz Hz(σn)−Kτ σnσn+1 (3) spin-flips. We avoid this effect by doing 1D cluster up-
CIM i i
dates instead of single-spin updates: In every iteration
n n,i
of the algorithm an initial flippable spin is selected and,
suchthatthepartitionfunctionsZQIM ≡Tre−βHˆQIM and in an intermediate phase, a 1D cluster is grown in the
ZCIM =(cid:80)σe−ECIM(σ) of the two models coincide (up to icmluasgteinraarsy-atimwehodleir.ecWtioenfubrethfoerredseucgregaessetinrgejetcotioflniprattheiss
a known constant A),
by taking account of changes in the number of flippable
Z =A·Z +O(∆β3) (4a) spinsduringtheclusterconstruction. SeeappendixCfor
QIM CIM
details.
with A=[sinh(2∆βt)/2]LN/2, (4b)
Besides computing in this way expectation values of
diagonal operators Oˆ = O({σˆz}), one can also evalu-
aswellasexpectationvaluesofobservablesOˆ =O({σˆz}) i
i ate expectation values of non-diagonal operators such as
that are diagonal in the in the {σˆz}-eigenbasis,
i certain correlation functions or, for example, the energy
(cid:104)Oˆ(cid:105) =(cid:104)O(cid:105) +O(∆β3), where (5a) expectation value as described in appendix D.
QIM CIM
1
(cid:104)Oˆ(cid:105) ≡ Tr(e−βHˆOˆ) and (5b)
QIM ZQIM IV. STUDIED OBSERVABLES
1 (cid:88)
(cid:104)O(cid:105) ≡ e−ECIM(σ)O(σn) ∀ . (5c)
CIM Z n In the next section, section V, we numerically charac-
CIM
σ
terize the phase diagram of the quantum dimer model
In these equations, σ = (σn|n = 1,...,N) is a vector using several observables: the magnetization of the as-
of classical spin configurations σn = (σn|i ∈ T) on the sociated Ising model, dimer densities, the ground-state
i
triangular lattice T for each of the imaginary-time slices energy, and the energy gap to the first excited state. Let
n=1,...,N, and L is the number of lattice sites i in T. us briefly describe them in the following.
5
A. Magnetization of this possible problem in analyzing the data. Specifi-
callyforthesublatticedimerdensities,duringthecourse
We compute the root mean square (RMS) magnetiza- of the Monte Carlo simulation, we have rearranged the
tion (cid:104)mˆ2(cid:105)1/2 := (cid:104)((cid:80) σˆz/L)2(cid:105)1/2 for the quantum Ising three sublattice labels according to the dimer occupan-
i i QIM cies, instead of keeping the ordering constant.
model (2) with L sites as an order parameter to distin-
guishthestarandplaquettephasesandtofacilitatecom-
parison to earlier work [12].
D. Ground state energy
B. Local and global dimer observables To study the phase diagram, it is certainly of high
interest to access the ground-state energy which directly
decideswhatphaseprevailsforgivenvaluesoftheHamil-
The simulations give access to dimer densities (cid:104)nˆ (cid:105),
i
tonian parameters. For sufficiently low temperatures
the average number of dimers on plaquette i. Two-
in the simulation, the expectation value (cid:104)Hˆ (cid:105) of the
dimensional (contrast) plots of these densities nicely il- QIM
quantum Ising model Hamiltonian corresponds to the
lustrate the ground-state structure, even when its peri-
ground-state energy. But Hˆ is not a diagonal op-
odicity has to be understood in a probabilistic way. QIM
erator, and hence Eq. (5) cannot be used. As detailed in
We also evaluate the normalized total numbers of j-
appendixD,itcanneverthelessbeevaluatedonthebasis
plaquettes ((cid:104)ρˆ (cid:105),(cid:104)ρˆ (cid:105),(cid:104)ρˆ (cid:105),(cid:104)ρˆ (cid:105)). Specifically, with j-
0 1 2 3 of imaginary-time correlators (cid:104)σnσn+1(cid:105) .
plaquettes being the plaquettes carrying j dimers and i i CIM
(cid:80)
L being the system size, ρˆ ≡ δ /L. As described
j i nˆi,j
inappendixA,theplaquettenumbers(cid:104)ρˆ (cid:105)obeythesum
j
rule E. Energy gap
(cid:104)ρˆ3(cid:105)−(cid:104)ρˆ1(cid:105)−2(cid:104)ρˆ0(cid:105)=0. (7) It is important to determine whether a given phase
has gapless excitations or not. As explained in ap-
Notice that (cid:104)ρˆ (cid:105) does not enter in the sum rule, while
2 pendix E, we can estimate the energy gap to the first
changes in the number of 3-plaquettes, which enter both
excited state by fitting imaginary-time correlation func-
thekineticandpotentialenergyterms,mustbecompen- tions (cid:104)Aˆ(0)Aˆ†(iτ)(cid:105). In the classical Ising model (3) they
sated by plaquettes with zero dimers or one dimer. (cid:104)ρˆ (cid:105)
2 correspond to inter-layer correlators with layer distance
is nevertheless constrained by the fact the total number
of plaquettes is of course constant, i.e., (cid:80)3 (cid:104)ρˆ (cid:105)=1.
j=0 j
C. Sublattice dimer densities Star Plaquette Stagerred
RK point
AsdescribedaboveandindicatedinFig.1,thehexago-
nalplaquettescanbeseparatedintothreesets(A,B,C), 0.35
eachformingatriangularlattice,suchthateveryhexagon
0.3
in a set shares a bond with three hexagons of the two
other sublattices each. The “prototype” states of the
0.25
star and the plaquette phases (Fig. 1) can be charac-
terized qualitatively in terms of dimer densities in the 1/2 0.2
three sublattices. To this purpose, we can analyze av- >
2
eraged dimer densities on each sublattice and call them m 0.15 L
(cid:104)nˆA,B,C(cid:105), i.e., nˆA ≡ L3 (cid:80)i∈Anˆi etc. < 14346
Itshouldbestressedthatthesystemsunderstudymay 0.1 324
576
have degenerate (or nearly degenerate) ground states.
1296
0.05
The star crystal (ground state for V/t = −∞) and the 3600
ideal plaquette state (not a ground state, see below) are Staggered
0
both 3-fold degenerate. For sufficiently large systems, -3 -2 -1 0 1 2
it is expected that this symmetry is kinetically broken
V/t
in the Monte Carlo simulation. However, one cannot
fully prevent the system from translating from one typ- Figure 3: The root-mean square magnetization (cid:104)mˆ2(cid:105)1/2, as
ical ground-state configuration to another (even at the defined in Section IVA, for the quantum dimer model. The
levelofmedium-sizepatches), smearingouttheinforma- different curves correspond to different system sizes L and
tion carried by these local parameters. This possibility are obtained from Monte-Carlo simulations for V/t≤1 with
β =19.2and∆β =0.02. ForallV/t>1,thestaggeredstate,
was minimized here by choosing large system sizes and
depictedinFig.1c,isthegroundstateandhence(cid:104)mˆ2(cid:105)1/2 =0.
low temperatures. We nevertheless carefully kept track
6
0.16
L
0.14 36 4.8
144 9.6
0.12 324 19.2
V/t
576 38.4
1/2 0.1 13269060 -0-.02.73
2> 0.08 6561 -0.24
m 11664 -0.23
< 0.06 -0.226
-0.22
0.04 -0.21
-0.15
0.02 -0.1
0
0
-0.3 -0.25 -0.2 -0.15 -0.1 -0.05 0 0 0.03 0.06 0.09 0.12 0.15 0.18 -0.3-0.25-0.2-0.15-0.1-0.05 0 0.05 0.1
V/t 1/l V/t
(a) (b) (c)
Figure4: Locatingthetransitionbetweenthestarandplaquettephases. (a)Theroot-meansquaremagnetization(cid:104)mˆ2(cid:105)1/2 for
different lattice sizes L as a function of V/t for β = 19.2 and ∆β = 0.02. (b) For the same temperature, (cid:104)mˆ2(cid:105)1/2 is plotted
for different values of V/t as a function of the inverse linear system size 1/(cid:96). (c) (cid:104)mˆ2(cid:105)1/2 as a function of V/t for different
temperatures, ∆β =0.02, and L=81×81.
∆n = τ/∆β. For sufficiently low temperatures, and τ third sublattice (C) are dimer-free.
and β−τ big compared to the gap to the second excited
(cid:79) (cid:79)
state,theleadingtermsinthecorrelationfunctionareof |ψstar(cid:105)= | i(cid:105) | j(cid:105) (8)
the form a+b·cosh((β/2−τ)∆E), allowing to fit the i∈A j∈B
upper bound ∆E of the gap.
Changing from the dimer to the Ising-spin represen-
tation, A and B carry spins of equal orientation, say
σ =+1, and all spins on sublattice C have the oppo-
A,B
V. SIMULATION RESULTS siteorientation(σ =−1)suchthattheRMSmagnetiza-
C
tion reaches its maximum possible value (cid:104)mˆ2(cid:105)1/2 = 1/3.
In the following, let us study in detail the phase dia- It will decrease as V/t is increased – a behavior which is
gram of the quantum dimer model, starting from large clearlyseeninFig.3. Noticethat,forV/t=−3,(cid:104)mˆ2(cid:105)1/2
negative V/t, i.e., in the star phase. The observables is still very close to the maximum value 1/3.
described in the previous section are evaluated in simu- TounderstandhowincreasingV/taffectstheidealstar
lations for patches of linear size (cid:96) with a 60◦ rhombus state, one can do perturbation theory in t/V. The cal-
shape,periodicboundaryconditions,andL=(cid:96)2 plaque- culation, done up to second order in t/V, is given in
ttes. In order to be able to separate the lattice into the appendix F. The result for the ground-state energy is
three sublattices A, B, and C (Section IVC), (cid:96) needs to plotted in Fig. 8. It compares well with the simulation
be a multiple of three. results up to V/t∼−1. The first correction to the ideal
star state amounts to mixing in configurations with one
flipped plaquette.
The ground state for small negative t/V is the ideal
A. The star phase (−∞<V/t<(V/t) )
C star state dressed with flipped plaquettes in both A and
B,and,atsomepointalsoinC,whenthreeormoreflips
This phase has previously been called the “columnar haveoccurredlocally. Theseschangesinthegroundsate
phase”, in analogy with a corresponding phase of the canbequantifiedbythenumbersofj-plaquettesasdone
square lattice quantum dimer model, where dimers are in Fig. 5. In the ideal star state (V/t → −∞), one has
aligned along columns. For the hexagonal lattice, this ((cid:104)ρˆ (cid:105),(cid:104)ρˆ (cid:105),(cid:104)ρˆ (cid:105),(cid:104)ρˆ (cid:105)) = (1/3,0,0,2/3). Say, sublattices
0 1 2 3
denomination is a bit misleading, and we follow Ref. [18] A and B contain the flippable plaquettes in this limit.
in calling it the “star phase” [27]. For large negative AfterflippingaplaquetteinA,thethreeneighboringpla-
V, the potential term dominates the kinetic term and quettes in sublattice B carry two instead of three dimers
the ground state is dominated by dimer configurations and the three neighboring plaquettes in C are no more
that maximize the number of flippable plaquettes. In dimer free, but carry one dimer each. The numbers of 0-
the limit V → −∞, there are three degenerate ground and 3-plaquettes are hence reduced by three and those
states given by ideal star states as depicted in Fig. 1a, of 1- and 2-plaquettes are increased by 3. This explains
where all plaquettes from two of the three sublattices, why the curves for (cid:104)ρˆ (cid:105) and (cid:104)ρˆ (cid:105) in Fig. 5 are almost
0 3
say A and B, are flippable, while all plaquettes of the parallel up to V/t∼−1 and why those for (cid:104)ρˆ (cid:105) and (cid:104)ρˆ (cid:105)
1 2
7
our simulations
The transition can also be observed in the dimer ob-
0
1 servables. The normalized j-plaquette numbers (cid:104)ρˆ (cid:105) all
1 j
2 show a small discontinuity at (V/t)C (Fig. 5). The dis-
3 continuityof(cid:104)ρˆ (cid:105)(seeFig7b)atteststhefirstorderchar-
3
acter of the transition, since (cid:104)ρˆ (cid:105) is the derivative of the
0.8 3
energy with respect to V. But the discontinuity is quite
small leading to a barely visible slope change for the en-
ergy (Fig. 7a).
0.6
> At least as spectacular as the magnetization drop is
i
< the sudden shift in sublattices dimer densities seen in
Fig. 6b. It nicely agrees with the qualitative properties
0.4 oftheidealplaquettestate,depictedinFig1b,wherethe
resonating 3-plaquettes are located on one of the three
sublattices.
0.2
C. The plaquette phase ((V/t) <V/t<1)
C
0
-3 -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 The plaquette phase is certainly more complex to de-
scribe than the star phase. The features of the ground
V/t
state in this phase are to some extent captured by the
“ideal”resonatingplaquettestate,asimpleproductstate
Figure 5: Normalized numbers of j-plaquettes, (cid:104)ρˆ (cid:105), for the such that all plaquettes of one of the three sublattices,
j
zero flux sector, system size L = 81×81, β = 19.2, and say A, are in a benzene-type resonance state such that
∆β =0.02. Around(V/t) ,afinergridofpointswasusedto
C √
(cid:79)
resolve the jumps in the densities at the transition. In that |ψ (cid:105) = (| (cid:105)+| (cid:105))/ 2, (9)
plaq dimer i i
region,datapointsarenotmarkedbysymbols. Althoughthe
i∈A
global ground state is not in the zero flux sector for V/t>1,
data obtained for the zero flux sector is also shown for that This would in principle lead to a three-fold degeneracy
region and is discussed in the text. (theresonating3-plaquettescouldjustaswellbelocated
onsublatticesBorC). Incontrasttothestarphasecase,
the ideal plaquette state is not an exact ground state for
increase correspondingly and are on top of each other.
any V/t. The actual ground states can be viewed as
The contrast plots of the dimer density (cid:104)nˆ (cid:105) in Fig. 6a
i the ideal plaquette state, dressed by a variable amount
andtheplotofthesublatticedimerdensities(cid:104)nˆ (cid:105)in
A,B,C of flip excitations in the other two plaquette sublattices.
Fig.6bshowthatthedifferencesbetweendimerdensities
For convenience, we use the Ising-spin notation. In this
on sublattices A and B on one hand and those on sub-
representation, the ideal plaquette state |ψ (cid:105) can be
plaq
lattice C on the other hand are reduced before reaching
denoted as
a critical point (V/t) .
C
(cid:79) (cid:79) (cid:79)
|ψ (cid:105)= |→(cid:105) |↑(cid:105) |↓(cid:105) , (10)
plaq i j k
i∈A j∈B k∈C
B. The star to plaquette phase transition at
√
(V/t)C =−0.228±0.002 where |→(cid:105) denotes the σˆx-eigenstate (|↑(cid:105) +|↓(cid:105) )/ 2.
i i i i
The spins in sublattices B and C must be anti-parallel
A first order transition occurring between the star with respect to each other. In accordance with the nu-
phase and the the so-called plaquette phase is found at mericalresults,theRMSmagnetization(cid:104)mˆ2(cid:105)1/2alsovan-
(V/t)C = −0.228±0.002. This critical value is consis- ishes for the ideal plaquette state |ψplaq(cid:105) in the thermo-
tent, but more precise than that given in Ref. [12]. At dynamical limit. As (cid:80) σˆz|ψ (cid:105)=(cid:80) σˆz|ψ (cid:105), we
(V/t)C, the RMS magnetization (cid:104)mˆ2(cid:105)1/2 suddenly drops have that i i plaq i∈A i plaq
to a much smaller value which goes to zero in the ther-
modynamic limit. Fig. 3 displays the RMS magnetiza- (cid:16)(cid:88) (cid:17)2
(cid:104)ψ |mˆ2|ψ (cid:105)=(cid:104)ψ | σˆz |ψ (cid:105)/L2
tion for the whole phase diagram, while Fig. 4a provides plaq plaq plaq i plaq
a zoom close to the transition. We determined (V/t) i∈A
C
by plotting (cid:104)mˆ2(cid:105)1/2 as a function of the inverse (linear) =(cid:88)(cid:104)ψ |(σˆz)2|ψ (cid:105)/L2 = 1 →0. (11)
size of system (Fig. 4b). The temperature dependence plaq i plaq 3L
i∈A
shown in Fig. 4c indicates that using a larger β should
not substantially modify the numerical results. Consid- Theenergydensityfor|ψ (cid:105)canbecomputedeasilyand
plaq
ering this, we set β = 19.2 and ∆β = 0.02 for most of yields an upper bound to the exact ground state energy.
8
3
2.5
> 2
C
n
<
>,
B 1.5
n
<
>,
A
n
< 1
0.5 L = 144
L = 3600, Sub. A
Sub. B
Sub. C
0
-1 -0.5 0 0.5 1
V/t
(a) (b)
Figure 6: (a) Local dimer density (cid:104)nˆ (cid:105) for different values of V/t with L=60×60 plaquettes, β =19.2, and ∆β =0.02. (b)
i
Sublatticedimerdensities(cid:104)nˆ (cid:105)asfunctionsofV/tforL=60×60(solidlines)andL=12×12plaquettes(dashedlines),
A,B,C
respectively.
See section VI. At V =0, it takes for example the value lowering the temperature, observables should converge,
−1/3 which is clearly above the numerically determined once the temperature is sufficiently below the gap. This
value of ≈−0.38 (Fig. 7a). One can improve |ψ (cid:105) as a is confirmed in Fig. 4c.
plaq
variational state by adding flip excitations in sublattices Incontrast, anearlieranalyticaltreatmentinRef.[13]
B and C. This is possible due the fact that 3-plaquettes suggests that the plaquette phase should be gapless. We
occur in B and C with density 1/8. believe that this is due to a mistake in that derivation.
In Ref. [13], the model for V = 0 and a hexagonal lat-
A finite energy gap for the plaquette phase was advo-
tice with mixed boundary conditions is mapped to a
cated in Ref. [4] with an indirect numerical confirmation
model of vertically fluctuating non-intersecting strings
based on the magnetization for three different tempera-
on a square lattice. First, one can obtain the ground
tures. As a matter of fact, it is possible to estimate ex-
state of a single string which corresponds to the ground
citation gaps on the basis of imaginary-time correlation
state of the XX chain (energy E(1) → −2(cid:96)t/π) and
functions. The computation is described in appendix E. 0
that of the quantum dimer model in a high-flux sec-
OurresultsarepresentedinFig9: Startingfromthestar
tor. One can now add further strings, each reducing the
phase, the gap estimate decreases in a marked behavior
flux by one. To construct an N-string ground state, in
around the first order phase transition at (V/t) . Then,
C
Ref. [13], the product of vertically shifted single-string
it increases again in the plaquette phase, and eventually
ground state wavefunctions is considered. To take ac-
goes to zero as we approach the RK point at V/t = 1.
count of the no-intersection constraint for the strings,
This is clear evidence for a finite gap in the plaque-
this wavefunction is anti-symmetrized with respect to
tte phase with a maximal value of roughly 0.7t around
the string positions, first with respect to all variables
V/t≈0.1. InthevicinityoftheRKpoint,thecurvesfor
gap estimates in Fig 9, computed at different β, are not y1(n), then with respect to all y2(n), etc., where (x,yx(n))
convergedwithrespecttothetemperatureanymore. The are the coordinates of string n. In analogy to the
reason is simply that, as the gap vanishes, temperatures anti-symmetrization for fermions, it is being assumed in
would have to be reduced more and more to obtain the Ref.[13]thattheresultingstatehasenergyNE(1) andis
0
actualgapfromtheimaginary-timecorrelators. Also,fit- hence the N-string ground state. Generalizing the pro-
ting the correlation functions becomes more difficult as cedure to excited states, gapless excitations are found
they ultimately change from an exponential to an alge- which simply correspond to gapless excitations of a sin-
braic decay. Further evidence for the finite gap is given gle string. The described anti-symmetrization, also em-
by the temperature dependence of observables. When ployed in Refs. [19, 20], appears to be flawed. Different
9
-0.38 0.46
L
36
144
-0.4 0.45 1296
11664
-0.42
0.44
> / L -0.44 >
M 3 0.43
QI <
<H -0.46
0.42
-0.48
-0.5 L = 36 sites 0.41
144 sites
1296 sites
11664 sites
-0.52 0.4
-0.3 -0.25 -0.2 -0.15 -0.1 -0.05 0 -0.3 -0.25 -0.2 -0.15 -0.1 -0.05 0
V/t V/t
(a) (b)
Figure 7: Zoom near the first order phase transition at (V/t) for different lattice sizes with β = 19.2 and ∆β = 0.02: (a)
C
While the ground-state energy density, (cid:104)Hˆ (cid:105)/L, is continuous near (V/t) , (b), the density of flippable 3-plaquettes, (cid:104)ρˆ (cid:105),
QDM C 3
displays a small but evident jump.
from the conventional anti-symmetrization for fermions, and sublattice dimer densities in Fig 6b show large fluc-
theresultingN-stringwavefunctionisnotasumofprod- tuations in the interval 0.7<(V/t) <1.
C
uct states but contains also entangled states. Hence, the
The most natural explanation for this behavior are fi-
resulting state is not an energy eigenstate. [28]
nitesizeeffects,andthevanishingofthegapinthevicin-
Let us look at further observables to better under-
ity of the RK point which leads to an enhancement of
stand the plaquette phase. The normalized j-plaquette
fluctuations, the divergence of dimer-dimer correlation
numbers (cid:104)ρˆ(cid:105) are shown in Fig 5. They appear to be
i lengths, and a critical slowing down of the Monte Carlo
much more sensitive to variations in V/t than the RMS
simulation. More precisely, the observed effects can be
magnetization. As V/t increases, (cid:104)ρˆ (cid:105) and (cid:104)ρˆ (cid:105) contin-
3 0 attributedtoacrystallineregimewithapproximateU(1)
uously decreases while (cid:104)ρˆ (cid:105) increases, and (cid:104)ρˆ (cid:105) stays al-
2 1 symmetryinthevicinityoftheRKpoint. Thecontinuum
most constant, assuming its maximal value in the phase
version of the height representation [10] of the quantum
diagram. Theconstantandmaximalvalueof(cid:104)ρˆ (cid:105)≈0.25
1 dimermodelhasU(1)symmetryandalgebraicallydecay-
seems to be a characteristic signature for the plaquette
ing correlations at the RK point V/t = 1. For V/t < 1,
phase. For the ideal plaquette state |ψ (cid:105), one ob-
plaq close to the RK point, there are two length scales, one
tains ((cid:104)ρˆ (cid:105),(cid:104)ρˆ (cid:105),(cid:104)ρˆ (cid:105),(cid:104)ρˆ (cid:105)) = (1/12,1/4,1/4,5/12). For
0 1 2 3 beyond which dimer-dimer correlators show exponential
no value of V/t do we find agreement with these values,
decay signaling crystalline order, and one beyond which
showing once again the difference between the ideal and
one can observe the breaking of the U(1) symmetry. A
realplaquettestates. Thediscussionabouttheapproach
linearsystemsizein-betweenthesetwolengthscalescor-
to the RK point is postponed to section VD.
responds to the crystalline U(1) regime [21].
D. From the plaquette phase to the RK point
The current belief is that, for bipartite lattices, there
occurs a continuous transition from the plaquette phase E. The Rokhsar-Kivelson point (V/t=1)
to the RK point, the latter being an isolated critical
point. Some of our measured parameters, like the dimer
densities (Fig 5), show indeed the expected smooth be- The Rokhsar-Kivelson point is the only point of the
havior. Nevertheless, the magnetization curves, dis- phase diagram where the system does not display local
playedinFig.3,showasmallbumpbeforetheRKpoint, order. At this point, the Hamiltonian Hˆ becomes a
QDM
10
0.25 1.2
E(0)Star 19.2
0 EEVa(2Er,) PS6ltaxa6qr 1 1537386...648
81x81
-0.25
0.8
-0.5
Lt)
>/(M -0.75 E/t 0.6
D
Q
H
<
-1
0.4
-1.25
0.2
-1.5
-1.75 0
-2.5 -2 -1.5 -1 -0.5 0 0.5 1 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1
V/t V/t
Figure 8: Numerically computed energy density (cid:104)Hˆ (cid:105)/L Figure 9: Estimates for the energy gap ∆E/t to the first
QDM
for β =19.2, ∆β =0.02, and L=81×81, compared to vari- excited state for different temperatures. The gaps were ob-
ationalandperturbativeestimatesasdescribedinsectionVI tained from fits of imaginary-time auto-correlation functions
and appendix F, respectively. (cid:104)σˆ (0)σˆ (iτ)(cid:105) , for a system with L = 36×36 sites. The
i i QIM
resultsshouldbeinterpretedasupperboundstotherealgap,
whichareclosetheactualgapafterconvergenceinβ. Theer-
rorbarsonlyindicatethequalityofthefit,notthestatistical
sum of projectors
Monte-Carlo error.
Hˆ =−V (cid:88)(| (cid:105)(cid:104) |+h.c.)
QDM,RK i i
i
(cid:88) sum of projectors with positive coefficients,
+V (| (cid:105)(cid:104) |+| (cid:105)(cid:104) |)
i i i i
i
(cid:88) Hˆ =t(cid:88)(| (cid:105)−| (cid:105))·((cid:104) |−(cid:104) |)
=V (| (cid:105)−| (cid:105))·((cid:104) |−(cid:104) |). QDM i i i i
i i i i
i
i
(cid:88)
+(V −t) (| (cid:105)(cid:104) |+| (cid:105)(cid:104) |), (12)
i i i i
Therefore, the ground-state energy vanishes. For each
i
topological sector, and for each flip-connected sub-space
in a topological sector, one can build a zero-energy state
the ground state energy is non-negative. The ground
as an equal-amplitude superposition of all dimer cover-
states are isolated states, corresponding to dimer cover-
ings in the corresponding basis set.
ings without any flippable plaquettes, for which the en-
AttheRKpoint,manyphysicalproperties,likedimer-
ergy vanishes and which are generally called staggered
dimer correlations, can be derived from the behavior of
configurations. These states belong to the topological
the classical dimer problem. See for instance Ref. [22],
sectors with highest flux, one of which being displayed
wheretherelationbetweenquantumdimermodelsatthe
in Fig. 1c. Notice that such states are zero-energy eigen-
RKpointandtheirclassicalcounterpartsisdiscussed. In states of Hˆ for all values of V/t and become ground
QDM
particular, diagonal operator expectation values amount
statesforV/t≥1. ThetransitionontherightoftheRK
to doing classical enumerations. We used such computa-
point is abrupt. At the RK point, all topological sectors
tions to benchmark the QMC simulations.
sectors contain (at least one) state of vanishing energy,
whileonlytheisolatedgroundstatesinthemaximalflux
sector persist for V/t>1.
F. Staggered phase (1<V/t<∞) In the zero flux sector, the RK point corresponds to
a first order transition to states with a large majority of
In the parameter region 1 < V/t < ∞, flippable pla- 2-plaquettes ((cid:104)ρˆ (cid:105) > 0.8), vanishing (cid:104)ρˆ (cid:105), and finite but
2 0
quettes are disfavored. As the Hamiltonian becomes a small values of (cid:104)ρˆ (cid:105)=(cid:104)ρˆ (cid:105). See Fig. 5.
1 3